]> git.donarmstrong.com Git - mothur.git/blobdiff - filterseqscommand.cpp
changed random forest output filename
[mothur.git] / filterseqscommand.cpp
index e96b7fd3762fa69c5930066729e867ff69705b3d..8b68988884acb408d0e72c5c3fa407a7b9f241b6 100644 (file)
 //**********************************************************************************************************************
 vector<string> FilterSeqsCommand::setParameters(){     
        try {
-               CommandParameter pfasta("fasta", "InputTypes", "", "", "none", "none", "none",false,true); parameters.push_back(pfasta);
-               CommandParameter phard("hard", "InputTypes", "", "", "none", "none", "none",false,false); parameters.push_back(phard);
-               CommandParameter ptrump("trump", "String", "", "*", "", "", "",false,false); parameters.push_back(ptrump);
-               CommandParameter psoft("soft", "Number", "", "0", "", "", "",false,false); parameters.push_back(psoft);
-               CommandParameter pvertical("vertical", "Boolean", "", "T", "", "", "",false,false); parameters.push_back(pvertical);
-               CommandParameter pprocessors("processors", "Number", "", "1", "", "", "",false,false); parameters.push_back(pprocessors);
-               CommandParameter pinputdir("inputdir", "String", "", "", "", "", "",false,false); parameters.push_back(pinputdir);
-               CommandParameter poutputdir("outputdir", "String", "", "", "", "", "",false,false); parameters.push_back(poutputdir);
+               CommandParameter pfasta("fasta", "InputTypes", "", "", "none", "none", "none","fasta-filter",false,true, true); parameters.push_back(pfasta);
+               CommandParameter phard("hard", "InputTypes", "", "", "none", "none", "none","",false,false); parameters.push_back(phard);
+               CommandParameter ptrump("trump", "String", "", "*", "", "", "","",false,false, true); parameters.push_back(ptrump);
+               CommandParameter psoft("soft", "Number", "", "0", "", "", "","",false,false); parameters.push_back(psoft);
+               CommandParameter pvertical("vertical", "Boolean", "", "T", "", "", "","",false,false, true); parameters.push_back(pvertical);
+               CommandParameter pprocessors("processors", "Number", "", "1", "", "", "","",false,false, true); parameters.push_back(pprocessors);
+               CommandParameter pinputdir("inputdir", "String", "", "", "", "", "","",false,false); parameters.push_back(pinputdir);
+               CommandParameter poutputdir("outputdir", "String", "", "", "", "", "","",false,false); parameters.push_back(poutputdir);
                
                vector<string> myArray;
                for (int i = 0; i < parameters.size(); i++) {   myArray.push_back(parameters[i].name);          }
@@ -57,6 +57,22 @@ string FilterSeqsCommand::getHelpString(){
        }
 }
 //**********************************************************************************************************************
+string FilterSeqsCommand::getOutputPattern(string type) {
+    try {
+        string pattern = "";
+        
+        if (type == "fasta") {  pattern = "[filename],filter.fasta"; } 
+        else if (type == "filter") {  pattern =  "[filename],filter"; }
+        else { m->mothurOut("[ERROR]: No definition for type " + type + " output pattern.\n"); m->control_pressed = true;  }
+        
+        return pattern;
+    }
+    catch(exception& e) {
+        m->errorOut(e, "FilterSeqsCommand", "getOutputPattern");
+        exit(1);
+    }
+}
+//**********************************************************************************************************************
 FilterSeqsCommand::FilterSeqsCommand(){        
        try {
                abort = true; calledHelp = true; 
@@ -125,7 +141,12 @@ FilterSeqsCommand::FilterSeqsCommand(string option)  {
                        fasta = validParameter.validFile(parameters, "fasta", false);
                        if (fasta == "not found") {                             
                                fasta = m->getFastaFile(); 
-                               if (fasta != "") { fastafileNames.push_back(fasta);  m->mothurOut("Using " + fasta + " as input file for the fasta parameter."); m->mothurOutEndLine(); }
+                               if (fasta != "") { 
+                    fastafileNames.push_back(fasta);  
+                    m->mothurOut("Using " + fasta + " as input file for the fasta parameter."); m->mothurOutEndLine();
+                    string simpleName = m->getSimpleName(fasta);
+                    filterFileName += simpleName.substr(0, simpleName.find_first_of('.'));
+                }
                                else {  m->mothurOut("You have no current fastafile and the fasta parameter is required."); m->mothurOutEndLine(); abort = true; }
                        }
                        else { 
@@ -222,7 +243,7 @@ FilterSeqsCommand::FilterSeqsCommand(string option)  {
                        
                        temp = validParameter.validFile(parameters, "processors", false);       if (temp == "not found"){       temp = m->getProcessors();      }
                        m->setProcessors(temp);
-                       convert(temp, processors); 
+                       m->mothurConvert(temp, processors); 
                        
                        vertical = validParameter.validFile(parameters, "vertical", false);             
                        if (vertical == "not found") { 
@@ -272,9 +293,10 @@ int FilterSeqsCommand::execute() {
                ofstream outFilter;
                
                //prevent giantic file name
-               string filterFile;
-               if (fastafileNames.size() > 3) { filterFile = outputDir + "merge.filter"; }
-               else {  filterFile = outputDir + filterFileName + ".filter";  }
+        map<string, string> variables; 
+        variables["[filename]"] = outputDir + filterFileName + ".";
+               if (fastafileNames.size() > 3) { variables["[filename]"] = outputDir + "merge."; }
+               string filterFile = getOutputFileName("filter", variables);  
                
                m->openOutputFile(filterFile, outFilter);
                outFilter << filter << endl;
@@ -337,7 +359,9 @@ int FilterSeqsCommand::filterSequences() {
                        
                                for (int i = 0; i < lines.size(); i++) {  delete lines[i];  }  lines.clear();
                                
-                               string filteredFasta = outputDir + m->getRootName(m->getSimpleName(fastafileNames[s])) + "filter.fasta";
+                map<string, string> variables; 
+                variables["[filename]"] = outputDir + m->getRootName(m->getSimpleName(fastafileNames[s]));
+                               string filteredFasta = getOutputFileName("fasta", variables);
 #ifdef USE_MPI 
                                int pid, numSeqsPerProcessor, num; 
                                int tag = 2001;
@@ -420,32 +444,56 @@ int FilterSeqsCommand::filterSequences() {
                                MPI_Barrier(MPI_COMM_WORLD); //make everyone wait - just in case
                                
 #else
-                       vector<unsigned long long> positions = m->divideFile(fastafileNames[s], processors);
-                               
+            
+            vector<unsigned long long> positions;
+            if (savedPositions.size() != 0) { positions = savedPositions[s]; }
+            else {
+#if defined (__APPLE__) || (__MACH__) || (linux) || (__linux) || (__linux__) || (__unix__) || (__unix)
+                               positions = m->divideFile(fastafileNames[s], processors);
+#else
+                if(processors != 1){
+                    int numFastaSeqs = 0;
+                    positions = m->setFilePosFasta(fastafileNames[s], numFastaSeqs); 
+                    if (positions.size() < processors) { processors = positions.size(); }
+                }
+#endif
+            }
+               #if defined (__APPLE__) || (__MACH__) || (linux) || (__linux) || (__linux__) || (__unix__) || (__unix)
+                       //vector<unsigned long long> positions = m->divideFile(fastafileNames[s], processors);
+                       
                        for (int i = 0; i < (positions.size()-1); i++) {
                                lines.push_back(new linePair(positions[i], positions[(i+1)]));
                        }       
-               #if defined (__APPLE__) || (__MACH__) || (linux) || (__linux)
+                       
                                if(processors == 1){
                                        int numFastaSeqs = driverRunFilter(filter, filteredFasta, fastafileNames[s], lines[0]);
                                        numSeqs += numFastaSeqs;
                                }else{
-                                       int numFastaSeqs = createProcessesRunFilter(filter, fastafileNames[s]); 
+                                       int numFastaSeqs = createProcessesRunFilter(filter, fastafileNames[s], filteredFasta); 
                                        numSeqs += numFastaSeqs;
-                               
-                                       rename((fastafileNames[s] + toString(processIDS[0]) + ".temp").c_str(), filteredFasta.c_str());
-                               
-                                       //append fasta files
-                                       for(int i=1;i<processors;i++){
-                                               m->appendFiles((fastafileNames[s] + toString(processIDS[i]) + ".temp"), filteredFasta);
-                                               m->mothurRemove((fastafileNames[s] + toString(processIDS[i]) + ".temp"));
-                                       }
                                }
                                
                                if (m->control_pressed) {  return 1; }
                #else
+            if(processors == 1){
+                lines.push_back(new linePair(0, 1000));
                                int numFastaSeqs = driverRunFilter(filter, filteredFasta, fastafileNames[s], lines[0]);
                                numSeqs += numFastaSeqs;
+            }else {
+                int numFastaSeqs = positions.size()-1;
+                //positions = m->setFilePosFasta(fastafileNames[s], numFastaSeqs); 
+                
+                //figure out how many sequences you have to process
+                int numSeqsPerProcessor = numFastaSeqs / processors;
+                for (int i = 0; i < processors; i++) {
+                    int startIndex =  i * numSeqsPerProcessor;
+                    if(i == (processors - 1)){ numSeqsPerProcessor = numFastaSeqs - i * numSeqsPerProcessor;   }
+                    lines.push_back(new linePair(positions[startIndex], numSeqsPerProcessor));
+                }
+                
+                numFastaSeqs = createProcessesRunFilter(filter, fastafileNames[s], filteredFasta); 
+                numSeqs += numFastaSeqs;
+            }
 
                                if (m->control_pressed) {  return 1; }
                #endif
@@ -510,7 +558,7 @@ int FilterSeqsCommand::driverMPIRun(int start, int num, MPI_File& inMPI, MPI_Fil
 
                        }
                        
-                       if((i+1) % 100 == 0){   cout << (i+1) << endl;   m->mothurOutJustToLog(toString(i+1) + "\n");   }
+                       if((i+1) % 100 == 0){   cout << (i+1) << endl;          }
                }
                
                if(outputString != ""){ //output to file 
@@ -524,7 +572,7 @@ int FilterSeqsCommand::driverMPIRun(int start, int num, MPI_File& inMPI, MPI_Fil
                        delete buf;
                }
                
-               if((num) % 100 != 0){   cout << (num) << endl;   m->mothurOutJustToLog(toString(num) + "\n");   }
+               if((num) % 100 != 0){   cout << (num) << endl;          }
                        
                return 0;
        }
@@ -567,7 +615,7 @@ int FilterSeqsCommand::driverRunFilter(string F, string outputFilename, string i
                                count++;
                        }
                        
-                       #if defined (__APPLE__) || (__MACH__) || (linux) || (__linux)
+                       #if defined (__APPLE__) || (__MACH__) || (linux) || (__linux) || (__linux__) || (__unix__) || (__unix)
                                unsigned long long pos = in.tellg();
                                if ((pos == -1) || (pos >= filePos->end)) { break; }
                        #else
@@ -575,10 +623,10 @@ int FilterSeqsCommand::driverRunFilter(string F, string outputFilename, string i
                        #endif
                        
                        //report progress
-                       if((count) % 100 == 0){ m->mothurOut(toString(count)); m->mothurOutEndLine();           }
+                       if((count) % 100 == 0){ m->mothurOutJustToScreen(toString(count)+"\n");         }
                }
                //report progress
-               if((count) % 100 != 0){ m->mothurOut(toString(count)); m->mothurOutEndLine();           }
+               if((count) % 100 != 0){ m->mothurOutJustToScreen(toString(count)+"\n");                 }
                
                
                out.close();
@@ -593,12 +641,15 @@ int FilterSeqsCommand::driverRunFilter(string F, string outputFilename, string i
 }
 /**************************************************************************************************/
 
-int FilterSeqsCommand::createProcessesRunFilter(string F, string filename) {
+int FilterSeqsCommand::createProcessesRunFilter(string F, string filename, string filteredFastaName) {
        try {
-#if defined (__APPLE__) || (__MACH__) || (linux) || (__linux)
-               int process = 0;
+        
+        int process = 1;
                int num = 0;
                processIDS.clear();
+        
+#if defined (__APPLE__) || (__MACH__) || (linux) || (__linux) || (__linux__) || (__unix__) || (__unix)
+               
                
                //loop through and create all the processes you want
                while (process != processors) {
@@ -626,8 +677,10 @@ int FilterSeqsCommand::createProcessesRunFilter(string F, string filename) {
                        }
                }
                
+        num = driverRunFilter(F, filteredFastaName, filename, lines[0]);
+        
                //force parent to wait until all the processes are done
-               for (int i=0;i<processors;i++) { 
+               for (int i=0;i<processIDS.size();i++) { 
                        int temp = processIDS[i];
                        wait(&temp);
                }       
@@ -638,11 +691,59 @@ int FilterSeqsCommand::createProcessesRunFilter(string F, string filename) {
                        m->openInputFile(tempFile, in);
                        if (!in.eof()) { int tempNum = 0; in >> tempNum; num += tempNum; }
                        in.close(); m->mothurRemove(tempFile);
+            
+            m->appendFiles((filename + toString(processIDS[i]) + ".temp"), filteredFastaName);
+            m->mothurRemove((filename + toString(processIDS[i]) + ".temp"));
                }
-
+               
+#else
+        
+        //////////////////////////////////////////////////////////////////////////////////////////////////////
+               //Windows version shared memory, so be careful when passing variables through the filterData struct. 
+               //Above fork() will clone, so memory is separate, but that's not the case with windows, 
+               //Taking advantage of shared memory to allow both threads to add info to F.
+               //////////////////////////////////////////////////////////////////////////////////////////////////////
+               
+               vector<filterRunData*> pDataArray; 
+               DWORD   dwThreadIdArray[processors-1];
+               HANDLE  hThreadArray[processors-1]; 
+               
+               //Create processor worker threads.
+               for( int i=0; i<processors-1; i++){
+                       
+            string extension = "";
+                       if (i != 0) { extension = toString(i) + ".temp"; }
+            
+                       filterRunData* tempFilter = new filterRunData(filter, filename, (filteredFastaName + extension), m, lines[i]->start, lines[i]->end, alignmentLength, i);
+                       pDataArray.push_back(tempFilter);
+                       processIDS.push_back(i);
+            
+                       hThreadArray[i] = CreateThread(NULL, 0, MyRunFilterThreadFunction, pDataArray[i], 0, &dwThreadIdArray[i]);   
+               }
+        
+        num = driverRunFilter(F, (filteredFastaName + toString(processors-1) + ".temp"), filename, lines[processors-1]);
+        
+               //Wait until all threads have terminated.
+               WaitForMultipleObjects(processors-1, hThreadArray, TRUE, INFINITE);
                
-               return num;
-#endif         
+               //Close all thread handles and free memory allocations.
+               for(int i=0; i < pDataArray.size(); i++){
+                       num += pDataArray[i]->count;
+            if (pDataArray[i]->count != pDataArray[i]->end) {
+                m->mothurOut("[ERROR]: process " + toString(i) + " only processed " + toString(pDataArray[i]->count) + " of " + toString(pDataArray[i]->end) + " sequences assigned to it, quitting. \n"); m->control_pressed = true; 
+            }
+            CloseHandle(hThreadArray[i]);
+                       delete pDataArray[i];
+               }
+        
+        for (int i = 1; i < processors; i++) {
+            m->appendFiles((filteredFastaName + toString(i) + ".temp"), filteredFastaName);
+            m->mothurRemove((filteredFastaName + toString(i) + ".temp"));
+               }
+#endif 
+        
+        return num;
+        
        }
        catch(exception& e) {
                m->errorOut(e, "FilterSeqsCommand", "createProcessesRunFilter");
@@ -736,12 +837,14 @@ string FilterSeqsCommand::createFilter() {
                                MPI_Barrier(MPI_COMM_WORLD); //make everyone wait - just in case
                                
 #else
-               vector<unsigned long long> positions = m->divideFile(fastafileNames[s], processors);
                                
-               for (int i = 0; i < (positions.size()-1); i++) {
-                       lines.push_back(new linePair(positions[i], positions[(i+1)]));
-               }       
-               #if defined (__APPLE__) || (__MACH__) || (linux) || (__linux)
+                vector<unsigned long long> positions;
+               #if defined (__APPLE__) || (__MACH__) || (linux) || (__linux) || (__linux__) || (__unix__) || (__unix)
+                               positions = m->divideFile(fastafileNames[s], processors);
+                               for (int i = 0; i < (positions.size()-1); i++) {
+                                       lines.push_back(new linePair(positions[i], positions[(i+1)]));
+                               }       
+                               
                                if(processors == 1){
                                        int numFastaSeqs = driverCreateFilter(F, fastafileNames[s], lines[0]);
                                        numSeqs += numFastaSeqs;
@@ -749,13 +852,32 @@ string FilterSeqsCommand::createFilter() {
                                        int numFastaSeqs = createProcessesCreateFilter(F, fastafileNames[s]); 
                                        numSeqs += numFastaSeqs;
                                }
-                               
-                               if (m->control_pressed) {  return filterString; }
                #else
-                               int numFastaSeqs = driverCreateFilter(F, fastafileNames[s], lines[0]);
-                               numSeqs += numFastaSeqs;
-                               if (m->control_pressed) {  return filterString; }
+                if(processors == 1){
+                    lines.push_back(new linePair(0, 1000));
+                    int numFastaSeqs = driverCreateFilter(F, fastafileNames[s], lines[0]);
+                    numSeqs += numFastaSeqs;
+                               }else {
+                    int numFastaSeqs = 0;
+                    positions = m->setFilePosFasta(fastafileNames[s], numFastaSeqs); 
+                    if (positions.size() < processors) { processors = positions.size(); }
+                    
+                    //figure out how many sequences you have to process
+                    int numSeqsPerProcessor = numFastaSeqs / processors;
+                    for (int i = 0; i < processors; i++) {
+                        int startIndex =  i * numSeqsPerProcessor;
+                        if(i == (processors - 1)){     numSeqsPerProcessor = numFastaSeqs - i * numSeqsPerProcessor;   }
+                        lines.push_back(new linePair(positions[startIndex], numSeqsPerProcessor));
+                    }
+                    
+                    numFastaSeqs = createProcessesCreateFilter(F, fastafileNames[s]); 
+                                       numSeqs += numFastaSeqs;
+                }
                #endif
+                //save the file positions so we can reuse them in the runFilter function
+                savedPositions[s] = positions;
+                
+                               if (m->control_pressed) {  return filterString; }
 #endif
                        
                        }
@@ -842,7 +964,7 @@ string FilterSeqsCommand::createFilter() {
        
        MPI_Barrier(MPI_COMM_WORLD);
 #endif
-                               
+            
                return filterString;
        }
        catch(exception& e) {
@@ -876,7 +998,7 @@ int FilterSeqsCommand::driverCreateFilter(Filters& F, string filename, linePair*
                                        count++;
                        }
                        
-                       #if defined (__APPLE__) || (__MACH__) || (linux) || (__linux)
+                       #if defined (__APPLE__) || (__MACH__) || (linux) || (__linux) || (__linux__) || (__unix__) || (__unix)
                                unsigned long long pos = in.tellg();
                                if ((pos == -1) || (pos >= filePos->end)) { break; }
                        #else
@@ -884,10 +1006,10 @@ int FilterSeqsCommand::driverCreateFilter(Filters& F, string filename, linePair*
                        #endif
                        
                        //report progress
-                       if((count) % 100 == 0){ m->mothurOut(toString(count)); m->mothurOutEndLine();           }
+                       if((count) % 100 == 0){ m->mothurOutJustToScreen(toString(count)+"\n");                 }
                }
                //report progress
-               if((count) % 100 != 0){ m->mothurOut(toString(count)); m->mothurOutEndLine();           }
+               if((count) % 100 != 0){ m->mothurOutJustToScreen(toString(count)+"\n");         }
                in.close();
                
                return count;
@@ -930,11 +1052,11 @@ int FilterSeqsCommand::MPICreateFilter(int start, int num, Filters& F, MPI_File&
                        cout.flush();
                                                
                        //report progress
-                       if((i+1) % 100 == 0){   cout << (i+1) << endl;   m->mothurOutJustToLog(toString(i+1) + "\n");   }
+                       if((i+1) % 100 == 0){   cout << (i+1) << endl;          }
                }
                
                //report progress
-               if((num) % 100 != 0){   cout << num << endl; m->mothurOutJustToLog(toString(num) + "\n");       }
+               if((num) % 100 != 0){   cout << num << endl;    }
                
                return 0;
        }
@@ -948,11 +1070,12 @@ int FilterSeqsCommand::MPICreateFilter(int start, int num, Filters& F, MPI_File&
 
 int FilterSeqsCommand::createProcessesCreateFilter(Filters& F, string filename) {
        try {
-#if defined (__APPLE__) || (__MACH__) || (linux) || (__linux)
-               int process = 1;
+        int process = 1;
                int num = 0;
                processIDS.clear();
-               
+
+#if defined (__APPLE__) || (__MACH__) || (linux) || (__linux) || (__linux__) || (__unix__) || (__unix)
+                               
                //loop through and create all the processes you want
                while (process != processors) {
                        int pid = fork();
@@ -1027,8 +1150,53 @@ int FilterSeqsCommand::createProcessesCreateFilter(Filters& F, string filename)
                        m->mothurRemove(tempFilename);
                }
                
-               return num;
-#endif         
+               
+#else
+        
+        //////////////////////////////////////////////////////////////////////////////////////////////////////
+               //Windows version shared memory, so be careful when passing variables through the filterData struct. 
+               //Above fork() will clone, so memory is separate, but that's not the case with windows, 
+               //Taking advantage of shared memory to allow both threads to add info to F.
+               //////////////////////////////////////////////////////////////////////////////////////////////////////
+               
+               vector<filterData*> pDataArray; 
+               DWORD   dwThreadIdArray[processors];
+               HANDLE  hThreadArray[processors]; 
+               
+               //Create processor worker threads.
+               for( int i=0; i<processors; i++ ){
+                       
+                       filterData* tempFilter = new filterData(filename, m, lines[i]->start, lines[i]->end, alignmentLength, trump, vertical, soft, hard, i);
+                       pDataArray.push_back(tempFilter);
+                       processIDS.push_back(i);
+            
+                       hThreadArray[i] = CreateThread(NULL, 0, MyCreateFilterThreadFunction, pDataArray[i], 0, &dwThreadIdArray[i]);   
+               }
+        
+               //Wait until all threads have terminated.
+               WaitForMultipleObjects(processors, hThreadArray, TRUE, INFINITE);
+               
+               //Close all thread handles and free memory allocations.
+               for(int i=0; i < pDataArray.size(); i++){
+                       num += pDataArray[i]->count;
+            if (pDataArray[i]->count != pDataArray[i]->end) {
+                m->mothurOut("[ERROR]: process " + toString(i) + " only processed " + toString(pDataArray[i]->count) + " of " + toString(pDataArray[i]->end) + " sequences assigned to it, quitting. \n"); m->control_pressed = true; 
+            }
+            F.mergeFilter(pDataArray[i]->F.getFilter());
+            
+                       for (int k = 0; k < alignmentLength; k++) {      F.a[k] += pDataArray[i]->F.a[k];       }
+                       for (int k = 0; k < alignmentLength; k++) {      F.t[k] += pDataArray[i]->F.t[k];       }
+                       for (int k = 0; k < alignmentLength; k++) {      F.g[k] += pDataArray[i]->F.g[k];       }
+                       for (int k = 0; k < alignmentLength; k++) {      F.c[k] += pDataArray[i]->F.c[k];       }
+                       for (int k = 0; k < alignmentLength; k++) {      F.gap[k] += pDataArray[i]->F.gap[k];   }
+
+                       CloseHandle(hThreadArray[i]);
+                       delete pDataArray[i];
+               }
+               
+#endif 
+        return num;
+        
        }
        catch(exception& e) {
                m->errorOut(e, "FilterSeqsCommand", "createProcessesCreateFilter");