+ fixed if/else cases bug in the increaseFeatureCount function of result_align
authorfabio <fabio@e1793c9e-67f9-0310-80fc-b846ff1f7b36>
Tue, 15 Jan 2008 11:10:22 +0000 (11:10 +0000)
committerfabio <fabio@e1793c9e-67f9-0310-80fc-b846ff1f7b36>
Tue, 15 Jan 2008 11:10:22 +0000 (11:10 +0000)
TODO
+ fix Swig/python double free issue

git-svn-id: http://svn.tuebingen.mpg.de/ag-raetsch/projects/QPalma@7356 e1793c9e-67f9-0310-80fc-b846ff1f7b36

QPalmaDP/QPalmaDP.i
QPalmaDP/debug_tools.h [new file with mode: 0644]
QPalmaDP/qpalma_dp.cpp
QPalmaDP/result_align.cpp
python/Configuration.py
python/qpalma.py

index 1e7c4e8..66ed219 100644 (file)
@@ -1,5 +1,6 @@
 %module QPalmaDP
 %{
+#define SWIG_FILE_WITH_INIT
 #include "common.h"
 #include "Mathmatics.h"
 #include "penalty_info.h"
diff --git a/QPalmaDP/debug_tools.h b/QPalmaDP/debug_tools.h
new file mode 100644 (file)
index 0000000..d7039c1
--- /dev/null
@@ -0,0 +1,17 @@
+#ifndef _DEBUG_TOOLS_
+#define _DEBUG_TOOLS_
+
+#include <cstdio>
+
+void fassert(bool exp,int line, char* file) {
+   if (! exp) {
+      printf("invalid fassert at line %d in file %s!\n",line,file);
+      exit(EXIT_FAILURE);
+   }
+}
+
+
+#define FA(expr) (fassert(expr,__LINE__,__FILE__))
+
+#endif // _DEBUG_TOOLS_
+
index c2af278..be35def 100644 (file)
@@ -1,4 +1,5 @@
 #include "qpalma_dp.h"
+#include "debug_tools.h"
 #include <cstring>
 using namespace std;
 
@@ -26,6 +27,9 @@ Alignment::Alignment(int numQPlifs, int numq) {
 
       numQualSuppPoints = numq;
       totalNumPlifs = numQPlifs;
+
+      FA( numQualSuppPoints > 0 );
+      FA( totalNumPlifs > 0 );
 }
 
 void Alignment::getDNAEST(){}
@@ -66,6 +70,7 @@ void Alignment::myalign(int nr_paths_p, char* dna, int dna_len_p, char* est,
     }  
   }
 
+
   int* max_score_positions = new int[nr_paths*2];
 
   Pre_score** matrices = new Pre_score*[nr_paths];
@@ -80,7 +85,7 @@ void Alignment::myalign(int nr_paths_p, char* dna, int dna_len_p, char* est,
 
   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);
   
-  // printf("after call to fill_matrix...\n");
+  //printf("after call to fill_matrix...\n");
   /***************************************************************************/ 
   // return arguments etc.
   /***************************************************************************/
