modified bin.seqs and get.oturep commands to include use of a groupfile if provided...

[mothur.git] / rarefact.cpp

diff --git a/rarefact.cpp b/rarefact.cpp
index c064125c37e846c538075f04d652eca6ec336708..67662496618d427c12363ea5ec120f21578e3f12 100644 (file)
--- a/rarefact.cpp
+++ b/rarefact.cpp
@@ -8,7 +8,7 @@
  */
 
 #include "rarefact.h"
-#include "ordervector.hpp"
+//#include "ordervector.hpp"
 
 /***********************************************************************/
 
@@ -52,11 +52,15 @@ void Rarefact::getCurve(int increment = 1, int nIters = 1000){
                        for(int i=0;i<displays.size();i++){
                                displays[i]->reset();
                        }
+                       
+                       delete lookup;
+                       delete rank;
                }
 
                for(int i=0;i<displays.size();i++){
                        displays[i]->close();
                }
+               delete rcd;
        }
        catch(exception& e) {
                cout << "Standard Error: " << e.what() << " has occurred in the Rarefact class Function getCurve. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";
@@ -72,65 +76,56 @@ void Rarefact::getCurve(int increment = 1, int nIters = 1000){
 
 void Rarefact::getSharedCurve(int increment = 1, int nIters = 1000){
 try {
-               globaldata = GlobalData::getInstance();
                SharedRarefactionCurveData* rcd = new SharedRarefactionCurveData();
                
+               label = lookup[0]->getLabel();
+               
                //register the displays
                for(int i=0;i<displays.size();i++){
                        rcd->registerDisplay(displays[i]);
                }
-
+               
                for(int iter=0;iter<nIters;iter++){
-                       //clear out old values for new iter
-                       lookup.clear();
                
-                       //create and initialize vector of sharedvectors, one for each group
-                       for (int i = 0; i < globaldata->Groups.size(); i++) { 
-                               SharedRAbundVector* temp = new SharedRAbundVector(sharedorder->getNumBins());
-                               temp->setLabel(sharedorder->getLabel());
-                               temp->setGroup(globaldata->Groups[i]);
-                               lookup.push_back(temp);
-                       }
-                       
                        for(int i=0;i<displays.size();i++){
                                displays[i]->init(label);                 
                        }
-               
-                       //sample all the members
-                       for(int i=0;i<numSeqs;i++){
-                               //get first sample
-                               individual chosen = sharedorder->get(i);
-                               int abundance; 
-                                       
-                               //set info for sharedvector in chosens group
-                               for (int j = 0; j < lookup.size(); j++) { 
-                                       if (chosen.group == lookup[j]->getGroup()) {
-                                               abundance = lookup[j]->getAbundance(chosen.bin);
-                                               lookup[j]->set(chosen.bin, (abundance + 1), chosen.group);
-                                               break;
-                                       }
-                               }
-                       }
-                       
+
                        //randomize the groups
                        random_shuffle(lookup.begin(), lookup.end());
-               
+                       
+                       //make merge the size of lookup[0]
+                       SharedRAbundVector* merge = new SharedRAbundVector(lookup[0]->size());
+                       
+                       //make copy of lookup zero
+                       for(int i = 0; i<lookup[0]->size(); i++) {
+                               merge->set(i, lookup[0]->getAbundance(i), "merge");
+                       }
+                       
+                       vector<SharedRAbundVector*> subset;
                        //send each group one at a time
                        for (int k = 0; k < lookup.size(); k++) { 
-                               rcd->updateSharedData(lookup[0], lookup[k], k+1, numGroupComb);
-                               mergeVectors(lookup[0], lookup[k]);
+                               subset.clear(); //clears out old pair of sharedrabunds
+                               //add in new pair of sharedrabunds
+                               subset.push_back(merge); subset.push_back(lookup[k]);
+                               
+                               rcd->updateSharedData(subset, k+1, numGroupComb);
+                               mergeVectors(merge, lookup[k]);
                        }
 
                        //resets output files
                        for(int i=0;i<displays.size();i++){
                                displays[i]->reset();
                        }
+                       
+                       delete merge;
                }
                
                for(int i=0;i<displays.size();i++){
                        displays[i]->close();
                }
-
+               
+               delete rcd;
        }
        catch(exception& e) {
                cout << "Standard Error: " << e.what() << " has occurred in the Rarefact class Function getSharedCurve. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";
@@ -148,9 +143,7 @@ void Rarefact::mergeVectors(SharedRAbundVector* shared1, SharedRAbundVector* sha
        try{
                for (int k = 0; k < shared1->size(); k++) {
                        //merge new species into shared1
-                       if ((shared1->getAbundance(k) == 0) && (shared2->getAbundance(k) != 0)) {
-                               shared1->set(k, shared2->getAbundance(k), "combo");  //set to 'combo' since this vector now contains multiple groups
-                       }
+                       shared1->set(k, (shared1->getAbundance(k) + shared2->getAbundance(k)), "combo");  //set to 'combo' since this vector now contains multiple groups
                }
        }
        catch(exception& e) {