hit_set.cpp

   1 /*
   2  * hit_set.cpp
   3  *
   4  *  Created on: Jul 31, 2009
   5  *      Author: Ben Langmead
   6  */
   7
   8 #include <iostream>
   9 #include <vector>
  10 #include <seqan/sequence.h>
  11 #include "alphabet.h"
  12 #include "hit_set.h"
  13
  14 using namespace std;
  15 using namespace seqan;
  16
  17 /**
  18  * Report up to 'khits' hits from this HitSet.
  19  */
  20 void HitSet::reportUpTo(ostream& os, int khits) {
  21         khits = min(khits, (int)size());
  22         String<Dna5> seqrc;
  23         String<char> qualr;
  24         for(int i = 0; i < khits; i++) {
  25                 const HitSetEnt& h = ents[i];
  26                 if(!h.fw && seqan::empty(seqrc)) {
  27                         // Lazily initialize seqrc and qualr
  28                         seqrc = seq;
  29                         reverseComplementInPlace(seqrc, color);
  30                         assert_eq(seqan::length(seqrc), seqan::length(seq));
  31                         qualr = qual;
  32                         reverseInPlace(qualr);
  33                         assert_eq(seqan::length(qualr), seqan::length(qual));
  34                 }
  35                 os << name << '\t'
  36                    << (h.fw ? '+' : '-') << '\t'
  37                    << h.h.first << '\t'
  38                    << h.h.second << '\t'
  39                    << (h.fw ? seq : seqrc) << '\t'
  40                    << (h.fw ? qual : qualr) << '\t'
  41                    << h.oms << '\t';
  42                 for(size_t i = 0; i < h.edits.size(); i++) {
  43                         const Edit& e = h.edits[i];
  44                         os << e.pos;
  45                         if(e.type == EDIT_TYPE_SNP) os << "S";
  46                         os << ":" << (char)e.chr << ">" << (e.qchr != 0 ? (char)e.qchr : (char)seq[e.pos]);
  47                         if(i < h.edits.size()-1 || !h.cedits.empty()) os << ",";
  48                 }
  49                 for(size_t i = 0; i < h.cedits.size(); i++) {
  50                         const Edit& e = h.cedits[i];
  51                         os << e.pos;
  52                         if(e.type == EDIT_TYPE_SNP) os << "S";
  53                         os << ":" << (char)e.chr << ">" << (e.qchr != 0 ? (char)e.qchr : (char)seq[e.pos]);
  54                         if(i < h.cedits.size()-1) os << ",";
  55                 }
  56                 os << endl;
  57         }
  58 }
  59
  60 ostream& operator << (ostream& os, const HitSetEnt& hs) {
  61         os << "\t" << hs.h.first << ":" << hs.h.second;
  62         return os;
  63 }
  64
  65 ostream& operator << (ostream& os, const HitSet& hs) {
  66         os << hs.name << ":" << hs.seq << ":" << hs.qual << endl;
  67         vector<HitSetEnt>::const_iterator it;
  68         for(it = hs.ents.begin(); it != hs.ents.end(); it++) {
  69                 os << (*it);
  70         }
  71         return os;
  72 }