arc_matches.hh

#ifndef LOCARNA_ARC_MATCHES_HH
#define LOCARNA_ARC_MATCHES_HH
 
#ifdef HAVE_CONFIG_H
#include <config.h>
#endif
 
#include <algorithm>
#include <vector>
#include <unordered_map>
 
#include "scoring_fwd.hh"
#include "aux.hh"
#include "matrix.hh"
#include "basepairs.hh"
 
#include <assert.h>
 
namespace LocARNA {
    class Scoring;
    class Sequence;
    class RnaData;
    class AnchorConstraints;
    class TraceController;
    class MatchController;
 
    class ArcMatch {
    public:
        typedef std::vector<int>::size_type size_type; 
        typedef size_type idx_type;                    
        typedef BasePairs__Arc Arc;                    
    private:
        const Arc *arcA_; 
        const Arc *arcB_; 
        idx_type idx_;    
 
    public:
        ArcMatch(const Arc *arcA, const Arc *arcB, idx_type idx)
            : arcA_(arcA), arcB_(arcB), idx_(idx) {}
 
        const Arc &
        arcA() const {
            return *arcA_;
        }
 
        const Arc &
        arcB() const {
            return *arcB_;
        }
 
        idx_type
        idx() const {
            return idx_;
        }
    };
 
    typedef std::vector<ArcMatch> ArcMatchVec;
 
    typedef std::vector<ArcMatch::idx_type> ArcMatchIdxVec;
 
    class ArcMatches {
    public:
        typedef std::vector<int>::size_type size_type; 
        typedef BasePairs__Arc Arc;                    
    protected:
        size_type lenA; 
        size_type lenB; 
 
        BasePairs *bpsA; 
        BasePairs *bpsB; 
 
        /* Constraints and Heuristics */
 
        size_type max_length_diff; 
 
        size_type max_diff_at_am; 
 
        const MatchController &match_controller; 
 
        const AnchorConstraints &constraints; 
 
        bool
        is_valid_arcmatch(const Arc &arcA, const Arc &arcB) const;
 
        /* END constraints and heuristics */
 
        bool maintain_explicit_scores; 
 
        ArcMatchVec arc_matches_vec;
 
        size_type number_of_arcmatches;
 
        std::vector<score_t> scores;
 
        Matrix<ArcMatchIdxVec> common_right_end_lists;
 
        Matrix<ArcMatchIdxVec> common_left_end_lists;
 
        ArcMatchIdxVec inner_arcmatch_idxs;
 
        void
        init_inner_arc_matchs();
 
        class lex_greater_left_ends {
            const ArcMatches &arc_matches;
 
        public:
            explicit lex_greater_left_ends(const ArcMatches &arc_matches_)
                : arc_matches(arc_matches_) {}
 
            bool
            operator()(const ArcMatch::idx_type &i,
                       const ArcMatch::idx_type &j) const {
                size_type ali = arc_matches.arcmatch(i).arcA().left();
                size_type bli = arc_matches.arcmatch(i).arcB().left();
                size_type alj = arc_matches.arcmatch(j).arcA().left();
                size_type blj = arc_matches.arcmatch(j).arcB().left();
 
                return (ali > alj) || (ali == alj && bli > blj);
            }
        };
 
        class tuple5 {
        public:
            typedef std::vector<int>::size_type size_type; 
 
            size_type i;   
            size_type j;   
            size_type k;   
            size_type l;   
            score_t score; 
 
            tuple5(size_type i_,
                   size_type j_,
                   size_type k_,
                   size_type l_,
                   score_t score_)
                : i(i_), j(j_), k(k_), l(l_), score(score_) {}
        };
 
    public:
        ArcMatches(const Sequence &seqA_,
                   const Sequence &seqB_,
                   const std::string &arcmatch_scores_file,
                   int probability_scale,
                   size_type max_length_diff,
                   size_type max_diff_at_am,
                   const MatchController &trace_controller,
                   const AnchorConstraints &constraints);
 
        ArcMatches(const RnaData &rnadataA,
                   const RnaData &rnadataB,
                   double min_prob,
                   size_type max_length_diff,
                   size_type max_diff_at_am,
                   const MatchController &trace_controller,
                   const AnchorConstraints &constraints);
 
        ~ArcMatches();
 
        // for the mea probabilistic consistency transformation, support to read
        // and write the arcmatch scores
        // this allows in general to have user defined arc-match scores
 
        void
        read_arcmatch_scores(const std::string &arcmatch_scores_file,
                             int probability_scale);
 
        void
        write_arcmatch_scores(const std::string &arcmatch_scores_file,
                              const Scoring &scoring) const;
 
