]> git.donarmstrong.com Git - mothur.git/blob - distancecommand.cpp
paralellized screen.seqs and added mpi code to it. fixed bug with all mpi commands...
[mothur.git] / distancecommand.cpp
1 /*
2  *  distancecommand.cpp
3  *  Mothur
4  *
5  *  Created by Sarah Westcott on 5/7/09.
6  *  Copyright 2009 Schloss Lab UMASS Amherst. All rights reserved.
7  *
8  */
9
10 #include "distancecommand.h"
11 #include "ignoregaps.h"
12 #include "eachgapdist.h"
13 #include "eachgapignore.h"
14 #include "onegapdist.h"
15 #include "onegapignore.h"
16
17 //**********************************************************************************************************************
18
19 DistanceCommand::DistanceCommand(string option) {
20         try {
21                 abort = false;
22                 Estimators.clear();
23                                 
24                 //allow user to run help
25                 if(option == "help") { help(); abort = true; }
26                 
27                 else {
28                         //valid paramters for this command
29                         string Array[] =  {"fasta", "output", "calc", "countends", "cutoff", "processors", "outputdir","inputdir"};
30                         vector<string> myArray (Array, Array+(sizeof(Array)/sizeof(string)));
31                         
32                         OptionParser parser(option);
33                         map<string, string> parameters = parser.getParameters();
34                         
35                         ValidParameters validParameter;
36                         map<string, string>::iterator it2;
37                 
38                         //check to make sure all parameters are valid for command
39                         for (it2 = parameters.begin(); it2 != parameters.end(); it2++) { 
40                                 if (validParameter.isValidParameter(it2->first, myArray, it2->second) != true) {  abort = true;  }
41                         }
42                         
43                         //if the user changes the input directory command factory will send this info to us in the output parameter 
44                         string inputDir = validParameter.validFile(parameters, "inputdir", false);              
45                         if (inputDir == "not found"){   inputDir = "";          }
46                         else {
47                                 string path;
48                                 it2 = parameters.find("fasta");
49                                 //user has given a template file
50                                 if(it2 != parameters.end()){ 
51                                         path = hasPath(it2->second);
52                                         //if the user has not given a path then, add inputdir. else leave path alone.
53                                         if (path == "") {       parameters["fasta"] = inputDir + it2->second;           }
54                                 }
55                         }
56
57                         //check for required parameters
58                         fastafile = validParameter.validFile(parameters, "fasta", true);
59                         if (fastafile == "not found") { m->mothurOut("fasta is a required parameter for the dist.seqs command."); m->mothurOutEndLine(); abort = true; }
60                         else if (fastafile == "not open") { abort = true; }     
61                         else{
62                                 ifstream inFASTA;
63                                 openInputFile(fastafile, inFASTA);
64                                 alignDB = SequenceDB(inFASTA); 
65                                 inFASTA.close();
66                         }
67                         
68                         //if the user changes the output directory command factory will send this info to us in the output parameter 
69                         outputDir = validParameter.validFile(parameters, "outputdir", false);           if (outputDir == "not found"){  
70                                 outputDir = ""; 
71                                 outputDir += hasPath(fastafile); //if user entered a file with a path then preserve it  
72                         }
73
74                         //check for optional parameter and set defaults
75                         // ...at some point should added some additional type checking...
76                         calc = validParameter.validFile(parameters, "calc", false);                     
77                         if (calc == "not found") { calc = "onegap";  }
78                         else { 
79                                  if (calc == "default")  {  calc = "onegap";  }
80                         }
81                         splitAtDash(calc, Estimators);
82
83                         string temp;
84                         temp = validParameter.validFile(parameters, "countends", false);        if(temp == "not found"){        temp = "T";     }
85                         convert(temp, countends); 
86                         
87                         temp = validParameter.validFile(parameters, "cutoff", false);           if(temp == "not found"){        temp = "1.0"; }
88                         convert(temp, cutoff); 
89                         
90                         temp = validParameter.validFile(parameters, "processors", false);       if(temp == "not found"){        temp = "1"; }
91                         convert(temp, processors); 
92                         
93                         output = validParameter.validFile(parameters, "output", false);         if(output == "not found"){      output = "column"; }
94                         
95                         if ((output != "column") && (output != "lt") && (output != "square")) { m->mothurOut(output + " is not a valid output form. Options are column, lt and square. I will use column."); m->mothurOutEndLine(); output = "column"; }
96                         
97                         ValidCalculators validCalculator;
98                         
99                         if (isTrue(countends) == true) {
100                                 for (int i=0; i<Estimators.size(); i++) {
101                                         if (validCalculator.isValidCalculator("distance", Estimators[i]) == true) { 
102                                                 if (Estimators[i] == "nogaps")                  {       distCalculator = new ignoreGaps();      }
103                                                 else if (Estimators[i] == "eachgap")    {       distCalculator = new eachGapDist();     }
104                                                 else if (Estimators[i] == "onegap")             {       distCalculator = new oneGapDist();      }
105                                         }
106                                 }
107                         }else {
108                                 for (int i=0; i<Estimators.size(); i++) {
109                                         if (validCalculator.isValidCalculator("distance", Estimators[i]) == true) { 
110                                                 if (Estimators[i] == "nogaps")          {       distCalculator = new ignoreGaps();                                      }
111                                                 else if (Estimators[i] == "eachgap"){   distCalculator = new eachGapIgnoreTermGapDist();        }
112                                                 else if (Estimators[i] == "onegap")     {       distCalculator = new oneGapIgnoreTermGapDist();         }
113                                         }
114                                 }
115                         }
116
117                 }
118                                 
119         }
120         catch(exception& e) {
121                 m->errorOut(e, "DistanceCommand", "DistanceCommand");
122                 exit(1);
123         }
124 }
125
126 //**********************************************************************************************************************
127
128 DistanceCommand::~DistanceCommand(){
129         
130         for(int i=0;i<lines.size();i++){
131                 delete lines[i];
132         }
133         
134 }
135         
136 //**********************************************************************************************************************
137
138 void DistanceCommand::help(){
139         try {
140                 m->mothurOut("The dist.seqs command reads a file containing sequences and creates a distance file.\n");
141                 m->mothurOut("The dist.seqs command parameters are fasta, calc, countends, output, cutoff and processors.  \n");
142                 m->mothurOut("The fasta parameter is required.\n");
143                 m->mothurOut("The calc parameter allows you to specify the method of calculating the distances.  Your options are: nogaps, onegap or eachgap. The default is onegap.\n");
144                 m->mothurOut("The countends parameter allows you to specify whether to include terminal gaps in distance.  Your options are: T or F. The default is T.\n");
145                 m->mothurOut("The cutoff parameter allows you to specify maximum distance to keep. The default is 1.0.\n");
146                 m->mothurOut("The output parameter allows you to specify format of your distance matrix. Options are column, lt, and square. The default is column.\n");
147                 m->mothurOut("The processors parameter allows you to specify number of processors to use.  The default is 1.\n");
148                 m->mothurOut("The dist.seqs command should be in the following format: \n");
149                 m->mothurOut("dist.seqs(fasta=yourFastaFile, calc=yourCalc, countends=yourEnds, cutoff= yourCutOff, processors=yourProcessors) \n");
150                 m->mothurOut("Example dist.seqs(fasta=amazon.fasta, calc=eachgap, countends=F, cutoff= 2.0, processors=3).\n");
151                 m->mothurOut("Note: No spaces between parameter labels (i.e. calc), '=' and parameters (i.e.yourCalc).\n\n");
152         }
153         catch(exception& e) {
154                 m->errorOut(e, "DistanceCommand", "help");
155                 exit(1);
156         }
157 }
158 //**********************************************************************************************************************
159
160 int DistanceCommand::execute(){
161         try {
162                 
163                 if (abort == true) { return 0; }
164                 
165                 int startTime = time(NULL);
166                 
167                 int numSeqs = alignDB.getNumSeqs();
168                 cutoff += 0.005;
169                 
170                 string outputFile;
171                                 
172                 if (output == "lt") { //does the user want lower triangle phylip formatted file 
173                         outputFile = outputDir + getRootName(getSimpleName(fastafile)) + "phylip.dist";
174                         remove(outputFile.c_str());
175                         
176                         //output numSeqs to phylip formatted dist file
177                 }else if (output == "column") { //user wants column format
178                         outputFile = outputDir + getRootName(getSimpleName(fastafile)) + "dist";
179                         remove(outputFile.c_str());
180                 }else { //assume square
181                         outputFile = outputDir + getRootName(getSimpleName(fastafile)) + "square.dist";
182                         remove(outputFile.c_str());
183                 }
184                 
185
186 #ifdef USE_MPI
187                 
188                 int pid, start, end; 
189                 int tag = 2001;
190                                 
191                 MPI_Status status; 
192                 MPI_Comm_size(MPI_COMM_WORLD, &processors); //set processors to the number of mpi processes running
193                 MPI_Comm_rank(MPI_COMM_WORLD, &pid); //find out who we are
194                 
195                 //each process gets where it should start and stop in the file
196                 start = int (sqrt(float(pid)/float(processors)) * numSeqs);
197                 end = int (sqrt(float(pid+1)/float(processors)) * numSeqs);
198                 
199                 if (output != "lt") {
200                         MPI_File outMPI;
201                         int amode=MPI_MODE_CREATE|MPI_MODE_WRONLY; 
202
203                         //char* filename = new char[outputFile.length()];
204                         //memcpy(filename, outputFile.c_str(), outputFile.length());
205                         
206                         char filename[1024];
207                         strcpy(filename, outputFile.c_str());
208                         
209                         MPI_File_open(MPI_COMM_WORLD, filename, amode, MPI_INFO_NULL, &outMPI);
210                         //delete filename;
211
212                         if (pid == 0) { //you are the root process 
213                         
214                                 //do your part
215                                 string outputMyPart;
216                                 driverMPI(start, end, outMPI, cutoff);
217                                 
218                                 if (m->control_pressed) { MPI_File_close(&outMPI);  delete distCalculator;  return 0; }
219                         
220                                 //wait on chidren
221                                 for(int i = 1; i < processors; i++) { 
222                                         if (m->control_pressed) { MPI_File_close(&outMPI);  delete distCalculator;  return 0; }
223                                         
224                                         char buf[4];
225                                         MPI_Recv(buf, 4, MPI_CHAR, i, tag, MPI_COMM_WORLD, &status); 
226                                 }
227                         }else { //you are a child process
228                                 //do your part
229                                 driverMPI(start, end, outMPI, cutoff);
230                                 
231                                 if (m->control_pressed) { MPI_File_close(&outMPI);  delete distCalculator;  return 0; }
232                         
233                                 char buf[4];
234                                 strcpy(buf, "done"); 
235                                 //tell parent you are done.
236                                 MPI_Send(buf, 4, MPI_CHAR, 0, tag, MPI_COMM_WORLD);
237                         }
238                         
239                         MPI_File_close(&outMPI);
240                         
241                 }else { //lower triangle format
242                         if (pid == 0) { //you are the root process 
243                         
244                                 //do your part
245                                 string outputMyPart;
246                                 long mySize;
247                                 driverMPI(start, end, outputFile, mySize);
248         
249                                 if (m->control_pressed) {  delete distCalculator;  return 0; }
250                                 
251                                 int amode=MPI_MODE_APPEND|MPI_MODE_WRONLY|MPI_MODE_CREATE; //
252                                 MPI_File outMPI;
253                                 MPI_File inMPI;
254
255                                 //char* filename = new char[outputFile.length()];
256                                 //memcpy(filename, outputFile.c_str(), outputFile.length());
257                                 
258                                 char filename[1024];
259                                 strcpy(filename, outputFile.c_str());
260
261                                 MPI_File_open(MPI_COMM_SELF, filename, amode, MPI_INFO_NULL, &outMPI);
262                                 //delete filename;
263
264                                 //wait on chidren
265                                 for(int b = 1; b < processors; b++) { 
266                                         long fileSize;
267                                         
268                                         if (m->control_pressed) { MPI_File_close(&outMPI);  delete distCalculator;  return 0; }
269                                         
270                                         MPI_Recv(&fileSize, 1, MPI_LONG, b, tag, MPI_COMM_WORLD, &status); 
271                                         
272                                         string outTemp = outputFile + toString(b) + ".temp";
273
274                                         char* buf = new char[outTemp.length()];
275                                         memcpy(buf, outTemp.c_str(), outTemp.length());
276                                         
277                                         MPI_File_open(MPI_COMM_SELF, buf, MPI_MODE_DELETE_ON_CLOSE|MPI_MODE_RDONLY, MPI_INFO_NULL, &inMPI);
278                                         delete buf;
279
280                                         int count = 0;
281                                         while (count < fileSize) { 
282                                                 char buf2[1];
283                                                 MPI_File_read(inMPI, buf2, 1, MPI_CHAR, &status);
284                                                 MPI_File_write(outMPI, buf2, 1, MPI_CHAR, &status);
285                                                 count += 1;
286                                         }
287                                         
288                                         MPI_File_close(&inMPI); //deleted on close
289                                 }
290                                 
291                                 MPI_File_close(&outMPI);
292                         }else { //you are a child process
293                                 //do your part
294                                 long size;
295                                 driverMPI(start, end, (outputFile + toString(pid) + ".temp"), size);
296                                 
297                                 if (m->control_pressed) { delete distCalculator;  return 0; }
298                         
299                                 //tell parent you are done.
300                                 MPI_Send(&size, 1, MPI_LONG, 0, tag, MPI_COMM_WORLD);
301                         }
302                 }
303 #else           
304                                 
305         #if defined (__APPLE__) || (__MACH__) || (linux) || (__linux)
306                 //if you don't need to fork anything
307                 if(processors == 1){
308                         driver(0, numSeqs, outputFile, cutoff);
309                 }else{ //you have multiple processors
310                         
311                         for (int i = 0; i < processors; i++) {
312                                 lines.push_back(new linePair());
313                                 lines[i]->start = int (sqrt(float(i)/float(processors)) * numSeqs);
314                                 lines[i]->end = int (sqrt(float(i+1)/float(processors)) * numSeqs);
315                         }
316
317                         createProcesses(outputFile); 
318                 
319                         map<int, int>::iterator it = processIDS.begin();
320                         rename((outputFile + toString(it->second) + ".temp").c_str(), outputFile.c_str());
321                         it++;
322                         
323                         //append and remove temp files
324                         for (; it != processIDS.end(); it++) {
325                                 appendFiles((outputFile + toString(it->second) + ".temp"), outputFile);
326                                 remove((outputFile + toString(it->second) + ".temp").c_str());
327                         }
328                 }
329         #else
330                 ifstream inFASTA;
331                 driver(0, numSeqs, outputFile, cutoff);
332         #endif
333         
334 #endif
335                 if (m->control_pressed) { delete distCalculator; remove(outputFile.c_str()); return 0; }
336                 
337                 #ifdef USE_MPI
338                         MPI_Comm_rank(MPI_COMM_WORLD, &pid); 
339                                         
340                         if (pid == 0) { //only one process should output to screen
341                 #endif
342                 
343                 if (output == "square") {  convertMatrix(outputFile); }
344                 
345                 ifstream fileHandle;
346                 fileHandle.open(outputFile.c_str());
347                 if(fileHandle) {
348                         gobble(fileHandle);
349                         if (fileHandle.eof()) { m->mothurOut(outputFile + " is blank. This can result if there are no distances below your cutoff.");  m->mothurOutEndLine(); }
350                 }
351
352                 
353                 #ifdef USE_MPI
354                         }
355                 #endif
356                 
357                 if (m->control_pressed) { delete distCalculator; remove(outputFile.c_str()); return 0; }
358                 
359                 delete distCalculator;
360                 
361                 m->mothurOutEndLine();
362                 m->mothurOut("Output File Name: "); m->mothurOutEndLine();
363                 m->mothurOut(outputFile); m->mothurOutEndLine();
364                 m->mothurOutEndLine();
365                 m->mothurOut("It took " + toString(time(NULL) - startTime) + " to calculate the distances for " + toString(numSeqs) + " sequences."); m->mothurOutEndLine();
366                 return 0;
367                 
368         }
369         catch(exception& e) {
370                 m->errorOut(e, "DistanceCommand", "execute");
371                 exit(1);
372         }
373 }
374 /**************************************************************************************************/
375 void DistanceCommand::createProcesses(string filename) {
376         try {
377 #if defined (__APPLE__) || (__MACH__) || (linux) || (__linux)
378                 int process = 0;
379                 processIDS.clear();
380                 
381                 //loop through and create all the processes you want
382                 while (process != processors) {
383                         int pid = fork();
384                         
385                         if (pid > 0) {
386                                 processIDS[lines[process]->end] = pid;  //create map from line number to pid so you can append files in correct order later
387                                 process++;
388                         }else if (pid == 0){
389                                 driver(lines[process]->start, lines[process]->end, filename + toString(getpid()) + ".temp", cutoff);
390                                 exit(0);
391                         }else { m->mothurOut("unable to spawn the necessary processes."); m->mothurOutEndLine(); exit(0); }
392                 }
393         
394                 //force parent to wait until all the processes are done
395                 for (map<int, int>::iterator it = processIDS.begin(); it != processIDS.end(); it++) { 
396                         int temp = it->second;
397                         wait(&temp);
398                 }
399 #endif
400         }
401         catch(exception& e) {
402                 m->errorOut(e, "DistanceCommand", "createProcesses");
403                 exit(1);
404         }
405 }
406
407 /**************************************************************************************************/
408 /////// need to fix to work with calcs and sequencedb
409 int DistanceCommand::driver(int startLine, int endLine, string dFileName, float cutoff){
410         try {
411
412                 int startTime = time(NULL);
413                 
414                 //column file
415                 ofstream outFile(dFileName.c_str(), ios::trunc);
416                 outFile.setf(ios::fixed, ios::showpoint);
417                 outFile << setprecision(4);
418                 
419                 if((output == "lt") && startLine == 0){ outFile << alignDB.getNumSeqs() << endl;        }
420                 
421                 for(int i=startLine;i<endLine;i++){
422                         if(output == "lt")      {       
423                                 string name = alignDB.get(i).getName();
424                                 if (name.length() < 10) { //pad with spaces to make compatible
425                                         while (name.length() < 10) {  name += " ";  }
426                                 }
427                                 outFile << name << '\t';        
428                         }
429                         for(int j=0;j<i;j++){
430                                 
431                                 if (m->control_pressed) { outFile.close(); return 0;  }
432                                 
433                                 distCalculator->calcDist(alignDB.get(i), alignDB.get(j));
434                                 double dist = distCalculator->getDist();
435                                 
436                                 if(dist <= cutoff){
437                                         if (output == "column") { outFile << alignDB.get(i).getName() << ' ' << alignDB.get(j).getName() << ' ' << dist << endl; }
438                                 }
439                                 if (output == "lt") {  outFile << dist << '\t'; }
440                                 
441                                 if (output == "square") { //make a square column you can convert to square phylip
442                                         outFile << alignDB.get(i).getName() << '\t' << alignDB.get(j).getName() << '\t' << dist << endl;
443                                         outFile << alignDB.get(j).getName() << '\t' << alignDB.get(i).getName() << '\t' << dist << endl;
444                                 }
445
446                         }
447                         
448                         if (output == "lt") { outFile << endl; }
449                         
450                         if(i % 100 == 0){
451                                 m->mothurOut(toString(i) + "\t" + toString(time(NULL) - startTime)); m->mothurOutEndLine();
452                         }
453                         
454                 }
455                 m->mothurOut(toString(endLine-1) + "\t" + toString(time(NULL) - startTime)); m->mothurOutEndLine();
456                 
457                 outFile.close();
458                 
459                 return 1;
460         }
461         catch(exception& e) {
462                 m->errorOut(e, "DistanceCommand", "driver");
463                 exit(1);
464         }
465 }
466 #ifdef USE_MPI
467 /**************************************************************************************************/
468 /////// need to fix to work with calcs and sequencedb
469 int DistanceCommand::driverMPI(int startLine, int endLine, MPI_File& outMPI, float cutoff){
470         try {
471                 MPI_Status status;
472                 int startTime = time(NULL);
473                 
474                 string outputString = "";
475                 
476                 for(int i=startLine;i<endLine;i++){
477         
478                         for(int j=0;j<i;j++){
479                                 
480                                 if (m->control_pressed) {  return 0;  }
481                                 
482                                 distCalculator->calcDist(alignDB.get(i), alignDB.get(j));
483                                 double dist = distCalculator->getDist();
484                                 
485                                 if(dist <= cutoff){
486                                         if (output == "column") { outputString += (alignDB.get(i).getName() + ' ' + alignDB.get(j).getName() + ' ' + toString(dist) + '\n'); }
487                                 }
488                                 
489                                 if (output == "square") { //make a square column you can convert to square phylip
490                                         outputString += (alignDB.get(i).getName() + ' ' + alignDB.get(j).getName() + ' ' + toString(dist) + '\n');
491                                         outputString += (alignDB.get(j).getName() + ' ' + alignDB.get(i).getName() + ' ' + toString(dist) + '\n');
492                                 }
493                         }
494                         
495                         if(i % 100 == 0){
496                                 //m->mothurOut(toString(i) + "\t" + toString(time(NULL) - startTime)); m->mothurOutEndLine();
497                                 cout << i << '\t' << (time(NULL) - startTime) << endl;
498                         }
499                         
500                          
501                         //send results to parent
502                         int length = outputString.length();
503
504                         char* buf = new char[length];
505                         memcpy(buf, outputString.c_str(), length);
506                         
507                         MPI_File_write_shared(outMPI, buf, length, MPI_CHAR, &status);
508                         outputString = "";
509                         delete buf;
510                         
511                 }
512                 
513                 //m->mothurOut(toString(endLine-1) + "\t" + toString(time(NULL) - startTime)); m->mothurOutEndLine();
514                 cout << (endLine-1) << '\t' << (time(NULL) - startTime) << endl;                
515                 return 1;
516         }
517         catch(exception& e) {
518                 m->errorOut(e, "DistanceCommand", "driverMPI");
519                 exit(1);
520         }
521 }
522 /**************************************************************************************************/
523 /////// need to fix to work with calcs and sequencedb
524 int DistanceCommand::driverMPI(int startLine, int endLine, string file, long& size){
525         try {
526                 MPI_Status status;
527                 
528                 MPI_File outMPI;
529                 int amode=MPI_MODE_CREATE|MPI_MODE_WRONLY; 
530
531                 //char* filename = new char[file.length()];
532                 //memcpy(filename, file.c_str(), file.length());
533                 
534                 char filename[1024];
535                 strcpy(filename, file.c_str());
536
537                 MPI_File_open(MPI_COMM_SELF, filename, amode, MPI_INFO_NULL, &outMPI);
538                 //delete filename;
539
540                 int startTime = time(NULL);
541                 
542                 string outputString = "";
543                 size = 0;
544                 
545                 if((output == "lt") && startLine == 0){ outputString += toString(alignDB.getNumSeqs()) + "\n";  }
546                 
547                 for(int i=startLine;i<endLine;i++){
548                         if(output == "lt")      {       
549                                 string name = alignDB.get(i).getName();
550                                 if (name.length() < 10) { //pad with spaces to make compatible
551                                         while (name.length() < 10) {  name += " ";  }
552                                 }
553                                 outputString += name + "\t";    
554                         }
555                         for(int j=0;j<i;j++){
556                                 
557                                 if (m->control_pressed) {  return 0;  }
558                                 
559                                 distCalculator->calcDist(alignDB.get(i), alignDB.get(j));
560                                 double dist = distCalculator->getDist();
561                                 
562                                 if (output == "lt") {  outputString += toString(dist) + "\t"; }
563                         }
564                         
565                         if (output == "lt") { outputString += "\n"; }
566
567                 
568                         if(i % 100 == 0){
569                                 //m->mothurOut(toString(i) + "\t" + toString(time(NULL) - startTime)); m->mothurOutEndLine();
570                                 cout << i << '\t' << (time(NULL) - startTime) << endl;
571                         }
572                         
573                         
574                         //send results to parent
575                         int length = outputString.length();
576                         char* buf = new char[length];
577                         memcpy(buf, outputString.c_str(), length);
578                         
579                         MPI_File_write(outMPI, buf, length, MPI_CHAR, &status);
580                         size += outputString.length();
581                         outputString = "";
582                         delete buf;
583                 }
584                 
585                 //m->mothurOut(toString(endLine-1) + "\t" + toString(time(NULL) - startTime)); m->mothurOutEndLine();
586                 cout << (endLine-1) << '\t' << (time(NULL) - startTime) << endl;
587                 MPI_File_close(&outMPI);
588                 
589                 return 1;
590         }
591         catch(exception& e) {
592                 m->errorOut(e, "DistanceCommand", "driverMPI");
593                 exit(1);
594         }
595 }
596 #endif
597 /**************************************************************************************************/
598 int DistanceCommand::convertMatrix(string outputFile) {
599         try{
600
601                 //sort file by first column so the distances for each row are together
602                 string outfile = getRootName(outputFile) + "sorted.dist.temp";
603                 
604                 //use the unix sort 
605                 #if defined (__APPLE__) || (__MACH__) || (linux) || (__linux)
606                         string command = "sort -n " + outputFile + " -o " + outfile;
607                         system(command.c_str());
608                 #else //sort using windows sort
609                         string command = "sort " + outputFile + " /O " + outfile;
610                         system(command.c_str());
611                 #endif
612                 
613
614                 //output to new file distance for each row and save positions in file where new row begins
615                 ifstream in;
616                 openInputFile(outfile, in);
617                 
618                 ofstream out;
619                 openOutputFile(outputFile, out);
620                 
621                 out.setf(ios::fixed, ios::floatfield); out.setf(ios::showpoint);
622
623                 out << alignDB.getNumSeqs() << endl;
624                 
625                 //get first currentRow
626                 string first, currentRow, second;
627                 float dist;
628                 map<string, float> rowDists; //take advantage of the fact that maps are already sorted by key 
629                 map<string, float>::iterator it;
630                 
631                 in >> first;
632                 currentRow = first;
633                 
634                 rowDists[first] = 0.00; //distance to yourself is 0.0
635                 
636                 in.seekg(0);
637                 //openInputFile(outfile, in);
638                 
639                 while(!in.eof()) {
640                         if (m->control_pressed) { in.close(); remove(outfile.c_str()); out.close(); return 0; }
641                         
642                         in >> first >> second >> dist; gobble(in);
643                                 
644                         if (first != currentRow) {
645                                 //print out last row
646                                 out << currentRow << '\t'; //print name
647
648                                 //print dists
649                                 for (it = rowDists.begin(); it != rowDists.end(); it++) {
650                                         out << it->second << '\t';
651                                 }
652                                 out << endl;
653                                 
654                                 //start new row
655                                 currentRow = first;
656                                 rowDists.clear();
657                                 rowDists[first] = 0.00;
658                                 rowDists[second] = dist;
659                         }else{
660                                 rowDists[second] = dist;
661                         }
662                 }
663                 //print out last row
664                 out << currentRow << '\t'; //print name
665                                 
666                 //print dists
667                 for (it = rowDists.begin(); it != rowDists.end(); it++) {
668                         out << it->second << '\t';
669                 }
670                 out << endl;
671                 
672                 in.close();
673                 out.close();
674                 
675                 remove(outfile.c_str());
676                 
677                 return 1;
678                 
679         }
680         catch(exception& e) {
681                 m->errorOut(e, "DistanceCommand", "convertMatrix");
682                 exit(1);
683         }
684 }
685 /**************************************************************************************************/
686 int DistanceCommand::convertToLowerTriangle(string outputFile) {
687         try{
688
689                 //sort file by first column so the distances for each row are together
690                 string outfile = getRootName(outputFile) + "sorted.dist.temp";
691                 
692                 //use the unix sort 
693                 #if defined (__APPLE__) || (__MACH__) || (linux) || (__linux)
694                         string command = "sort -n " + outputFile + " -o " + outfile;
695                         system(command.c_str());
696                 #else //sort using windows sort
697                         string command = "sort " + outputFile + " /O " + outfile;
698                         system(command.c_str());
699                 #endif
700                 
701
702                 //output to new file distance for each row and save positions in file where new row begins
703                 ifstream in;
704                 openInputFile(outfile, in);
705                 
706                 ofstream out;
707                 openOutputFile(outputFile, out);
708                 
709                 out.setf(ios::fixed, ios::floatfield); out.setf(ios::showpoint);
710
711                 out << alignDB.getNumSeqs() << endl;
712                 
713                 //get first currentRow
714                 string first, currentRow, second;
715                 float dist;
716                 int i, j;
717                 i = 0; j = 0;
718                 map<string, float> rowDists; //take advantage of the fact that maps are already sorted by key 
719                 map<string, float>::iterator it;
720                 
721                 in >> first;
722                 currentRow = first;
723                 
724                 rowDists[first] = 0.00; //distance to yourself is 0.0
725                 
726                 in.seekg(0);
727                 //openInputFile(outfile, in);
728                 
729                 while(!in.eof()) {
730                         if (m->control_pressed) { in.close(); remove(outfile.c_str()); out.close(); return 0; }
731                         
732                         in >> first >> second >> dist; gobble(in);
733                                 
734                         if (first != currentRow) {
735                                 //print out last row
736                                 out << currentRow << '\t'; //print name
737
738                                 //print dists
739                                 for (it = rowDists.begin(); it != rowDists.end(); it++) {
740                                         if (j >= i) { break; }
741                                         out << it->second << '\t';
742                                         j++;
743                                 }
744                                 out << endl;
745                                 
746                                 //start new row
747                                 currentRow = first;
748                                 rowDists.clear();
749                                 rowDists[first] = 0.00;
750                                 rowDists[second] = dist;
751                                 j = 0;
752                                 i++;
753                         }else{
754                                 rowDists[second] = dist;
755                         }
756                 }
757                 //print out last row
758                 out << currentRow << '\t'; //print name
759                                 
760                 //print dists
761                 for (it = rowDists.begin(); it != rowDists.end(); it++) {
762                         out << it->second << '\t';
763                 }
764                 out << endl;
765                 
766                 in.close();
767                 out.close();
768                 
769                 remove(outfile.c_str());
770                 
771                 return 1;
772                 
773         }
774         catch(exception& e) {
775                 m->errorOut(e, "DistanceCommand", "convertToLowerTriangle");
776                 exit(1);
777         }
778 }
779 /**************************************************************************************************
780 void DistanceCommand::appendFiles(string temp, string filename) {
781         try{
782                 ofstream output;
783                 ifstream input;
784         
785                 //open output file in append mode
786                 openOutputFileAppend(filename, output);
787                 openInputFile(temp, input);
788                 
789                 while(char c = input.get()){
790                         if(input.eof())         {       break;                  }
791                         else                            {       output << c;    }
792                 }
793                 
794                 input.close();
795                 output.close();
796         }
797         catch(exception& e) {
798                 m->errorOut(e, "DistanceCommand", "appendFiles");
799                 exit(1);
800         }
801 }
802 /**************************************************************************************************/
803
804