]> git.donarmstrong.com Git - mothur.git/blobdiff - trimseqscommand.cpp
added check for blank files in openInputFile, add check for groupfiles that contain...
[mothur.git] / trimseqscommand.cpp
index 937b130deb2c4b2dd6d3e72d191c090b34ee51e5..455cdc29aa47ecb1ca26a7f9b4021bef151d6314 100644 (file)
@@ -160,7 +160,7 @@ int TrimSeqsCommand::execute(){
                if(qFileName != "")     {       openInputFile(qFileName, qFile);        }
                
                bool success;
-               
+                       
                while(!inFASTA.eof()){
                        Sequence currSeq(inFASTA);
                        string origSeq = currSeq.getUnaligned();
@@ -171,9 +171,9 @@ int TrimSeqsCommand::execute(){
                                if(qThreshold != 0)             {       success = stripQualThreshold(currSeq, qFile);   }
                                else if(qAverage != 0)  {       success = cullQualAverage(currSeq, qFile);              }
                                if(!success)                    {       trashCode += 'q';                                                               }
-                               qFile.close();
                        }
                        if(barcodes.size() != 0){
+       
                                success = stripBarcode(currSeq, group);
                                if(!success){   trashCode += 'b';       }
                        }
@@ -187,7 +187,8 @@ int TrimSeqsCommand::execute(){
                        }
                        if(minLength > 0 || maxLength > 0){
                                success = cullLength(currSeq);
-                               if(!success){   trashCode += 'l';       }
+                       if ((currSeq.getUnaligned().length() > 300) && (success)) {  cout << "too long " << currSeq.getUnaligned().length() << endl;  }
+                               if(!success){   trashCode += 'l'; }
                        }
                        if(maxHomoP > 0){
                                success = cullHomoP(currSeq);
@@ -201,6 +202,7 @@ int TrimSeqsCommand::execute(){
                        if(flip){       currSeq.reverseComplement();    }               // should go last                       
                        
                        if(trashCode.length() == 0){
+                               currSeq.setAligned(currSeq.getUnaligned());  //this is because of a modification we made to the sequence class to fix a bug.  all seqs have an aligned version, which is the version that gets printed.
                                currSeq.printSequence(outFASTA);
                                if(barcodes.size() != 0){
                                        outGroups << currSeq.getName() << '\t' << groupVector[group] << endl;
@@ -221,6 +223,7 @@ int TrimSeqsCommand::execute(){
                outFASTA.close();
                scrapFASTA.close();
                outGroups.close();
+               if(qFileName != "")     {       qFile.close();  }
                
                for(int i=0;i<fastaFileNames.size();i++){
                        fastaFileNames[i]->close();
@@ -548,7 +551,8 @@ bool TrimSeqsCommand::cullQualAverage(Sequence& seq, ifstream& qFile){
                string name;
                
                qFile >> name;
-               if (name.length() != 0) {  if(name.substr(1) != seq.getName())  {       mothurOut("sequence name mismatch btwn fasta and qual file"); mothurOutEndLine();       } }
+               if (name[0] == '>') {  if(name.substr(1) != seq.getName())      {       mothurOut("sequence name mismatch btwn fasta: " + seq.getName() + " and qual file: " + name); mothurOutEndLine();       } }
+               
                while (!qFile.eof())    {       char c = qFile.get(); if (c == 10 || c == 13){  break;  }       }
                
                float score;    
@@ -559,7 +563,7 @@ bool TrimSeqsCommand::cullQualAverage(Sequence& seq, ifstream& qFile){
                        average += score;
                }
                average /= seqLength;
-               
+
                if(average >= qAverage) {       success = 1;    }
                else                                    {       success = 0;    }
                
@@ -572,5 +576,3 @@ bool TrimSeqsCommand::cullQualAverage(Sequence& seq, ifstream& qFile){
 }
 
 //***************************************************************************************************************
-
-