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