update .gitignore

[mothur.git] / treegroupscommand.cpp

diff --git a/treegroupscommand.cpp b/treegroupscommand.cpp
index 0df53d5649881c924be60140338e5fff7c2bc840..4b4f21c86366000b8efb195bf94232059a78d03c 100644 (file)
--- a/treegroupscommand.cpp
+++ b/treegroupscommand.cpp
@@ -25,7 +25,7 @@ vector<string> TreeGroupCommand::setParameters(){
                CommandParameter pprecision("precision", "Number", "", "100", "", "", "","",false,false); parameters.push_back(pprecision);             
                CommandParameter plabel("label", "String", "", "", "", "", "","",false,false); parameters.push_back(plabel);
                CommandParameter pgroups("groups", "String", "", "", "", "", "","",false,false); parameters.push_back(pgroups);
-               CommandParameter pcalc("calc", "Multiple", "sharedsobs-sharedchao-sharedace-jabund-sorabund-jclass-sorclass-jest-sorest-thetayc-thetan-kstest-sharednseqs-ochiai-anderberg-kulczynski-kulczynskicody-lennon-morisitahorn-braycurtis-whittaker-odum-canberra-structeuclidean-structchord-hellinger-manhattan-structpearson-soergel-spearman-structkulczynski-speciesprofile-hamming-structchi2-gower-memchi2-memchord-memeuclidean-mempearson", "jclass-thetayc", "", "", "","",true,false,true); parameters.push_back(pcalc);
+               CommandParameter pcalc("calc", "Multiple", "sharedsobs-sharedchao-sharedace-jabund-sorabund-jclass-sorclass-jest-sorest-thetayc-thetan-kstest-sharednseqs-ochiai-anderberg-kulczynski-kulczynskicody-lennon-morisitahorn-braycurtis-whittaker-odum-canberra-structeuclidean-structchord-hellinger-manhattan-structpearson-soergel-spearman-structkulczynski-speciesprofile-hamming-structchi2-gower-memchi2-memchord-memeuclidean-mempearson-jsd-rjsd", "jclass-thetayc", "", "", "","",true,false,true); parameters.push_back(pcalc);
                
         CommandParameter pprocessors("processors", "Number", "", "1", "", "", "","",false,false,true); parameters.push_back(pprocessors);
 //CommandParameter poutput("output", "Multiple", "lt-square", "lt", "", "", "",false,false); parameters.push_back(poutput);
@@ -397,7 +397,12 @@ int TreeGroupCommand::execute(){
                                                treeCalculators.push_back(new MemEuclidean());
                                        }else if (Estimators[i] == "mempearson") { 
                                                treeCalculators.push_back(new MemPearson());
-                                       }
+                    }else if (Estimators[i] == "jsd") {
+                        treeCalculators.push_back(new JSD());
+                    }else if (Estimators[i] == "rjsd") {
+                        treeCalculators.push_back(new RJSD());
+                    }
+
                                }
                        }
                        
