X-Git-Url: https://git.donarmstrong.com/?p=mothur.git;a=blobdiff_plain;f=rarefact.cpp;h=0454d9e253473b368af34f35021dabbc4a50ae1b;hp=2f5f5067a388f1c9c6026a32cff81f003fec1106;hb=cf9987b67aa49777a4c91c2d21f96e58bf17aa82;hpb=c5c7502f435e1413c19e373dab1dfebcaa67588d

diff --git a/rarefact.cpp b/rarefact.cpp
index 2f5f506..0454d9e 100644
--- a/rarefact.cpp
+++ b/rarefact.cpp
@@ -12,13 +12,57 @@
 
 /***********************************************************************/
 
-void Rarefact::getCurve(int increment = 1, int nIters = 1000){
+int Rarefact::getCurve(float percentFreq = 0.01, int nIters = 1000){
 	try {
 		RarefactionCurveData* rcd = new RarefactionCurveData();
 		for(int i=0;i<displays.size();i++){
 			rcd->registerDisplay(displays[i]);
 		}
-	
+		
+		//convert freq percentage to number
+		int increment = 1;
+		if (percentFreq < 1.0) {  increment = numSeqs * percentFreq;  }
+		else { increment = percentFreq;  }	
+		
+		#if defined (__APPLE__) || (__MACH__) || (linux) || (__linux) || (__linux__) || (__unix__) || (__unix)
+				if(processors == 1){
+					driver(rcd, increment, nIters);	
+				}else{
+					vector<int> procIters;
+					
+					int numItersPerProcessor = nIters / processors;
+					
+					//divide iters between processes
+					for (int i = 0; i < processors; i++) {
+						if(i == processors - 1){
+							numItersPerProcessor = nIters - i * numItersPerProcessor;
+						}
+						procIters.push_back(numItersPerProcessor);
+					}
+					
+					createProcesses(procIters, rcd, increment); 
+				}
+
+		#else
+			driver(rcd, increment, nIters);	
+		#endif
+
+		for(int i=0;i<displays.size();i++){
+			displays[i]->close();
+		}
+		
+		delete rcd;
+		return 0;
+	}
+	catch(exception& e) {
+		m->errorOut(e, "Rarefact", "getCurve");
+		exit(1);
+	}
+}
+/***********************************************************************/
+int Rarefact::driver(RarefactionCurveData* rcd, int increment, int nIters = 1000){
+	try {
+			
 		for(int iter=0;iter<nIters;iter++){
 		
 			for(int i=0;i<displays.size();i++){
@@ -31,6 +75,8 @@ void Rarefact::getCurve(int increment = 1, int nIters = 1000){
 		
 			for(int i=0;i<numSeqs;i++){
 			
+				if (m->control_pressed) { delete lookup; delete rank; delete rcd; return 0;  }
+			
 				int binNumber = order->get(i);
 				int abundance = lookup->get(binNumber);
 			
@@ -40,107 +86,161 @@ void Rarefact::getCurve(int increment = 1, int nIters = 1000){
 				lookup->set(binNumber, abundance);
 				rank->set(abundance, rank->get(abundance)+1);
 
-				if((i == 0) || (i+1) % increment == 0){
+				if((i == 0) || ((i+1) % increment == 0) || (ends.count(i+1) != 0)){
 					rcd->updateRankData(rank);
 				}
 			}
 	
-			if(numSeqs % increment != 0){
+			if((numSeqs % increment != 0) || (ends.count(numSeqs) != 0)){
 				rcd->updateRankData(rank);
 			}
 
 			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();
-		}
+		return 0;
 	}
 	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";
+		m->errorOut(e, "Rarefact", "driver");
 		exit(1);
 	}
-	catch(...) {
-		cout << "An unknown error has occurred in the Rarefact class function getCurve. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";
+}
+/**************************************************************************************************/
+
+int Rarefact::createProcesses(vector<int>& procIters, RarefactionCurveData* rcd, int increment) {
+	try {
+#if defined (__APPLE__) || (__MACH__) || (linux) || (__linux) || (__linux__) || (__unix__) || (__unix)
+		int process = 1;
+		
+		vector<int> processIDS;
+		
+		EstOutput results;
+		
+		//loop through and create all the processes you want
+		while (process != processors) {
+			int pid = fork();
+			
+			if (pid > 0) {
+				processIDS.push_back(pid);  //create map from line number to pid so you can append files in correct order later
+				process++;
+			}else if (pid == 0){
+				driver(rcd, increment, procIters[process]);
+			
+				//pass numSeqs to parent
+				for(int i=0;i<displays.size();i++){
+					string tempFile = toString(getpid()) + toString(i) + ".rarefact.temp";
+					displays[i]->outputTempFiles(tempFile);
+				}
+				exit(0);
+			}else { 
+				m->mothurOut("[ERROR]: unable to spawn the necessary processes."); m->mothurOutEndLine(); 
+				for (int i = 0; i < processIDS.size(); i++) { kill (processIDS[i], SIGINT); }
+				exit(0);
+			}
+		}
+		
+		driver(rcd, increment, procIters[0]);
+		
+		//force parent to wait until all the processes are done
+		for (int i=0;i<(processors-1);i++) { 
+			int temp = processIDS[i];
+			wait(&temp);
+		}
+		
+		//get data created by processes
+		for (int i=0;i<(processors-1);i++) { 
+			for(int j=0;j<displays.size();j++){
+				string s = toString(processIDS[i]) + toString(j) + ".rarefact.temp";
+				displays[j]->inputTempFiles(s);
+				m->mothurRemove(s);
+			}
+		}
+		
+		return 0;
+#endif		
+	}
+	catch(exception& e) {
+		m->errorOut(e, "Rarefact", "createProcesses");
 		exit(1);
 	}
 }
-
 /***********************************************************************/
 
-void Rarefact::getSharedCurve(int increment = 1, int nIters = 1000){
+int Rarefact::getSharedCurve(float percentFreq = 0.01, 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]);
 		}
-
+		
+		//if jumble is false all iters will be the same
+		if (m->jumble == false)  {  nIters = 1;  }
+		
+		//convert freq percentage to number
+		int increment = 1;
+		if (percentFreq < 1.0) {  increment = numSeqs * percentFreq;  }
+		else { increment = percentFreq;  }
+		
 		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;
-					}
-				}
+			
+			if (m->jumble == true)  {
+				//randomize the groups
+				random_shuffle(lookup.begin(), lookup.end());
 			}
 			
-			//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]);
+				if (m->control_pressed) {  delete merge; delete rcd; return 0;  }
+				
+				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;
+		return 0;
 	}
 	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";
+		m->errorOut(e, "Rarefact", "getSharedCurve");
 		exit(1);
 	}
-	catch(...) {
-		cout << "An unknown error has occurred in the Rarefact class function getSharedCurve. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";
-		exit(1);
-	}
-
 }
 
 /**************************************************************************************/
@@ -148,18 +248,12 @@ 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) {
-		cout << "Standard Error: " << e.what() << " has occurred in the Rarefact class Function mergeVectors. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";
+		m->errorOut(e, "Rarefact", "mergeVectors");
 		exit(1);
 	}
-	catch(...) {
-		cout << "An unknown error has occurred in the Rarefact class function mergeVectors. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";
-		exit(1);
-	}	
 }