sequence_map_full.hpp

Go to the documentation of this file.

#ifndef SEQUENCE_MAP_FULL_HPP_
#define SEQUENCE_MAP_FULL_HPP_
#include <vector>
#include "boost/foreach.hpp"
#include "boost/utility/string_ref.hpp"
#include "vdj_pipe/best_match_pair.hpp"
#include "vdj_pipe/detail/id_bimap.hpp"
#include "vdj_pipe/detail/id_iterator.hpp"
#include "vdj_pipe/detail/id_map.hpp"
#include "vdj_pipe/detail/unused_variable.hpp"
#include "vdj_pipe/nucleotide_index.hpp"

namespace vdj_pipe{

class Seq_map_full {
public:
   typedef std::vector<Nucleotide> seq_type;
   typedef Best_match_pair<Mid_id, int> match_type;

private:
   typedef detail::Id_bimap<Mid_id, seq_type> seq_map;
   typedef detail::Id_map<Mid_id, std::string> name_map;

   static seq_type to_vector(const boost::string_ref s) {
      seq_type sv(s.size());
      for(std::size_t n = 0; n != sv.size(); ++n) {
         sv[n] = nucleotide_index(s[n]);
      }
      return sv;
   }

public:
   typedef Id_iterator<Mid_id> iterator;
   typedef iterator const_iterator;
   struct Err : public base_exception{};

   static std::string to_string(seq_type const& sv) {
      std::string s(sv.size(), 'X');
      for(std::size_t n = 0; n != sv.size(); ++n) {
         s[n] = to_capital(sv[n]);
      }
      return s;
   }

   explicit Seq_map_full(
      const bool allow_duplicate_seq = false
   )
   : seqs_(Mid_id(1)),
     names_(Mid_id(1)),
     allow_duplicate_seq_(allow_duplicate_seq),
     min_len_(std::numeric_limits<std::size_t>::max()),
     max_len_(0)
   {}

   iterator begin() const {return iterator(seqs_.min_id());}
   iterator end() const {return ++iterator(seqs_.max_id());}
   seq_type seq(const Mid_id id) const {return seqs_[id];}
   std::string seq_string(const Mid_id id) const {return to_string(seqs_[id]);}
   std::string const& name(const Mid_id id) const {return names_[id];}
   Mid_id max_id() const {return seqs_.max_id();}
   std::size_t min_len() const {return min_len_;}
   std::size_t max_len() const {return max_len_;}
   Mid_id const* find_exact(seq_type const& seq) const {return seqs_.find(seq);}

   Mid_id const* find_exact(const boost::string_ref seq) const {
      return find_exact(to_vector(seq));
   }

   template<class Seq> bool find_closest(
            const Seq s,
            match_type& m,
            scoring_matrix_t const& sm
   ) const {
      if( Mid_id const* id = find_exact(s) ) {
         return m.combine(*id, identity(seq(*id), s, sm));
      }

      bool found = false;
      BOOST_FOREACH(const Mid_id id, *this) {
         const int n = identity(seq(id), s, sm);
         if( m.combine(id, n) ) found = true;
      }
      return found;
   }

   Mid_id insert(std::string name, std::string const& seq) {
      if( name.empty() ) name = seq;

      BOOST_ASSERT(seqs_.size() == names_.size());
      const seq_type s = to_vector(seq);
      if( Mid_id const* id = seqs_.find(s) ) {
         if( allow_duplicate_seq_ ) {
            if( names_[*id] != name ) BOOST_THROW_EXCEPTION(
               Err()
               << Err::msg_t("same sequence, different name")
               << Err::str1_t(sanitize(seq))
               << Err::str2_t(sanitize(names_[*id]))
               << Err::str3_t(sanitize(name))
            );
            return *id;
         }
         BOOST_THROW_EXCEPTION(
                  Err()
                  << Err::msg_t("duplicate match element sequence")
                  << Err::str1_t(sanitize(seq))
                  << Err::str2_t(sanitize(name))
         );
      }

      const std::pair<Mid_id,bool> p1 = seqs_.insert(s);
      const Mid_id id = names_.insert(name);
      unused_variable(id);
      BOOST_ASSERT(id);
      BOOST_ASSERT(p1.first == id);
      if( min_len_ > seq.size() ) min_len_ = seq.size();
      if( max_len_ < seq.size() ) max_len_ = seq.size();
      return p1.first;
   }

private:
   seq_map seqs_;
   name_map names_;
   bool allow_duplicate_seq_;
   std::size_t min_len_, max_len_;
};

}//namespace vdj_pipe
#endif /* SEQUENCE_MAP_FULL_HPP_ */