~ubuntu-branches/ubuntu/saucy/ncbi-tools6/saucy-proposed

« back to all changes in this revision

Viewing changes to algo/blast/core/blast_traceback.c

  • Committer: Bazaar Package Importer
  • Author(s): Aaron M. Ucko
  • Date: 2009-08-11 22:03:47 UTC
  • mfrom: (1.4.1 upstream)
  • mto: This revision was merged to the branch mainline in revision 10.
  • Revision ID: james.westby@ubuntu.com-20090811220347-g4b6lzdvphvvbpiu
* New upstream release.
* debian/libncbi6.symbols: update accordingly.
* debian/control: clean up obsolete or redundant relationship declarations.

Show diffs side-by-side

added added

removed removed

Lines of Context:
1
 
/* $Id: blast_traceback.c,v 1.209 2009/01/12 14:19:33 kazimird Exp $
 
1
/* $Id: blast_traceback.c,v 1.212 2009/07/09 15:49:31 kazimird Exp $
2
2
 * ===========================================================================
3
3
 *
4
4
 *                            PUBLIC DOMAIN NOTICE
50
50
 
51
51
#ifndef SKIP_DOXYGEN_PROCESSING
52
52
static char const rcsid[] = 
53
 
    "$Id: blast_traceback.c,v 1.209 2009/01/12 14:19:33 kazimird Exp $";
 
53
    "$Id: blast_traceback.c,v 1.212 2009/07/09 15:49:31 kazimird Exp $";
54
54
#endif /* SKIP_DOXYGEN_PROCESSING */
55
55
 
56
56
#include <algo/blast/core/blast_traceback.h>
359
359
   const Boolean kGreedyTraceback = (ext_options->eTbackExt == eGreedyTbck);
360
360
   const Boolean kTranslateSubject = 
361
361
        (Blast_SubjectIsTranslated(program_number) || program_number == eBlastTypeRpsTblastn);
 
362
   const Boolean kFullTranslation = (fence_hit && *fence_hit);
362
363
   const Boolean kSmithWaterman = (ext_options->eTbackExt == 
363
364
                                   eSmithWatermanTbckFull);
364
365
   BlastQueryInfo* query_info = query_info_in;
375
376
   
376
377
   if (fence_hit) {
377
378
       orig_hsplist = BlastHSPListDup(hsp_list);
 
379
       *fence_hit = FALSE;  /* reset fence_hit */
378
380
   }
379
381
   
380
382
   hsp_array = hsp_list->hsp_array;
479
481
                 subject = subject_blk->oof_sequence + CODON_LENGTH;
480
482
                 if (hsp->subject.frame < 0)
481
483
                      subject += subject_length + 1;
 
484
            } else if (kFullTranslation) {
 
485
                 /* previous traceback has hit the fence, will do full translation */
 
486
                 BlastHSP* hsp_new = Blast_HSPNew();
 
487
                 hsp_new->subject.frame = hsp->subject.frame;
 
488
                 hsp_new->subject.offset = -1;
 
489
                 subject = Blast_HSPGetTargetTranslation(target_t, hsp_new, &subject_length);
 
490
                 hsp_new = Blast_HSPFree(hsp_new);
482
491
            } else
483
 
                subject = Blast_HSPGetTargetTranslation(target_t, hsp, &subject_length);
 
492
                 subject = Blast_HSPGetTargetTranslation(target_t, hsp, &subject_length);
484
493
         }
485
494
 
486
495
         if (!kIsOutOfFrame && (((hsp->query.gapped_start == 0 && 
560
569
                 adjusted_subject, gap_align, score_params, q_start, s_start, 
561
570
                 query_length, adjusted_s_length,
562
571
                 fence_hit);
 
572
             ASSERT(!(kFullTranslation && *fence_hit));
563
573
         }
564
574
         
565
575
         fence_error = (fence_hit && *fence_hit);
1078
1088
                                    hit_params->options->hitlist_size);
1079
1089
   }
1080
1090
 
 
1091
   /* post-traceback pipes */
 
1092
   BlastHSPStreamTBackClose(hsp_stream, results);
 
1093
 
1081
1094
   if (make_up_kbp)
1082
1095
       sbp->kbp_gap[0] = NULL;
1083
1096
 
1150
1163
                                 interrupt_search, progress_info);
