1 // ==========================================================================
2 // SeqAn - The Library for Sequence Analysis
3 // ==========================================================================
4 // Copyright (c) 2006-2010, Knut Reinert, FU Berlin
5 // All rights reserved.
6 //
7 // Redistribution and use in source and binary forms, with or without
8 // modification, are permitted provided that the following conditions are met:
9 //
10 // * Redistributions of source code must retain the above copyright
11 // notice, this list of conditions and the following disclaimer.
12 // * Redistributions in binary form must reproduce the above copyright
13 // notice, this list of conditions and the following disclaimer in the
14 // documentation and/or other materials provided with the distribution.
15 // * Neither the name of Knut Reinert or the FU Berlin nor the names of
16 // its contributors may be used to endorse or promote products derived
17 // from this software without specific prior written permission.
18 //
19 // THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS"
20 // AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
21 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
22 // ARE DISCLAIMED. IN NO EVENT SHALL KNUT REINERT OR THE FU BERLIN BE LIABLE
23 // FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL
24 // DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR
25 // SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER
26 // CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT
27 // LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY
28 // OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH
29 // DAMAGE.
30 //
31 // ==========================================================================
32
33 #ifndef SEQAN_HEADER_BLAST_HSP_H
34 #define SEQAN_HEADER_BLAST_HSP_H
35
36
37 namespace SEQAN_NAMESPACE_MAIN
38 {
39
40
41
42
43
44
45
46 template<typename TBlastSpec = NucleotideBlast<>, typename TInfoSpec = BasicInfo>
47 class BlastHsp;
48
49
50
51
52 /////////////////////////////////////////////////////////////////////////////
53 //////////////////////// FullInfo Nucleotide Spec //////////////////////////////
54 /////////////////////////////////////////////////////////////////////////////
55
56 //Info types
57 /**
58 .Spec.FullInfo:
59 ..cat:Blast
60 ..general:Class.BlastHsp
61 ..summary:Stores all pieces of information delivered with an alignment in a Blast report.
62 ..signature:FullInfo
63 ..include:blast.h
64 */
65
66 /**
67 .Spec.BasicInfo:
68 ..cat:Blast
69 ..general:Class.BlastHsp
70 ..summary:Stores only the basic pieces of information delivered with an alignment in a Blast report.
71 ..signature:BasicInfo
72 ..include:blast.h
73 */
74
75
76
77 /**
78 .Class.BlastHsp:
79 ..cat:Blast
80 ..summary:Object for storing Blast HSPs.
81 ..signature:BlastHsp<TBlastSpec, TInfoSpec>
82 ..param.TBlastSpec:The type of Blast report to be parsed.
83 ...type:Spec.BlastN
84 ...type:Spec.BlastP
85 ...default:BlastN
86 ..param.TInfoSpec:The specializing type determining the amount of information to be stored.
87 ...type:Spec.BasicInfo
88 ...type:Spec.FullInfo
89 ...default:BasicInfo
90 ...remarks:BasicInfo stores begin and end positions on query and database sequence, as well as the alignment. FullInfo stores additional information such as score, e-value...
91 ..include:blast.h
92 */
93 template<typename TSpec>
94 class BlastHsp<NucleotideBlast<TSpec>, FullInfo>
95 {
96
97 public:
98 float score;
99 float bits;
100 double expect;
101 unsigned int identity;
102 unsigned int gaps;
103 unsigned int abs_gaps;
104 bool query_strand; //false=>minus, true=>plus
105 bool db_strand;
106 unsigned int query_begin;
107 unsigned int db_begin;
108 unsigned int query_end; //included in the alignment
109 unsigned int db_end; //included in the alignment
110
111 String<char> query_string;
112 String<char> db_string;
113
114
BlastHsp()115 BlastHsp()
116 {
117 SEQAN_CHECKPOINT
118 //defaults for those values that are not necessarily present in each HSP in the Blast Report
119 score = 0.0;
120 bits = 0.0;
121 expect = 0.0;
122 identity = 0;
123 gaps = 0;
124 abs_gaps = 0;
125 query_strand = true; //false=>minus, true=>plus
126 db_strand = true;
127 }
128
BlastHsp(BlastHsp const & other)129 BlastHsp(BlastHsp const& other)
130 {
131 SEQAN_CHECKPOINT
132 score = other.score;
133 bits = other.bits;
134 expect = other.expect;
135 identity = other.identity;
136 gaps = other.gaps;
137 abs_gaps = other.abs_gaps;
138 query_strand = other.query_strand;
139 db_strand = other.db_strand;
140 query_begin = other.query_begin;
141 db_begin = other.db_begin;
142 query_end = other.query_end;
143 db_end = other.db_end;
144 query_string = other.query_string;
145 db_string = other.db_string;
146
147 }
148
149 BlastHsp & operator = (BlastHsp const & other)
150 {
151 SEQAN_CHECKPOINT
152 score = other.score;
153 bits = other.bits;
154 expect = other.expect;
155 identity = other.identity;
156 gaps = other.gaps;
157 abs_gaps = other.abs_gaps;
158 query_strand = other.query_strand;
159 db_strand = other.db_strand;
160 query_begin = other.query_begin;
161 db_begin = other.db_begin;
162 query_end = other.query_end;
163 db_end = other.db_end;
164 query_string = other.query_string;
165 db_string = other.db_string;
166
167 return *this;
168 }
169
170
~BlastHsp()171 ~BlastHsp()
172 {
173 SEQAN_CHECKPOINT
174 }
175 };
176
177
178 //parse BlastHsp
179 template<typename TFile, typename TChar, typename TSpec>
180 inline typename Position<TFile>::Type
_parseBlastHsp(TFile & file,TChar & c,BlastHsp<NucleotideBlast<TSpec>,FullInfo> & hsp)181 _parseBlastHsp(TFile & file,
182 TChar & c,
183 BlastHsp<NucleotideBlast<TSpec>,FullInfo > & hsp)
184 {
185 SEQAN_CHECKPOINT
186 typedef typename Position<TFile>::Type TPosition;
187
188 clear(hsp);
189 String<char> pword;
190 int pint;
191 float pfloat;
192 double pdouble;
193 TPosition start_pos,act_pos;
194 act_pos = _streamTellG(file);
195
196 String<char> query_string;
197 String<char> db_string;
198
199 _parseSkipWhitespace(file,c);
200 c = _streamGet(file);
201 _parseSkipWhitespace(file,c);
202 pfloat = (float)_parseReadEValue(file, c);
203 hsp.bits = pfloat;
204 if(_parseLineUntil(file,c,'('))
205 {
206 c = _streamGet(file);
207 pfloat = _parseReadFloat(file, c);
208 hsp.score = pfloat;
209 }
210
211 _parseSkipWhitespace(file,c);
212 String<char> query = "Expect";
213 if(_parseLineUntil(file,c,query,6))
214 {
215 c = _streamGet(file);
216 _parseSkipWhitespace(file,c);
217 if(c == '=')
218 c = _streamGet(file);
219 _parseSkipWhitespace(file,c);
220 pdouble = _parseReadEValue(file, c);
221 hsp.expect = pdouble;
222 }
223
224 query = "Identities";
225 if(_parseUntilBeginLine(file,c,query,10,3))
226 {
227 _parseLineUntil(file,c,'(');
228 c = _streamGet(file);
229 pint = _parseReadNumber(file, c);
230 hsp.identity = pint;
231 }
232
233 query = "Gaps";
234 if(_parseLineUntil(file,c,query,4))
235 {
236 _parseLineUntil(file,c,'=');
237 c = _streamGet(file);
238 pint = _parseReadNumber(file, c);
239 hsp.abs_gaps = pint;
240 _parseLineUntil(file,c,'(');
241 c = _streamGet(file);
242 pint = _parseReadNumber(file, c);
243 hsp.gaps = pint;
244 }
245 //else
246 //{
247 // hsp.abs_gaps = 0;
248 // hsp.gaps = 0;
249 //}
250
251 query = "Strand";
252 if(_parseUntilBeginLine(file,c,query,6,3))
253 {
254 c = _streamGet(file);
255 _parseSkipWhitespace(file,c);
256 if(c == '=')
257 c = _streamGet(file);
258 _parseSkipWhitespace(file,c);
259 pword = _parseReadWord(file, c);
260 if(pword == "Plus")
261 hsp.query_strand = true;
262 else
263 hsp.query_strand = false;
264 c = _streamGet(file);
265 _parseSkipWhitespace(file,c);
266 if(c == '/')
267 c = _streamGet(file);
268 _parseSkipWhitespace(file,c);
269 pword = _parseReadWord(file, c);
270 if(pword == "Plus")
271 hsp.db_strand = true;
272 else
273 hsp.db_strand = false;
274 }
275
276 bool first = true;
277 query = "Query";
278 if(_parseUntilBeginLine(file,c,query,5))
279 {
280 bool in_hsp = true;
281 while(in_hsp)
282 {
283 int end_query,end_db = -1;
284 _parseSkipWhitespace(file,c);
285 if(c == ':')
286 c = _streamGet(file);
287 _parseSkipWhitespace(file,c);
288 pint = _parseReadNumber(file,c);
289 _parseSkipWhitespace(file,c);
290 //TPosition act_posQ = _streamTellG(file);
291 if(first)
292 hsp.query_begin = pint;
293 query_string += _parseReadAlignmentString(file,c);
294 _parseSkipWhitespace(file,c);
295 end_query = _parseReadNumber(file,c);
296 query = "Sbjct";
297 _parseUntilBeginLine(file,c,query,5);
298 _parseSkipWhitespace(file,c);
299 if(c == ':')
300 c = _streamGet(file);
301 _parseSkipWhitespace(file,c);
302 pint = _parseReadNumber(file,c);
303 _parseSkipWhitespace(file,c);
304 //TPosition act_posS = _streamTellG(file);
305 if(first)
306 hsp.db_begin = pint;
307 db_string += _parseReadAlignmentString(file,c);
308 _parseSkipWhitespace(file,c);
309 end_db = _parseReadNumber(file,c);
310 first = false;
311 String<TChar> delim = "QS>R";
312 if(_parseUntilBeginLineOneOf(file,c,delim,4))
313 {
314 TPosition pos = _streamTellG(file);
315 pword = _parseReadWord(file,c);
316
317 //still in the same HSP
318 if(pword == "Query")
319 in_hsp = true;
320 else
321 {
322 hsp.query_end = end_query;
323 hsp.db_end = end_db;
324 hsp.query_string = query_string;
325 hsp.db_string = db_string;
326 //next HSP
327 if(pword == "Score")
328 return _streamTellG(file);
329 if(pword[0] == '>')
330 return act_pos;
331 if(pword[0] == 'R')
332 return (TPosition) 0;
333 in_hsp = false;
334 }
335 }
336 else
337 return act_pos;
338 }
339 }
340
341 return act_pos;
342
343 }
344
345
346 template<typename TSpec>
347 inline void
clear(BlastHsp<NucleotideBlast<TSpec>,FullInfo> & blastHsp)348 clear(BlastHsp<NucleotideBlast<TSpec>, FullInfo>& blastHsp)
349 {
350 SEQAN_CHECKPOINT
351
352 blastHsp.identity = 0;
353 blastHsp.score = 0;
354 blastHsp.bits = 0.0;
355 blastHsp.expect = 0.0;
356 blastHsp.gaps = 0;
357 blastHsp.abs_gaps = 0;
358 blastHsp.query_begin = 0;
359 blastHsp.db_end = 0;
360 blastHsp.query_end = 0;
361 blastHsp.db_begin = 0;
362 blastHsp.query_strand = true;
363 blastHsp.db_strand = true;
364 resize(blastHsp.query_string,0);
365 resize(blastHsp.db_string,0);
366 }
367
368 //////////////////////////////////////////////////////////////////////////
369
370
371 /////////////////////////////////////////////////////////////////////////////
372 //////////////////////// FullInfo Protein Spec //////////////////////////////
373 /////////////////////////////////////////////////////////////////////////////
374
375 //Protein BlastHsp
376 template<typename TSpec>
377 class BlastHsp<ProteinBlast<TSpec>, FullInfo>
378 {
379
380 public:
381 float score;
382 float bits;
383 double expect;
384 unsigned int identity;
385 unsigned int positives;
386 unsigned int gaps;
387 unsigned int abs_gaps;
388 bool query_strand; //false=>minus, true=>plus
389 bool db_strand;
390 unsigned int query_begin;
391 unsigned int db_begin;
392 unsigned int query_end; //included in the alignment
393 unsigned int db_end; //included in the alignment
394 int query_frame;
395 int db_frame;
396 String<char> query_string;
397 String<char> db_string;
398
399
BlastHsp()400 BlastHsp()
401 {
402 SEQAN_CHECKPOINT
403 score = 0.0;
404 bits = 0.0 ;
405 expect = 0.0;
406 identity = 0;
407 positives = 0;
408 gaps = 0;
409 abs_gaps = 0;
410 query_strand = true; //false=>minus, true=>plus
411 db_strand = true;
412 query_frame = 1;
413 db_frame = 1;
414 }
415
BlastHsp(BlastHsp const & other)416 BlastHsp(BlastHsp const& other)
417 {
418 SEQAN_CHECKPOINT
419 score = other.score;
420 bits = other.bits;
421 expect = other.expect;
422 identity = other.identity;
423 gaps = other.gaps;
424 abs_gaps = other.abs_gaps;
425 query_strand = other.query_strand;
426 db_strand = other.db_strand;
427 query_begin = other.query_begin;
428 db_begin = other.db_begin;
429 query_end = other.query_end;
430 db_end = other.db_end;
431 query_string = other.query_string;
432 db_string = other.db_string;
433 positives = other.positives;
434 query_frame = other.query_frame;
435 db_frame = other.db_frame;
436 }
437
438 BlastHsp & operator = (BlastHsp const & other)
439 {
440 SEQAN_CHECKPOINT
441 score = other.score;
442 bits = other.bits;
443 expect = other.expect;
444 identity = other.identity;
445 gaps = other.gaps;
446 abs_gaps = other.abs_gaps;
447 query_strand = other.query_strand;
448 db_strand = other.db_strand;
449 query_begin = other.query_begin;
450 db_begin = other.db_begin;
451 query_end = other.query_end;
452 db_end = other.db_end;
453 query_string = other.query_string;
454 db_string = other.db_string;
455 positives = other.positives;
456 query_frame = other.query_frame;
457 db_frame = other.db_frame;
458 return *this;
459 }
460
~BlastHsp()461 ~BlastHsp()
462 {
463 SEQAN_CHECKPOINT
464 }
465 };
466
467
468
469 //parse BlastHsp
470 template<typename TFile, typename TChar, typename TSpec>
471 inline typename Position<TFile>::Type
_parseBlastHsp(TFile & file,TChar & c,BlastHsp<ProteinBlast<TSpec>,FullInfo> & hsp)472 _parseBlastHsp(TFile & file,
473 TChar & c,
474 BlastHsp<ProteinBlast<TSpec>,FullInfo> & hsp)
475 {
476 SEQAN_CHECKPOINT
477 typedef typename Position<TFile>::Type TPosition;
478
479 clear(hsp);
480 String<char> pword;
481 int pint;
482 float pfloat;
483 double pdouble;
484 TPosition start_pos,act_pos;
485 act_pos = _streamTellG(file);
486
487 //String<char> query = "Score";
488 //if(_parseUntilBeginLine(file,c,query,6))
489 //{
490 String<char> query_string;
491 String<char> db_string;
492
493 _parseSkipWhitespace(file,c);
494 if(c == '=')
495 c = _streamGet(file);
496 _parseSkipWhitespace(file,c);
497 pfloat = (float)_parseReadEValue(file, c);
498 hsp.bits = pfloat;
499 if(_parseLineUntil(file,c,'('))
500 {
501 c = _streamGet(file);
502 pfloat = _parseReadFloat(file, c);
503 hsp.score = pfloat;
504 }
505
506 _parseSkipWhitespace(file,c);
507 String<char> query = "Expect";
508 if(_parseLineUntil(file,c,query,6))
509 {
510 c = _streamGet(file);
511 _parseSkipWhitespace(file,c);
512 if(c == '=')
513 c = _streamGet(file);
514 _parseSkipWhitespace(file,c);
515 pdouble = _parseReadEValue(file, c);
516 hsp.expect = pdouble;
517 }
518
519 query = "Identities";
520 if(_parseUntilBeginLine(file,c,query,10,3))
521 {
522 _parseLineUntil(file,c,'(');
523 c = _streamGet(file);
524 pint = _parseReadNumber(file, c);
525 hsp.identity = pint;
526 }
527
528 query = "Positives";
529 if(_parseLineUntil(file,c,query,9))
530 {
531 _parseUntil(file,c,'(');
532 c = _streamGet(file);
533 pint = _parseReadNumber(file, c);
534 hsp.positives = pint;
535 }
536
537 query = "Gaps";
538 if(_parseLineUntil(file,c,query,4))
539 {
540 _parseLineUntil(file,c,'=');
541 c = _streamGet(file);
542 pint = _parseReadNumber(file, c);
543 hsp.abs_gaps = pint;
544 _parseLineUntil(file,c,'(');
545 c = _streamGet(file);
546 pint = _parseReadNumber(file, c);
547 hsp.gaps = pint;
548 }
549
550 query = "Strand";
551 if(_parseUntilBeginLine(file,c,query,6,3))
552 {
553 c = _streamGet(file);
554 _parseSkipWhitespace(file,c);
555 if(c == '=')
556 c = _streamGet(file);
557 _parseSkipWhitespace(file,c);
558 pword = _parseReadWord(file, c);
559 if(pword == "Plus")
560 hsp.query_strand = true;
561 else
562 hsp.query_strand = false;
563 c = _streamGet(file);
564 _parseSkipWhitespace(file,c);
565 if(c == '/')
566 c = _streamGet(file);
567 _parseSkipWhitespace(file,c);
568 pword = _parseReadWord(file, c);
569 if(pword == "Plus")
570 hsp.db_strand = true;
571 else
572 hsp.db_strand = false;
573 }
574
575 bool first = true;
576 query = "Query";
577 if(_parseUntilBeginLine(file,c,query,5))
578 {
579 bool in_hsp = true;
580 while(in_hsp)
581 {
582 int end_query,end_db = -1;
583 _parseSkipWhitespace(file,c);
584 if(c == ':')
585 c = _streamGet(file);
586 _parseSkipWhitespace(file,c);
587 pint = _parseReadNumber(file,c);
588 _parseSkipWhitespace(file,c);
589 //TPosition act_posQ = _streamTellG(file);
590 if(first)
591 hsp.query_begin = pint;
592 query_string += _parseReadAlignmentString(file,c);
593 _parseSkipWhitespace(file,c);
594 end_query = _parseReadNumber(file,c);
595 query = "Sbjct";
596 _parseUntilBeginLine(file,c,query,5);
597 _parseSkipWhitespace(file,c);
598 if(c == ':')
599 c = _streamGet(file);
600 _parseSkipWhitespace(file,c);
601 pint = _parseReadNumber(file,c);
602 _parseSkipWhitespace(file,c);
603 //TPosition act_posS = _streamTellG(file);
604 if(first)
605 hsp.db_begin = pint;
606 db_string += _parseReadAlignmentString(file,c);
607 _parseSkipWhitespace(file,c);
608 end_db = _parseReadNumber(file,c);
609 first = false;
610 String<TChar> delim = "QS>R";
611 if(_parseUntilBeginLineOneOf(file,c,delim,4))
612 {
613 TPosition pos = _streamTellG(file);
614 pword = _parseReadWord(file,c);
615 // _streamSeekG(file,pos);
616
617 //still in the same HSP
618 if(pword == "Query")
619 in_hsp = true;
620 else
621 {
622 hsp.query_end = end_query;
623 hsp.db_end = end_db;
624 hsp.query_string = query_string;
625 hsp.db_string = db_string;
626 //next HSP
627 if(pword == "Score")
628 return _streamTellG(file);
629 if(pword[0] == '>')
630 return act_pos;
631 if(pword[0] == 'R')
632 return (TPosition) 0;
633 in_hsp = false;
634 }
635 }
636 else
637 return act_pos;
638 }
639 }
640 return act_pos;
641
642 }
643
644
645 template<typename TSpec>
646 inline void
clear(BlastHsp<ProteinBlast<TSpec>,FullInfo> & blastHsp)647 clear(BlastHsp<ProteinBlast<TSpec>, FullInfo>& blastHsp)
648 {
649 SEQAN_CHECKPOINT
650 blastHsp.score = 0;
651 blastHsp.bits = 0.0;
652 blastHsp.expect = 0.0;
653 blastHsp.abs_gaps = 0;
654 blastHsp.query_begin = 0;
655 blastHsp.db_end = 0;
656 blastHsp.query_end = 0;
657 blastHsp.db_begin = 0;
658 blastHsp.positives = 0;
659 blastHsp.query_frame = 1;
660 blastHsp.db_frame = 1;
661 blastHsp.identity = 0;
662 blastHsp.gaps = 0;
663 blastHsp.query_strand = true;
664 blastHsp.db_strand = true;
665 resize(blastHsp.query_string,0);
666 resize(blastHsp.db_string,0);
667 }
668
669
670
671
672
673
674 /////////////////////////////////////////////////////////////////////////////
675 /////////////////////////// BasicInfo Spec //////////////////////////////////
676 /////////////////////////////////////////////////////////////////////////////
677
678
679
680 template<typename TBlastSpec>
681 class BlastHsp<TBlastSpec, BasicInfo>
682 {
683
684 public:
685 double expect;
686 unsigned int query_begin;
687 unsigned int db_begin;
688 unsigned int query_end; //included in the alignment
689 unsigned int db_end; //included in the alignment
690
691 String<char> query_string;
692 String<char> db_string;
693
694
BlastHsp()695 BlastHsp()
696 {
697 SEQAN_CHECKPOINT
698 expect = 0;
699 }
700
BlastHsp(BlastHsp const & other)701 BlastHsp(BlastHsp const& other)
702 {
703 SEQAN_CHECKPOINT
704 expect = other.expect;
705 query_begin = other.query_begin;
706 db_begin = other.db_begin;
707 query_end = other.query_end;
708 db_end = other.db_end;
709 query_string = other.query_string;
710 db_string = other.db_string;
711 }
712
713 BlastHsp & operator = (BlastHsp const & other)
714 {
715 SEQAN_CHECKPOINT
716 expect = other.expect;
717 query_begin = other.query_begin;
718 db_begin = other.db_begin;
719 query_end = other.query_end;
720 db_end = other.db_end;
721 query_string = other.query_string;
722 db_string = other.db_string;
723 return *this;
724 }
725
~BlastHsp()726 ~BlastHsp()
727 {
728 SEQAN_CHECKPOINT
729 }
730 };
731
732
733 template<typename TBlastSpec>
734 inline void
clear(BlastHsp<TBlastSpec,BasicInfo> & blastHsp)735 clear(BlastHsp<TBlastSpec, BasicInfo>& blastHsp)
736 {
737 SEQAN_CHECKPOINT
738
739 blastHsp.expect = 0.0;
740 blastHsp.query_begin = 0;
741 blastHsp.db_end = 0;
742 blastHsp.query_end = 0;
743 blastHsp.db_begin = 0;
744 resize(blastHsp.query_string,0);
745 resize(blastHsp.db_string,0);
746 }
747
748
749
750 //parse BlastHsp
751 template<typename TFile, typename TChar, typename TBlastSpec, typename TInfoSpec>
752 inline typename Position<TFile>::Type
_parseBlastHsp(TFile & file,TChar & c,BlastHsp<TBlastSpec,TInfoSpec> & hsp)753 _parseBlastHsp(TFile & file,
754 TChar & c,
755 BlastHsp<TBlastSpec,TInfoSpec> & hsp)
756 {
757 SEQAN_CHECKPOINT
758 typedef typename Position<TFile>::Type TPosition;
759
760 clear(hsp);
761
762 String<char> pword;
763 int pint;
764 double pdouble;
765 TPosition start_pos,act_pos;
766 act_pos = _streamTellG(file);
767
768 String<char> query_string;
769 String<char> db_string;
770
771 _parseSkipWhitespace(file,c);
772 String<char> query = "Expect";
773 if(_parseLineUntil(file,c,query,6))
774 {
775 c = _streamGet(file);
776 _parseSkipWhitespace(file,c);
777 c = _streamGet(file); // =
778 _parseSkipWhitespace(file,c);
779 pdouble = _parseReadEValue(file, c);
780 hsp.expect = pdouble;
781 }
782
783 bool first = true;
784 query = "Query";
785 if(_parseUntilBeginLine(file,c,query,5))
786 {
787 bool in_hsp = true;
788 while(in_hsp)
789 {
790 int end_query,end_db = -1;
791 _parseSkipWhitespace(file,c);
792 c = _streamGet(file);
793 _parseSkipWhitespace(file,c);
794 pint = _parseReadNumber(file,c);
795 _parseSkipWhitespace(file,c);
796 if(first)
797 hsp.query_begin = pint;
798 query_string += _parseReadAlignmentString(file,c);
799 _parseSkipWhitespace(file,c);
800 end_query = _parseReadNumber(file,c);
801 query = "Sbjct";
802 _parseUntilBeginLine(file,c,query,5);
803 _parseSkipWhitespace(file,c);
804 c = _streamGet(file);
805 _parseSkipWhitespace(file,c);
806 pint = _parseReadNumber(file,c);
807 _parseSkipWhitespace(file,c);
808 if(first)
809 hsp.db_begin = pint;
810 db_string += _parseReadAlignmentString(file,c);
811 _parseSkipWhitespace(file,c);
812 end_db = _parseReadNumber(file,c);
813 first = false;
814 String<TChar> delim = "QS>R";
815 if(_parseUntilBeginLineOneOf(file,c,delim,4))
816 {
817 TPosition pos = _streamTellG(file);
818 pword = _parseReadWord(file,c);
819
820 //still in the same HSP
821 if(pword == "Query")
822 in_hsp = true;
823 else
824 {
825 hsp.query_end = end_query;
826 hsp.db_end = end_db;
827 hsp.query_string = query_string;
828 hsp.db_string = db_string;
829 //next HSP
830 if(pword == "Score")
831 return _streamTellG(file);
832 if(pword[0] == '>')
833 return act_pos;
834 if(pword[0] == 'R')
835 return (TPosition) 0;
836 in_hsp = false;
837 }
838 }
839 else
840 return act_pos;
841 }
842 }
843 return act_pos;
844
845 }
846
847
848
849
850 /////////////////////////////////////////////////////////////
851 /////////////////////// get functions ///////////////////////
852 /////////////////////////////////////////////////////////////
853
854 struct TagUnknownSource_;
855 typedef Tag<TagUnknownSource_> const UnknownSource;
856
857 struct TagKnownSource_;
858 typedef Tag<TagKnownSource_> const KnownSource;
859
860 /////////////////////////////////////////////////////////////
861 // get Alignment for Align<TSource,TSpec>
862 template<typename TBlastHsp, typename TSpec, typename TSource>
863 inline unsigned int
getAlignment(TBlastHsp & hsp,Align<TSource,TSpec> & ali,UnknownSource)864 getAlignment(TBlastHsp & hsp,
865 Align<TSource,TSpec> & ali,
866 UnknownSource)
867 {
868 SEQAN_CHECKPOINT
869
870 typedef Align<TSource, TSpec> TAlign;
871 typedef typename Row<TAlign>::Type TRow;
872
873 typedef Pair<unsigned int, unsigned int> TPosLen;
874
875 //for storing node info
876 String<TPosLen> seq0info;
877 String<TPosLen> seq1info;
878
879 TSource str0,str1;
880
881 typedef typename Iterator<String<char>, Standard >::Type TStringIterator;
882 TStringIterator it_0 = begin(queryAlignmentString(hsp),Standard());
883 TStringIterator it_0_end = end(queryAlignmentString(hsp),Standard());
884 TStringIterator it_1 = begin(databaseAlignmentString(hsp),Standard());
885 TStringIterator it_1_end = end(databaseAlignmentString(hsp),Standard());
886
887 unsigned int act0_pos = 0;
888 unsigned int act1_pos = 0;
889
890 while ((it_0 != it_0_end) && (it_1 != it_1_end))
891 {
892 if(*it_0 == '-')
893 {
894 unsigned int gap_len = 0;
895 while(it_0 != it_0_end && *it_0 == '-')
896 {
897 append(str1,getValue(it_1));
898 ++it_0;
899 ++it_1;
900 ++act1_pos;
901 ++gap_len;
902 }
903 append(seq0info,TPosLen(act0_pos,gap_len));
904 act0_pos += gap_len;
905 continue;
906 }
907 if(*it_1 == '-')
908 {
909 unsigned int gap_len = 0;
910 while(it_1 != it_1_end && *it_1 == '-')
911 {
912 append(str0,getValue(it_0));
913 ++it_1;
914 ++it_0;
915 ++act0_pos;
916 ++gap_len;
917 }
918 append(seq1info,TPosLen(act1_pos,gap_len));
919 act1_pos += gap_len;
920 continue;
921 }
922 append(str0,getValue(it_0));
923 append(str1,getValue(it_1));
924 ++it_0;
925 ++it_1;
926 ++act0_pos;
927 ++act1_pos;
928
929 }
930
931 resize(rows(ali),2);
932 setSource(row(ali,0),str0);
933 setSource(row(ali,1),str1);
934
935 typedef typename Iterator<String<TPosLen>, Standard >::Type TNodeIterator;
936 TNodeIterator nodes0it = begin(seq0info,Standard());
937 TNodeIterator nodes0end = end(seq0info,Standard());
938
939 while(nodes0it != nodes0end)
940 {
941 insertGaps(row(ali,0),nodes0it->i1,nodes0it->i2);
942 ++nodes0it;
943 }
944
945
946 typedef typename Iterator<String<TPosLen>, Standard >::Type TNodeIterator;
947 TNodeIterator nodes1it = begin(seq1info,Standard());
948 TNodeIterator nodes1end = end(seq1info,Standard());
949
950 while(nodes1it != nodes1end)
951 {
952 insertGaps(row(ali,1),nodes1it->i1,nodes1it->i2);
953 ++nodes1it;
954 }
955 detach(ali);
956
957 return length(queryAlignmentString(hsp));
958
959 }
960
961 /**
962 .Function.blast#getAlignment
963 ..cat:Blast
964 ..summary:Turns a HSP from a Blast search into an Alignment object.
965 ..signature:getAlignment(hsp,alignment)
966 ..param.hsp:A Blast HSP.
967 ...type:Class.BlastHsp
968 ..param.alignment:An Alignment object to be filled.
969 ...type:Class.Align
970 ...type:Spec.Alignment Graph
971 ..include:seqan/blast.h
972 */
973 template<typename TBlastHsp, typename TSpec, typename TSource>
974 inline unsigned int
getAlignment(TBlastHsp & hsp,Align<TSource,TSpec> & ali)975 getAlignment(TBlastHsp & hsp,
976 Align<TSource,TSpec> & ali)
977 {
978 SEQAN_CHECKPOINT
979
980 if(length(rows(ali))>1)
981 if(!emptySource(row(ali,0)) && !!emptySource(row(ali,1)))
982 return getAlignment(hsp,ali,KnownSource());
983
984 return getAlignment(hsp,ali,UnknownSource());
985
986 }
987
988
989
990 /////////////////////////////////////////////////////////////
991 // get Alignment for Align<TSource,TSpec>
992 template<typename TBlastHsp, typename TSpec, typename TSource>
993 inline unsigned int
getAlignment(TBlastHsp & hsp,Align<TSource,TSpec> & ali,KnownSource)994 getAlignment(TBlastHsp & hsp,
995 Align<TSource,TSpec> & ali,
996 KnownSource)
997 {
998 SEQAN_CHECKPOINT
999
1000 typedef Align<TSource, TSpec> TAlign;
1001 typedef typename Row<TAlign>::Type TRow;
1002
1003 typedef typename Iterator<TRow, Standard>::Type TAliIterator;
1004 TAliIterator ali_0 = begin(row(ali, 0), Standard());
1005 TAliIterator ali_1 = begin(row(ali, 1), Standard());
1006
1007 typedef typename Iterator<String<char>, Standard >::Type TStringIterator;
1008 TStringIterator it_0 = begin(queryAlignmentString(hsp),Standard());
1009 TStringIterator it_0_end = end(queryAlignmentString(hsp),Standard());
1010
1011 TStringIterator it_1 = begin(databaseAlignmentString(hsp),Standard());
1012 TStringIterator it_1_end = end(databaseAlignmentString(hsp),Standard());
1013
1014 TSource str0,str1;
1015
1016 while ((it_0 != it_0_end) && (it_1 != it_1_end))
1017 {
1018 if (*it_0 == '-')
1019 {
1020 insertGap(ali_0);
1021 }
1022 else
1023 {
1024 append(str0,getValue(it_0));
1025 }
1026 if (*it_1 == '-')
1027 {
1028 insertGap(ali_1);
1029 }
1030 else
1031 {
1032 append(str1,getValue(it_1));
1033 }
1034
1035 ++it_0;
1036 ++it_1;
1037 ++ali_0;
1038 ++ali_1;
1039 }
1040
1041 setClippedBeginPosition(row(ali,0),queryBegin(hsp)-1);
1042 setClippedBeginPosition(row(ali,1),databaseBegin(hsp)-1);
1043 setClippedEndPosition(row(ali,0),queryEnd(hsp));
1044 setClippedEndPosition(row(ali,1),databaseEnd(hsp));
1045
1046 return length(databaseAlignmentString(hsp));
1047
1048 }
1049
1050
1051
1052
1053
1054 /**
1055 .Function.blast#getAlignment
1056 ..cat:Blast
1057 ..signature:getAlignment(hsp,alignment,id0,id1)
1058 ..param.id0:The Id of the query sequence in the StringSet of the Alignment Graph.
1059 ..param.id1:The Id of the hit sequence in the StringSet of the Alignment Graph.
1060 ..include:seqan/blast.h
1061 */
1062 /////////////////////////////////////////////////////////////
1063 // get Alignment for Graph<TAlign>
1064 template<typename TBlastHsp, typename TAlign, typename TId>
1065 inline unsigned int
getAlignment(TBlastHsp & hsp,Graph<TAlign> & ali,TId id0,TId id1)1066 getAlignment(TBlastHsp & hsp,
1067 Graph<TAlign> & ali,
1068 TId id0, //query ID
1069 TId id1) //hit ID
1070 {
1071 SEQAN_CHECKPOINT
1072
1073 typedef Graph<TAlign> TAliGraph;
1074 typedef typename Iterator<String<char>, Standard >::Type TStringIterator;
1075 typedef typename VertexDescriptor<TAliGraph>::Type TVertexDescriptor;
1076
1077 TStringIterator it_0 = begin(queryAlignmentString(hsp),Standard());
1078 TStringIterator it_0_end = end(queryAlignmentString(hsp),Standard());
1079
1080 TStringIterator it_1 = begin(databaseAlignmentString(hsp),Standard());
1081 TStringIterator it_1_end = end(databaseAlignmentString(hsp),Standard());
1082 unsigned int act0_pos = getQueryBegin(hsp)-1;
1083 unsigned int act1_pos = getDatabaseBegin(hsp)-1;
1084
1085 unsigned int act_len = 0;
1086 while ((it_0 != it_0_end) && (it_1 != it_1_end))
1087 {
1088 if(*it_0 == '-')
1089 {
1090 if(act_len > 0)
1091 {
1092 TVertexDescriptor vd0 = addVertex(ali,id0,act0_pos,act_len);
1093 TVertexDescriptor vd1 = addVertex(ali,id1,act1_pos,act_len);
1094 addEdge(ali,vd0,vd1);
1095 act0_pos += act_len;
1096 act1_pos += act_len;
1097 act_len = 0;
1098 }
1099 unsigned int gap_len = 0;
1100 while(it_0 != it_0_end && *it_0 == '-')
1101 {
1102 ++it_0;
1103 ++it_1;
1104 ++gap_len;
1105 }
1106 addVertex(ali,id1,act1_pos,gap_len);
1107 act1_pos += gap_len;
1108 continue;
1109 }
1110 if(*it_1 == '-')
1111 {
1112 if(act_len > 0)
1113 {
1114 TVertexDescriptor vd0 = addVertex(ali,id0,act0_pos,act_len);
1115 TVertexDescriptor vd1 = addVertex(ali,id1,act1_pos,act_len);
1116 addEdge(ali,vd0,vd1);
1117 act0_pos += act_len;
1118 act1_pos += act_len;
1119 act_len = 0;
1120 }
1121 unsigned int gap_len = 0;
1122 while(it_1 != it_1_end && *it_1 == '-')
1123 {
1124 ++it_1;
1125 ++it_0;
1126 ++gap_len;
1127 }
1128 addVertex(ali,id0,act0_pos,gap_len);
1129 act0_pos += gap_len;
1130 continue;
1131 }
1132 ++it_0;
1133 ++it_1;
1134 ++act_len;
1135
1136 }
1137 SEQAN_ASSERT_TRUE(act0_pos+act_len == getQueryEnd(hsp));
1138 SEQAN_ASSERT_TRUE(act1_pos+act_len == getDatabaseEnd(hsp));
1139 if(act_len>0)
1140 {
1141 TVertexDescriptor vd0 = addVertex(ali,id0,act0_pos,act_len);
1142 TVertexDescriptor vd1 = addVertex(ali,id1,act1_pos,act_len);
1143 addEdge(ali,vd0,vd1);
1144 }
1145
1146 return length(queryAlignmentString(hsp));
1147
1148 }
1149
1150
1151
1152 /////////////////////////////////////////////////////////////
1153 // get Alignment for Graph<TAlign> (stringSet and sequence ids not given)
1154 template<typename TBlastHsp, typename TStringSet, typename TCargo, typename TSpec>
1155 inline unsigned int //returns the length of the alignment
getAlignment(TBlastHsp & hsp,Graph<Alignment<TStringSet,TCargo,TSpec>> & ali)1156 getAlignment(TBlastHsp & hsp,
1157 Graph<Alignment<TStringSet,TCargo,TSpec> > & ali) //hit ID
1158 {
1159 SEQAN_CHECKPOINT
1160
1161 typedef Graph<Alignment<TStringSet,TCargo,TSpec> > TAliGraph;
1162 typedef typename Iterator<String<char>, Standard >::Type TStringIterator;
1163 typedef typename VertexDescriptor<TAliGraph>::Type TVertexDescriptor;
1164 typedef Pair<unsigned int, unsigned int> TPosLen;
1165 typedef Triple<unsigned int, unsigned int, int> TPosLenLink;
1166 typedef typename Value<TStringSet>::Type TString;
1167
1168 //for storing node info
1169 String<TPosLen> seq0nodes;
1170 //for storing node + edge info
1171 String<TPosLenLink> seq1nodes;
1172
1173 //for storing the sequences
1174 TString str0,str1;
1175
1176 TStringIterator it_0 = begin(queryAlignmentString(hsp),Standard());
1177 TStringIterator it_0_end = end(queryAlignmentString(hsp),Standard());
1178
1179 TStringIterator it_1 = begin(databaseAlignmentString(hsp),Standard());
1180 TStringIterator it_1_end = end(databaseAlignmentString(hsp),Standard());
1181 unsigned int act0_pos = 0;
1182 unsigned int act1_pos = 0;
1183 unsigned int act_len = 0;
1184
1185 //walk over (gapped) sequences, retrieve ungapped sequences and store node+edge info for alignment graph
1186 while ((it_0 != it_0_end) && (it_1 != it_1_end))
1187 {
1188 if(*it_0 == '-')
1189 {
1190 if(act_len > 0)
1191 {
1192 append(seq0nodes,TPosLen(act0_pos,act_len));
1193 append(seq1nodes,TPosLenLink(act1_pos,act_len,act0_pos));
1194 act0_pos += act_len;
1195 act1_pos += act_len;
1196 act_len = 0;
1197 }
1198 unsigned int gap_len = 0;
1199 while(it_0 != it_0_end && *it_0 == '-')
1200 {
1201 append(str1,getValue(it_1));
1202 ++it_0;
1203 ++it_1;
1204 ++gap_len;
1205 }
1206 append(seq1nodes,TPosLenLink(act1_pos,gap_len,-1));
1207 act1_pos += gap_len;
1208 continue;
1209 }
1210 if(*it_1 == '-')
1211 {
1212 if(act_len > 0)
1213 {
1214 append(seq0nodes,TPosLen(act0_pos,act_len));
1215 append(seq1nodes,TPosLenLink(act1_pos,act_len,act0_pos));
1216 act0_pos += act_len;
1217 act1_pos += act_len;
1218 act_len = 0;
1219 }
1220 unsigned int gap_len = 0;
1221 while(it_1 != it_1_end && *it_1 == '-')
1222 {
1223 append(str0,getValue(it_0));
1224 ++it_1;
1225 ++it_0;
1226 ++gap_len;
1227 }
1228 append(seq0nodes,TPosLen(act0_pos,gap_len));
1229 act0_pos += gap_len;
1230 continue;
1231 }
1232 append(str0,getValue(it_0));
1233 append(str1,getValue(it_1));
1234 ++it_0;
1235 ++it_1;
1236 ++act_len;
1237
1238 }
1239
1240 if(act_len>0)
1241 {
1242 append(seq0nodes,TPosLen(act0_pos,act_len));
1243 append(seq1nodes,TPosLenLink(act1_pos,act_len,act0_pos));
1244 }
1245
1246 TStringSet str;
1247 unsigned int id0 = assignValueById(str,str0);
1248 unsigned int id1 = assignValueById(str,str1);
1249
1250 assignStringSet(ali,str);
1251
1252 typedef typename Iterator<String<TPosLen>, Standard >::Type TNodeIterator;
1253 TNodeIterator nodes0it = begin(seq0nodes,Standard());
1254 TNodeIterator nodes0end = end(seq0nodes,Standard());
1255
1256 while(nodes0it != nodes0end)
1257 {
1258 addVertex(ali,id0,nodes0it->i1,nodes0it->i2);
1259 ++nodes0it;
1260 }
1261
1262 typedef typename Iterator<String<TPosLenLink>, Standard >::Type TNodeEdgeIterator;
1263 TNodeEdgeIterator nodes1it = begin(seq1nodes,Standard());
1264 TNodeEdgeIterator nodes1end = end(seq1nodes,Standard());
1265
1266 while(nodes1it != nodes1end)
1267 {
1268 TVertexDescriptor vd1 = addVertex(ali,id1,nodes1it->i1,nodes1it->i2);
1269 if(nodes1it->i3 >= 0)
1270 {
1271 TVertexDescriptor vd2 = findVertex(ali,id0,nodes1it->i3);
1272 addEdge(ali,vd1,vd2);
1273 }
1274 ++nodes1it;
1275 }
1276
1277 return length(databaseAlignmentString(hsp));
1278
1279 }
1280
1281
1282
1283 // general (for all specs)
1284
1285 template<typename TBlastSpec, typename TInfoSpec>
1286 inline unsigned int &
queryBegin(BlastHsp<TBlastSpec,TInfoSpec> & blastHsp)1287 queryBegin(BlastHsp<TBlastSpec, TInfoSpec>& blastHsp)
1288 {
1289 SEQAN_CHECKPOINT
1290 return blastHsp.query_begin;
1291 }
1292
1293
1294
1295
1296 /**
1297 .Function.getQueryBegin:
1298 ..cat:Blast
1299 ..summary:The begin position of the HSP on the query sequence.
1300 ..signature:getQueryBegin(object);
1301 ..param.object:A Blast HSP object.
1302 ...type:Class.BlastHsp
1303 ..returns:The begin position.
1304 ...type:nolink:unsigned
1305 ..include:seqan/blast.h
1306 */
1307 template<typename TBlastSpec, typename TInfoSpec>
1308 inline unsigned int
getQueryBegin(BlastHsp<TBlastSpec,TInfoSpec> & blastHsp)1309 getQueryBegin(BlastHsp<TBlastSpec, TInfoSpec>& blastHsp)
1310 {
1311 SEQAN_CHECKPOINT
1312 return blastHsp.query_begin;
1313 }
1314
1315
1316
1317 template<typename TBlastSpec, typename TInfoSpec>
1318 inline unsigned int &
databaseBegin(BlastHsp<TBlastSpec,TInfoSpec> & blastHsp)1319 databaseBegin(BlastHsp<TBlastSpec, TInfoSpec>& blastHsp)
1320 {
1321 SEQAN_CHECKPOINT
1322 return blastHsp.db_begin;
1323 }
1324
1325 /**
1326 .Function.getDatabaseBegin:
1327 ..cat:Blast
1328 ..summary:The begin position of the HSP on the database sequence.
1329 ..signature:getDatabaseBegin(object);
1330 ..param.object:A Blast HSP object.
1331 ...type:Class.BlastHsp
1332 ..returns:The begin position.
1333 ...type:nolink:unsigned
1334 ..include:seqan/blast.h
1335 */
1336 template<typename TBlastSpec, typename TInfoSpec>
1337 inline unsigned int
getDatabaseBegin(BlastHsp<TBlastSpec,TInfoSpec> & blastHsp)1338 getDatabaseBegin(BlastHsp<TBlastSpec, TInfoSpec>& blastHsp)
1339 {
1340 SEQAN_CHECKPOINT
1341 return blastHsp.db_begin;
1342 }
1343
1344
1345
1346 template<typename TBlastSpec, typename TInfoSpec>
1347 inline unsigned int &
queryEnd(BlastHsp<TBlastSpec,TInfoSpec> & blastHsp)1348 queryEnd(BlastHsp<TBlastSpec, TInfoSpec>& blastHsp)
1349 {
1350 SEQAN_CHECKPOINT
1351 return blastHsp.query_end;
1352 }
1353
1354 /**
1355 .Function.getQueryEnd:
1356 ..cat:Blast
1357 ..summary:The end position of the HSP on the query sequence.
1358 ..signature:getQueryEnd(object);
1359 ..param.object:A Blast HSP object.
1360 ...type:Class.BlastHsp
1361 ..returns:The end position.
1362 ...type:nolink:unsigned
1363 ..include:seqan/blast.h
1364 */
1365 template<typename TBlastSpec, typename TInfoSpec>
1366 inline unsigned int
getQueryEnd(BlastHsp<TBlastSpec,TInfoSpec> & blastHsp)1367 getQueryEnd(BlastHsp<TBlastSpec, TInfoSpec>& blastHsp)
1368 {
1369 SEQAN_CHECKPOINT
1370 return blastHsp.query_end;
1371 }
1372
1373
1374
1375 /**
1376 .Function.getDatabaseEnd:
1377 ..cat:Blast
1378 ..summary:The end position of the HSP on the database sequence.
1379 ..signature:getDatabaseEnd(object);
1380 ..param.object:A Blast HSP object.
1381 ...type:Class.BlastHsp
1382 ..returns:The end position.
1383 ...type:nolink:unsigned
1384 ..include:seqan/blast.h
1385 */
1386 template<typename TBlastSpec, typename TInfoSpec>
1387 inline unsigned int &
databaseEnd(BlastHsp<TBlastSpec,TInfoSpec> & blastHsp)1388 databaseEnd(BlastHsp<TBlastSpec, TInfoSpec>& blastHsp)
1389 {
1390 SEQAN_CHECKPOINT
1391 return blastHsp.db_end;
1392 }
1393
1394 template<typename TBlastSpec, typename TInfoSpec>
1395 inline unsigned int
getDatabaseEnd(BlastHsp<TBlastSpec,TInfoSpec> & blastHsp)1396 getDatabaseEnd(BlastHsp<TBlastSpec, TInfoSpec>& blastHsp)
1397 {
1398 SEQAN_CHECKPOINT
1399 return blastHsp.db_end;
1400 }
1401
1402
1403 template<typename TBlastSpec, typename TInfoSpec>
1404 inline String<char> &
queryAlignmentString(BlastHsp<TBlastSpec,TInfoSpec> & blastHsp)1405 queryAlignmentString(BlastHsp<TBlastSpec, TInfoSpec>& blastHsp)
1406 {
1407 SEQAN_CHECKPOINT
1408 return blastHsp.query_string;
1409 }
1410
1411
1412 template<typename TBlastSpec, typename TInfoSpec>
1413 inline String<char>
getQueryAlignmentString(BlastHsp<TBlastSpec,TInfoSpec> & blastHsp)1414 getQueryAlignmentString(BlastHsp<TBlastSpec, TInfoSpec>& blastHsp)
1415 {
1416 SEQAN_CHECKPOINT
1417 return blastHsp.query_string;
1418 }
1419
1420
1421
1422 template<typename TBlastSpec, typename TInfoSpec>
1423 inline String<char> &
databaseAlignmentString(BlastHsp<TBlastSpec,TInfoSpec> & blastHsp)1424 databaseAlignmentString(BlastHsp<TBlastSpec, TInfoSpec>& blastHsp)
1425 {
1426 SEQAN_CHECKPOINT
1427 return blastHsp.db_string;
1428 }
1429
1430 template<typename TBlastSpec, typename TInfoSpec>
1431 inline String<char>
getDatabaseAlignmentString(BlastHsp<TBlastSpec,TInfoSpec> & blastHsp)1432 getDatabaseAlignmentString(BlastHsp<TBlastSpec, TInfoSpec>& blastHsp)
1433 {
1434 SEQAN_CHECKPOINT
1435 return blastHsp.db_string;
1436 }
1437
1438
1439 template<typename TBlastSpec, typename TInfoSpec>
1440 inline double &
eValue(BlastHsp<TBlastSpec,TInfoSpec> & blastHsp)1441 eValue(BlastHsp<TBlastSpec, TInfoSpec>& blastHsp)
1442 {
1443 SEQAN_CHECKPOINT
1444 return blastHsp.expect;
1445 }
1446
1447 /**
1448 .Function.getEValue:
1449 ..cat:Blast
1450 ..summary:The e-value associated with a Blast HSP.
1451 ..signature:getEValue(object);
1452 ..param.object:A Blast HSP object.
1453 ...type:Class.BlastHsp
1454 ..returns:The e-value.
1455 ...type:nolink:double
1456 ..include:seqan/blast.h
1457 */
1458 template<typename TBlastSpec, typename TInfoSpec>
1459 inline double
getEValue(BlastHsp<TBlastSpec,TInfoSpec> & blastHsp)1460 getEValue(BlastHsp<TBlastSpec, TInfoSpec>& blastHsp)
1461 {
1462 SEQAN_CHECKPOINT
1463 return blastHsp.expect;
1464 }
1465
1466
1467 // for FullInfo specs
1468
1469 template<typename TBlastSpec>
1470 inline float &
score(BlastHsp<TBlastSpec,FullInfo> & blastHsp)1471 score(BlastHsp<TBlastSpec, FullInfo>& blastHsp)
1472 {
1473 SEQAN_CHECKPOINT
1474 return blastHsp.score;
1475 }
1476
1477 /**
1478 .Function.Blast#getBlastMatchScore:
1479 ..cat:Blast
1480 ..summary:The Smith-Waterman score associated with a Blast HSP.
1481 ..signature:getBlastMatchScore(object);
1482 ..param.object:A Blast HSP object.
1483 ...type:Spec.FullInfo
1484 ..returns:The score.
1485 ...type:nolink:float
1486 ..include:seqan/blast.h
1487 */
1488 template<typename TBlastSpec>
1489 inline float
getBlastMatchScore(BlastHsp<TBlastSpec,FullInfo> & blastHsp)1490 getBlastMatchScore(BlastHsp<TBlastSpec, FullInfo>& blastHsp)
1491 {
1492 SEQAN_CHECKPOINT
1493 return blastHsp.score;
1494 }
1495
1496 template<typename TBlastSpec>
1497 inline float &
bitScore(BlastHsp<TBlastSpec,FullInfo> & blastHsp)1498 bitScore(BlastHsp<TBlastSpec, FullInfo>& blastHsp)
1499 {
1500 SEQAN_CHECKPOINT
1501 return blastHsp.bits;
1502 }
1503
1504
1505 /**
1506 .Function.getBitScore:
1507 ..cat:Blast
1508 ..summary:The bit score associated with a Blast HSP.
1509 ..signature:getBitScore(object);
1510 ..param.object:A Blast HSP object.
1511 ...type:Spec.FullInfo
1512 ..returns:The bit score.
1513 ...type:nolink:float
1514 ..include:seqan/blast.h
1515 */
1516 template<typename TBlastSpec>
1517 inline float
getBitScore(BlastHsp<TBlastSpec,FullInfo> & blastHsp)1518 getBitScore(BlastHsp<TBlastSpec, FullInfo>& blastHsp)
1519 {
1520 SEQAN_CHECKPOINT
1521 return blastHsp.bits;
1522 }
1523
1524
1525 template<typename TBlastSpec>
1526 inline unsigned int &
percentIdentity(BlastHsp<TBlastSpec,FullInfo> & blastHsp)1527 percentIdentity(BlastHsp<TBlastSpec, FullInfo>& blastHsp)
1528 {
1529 SEQAN_CHECKPOINT
1530 return blastHsp.identity;
1531 }
1532
1533
1534
1535 template<typename TBlastSpec>
1536 inline unsigned int
getPercentIdentity(BlastHsp<TBlastSpec,FullInfo> & blastHsp)1537 getPercentIdentity(BlastHsp<TBlastSpec, FullInfo>& blastHsp)
1538 {
1539 SEQAN_CHECKPOINT
1540 return blastHsp.identity;
1541 }
1542
1543 template<typename TBlastSpec>
1544 inline unsigned int &
percentGaps(BlastHsp<TBlastSpec,FullInfo> & blastHsp)1545 percentGaps(BlastHsp<TBlastSpec, FullInfo>& blastHsp)
1546 {
1547 SEQAN_CHECKPOINT
1548 return blastHsp.gaps;
1549 }
1550
1551 template<typename TBlastSpec>
1552 inline unsigned int
getPercentGaps(BlastHsp<TBlastSpec,FullInfo> & blastHsp)1553 getPercentGaps(BlastHsp<TBlastSpec, FullInfo>& blastHsp)
1554 {
1555 SEQAN_CHECKPOINT
1556 return blastHsp.gaps;
1557 }
1558
1559 template<typename TBlastSpec>
1560 inline unsigned int &
numGaps(BlastHsp<TBlastSpec,FullInfo> & blastHsp)1561 numGaps(BlastHsp<TBlastSpec, FullInfo>& blastHsp)
1562 {
1563 SEQAN_CHECKPOINT
1564 return blastHsp.abs_gaps;
1565 }
1566
1567
1568
1569 /**
1570 .Function.getNumGaps:
1571 ..cat:Blast
1572 ..summary:The number of gaps within a Blast HSP alignment.
1573 ..signature:getNumGaps(object);
1574 ..param.object:A Blast HSP object.
1575 ...type:Spec.FullInfo
1576 ..returns:The number of gaps.
1577 ...type:nolink:unsigned
1578 ..include:seqan/blast.h
1579 */
1580 template<typename TBlastSpec>
1581 inline unsigned int
getNumGaps(BlastHsp<TBlastSpec,FullInfo> & blastHsp)1582 getNumGaps(BlastHsp<TBlastSpec, FullInfo>& blastHsp)
1583 {
1584 SEQAN_CHECKPOINT
1585 return blastHsp.abs_gaps;
1586 }
1587
1588
1589
1590 /**
1591 .Function.queryOrientationPlus:
1592 ..cat:Blast
1593 ..summary:Orientation of the query sequence within a Blast HSP alignment.
1594 ..signature:queryOrientationPlus(object);
1595 ..param.object:A Blast HSP object.
1596 ...type:Spec.FullInfo
1597 ..returns:True if the query is in forward orientation.
1598 ...type:nolink:bool
1599 ..include:seqan/blast.h
1600 */
1601 template<typename TBlastSpec>
1602 inline bool
queryOrientationPlus(BlastHsp<TBlastSpec,FullInfo> & blastHsp)1603 queryOrientationPlus(BlastHsp<TBlastSpec, FullInfo>& blastHsp)
1604 {
1605 SEQAN_CHECKPOINT
1606 return blastHsp.query_strand;
1607 }
1608
1609 /**
1610 .Function.databaseOrientationPlus:
1611 ..cat:Blast
1612 ..summary:Orientation of the database sequence within a Blast HSP alignment.
1613 ..signature:databaseOrientationPlus(object);
1614 ..param.object:A Blast HSP object.
1615 ...type:Spec.FullInfo
1616 ..returns:True if the database sequence is in forward orientation.
1617 ...type:nolink:bool
1618 ..include:seqan/blast.h
1619 */
1620 template<typename TBlastSpec>
1621 inline bool
databaseOrientationPlus(BlastHsp<TBlastSpec,FullInfo> & blastHsp)1622 databaseOrientationPlus(BlastHsp<TBlastSpec, FullInfo>& blastHsp)
1623 {
1624 SEQAN_CHECKPOINT
1625 return blastHsp.db_strand;
1626 }
1627
1628
1629
1630 // only for fullinfo protein spec
1631 template<typename TSpec>
1632 inline int &
queryFrame(BlastHsp<ProteinBlast<TSpec>,FullInfo> & blastHsp)1633 queryFrame(BlastHsp<ProteinBlast<TSpec>, FullInfo>& blastHsp)
1634 {
1635 SEQAN_CHECKPOINT
1636 return blastHsp.query_frame;
1637 }
1638
1639 template<typename TSpec>
1640 inline int
getQueryFrame(BlastHsp<ProteinBlast<TSpec>,FullInfo> & blastHsp)1641 getQueryFrame(BlastHsp<ProteinBlast<TSpec>, FullInfo>& blastHsp)
1642 {
1643 SEQAN_CHECKPOINT
1644 return blastHsp.query_frame;
1645 }
1646
1647
1648 template<typename TSpec>
1649 inline int &
databaseFrame(BlastHsp<ProteinBlast<TSpec>,FullInfo> & blastHsp)1650 databaseFrame(BlastHsp<ProteinBlast<TSpec>, FullInfo>& blastHsp)
1651 {
1652 SEQAN_CHECKPOINT
1653 return blastHsp.db_frame;
1654 }
1655
1656 template<typename TSpec>
1657 inline int
getDatabaseFrame(BlastHsp<ProteinBlast<TSpec>,FullInfo> & blastHsp)1658 getDatabaseFrame(BlastHsp<ProteinBlast<TSpec>, FullInfo>& blastHsp)
1659 {
1660 SEQAN_CHECKPOINT
1661 return blastHsp.db_frame;
1662 }
1663
1664 template<typename TSpec>
1665 inline unsigned int &
percentPositives(BlastHsp<ProteinBlast<TSpec>,FullInfo> & blastHsp)1666 percentPositives(BlastHsp<ProteinBlast<TSpec>, FullInfo>& blastHsp)
1667 {
1668 SEQAN_CHECKPOINT
1669 return blastHsp.positives;
1670 }
1671
1672 template<typename TSpec>
1673 inline unsigned int
getPercentPositives(BlastHsp<ProteinBlast<TSpec>,FullInfo> & blastHsp)1674 getPercentPositives(BlastHsp<ProteinBlast<TSpec>, FullInfo>& blastHsp)
1675 {
1676 SEQAN_CHECKPOINT
1677 return blastHsp.positives;
1678 }
1679
1680 /**
1681 .Function.length:
1682 ..cat:Blast
1683 ..param.object:
1684 ...type:Class.BlastHsp
1685 ..include:seqan/blast.h
1686 */
1687 template<typename TBlast, typename TSpec>
1688 inline unsigned int
length(BlastHsp<TBlast,TSpec> & blastHsp)1689 length(BlastHsp<TBlast,TSpec >& blastHsp)
1690 {
1691 return length(databaseAlignmentString(blastHsp));
1692 }
1693
1694 }// namespace SEQAN_NAMESPACE_MAIN
1695
1696 #endif //#ifndef SEQAN_HEADER_...
1697