X-Git-Url: https://git.donarmstrong.com/?a=blobdiff_plain;f=shhhercommand.h;h=b3a0071079bdbcbdb6a1069cf6bc1d2f53f8fad9;hb=2bb9267aa4b4ecdf8488b06605cc9f3f36fa4332;hp=8d29d86e0ec606379e471a363dbe33af34336d8e;hpb=ca9ac1d80c62f57270b0dcd49410ebe08a8aecd6;p=mothur.git

diff --git a/shhhercommand.h b/shhhercommand.h
index 8d29d86..b3a0071 100644
--- a/shhhercommand.h
+++ b/shhhercommand.h
@@ -13,6 +13,15 @@
 #include "mothur.h"
 #include "command.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 {
 	
@@ -22,9 +31,12 @@ public:
 	~ShhherCommand() {}
 	
 	vector<string> setParameters();
-	string getCommandName()			{ return "shhh.seqs";	}
-	string getCommandCategory()		{ return "Hidden";		}
+	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()); }		
@@ -32,7 +44,7 @@ private:
 	
 	int abort;
 	
-	string outputDir, flowFileName, flowFilesFileName, lookupFileName, compositeFASTAFileName;
+	string outputDir, flowFileName, flowFilesFileName, lookupFileName, compositeFASTAFileName, compositeNamesFileName;
 
 	int processors, maxIters;
 	float cutoff, sigma, minDelta;
@@ -65,7 +77,9 @@ private:
 	vector<int> mapSeqToUnique;
 	vector<int> mapUniqueToSeq;
 	vector<int> uniqueLengths;
-	
+
+	vector<string> outputNames;
+
 	int numSeqs, numUniques, numOTUs, numFlowCells;
 	
 	void getSingleLookUp();
@@ -111,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