sffinfo bug with flow grams right index when clipQualRight=0

[mothur.git] / raredisplay.cpp

diff --git a/raredisplay.cpp b/raredisplay.cpp
index b86312bc4411e8cdd2faca9934aab3abf5c96c52..b82127c51e6a64b50480a0f9b359cda3cd8e75f0 100644 (file)
--- a/raredisplay.cpp
+++ b/raredisplay.cpp
@@ -14,17 +14,9 @@
 void RareDisplay::init(string label){
        try {
                this->label = label;
-               if(nIters != 1){
-                       tempInFile.clear();
-                       openOutputFile(tempOutName, tempOutFile);
-                       openInputFile(tempInName, tempInFile);
-               }
-               else{
-                       openOutputFile(tempOutName, tempOutFile);
-               }
        }
        catch(exception& e) {
-               errorOut(e, "RareDisplay", "init");
+               m->errorOut(e, "RareDisplay", "init");
                exit(1);
        }
 }
@@ -38,24 +30,26 @@ void RareDisplay::update(SAbundVector* rank){
 
                if(nIters != 1){
 
-                       double oldNSeqs, oldMean, oldVar;
-               
-                       tempInFile >> oldNSeqs >> oldMean >> oldVar;
-               
-                       double oldS = oldVar * ( nIters - 2 );
-                       double delta = data[0] - oldMean;
-                       double newMean = oldMean + delta / nIters;
+                       double oldS = var[index] * ( nIters - 2 );
+                       double delta = data[0] - results[index];
+                       double newMean = results[index] + delta / nIters;
                        double newS = oldS + delta * ( data[0] - newMean );
                        double newVar = newS / ( nIters - 1 );
 
-                       tempOutFile << newNSeqs << '\t' << newMean << '\t' << newVar << endl;
+                       seqs[index] = newNSeqs;
+                       results[index] = newMean; 
+                       var[index] = newVar;
+                       
+                       index++;  
                }
                else{
-                       tempOutFile << newNSeqs << '\t' << data[0] << '\t' << 0 << endl;
+                       seqs.push_back(newNSeqs); 
+                       results.push_back(data[0]);
+                       var.push_back(0.0);
                }
        }
        catch(exception& e) {
-               errorOut(e, "RareDisplay", "update");
+               m->errorOut(e, "RareDisplay", "update");
                exit(1);
        }
 }
@@ -67,24 +61,28 @@ void RareDisplay::update(vector<SharedRAbundVector*> shared, int numSeqs, int nu
                double newNSeqs = data[0];
                
                if(nIters != 1){
-                       double oldNSeqs, oldMean, oldVar;
-               
-                       tempInFile >> oldNSeqs >> oldMean >> oldVar;
                
-                       double oldS = oldVar * ( nIters - 2 );
-                       double delta = data[0] - oldMean;
-                       double newMean = oldMean + delta / nIters;
+                       double oldS = var[index] * ( nIters - 2 );
+                       double delta = data[0] - results[index];
+                       double newMean = results[index] + delta / nIters;
                        double newS = oldS + delta * ( data[0] - newMean );
                        double newVar = newS / ( nIters - 1 );
-
-                       tempOutFile << newNSeqs << '\t' << newMean << '\t' << newVar << endl;
+                       seqs[index] = newNSeqs;
+                       results[index] = newMean; 
+                       var[index] = newVar;
+                       
+                       index++;  
                }
                else{
-                       tempOutFile << newNSeqs << '\t' << data[0] << '\t' << 0 << endl;
+                       
+                       seqs.push_back(newNSeqs); 
+                       results.push_back(data[0]);
+                       var.push_back(0.0);
+
                }
        }
        catch(exception& e) {
-               errorOut(e, "RareDisplay", "update");
+               m->errorOut(e, "RareDisplay", "update");
                exit(1);
        }
 }
@@ -93,24 +91,11 @@ void RareDisplay::update(vector<SharedRAbundVector*> shared, int numSeqs, int nu
 
 void RareDisplay::reset(){
        try {
-               if(nIters != 1){
-                       tempOutFile.close();
-                       tempInFile.close();
-               }
-               else{
-                       tempOutFile.close();
-               }
-       
                nIters++;
-       
-               remove(tempInName.c_str());
-               renameOk = rename(tempOutName.c_str(), tempInName.c_str());     
-               
-               //checks to make sure user was able to rename and remove successfully
-               if (renameOk != 0) { mothurOut("Unable to rename the necessary temp files."); mothurOutEndLine(); }
+               index = 0;
        }
        catch(exception& e) {
-               errorOut(e, "RareDisplay", "reset");
+               m->errorOut(e, "RareDisplay", "reset");
                exit(1);
        }
 }
@@ -122,39 +107,83 @@ void RareDisplay::close(){
                
                output->initFile(label);
        
-               openInputFile(tempInName, tempInFile);
-       
-               while(!tempInFile.eof()){
-                       int nSeqs;
-                       tempInFile >> nSeqs;
+               for (int i = 0; i < seqs.size(); i++) {
                
                        vector<double> data(3,0);
-                       double variance = 0;
-               
-                       tempInFile >> data[0];
-                       tempInFile >> variance;
-               
+                       double variance = var[i];
+                       
+                       data[0] = results[i];
+                       
                        double ci = 1.96 * pow(variance, 0.5);
                        data[1] = data[0] - ci;
                        data[2] = data[0] + ci;
                
-                       output->output(nSeqs, data);
-               
-                       gobble(tempInFile);
+                       output->output(seqs[i], data);
                }
-               tempInFile.close();
-       
-               remove(tempInName.c_str());
-               remove(tempOutName.c_str());
                
                nIters = 1;
+               index = 0;
+               
+               seqs.clear();
+               results.clear();
+               var.clear();
+               
                output->resetFile();
        }
        catch(exception& e) {
-               errorOut(e, "RareDisplay", "close");
+               m->errorOut(e, "RareDisplay", "close");
                exit(1);
        }
 }
+/***********************************************************************/
+
+void RareDisplay::inputTempFiles(string filename){
+       try {
+               ifstream in;
+               m->openInputFile(filename, in);
+               
+               int thisIters;
+               in >> thisIters; m->gobble(in);
+               
+               for (int i = 0; i < seqs.size(); i++) {
+                       double tempresult, tempvar;
+                       in >> tempresult >> tempvar; m->gobble(in);
+                       
+                       //find weighted result
+                       results[i] = ((nIters * results[i]) + (thisIters * tempresult)) / (float)(nIters + thisIters);
+                       
+                       var[i] = ((nIters * var[i]) + (thisIters * tempvar)) / (float)(nIters + thisIters);
+               }
+               
+               in.close();
+       }
+       catch(exception& e) {
+               m->errorOut(e, "RareDisplay", "inputTempFiles");
+               exit(1);
+       }
+}
+
+/***********************************************************************/
+
+void RareDisplay::outputTempFiles(string filename){
+       try {
+               ofstream out;
+               m->openOutputFile(filename, out);
+               
+               out << nIters << endl;
+               
+               for (int i = 0; i < seqs.size(); i++) {
+                       out << results[i] << '\t' << var[i] << endl;
+               }
+               
+               out.close();
+       }
+       catch(exception& e) {
+               m->errorOut(e, "RareDisplay", "outputTempFiles");
+               exit(1);
+       }
+}
+
 
 /***********************************************************************/