@@ -458,7 +463,7 @@ int TreeGroupCommand::execute(){
                 readMatrix->read(nameMap);
             }else if (countfile != "") {
                 ct = new CountTable();
-                ct->readTable(countfile);
+                ct->readTable(countfile, true, false);
                 readMatrix->read(ct);
             }else {
                 readMatrix->read(nameMap);
@@ -644,6 +649,7 @@ int TreeGroupCommand::makeSimsShared() {
             }else {
                 m->clearGroups();
                 Groups.clear();
+                m->Treenames.clear();
                 vector<SharedRAbundVector*> temp;
                 for (int i = 0; i < lookup.size(); i++) {
                     if (lookup[i]->getNumSeqs() < subsampleSize) { 
@@ -652,6 +658,7 @@ int TreeGroupCommand::makeSimsShared() {
                     }else { 
                         Groups.push_back(lookup[i]->getGroup()); 
                         temp.push_back(lookup[i]);
+                        m->Treenames.push_back(lookup[i]->getGroup());
                     }
                 } 
                 lookup = temp;
@@ -660,8 +667,18 @@ int TreeGroupCommand::makeSimsShared() {
             
             if (lookup.size() < 2) { m->mothurOut("You have not provided enough valid groups.  I cannot run the command."); m->mothurOutEndLine(); m->control_pressed = true; return 0; }
         }
-        
         numGroups = lookup.size();
+        
+        //sanity check to make sure processors < numComparisions
+        int numDists = 0;
+        for(int i=0;i<numGroups;i++){
+            for(int j=0;j<i;j++){
+                numDists++;
+                if (numDists > processors) { break; }
+            }
+        }
+        if (numDists < processors) { processors = numDists; }
+        
                lines.resize(processors);
                for (int i = 0; i < processors; i++) {
                        lines[i].start = int (sqrt(float(i)/float(processors)) * numGroups);
@@ -783,7 +800,7 @@ int TreeGroupCommand::process(vector<SharedRAbundVector*> thisLookup) {
 #if defined (__APPLE__) || (__MACH__) || (linux) || (__linux) || (__linux__) || (__unix__) || (__unix)
                 //loop through and create all the processes you want
                 while (process != processors) {
-                    int pid = fork();
+                    pid_t pid = fork();
                     
                     if (pid > 0) {
                         processIDS.push_back(pid); 
@@ -792,7 +809,7 @@ int TreeGroupCommand::process(vector<SharedRAbundVector*> thisLookup) {
                         
                         driver(thisItersLookup, lines[process].start, lines[process].end, calcDists);   
                         
-                        string tempdistFileName = m->getRootName(m->getSimpleName(sharedfile)) + toString(getpid()) + ".dist";
+                        string tempdistFileName = m->getRootName(m->getSimpleName(sharedfile)) + m->mothurGetpid(process) + ".dist";
                         ofstream outtemp;
                         m->openOutputFile(tempdistFileName, outtemp);
                         
@@ -916,12 +933,18 @@ int TreeGroupCommand::process(vector<SharedRAbundVector*> thisLookup) {
                 thisItersLookup.clear();
                 for (int i = 0; i < calcDists.size(); i++) {  calcDists[i].clear(); }
             }
+            
+            if (m->debug) {  m->mothurOut("[DEBUG]: iter = " + toString(thisIter) + ".\n"); }
                }
-               
+        
+               if (m->debug) {  m->mothurOut("[DEBUG]: done with iters.\n"); }
+            
         if (iters != 1) {
             //we need to find the average distance and standard deviation for each groups distance
             vector< vector<seqDist>  > calcAverages = m->getAverages(calcDistsTotals);  
             
+            if (m->debug) {  m->mothurOut("[DEBUG]: found averages.\n"); }
+            
             //create average tree for each calc
             for (int i = 0; i < calcDists.size(); i++) {
                 vector< vector<double> > matrix; //square matrix to represent the distance
@@ -951,6 +974,8 @@ int TreeGroupCommand::process(vector<SharedRAbundVector*> thisLookup) {
                 if (newTree != NULL) { writeTree(outputFile, newTree); }                
             }
             
+            if (m->debug) {  m->mothurOut("[DEBUG]: done averages trees.\n"); }
+            
             //create all trees for each calc and find their consensus tree
             for (int i = 0; i < calcDists.size(); i++) {
                 if (m->control_pressed) { break; }
@@ -982,7 +1007,7 @@ int TreeGroupCommand::process(vector<SharedRAbundVector*> thisLookup) {
                         int row = calcDistsTotals[myIter][i][j].seq1;
                         int column = calcDistsTotals[myIter][i][j].seq2;
                         double dist = calcDistsTotals[myIter][i][j].dist;
-                        
+                       
                         matrix[row][column] = dist;
                         matrix[column][row] = dist;
                     }
@@ -997,11 +1022,15 @@ int TreeGroupCommand::process(vector<SharedRAbundVector*> thisLookup) {
                 outAll.close();
                 if (m->control_pressed) { for (int k = 0; k < trees.size(); k++) { delete trees[k]; } }
                 
+                if (m->debug) {  m->mothurOut("[DEBUG]: done all trees.\n"); }
+                
                 Consensus consensus;
                 //clear old tree names if any
                 m->Treenames.clear(); m->Treenames = m->getGroups(); //may have changed if subsample eliminated groups
                 Tree* conTree = consensus.getTree(trees);
                 
+                if (m->debug) {  m->mothurOut("[DEBUG]: done cons tree.\n"); }
+                
                 //create a new filename
                 variables["[tag]"] = "cons";
                 string conFile = getOutputFileName("tree",variables);