aa_ungapped.c

ncbi源码

            /* if the hsp meets the score threshold, report it */
            if (score >= cutoff) {
               BlastSaveInitHsp(ungapped_hsps, hsp_q, hsp_s, 
                  query_offsets[i], subject_offsets[i], hsp_len, score);
            }
            diag_array[diag_coord].diag_level = 
                  s_last_off - (lookup->wordsize - 1) + diag_offset;
         }
      } /* end for */
   } /* end while */

   /* increment the offset in the diagonal array (no windows used) */
   DiagUpdate(diag, subject->length); 
 
   Blast_UngappedStatsUpdate(ungapped_stats, totalhits, hits_extended, 
                             ungapped_hsps->total);
   return 0;
}

Int4 BlastAaExtendRight(Int4 ** matrix,
			const BLAST_SequenceBlk* subject,
			const BLAST_SequenceBlk* query,
			Int4 s_off,
			Int4 q_off,
			Int4 dropoff,
			Int4* displacement,
                        Int4 maxscore,
                        Int4* s_last_off)
{
  Int4 i, n, best_i = -1;
  Int4 score = maxscore;

  Uint1* s,* q;
  n = MIN( subject->length - s_off , query->length - q_off );

  s=subject->sequence + s_off;
  q=query->sequence + q_off;

  for(i = 0; i < n; i++) {
     score += matrix[ q[i] ] [ s[i] ];

     if (score > maxscore) {
        maxscore = score;
        best_i = i;
     }

     if (score <= 0 || (maxscore - score) >= dropoff)
        break;
  }

  *displacement = best_i;
  *s_last_off = s_off + i;
  return maxscore;
}

Int4 BlastAaExtendLeft(Int4 ** matrix,
		       const BLAST_SequenceBlk* subject,
		       const BLAST_SequenceBlk* query,
		       Int4 s_off,
		       Int4 q_off,
		       Int4 dropoff,
		       Int4* displacement)
{
   Int4 i, n, best_i;
   Int4 score = 0;
   Int4 maxscore = 0;
   
   Uint1* s,* q;
   
   n = MIN( s_off , q_off );
   best_i = n + 1;

   s = subject->sequence + s_off - n;
   q = query->sequence + q_off - n;
   
   for(i = n; i >= 0; i--) {
      score += matrix[ q[i] ] [ s[i] ];
      
      if (score > maxscore) {
         maxscore = score;
         best_i = i;
      }
      if ((maxscore - score) >= dropoff)
         break;
   }
   
   *displacement = n - best_i;
   return maxscore;
}

Int4 BlastPSSMExtendRight(Int4 ** matrix,
			const BLAST_SequenceBlk* subject,
			Int4 query_size,
			Int4 s_off,
			Int4 q_off,
			Int4 dropoff,
			Int4* displacement,
                        Int4 maxscore,
			Int4* s_last_off)
{
  Int4 i, n, best_i = -1;
  Int4 score = maxscore;
  Uint1* s;

  n = MIN( subject->length - s_off , query_size - q_off );
  s = subject->sequence + s_off;

  for(i = 0; i < n; i++) {
     score += matrix[ q_off+i ] [ s[i] ];

     if (score > maxscore) {
        maxscore = score;
        best_i = i;
     }

     if (score <= 0 || (maxscore - score) >= dropoff)
        break;
  }

  *displacement = best_i;
  *s_last_off = s_off + i;
  return maxscore;
}

Int4 BlastPSSMExtendLeft(Int4 ** matrix,
		       const BLAST_SequenceBlk* subject,
		       Int4 query_size,
		       Int4 s_off,
		       Int4 q_off,
		       Int4 dropoff,
		       Int4* displacement)
{
   Int4 i, n, best_i;
   Int4 score = 0;
   Int4 maxscore = 0;
   Uint1* s;
   
   n = MIN( s_off , q_off );
   best_i = n + 1;
   s = subject->sequence + s_off - n;
   
   for(i = n; i >= 0; i--) {
      score += matrix[ q_off-n+i ] [ s[i] ];

      if (score > maxscore) {
         maxscore = score;
         best_i = i;
      }
      if ((maxscore - score) >= dropoff)
         break;
   }
   
   *displacement = n - best_i;
   return maxscore;
}