        const BasePairs &
        get_base_pairsA() const {
            return *bpsA;
        }
 
        const BasePairs &
        get_base_pairsB() const {
            return *bpsB;
        }
 
        bool
        explicit_scores() const {
            return maintain_explicit_scores;
        }
 
        void
        make_scores_explicit(const Scoring &scoring);
 
        score_t
        get_score(const ArcMatch &am) const {
            assert(maintain_explicit_scores);
            return scores[am.idx()];
        }
 
        size_type
        num_arc_matches() const {
            return number_of_arcmatches;
        }
 
        const ArcMatch &
        arcmatch(size_type idx) const {
            assert(idx < number_of_arcmatches);
            return arc_matches_vec[idx];
        }
 
        // ============================================================
        // Iteration over arc matches
        //
 
        const ArcMatchIdxVec &
        common_right_end_list(size_type i, size_type j) const {
            return common_right_end_lists(i, j);
        }
 
        const ArcMatchIdxVec &
        common_left_end_list(size_type i, size_type j) const {
            return common_left_end_lists(i, j);
        }
 
        // ============================================================
 
        void
        get_max_right_ends(size_type al,
                           size_type bl,
                           size_type *max_ar,
                           size_type *max_br,
                           bool no_lonely_pairs) const;
 
        void
        get_min_right_ends(size_type al,
                           size_type bl,
                           size_type *min_ar,
                           size_type *min_br) const;
 
        // ------------------------------------------------------------
        // inner arc matches
 
        bool
        exists_inner_arc_match(const ArcMatch &am) const {
            return inner_arcmatch_idxs[am.idx()] < num_arc_matches();
        }
 
        const ArcMatch &
        inner_arc_match(const ArcMatch &am) const {
            return arcmatch(inner_arcmatch_idxs[am.idx()]);
        }
 
        void
        sort_right_adjacency_lists();
 
        // ------------------------------------------------------------
        // iteration (in no specific order)
 
        typedef ArcMatchVec::const_iterator const_iterator;
 
        const_iterator
        begin() const {
            return arc_matches_vec.begin();
        }
 
        const_iterator
        end() const {
            return arc_matches_vec.begin() + number_of_arcmatches;
        }
    };
 
    class ArcMatchesIndexed : public ArcMatches {
    public:
        ArcMatchesIndexed(const Sequence &seqA_,
                          const Sequence &seqB_,
                          const std::string &arcmatch_scores_file,
                          int probability_scale,
                          size_type max_length_diff,
                          size_type max_diff_at_am,
                          const MatchController &trace_controller,
                          const AnchorConstraints &constraints)
            : ArcMatches(seqA_,
                         seqB_,
                         arcmatch_scores_file,
                         probability_scale,
                         max_length_diff,
                         max_diff_at_am,
                         trace_controller,
                         constraints),
              am_index_() {
            build_arcmatch_index();
        }
 
        ArcMatchesIndexed(const RnaData &rnadataA,
                          const RnaData &rnadataB,
                          double min_prob,
                          size_type max_length_diff,
                          size_type max_diff_at_am,
                          const MatchController &trace_controller,
                          const AnchorConstraints &constraints)
            : ArcMatches(rnadataA,
                         rnadataB,
                         min_prob,
                         max_length_diff,
                         max_diff_at_am,
                         trace_controller,
                         constraints),
              am_index_() {
            build_arcmatch_index();
        }
 
    private:
        typedef std::pair<size_type, size_type> idx_pair_t;
 
        typedef std::unordered_map<idx_pair_t, ArcMatch::idx_type>
            am_index_type;
 
        am_index_type am_index_;
 
        void
        build_arcmatch_index();
 
    public:
        const ArcMatch::idx_type
        invalid_am_index() const {
            // The invalid arc match index is implemented to be the
            // maximum valid index + 1.
            // This allows an efficient implementation, where we push
            // an invalid arc match to the end of vector arc_matches_vec.
            // Thus, we return size-1!
 
            return number_of_arcmatches;
        }
 
        const ArcMatch::idx_type
        am_index(const size_type &arcAIdx, const size_type &arcBIdx) const {
            am_index_type::const_iterator it =
                am_index_.find(idx_pair_t(arcAIdx, arcBIdx));
            if (am_index_.end() != it) {
                return it->second;
            } else {
                return invalid_am_index();
            }
        }
 
        const ArcMatch &
        am_index(const Arc &arcA, const Arc &arcB) const {
            return arc_matches_vec[am_index(arcA.idx(), arcB.idx())];
        }
    };
 
} // end namespace LocARNA
 
#endif // LOCARNA_ARC_MATCHES_HH