+ fixed some index bugs in the evaluation
[qpalma.git] / dyn_prog / qpalma_dp.cpp
index e2bdb35..a20af18 100644 (file)
@@ -1,3 +1,4 @@
+
 #include "qpalma_dp.h"
 #include <cstring>
 
@@ -24,7 +25,7 @@ Alignment::Alignment(int numQPlifs, int numq, bool use_qscores) {
       est_align            = 0;
       mmatrix_param        = 0;
       alignmentscores      = 0;
-      qualityScoresAllPaths = 0;
+      qualityFeaturesAllPaths = 0;
       mlen = 6; // score matrix: length of 6 for "- A C G T N"
 
       numQualSuppPoints = numq;
@@ -62,6 +63,8 @@ void Alignment::myalign(int nr_paths_p, char* dna, int dna_len_p, char* est,
     // dnaest
    DNA_ARRAY = 0;
    EST_ARRAY = 0;
+
+  //printf("[qpalma_dp] read is: %s\n",est);
   
   //possible donor positions
   int nr_donor_sites = 0 ;
@@ -96,14 +99,13 @@ void Alignment::myalign(int nr_paths_p, char* dna, int dna_len_p, char* est,
    //      printf("%f ",pidx,p.penalties[pidx]);
    //   printf("\n");
    //}
-   //
 
    //int h_idx;
    //for(h_idx=0;h_idx<h.len;h_idx++)
    //   printf("%d %f %f\n",h_idx,h.limits[h_idx],h.penalties[h_idx]);
   
   //printf("calling fill_matrix...\n");
-  fill_matrix(nr_paths, matrices, est_len, dna_len, est, dna, prb, &h, matchmatrix, qualityScores, donor, acceptor, remove_duplicate_scores, nr_donor_sites, donor_sites, max_score_positions,currentMode);
+  fill_matrix(nr_paths, matrices, est_len, dna_len, est, dna, prb, &h, matchmatrix, qualityScores, donor, acceptor, remove_duplicate_scores, nr_donor_sites, donor_sites, max_score_positions);
   //printf("after call to fill_matrix...\n");
 
   /***************************************************************************/ 
@@ -142,7 +144,7 @@ void Alignment::myalign(int nr_paths_p, char* dna, int dna_len_p, char* est,
   memset(alignmentscores, -1, nr_paths*sizeof(double)); //fills alignmentscores with zeros
   //printf("after memset...\n");
   
-  qualityScoresAllPaths= new penalty_struct*[nr_paths];
+  qualityFeaturesAllPaths= new penalty_struct*[nr_paths];
 
   for (int z=0; z<nr_paths; z++) {
     result_length = 0 ;
@@ -151,37 +153,43 @@ void Alignment::myalign(int nr_paths_p, char* dna, int dna_len_p, char* est,
     int* e_align = est_align + (est_len-1)*z ;    //pointer
     int* mparam  = mmatrix_param + mmatrix_size*z; //pointer
 
-   qualityScoresAllPaths[z] = new penalty_struct[numPlifs];
+   qualityFeaturesAllPaths[z] = new penalty_struct[numPlifs];
 
    int qidx, pidx;
    for(qidx=0;qidx<numPlifs;qidx++) {
       penalty_struct p;
+      init_penalty_struct(p);
       p.len = numQualSuppPoints;
-      p.limits = (REAL*) calloc(p.len,sizeof(REAL));
+      p.limits = (double*) calloc(p.len,sizeof(double));
 
       for(pidx=0;pidx<p.len;pidx++)
-         p.limits[pidx] = qualityScores[qidx%numPlifs].limits[pidx];
+         p.limits[pidx] = qualityScores[qidx].limits[pidx];
 
-      p.penalties = (REAL*) calloc(p.len,sizeof(REAL));
-      qualityScoresAllPaths[z][qidx] = p;
+      p.penalties = (double*) calloc(p.len,sizeof(double));
+      qualityFeaturesAllPaths[z][qidx] = p;
    }
 
    //penalty_struct* qparam = qualityScoresAllPaths + (numPlifs*z);
 
     //printf("before call to result_align...\n");
-    bool no_more_path = result_align(matrices, z, est_len, dna_len, &result_length, est, dna, prb, chastity, s_align, e_align, mparam, alignmentscores, max_score_positions, qualityScoresAllPaths[z] , currentMode);
+    bool no_more_path = result_align(matrices, z, est_len, dna_len, &result_length, est, dna, prb, chastity, s_align, e_align, mparam, alignmentscores, max_score_positions, qualityFeaturesAllPaths[z] , currentMode);
     //printf("after call to result_align...\n");
 
    //printf("z is %d\n",z);
+   //int len;
    //for(qidx=0;qidx<numPlifs;qidx++) {
    //   penalty_struct p;
-   //   p = qualityScoresAllPaths[qidx+numPlifs*z];
-   //   printf("%d ",qidx);
+   //   p = qualityScoresAllPaths[z][qidx];
+   //   printf("%d: ",qidx);
    //   for(pidx=0;pidx<p.len;pidx++)
-   //      printf("%f ",pidx,p.penalties[pidx]);
+   //      printf("%f ",p.limits[pidx]);
+   //   printf("\n");
+
+   //   for(pidx=0;pidx<p.len;pidx++)
+   //      printf("%f ",p.penalties[pidx]);
    //   printf("\n");
    //}
-   //
+
    if(z==0) {
 
    if(DNA_ARRAY != 0) {
@@ -260,25 +268,20 @@ void Alignment::getAlignmentResults(int* s_align, int* e_align,
       for (int z=0; z<nr_paths; z++) {
          for(int estChar=1;estChar<6;estChar++) {
             for(int dnaChar=0;dnaChar<6;dnaChar++) {
+
                int currentPos = (estChar-1)*6+dnaChar;
-               currentPlif = qualityScoresAllPaths[z][currentPos];
+               currentPlif = qualityFeaturesAllPaths[z][currentPos];
+
                for(int pidx=0; pidx<currentPlif.len; pidx++) {
                   qScores[ctr] = currentPlif.penalties[pidx];
+                  //printf("%f ",qScores[ctr]);
                   ctr++;
                }
+               //printf("\n");
       }}}
-   }
 
-   //for(idx=0; idx<numPathsPlifs; idx++) {
-   //   currentPlif = qualityScoresAllPaths[idx];
-   //   //printf("Size is %dx%d\n",numPathsPlifs,currentPlif.len);
-   //   for(pidx=0; pidx<currentPlif.len; pidx++) {
-   //      qScores[ctr] = currentPlif.penalties[pidx];
-   //      ctr++;
-   //   }
-   //}
-
-   //printf("\nctr is %d\n",ctr);
+      //printf("\nctr is %d\n",ctr);
+   }
    //printf("Leaving getAlignmentResults...\n");
 }