X-Git-Url: https://git.donarmstrong.com/?a=blobdiff_plain;ds=sidebyside;f=shhhercommand.h;h=b3a0071079bdbcbdb6a1069cf6bc1d2f53f8fad9;hb=8dd3c225255d7084e3aff8740aa4f1f1cabb367a;hp=9e50824b05f9c70a151a903b8119ba09c447b963;hpb=2a99ff8ba22dda9d0c1b125e5c9ff91e23ed6879;p=mothur.git

diff --git a/shhhercommand.h b/shhhercommand.h
index 9e50824..b3a0071 100644
--- a/shhhercommand.h
+++ b/shhhercommand.h
@@ -12,47 +12,60 @@
 
 #include "mothur.h"
 #include "command.hpp"
-#include "globaldata.hpp"
+
+//**********************************************************************************************************************
+
+#define NUMBINS 1000
+#define HOMOPS 10
+#define MIN_COUNT 0.1
+#define MIN_WEIGHT 0.1
+#define MIN_TAU 0.0001
+#define MIN_ITER 10
+//**********************************************************************************************************************
 
 class ShhherCommand : public Command {
 	
 public:
 	ShhherCommand(string);
 	ShhherCommand();
-	~ShhherCommand();
-	vector<string> getRequiredParameters();
-	vector<string> getValidParameters();
-	vector<string> getRequiredFiles();
-	map<string, vector<string> > getOutputFiles() { return outputTypes; }
-	int execute();
-	void help();
+	~ShhherCommand() {}
 	
+	vector<string> setParameters();
+	string getCommandName()			{ return "shhh.flows";	}
+	string getCommandCategory()		{ return "Sequence Processing";		}
+	string getHelpString();	
+	string getCitation() { return "http://www.mothur.org/wiki/Shhh.flows"; }
+	string getDescription()		{ return "shhh.flows"; }
+
+	
+	int execute(); 
+	void help() { m->mothurOut(getHelpString()); }		
 private:
-	GlobalData* globaldata;
 	
 	int abort;
-	map<string, vector<string> > outputTypes;
 	
-	string outputDir, flowFileName, flowFilesFileName, lookupFileName, compositeFASTAFileName;
+	string outputDir, flowFileName, flowFilesFileName, lookupFileName, compositeFASTAFileName, compositeNamesFileName;
 
 	int processors, maxIters;
 	float cutoff, sigma, minDelta;
+	string flowOrder;
 	
 	vector<int> nSeqsBreaks;
 	vector<int> nOTUsBreaks;
-	vector<double> flowDataPrI;
-	vector<short> flowDataIntI;
-	vector<int> lengths;
-	vector<string> seqNameVector;
 	vector<double> singleLookUp;
 	vector<double> jointLookUp;
+	
+	vector<string> seqNameVector;
+	vector<int> lengths;
+	vector<short> flowDataIntI;
+	vector<double> flowDataPrI;
 	map<string, int> nameMap;
 	vector<int> otuData;
 	vector<int> cumNumSeqs;
 	vector<int> nSeqsPerOTU;
 	vector<vector<int> > aaP;	//tMaster->aanP:	each row is a different otu / each col contains the sequence indices
-	vector<int> seqNumber;		//tMaster->anP:		the sequence id number sorted by OTU
 	vector<vector<int> > aaI;	//tMaster->aanI:	that are in each otu - can't differentiate between aaP and aaI 
+	vector<int> seqNumber;		//tMaster->anP:		the sequence id number sorted by OTU
 	vector<int> seqIndex;		//tMaster->anI;		the index that corresponds to seqNumber
 	vector<double> dist;		//adDist - distance of sequences to centroids
 	vector<short> change;		//did the centroid sequence change? 0 = no; 1 = yes
@@ -61,10 +74,12 @@ private:
 	vector<double> singleTau;	//tMaster->adTau:	1-D Tau vector (1xnumSeqs)
 	vector<short> uniqueFlowgrams;
 	vector<int> uniqueCount;
-	vector<int> uniqueLengths;
 	vector<int> mapSeqToUnique;
 	vector<int> mapUniqueToSeq;
-	
+	vector<int> uniqueLengths;
+
+	vector<string> outputNames;
+
 	int numSeqs, numUniques, numOTUs, numFlowCells;
 	
 	void getSingleLookUp();
@@ -110,6 +125,125 @@ private:
 	
 };
 
