search_seeded_phase4.c

   1 /*
   2  * This is a fragment, included from multiple places in ebwt_search.cpp.
   3  * It implements the logic of the fourth phase of the seeded, quality-
   4  * aware search routine.  It is implemented as a code fragment so that
   5  * it can be reused in both the half-index-in-memory and full-index-in-
   6  * memory situations.
   7  */
   8 {
   9         if(!nofw) {
  10                 params.setFw(true);
  11                 btf4.setReportExacts(true);
  12                 btf4.setQuery(patsrc->bufa());
  13                 // Get all partial alignments for this read's reverse
  14                 // complement
  15                 pals.clear();
  16                 if(pamFw != NULL && pamFw->size() > 0) {
  17                         // We can get away with an unsynchronized call because there
  18                         // are no writers for pamFw in this phase
  19                         pamFw->getPartialsUnsync(patid, pals);
  20                         pamFw->clear(patid);
  21                         assert_eq(0, pamFw->size());
  22                 }
  23                 bool done = false;
  24                 if(pals.size() > 0) {
  25                         // Partial alignments exist - extend them
  26                         // Set up seed bounds
  27                         if(qs < s) {
  28                                 btf4.setOffs(0, 0, qs, qs, qs, qs);
  29                         } else {
  30                                 btf4.setOffs(0, 0, s, s, s, s);
  31                         }
  32                         for(size_t i = 0; i < pals.size(); i++) {
  33                                 seqan::clear(muts);
  34                                 uint8_t oldQuals =
  35                                         PartialAlignmentManager::toMutsString(
  36                                                 pals[i], patFwRev, qualRev, muts, !noMaqRound);
  37                                 assert_leq(oldQuals, qualThresh);
  38                                 // Set the backtracking thresholds appropriately
  39                                 // Now begin the backtracking, treating the first
  40                                 // 24 bases as unrevisitable
  41                                 ASSERT_ONLY(String<Dna5> tmp = patFwRev);
  42                                 btf4.setMuts(&muts);
  43                                 done = btf4.backtrack(oldQuals);
  44                                 btf4.setMuts(NULL);
  45                                 assert_eq(tmp, patFwRev); // assert mutations were undone
  46                                 if(done) {
  47                                         // Got a hit; stop processing partial
  48                                         // alignments
  49                                         break;
  50                                 }
  51                         } // Loop over partial alignments
  52                 }
  53                 seqan::clear(muts);
  54
  55                 // Case 4F yielded a hit; continue to next pattern
  56                 if(done) continue;
  57
  58                 if(sink->finishedWithStratum(1)) { // no more exact hits are possible
  59                         continue;
  60                 }
  61
  62                 // If we're in two-mismatch mode, then now is the time to
  63                 // try the final case that might apply to the forward
  64                 // pattern: 1 mismatch in each of the 3' and 5' halves of
  65                 // the seed.
  66                 bool gaveUp = false;
  67                 if(seedMms >= 2) {
  68                         btf24.setQuery(patsrc->bufa());
  69                         // Set up seed bounds
  70                         if(qs < s) {
  71                                 btf24.setOffs(qs5, qs,
  72                                                          0,                        // unrevOff
  73                                                          (seedMms <= 2)? qs5 : 0,  // 1revOff
  74                                                          (seedMms < 3)? qs : qs5,  // 2revOff
  75                                                          qs);                      // 3revOff
  76                         } else {
  77                                 btf24.setOffs(s5, s,
  78                                                           0,                       // unrevOff
  79                                                           (seedMms <= 2)? s5 : 0,  // 1revOff
  80                                                           (seedMms < 3)?  s : s5,  // 2revOff
  81                                                           s);                      // 3revOff
  82                         }
  83                         done = btf24.backtrack();
  84                         if(btf24.numBacktracks() == btf24.maxBacktracks()) {
  85                                 gaveUp = true;
  86                         }
  87                         if(done) {
  88                                 btf24.resetNumBacktracks();
  89                                 continue;
  90                         }
  91                         btf24.resetNumBacktracks();
  92                 }
  93         }
  94 }