alg/seq_span.cpp

   1 #include <stdexcept>
   2 #include <set>
   3 #include <algorithm>
   4
   5 #include "seq_span.hpp"
   6 #include "mussa_exceptions.hpp"
   7
   8 SeqSpan::SeqSpan(const SeqSpan &o)
   9   : seq(o.seq),
  10     seq_start(o.seq_start),
  11     seq_count(o.seq_count),
  12     parent(o.parent)
  13 {
  14 }
  15
  16 SeqSpan::SeqSpan(const SeqSpan *p)
  17   : seq(p->seq),
  18     seq_start(p->seq_start),
  19     seq_count(p->seq_count),
  20     parent(p->parent)
  21 {
  22 }
  23
  24 SeqSpan::SeqSpan(const std::string &seq_,
  25                  AlphabetRef a)
  26   : seq(new SeqString(seq_, a)),
  27     seq_start(0),
  28     seq_count(seq_.length()),
  29     parent()
  30 {
  31 }
  32
  33 SeqSpan::SeqSpan(const SeqSpanRef parent_, size_type start_, size_type count_)
  34   : seq(parent_->seq),
  35     seq_start(parent_->seq_start + start_),
  36     parent(parent_)
  37 {
  38   if (count_ == npos)
  39     seq_count = parent_->seq_count;
  40   else
  41     seq_count = count_;
  42 }
  43
  44 //////
  45 // Helpful operators
  46 SeqSpan &SeqSpan::operator=(const SeqSpan& s)
  47 {
  48   if (this != &s) {
  49     seq = s.seq;
  50     seq_start = s.seq_start;
  51     seq_count = s.seq_count;
  52     parent = s.parent;
  53   }
  54   return *this;
  55 }
  56
  57 std::ostream& operator<<(std::ostream& out, const SeqSpan& s)
  58 {
  59   out << s.sequence();
  60 }
  61
  62 /* Not implemented yet
  63 //! compare two spans
  64 //! \throws sequence_invalid_comparison
  65 friend bool operator<(const SeqSpan&, const SeqSpan&);
  66 bool operator<(const SeqSpan& a, const SeqSpan& b)
  67 {
  68   // are we subcomponents of the same sequence region?
  69   if (a.seq.get() == b.seq.get()) {
  70     if (a.seq_start < b.seq_start)
  71       return true;
  72     else
  73       return false;
  74   } else {
  75     throw mussa_error("can only compare two spans from the same sequence");
  76   }
  77 }
  78 */
  79 #include <iostream>
  80 bool operator==(const SeqSpan& a, const SeqSpan& b)
  81 {
  82   if (SeqSpan::isFamily(a, b)) {
  83     // can do fast comparison
  84     if (a.seq_start == b.seq_start and a.seq_count == b.seq_count) {
  85       return true;
  86     } else {
  87       return false;
  88     }
  89   }
  90   return false;
  91 }
  92
  93 bool operator!=(const SeqSpan& a, const SeqSpan& b)
  94 {
  95   return not (a == b);
  96 }
  97
  98 SeqSpan::const_reference SeqSpan::operator[](SeqSpan::size_type i) const
  99 {
 100   return at(i);
 101 }
 102
 103 SeqSpan::const_reference SeqSpan::at(SeqSpan::size_type i) const
 104 {
 105   if (!seq) throw std::out_of_range("empty sequence");
 106   return seq->at(i+seq_start);
 107 }
 108
 109 const char *SeqSpan::data() const
 110 {
 111   if (seq) {
 112     return seq->c_str()+seq_start;
 113   } else
 114     return 0;
 115 }
 116
 117 SeqSpan::const_iterator SeqSpan::begin() const
 118 {
 119   if (seq and seq_count != 0)
 120     return seq->begin()+seq_start;
 121   else
 122     return SeqSpan::const_iterator(0);
 123 }
 124
 125 SeqSpan::const_iterator SeqSpan::end() const
 126 {
 127   if (seq and seq_count != 0) {
 128     return seq->begin() + seq_start + seq_count;
 129   } else {
 130     return SeqSpan::const_iterator(0);
 131   }
 132 }
 133
 134 SeqSpan::const_reverse_iterator SeqSpan::rbegin() const
 135 {
 136   if (seq and seq_count != 0)
 137     return seq->rbegin()+(seq->size()-(seq_start+seq_count));
 138   else
 139     return SeqSpan::const_reverse_iterator();
 140 }
 141
 142 SeqSpan::const_reverse_iterator SeqSpan::rend() const
 143 {
 144   if (seq and seq_count != 0) {
 145     return rbegin() + seq_count;
 146   } else {
 147     return SeqSpan::const_reverse_iterator();
 148   }
 149 }
 150
 151 bool SeqSpan::empty() const
 152 {
 153   return (seq_count == 0) ? true : false;
 154 }
 155
 156 SeqSpan::size_type SeqSpan::find_first_not_of(
 157   const std::string& query,
 158   SeqSpan::size_type index) const
 159 {
 160   typedef std::set<std::string::value_type> sequence_set;
 161   sequence_set match_set;
 162
 163   for(const_iterator query_item = query.begin();
 164       query_item != query.end();
 165       ++query_item)
 166   {
 167     match_set.insert(*query_item);
 168   }
 169   for(const_iterator base = begin();
 170       base != end();
 171       ++base)
 172   {
 173     if(match_set.find(*base) == match_set.end()) {
 174       return base-begin();
 175     }
 176   }
 177   return SeqSpan::npos;
 178 }
 179
 180 void SeqSpan::setStart(SeqSpan::size_type v)
 181 {
 182   if (v > stop()) {
 183     // cry
 184     throw mussa_error("can't set Start > Stop");
 185   }
 186   seq_count += seq_start - v;
 187   seq_start = v;
 188 }
 189
 190 void SeqSpan::setStop(SeqSpan::size_type v)
 191 {
 192   if ( v < start() ) {
 193     // negative sized sequences are bad
 194     throw mussa_error("can't set Stop < Start");
 195   }
 196   seq_count = std::min<size_type>(v - seq_start, parentSize()-parentStart());
 197 }
 198
 199 SeqSpan::size_type SeqSpan::parentStart() const
 200 {
 201   if (!parent) {
 202     // no parent
 203     return start();
 204   } else {
 205     return start() - parent->start();
 206   }
 207 }
 208
 209 void SeqSpan::setParentStart(SeqSpan::size_type v)
 210 {
 211   setStart(parent->start() + v);
 212 }
 213
 214 SeqSpan::size_type SeqSpan::parentStop() const
 215 {
 216   if (!parent) {
 217     // no parent
 218     return stop();
 219   } else {
 220     return stop() - parent->start();
 221   }
 222 }
 223
 224 void SeqSpan::setParentStop(SeqSpan::size_type v)
 225 {
 226   setStop(parent->start() + v);
 227 }
 228
 229 SeqSpanRef SeqSpan::subseq(size_type start, size_type count)
 230 {
 231   count = std::min<size_type>(count, seq_count - start);
 232
 233   SeqSpanRef new_span(new SeqSpan(this->shared_from_this(), start, count));
 234   return new_span;
 235 }
 236
 237 std::string SeqSpan::sequence() const
 238 {
 239   if (seq) {
 240     return seq->substr(seq_start, seq_count);
 241   } else {
 242     return std::string();
 243   }
 244 }
 245
 246 bool SeqSpan::isFamily(const SeqSpan& a, const SeqSpan& b)
 247 {
 248   return a.seq.get() == b.seq.get();
 249 }