@@ -129,7 +134,7 @@ void Alignment::myalign(int nr_paths_p, char* dna, int dna_len_p, char* est,
   memset((char*)est_align, -1, (est_len-1)*nr_paths*sizeof(int)); // fills est_align with zeros
   memset((char*)mmatrix_param, 0, mlen*mlen*nr_paths*sizeof(int)); //fills mmatrix_param with zeros
   memset(alignmentscores, -1, nr_paths*sizeof(double)); //fills alignmentscores with zeros
-  // printf("after memset...\n");
+  //printf("after memset...\n");
   
 
   for (int z=0; z<nr_paths; z++) {
@@ -144,15 +149,13 @@ void Alignment::myalign(int nr_paths_p, char* dna, int dna_len_p, char* est,
     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, qparam, currentMode);
     printf("after call to result_align...\n");
 
-  for(qidx=0;qidx<totalNumPlifs*nr_paths;qidx++) {
-   penalty_struct p;
-   p = qualityScoresAllPaths[qidx];
-   for(pidx=0;pidx<p.len;pidx++)
-      printf("%f ",p.penalties[pidx]);
-   printf("\n");
-  }
-
- printf("after printf...\n");
+  //for(qidx=0;qidx<totalNumPlifs*nr_paths;qidx++) {
+  // penalty_struct p;
+  // p = qualityScoresAllPaths[qidx];
+  // for(pidx=0;pidx<p.len;pidx++)
+  //    printf("%f ",p.penalties[pidx]);
+  // printf("\n");
+  //}
 
    //if(DNA_ARRAY != 0) {
    //    delete[] DNA_ARRAY;
@@ -234,6 +237,7 @@ void Alignment::getAlignmentResults(int* s_align, int* e_align,
    int ctr=0;
    for(idx=0; idx<qScores_size; idx++) {
       currentPlif = qualityScoresAllPaths[idx];
+      //printf("Size is %dx%d\n",qScores_size,currentPlif.len);
       for(pidx=0; pidx<currentPlif.len; pidx++) {
          qScores[ctr] = currentPlif.penalties[pidx];
          ctr++;
index 20dd4eb..c98b36b 100644 (file)
@@ -5,20 +5,20 @@ using namespace std;
 
 void increaseFeatureCount(penalty_struct* qparam, int dnanum, int estnum, double estprb) {
 
+   //printf("Current index %d dna/est: %d %d\n",(estnum-1)*6+(dnanum),dnanum,estnum);
    penalty_struct currentStruct = qparam[(estnum-1)*6+(dnanum)];
-   printf("Current index %d dna/est: %d %d\n",(estnum-1)*6+(dnanum),dnanum,estnum);
 
-      printf("before\n");
-      int p_idx;
-      for(p_idx=0;p_idx<currentStruct.len;p_idx++) {
-         printf("%f ",currentStruct.limits[p_idx]);
-      }
-      printf("\n");
+      //printf("before\n");
+      //int p_idx;
+      //for(p_idx=0;p_idx<currentStruct.len;p_idx++) {
+      //   printf("%f ",currentStruct.limits[p_idx]);
+      //}
+      //printf("\n");
    
-      for(p_idx=0;p_idx<currentStruct.len;p_idx++) {
-         printf("%f ",currentStruct.penalties[p_idx]);
-      }
-      printf("\n");
+      //for(p_idx=0;p_idx<currentStruct.len;p_idx++) {
+      //   printf("%f ",currentStruct.penalties[p_idx]);
+      //}
+      //printf("\n");
 
    double value = estprb;
    int Lower = 0;
@@ -28,11 +28,15 @@ void increaseFeatureCount(penalty_struct* qparam, int dnanum, int estnum, double
          Lower++;
    }
 
-   if (Lower == 0)
+   if (Lower == 0) {
          currentStruct.penalties[0] += 1;
+         return;
+   }
 
-   if (Lower == currentStruct.len)
+   if (Lower == currentStruct.len) {
          currentStruct.penalties[currentStruct.len-1] += 1;
+         return;
+   }
 
    Lower -= 1;
    int Upper = Lower+1; // x-werte bleiben fest
@@ -46,14 +50,14 @@ void increaseFeatureCount(penalty_struct* qparam, int dnanum, int estnum, double
 
    qparam[(estnum-1)*5+(dnanum-1)] = currentStruct;
 
-   printf("after\n");
-   for(p_idx=0;p_idx<currentStruct.len;p_idx++) {
-      printf("%f ",currentStruct.limits[p_idx]);
-   }
-   printf("\n");
-   for(p_idx=0;p_idx<currentStruct.len;p_idx++)
-      printf("%f ",currentStruct.penalties[p_idx]);
-   printf("\n");
+   //printf("after\n");
+   //for(p_idx=0;p_idx<currentStruct.len;p_idx++) {
+   //   printf("%f ",currentStruct.limits[p_idx]);
+   //}
+   //printf("\n");
+   //for(p_idx=0;p_idx<currentStruct.len;p_idx++)
+   //   printf("%f ",currentStruct.penalties[p_idx]);
+   //printf("\n");
 }
 
 bool result_align(Pre_score* matrices[], int z, int est_len, int dna_len, int* result_length_ptr, char* est, char* dna, double* prb, double* chastity, int* s_align, int* e_align, int* mparam, double* alignmentscores, int* max_score_positions, penalty_struct* qparam, mode currentMode)
@@ -70,6 +74,7 @@ bool result_align(Pre_score* matrices[], int z, int est_len, int dna_len, int* r
 
   double prbnum ;
   double chastitynum ;
+
   
   int dna_pos ;
   int est_pos ;
@@ -111,6 +116,7 @@ bool result_align(Pre_score* matrices[], int z, int est_len, int dna_len, int* r
   splice_state = 0 ; //exon
   est_state = 1 ; //[112222111 ... for exons of length 2,4,3 ...
 
+
   //compute length of alignment (between max(m,n) and m+n)     
   //local_alignment: backtracking to first Zero: look at actual value
   while(!(((Pre_score*)matrices[act_m_no] + i*dna_len +j)->value)==0.0) {
index 554ec24..f486b66 100644 (file)
@@ -114,7 +114,7 @@ if mode == 'normal':
    numQualPlifs = 0
 elif mode == 'using_quality_scores':
    sizeMatchmatrix   = (6,6)
-   numQualPlifs = 5*5 #5*6
+   numQualPlifs = 5*6
 else:
    assert False, 'Wrong operation mode specified'
 
index ff4fc31..fc7bb3d 100644 (file)
@@ -206,6 +206,7 @@ class QPalma:
             acceptor = QPalmaDP.createDoubleArrayFromList(acceptor)
 
             currentAlignment = QPalmaDP.Alignment(Configuration.numQualPlifs,Configuration.numQualSuppPoints)
+
             c_qualityPlifs = QPalmaDP.createPenaltyArrayFromList([elem.convert2SWIG() for elem in qualityPlifs])
 
             # calculates SpliceAlign, EstAlign, weightMatch, Gesamtscores, dnaest
@@ -219,13 +220,11 @@ class QPalma:
             c_EstAlign          = QPalmaDP.createIntArrayFromList([0]*(est_len*num_path[exampleIdx]))
             c_WeightMatch       = QPalmaDP.createIntArrayFromList([0]*(mm_len*num_path[exampleIdx]))
             c_AlignmentScores   = QPalmaDP.createDoubleArrayFromList([.0]*num_path[exampleIdx])
-           
-            emptyPlif = Plf(Configuration.numQualSuppPoints)
-            emptyPlif = emptyPlif.convert2SWIG()
+
             c_qualityPlifsFeatures = QPalmaDP.createDoubleArrayFromList([.0]*(Configuration.numQualPlifs*num_path[exampleIdx]))
 
-            currentAlignment.getAlignmentResults(c_SpliceAlign, c_EstAlign,\
-            c_WeightMatch, c_AlignmentScores, c_qualityPlifsFeatures)
+            #currentAlignment.getAlignmentResults(c_SpliceAlign, c_EstAlign,\
+            #c_WeightMatch, c_AlignmentScores, c_qualityPlifsFeatures)
 
             print 'Python: after getAlignmentResults...'
 
@@ -253,12 +252,11 @@ class QPalma:
                AlignmentScores[i+1] = c_AlignmentScores[i]
 
             print "Calling destructors"
-
             del c_SpliceAlign
             del c_EstAlign
             del c_WeightMatch
             del c_AlignmentScores
-            del c_qualityPlifs
+            del c_qualityPlifsFeatures
             del currentAlignment
 
             newSpliceAlign = newSpliceAlign.reshape(num_path[exampleIdx],dna_len)
@@ -274,13 +272,15 @@ class QPalma:
             true_map[0] = 1
             path_loss = [0]*(num_path[exampleIdx]+1)
 
+            print 'Calculating decoded features'
+
             for pathNr in range(num_path[exampleIdx]):
                weightDon, weightAcc, weightIntron = computeSpliceWeights(d, a, h, newSpliceAlign[pathNr,:].flatten().tolist()[0], don_supp, acc_supp)
 
                decodedQualityFeatures = zeros((Configuration.totalQualSuppPoints,1))
 
-               for qidx in range(Configuration.numQualPlifs*pathNr,Configuration.numQualPlifs*(pathNr+1)):
-                decodedQualityFeatures[qidx%Configuration.numQualPlifs] = c_qualityPlifsFeatures[qidx]
+               #for qidx in range(Configuration.numQualPlifs*pathNr,Configuration.numQualPlifs*(pathNr+1)):
+               #   decodedQualityFeatures[qidx%Configuration.numQualPlifs] = c_qualityPlifsFeatures[1]
 
                # sum up positionwise loss between alignments
                for alignPosIdx in range(len(newSpliceAlign[pathNr,:])):
@@ -334,7 +334,6 @@ class QPalma:
                #
                # end of one example processing 
                #
-            del c_qualityPlifsFeatures
 
             # call solver every nth step
             if exampleIdx != 0 and exampleIdx % 3 == 0 and not __debug__:
@@ -363,7 +362,7 @@ class QPalma:
       # end of optimization 
       #  
       print 'Training completed'
-      pdb.set_trace()
+      #pdb.set_trace()
       export_param('elegans.param',h,d,a,mmatrix)
       self.logfh.close()
 
@@ -390,5 +389,3 @@ def fetch_genome_info(ginfo_filename):
 if __name__ == '__main__':
    qpalma = QPalma()
    qpalma.run()
-   pdb.set_trace()
-   l = [elem+1 for elem in range(44)]