working on pca

[mothur.git] / seqerrorcommand.cpp

diff --git a/seqerrorcommand.cpp b/seqerrorcommand.cpp
index 2e20c57e6245cc66f4b4f9483f643aacd6c1201a..5a0cc02fab59de9960e6626b4cbcc87e8ce82d1f 100644 (file)
--- a/seqerrorcommand.cpp
+++ b/seqerrorcommand.cpp
@@ -15,7 +15,7 @@
 //**********************************************************************************************************************
 vector<string> SeqErrorCommand::getValidParameters(){  
        try {
-               string Array[] =  {"query", "reference", "name", "qfile", "report", "threshold", "inputdir", "outputdir"};
+               string Array[] =  {"query", "reference", "name", "qfile", "report", "threshold", "ignorechimeras", "inputdir", "outputdir"};
                vector<string> myArray (Array, Array+(sizeof(Array)/sizeof(string)));
                return myArray;
        }
@@ -230,6 +230,8 @@ int SeqErrorCommand::execute(){
        try{
                if (abort == true) { if (calledHelp) { return 0; }  return 2;   }
 
+               maxLength = 2000;
+               
                string errorSummaryFileName = queryFileName.substr(0,queryFileName.find_last_of('.')) + ".error.summary";
                m->openOutputFile(errorSummaryFileName, errorSummaryFile);
                outputNames.push_back(errorSummaryFileName); outputTypes["error.summary"].push_back(errorSummaryFileName);
@@ -259,11 +261,11 @@ int SeqErrorCommand::execute(){
                        m->openInputFile(qualFileName, qualFile);
                        report = ReportFile(reportFile, reportFileName);
                        
-                       qualForwardMap.resize(1000);
-                       qualReverseMap.resize(1000);
-                       for(int i=0;i<1000;i++){
-                               qualForwardMap[i].assign(100,0);
-                               qualReverseMap[i].assign(100,0);
+                       qualForwardMap.resize(maxLength);
+                       qualReverseMap.resize(maxLength);
+                       for(int i=0;i<maxLength;i++){
+                               qualForwardMap[i].assign(41,0);
+                               qualReverseMap[i].assign(41,0);
                        }                               
                }
                
@@ -282,18 +284,18 @@ int SeqErrorCommand::execute(){
                qScoreErrorMap['a'].assign(41, 0);
                
                map<char, vector<int> > errorForward;
-               errorForward['m'].assign(1000,0);
-               errorForward['s'].assign(1000,0);
-               errorForward['i'].assign(1000,0);
-               errorForward['d'].assign(1000,0);
-               errorForward['a'].assign(1000,0);
+               errorForward['m'].assign(maxLength,0);
+               errorForward['s'].assign(maxLength,0);
+               errorForward['i'].assign(maxLength,0);
+               errorForward['d'].assign(maxLength,0);
+               errorForward['a'].assign(maxLength,0);
                
                map<char, vector<int> > errorReverse;
-               errorReverse['m'].assign(1000,0);
-               errorReverse['s'].assign(1000,0);
-               errorReverse['i'].assign(1000,0);
-               errorReverse['d'].assign(1000,0);
-               errorReverse['a'].assign(1000,0);       
+               errorReverse['m'].assign(maxLength,0);
+               errorReverse['s'].assign(maxLength,0);
+               errorReverse['i'].assign(maxLength,0);
+               errorReverse['d'].assign(maxLength,0);
+               errorReverse['a'].assign(maxLength,0);  
                
 
                string errorChimeraFileName = queryFileName.substr(0,queryFileName.find_last_of('.')) + ".error.chimera";
@@ -367,7 +369,8 @@ int SeqErrorCommand::execute(){
                        }
 
                        index++;
-                       if(index % 1000 == 0){  cout << index << endl;  }
+                       
+                       if(index % 1000 == 0){  m->mothurOut(toString(index) + '\n');   }
                }
                queryFile.close();
                errorSummaryFile.close();       
@@ -715,7 +718,7 @@ void SeqErrorCommand::printErrorFRFile(map<char, vector<int> > errorForward, map
                outputNames.push_back(errorForwardFileName);  outputTypes["error.forward"].push_back(errorForwardFileName);
 
                errorForwardFile << "position\ttotalseqs\tmatch\tsubstitution\tinsertion\tdeletion\tambiguous" << endl;
-               for(int i=0;i<1000;i++){
+               for(int i=0;i<maxLength;i++){
                        float match = (float)errorForward['m'][i];
                        float subst = (float)errorForward['s'][i];
                        float insert = (float)errorForward['i'][i];
@@ -733,7 +736,7 @@ void SeqErrorCommand::printErrorFRFile(map<char, vector<int> > errorForward, map
                outputNames.push_back(errorReverseFileName);  outputTypes["error.reverse"].push_back(errorReverseFileName);
 
                errorReverseFile << "position\ttotalseqs\tmatch\tsubstitution\tinsertion\tdeletion\tambiguous" << endl;
-               for(int i=0;i<1000;i++){
+               for(int i=0;i<maxLength;i++){
                        float match = (float)errorReverse['m'][i];
                        float subst = (float)errorReverse['s'][i];
                        float insert = (float)errorReverse['i'][i];
@@ -777,17 +780,15 @@ void SeqErrorCommand::printErrorQuality(map<char, vector<int> > qScoreErrorMap){
 //***************************************************************************************************************
 
 void SeqErrorCommand::printQualityFR(vector<vector<int> > qualForwardMap, vector<vector<int> > qualReverseMap){
-       try{
-
 
-               int lastRow = 0;
-               int lastColumn = 0;
+       try{
+               int numRows = 0;
+               int numColumns = qualForwardMap[0].size();
 
                for(int i=0;i<qualForwardMap.size();i++){
-                       for(int j=0;j<qualForwardMap[i].size();j++){
+                       for(int j=0;j<numColumns;j++){
                                if(qualForwardMap[i][j] != 0){
-                                       if(lastRow < i)         {       lastRow = i+2;          }
-                                       if(lastColumn < j)      {       lastColumn = j+2;       }
+                                       if(numRows < i)         {       numRows = i+20;         }
                                }
                        }
                }
@@ -797,11 +798,11 @@ void SeqErrorCommand::printQualityFR(vector<vector<int> > qualForwardMap, vector
                m->openOutputFile(qualityForwardFileName, qualityForwardFile);
                outputNames.push_back(qualityForwardFileName);  outputTypes["error.qual.forward"].push_back(qualityForwardFileName);
 
-               for(int i=0;i<lastColumn;i++){  qualityForwardFile << '\t' << i;        }       qualityForwardFile << endl;
+               for(int i=0;i<numColumns;i++){  qualityForwardFile << '\t' << i;        }       qualityForwardFile << endl;
 
-               for(int i=0;i<lastRow;i++){
+               for(int i=0;i<numRows;i++){
                        qualityForwardFile << i+1;
-                       for(int j=0;j<lastColumn;j++){
+                       for(int j=0;j<numColumns;j++){
                                qualityForwardFile << '\t' << qualForwardMap[i][j];
                        }
 
@@ -815,11 +816,11 @@ void SeqErrorCommand::printQualityFR(vector<vector<int> > qualForwardMap, vector
                m->openOutputFile(qualityReverseFileName, qualityReverseFile);
                outputNames.push_back(qualityReverseFileName);  outputTypes["error.qual.reverse"].push_back(qualityReverseFileName);
                
-               for(int i=0;i<lastColumn;i++){  qualityReverseFile << '\t' << i;        }       qualityReverseFile << endl;
-               for(int i=0;i<lastRow;i++){
+               for(int i=0;i<numColumns;i++){  qualityReverseFile << '\t' << i;        }       qualityReverseFile << endl;
+               for(int i=0;i<numRows;i++){
                        
                        qualityReverseFile << i+1;
-                       for(int j=0;j<lastColumn;j++){
+                       for(int j=0;j<numColumns;j++){
                                qualityReverseFile << '\t' << qualReverseMap[i][j];
                        }
                        qualityReverseFile << endl;
@@ -830,7 +831,7 @@ void SeqErrorCommand::printQualityFR(vector<vector<int> > qualForwardMap, vector
                m->errorOut(e, "SeqErrorCommand", "printErrorFRFile");
                exit(1);
        }
+       
 }
 
-
 //***************************************************************************************************************