sequence_store.hpp

Go to the documentation of this file.

#ifndef SEQUENCE_STORE_HPP_
#define SEQUENCE_STORE_HPP_
#include "boost/assert.hpp"
#include "boost/foreach.hpp"
#include "boost/multi_index_container.hpp"
#include "boost/multi_index/hashed_index.hpp"
#include "boost/multi_index/ordered_index.hpp"
#include "boost/range.hpp"
#include "boost/shared_ptr.hpp"

#include "vdj_pipe/detail/get_by_id.hpp"
#include "vdj_pipe/detail/id_map.hpp"
#include "vdj_pipe/detail/string_ref.hpp"
#include "vdj_pipe/detail/vector_set.hpp"
#include "vdj_pipe/exception.hpp"
#include "vdj_pipe/object_ids.hpp"
#include "vdj_pipe/sequence_entry.hpp"
#include "vdj_pipe/sequence_position.hpp"

namespace vdj_pipe{

class Seq_store {
public:
   typedef detail::Id_map<Seq_id, Seq_entry> seq_map;
   typedef boost::shared_ptr<seq_map> seq_map_ptr;

private:
   typedef detail::Getter<
            seq_map, Seq_id, Seq_entry,
            std::string const&, &Seq_entry::sequence
            > get_seq;

   typedef detail::Getter<
            seq_map, Seq_id, Seq_entry,
            unsigned, &Seq_entry::size
            > get_size;

   typedef boost::multi_index_container<
            Seq_id,
            boost::multi_index::indexed_by<
               boost::multi_index::hashed_unique<
                  boost::multi_index::tag<struct seq_tag>,
                  get_seq
               >,
               boost::multi_index::ordered_non_unique<
                  boost::multi_index::tag<struct size_tag>,
                  get_size
               >
            >
   > seq_mi_t;
   typedef seq_mi_t::index<seq_tag>::type seq_index;
   typedef seq_mi_t::index<size_tag>::type size_index;

   static seq_mi_t::ctor_args_list seq_index_init(seq_map const& sm) {
      return boost::make_tuple(
               boost::make_tuple(
                        0,
                        get_seq(sm),
                        boost::hash<std::string>(),
                        std::equal_to<std::string>()
               ),
               boost::make_tuple(get_size(sm), std::less<unsigned>())
      );
   }

public:
   typedef detail::Vector_set<sub_seq> sub_seq_set;
   typedef detail::Vector_set<super_seq> super_seq_set;

private:
   typedef detail::Id_map<Seq_id, sub_seq_set> seq_to_info;
   typedef detail::Id_map<Read_id, super_seq_set> info_to_seq;
   typedef size_index::const_iterator size_iterator;
   typedef boost::iterator_range<size_iterator> size_range;

public:
   typedef seq_index::iterator iterator;
   typedef seq_index::const_iterator const_iterator;

   Seq_store()
   : sm_(new seq_map(Seq_id(1))),
     si_(seq_index_init(*sm_)),
     s2i_(Seq_id(1)),
     i2s_(Read_id(1))
   {}

   std::size_t size() const {return si_.size();}
   const_iterator begin() const {return si_.begin();}
   const_iterator end() const {return si_.end();}
   bool empty() const {return si_.empty();}
   Seq_entry const& operator[](const Seq_id sid) const {return (*sm_)[sid];}
   boost::shared_ptr<seq_map> sequence_map() {return sm_;}

   size_range by_size() const {return boost::make_iterator_range(si_.get<size_tag>());}

   size_range by_size(const unsigned size) const {
      return boost::make_iterator_range(
               si_.get<size_tag>().equal_range(size)
      );
   }

   size_range by_size_range(
            const unsigned from,
            const unsigned to = std::numeric_limits<unsigned>::max()
   ) const {
      return boost::make_iterator_range(
               si_.get<size_tag>().lower_bound(from),
               si_.get<size_tag>().upper_bound(to)
      );
   }

   Seq_id insert(const Read_id rid, const boost::string_ref seq) {
      seq_index& ind = si_.get<seq_tag>();
      const const_iterator i = ind.find(
               seq,
               boost::hash<boost::string_ref>(),
               detail::Equal_string_ref()
      );

      Seq_id sid;
      if( i == ind.end() ) {
         sid = sm_->insert(Seq_entry(seq, rid));
         si_.insert(sid);
      } else {
         sid = *i;
      }
      s2i_.insert(sid).insert(sub_seq(rid, 0));
      i2s_.insert(rid).insert(super_seq(sid, 0));
      return sid;
   }

   template<class Super_seqs> void
   remove_subsequence(const Seq_id sid, Super_seqs const& ss) {
      //re-map all Read_id-s to sequences in ss and back
      BOOST_FOREACH(sub_seq const& subs, s2i_[sid]) {
         const Read_id rid = subs.id_;
         super_seq_set& super_ss = i2s_[rid];
         BOOST_ASSERT(super_ss.find(sid));
         super_ss.erase(sid);
         BOOST_FOREACH(super_seq const& supers, ss) {
            const unsigned pos = subs.pos_ + supers.pos_;
            super_ss.insert(super_seq(supers.id_, pos));
            s2i_[supers.id_].insert(sub_seq(rid, pos));
         }
      }
      //remove sid
      s2i_.erase(sid);
      si_.get<seq_tag>().erase((*sm_)[sid].sequence());
      sm_->erase(sid);
   }

   super_seq_set const& maps_to(const Read_id iid) const {return i2s_[iid];}
   sub_seq_set const& maps_from(const Seq_id sid) const {return s2i_[sid];}

private:
   seq_map_ptr sm_;
   seq_mi_t si_;
   seq_to_info s2i_;
   info_to_seq i2s_;
};

}//namespace vdj_pipe
#endif /* SEQUENCE_STORE_HPP_ */