more chimera.slayer debugging

diff --git a/chimeraslayer.cpp b/chimeraslayer.cpp
index efd29d64035e390c7b48aba87f3b93fdd8c1eb68..423f7c0e713eef53a8904c2dcaff49a35844f35f 100644 (file)
--- a/chimeraslayer.cpp
+++ b/chimeraslayer.cpp
@@ -785,6 +785,7 @@ int ChimeraSlayer::getChimeras(Sequence* query) {
 
                        }
                        
+//                     cout << query->getAligned() << endl;
                        //get sequence that were given from maligner results
                        vector<SeqDist> seqs;
                        map<string, float> removeDups;
@@ -792,9 +793,12 @@ int ChimeraSlayer::getChimeras(Sequence* query) {
                        map<string, string> parentNameSeq;
                        map<string, string>::iterator itSeq;
                        for (int j = 0; j < Results.size(); j++) {
+
                                float dist = (Results[j].regionEnd - Results[j].regionStart + 1) * Results[j].queryToParentLocal;
                                //only add if you are not a duplicate
-
+//                             cout << Results[j].parent << '\t' << Results[j].regionEnd << '\t' << Results[j].regionStart << '\t' << Results[j].regionEnd - Results[j].regionStart +1 << '\t' << Results[j].queryToParentLocal << '\t' << dist << endl;
+                               
+                               
                                if(Results[j].queryToParentLocal >= 90){        //local match has to be over 90% similarity
                                
                                        itDup = removeDups.find(Results[j].parent);
@@ -840,7 +844,7 @@ int ChimeraSlayer::getChimeras(Sequence* query) {
                        for (int k = 0; k < seqs.size(); k++) {  
 //                             cout << seqs[k].seq->getAligned() << endl;
                                seqsForSlayer.push_back(seqs[k].seq);   
-                       
+//                             cout << seqs[k].seq->getName() << endl;
                        }
                        
                        if (m->control_pressed) {  for (int k = 0; k < seqs.size(); k++) {  delete seqs[k].seq;   }  return 0;  }

diff --git a/maligner.cpp b/maligner.cpp
index 1205cef5e36b8bb4381f5ae8bb516c30e934318d..0cd8336c6de865d85d0fea9d607199599f10ab9e 100644 (file)
--- a/maligner.cpp
+++ b/maligner.cpp
@@ -144,8 +144,8 @@ string Maligner::chimeraMaligner(int chimeraPenalty, DeCalculator* decalc) {
                        temp.parentAligned = db[seqIndex]->getAligned();
                        temp.nastRegionStart = spotMap[regionStart];
                        temp.nastRegionEnd = spotMap[regionEnd];
-                       temp.regionStart = regionStart;
-                       temp.regionEnd = regionEnd;
+                       temp.regionStart = unalignedMap[regionStart];
+                       temp.regionEnd = unalignedMap[regionEnd];
                        
                        string parentInRange = refSeqs[seqIndex]->getAligned();
                        parentInRange = parentInRange.substr(traceStart, (traceEnd-traceStart+1));
@@ -280,7 +280,22 @@ void Maligner::verticalFilter(vector<Sequence*> seqs) {
                        
                        seqs[i]->setAligned(newAligned);
                }
-
+               
+               string query = seqs[seqs.size()-1]->getAligned();
+               int queryLength = query.length();
+               
+               unalignedMap.resize(queryLength, 0);
+               
+               
+               for(int i=1;i<queryLength;i++){
+                       if(query[i] != '.' && query[i] != '-'){
+                               unalignedMap[i] = unalignedMap[i-1] + 1;
+                       }
+                       else{
+                               unalignedMap[i] = unalignedMap[i-1];
+                       }
+               }
+               
                spotMap = newMap;
        }
        catch(exception& e) {
@@ -436,12 +451,13 @@ vector<score_struct> Maligner::extractHighestPath(vector<vector<score_struct> >
                        }
                }
                
-//             cout << highestScore << endl;
                vector<score_struct> path;
                
                int rowIndex = highestStruct.row;
                int pos = highestStruct.col;
                int score = highestStruct.score;
+
+//             cout << rowIndex << '\t' << pos << '\t' << score << endl;
                
                while (pos >= 0 && score > 0) {
                        score_struct temp = ms[rowIndex][pos];

diff --git a/maligner.h b/maligner.h
index da270d6ac10c2b43a5ee6100dd9e8c33ea5c71ae..4f5c949a45112ebed570f4f6df46cd60f1aeb23e 100644 (file)
--- a/maligner.h
+++ b/maligner.h
@@ -36,7 +36,7 @@ class Maligner {
                float minDivR, percentIdenticalQueryChimera;
                vector<results> outputResults;
                map<int, int> spotMap;
-               
+               vector<int> unalignedMap;
                vector<Sequence*> minCoverageFilter(vector<Sequence*>);  //removes top matches that do not have minimum coverage with query.
                int computeChimeraPenalty();
                void verticalFilter(vector<Sequence*>);

diff --git a/slayer.cpp b/slayer.cpp
index 8e17f87b109500f34abb0aa573b263627011c19a..92b4673979197d4591b771a1011470d46f22b973 100644 (file)
--- a/slayer.cpp
+++ b/slayer.cpp
@@ -34,7 +34,7 @@ string Slayer::getResults(Sequence* query, vector<Sequence*> refSeqs) {
                                vector<data_struct> divs = runBellerophon(q, leftParent, rightParent, spots);
        
                                if (m->control_pressed) { delete q; delete leftParent; delete rightParent; return "no"; }
-                                       
+//                             cout << "examining:\t" << refSeqs[i]->getName() << '\t' << refSeqs[j]->getName() << endl;
                                vector<data_struct> selectedDivs;
                                for (int k = 0; k < divs.size(); k++) {
                                        
@@ -46,6 +46,7 @@ string Slayer::getResults(Sequence* query, vector<Sequence*> refSeqs) {
                                        int numSNPSLeft = snpsLeft.size();
                                        int numSNPSRight = snpsRight.size();
                                        
+//                                     cout << numSNPSLeft << '\t' << numSNPSRight << endl;
                                        //require at least 4 SNPs on each side of the break
                                        if ((numSNPSLeft >= 4) && (numSNPSRight >= 4)) {