smith_waterman_altivec.c

序列对齐 Compare a protein sequence to a protein sequence database or a DNA sequence to a DNA sequenc

            // save old values of F and H to use on next row
            vec_st(F, 0,  p);
            vec_st(H, 16, p);
            p += 16; // move ahead 32 bytes
            
            // shift into place so we have complete F and H vectors
            // that refer to the values one unit up from each cell
            // that we are currently working on.
            Fup    = vec_sld(Fup,F,14);
            Hup2   = vec_sld(Hup2,H,14);            
            
            // do the dynamic programming 
            
            // update E value
            E   = vec_subs(E,v_gapextend);
            tmp = vec_subs(H,v_gapopen);
            E   = vec_max(E,tmp);
            
            // update F value
            F   = vec_subs(Fup,v_gapextend);
            tmp = vec_subs(Hup2,v_gapopen);
            F   = vec_max(F,tmp);
            
            // add score to H
            H   = vec_adds(Hup1,v_score);
            H   = vec_subs(H,v_bias);
            
            // set H to max of H,E,F
            H   = vec_max(H,E);
            H   = vec_max(H,F); 

            
            // Update highest score encountered this far
            v_maxscore = vec_max(v_maxscore,H);
            


            // STEP 5
            
            // prefetch next residue
            k          = *p_dbseq++;
            
            // Create the actual diagonal score vector
            // and update the queue of incomplete score vectors
            
            v_score    = vec_perm(v_score_q1, v_score_load, queue1_to_score);
            v_score_q1 = vec_perm(v_score_q2, v_score_load, queue2_to_queue1);
            v_score_q2 = vec_perm(v_score_q3, v_score_load, queue3_to_queue2);
            v_score_q3 = vec_perm(v_score_q3, v_score_load, queue3_with_load);
            
            // prefetch score for next step
            v_score_load = vec_ld(16*k,query_profile_word);
            
            // load values of F and H from previous row (one unit up)
            Fup    = vec_ld(256, p);
            Hup1   = vec_ld(272, p);
            
            // save old values of F and H to use on next row
            vec_st(F, 0,  p);
            vec_st(H, 16, p);
            p += 16; // move ahead 32 bytes
            
            // shift into place so we have complete F and H vectors
            // that refer to the values one unit up from each cell
            // that we are currently working on.
            Fup    = vec_sld(Fup,F,14);
            Hup1   = vec_sld(Hup1,H,14);            
            
            // do the dynamic programming 
            
            // update E value
            E   = vec_subs(E,v_gapextend);
            tmp = vec_subs(H,v_gapopen);
            E   = vec_max(E,tmp);
            
            // update F value
            F   = vec_subs(Fup,v_gapextend);
            tmp = vec_subs(Hup1,v_gapopen);
            F   = vec_max(F,tmp);
            
            // add score to H
            H   = vec_adds(Hup2,v_score);
            H   = vec_subs(H,v_bias);
            
            // set H to max of H,E,F
            H   = vec_max(H,E);
            H   = vec_max(H,F); 
            
            
            // Update highest score encountered this far
            v_maxscore = vec_max(v_maxscore,H);
            


            // STEP 6
            
            // prefetch next residue
            k          = *p_dbseq++;
            
            // Create the actual diagonal score vector
            // and update the queue of incomplete score vectors
            
            v_score    = vec_perm(v_score_q1, v_score_load, queue1_to_score);
            v_score_q1 = vec_perm(v_score_q2, v_score_load, queue2_to_queue1);
            v_score_q2 = vec_perm(v_score_q3, v_score_load, queue3_to_queue2);
            v_score_q3 = vec_perm(v_score_q3, v_score_load, queue3_with_load);
            
            // prefetch score for next step
            v_score_load = vec_ld(16*k,query_profile_word);
            
            // load values of F and H from previous row (one unit up)
            Fup    = vec_ld(256, p);
            Hup2   = vec_ld(272, p);
            
            // save old values of F and H to use on next row
            vec_st(F, 0,  p);
            vec_st(H, 16, p);
            p += 16; // move ahead 32 bytes
            
            // shift into place so we have complete F and H vectors
            // that refer to the values one unit up from each cell
            // that we are currently working on.
            Fup    = vec_sld(Fup,F,14);
            Hup2   = vec_sld(Hup2,H,14);            
            
            // do the dynamic programming 
            
            // update E value
            E   = vec_subs(E,v_gapextend);
            tmp = vec_subs(H,v_gapopen);
            E   = vec_max(E,tmp);
            
            // update F value
            F   = vec_subs(Fup,v_gapextend);
            tmp = vec_subs(Hup2,v_gapopen);
            F   = vec_max(F,tmp);
            
            // add score to H
            H   = vec_adds(Hup1,v_score);
            H   = vec_subs(H,v_bias);
            
            // set H to max of H,E,F
            H   = vec_max(H,E);
            H   = vec_max(H,F); 
            

            
            // Update highest score encountered this far
            v_maxscore = vec_max(v_maxscore,H);
            

            
            // STEP 7
            
            // prefetch next residue
            k          = *p_dbseq++;
            
            // Create the actual diagonal score vector
            // and update the queue of incomplete score vectors
            
            v_score    = vec_perm(v_score_q1, v_score_load, queue1_to_score);
            v_score_q1 = vec_perm(v_score_q2, v_score_load, queue2_to_queue1);
            v_score_q2 = vec_perm(v_score_q3, v_score_load, queue3_to_queue2);
            v_score_q3 = vec_perm(v_score_q3, v_score_load, queue3_with_load);
            
            // prefetch score for next step
            v_score_load = vec_ld(16*k,query_profile_word);
            
            // load values of F and H from previous row (one unit up)
            Fup    = vec_ld(256, p);
            Hup1   = vec_ld(272, p);
            
            // save old values of F and H to use on next row
            vec_st(F, 0,  p);
            vec_st(H, 16, p);
            p += 16; // move ahead 32 bytes
            
            // shift into place so we have complete F and H vectors
            // that refer to the values one unit up from each cell
            // that we are currently working on.
            Fup    = vec_sld(Fup,F,14);
            Hup1   = vec_sld(Hup1,H,14);            
            
            // do the dynamic programming 
            
            // update E value
            E   = vec_subs(E,v_gapextend);
            tmp = vec_subs(H,v_gapopen);
            E   = vec_max(E,tmp);
            
            // update F value
            F   = vec_subs(Fup,v_gapextend);
            tmp = vec_subs(Hup1,v_gapopen);
            F   = vec_max(F,tmp);
            
            // add score to H
            H   = vec_adds(Hup2,v_score);
            H   = vec_subs(H,v_bias);
            
            // set H to max of H,E,F
            H   = vec_max(H,E);
            H   = vec_max(H,F); 
            

            
            // Update highest score encountered this far
            v_maxscore = vec_max(v_maxscore,H);
            


            // STEP 8
            
            // prefetch next residue
            k          = *p_dbseq++;
            
            // Create the actual diagonal score vector
            // and update the queue of incomplete score vectors
            
            v_score    = vec_perm(v_score_q1, v_score_load, queue1_to_score);
            v_score_q1 = vec_perm(v_score_q2, v_score_load, queue2_to_queue1);
            v_score_q2 = vec_perm(v_score_q3, v_score_load, queue3_to_queue2);
            v_score_q3 = vec_perm(v_score_q3, v_score_load, queue3_with_load);
            
            // prefetch score for next step
            v_score_load = vec_ld(16*k,query_profile_word);
            
            // load values of F and H from previous row (one unit up)
            Fup    = vec_ld(256, p);
            Hup2   = vec_ld(272, p);
            
            // save old values of F and H to use on next row
            vec_st(F, 0,  p);
            vec_st(H, 16, p);
            p += 16; // move ahead 32 bytes
            
            // shift into place so we have complete F and H vectors
            // that refer to the values one unit up from each cell
            // that we are currently working on.
            Fup    = vec_sld(Fup,F,14);
            Hup2   = vec_sld(Hup2,H,14);            
            
            // do the dynamic programming 
            
            // update E value
            E   = vec_subs(E,v_gapextend);
            tmp = vec_subs(H,v_gapopen);
            E   = vec_max(E,tmp);
            
            // update F value
            F   = vec_subs(Fup,v_gapextend);
            tmp = vec_subs(Hup2,v_gapopen);
            F   = vec_max(F,tmp);
            
            // add score to H
            H   = vec_adds(Hup1,v_score);
            H   = vec_subs(H,v_bias);
            
            // set H to max of H,E,F
            H   = vec_max(H,E);
            H   = vec_max(H,F); 
            
            
            // Update highest score encountered this far
            v_maxscore = vec_max(v_maxscore,H);
        }
        
        v_score_load = vec_splat_u16(0);
        
        for(;j<db_length+7;j++)
        {
            // Create the actual diagonal score vector
            // and update the queue of incomplete score vectors
            //
            // This could of course be done with only vec_perm or vec_sel,
            // but since they use different execution units we have found
            // it to be slightly faster to mix them.
            v_score    = vec_perm(v_score_q1, v_score_load, queue1_to_score);
            v_score_q1 = vec_perm(v_score_q2, v_score_load, queue2_to_queue1);
            v_score_q2 = vec_perm(v_score_q3, v_score_load, queue3_to_queue2);
            v_score_q3 = vec_perm(v_score_q3, v_score_load, queue3_with_load);
            
            // save old values of F and H to use on next row
            vec_st(F, 0,  p);
            vec_st(H, 16, p);
            p += 16; // move ahead 32 bytes
            
            // v_score_load contains all zeros
            Fup    = vec_sld(v_score_load,F,14);
            Hup1   = vec_sld(v_score_load,H,14);            
            
            // do the dynamic programming 
            
            // update E value
            E   = vec_subs(E,v_gapextend);
            tmp = vec_subs(H,v_gapopen);
            E   = vec_max(E,tmp);
            
            // update F value
            F   = vec_subs(Fup,v_gapextend);
            tmp = vec_subs(Hup1,v_gapopen);
            F   = vec_max(F,tmp);
            
            // add score to H
            H   = vec_adds(Hup2,v_score);
            H   = vec_subs(H,v_bias);
            
            // set H to max of H,E,F
            H   = vec_max(H,E);
            H   = vec_max(H,F);
            
            // Save value to use for next diagonal H 
            Hup2 = Hup1;
            
            // Update highest score encountered this far
            v_maxscore = vec_max(v_maxscore,H);
        }
        vec_st(F, 0,  p);
        vec_st(H, 16, p);

        query_profile_word += 8*alphabet_size;
    }

    // find largest score in the v_maxscore vector
    tmp = vec_sld(v_maxscore,v_maxscore,8);
    v_maxscore = vec_max(v_maxscore,tmp);
    tmp = vec_sld(v_maxscore,v_maxscore,4);
    v_maxscore = vec_max(v_maxscore,tmp);
    tmp = vec_sld(v_maxscore,v_maxscore,2);
    v_maxscore = vec_max(v_maxscore,tmp);

    // store in temporary variable
    vec_ste(v_maxscore,0,&score);
    
    // return largest score
    return score;
}

int
smith_waterman_altivec_byte(unsigned char *     query_sequence,
                            unsigned char *     query_profile_byte,
                            int                 query_length,
                            unsigned char *     db_sequence,
                            int                 db_length,
                            unsigned char       bias,
                            unsigned char       gap_open,
                            unsigned char       gap_extend,
                            struct f_struct *   f_str)
{
    int                     i,j,k,k8;
    int                     overflow;
    unsigned char *         p;
    unsigned char           score;   
    int                     alphabet_size = f_str->alphabet_size;
    unsigned char *         workspace     = (unsigned char *)f_str->workspace;
    
    vector unsigned char    Fup,Hup1,Hup2,E,F,H,tmp;
    vector unsigned char    perm;
    vector unsigned char    v_maxscore;
    vector unsigned char    v_bias,v_gapopen,v_gapextend;
    vector unsigned char    v_score;
    vector unsigned char    v_score_q1;
    vector unsigned char    v_score_q2;