v 19.3

[mothur.git] / chimeraslayer.cpp

diff --git a/chimeraslayer.cpp b/chimeraslayer.cpp
index a182b846cc4f84288e029188f53c25dcf4ea93cc..70208607a00373eac0e13cfa4eb639d9025be99d 100644 (file)
--- a/chimeraslayer.cpp
+++ b/chimeraslayer.cpp
@@ -773,7 +773,7 @@ int ChimeraSlayer::getChimeras(Sequence* query) {
                for (int i = 0; i < refSeqs.size(); i++) {  delete refSeqs[i];  }
                
                if (chimeraFlag == "yes") {
-
+                       
                        if (realign) {
                                vector<string> parents;
                                for (int i = 0; i < Results.size(); i++) {
@@ -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;  }
@@ -1109,21 +1113,24 @@ vector<Sequence*> ChimeraSlayer::getBlastSeqs(Sequence* q, vector<Sequence*>& db
                                seen[tempIndexesRight[i]] = tempIndexesRight[i];
                        }
                }
-               
+               //string qname = q->getName().substr(0, q->getName().find_last_of('_'));        
+               //cout << qname << endl;        
                
                for (int i = 0; i < mergedResults.size(); i++) {
-                       //cout << mergedResults[i]  << '\t' << db[mergedResults[i]]->getName() << endl; 
+                       //cout << q->getName() << mergedResults[i]  << '\t' << db[mergedResults[i]]->getName() << endl; 
                        if (db[mergedResults[i]]->getName() != q->getName()) { 
                                Sequence* temp = new Sequence(db[mergedResults[i]]->getName(), db[mergedResults[i]]->getAligned());
                                refResults.push_back(temp);
                                
                        }
                }
-               
+               //cout << endl << endl;
 
                delete queryRight;
                delete queryLeft;
                
+               if (refResults.size() == 0) { m->mothurOut("[WARNING]: mothur found 0 potential parents, so we are not able to check " + q->getName() + ". This could be due to formatdb.exe not being setup properly, please check formatdb.log for errors."); m->mothurOutEndLine(); }
+               
                return refResults;
        }
        catch(exception& e) {
@@ -1150,17 +1157,44 @@ vector<Sequence*> ChimeraSlayer::getKmerSeqs(Sequence* q, vector<Sequence*>& db,
                //merge results         
                map<int, int> seen;
                map<int, int>::iterator it;
-                       vector<int> mergedResults;
-               for (int i = 0; i < tempIndexesLeft.size(); i++) {
+               vector<int> mergedResults;
+               
+               int index = 0;
+               //              for (int i = 0; i < smaller.size(); i++) {
+               while(index < tempIndexesLeft.size() && index < tempIndexesRight.size()){
                        
                        if (m->control_pressed) { delete queryRight; delete queryLeft; return refResults; }
                        
                        //add left if you havent already
+                       it = seen.find(tempIndexesLeft[index]);
+                       if (it == seen.end()) {  
+                               mergedResults.push_back(tempIndexesLeft[index]);
+                               seen[tempIndexesLeft[index]] = tempIndexesLeft[index];
+                       }
+                       
+                       //add right if you havent already
+                       it = seen.find(tempIndexesRight[index]);
+                       if (it == seen.end()) {  
+                               mergedResults.push_back(tempIndexesRight[index]);
+                               seen[tempIndexesRight[index]] = tempIndexesRight[index];
+                       }
+                       index++;
+               }
+               
+               
+               for (int i = index; i < tempIndexesLeft.size(); i++) {
+                       if (m->control_pressed) { delete queryRight; delete queryLeft; return refResults; }
+                       
+                       //add right if you havent already
                        it = seen.find(tempIndexesLeft[i]);
                        if (it == seen.end()) {  
                                mergedResults.push_back(tempIndexesLeft[i]);
                                seen[tempIndexesLeft[i]] = tempIndexesLeft[i];
                        }
+               }
+               
+               for (int i = index; i < tempIndexesRight.size(); i++) {
+                       if (m->control_pressed) { delete queryRight; delete queryLeft; return refResults; }
                        
                        //add right if you havent already
                        it = seen.find(tempIndexesRight[i]);
@@ -1170,17 +1204,15 @@ vector<Sequence*> ChimeraSlayer::getKmerSeqs(Sequence* q, vector<Sequence*>& db,
                        }
                }
                
-               //numWanted = mergedResults.size();
-                       
-               //cout << q->getName() << endl;         
-               
                for (int i = 0; i < mergedResults.size(); i++) {
-                       //cout << db[mergedResults[i]]->getName() << endl;      
+                       //cout << mergedResults[i]  << '\t' << db[mergedResults[i]]->getName() << endl; 
                        if (db[mergedResults[i]]->getName() != q->getName()) { 
                                Sequence* temp = new Sequence(db[mergedResults[i]]->getName(), db[mergedResults[i]]->getAligned());
                                refResults.push_back(temp);
+                               
                        }
                }
+
                //cout << endl;         
                delete queryRight;
                delete queryLeft;