+/**************************************************************************************************/
+//custom data structure for threads to use.
+// This is passed by void pointer so it can be any data type
+// that can be passed using a single void pointer (LPVOID).
+struct flowDistParentForkData {
+	string distFileName; 
+	vector<int> mapUniqueToSeq;
+	vector<int> mapSeqToUnique;
+	vector<int> lengths;
+	vector<short> flowDataIntI;
+	vector<double> flowDataPrI;
+	vector<double> jointLookUp;
+	MothurOut* m;
+	int threadID, startSeq, stopSeq, numFlowCells;
+	float cutoff;
+	
+	flowDistParentForkData(){}
+	flowDistParentForkData(string d, vector<int> mapU, vector<int> mapS, vector<int> l, vector<short> flowD, vector<double> flowDa, vector<double> j, MothurOut* mout, int st, int sp, int n, float cut, int tid) {
+		distFileName = d;
+		mapUniqueToSeq = mapU;
+		mapSeqToUnique = mapS;
+		lengths = l;
+		flowDataIntI = flowD;
+		flowDataPrI = flowDa;
+		jointLookUp = j;
+		m = mout;
+		startSeq = st;
+		stopSeq = sp;
+		numFlowCells = n;
+		cutoff= cut;
+		threadID = tid;
+	}
+};
+
+/**************************************************************************************************/
+#if defined (__APPLE__) || (__MACH__) || (linux) || (__linux)
+#else
+static DWORD WINAPI MyflowDistParentForkThreadFunction(LPVOID lpParam){ 
+	flowDistParentForkData* pDataArray;
+	pDataArray = (flowDistParentForkData*)lpParam;
+	
+	try {
+		ostringstream outStream;
+		outStream.setf(ios::fixed, ios::floatfield);
+		outStream.setf(ios::dec, ios::basefield);
+		outStream.setf(ios::showpoint);
+		outStream.precision(6);
+		
+		int begTime = time(NULL);
+		double begClock = clock();
+		string tempOut = "start and end = " + toString(pDataArray->startSeq) +'\t' + toString(pDataArray->stopSeq) + "-";
+		cout << tempOut << endl;
+	
+		for(int i=pDataArray->startSeq;i<pDataArray->stopSeq;i++){
+			
+			if (pDataArray->m->control_pressed) { break; }
+			cout << "thread i = " << i << endl;
+			for(int j=0;j<i;j++){
+				
+				cout << "thread j = " << j << endl;
+				float flowDistance = 0.0;
+				////////////////// calcPairwiseDist ///////////////////
+				//needed because this is a static global function that can't see the classes internal functions
+				
+				int minLength = pDataArray->lengths[pDataArray->mapSeqToUnique[pDataArray->mapUniqueToSeq[i]]];
+				if(pDataArray->lengths[pDataArray->mapUniqueToSeq[j]] < minLength){	minLength = pDataArray->lengths[pDataArray->mapSeqToUnique[pDataArray->mapUniqueToSeq[j]]];	}
+				
+				int ANumFlowCells = pDataArray->mapUniqueToSeq[i] * pDataArray->numFlowCells;
+				int BNumFlowCells = pDataArray->mapUniqueToSeq[j] * pDataArray->numFlowCells;
+				
+				for(int k=0;k<minLength;k++){
+					
+					if (pDataArray->m->control_pressed) { break; }
+					
+					int flowAIntI = pDataArray->flowDataIntI[ANumFlowCells + k];
+					float flowAPrI = pDataArray->flowDataPrI[ANumFlowCells + k];
+					
+					int flowBIntI = pDataArray->flowDataIntI[BNumFlowCells + k];
+					float flowBPrI = pDataArray->flowDataPrI[BNumFlowCells + k];
+					flowDistance += pDataArray->jointLookUp[flowAIntI * NUMBINS + flowBIntI] - flowAPrI - flowBPrI;
+				}
+				
+				flowDistance /= (float) minLength;
+				//cout << flowDistance << endl;
+				////////////////// end of calcPairwiseDist ///////////////////
+								
+				if(flowDistance < 1e-6){
+					outStream << pDataArray->mapUniqueToSeq[i] << '\t' << pDataArray->mapUniqueToSeq[j] << '\t' << 0.000000 << endl;
+				}
+				else if(flowDistance <= pDataArray->cutoff){
+					outStream << pDataArray->mapUniqueToSeq[i] << '\t' << pDataArray->mapUniqueToSeq[j] << '\t' << flowDistance << endl;
+				}
+			}
+			if(i % 100 == 0){
+				pDataArray->m->mothurOut(toString(i) + "\t" + toString(time(NULL) - begTime));
+				pDataArray->m->mothurOut("\t" + toString((clock()-begClock)/CLOCKS_PER_SEC));
+				pDataArray->m->mothurOutEndLine();
+			}
+		}
+		
+		ofstream distFile(pDataArray->distFileName.c_str());
+		distFile << outStream.str();		
+		distFile.close();
+		
+		if (pDataArray->m->control_pressed) {}
+		else {
+			pDataArray->m->mothurOut(toString(pDataArray->stopSeq-1) + "\t" + toString(time(NULL) - begTime));
+			pDataArray->m->mothurOut("\t" + toString((clock()-begClock)/CLOCKS_PER_SEC));
+			pDataArray->m->mothurOutEndLine();
+		}		
+		
+	}
+	catch(exception& e) {
+		pDataArray->m->errorOut(e, "ShhherCommand", "MyflowDistParentForkThreadFunction");
+		exit(1);
+	}
+} 
+#endif
+
 
 #endif