Int4 BlastAaExtendOneHit(Int4 ** matrix,
                         const BLAST_SequenceBlk* subject,
                         const BLAST_SequenceBlk* query,
                         Int4 s_off,
                         Int4 q_off,
                         Int4 dropoff,
			 Int4* hsp_q,
			 Int4* hsp_s,
			 Int4* hsp_len,
                         Boolean use_pssm,
                         Int4* s_last_off)
{
  Int4 left_score, right_score;
  Int4 left_disp, right_disp;

  if (use_pssm) {
     left_score = BlastPSSMExtendLeft(matrix, subject, query->length, 
                           s_off, q_off, dropoff, &left_disp);
     right_score = BlastPSSMExtendRight(matrix, subject, query->length, 
                           s_off+1, q_off+1, dropoff, &right_disp, 
                           left_score, s_last_off);
  }
  else {
     left_score = BlastAaExtendLeft(matrix, subject, query, 
                            s_off, q_off, dropoff, &left_disp);
     right_score = BlastAaExtendRight(matrix, subject, query, 
                            s_off+1, q_off+1, dropoff, &right_disp, 
                            left_score, s_last_off);
  }
  
  *hsp_q   = q_off - left_disp;
  *hsp_s   = s_off - left_disp;
  *hsp_len = left_disp + right_disp + 2;
  
  return right_score;
}

Int4 
BlastAaExtendTwoHit(Int4 ** matrix,
                    const BLAST_SequenceBlk* subject,
                    const BLAST_SequenceBlk* query,
                    Int4 s_left_off,
                    Int4 s_right_off,
                    Int4 q_right_off,
                    Int4 dropoff,
                    Int4* hsp_q,
                    Int4* hsp_s,
                    Int4* hsp_len,
                    Boolean use_pssm,
                    Int4 word_size,
                    Boolean *right_extend,
                    Int4 *s_last_off)
{
   Int4 left_d = 0, right_d = 0; /* left and right displacements */
   Int4 left_score = 0, right_score = 0; /* left and right scores */
   Int4 i, score = 0;
   Uint1 *s = subject->sequence;
   Uint1 *q = query->sequence;

   /* find the position (up to word_size-1 letters to the
      right) that gives the largest starting score */
   for (i = 0; i < word_size; i++) {
      if (use_pssm)
         score += matrix[ q_right_off+i ][ s[s_right_off+i] ];
      else
         score += matrix[ q[q_right_off+i] ][ s[s_right_off+i] ];

      if (score > left_score) {
         left_score = score;
         right_d = i;
      }
   }
   q_right_off += right_d;
   s_right_off += right_d;

   right_d = -1;
   *right_extend = FALSE;
   *s_last_off = s_right_off;

   /* first, try to extend left, from the second hit to the first hit. */
   if (use_pssm)
      left_score = BlastPSSMExtendLeft(matrix, subject, query->length,
    			     s_right_off, q_right_off, dropoff, &left_d);
   else
      left_score = BlastAaExtendLeft(matrix, subject, query,
    			     s_right_off, q_right_off, dropoff, &left_d);

   /* Extend to the right only if left extension reached the first hit. */
   if (left_d >= (s_right_off - s_left_off)) {

      *right_extend = TRUE;
      if (use_pssm)
         right_score = BlastPSSMExtendRight(matrix, subject, query->length,
                                     s_right_off + 1, q_right_off + 1,
                                     dropoff, &right_d, left_score, s_last_off);
      else 
         right_score = BlastAaExtendRight(matrix, subject, query,
                                     s_right_off + 1, q_right_off + 1,
                                     dropoff, &right_d, left_score, s_last_off);
   }

   *hsp_q = q_right_off - left_d;
   *hsp_s = s_right_off - left_d;
   *hsp_len = left_d + right_d + 2;
   return MAX(left_score,right_score);
}

Int4 DiagUpdate(BLAST_DiagTable* diag, Int4 length)
{
  if (diag == NULL)
    return 0;

  if (diag->offset >= INT4_MAX/2)
    {
      DiagClear(diag);
    }
  else
    {
      diag->offset += length;
    }
  return 0;
}


Int4 DiagClear(BLAST_DiagTable* diag)
{
  Int4 i,n;

  if (diag==NULL)
    return 0;

  n=diag->diag_array_length;

  for(i=0;i<n;i++)
    {
      diag->diag_array[i].diag_level = 0;
      diag->diag_array[i].last_hit = - diag->window;
    }
  return 0;
}