]> git.donarmstrong.com Git - mothur.git/commitdiff
finished with ccode, returned bellerophon to last save before move, cleaned up pintai...
authorwestcott <westcott>
Tue, 8 Sep 2009 13:53:16 +0000 (13:53 +0000)
committerwestcott <westcott>
Tue, 8 Sep 2009 13:53:16 +0000 (13:53 +0000)
bellerophon.cpp
ccode.cpp
ccode.h
chimera.cpp
decalc.cpp
pintail.cpp

index 2ce395189520a83837805a8322f2619a5f227508..9aa4f7ac7aa566f39a9b68c7b9e0103f516e0ba4 100644 (file)
@@ -35,7 +35,7 @@ void Bellerophon::print(ostream& out) {
                        out << pref[i].name << '\t' << setprecision(3) << pref[i].score[0] << '\t' << pref[i].leftParent[0] << '\t' << pref[i].rightParent[0] << endl;
                        
                        //calc # of seqs with preference above 1.0
-                       if (pref[i].score[0] > 1.1) { 
+                       if (pref[i].score[0] > 1.0) { 
                                above1++; 
                                mothurOut(pref[i].name + " is a suspected chimera at breakpoint " + toString(pref[i].midpoint)); mothurOutEndLine();
                                mothurOut("It's score is " + toString(pref[i].score[0]) + " with suspected left parent " + pref[i].leftParent[0] + " and right parent " + pref[i].rightParent[0]); mothurOutEndLine();
@@ -44,7 +44,7 @@ void Bellerophon::print(ostream& out) {
                
                //output results to screen
                mothurOutEndLine();
-               mothurOut("Sequence with preference score above 1.1: " + toString(above1)); mothurOutEndLine();
+               mothurOut("Sequence with preference score above 1.0: " + toString(above1)); mothurOutEndLine();
                int spot;
                spot = pref.size()-1;
                mothurOut("Minimum:\t" + toString(pref[spot].score[0])); mothurOutEndLine();
@@ -203,7 +203,7 @@ void Bellerophon::getChimeras() {
                delete distCalculator;
                
                //rank preference score to eachother
-               /*float dme = 0.0;
+               float dme = 0.0;
                float expectedPercent = 1 / (float) (pref.size());
                
                for (int i = 0; i < pref.size(); i++) {  dme += pref[i].score[0];  }
@@ -216,7 +216,7 @@ void Bellerophon::getChimeras() {
                        //how much higher or lower is this than expected
                        pref[i].score[0] = pref[i].score[0] / expectedPercent;
                
-               }*/
+               }
                
                //sort Preferences highest to lowest
                sort(pref.begin(), pref.end(), comparePref);
@@ -335,24 +335,21 @@ void Bellerophon::generatePreferences(vector<SeqMap> left, vector<SeqMap> right,
                  
                //calculate the dme
                int count0 = 0;
-               for (int i = 0; i < pref.size(); i++) {  dme += pref[i].score[1];  if (pref[i].score[1] == 0.0) {  count0++; }  }
+               for (int i = 0; i < pref.size(); i++) {  dme += pref[i].score[1];  if (pref[i].score[1] == 0.0) { count0++; }  }
                
-               //float expectedPercent = 1 / (float) (pref.size() - count0);
+               float expectedPercent = 1 / (float) (pref.size() - count0);
 //cout << endl << "dme = " << dme << endl;
                //recalculate prefernences based on dme
                for (int i = 0; i < pref.size(); i++) {
-//cout << "unadjusted col[i] " << i << " = " << pref[i].score[1] << endl;      
-//cout << i << '\t' <<  (dme / (dme - 2 * pref[i].score[1])) << endl;
-                       
+//cout << "unadjusted pref " << i << " = " << pref[i].score[1] << endl;        
                        // gives the actual percentage of the dme this seq adds
-                       //pref[i].score[1] = pref[i].score[1] / dme;
+                       pref[i].score[1] = pref[i].score[1] / dme;
                        
                        //how much higher or lower is this than expected
-                       //pref[i].score[1] = pref[i].score[1] / expectedPercent;
+                       pref[i].score[1] = pref[i].score[1] / expectedPercent;
+                       
+                       //pref[i].score[1] = dme / (dme - 2 * pref[i].score[1]);
                        
-                       //not 2 * pref[i].score[1] because i only calulate the pref scores once.
-                       pref[i].score[1] = dme / (dme - pref[i].score[1]);
-//cout << i << '\t' << pref[i].score[1] << endl;
                        //so a non chimeric sequence would be around 1, and a chimeric would be signifigantly higher.
 //cout << "adjusted pref " << i << " = " << pref[i].score[1] << endl;                                  
                }
index 7e6154b6c5af1d69259223025d9da9f7f3f5d4df..ec2e355f882d86a60aba6f520604ef06ead6e6c8 100644 (file)
--- a/ccode.cpp
+++ b/ccode.cpp
@@ -20,7 +20,6 @@ Ccode::~Ccode() {
        try {
                for (int i = 0; i < querySeqs.size(); i++)              {  delete querySeqs[i];         }
                for (int i = 0; i < templateSeqs.size(); i++)   {  delete templateSeqs[i];      }
-               delete distCalc;
        }
        catch(exception& e) {
                errorOut(e, "Ccode", "~Ccode");
@@ -202,21 +201,7 @@ void Ccode::getChimeras() {
                        mothurOut("Done."); mothurOutEndLine();
                }else {         createProcessesClosest();               }
 
-               
-for (int i = 0; i < closest.size(); i++) {
-       //cout << querySeqs[i]->getName() << ": ";
-       string file = querySeqs[i]->getName() + ".input";
-       ofstream o;
-       openOutputFile(file, o);
-       
-       querySeqs[i]->printSequence(o);
-       for (int j = 0; j < closest[i].size(); j++) {
-               //cout << closest[i][j].seq->getName() << '\t';
-               closest[i][j].seq->printSequence(o);
-       }
-       //cout << endl;
-       o.close();
-}      
+               //initialize spotMap
                for (int j = 0; j < numSeqs; j++) {
                        for (int i = 0; i < querySeqs[0]->getAligned().length(); i++) {
                                spotMap[j][i] = i;
@@ -280,39 +265,49 @@ for (int i = 0; i < closest.size(); i++) {
                        windows[i] = findWindows(i);  
                }
 
-               //remove sequences that are more than 20% different and less than 0.5% different - may want to allow user to specify this later - should be paralellized
-               for (int i = 0; i < closest.size(); i++) {
-                       removeBadReferenceSeqs(closest[i], i);
-               }
+               //remove sequences that are more than 20% different and less than 0.5% different - may want to allow user to specify this later 
+               if (processors == 1) { 
+                       for (int i = 0; i < closest.size(); i++) {
+                               removeBadReferenceSeqs(closest[i], i);
+                       }
+               }else {         createProcessesRemoveBad();             }
+
                
-       
-               //find the averages for each querys references - should be paralellized
-               for (int i = 0; i < numSeqs; i++) {
-                       getAverageRef(closest[i], i);  //fills sumRef[i], averageRef[i], sumSquaredRef[i] and refCombo[i].
-               }
-       
-               //find the averages for the query - should be paralellized
-               for (int i = 0; i < numSeqs; i++) {
-                       getAverageQuery(closest[i], i);  //fills sumQuery[i], averageQuery[i], sumSquaredQuery[i].
-               }
+               if (processors == 1) { 
+                       //find the averages for each querys references
+                       for (int i = 0; i < numSeqs; i++) {
+                               getAverageRef(closest[i], i);  //fills sumRef[i], averageRef[i], sumSquaredRef[i] and refCombo[i].
+                       }
+                       
+                       //find the averages for the query 
+                       for (int i = 0; i < numSeqs; i++) {
+                               getAverageQuery(closest[i], i);  //fills sumQuery[i], averageQuery[i], sumSquaredQuery[i].
+                       }
+               }else {         createProcessesAverages();              }
                
-               //find the averages for each querys references - should be paralellized
-               for (int i = 0; i < numSeqs; i++) {
-                       findVarianceRef(i);  //fills varRef[i] and sdRef[i] also sets minimum error rate to 0.001 to avoid divide by 0.
-               }
+               if (processors == 1) { 
+                       //find the averages for each querys references 
+                       for (int i = 0; i < numSeqs; i++) {
+                               findVarianceRef(i);  //fills varRef[i] and sdRef[i] also sets minimum error rate to 0.001 to avoid divide by 0.
+                       }
+                       
+                       //find the averages for the query 
+                       for (int i = 0; i < numSeqs; i++) {
+                               findVarianceQuery(i);  //fills varQuery[i] and sdQuery[i] also sets minimum error rate to 0.001 to avoid divide by 0.
+                       }
+               }else {         createProcessesVariances();             }
+               
+               if (processors == 1) { 
+                       for (int i = 0; i < numSeqs; i++) {
+                               determineChimeras(i);  //fills anova, isChimericConfidence[i], isChimericTStudent[i] and isChimericANOVA[i]. 
+                       }
+               }else {         createProcessesDetermine();             }
                
-               //find the averages for the query - should be paralellized
-               for (int i = 0; i < numSeqs; i++) {
-                       findVarianceQuery(i);  //fills varQuery[i] and sdQuery[i] also sets minimum error rate to 0.001 to avoid divide by 0.
-               }
-
-               for (int i = 0; i < numSeqs; i++) {
-                       determineChimeras(i);  //fills anova, isChimericConfidence[i], isChimericTStudent[i] and isChimericANOVA[i]. - should be paralellized
-               }
-
                //free memory
                for (int i = 0; i < lines.size(); i++)                                  {       delete lines[i];                                }
                for (int i = 0; i < templateLines.size(); i++)                  {       delete templateLines[i];                }
+               delete distCalc;
+               delete decalc;
                        
        }
        catch(exception& e) {
@@ -908,6 +903,8 @@ void Ccode::determineChimeras (int query) {
                       isChimericANOVA[query][i] = true;   /* significant P<0.05 */
                }
                        
+                       if (isnan(anovaValue) || isinf(anovaValue)) { anovaValue = 0.0; }
+                       
                        anova[query][i] = anovaValue;
                }
                
@@ -1045,7 +1042,7 @@ void Ccode::createProcessesClosest() {
                                ofstream out;
                                string s = toString(getpid()) + ".temp";
                                openOutputFile(s, out);
-                               
+                                                       
                                //output pairs
                                for (int i = lines[process]->start; i < lines[process]->end; i++) {
                                         for (int j = 0; j < closest[i].size(); j++) {
@@ -1060,9 +1057,7 @@ void Ccode::createProcessesClosest() {
                                         }
                                         out << endl;
                                }
-
                                out.close();
-                               
                                exit(0);
                        }else { mothurOut("unable to spawn the necessary processes."); mothurOutEndLine(); exit(0); }
                }
@@ -1083,7 +1078,7 @@ void Ccode::createProcessesClosest() {
                        //get pairs
                        for (int k = lines[i]->start; k < lines[i]->end; k++) {
                                vector<Sequence*> tempVector;
-                               
+       
                                for (int j = 0; j < numWanted; j++) {
                                
                                        Sequence* temp = new Sequence(in);
@@ -1097,21 +1092,23 @@ void Ccode::createProcessesClosest() {
                        
                        string junk;
                        in >> junk;  gobble(in);  // to get ">"
-                       
+               
                        vector< vector<float> > dists; dists.resize(querySeqs.size());
                        
-                       for (int i = lines[process]->start; i < lines[process]->end; i++) {
-                               dists[i].resize(closest[i].size());
-                               for (int j = 0; j < closest[i].size(); j++) {
-                                       in >> dists[i][j];
+                       for (int k = lines[i]->start; k < lines[i]->end; k++) {
+                               dists[k].resize(numWanted);
+                               for (int j = 0; j < numWanted; j++) {
+                                       in >> dists[k][j];
                                }
                                gobble(in);
-                       } 
                        
-                       for (int i = lines[process]->start; i < lines[process]->end; i++) {
-                               for (int j = 0; j < closest[i].size(); j++) {
-                                       closest[i][j].seq = tempClosest[i][j];
-                                       closest[i][j].dist = dists[i][j]; 
+                       } 
+
+                       for (int k = lines[i]->start; k < lines[i]->end; k++) {
+                               closest[k].resize(numWanted);
+                               for (int j = 0; j < closest[k].size(); j++) {
+                                       closest[k][j].seq = tempClosest[k][j];
+                                       closest[k][j].dist = dists[k][j]; 
                                }
                        } 
 
@@ -1132,4 +1129,486 @@ void Ccode::createProcessesClosest() {
 }
 
 //***************************************************************************************************************
+void Ccode::createProcessesRemoveBad() {
+       try {
+#if defined (__APPLE__) || (__MACH__) || (linux) || (__linux)
+               int process = 0;
+               vector<int> processIDS;
+               
+               //loop through and create all the processes you want
+               while (process != processors) {
+                       int pid = fork();
+                       
+                       if (pid > 0) {
+                               processIDS.push_back(pid);  
+                               process++;
+                       }else if (pid == 0){
+                                                       
+                               for (int i = lines[process]->start; i < lines[process]->end; i++) {
+                                       removeBadReferenceSeqs(closest[i], i);
+                               }
+                               
+                               //write out data to file so parent can read it
+                               ofstream out;
+                               string s = toString(getpid()) + ".temp";
+                               openOutputFile(s, out);
+                               
+                               //output pairs
+                               for (int i = lines[process]->start; i < lines[process]->end; i++) {
+                                       out << closest[i].size() << endl;
+                                        for (int j = 0; j < closest[i].size(); j++) {
+                                               closest[i][j].seq->printSequence(out);
+                                        }
+                                        out << ">" << endl; //to stop sequence read
+                               }
+                               
+                               for (int i = lines[process]->start; i < lines[process]->end; i++) {
+                                        for (int j = 0; j < closest[i].size(); j++) {
+                                               out << closest[i][j].dist << '\t';
+                                        }
+                                        out << endl;
+                               }
+
+                               out.close();
+                               
+                               exit(0);
+                       }else { mothurOut("unable to spawn the necessary processes."); mothurOutEndLine(); exit(0); }
+               }
+               
+               //force parent to wait until all the processes are done
+               for (int i=0;i<processors;i++) { 
+                       int temp = processIDS[i];
+                       wait(&temp);
+               }
+       
+               //get data created by processes
+               for (int i=0;i<processors;i++) { 
+                       ifstream in;
+                       string s = toString(processIDS[i]) + ".temp";
+                       openInputFile(s, in);
+               
+                       vector< vector<Sequence*> > tempClosest;  tempClosest.resize(querySeqs.size());
+                       vector<int> sizes; 
+                       //get pairs
+                       for (int k = lines[i]->start; k < lines[i]->end; k++) {
+               
+                               int num;
+                               in >> num; 
+                               sizes.push_back(num);
+                               gobble(in);
+                               
+                               vector<Sequence*> tempVector;
+                               
+                               for (int j = 0; j < num; j++) {
+                               
+                                       Sequence* temp = new Sequence(in);
+                                       gobble(in);
+                                               
+                                       tempVector.push_back(temp);
+                               }
+                               string junk;
+                               in >> junk;  gobble(in);  // to get ">"
+
+                               tempClosest[k] = tempVector;
+                       }
+                       
+                       vector< vector<float> > dists; dists.resize(querySeqs.size());
+                       int count = 0;
+                       for (int k = lines[i]->start; k < lines[i]->end; k++) {
+                               dists[k].resize(sizes[count]);
+                               for (int j = 0; j < sizes[count]; j++) {
+                                       in >> dists[k][j];
+                               }
+                               gobble(in);
+                               count++;
+                       } 
+                       
+                       count = 0;
+                       for (int k = lines[i]->start; k < lines[i]->end; k++) {
+                               for (int j = 0; j < sizes[count]; j++) {
+                                       closest[k][j].seq = tempClosest[k][j];
+                                       closest[k][j].dist = dists[k][j]; 
+                               }
+                               count++;
+                       } 
+
+                       in.close();
+                       remove(s.c_str());
+               }
+#else
+               for (int i = 0; i < closest.size(); i++) {
+                       removeBadReferenceSeqs(closest[i], i);
+               }
+#endif         
+       
+       }
+       catch(exception& e) {
+               errorOut(e, "Ccode", "createProcessesRemoveBad");
+               exit(1);
+       }
+}
+//***************************************************************************************************************
+void Ccode::createProcessesAverages() {
+       try {
+       #if defined (__APPLE__) || (__MACH__) || (linux) || (__linux)
+               int process = 0;
+               vector<int> processIDS;
+               
+               //loop through and create all the processes you want
+               while (process != processors) {
+                       int pid = fork();
+                       
+                       if (pid > 0) {
+                               processIDS.push_back(pid);  
+                               process++;
+                       }else if (pid == 0){
+                                                       
+                               //find the averages for each querys references 
+                               for (int i = lines[process]->start; i < lines[process]->end; i++) {
+                                       getAverageRef(closest[i], i);  //fills sumRef[i], averageRef[i], sumSquaredRef[i] and refCombo[i].
+                               }
+                               
+                               //find the averages for the query 
+                               for (int i = lines[process]->start; i < lines[process]->end; i++) {
+                                       getAverageQuery(closest[i], i);  //fills sumQuery[i], averageQuery[i], sumSquaredQuery[i].
+                               }
+                               
+                               
+                               //write out data to file so parent can read it
+                               ofstream out;
+                               string s = toString(getpid()) + ".temp";
+                               openOutputFile(s, out);
+                               
+                               //output pairs
+                               for (int i = lines[process]->start; i < lines[process]->end; i++) {
+                                       for (int j = 0; j < windows[i].size(); j++) {
+                                               out << sumRef[i][j] << '\t';
+                                       }
+                                       out << endl;
+                                       for (int j = 0; j < windows[i].size(); j++) {
+                                               out << averageRef[i][j] << '\t';
+                                       }
+                                       out << endl;
+                                       for (int j = 0; j < windows[i].size(); j++) {
+                                               out << sumSquaredRef[i][j] << '\t';
+                                       }
+                                       out << endl;
+                                       for (int j = 0; j < windows[i].size(); j++) {
+                                               out << sumQuery[i][j] << '\t';
+                                       }
+                                       out << endl;
+                                       for (int j = 0; j < windows[i].size(); j++) {
+                                               out << averageQuery[i][j] << '\t';
+                                       }
+                                       out << endl;
+                                       for (int j = 0; j < windows[i].size(); j++) {
+                                               out << sumSquaredQuery[i][j] << '\t';
+                                       }
+                                       out << endl;
+                                       out << refCombo[i] << endl;
+                               }
+
+                               out.close();
+                               
+                               exit(0);
+                       }else { mothurOut("unable to spawn the necessary processes."); mothurOutEndLine(); exit(0); }
+               }
+               
+               //force parent to wait until all the processes are done
+               for (int i=0;i<processors;i++) { 
+                       int temp = processIDS[i];
+                       wait(&temp);
+               }
+       
+               //get data created by processes
+               for (int i=0;i<processors;i++) { 
+                       ifstream in;
+                       string s = toString(processIDS[i]) + ".temp";
+                       openInputFile(s, in);
+               
+                       //get pairs
+                       for (int k = lines[i]->start; k < lines[i]->end; k++) {
+                               sumRef[k].resize(windows[k].size());
+                               averageRef[k].resize(windows[k].size());
+                               sumSquaredRef[k].resize(windows[k].size());
+                               averageQuery[k].resize(windows[k].size());
+                               sumQuery[k].resize(windows[k].size());
+                               sumSquaredQuery[k].resize(windows[k].size());
+                               
+                               for (int j = 0; j < windows[k].size(); j++) {
+                                       in >> sumRef[k][j];
+                               }
+                               gobble(in);
+                               for (int j = 0; j < windows[k].size(); j++) {
+                                       in >> averageRef[k][j];
+                               }
+                               gobble(in);
+                               for (int j = 0; j < windows[k].size(); j++) {
+                                       in >> sumSquaredRef[k][j];
+                               }
+                               gobble(in);
+                               for (int j = 0; j < windows[k].size(); j++) {
+                                       in >> sumQuery[k][j];
+                               }
+                               gobble(in);
+                               for (int j = 0; j < windows[k].size(); j++) {
+                                       in >> averageQuery[k][j];
+                               }
+                               gobble(in);
+                               for (int j = 0; j < windows[k].size(); j++) {
+                                       in >> sumSquaredQuery[k][j];
+                               }
+                               gobble(in);
+                               in >> refCombo[k];
+                               gobble(in);
+                       }
+                       
+                       in.close();
+                       remove(s.c_str());
+               }
+#else
+               //find the averages for each querys references 
+               for (int i = 0; i < querySeqs.size(); i++) {
+                       getAverageRef(closest[i], i);  //fills sumRef[i], averageRef[i], sumSquaredRef[i] and refCombo[i].
+               }
+       
+               //find the averages for the query 
+               for (int i = 0; i < querySeqs.size(); i++) {
+                       getAverageQuery(closest[i], i);  //fills sumQuery[i], averageQuery[i], sumSquaredQuery[i].
+               }
+
+#endif         
+       
+       }
+       catch(exception& e) {
+               errorOut(e, "Ccode", "createProcessesAverages");
+               exit(1);
+       }
+}
+//***************************************************************************************************************
+void Ccode::createProcessesVariances() {
+       try {
+       #if defined (__APPLE__) || (__MACH__) || (linux) || (__linux)
+               int process = 0;
+               vector<int> processIDS;
+               
+               //loop through and create all the processes you want
+               while (process != processors) {
+                       int pid = fork();
+                       
+                       if (pid > 0) {
+                               processIDS.push_back(pid);  
+                               process++;
+                       }else if (pid == 0){
+                                                       
+                               //find the averages for each querys references 
+                               for (int i = lines[process]->start; i < lines[process]->end; i++) {
+                                       findVarianceRef(i);  //fills varRef[i] and sdRef[i] also sets minimum error rate to 0.001 to avoid divide by 0.
+                               }
+                               
+                               //find the averages for the query 
+                               for (int i = lines[process]->start; i < lines[process]->end; i++) {
+                                       findVarianceQuery(i);  //fills varQuery[i] and sdQuery[i] also sets minimum error rate to 0.001 to avoid divide by 0.
+                               }
+                               
+                               
+                               //write out data to file so parent can read it
+                               ofstream out;
+                               string s = toString(getpid()) + ".temp";
+                               openOutputFile(s, out);
+                               
+                               //output pairs
+                               for (int i = lines[process]->start; i < lines[process]->end; i++) {
+                                       for (int j = 0; j < windows[i].size(); j++) {
+                                               out << varRef[i][j] << '\t';
+                                       }
+                                       out << endl;
+                                       for (int j = 0; j < windows[i].size(); j++) {
+                                               out << sdRef[i][j] << '\t';
+                                       }
+                                       out << endl;
+                                       for (int j = 0; j < windows[i].size(); j++) {
+                                               out << varQuery[i][j] << '\t';
+                                       }
+                                       out << endl;
+                                       for (int j = 0; j < windows[i].size(); j++) {
+                                               out << sdQuery[i][j] << '\t';
+                                       }
+                                       out << endl;
+                               }
+
+                               out.close();
+                               
+                               exit(0);
+                       }else { mothurOut("unable to spawn the necessary processes."); mothurOutEndLine(); exit(0); }
+               }
+               
+               //force parent to wait until all the processes are done
+               for (int i=0;i<processors;i++) { 
+                       int temp = processIDS[i];
+                       wait(&temp);
+               }
+       
+               //get data created by processes
+               for (int i=0;i<processors;i++) { 
+                       ifstream in;
+                       string s = toString(processIDS[i]) + ".temp";
+                       openInputFile(s, in);
+               
+                       //get pairs
+                       for (int k = lines[i]->start; k < lines[i]->end; k++) {
+                               varRef[k].resize(windows[k].size());
+                               sdRef[k].resize(windows[k].size());
+                               varQuery[k].resize(windows[k].size());
+                               sdQuery[k].resize(windows[k].size());
+                                                               
+                               for (int j = 0; j < windows[k].size(); j++) {
+                                       in >> varRef[k][j];
+                               }
+                               gobble(in);
+                               for (int j = 0; j < windows[k].size(); j++) {
+                                       in >> sdRef[k][j];
+                               }
+                               gobble(in);
+                               for (int j = 0; j < windows[k].size(); j++) {
+                                       in >> varQuery[k][j];
+                               }
+                               gobble(in);
+                               for (int j = 0; j < windows[k].size(); j++) {
+                                       in >> sdQuery[k][j];
+                               }
+                               gobble(in);
+                       }
+                       
+                       in.close();
+                       remove(s.c_str());
+               }
+#else
+                       //find the averages for each querys references 
+                       for (int i = 0; i < querySeqs.size(); i++) {
+                               findVarianceRef(i);  //fills varRef[i] and sdRef[i] also sets minimum error rate to 0.001 to avoid divide by 0.
+                       }
+                       
+                       //find the averages for the query 
+                       for (int i = 0; i < querySeqs.size(); i++) {
+                               findVarianceQuery(i);  //fills varQuery[i] and sdQuery[i] also sets minimum error rate to 0.001 to avoid divide by 0.
+                       }
+#endif         
+       }
+       catch(exception& e) {
+               errorOut(e, "Ccode", "createProcessesVariances");
+               exit(1);
+       }
+}
+//***************************************************************************************************************
+void Ccode::createProcessesDetermine() {
+       try {
+       #if defined (__APPLE__) || (__MACH__) || (linux) || (__linux)
+               int process = 0;
+               vector<int> processIDS;
+               
+               //loop through and create all the processes you want
+               while (process != processors) {
+                       int pid = fork();
+                       
+                       if (pid > 0) {
+                               processIDS.push_back(pid);  
+                               process++;
+                       }else if (pid == 0){
+                                                       
+                               //find the averages for each querys references 
+                               for (int i = lines[process]->start; i < lines[process]->end; i++) {
+                                       determineChimeras(i);  //fills anova, isChimericConfidence[i], isChimericTStudent[i] and isChimericANOVA[i]. 
+                               }
+
+                               //write out data to file so parent can read it
+                               ofstream out;
+                               string s = toString(getpid()) + ".temp";
+                               openOutputFile(s, out);
+                               
+                               //output pairs
+                               for (int i = lines[process]->start; i < lines[process]->end; i++) {
+                                       for (int j = 0; j < windows[i].size(); j++) {
+                                               out << anova[i][j] << '\t';
+                                       }
+                                       out << endl;
+                                       for (int j = 0; j < windows[i].size(); j++) {
+                                               out << isChimericConfidence[i][j] << '\t';
+                                       }
+                                       out << endl;
+                                       for (int j = 0; j < windows[i].size(); j++) {
+                                               out << isChimericTStudent[i][j] << '\t';
+                                       }
+                                       out << endl;
+                                       for (int j = 0; j < windows[i].size(); j++) {
+                                               out << isChimericANOVA[i][j] << '\t';
+                                       }
+                                       out << endl;
+                               }
+
+                               out.close();
+                               
+                               exit(0);
+                       }else { mothurOut("unable to spawn the necessary processes."); mothurOutEndLine(); exit(0); }
+               }
+               
+               //force parent to wait until all the processes are done
+               for (int i=0;i<processors;i++) { 
+                       int temp = processIDS[i];
+                       wait(&temp);
+               }
+       
+               //get data created by processes
+               for (int i=0;i<processors;i++) { 
+                       ifstream in;
+                       string s = toString(processIDS[i]) + ".temp";
+                       openInputFile(s, in);
+               
+                       //get pairs
+                       for (int k = lines[i]->start; k < lines[i]->end; k++) {
+                               anova[k].resize(windows[k].size());
+                               isChimericConfidence[k].resize(windows[k].size(), false);
+                               isChimericTStudent[k].resize(windows[k].size(), false);
+                               isChimericANOVA[k].resize(windows[k].size(), false);
+                               int tempBool;
+                                                               
+                               for (int j = 0; j < windows[k].size(); j++) {
+                                       in >> anova[k][j];
+                               }
+                               gobble(in);
+                               for (int j = 0; j < windows[k].size(); j++) {
+                                       in >> tempBool;
+                                       if (tempBool == 1) {  isChimericConfidence[k][j] = true;  }
+                               }
+                               gobble(in);
+                               for (int j = 0; j < windows[k].size(); j++) {
+                                       in >> tempBool;
+                                       if (tempBool == 1) {  isChimericTStudent[k][j] = true;  }
+                               }
+                               gobble(in);
+                               for (int j = 0; j < windows[k].size(); j++) {
+                                       in >> tempBool;
+                                       if (tempBool == 1) {  isChimericANOVA[k][j] = true;  }
+                               }
+                               gobble(in);
+                       }
+                       
+                       in.close();
+                       remove(s.c_str());
+               }
+       #else
+                       for (int i = 0; i < querySeqs.size(); i++) {
+                               determineChimeras(i);  //fills anova, isChimericConfidence[i], isChimericTStudent[i] and isChimericANOVA[i].
+                       }
+       #endif          
+
+       }
+       catch(exception& e) {
+               errorOut(e, "Ccode", "createProcessesDetermine");
+               exit(1);
+       }
+}
+//***************************************************************************************************************
+
+
 
diff --git a/ccode.h b/ccode.h
index 6322292c12534fc44838255ed6490da4292bbf46..2888aea3b7cc8bf2e063c078b7abbaae43f1d62c 100644 (file)
--- a/ccode.h
+++ b/ccode.h
@@ -83,7 +83,11 @@ class Ccode : public Chimera {
                float getF(int); 
                
                void createProcessesClosest();
-               
+               void createProcessesRemoveBad();
+               void createProcessesAverages();
+               void createProcessesVariances();
+               void createProcessesDetermine();
+                               
 };
 
 /***********************************************************/
index 75e81c076a9ce113ff26e1d62027b86b16571527..7e997b9824439e4e42f495f95a2d101b58a89600 100644 (file)
@@ -142,6 +142,10 @@ vector< vector<float> > Chimera::readQuantiles() {
                openInputFile(quanfile, in);
                
                vector< vector<float> > quan;
+               vector <float> temp;
+               
+               //to fill 0
+               quan.push_back(temp); 
        
                int num; float ten, twentyfive, fifty, seventyfive, ninetyfive, ninetynine; 
                
@@ -149,7 +153,7 @@ vector< vector<float> > Chimera::readQuantiles() {
                        
                        in >> num >> ten >> twentyfive >> fifty >> seventyfive >> ninetyfive >> ninetynine; 
                        
-                       vector <float> temp;
+                       temp.clear();
                        
                        temp.push_back(ten); 
                        temp.push_back(twentyfive);
index b2e8be61f8df0aff52753058532a3ef42c6b6f0f..bfaae000386c872584952b27b5dcc64f52a785ea 100644 (file)
@@ -426,8 +426,7 @@ vector< vector<quanMember> > DeCalculator::getQuantiles(vector<Sequence*> seqs,
                                
                                quanMember newScore(de, i, j);
                                
-                               //dist-1 because vector indexes start at 0.
-                               quan[dist-1].push_back(newScore);
+                               quan[dist].push_back(newScore);
                                
                                delete subject;
                        }
index 04dfb0a76467ab6b030cfa60b8462dd7d30afd73..7b99cab9e7eed714a0cefeee4b45ed8057429c6f 100644 (file)
@@ -41,13 +41,15 @@ void Pintail::print(ostream& out) {
                for (int i = 0; i < querySeqs.size(); i++) {
                        
                        int index = ceil(deviation[i]);
-                       
-                       if (index == 0) { index=1;  }
                                                
                        //is your DE value higher than the 95%
                        string chimera;
-                       if (DE[i] > quantiles[index-1][4])              {       chimera = "Yes";        }
-                       else                                                                    {       chimera = "No";         }
+                       if (quantiles[index][4] == 0.0) {
+                               chimera = "Your template does not include sequences that provide quantile values at distance " + toString(index);
+                       }else {
+                               if (DE[i] > quantiles[index][4])                {       chimera = "Yes";        }
+                               else                                                                    {       chimera = "No";         }
+                       }
                        
                        out << querySeqs[i]->getName() << '\t' << "div: " << deviation[i] << "\tstDev: " << DE[i] << "\tchimera flag: " << chimera << endl;
                        if (chimera == "Yes") {
@@ -279,10 +281,19 @@ out7.close();/*/
                        ofstream out4, out5;
                        string noOutliers, outliers;
                        
-                       noOutliers = getRootName(templateFile) + "pintail.quanNOOUTLIERS";
-                       outliers = getRootName(templateFile) + "pintail.quanYESOUTLIERS";
+                       if ((!filter) && (seqMask == "")) {
+                               noOutliers = getRootName(templateFile) + "pintail.quan";
+                       }else if ((filter) && (seqMask == "")) { 
+                               noOutliers = getRootName(templateFile) + "pintail.filtered.quan";
+                       }else if ((!filter) && (seqMask != "")) { 
+                               noOutliers = getRootName(templateFile) + "pintail.masked.quan";
+                       }else if ((filter) && (seqMask != "")) { 
+                               noOutliers = getRootName(templateFile) + "pintail.filtered.masked.quan";
+                       }
+
+                       //outliers = getRootName(templateFile) + "pintail.quanYESOUTLIERS";
                        
-                       openOutputFile(outliers, out4);
+                       /*openOutputFile(outliers, out4);
                        
                        //adjust quantiles
                        for (int i = 0; i < quantilesMembers.size(); i++) {
@@ -321,7 +332,7 @@ out7.close();/*/
                                
                        }
                        
-                       out4.close();
+                       out4.close();*/
                        
                        decalc->removeObviousOutliers(quantilesMembers, templateSeqs.size());