1151
1164
   } else if (ext_params->options->compositionBasedStats > 0 ||
1152
1165
              ext_params->options->eTbackExt == eSmithWatermanTbck) {
 
1166
      /* FIXME partial sequence fetching/translation could lead to fence hit
 
1167
         and seg fault */
1153
1168
      status =
1154
1169
          Blast_RedoAlignmentCore(program_number, query, query_info, sbp,
1155
1170
                                  hsp_stream, seq_src, default_db_genetic_code,
1163
1178
      const Boolean kPhiBlast = Blast_ProgramIsPhiBlast(program_number);
1164
1179
      BlastHSPStreamResultBatch *batch = 
1165
1180
                      Blast_HSPStreamResultBatchInit(query_info->num_queries);
1166
 
      const Boolean kSubjectRanges = !(hit_params->restricted_align ||
1167
 
                                     score_params->options->is_ooframe ||
1168
 
                                     ext_params->options->eTbackExt == 
1169
 
                                                     eSmithWatermanTbck);
1170
1181
 
1171
1182
      memset((void*) &seq_arg, 0, sizeof(seq_arg));
1172
1183
 
1187
1198
         if (perform_traceback) {
1188
1199
            seq_arg.oid = batch->hsplist_array[0]->oid;
1189
1200
            seq_arg.encoding = encoding;
1190
 
            seq_arg.enable_ranges = kSubjectRanges;
1191
1201
            seq_arg.check_oid_exclusion = TRUE;
 
1202
            seq_arg.reset_ranges = FALSE;
1192
1203
            
1193
1204
            BlastSequenceBlkClean(seq_arg.seq);
1194
1205
            if (BlastSeqSrcGetSequence(seq_src, &seq_arg) < 0) {
1242
1253
                                            query_info, pattern_blk);
1243
1254
               } else {
1244
1255
                  Boolean fence_hit = FALSE;
1245
 
                  Boolean * fence_hit_ptr = NULL;
1246
 
                    
1247
 
                  /* prepare to deal with partial subject sequence
1248
 
                     ranges, if these are configured and a previous
1249
 
                     HSP list has not turned them off */
1250
 
                  if (seq_arg.enable_ranges == TRUE)
1251
 
                     fence_hit_ptr = &fence_hit;
1252
 
                    
1253
1256
                  Blast_TracebackFromHSPList(program_number, hsp_list, query,
1254
1257
                                             seq_arg.seq, query_info, 
1255
1258
                                             gap_align, sbp, score_params,
1256
1259
                                             ext_params->options, hit_params,
1257
1260
                                             seq_arg.seq->gen_code_string,
1258
 
                                             fence_hit_ptr);
 
1261
                                             &fence_hit);
1259
1262
                    
1260
1263
                  if (fence_hit) {
1261
1264
                     /* Disable range support and refetch the 
1262
1265
                        (whole) subject sequence */
1263
1266
                        
1264
 
                     seq_arg.enable_ranges = FALSE;
 
1267
                     seq_arg.reset_ranges = TRUE;
1265
1268
                     BlastSeqSrcReleaseSequence(seq_src, &seq_arg);
1266
1269
                     BlastSeqSrcGetSequence(seq_src, &seq_arg);
1267
1270
                        
1268
 
                     /* Retry the alignment */
1269
 
                       
 
1271
                     /* The C toolkit will erase genetic_code, so do it again */
 
1272
                     if (Blast_SubjectIsTranslated(program_number) && 
 
1273
                         seq_arg.seq->gen_code_string == NULL) {
 
1274
                         seq_arg.seq->gen_code_string = 
 
1275
                             GenCodeSingletonFind(default_db_genetic_code);
 
1276
                         ASSERT(seq_arg.seq->gen_code_string);
 
1277
                     }
 
1278
            
 
1279
                     /* Retry the alignment with fence_hit set*/
1270
1280
                     Blast_TracebackFromHSPList(program_number, hsp_list, 
1271
1281
                                                query, seq_arg.seq, 
1272
1282
                                                query_info, gap_align,
1274
1284
                                                ext_params->options, 
1275
1285
                                                hit_params, 
1276
1286
                                                seq_arg.seq->gen_code_string,
1277
 
                                                NULL);
 
1287
                                                &fence_hit);
 
1288
                     ASSERT(fence_hit == FALSE);
1278
1289
                  } /* fence_hit */
1279
1290
               }    /* !phi_blast */
1280
1291
 
1305
1316
      }         /* loop over all batches */
1306
1317
 
1307
1318
      batch = Blast_HSPStreamResultBatchFree(batch);
 
1319
      /* post-traceback pipes */
 
1320
      BlastHSPStreamTBackClose(hsp_stream, results);
 
1321
 
1308
1322
      BlastSequenceBlkFree(seq_arg.seq);
1309
1323
   }
1310
1324
 
1311
 
   if (hit_params->options->culling_limit > 0)
1312
 
      Blast_HSPResultsPerformCulling(results, query_info, 
1313
 
                                    hit_params->options->culling_limit,
1314
 
                                    query->length);
1315
 
 
1316
1325
   /* Re-sort the hit lists according to their best e-values, because
1317
1326
      they could have changed. Only do this for a database search. */
1318
1327
   if (BlastSeqSrcGetTotLen(seq_src) > 0)