ppALIGN API documentation
src/align.hpp
00001 /****************************************************************************** 00002 Copyright 2009 Stefan Wolfsheimer & Gregory Nuel. 00003 00004 This file is part of ppALIGN 00005 00006 ppALIGN is free software; you can redistribute it and/or modify 00007 it under the terms of the GNU General Public License as published by 00008 the Free Software Foundation; either version 2 of the License, or 00009 (at your option) any later version. 00010 00011 ppALIGN is distributed in the hope that it will be useful, 00012 but WITHOUT ANY WARRANTY; without even the implied warranty of 00013 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 00014 GNU General Public License for more details. 00015 00016 You should have received a copy of the GNU General Public License 00017 along with ppALIGN; if not, write to the Free Software 00018 Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA 00019 *******************************************************************************/ 00020 #ifndef _ALIGN_HPP_ 00021 #define _ALIGN_HPP_ 00022 00023 #include "sequence.hpp" 00024 #include "scorematrix.hpp" 00025 00026 #include <iostream> 00027 00030 class Pair 00031 { 00032 public: 00035 long i; 00036 00039 long j; 00040 00044 Pair() : i(0),j(0) {} 00045 00055 Pair(long _i, long _j) : i(_i),j(_j) {} 00056 00060 inline Pair & operator+=(const Pair & p) 00061 { i+= p.i; j+= p.j; return *this; } 00062 00063 00064 inline Pair & operator-=(const Pair & p) 00065 { i-= p.i; j-= p.j; return *this; } 00066 00067 friend inline Pair operator+(const Pair & p1, const Pair & p2) 00068 { Pair tmp(p1); tmp+= p2; return tmp; } 00069 00070 friend inline Pair operator-(const Pair & p1, const Pair & p2) 00071 { Pair tmp(p1); tmp-= p2; return tmp; } 00073 00074 00078 00079 friend inline bool operator < (const Pair & p1, const Pair & p2) 00080 { if(p1.i < p2.i) return true; 00081 else if(p1.i == p2.i) return p1.j < p2.j; 00082 else return false; } 00083 00084 friend inline bool operator <= (const Pair & p1, const Pair & p2) 00085 { 00086 if(p1.i < p2.i) return true; 00087 else if(p1.i == p2.i) return p1.j <= p2.j; 00088 else return false; 00089 } 00090 00091 friend inline bool operator > (const Pair & p1, const Pair & p2) 00092 { 00093 if(p1.i > p2.i) return true; 00094 else if(p1.i == p2.i) return p1.j > p2.j; 00095 else return false; 00096 } 00097 00098 friend inline bool operator >= (const Pair & p1, const Pair & p2) 00099 { 00100 if(p1.i > p2.i) return true; 00101 else if(p1.i == p2.i) return p1.j >= p2.j; 00102 else return false; 00103 } 00104 00105 friend inline bool operator == (const Pair & p1, const Pair & p2) 00106 { 00107 return (p1.i == p2.i) && (p1.j == p2.j); 00108 } 00109 00110 friend inline bool operator != (const Pair & p1, const Pair & p2) 00111 { return (p1.i != p2.i) || (p1.j != p2.j); } 00112 00114 00120 friend std::ostream & operator<<(std::ostream & ost, const Pair & p); 00121 00122 // friend std::ostream & operator<<(std::ostream & ost, const std::set<Pair> & s); 00123 }; 00124 00126 // 00127 // AlignState 00128 // 00130 00132 enum AlignState 00133 { 00134 StateMatch, /*<- Match or mismatch, corresponds to Pair(i,j) with i,j >= 0 */ 00135 StateInsert, /*<- Insertion corresponds to Pair(i,-1) with i >= 0 */ 00136 StateDelete, /*<- Insertion corresponds to Pair(-1,j) with j >= 0 */ 00137 StateStop /*<- begin or end of alignment */ 00138 00139 }; 00140 00143 std::ostream & operator << (std::ostream & ost, const AlignState & s); 00144 00148 class AlignEntry 00149 { 00150 public: 00156 Pair pair; 00157 00160 double posterior; 00161 00167 AlignEntry(long i, long j) : pair(i,j),posterior(0) { } 00168 00172 AlignEntry() : pair(-1,-1),posterior(0) {} 00173 00174 00178 inline operator const AlignState () const 00179 { if(pair.j == -1) return StateInsert; 00180 else if(pair.i ==-1) return StateDelete; 00181 else return StateMatch; 00182 } 00183 }; 00184 00185 00186 00189 class Align : public std::vector<AlignEntry> 00190 { 00191 00192 public: 00193 00196 typedef std::vector<AlignEntry> parent_t; 00197 00198 00202 class GapsAlignedError : public ExceptionBase 00203 { 00204 public: 00209 GapsAlignedError(size_t _pos); 00212 const size_t pos; 00213 }; 00214 00215 00216 00220 Align() {} 00221 00226 Align(const std::string & str); 00227 00228 00244 Align(const Sequence & s1, 00245 const Sequence & s2, 00246 Pair startPos = Pair(0,0)); 00247 00248 00254 void SetStartPos(Pair p); 00255 00256 00263 inline void push_back(const AlignEntry & v) 00264 { parent_t::push_back(v); } 00265 00266 00276 inline void push_back(long i, long j) 00277 { parent_t::push_back(AlignEntry(i,j)); } 00278 00283 inline Pair FirstPair() const 00284 { 00285 const_iterator itr; 00286 const_iterator end = parent_t::end(); 00287 for(itr = parent_t::begin(); itr != end && *itr != StateMatch; ++itr) ; 00288 if(itr != end) return itr->pair; 00289 else return Pair(-1,-1); 00290 } 00291 00296 inline Pair LastPair() const 00297 { 00298 const_reverse_iterator itr; 00299 const_reverse_iterator end = parent_t::rend(); 00300 Pair ret(-1,-1); 00301 for(itr = parent_t::rbegin(); itr != end && *itr != StateMatch; ++itr) ; 00302 if(itr != end) return itr->pair; 00303 else return Pair(-1,-1); 00304 } 00305 00306 00309 inline Pair SeqOffset() const 00310 { const_iterator itr; 00311 const_iterator end = parent_t::end(); 00312 Pair ret; 00313 for(itr = parent_t::begin(); itr != end && itr->pair.i == -1; ++itr) ; 00314 if(itr == end) ret.i = 0; 00315 else ret.i = itr->pair.i; 00316 for(itr = parent_t::begin(); itr != end && itr->pair.j == -1; ++itr); 00317 if(itr == end) ret.j = 0; 00318 else ret.j = itr->pair.j; 00319 return ret; 00320 } 00321 00328 inline Pair SeqEnd() const 00329 { const_reverse_iterator itr; 00330 const_reverse_iterator rend = parent_t::rend(); 00331 Pair ret; 00332 for(itr = parent_t::rbegin(); itr != rend && itr->pair.i == -1; ++itr) ; 00333 if(itr == rend) ret.i = 0; 00334 else ret.i = itr->pair.i; 00335 for(itr = parent_t::rbegin(); itr != rend && itr->pair.j == -1; ++itr); 00336 if(itr == rend) ret.j = 0; 00337 else ret.j = itr->pair.j; 00338 return ret + Pair(1,1); 00339 } 00340 00341 // @} 00342 00344 friend std::ostream & operator<<(std::ostream & ost, const Align & a); 00345 00346 00347 00350 00360 int Score(const Sequence::const_iterator & first_begin, 00361 const Sequence::const_iterator & second_begin, 00362 const ScoreMatrix & scoreMatrix, 00363 int gap_open, 00364 int gap_extension, 00365 bool local=false) const; 00366 00376 inline int Score(const Sequence & s1, 00377 const Sequence & s2, 00378 const ScoreMatrix & scoreMatrix, 00379 int gap_open, 00380 int gap_extension, 00381 bool local=false) const 00382 { 00383 return Score(s1.begin(), 00384 s2.begin(), 00385 scoreMatrix, 00386 gap_open, 00387 gap_extension, 00388 local); } 00389 00393 inline int Gaps() const 00394 { int ret=0; 00395 for(const_iterator itr = begin(); 00396 itr!= end(); 00397 ++itr) 00398 if(*itr != StateMatch) ret++; 00399 return ret; 00400 } 00401 00402 inline int Identity(const Sequence & s1, 00403 const Sequence & s2) const 00404 { 00405 int ret = 0; 00406 for(const_iterator itr = begin(); 00407 itr!= end(); 00408 ++itr) 00409 if(*itr == StateMatch && s1[itr->pair.i]==s2[itr->pair.j]) ret++; 00410 return ret; 00411 } 00412 00416 inline int Positive(const Sequence & s1, 00417 const Sequence & s2, 00418 const ScoreMatrix & s) const 00419 { 00420 int ret = 0; 00421 for(const_iterator itr = begin(); 00422 itr!= end(); 00423 ++itr) 00424 if(*itr == StateMatch && s[s1[itr->pair.i]][s2[itr->pair.j]] > 0) ret++; 00425 return ret; 00426 } 00427 00437 Pair MakeLocal(Sequence & seq1, 00438 Sequence & seq2); 00439 00452 Pair MakeLocal(Align & align, 00453 Sequence & seq1, 00454 Sequence & seq2, 00455 long offset=0) const; 00456 00457 00458 }; 00459 00460 00461 #endif