]> git.donarmstrong.com Git - mothur.git/blobdiff - distancecommand.cpp
finished mpi for filter.seqs
[mothur.git] / distancecommand.cpp
index 13f56eb90a422d43e821a9fa378f582691468cb8..4720df3cf0d968a84ba49c4af449b226a304c1cd 100644 (file)
@@ -166,7 +166,7 @@ int DistanceCommand::execute(){
                cutoff += 0.005;
                
                string outputFile;
-               
+                               
                if (output == "lt") { //does the user want lower triangle phylip formatted file 
                        outputFile = outputDir + getRootName(getSimpleName(fastafile)) + "phylip.dist";
                        remove(outputFile.c_str());
@@ -179,8 +179,60 @@ int DistanceCommand::execute(){
                        outputFile = outputDir + getRootName(getSimpleName(fastafile)) + "square.dist";
                        remove(outputFile.c_str());
                }
+               
+
+#ifdef USE_MPI
+               
+               int pid, start, end; 
+               int tag = 2001;
                                
-#if defined (__APPLE__) || (__MACH__) || (linux) || (__linux)
+               MPI_Status status; 
+               MPI_Comm_size(MPI_COMM_WORLD, &processors); //set processors to the number of mpi processes running
+               MPI_Comm_rank(MPI_COMM_WORLD, &pid); //find out who we are
+               
+               //each process gets where it should start and stop in the file
+               start = int (sqrt(float(pid)/float(processors)) * numSeqs);
+               end = int (sqrt(float(pid+1)/float(processors)) * numSeqs);
+       
+               MPI_File outMPI;
+               int amode=MPI_MODE_CREATE|MPI_MODE_WRONLY; 
+               
+               char filename[outputFile.length()];
+               strcpy(filename, outputFile.c_str());
+               
+               MPI_File_open(MPI_COMM_WORLD, filename, amode, MPI_INFO_NULL, &outMPI);
+               
+               if (pid == 0) { //you are the root process 
+               
+                       //do your part
+                       string outputMyPart;
+                       driverMPI(start, end, outMPI, cutoff);
+                       
+                       //wait on chidren
+                       for(int i = 1; i < processors; i++) { 
+                               char buf[4];
+                               MPI_Recv(buf, 4, MPI_CHAR, i, tag, MPI_COMM_WORLD, &status); 
+                       }
+                       
+                       if (output == "lt") {
+                               convertToLowerTriangle(outputFile);
+                       }
+                       
+               }else { //you are a child process
+                       //do your part
+                       driverMPI(start, end, outMPI, cutoff);
+               
+                       char buf[4];
+                       strcpy(buf, "done"); 
+                       
+                       //tell parent you are done.
+                       MPI_Send(buf, 4, MPI_CHAR, 0, tag, MPI_COMM_WORLD);
+               }
+               
+               MPI_File_close(&outMPI);
+#else          
+                               
+       #if defined (__APPLE__) || (__MACH__) || (linux) || (__linux)
                //if you don't need to fork anything
                if(processors == 1){
                        driver(0, numSeqs, outputFile, cutoff);
@@ -204,13 +256,28 @@ int DistanceCommand::execute(){
                                remove((outputFile + toString(it->second) + ".temp").c_str());
                        }
                }
-#else
+       #else
                ifstream inFASTA;
                driver(0, numSeqs, outputFile, cutoff);
+       #endif
+       
 #endif
+               if (m->control_pressed) { delete distCalculator; remove(outputFile.c_str()); return 0; }
+               
+               #ifdef USE_MPI
+                       MPI_Comm_rank(MPI_COMM_WORLD, &pid); 
+                                       
+                       if (pid == 0) { //only one process should output to screen
+               #endif
                
                if (output == "square") {  convertMatrix(outputFile); }
                
+               #ifdef USE_MPI
+                       }
+               #endif
+               
+               if (m->control_pressed) { delete distCalculator; remove(outputFile.c_str()); return 0; }
+               
                delete distCalculator;
                
                m->mothurOutEndLine();
@@ -282,6 +349,9 @@ int DistanceCommand::driver(int startLine, int endLine, string dFileName, float
                                outFile << name << '\t';        
                        }
                        for(int j=0;j<i;j++){
+                               
+                               if (m->control_pressed) { outFile.close(); return 0;  }
+                               
                                distCalculator->calcDist(alignDB.get(i), alignDB.get(j));
                                double dist = distCalculator->getDist();
                                
@@ -315,8 +385,73 @@ int DistanceCommand::driver(int startLine, int endLine, string dFileName, float
                exit(1);
        }
 }
+#ifdef USE_MPI
+/**************************************************************************************************/
+/////// need to fix to work with calcs and sequencedb
+int DistanceCommand::driverMPI(int startLine, int endLine, MPI_File& outMPI, float cutoff){
+       try {
+               MPI_Status status;
+               int startTime = time(NULL);
+               
+               string outputString = "";
+               
+               for(int i=startLine;i<endLine;i++){
+       
+                       for(int j=0;j<i;j++){
+                               
+                               if (m->control_pressed) {  return 0;  }
+                               
+                               distCalculator->calcDist(alignDB.get(i), alignDB.get(j));
+                               double dist = distCalculator->getDist();
+                               
+                               if(dist <= cutoff){
+                                       if (output == "column") { outputString += (alignDB.get(i).getName() + ' ' + alignDB.get(j).getName() + ' ' + toString(dist) + '\n'); }
+                               }
+                               
+                               if ((output == "square") || (output == "lt")){ //make a square column you can convert to square phylip
+                                       outputString += (alignDB.get(i).getName() + ' ' + alignDB.get(j).getName() + ' ' + toString(dist) + '\n');
+                                       outputString += (alignDB.get(j).getName() + ' ' + alignDB.get(i).getName() + ' ' + toString(dist) + '\n');
+                               }
+
+                       }
+                       
+                       if(i % 100 == 0){
+                               m->mothurOut(toString(i) + "\t" + toString(time(NULL) - startTime)); m->mothurOutEndLine();
+                       }
+                       
+                       if(i % 10 == 0){ //output to file 
+                               //send results to parent
+                               int length = outputString.length();
+                               char buf[length];
+                               strcpy(buf, outputString.c_str()); 
+                               
+                               MPI_File_write_shared(outMPI, buf, length, MPI_CHAR, &status);
+                               outputString = "";
+                       }
+                       
+               }
+               
+               m->mothurOut(toString(endLine-1) + "\t" + toString(time(NULL) - startTime)); m->mothurOutEndLine();
+               if(outputString != ""){ //output to file 
+                               //send results to parent
+                               int length = outputString.length();
+                               char buf[length];
+                               strcpy(buf, outputString.c_str()); 
+                               
+                               MPI_File_write_shared(outMPI, buf, length, MPI_CHAR, &status);
+                               outputString = "";
+               }
+               
+               return 1;
+       }
+       catch(exception& e) {
+               m->errorOut(e, "DistanceCommand", "driverMPI");
+               exit(1);
+       }
+}
+#endif
 /**************************************************************************************************/
-void DistanceCommand::convertMatrix(string outputFile) {
+int DistanceCommand::convertMatrix(string outputFile) {
        try{
 
                //sort file by first column so the distances for each row are together
@@ -358,6 +493,8 @@ void DistanceCommand::convertMatrix(string outputFile) {
                //openInputFile(outfile, in);
                
                while(!in.eof()) {
+                       if (m->control_pressed) { in.close(); remove(outfile.c_str()); out.close(); return 0; }
+                       
                        in >> first >> second >> dist; gobble(in);
                                
                        if (first != currentRow) {
@@ -393,12 +530,108 @@ void DistanceCommand::convertMatrix(string outputFile) {
                
                remove(outfile.c_str());
                
+               return 1;
+               
        }
        catch(exception& e) {
                m->errorOut(e, "DistanceCommand", "convertMatrix");
                exit(1);
        }
 }
+/**************************************************************************************************/
+int DistanceCommand::convertToLowerTriangle(string outputFile) {
+       try{
+
+               //sort file by first column so the distances for each row are together
+               string outfile = getRootName(outputFile) + "sorted.dist.temp";
+               
+               //use the unix sort 
+               #if defined (__APPLE__) || (__MACH__) || (linux) || (__linux)
+                       string command = "sort -n " + outputFile + " -o " + outfile;
+                       system(command.c_str());
+               #else //sort using windows sort
+                       string command = "sort " + outputFile + " /O " + outfile;
+                       system(command.c_str());
+               #endif
+               
+
+               //output to new file distance for each row and save positions in file where new row begins
+               ifstream in;
+               openInputFile(outfile, in);
+               
+               ofstream out;
+               openOutputFile(outputFile, out);
+               
+               out.setf(ios::fixed, ios::floatfield); out.setf(ios::showpoint);
+
+               out << alignDB.getNumSeqs() << endl;
+               
+               //get first currentRow
+               string first, currentRow, second;
+               float dist;
+               int i, j;
+               i = 0; j = 0;
+               map<string, float> rowDists; //take advantage of the fact that maps are already sorted by key 
+               map<string, float>::iterator it;
+               
+               in >> first;
+               currentRow = first;
+               
+               rowDists[first] = 0.00; //distance to yourself is 0.0
+               
+               in.seekg(0);
+               //openInputFile(outfile, in);
+               
+               while(!in.eof()) {
+                       if (m->control_pressed) { in.close(); remove(outfile.c_str()); out.close(); return 0; }
+                       
+                       in >> first >> second >> dist; gobble(in);
+                               
+                       if (first != currentRow) {
+                               //print out last row
+                               out << currentRow << '\t'; //print name
+
+                               //print dists
+                               for (it = rowDists.begin(); it != rowDists.end(); it++) {
+                                       if (j >= i) { break; }
+                                       out << it->second << '\t';
+                                       j++;
+                               }
+                               out << endl;
+                               
+                               //start new row
+                               currentRow = first;
+                               rowDists.clear();
+                               rowDists[first] = 0.00;
+                               rowDists[second] = dist;
+                               j = 0;
+                               i++;
+                       }else{
+                               rowDists[second] = dist;
+                       }
+               }
+               //print out last row
+               out << currentRow << '\t'; //print name
+                               
+               //print dists
+               for (it = rowDists.begin(); it != rowDists.end(); it++) {
+                       out << it->second << '\t';
+               }
+               out << endl;
+               
+               in.close();
+               out.close();
+               
+               remove(outfile.c_str());
+               
+               return 1;
+               
+       }
+       catch(exception& e) {
+               m->errorOut(e, "DistanceCommand", "convertToLowerTriangle");
+               exit(1);
+       }
+}
 /**************************************************************************************************
 void DistanceCommand::appendFiles(string temp, string filename) {
        try{