src/LocARNA/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 "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 unordered_map<idx_pair_t,
                              ArcMatch::idx_type,
                              pair_of_size_t_hash>::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