From: westcott Date: Thu, 19 Aug 2010 10:41:32 +0000 (+0000) Subject: modified calculators to use doubles, added numotu and fontsize parameters to heatmap... X-Git-Url: https://git.donarmstrong.com/?p=mothur.git;a=commitdiff_plain;h=e4c80376cc4533f66c8dfc18f3e1a86a60ac17fe modified calculators to use doubles, added numotu and fontsize parameters to heatmap as well as adding more sorting options. fixed bug in mothur.h setfilepositions that caused accession names to be truncated when windows line endings were used. --- diff --git a/Mothur.xcodeproj/project.pbxproj b/Mothur.xcodeproj/project.pbxproj index 919bc12..b081984 100644 --- a/Mothur.xcodeproj/project.pbxproj +++ b/Mothur.xcodeproj/project.pbxproj @@ -50,6 +50,8 @@ A7825503116519F70002E2DD /* chimerabellerophoncommand.cpp */ = {isa = PBXFileReference; fileEncoding = 4; lastKnownFileType = sourcecode.cpp.cpp; path = chimerabellerophoncommand.cpp; sourceTree = ""; }; A78434881162224F00100BE0 /* chimeraccodecommand.h */ = {isa = PBXFileReference; fileEncoding = 4; lastKnownFileType = sourcecode.c.h; path = chimeraccodecommand.h; sourceTree = ""; }; A78434891162224F00100BE0 /* chimeraccodecommand.cpp */ = {isa = PBXFileReference; fileEncoding = 4; lastKnownFileType = sourcecode.cpp.cpp; path = chimeraccodecommand.cpp; sourceTree = ""; }; + A7AE6302121C3408001DE6FC /* sharedrabundfloatvector.h */ = {isa = PBXFileReference; fileEncoding = 4; lastKnownFileType = sourcecode.c.h; path = sharedrabundfloatvector.h; sourceTree = ""; }; + A7AE6303121C3408001DE6FC /* sharedrabundfloatvector.cpp */ = {isa = PBXFileReference; fileEncoding = 4; lastKnownFileType = sourcecode.cpp.cpp; path = sharedrabundfloatvector.cpp; sourceTree = ""; }; A7BBDA7B11B5694E006E6551 /* classifyotucommand.h */ = {isa = PBXFileReference; fileEncoding = 4; lastKnownFileType = sourcecode.c.h; path = classifyotucommand.h; sourceTree = ""; }; A7BBDA7C11B5694E006E6551 /* classifyotucommand.cpp */ = {isa = PBXFileReference; fileEncoding = 4; lastKnownFileType = sourcecode.cpp.cpp; path = classifyotucommand.cpp; sourceTree = ""; }; A7D215C811996C6E00F13F13 /* clearcutcommand.h */ = {isa = PBXFileReference; fileEncoding = 4; lastKnownFileType = sourcecode.c.h; path = clearcutcommand.h; sourceTree = ""; }; @@ -870,8 +872,10 @@ A7DA212B113FECD400BF472F /* sharedordervector.h */, A7DA212C113FECD400BF472F /* sharedrabundvector.cpp */, A7DA212D113FECD400BF472F /* sharedrabundvector.h */, - A7DA212E113FECD400BF472F /* sharedsabundvector.cpp */, + A7AE6302121C3408001DE6FC /* sharedrabundfloatvector.h */, + A7AE6303121C3408001DE6FC /* sharedrabundfloatvector.cpp */, A7DA212F113FECD400BF472F /* sharedsabundvector.h */, + A7DA212E113FECD400BF472F /* sharedsabundvector.cpp */, A7DA214C113FECD400BF472F /* sparsematrix.cpp */, A7DA214D113FECD400BF472F /* sparsematrix.hpp */, A7DA214E113FECD400BF472F /* suffixdb.cpp */, diff --git a/ace.cpp b/ace.cpp index 1fde511..92c023f 100644 --- a/ace.cpp +++ b/ace.cpp @@ -16,12 +16,12 @@ EstOutput Ace::getValues(SAbundVector* rank) { data.resize(3,0); double ace, acelci, acehci; - int nrare = 0; - int srare = 0; - int sabund = 0; + double nrare = 0; + double srare = 0; + double sabund = 0; double Cace, term1, gamace; - int numsum = 0; + double numsum = 0; double maxRank = (double)rank->getMaxRank(); @@ -33,7 +33,7 @@ EstOutput Ace::getValues(SAbundVector* rank) { } else if(i>abund) {sabund += rank->get(i);} } - int sobs = srare + sabund; + double sobs = srare + sabund; if (nrare == 0){ Cace = 0.0000; } else { Cace = 1.0000 -(double)rank->get(1)/(double)nrare; } @@ -62,7 +62,8 @@ EstOutput Ace::getValues(SAbundVector* rank) { I have also added the forumlae to calculate the 95% confidence intervals. */ - int j,D_s=0,nn=0,ww=0,Max_Index=rank->getMaxRank()+1; + double j,D_s=0,nn=0,ww=0; + int Max_Index=rank->getMaxRank()+1; double pp, temp1, temp2; vector Part_N_Part_F(Max_Index+1,0.0); diff --git a/aligncommand.cpp b/aligncommand.cpp index f19810e..4641ed7 100644 --- a/aligncommand.cpp +++ b/aligncommand.cpp @@ -98,14 +98,6 @@ AlignCommand::AlignCommand(string option) { int ableToOpen; ifstream in; - - #ifdef USE_MPI - int pid; - MPI_Comm_size(MPI_COMM_WORLD, &processors); //set processors to the number of mpi processes running - MPI_Comm_rank(MPI_COMM_WORLD, &pid); //find out who we are - - if (pid == 0) { - #endif ableToOpen = openInputFile(candidateFileNames[i], in, "noerror"); @@ -119,16 +111,6 @@ AlignCommand::AlignCommand(string option) { } } in.close(); - #ifdef USE_MPI - for (int j = 1; j < processors; j++) { - MPI_Send(&ableToOpen, 1, MPI_INT, j, 2001, MPI_COMM_WORLD); - } - }else{ - MPI_Status status; - MPI_Recv(&ableToOpen, 1, MPI_INT, 0, 2001, MPI_COMM_WORLD, &status); - } - - #endif if (ableToOpen == 1) { m->mothurOut("Unable to open " + candidateFileNames[i] + ". It will be disregarded."); m->mothurOutEndLine(); @@ -263,7 +245,7 @@ int AlignCommand::execute(){ int tag = 2001; vector MPIPos; MPIWroteAccnos = false; - + MPI_Status status; MPI_Comm_rank(MPI_COMM_WORLD, &pid); //find out who we are MPI_Comm_size(MPI_COMM_WORLD, &processors); @@ -580,7 +562,7 @@ int AlignCommand::driver(linePair* line, string alignFName, string reportFName, //if there is a possibility that this sequence should be reversed if (candidateSeq->getNumBases() < numBasesNeeded) { - string wasBetter = ""; + string wasBetter = ""; //if the user wants you to try the reverse if (flip) { //get reverse compliment @@ -602,8 +584,9 @@ int AlignCommand::driver(linePair* line, string alignFName, string reportFName, delete nast; nast = nast2; needToDeleteCopy = true; + wasBetter = "\treverse complement produced a better alignment, so mothur used the reverse complement."; }else{ - wasBetter = "\treverse complement did NOT produce a better alignment, please check sequence."; + wasBetter = "\treverse complement did NOT produce a better alignment so it was not used, please check sequence."; delete nast2; delete copy; } @@ -678,7 +661,7 @@ int AlignCommand::driverMPI(int start, int num, MPI_File& inMPI, MPI_File& align int length = MPIPos[start+i+1] - MPIPos[start+i]; char* buf4 = new char[length]; - memcpy(buf4, outputString.c_str(), length); + //memcpy(buf4, outputString.c_str(), length); MPI_File_read_at(inMPI, MPIPos[start+i], buf4, length, MPI_CHAR, &status); @@ -687,9 +670,11 @@ int AlignCommand::driverMPI(int start, int num, MPI_File& inMPI, MPI_File& align delete buf4; if (tempBuf.length() > length) { tempBuf = tempBuf.substr(0, length); } - istringstream iss (tempBuf,istringstream::in); + istringstream iss (tempBuf,istringstream::in); + Sequence* candidateSeq = new Sequence(iss); + int origNumBases = candidateSeq->getNumBases(); string originalUnaligned = candidateSeq->getUnaligned(); int numBasesNeeded = origNumBases * threshold; diff --git a/alignment.cpp b/alignment.cpp index d543dbd..235216e 100644 --- a/alignment.cpp +++ b/alignment.cpp @@ -21,6 +21,7 @@ Alignment::Alignment() { /* do nothing */ } Alignment::Alignment(int A) : nCols(A), nRows(A) { try { + m = MothurOut::getInstance(); alignment.resize(nRows); // For the Gotoh and Needleman-Wunsch we initialize the dynamic programming for(int i=0;ireadKmerDB(kmerFileTest); } - + search->setNumSeqs(numSeqs); } } diff --git a/bootstrap.cpp b/bootstrap.cpp index 833f8f4..ce261ce 100644 --- a/bootstrap.cpp +++ b/bootstrap.cpp @@ -16,8 +16,8 @@ EstOutput Bootstrap::getValues(SAbundVector* rank){ //vector bootData(3,0); data.resize(1,0); double maxRank = (double)rank->getMaxRank(); - int sampled = rank->getNumSeqs(); - int sobs = rank->getNumBins(); + double sampled = rank->getNumSeqs(); + double sobs = rank->getNumBins(); double boot = (double)sobs; diff --git a/bstick.cpp b/bstick.cpp index d064ebe..9c36cd9 100644 --- a/bstick.cpp +++ b/bstick.cpp @@ -18,7 +18,7 @@ double BStick::invSum(int index, double numSpec) sum += 1/(double)i; return sum; } - +/***********************************************************************/ RAbundVector BStick::getRAbundVector(SAbundVector* rank){ vector rData; int mr = 1; @@ -26,7 +26,7 @@ RAbundVector BStick::getRAbundVector(SAbundVector* rank){ int ns = 0; for(int i = rank->size()-1; i > 0; i--) { - int cur = rank->get(i); + double cur = rank->get(i); if(mr == 1 && cur > 0) mr = i; nb += cur; diff --git a/chimerabellerophoncommand.cpp b/chimerabellerophoncommand.cpp index 9199909..6d09616 100644 --- a/chimerabellerophoncommand.cpp +++ b/chimerabellerophoncommand.cpp @@ -54,14 +54,6 @@ ChimeraBellerophonCommand::ChimeraBellerophonCommand(string option) { int ableToOpen; ifstream in; - - #ifdef USE_MPI - int pid; - MPI_Comm_size(MPI_COMM_WORLD, &processors); //set processors to the number of mpi processes running - MPI_Comm_rank(MPI_COMM_WORLD, &pid); //find out who we are - - if (pid == 0) { - #endif ableToOpen = openInputFile(fastaFileNames[i], in, "noerror"); @@ -75,17 +67,6 @@ ChimeraBellerophonCommand::ChimeraBellerophonCommand(string option) { } } in.close(); - - #ifdef USE_MPI - for (int j = 1; j < processors; j++) { - MPI_Send(&ableToOpen, 1, MPI_INT, j, 2001, MPI_COMM_WORLD); - } - }else{ - MPI_Status status; - MPI_Recv(&ableToOpen, 1, MPI_INT, 0, 2001, MPI_COMM_WORLD, &status); - } - - #endif if (ableToOpen == 1) { m->mothurOut("Unable to open " + fastaFileNames[i] + ". It will be disregarded."); m->mothurOutEndLine(); diff --git a/chimeraccodecommand.cpp b/chimeraccodecommand.cpp index 962c15b..4c8848d 100644 --- a/chimeraccodecommand.cpp +++ b/chimeraccodecommand.cpp @@ -66,14 +66,6 @@ ChimeraCcodeCommand::ChimeraCcodeCommand(string option) { int ableToOpen; ifstream in; - #ifdef USE_MPI - int pid; - MPI_Comm_size(MPI_COMM_WORLD, &processors); //set processors to the number of mpi processes running - MPI_Comm_rank(MPI_COMM_WORLD, &pid); //find out who we are - - if (pid == 0) { - #endif - ableToOpen = openInputFile(fastaFileNames[i], in, "noerror"); //if you can't open it, try default location @@ -87,17 +79,6 @@ ChimeraCcodeCommand::ChimeraCcodeCommand(string option) { } in.close(); - #ifdef USE_MPI - for (int j = 1; j < processors; j++) { - MPI_Send(&ableToOpen, 1, MPI_INT, j, 2001, MPI_COMM_WORLD); - } - }else{ - MPI_Status status; - MPI_Recv(&ableToOpen, 1, MPI_INT, 0, 2001, MPI_COMM_WORLD, &status); - } - - #endif - if (ableToOpen == 1) { m->mothurOut("Unable to open " + fastaFileNames[i] + ". It will be disregarded."); m->mothurOutEndLine(); //erase from file list diff --git a/chimeracheckcommand.cpp b/chimeracheckcommand.cpp index 264eb4c..6cea5d5 100644 --- a/chimeracheckcommand.cpp +++ b/chimeracheckcommand.cpp @@ -64,14 +64,6 @@ ChimeraCheckCommand::ChimeraCheckCommand(string option) { int ableToOpen; ifstream in; - #ifdef USE_MPI - int pid; - MPI_Comm_size(MPI_COMM_WORLD, &processors); //set processors to the number of mpi processes running - MPI_Comm_rank(MPI_COMM_WORLD, &pid); //find out who we are - - if (pid == 0) { - #endif - ableToOpen = openInputFile(fastaFileNames[i], in, "noerror"); //if you can't open it, try default location @@ -85,17 +77,6 @@ ChimeraCheckCommand::ChimeraCheckCommand(string option) { } in.close(); - #ifdef USE_MPI - for (int j = 1; j < processors; j++) { - MPI_Send(&ableToOpen, 1, MPI_INT, j, 2001, MPI_COMM_WORLD); - } - }else{ - MPI_Status status; - MPI_Recv(&ableToOpen, 1, MPI_INT, 0, 2001, MPI_COMM_WORLD, &status); - } - - #endif - if (ableToOpen == 1) { m->mothurOut("Unable to open " + fastaFileNames[i] +". It will be disregarded."); m->mothurOutEndLine(); //erase from file list @@ -131,14 +112,6 @@ ChimeraCheckCommand::ChimeraCheckCommand(string option) { int ableToOpen; ifstream in; - #ifdef USE_MPI - int pid; - MPI_Comm_size(MPI_COMM_WORLD, &processors); //set processors to the number of mpi processes running - MPI_Comm_rank(MPI_COMM_WORLD, &pid); //find out who we are - - if (pid == 0) { - #endif - ableToOpen = openInputFile(nameFileNames[i], in, "noerror"); //if you can't open it, try default location @@ -152,17 +125,6 @@ ChimeraCheckCommand::ChimeraCheckCommand(string option) { } in.close(); - #ifdef USE_MPI - for (int j = 1; j < processors; j++) { - MPI_Send(&ableToOpen, 1, MPI_INT, j, 2001, MPI_COMM_WORLD); - } - }else{ - MPI_Status status; - MPI_Recv(&ableToOpen, 1, MPI_INT, 0, 2001, MPI_COMM_WORLD, &status); - } - - #endif - if (ableToOpen == 1) { m->mothurOut("Unable to open " + nameFileNames[i] + ". It will be disregarded."); m->mothurOutEndLine(); //erase from file list diff --git a/chimerapintailcommand.cpp b/chimerapintailcommand.cpp index 5818ab2..fd8724d 100644 --- a/chimerapintailcommand.cpp +++ b/chimerapintailcommand.cpp @@ -83,14 +83,6 @@ ChimeraPintailCommand::ChimeraPintailCommand(string option) { int ableToOpen; ifstream in; - #ifdef USE_MPI - int pid; - MPI_Comm_size(MPI_COMM_WORLD, &processors); //set processors to the number of mpi processes running - MPI_Comm_rank(MPI_COMM_WORLD, &pid); //find out who we are - - if (pid == 0) { - #endif - ableToOpen = openInputFile(fastaFileNames[i], in, "noerror"); //if you can't open it, try default location @@ -104,17 +96,6 @@ ChimeraPintailCommand::ChimeraPintailCommand(string option) { } in.close(); - #ifdef USE_MPI - for (int j = 1; j < processors; j++) { - MPI_Send(&ableToOpen, 1, MPI_INT, j, 2001, MPI_COMM_WORLD); - } - }else{ - MPI_Status status; - MPI_Recv(&ableToOpen, 1, MPI_INT, 0, 2001, MPI_COMM_WORLD, &status); - } - - #endif - if (ableToOpen == 1) { m->mothurOut("Unable to open " + fastaFileNames[i] + ". It will be disregarded."); m->mothurOutEndLine(); //erase from file list diff --git a/chimeraslayercommand.cpp b/chimeraslayercommand.cpp index 1dd3d42..e675c00 100644 --- a/chimeraslayercommand.cpp +++ b/chimeraslayercommand.cpp @@ -69,14 +69,6 @@ ChimeraSlayerCommand::ChimeraSlayerCommand(string option) { int ableToOpen; ifstream in; - #ifdef USE_MPI - int pid; - MPI_Comm_size(MPI_COMM_WORLD, &processors); //set processors to the number of mpi processes running - MPI_Comm_rank(MPI_COMM_WORLD, &pid); //find out who we are - - if (pid == 0) { - #endif - ableToOpen = openInputFile(fastaFileNames[i], in, "noerror"); //if you can't open it, try default location @@ -90,17 +82,6 @@ ChimeraSlayerCommand::ChimeraSlayerCommand(string option) { } in.close(); - #ifdef USE_MPI - for (int j = 1; j < processors; j++) { - MPI_Send(&ableToOpen, 1, MPI_INT, j, 2001, MPI_COMM_WORLD); - } - }else{ - MPI_Status status; - MPI_Recv(&ableToOpen, 1, MPI_INT, 0, 2001, MPI_COMM_WORLD, &status); - } - - #endif - if (ableToOpen == 1) { m->mothurOut("Unable to open " + fastaFileNames[i] + ". It will be disregarded."); m->mothurOutEndLine(); //erase from file list diff --git a/classifyseqscommand.cpp b/classifyseqscommand.cpp index 6c42418..128fd01 100644 --- a/classifyseqscommand.cpp +++ b/classifyseqscommand.cpp @@ -98,14 +98,6 @@ ClassifySeqsCommand::ClassifySeqsCommand(string option) { int ableToOpen; - #ifdef USE_MPI - int pid; - MPI_Comm_size(MPI_COMM_WORLD, &processors); //set processors to the number of mpi processes running - MPI_Comm_rank(MPI_COMM_WORLD, &pid); //find out who we are - - if (pid == 0) { - #endif - ifstream in; ableToOpen = openInputFile(fastaFileNames[i], in, "noerror"); @@ -120,17 +112,6 @@ ClassifySeqsCommand::ClassifySeqsCommand(string option) { } in.close(); - #ifdef USE_MPI - for (int j = 1; j < processors; j++) { - MPI_Send(&ableToOpen, 1, MPI_INT, j, 2001, MPI_COMM_WORLD); - } - }else{ - MPI_Status status; - MPI_Recv(&ableToOpen, 1, MPI_INT, 0, 2001, MPI_COMM_WORLD, &status); - } - - #endif - if (ableToOpen == 1) { m->mothurOut("Unable to open " + fastaFileNames[i] + ". It will be disregarded."); m->mothurOutEndLine(); //erase from file list @@ -169,14 +150,6 @@ ClassifySeqsCommand::ClassifySeqsCommand(string option) { } int ableToOpen; - #ifdef USE_MPI - int pid; - MPI_Comm_size(MPI_COMM_WORLD, &processors); //set processors to the number of mpi processes running - MPI_Comm_rank(MPI_COMM_WORLD, &pid); //find out who we are - - if (pid == 0) { - #endif - ifstream in; ableToOpen = openInputFile(namefileNames[i], in, "noerror"); @@ -191,17 +164,6 @@ ClassifySeqsCommand::ClassifySeqsCommand(string option) { } in.close(); - #ifdef USE_MPI - for (int j = 1; j < processors; j++) { - MPI_Send(&ableToOpen, 1, MPI_INT, j, 2001, MPI_COMM_WORLD); - } - }else{ - MPI_Status status; - MPI_Recv(&ableToOpen, 1, MPI_INT, 0, 2001, MPI_COMM_WORLD, &status); - } - - #endif - if (ableToOpen == 1) { m->mothurOut("Unable to open " + namefileNames[i] + ". It will be disregarded."); m->mothurOutEndLine(); abort = true; //erase from file list @@ -230,14 +192,6 @@ ClassifySeqsCommand::ClassifySeqsCommand(string option) { } int ableToOpen; - #ifdef USE_MPI - int pid; - MPI_Comm_size(MPI_COMM_WORLD, &processors); //set processors to the number of mpi processes running - MPI_Comm_rank(MPI_COMM_WORLD, &pid); //find out who we are - - if (pid == 0) { - #endif - ifstream in; ableToOpen = openInputFile(groupfileNames[i], in, "noerror"); @@ -252,17 +206,6 @@ ClassifySeqsCommand::ClassifySeqsCommand(string option) { } in.close(); - #ifdef USE_MPI - for (int j = 1; j < processors; j++) { - MPI_Send(&ableToOpen, 1, MPI_INT, j, 2001, MPI_COMM_WORLD); - } - }else{ - MPI_Status status; - MPI_Recv(&ableToOpen, 1, MPI_INT, 0, 2001, MPI_COMM_WORLD, &status); - } - - #endif - if (ableToOpen == 1) { m->mothurOut("Unable to open " + groupfileNames[i] + ". It will be disregarded."); m->mothurOutEndLine(); groupfileNames[i] = ""; //erase from file list diff --git a/clearcutcommand.cpp b/clearcutcommand.cpp index 7bea292..b1f40e0 100644 --- a/clearcutcommand.cpp +++ b/clearcutcommand.cpp @@ -92,7 +92,7 @@ ClearcutCommand::ClearcutCommand(string option) { temp = validParameter.validFile(parameters, "shuffle", false); if (temp == "not found"){ temp = "F"; } shuffle = isTrue(temp); - temp = validParameter.validFile(parameters, "neighbor", false); if (temp == "not found"){ temp = "F"; } + temp = validParameter.validFile(parameters, "neighbor", false); if (temp == "not found"){ temp = "T"; } neighbor = isTrue(temp); temp = validParameter.validFile(parameters, "DNA", false); if (temp == "not found"){ temp = "F"; } @@ -146,7 +146,7 @@ void ClearcutCommand::help(){ m->mothurOut("The seed parameter allows you to explicitly set the PRNG seed to a specific value. \n"); m->mothurOut("The norandom parameter allows you to attempt joins deterministically, default=F. \n"); m->mothurOut("The shuffle parameter allows you to randomly shuffle the distance matrix, default=F. \n"); - m->mothurOut("The neighbor parameter allows you to use traditional Neighbor-Joining algorithm, default=F. \n"); + m->mothurOut("The neighbor parameter allows you to use traditional Neighbor-Joining algorithm, default=T. \n"); m->mothurOut("The DNA parameter allows you to indicate your fasta file contains DNA sequences, default=F. \n"); m->mothurOut("The protein parameter allows you to indicate your fasta file contains protein sequences, default=F. \n"); diff --git a/clustercommand.cpp b/clustercommand.cpp index e9df82f..3f27560 100644 --- a/clustercommand.cpp +++ b/clustercommand.cpp @@ -232,7 +232,12 @@ int ClusterCommand::execute(){ rabundFile.close(); listFile.close(); - if (saveCutoff != cutoff) { m->mothurOut("changed cutoff to " + toString(cutoff)); m->mothurOutEndLine(); } + if (saveCutoff != cutoff) { + if (hard) { saveCutoff = ceilDist(saveCutoff, precision); } + else { saveCutoff = roundDist(saveCutoff, precision); } + + m->mothurOut("changed cutoff to " + toString(cutoff)); m->mothurOutEndLine(); + } m->mothurOutEndLine(); m->mothurOut("Output File Names: "); m->mothurOutEndLine(); diff --git a/clustersplitcommand.cpp b/clustersplitcommand.cpp index d6855cc..ced5a63 100644 --- a/clustersplitcommand.cpp +++ b/clustersplitcommand.cpp @@ -994,7 +994,12 @@ vector ClusterSplitCommand::cluster(vector< map > distNa remove(thisDistFile.c_str()); remove(thisNamefile.c_str()); - if (saveCutoff != cutoff) { m->mothurOut("Cutoff was " + toString(cutoff) + " changed cutoff to " + toString(saveCutoff)); m->mothurOutEndLine(); } + if (saveCutoff != cutoff) { + if (hard) { saveCutoff = ceilDist(saveCutoff, precision); } + else { saveCutoff = roundDist(saveCutoff, precision); } + + m->mothurOut("Cutoff was " + toString(cutoff) + " changed cutoff to " + toString(saveCutoff)); m->mothurOutEndLine(); + } if (saveCutoff < smallestCutoff) { smallestCutoff = saveCutoff; } } diff --git a/filterseqscommand.cpp b/filterseqscommand.cpp index 7ab00a7..e6e616a 100644 --- a/filterseqscommand.cpp +++ b/filterseqscommand.cpp @@ -263,7 +263,7 @@ int FilterSeqsCommand::filterSequences() { MPI_Status status; MPI_Comm_size(MPI_COMM_WORLD, &processors); //set processors to the number of mpi processes running MPI_Comm_rank(MPI_COMM_WORLD, &pid); //find out who we are - cout << pid << "is in create filter " << endl; + MPI_File outMPI; MPI_File tempMPI; MPI_File inMPI; @@ -494,7 +494,14 @@ int FilterSeqsCommand::driverRunFilter(string F, string outputFilename, string i out << '>' << seq.getName() << endl << filterSeq << endl; } gobble(in); + + //report progress + if((i+1) % 100 == 0){ m->mothurOut(toString(i+1)); m->mothurOutEndLine(); } } + + //report progress + if((line->num) % 100 != 0){ m->mothurOut(toString(line->num)); m->mothurOutEndLine(); } + out.close(); in.close(); @@ -853,13 +860,27 @@ int FilterSeqsCommand::createProcessesCreateFilter(Filters& F, string filename) //loop through and create all the processes you want while (process != processors) { - int pid = vfork(); + 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){ driverCreateFilter(F, filename, lines[process]); + + //write out filter counts to file + filename += toString(getpid()) + "filterValues.temp"; + ofstream out; + openOutputFile(filename, out); + + for (int k = 0; k < alignmentLength; k++) { out << F.a[k] << '\t'; } out << endl; + for (int k = 0; k < alignmentLength; k++) { out << F.t[k] << '\t'; } out << endl; + for (int k = 0; k < alignmentLength; k++) { out << F.g[k] << '\t'; } out << endl; + for (int k = 0; k < alignmentLength; k++) { out << F.c[k] << '\t'; } out << endl; + for (int k = 0; k < alignmentLength; k++) { out << F.gap[k] << '\t'; } out << endl; + + out.close(); + exit(0); }else { m->mothurOut("unable to spawn the necessary processes."); m->mothurOutEndLine(); exit(0); } } @@ -870,6 +891,23 @@ int FilterSeqsCommand::createProcessesCreateFilter(Filters& F, string filename) wait(&temp); } + //parent reads in and combine Filter info + for (int i = 0; i < processIDS.size(); i++) { + string tempFilename = filename + toString(processIDS[i]) + "filterValues.temp"; + ifstream in; + openInputFile(tempFilename, in); + + int temp; + for (int k = 0; k < alignmentLength; k++) { in >> temp; F.a[k] += temp; } gobble(in); + for (int k = 0; k < alignmentLength; k++) { in >> temp; F.t[k] += temp; } gobble(in); + for (int k = 0; k < alignmentLength; k++) { in >> temp; F.g[k] += temp; } gobble(in); + for (int k = 0; k < alignmentLength; k++) { in >> temp; F.c[k] += temp; } gobble(in); + for (int k = 0; k < alignmentLength; k++) { in >> temp; F.gap[k] += temp; } gobble(in); + + in.close(); + remove(tempFilename.c_str()); + } + return exitCommand; #endif } diff --git a/geom.cpp b/geom.cpp index 201f8f5..5bd150e 100644 --- a/geom.cpp +++ b/geom.cpp @@ -44,9 +44,9 @@ EstOutput Geom::getValues(SAbundVector* rank){ data.resize(3,0); rdata = getRAbundVector(rank); - int numInd = rdata.getNumSeqs(); - int numSpec = rdata.getNumBins(); - int min = rdata.get(rdata.size()-1); + double numInd = rdata.getNumSeqs(); + double numSpec = rdata.getNumBins(); + double min = rdata.get(rdata.size()-1); double k = .5; double step = .49999; diff --git a/heatmap.cpp b/heatmap.cpp index 3dca535..aef9f37 100644 --- a/heatmap.cpp +++ b/heatmap.cpp @@ -10,7 +10,7 @@ #include "heatmap.h" //********************************************************************************************************************** -HeatMap::HeatMap(string sort, string scale, string dir){ +HeatMap::HeatMap(string sort, string scale, int num, int fsize, string dir){ try { globaldata = GlobalData::getInstance(); m = MothurOut::getInstance(); @@ -18,6 +18,8 @@ HeatMap::HeatMap(string sort, string scale, string dir){ sorted = sort; scaler = scale; outputDir = dir; + numOTU = num; + fontSize = fsize; } catch(exception& e) { m->errorOut(e, "HeatMap", "HeatMap"); @@ -70,7 +72,7 @@ string HeatMap::getPic(RAbundVector* rabund) { //white backround outsvg << "getNumBins()*5 + 120)) + "\"/>"; - outsvg << "Heatmap at distance " + rabund->getLabel() + "\n"; + outsvg << "Heatmap at distance " + rabund->getLabel() + "\n"; //output legend and color labels string color; @@ -101,8 +103,15 @@ string HeatMap::getPic(RAbundVector* rabund) { string HeatMap::getPic(vector lookup) { try { + + int numBinsToDisplay = lookup[0]->size(); + + if (numOTU != 0) { //user want to display a portion of the otus + if (numOTU < numBinsToDisplay) { numBinsToDisplay = numOTU; } + } + //sort lookup so shared bins are on top - if (isTrue(sorted) == true) { sortSharedVectors(lookup); } + if (sorted != "none") { sortSharedVectors(lookup); } vector > scaleRelAbund; vector maxRelAbund(lookup.size(), 0.0); @@ -148,13 +157,13 @@ string HeatMap::getPic(vector lookup) { //white backround outsvg << "getNumBins()*5 + 120)) + "\"/>"; - outsvg << "Heatmap at distance " + lookup[0]->getLabel() + "\n"; + outsvg << "Heatmap at distance " + lookup[0]->getLabel() + "\n"; //column labels for (int h = 0; h < lookup.size(); h++) { - outsvg << "getGroup().length() / 2)) + "\" y=\"50\">" + lookup[h]->getGroup() + "\n"; + outsvg << "getGroup().length() / 2)) + "\" y=\"50\">" + lookup[h]->getGroup() + "\n"; } - + //output legend and color labels string color; int x = 0; @@ -162,7 +171,7 @@ string HeatMap::getPic(vector lookup) { printLegend(y, superMaxRelAbund); y = 70; - for (int i = 0; i < scaleRelAbund[0].size(); i++) { + for (int i = 0; i < numBinsToDisplay; i++) { for (int j = 0; j < scaleRelAbund.size(); j++) { if (m->control_pressed) { outsvg.close(); return "control"; } @@ -186,18 +195,24 @@ string HeatMap::getPic(vector lookup) { } //********************************************************************************************************************** -void HeatMap::sortSharedVectors(vector& lookup){ +int HeatMap::sortSharedVectors(vector& lookup){ try { - //copy lookup and then clear it to refill with sorted. - //loop though lookup and determine if they are shared - //if they are then insert in the front - //if not push to back - + vector looktemp; map place; //spot in lookup where you insert shared by, ie, 3 -> 2 if they are shared by 3 inset into location 2. map::iterator it; - int count; + /****************** find order of otus **********************/ + if (sorted == "shared") { + place = orderShared(lookup); + }else if (sorted == "topotu") { + place = orderTopOtu(lookup); + }else if (sorted == "topgroup") { + place = orderTopGroup(lookup); + }else { m->mothurOut("Error: invalid sort option."); m->mothurOutEndLine(); return 1; } + + + /******************* create copy of lookup *********************/ //create and initialize looktemp as a copy of lookup for (int i = 0; i < lookup.size(); i++) { SharedRAbundVector* temp = new SharedRAbundVector(lookup[i]->getNumBins()); @@ -209,62 +224,139 @@ void HeatMap::sortSharedVectors(vector& lookup){ } looktemp.push_back(temp); } + + /************************ fill lookup in order given by place *********************/ + //for each bin + for (int i = 0; i < looktemp[0]->size(); i++) { //place + //fill lookup // 2 -> 1 + for (int j = 0; j < looktemp.size(); j++) { // 3 -> 2 + int newAbund = looktemp[j]->getAbundance(i); // 1 -> 3 + lookup[j]->set(place[i], newAbund, looktemp[j]->getGroup()); //binNumber, abundance, group + } + } - //clear out lookup to create sorted lookup -- Sarah look at - this is causing segmentation faults - for (int j = 0; j < lookup.size(); j++) { -// delete lookup[j]; + //delete looktemp -- Sarah look at - this is causing segmentation faults + for (int j = 0; j < looktemp.size(); j++) { +// delete looktemp[j]; } - lookup.clear(); //doesn't this do the job? - - //create and initialize lookup to empty vectors - for (int i = 0; i < looktemp.size(); i++) { - SharedRAbundVector* temp = new SharedRAbundVector(); - temp->setLabel(looktemp[i]->getLabel()); - temp->setGroup(looktemp[i]->getGroup()); - lookup.push_back(temp); + + return 0; + + } + catch(exception& e) { + m->errorOut(e, "HeatMap", "sortSharedVectors"); + exit(1); + } +} +//********************************************************************************************************************** +map HeatMap::orderShared(vector& lookup){ + try { + + map place; //spot in lookup where you insert shared by, ie, 3 -> 2 if they are shared by 3 inset into location 2. + map::iterator it; + + vector sharedBins; + vector uniqueBins; + + //for each bin + for (int i = 0; i < lookup[0]->size(); i++) { + int count = 0; + + //is this bin shared + for (int j = 0; j < lookup.size(); j++) { if (lookup[j]->getAbundance(i) != 0) { count++; } } - //initialize place map - place[i] = 0; + if (count < 2) { uniqueBins.push_back(i); } + else { sharedBins.push_back(i); } } + //fill place + for (int i = 0; i < sharedBins.size(); i++) { place[sharedBins[i]] = i; } + for (int i = 0; i < uniqueBins.size(); i++) { place[uniqueBins[i]] = (sharedBins.size() + i); } + + return place; + + } + catch(exception& e) { + m->errorOut(e, "HeatMap", "orderShared"); + exit(1); + } +} +//********************************************************************************************************************** +map HeatMap::orderTopOtu(vector& lookup){ + try { + + map place; //spot in lookup where you insert shared by, ie, 3 -> 2 if they are shared by 3 inset into location 2. + map::iterator it; + + vector totals; //for each bin - for (int i = 0; i < looktemp[0]->size(); i++) { - count = 0; - bool updatePlace = false; - //for each group - for (int j = 0; j < looktemp.size(); j++) { - if (looktemp[j]->getAbundance(i) != 0) { count++; } - } + for (int i = 0; i < lookup[0]->size(); i++) { + int total = 0; - //fill lookup - for (int j = 0; j < looktemp.size(); j++) { - //if they are not shared then push to back, if they are not insert in front - if (count < 2) { lookup[j]->push_back(looktemp[j]->getAbundance(i), looktemp[j]->getGroup()); } - //they are shared by some - else { lookup[j]->insert(looktemp[j]->getAbundance(i), place[count], looktemp[j]->getGroup()); updatePlace = true; } - } + for (int j = 0; j < lookup.size(); j++) { total += lookup[j]->getAbundance(i); } - if (updatePlace == true) { - //move place holders below where you entered up to "make space" for you entry - for (it = place.begin(); it!= place.end(); it++) { - if (it->first < count) { it->second++; } - } + binCount temp(i, total); + + totals.push_back(temp); + } + + sort(totals.begin(), totals.end(), comparebinCounts); + + //fill place + for (int i = 0; i < totals.size(); i++) { place[totals[i].bin] = i; } + + return place; + + } + catch(exception& e) { + m->errorOut(e, "HeatMap", "orderTopOtu"); + exit(1); + } +} +//********************************************************************************************************************** +map HeatMap::orderTopGroup(vector& lookup){ + try { + + map place; //spot in lookup where you insert shared by, ie, 3 -> 2 if they are shared by 3 inset into location 2. + map::iterator it; + + vector < vector > totals; //totals[0] = bin totals for group 0, totals[1] = bin totals for group 1, ... + totals.resize(lookup.size()); + + //for each bin + for (int i = 0; i < lookup[0]->size(); i++) { + for (int j = 0; j < lookup.size(); j++) { + binCount temp(i, (lookup[j]->getAbundance(i))); + totals[j].push_back(temp); } } - //delete looktemp -- Sarah look at - this is causing segmentation faults - for (int j = 0; j < looktemp.size(); j++) { -// delete looktemp[j]; + for (int i = 0; i < totals.size(); i++) { sort(totals[i].begin(), totals[i].end(), comparebinCounts); } + + //fill place + //grab the top otu for each group adding it if its not already added + int count = 0; + for (int i = 0; i < totals[0].size(); i++) { + + for (int j = 0; j < totals.size(); j++) { + it = place.find(totals[j][i].bin); + + if (it == place.end()) { //not added yet + place[totals[j][i].bin] = count; + count++; + } + } } + + return place; } catch(exception& e) { - m->errorOut(e, "HeatMap", "sortSharedVectors"); + m->errorOut(e, "HeatMap", "orderTopGroup"); exit(1); } } - //********************************************************************************************************************** void HeatMap::printLegend(int y, float maxbin) { @@ -304,6 +396,8 @@ void HeatMap::printLegend(int y, float maxbin) { exit(1); } } +//********************************************************************************************************************** + diff --git a/heatmap.h b/heatmap.h index ffbe946..11c67ea 100644 --- a/heatmap.h +++ b/heatmap.h @@ -14,25 +14,41 @@ #include "datavector.hpp" #include "globaldata.hpp" +/***********************************************************************/ +struct binCount { + int bin; + int abund; + binCount(int i, int j) : bin(i), abund(j) {} +}; +/***********************************************************************/ +//sorts highest abund to lowest +inline bool comparebinCounts(binCount left, binCount right){ + return (left.abund > right.abund); +} /***********************************************************************/ class HeatMap { public: - HeatMap(string, string, string); + HeatMap(string, string, int, int, string); ~HeatMap(){}; string getPic(RAbundVector*); string getPic(vector); private: - void sortSharedVectors(vector& ); + int sortSharedVectors(vector& ); void printLegend(int, float); GlobalData* globaldata; string format, sorted, groupComb, scaler, outputDir; ofstream outsvg; MothurOut* m; + int numOTU, fontSize; + + map orderTopGroup(vector&); + map orderTopOtu(vector&); + map orderShared(vector&); }; diff --git a/heatmapcommand.cpp b/heatmapcommand.cpp index 7aa978d..b57f3a0 100644 --- a/heatmapcommand.cpp +++ b/heatmapcommand.cpp @@ -24,7 +24,7 @@ HeatMapCommand::HeatMapCommand(string option) { else { //valid paramters for this command - string AlignArray[] = {"groups","label","sorted","scale","outputdir","inputdir"}; + string AlignArray[] = {"groups","label","sorted","scale","fontsize","numotu","outputdir","inputdir"}; vector myArray (AlignArray, AlignArray+(sizeof(AlignArray)/sizeof(string))); OptionParser parser(option); @@ -70,10 +70,22 @@ HeatMapCommand::HeatMapCommand(string option) { globaldata->Groups = Groups; } - sorted = validParameter.validFile(parameters, "sorted", false); if (sorted == "not found") { sorted = "T"; } + string temp = validParameter.validFile(parameters, "numotu", false); if (temp == "not found") { temp = "0"; } + convert(temp, numOTU); + + temp = validParameter.validFile(parameters, "fontsize", false); if (temp == "not found") { temp = "24"; } + convert(temp, fontSize); + + sorted = validParameter.validFile(parameters, "sorted", false); + if (sorted == "not found") { + //if numOTU is used change default + if (numOTU != 0) { sorted = "topotu"; } + else { sorted = "shared"; } + } scale = validParameter.validFile(parameters, "scale", false); if (scale == "not found") { scale = "log10"; } + if ((sorted != "none") && (sorted != "shared") && (sorted != "topotu") && (sorted != "topgroup")) { m->mothurOut(sorted + " is not a valid sorting option. Sorted options are: none, shared, topotu, topgroup"); m->mothurOutEndLine(); abort=true; } } } @@ -88,15 +100,17 @@ HeatMapCommand::HeatMapCommand(string option) { void HeatMapCommand::help(){ try { m->mothurOut("The heatmap.bin command can only be executed after a successful read.otu command.\n"); - m->mothurOut("The heatmap.bin command parameters are groups, sorted, scale label. No parameters are required.\n"); + m->mothurOut("The heatmap.bin command parameters are groups, sorted, scale, numotu, fontsize and label. No parameters are required.\n"); m->mothurOut("The groups parameter allows you to specify which of the groups in your groupfile you would like included in your heatmap.\n"); - m->mothurOut("The sorted parameter allows you to choose to see the file with the shared otus at the top or the otus in the order they appear in your input file. \n"); + m->mothurOut("The sorted parameter allows you to order the otus displayed, default=shared, meaning display the shared otus first. Other options for sorted are none, meaning the exact representation of your otus, \n"); + m->mothurOut("topotu, meaning the otus with the greatest abundance when totaled across groups, topgroup, meaning the top otus for each group. \n"); m->mothurOut("The scale parameter allows you to choose the range of color your bin information will be displayed with.\n"); + m->mothurOut("The numotu parameter allows you to display only the top N otus, by default all the otus are displayed. You could choose to look at the top 10, by setting numotu=10. The default for sorted is topotu when numotu is used.\n"); m->mothurOut("The group names are separated by dashes. The label parameter allows you to select what distance levels you would like a heatmap created for, and are also separated by dashes.\n"); + m->mothurOut("The fontsize parameter allows you to adjust the font size of the picture created, default=24.\n"); m->mothurOut("The heatmap.bin command should be in the following format: heatmap.bin(groups=yourGroups, sorted=yourSorted, label=yourLabels).\n"); m->mothurOut("Example heatmap.bin(groups=A-B-C, sorted=F, scale=log10).\n"); m->mothurOut("The default value for groups is all the groups in your groupfile, and all labels in your inputfile will be used.\n"); - m->mothurOut("The default value for sorted is T meaning you want the shared otus on top, you may change it to F meaning the exact representation of your input file.\n"); m->mothurOut("The default value for scale is log10; your other options are log2 and linear.\n"); m->mothurOut("The heatmap.bin command outputs a .svg file for each label you specify.\n"); m->mothurOut("Note: No spaces between parameter labels (i.e. groups), '=' and parameters (i.e.yourGroups).\n\n"); @@ -120,7 +134,7 @@ int HeatMapCommand::execute(){ if (abort == true) { return 0; } - heatmap = new HeatMap(sorted, scale, outputDir); + heatmap = new HeatMap(sorted, scale, numOTU, fontSize, outputDir); format = globaldata->getFormat(); string lastLabel; diff --git a/heatmapcommand.h b/heatmapcommand.h index 47838d8..f711cc3 100644 --- a/heatmapcommand.h +++ b/heatmapcommand.h @@ -41,6 +41,7 @@ private: set labels; //holds labels to be used string format, groups, sorted, scale, label, outputDir; vector Groups, outputNames; + int numOTU, fontSize; }; diff --git a/jackknife.cpp b/jackknife.cpp index 7f61366..150b7cd 100644 --- a/jackknife.cpp +++ b/jackknife.cpp @@ -67,7 +67,7 @@ EstOutput Jackknife::getValues(SAbundVector* rank){ double jack, jacklci, jackhci; - double maxRank = (double)rank->getMaxRank(); + int maxRank = (double)rank->getMaxRank(); int S = rank->getNumBins(); double N[maxOrder+1]; diff --git a/logsd.cpp b/logsd.cpp index 009a404..6a79541 100644 --- a/logsd.cpp +++ b/logsd.cpp @@ -14,6 +14,7 @@ double LogSD::logS(double x){ return -(1-x)*log(1-x)/x; } +/***********************************************************************/ EstOutput LogSD::getValues(SAbundVector* rank){ try { @@ -29,8 +30,8 @@ EstOutput LogSD::getValues(SAbundVector* rank){ SAbundVector *rank = &rankw;*/ data.resize(3,0); - int numInd = rank->getNumSeqs(); - int numSpec = rank->getNumBins(); + double numInd = rank->getNumSeqs(); + double numSpec = rank->getNumBins(); double snRatio = (double)numSpec/numInd; double x = .5; double step = .4999999999; @@ -44,7 +45,7 @@ EstOutput LogSD::getValues(SAbundVector* rank){ } double alpha = numInd*(1-x)/x; - int oct = 1; + double oct = 1; double octSumObs = 0; double sumObs = 0; double octSumExp = 0; diff --git a/mothur.h b/mothur.h index 0c4429e..765c90a 100644 --- a/mothur.h +++ b/mothur.h @@ -222,6 +222,16 @@ inline void gobble(istream& f){ } /***********************************************************************/ +inline void gobble(istringstream& f){ + + char d; + while(isspace(d=f.get())) {;} + f.putback(d); + +} + +/***********************************************************************/ + inline string getline(istringstream& fileHandle) { try { @@ -1044,13 +1054,14 @@ inline vector setFilePosFasta(string filename, int& num) { vector positions; ifstream inFASTA; openInputFile(filename, inFASTA); - + string input; while(!inFASTA.eof()){ - input = getline(inFASTA); gobble(inFASTA); + input = getline(inFASTA); if (input.length() != 0) { if(input[0] == '>'){ unsigned long int pos = inFASTA.tellg(); positions.push_back(pos - input.length() - 1); } } + gobble(inFASTA); //has to be here since windows line endings are 2 characters and mess up the positions } inFASTA.close(); @@ -1094,12 +1105,13 @@ inline vector setFilePosEachLine(string filename, int& num) { string input; while(!in.eof()){ unsigned long int lastpos = in.tellg(); - input = getline(in); gobble(in); + input = getline(in); if (input.length() != 0) { unsigned long int pos = in.tellg(); if (pos != -1) { positions.push_back(pos - input.length() - 1); } else { positions.push_back(lastpos); } } + gobble(in); //has to be here since windows line endings are 2 characters and mess up the positions } in.close(); diff --git a/npshannon.cpp b/npshannon.cpp index f14f922..a855b31 100644 --- a/npshannon.cpp +++ b/npshannon.cpp @@ -18,7 +18,7 @@ EstOutput NPShannon::getValues(SAbundVector* rank){ float npShannon = 0.0000; double maxRank = (double)rank->getMaxRank(); - int sampled = rank->getNumSeqs(); + double sampled = rank->getNumSeqs(); double Chat = 1.0000 - (double)rank->get(1)/(double)sampled; diff --git a/qstat.cpp b/qstat.cpp index 44da363..e07cdb3 100644 --- a/qstat.cpp +++ b/qstat.cpp @@ -31,8 +31,8 @@ EstOutput QStat::getValues(SAbundVector* rank){ int r3 = -1; int r1Ind = 0; int r3Ind = 0; - int sumSpec = 0; - int iqSum = 0; + double sumSpec = 0; + double iqSum = 0; for(int i = 1; i < rank->size(); i++) { if(r1 != -1 && r3 != -1) i = rank->size(); diff --git a/sharedace.cpp b/sharedace.cpp index 24fe5ef..c8ba4fa 100644 --- a/sharedace.cpp +++ b/sharedace.cpp @@ -18,10 +18,10 @@ EstOutput SharedAce::getValues(vector shared) { string label; label = shared[0]->getLabel(); - int fARare, fBRare, S12Rare, S12Abund, S12, f11, tempA, tempB, t10, t01, t11, t21, t12, t22, C12Numerator; + double fARare, fBRare, S12Rare, S12Abund, S12, f11, tempA, tempB, t10, t01, t11, t21, t12, t22, C12Numerator; fARare = 0; fBRare = 0; S12Rare = 0; S12Abund = 0; S12 = 0; f11 = 0; t10 = 0; t01 = 0; t11= 0; t21= 0; t12= 0; t22= 0; C12Numerator = 0; - float Sharedace, C12, part1, part2, part3, part4, part5, Gamma1, Gamma2, Gamma3; + double Sharedace, C12, part1, part2, part3, part4, part5, Gamma1, Gamma2, Gamma3; /*fARare = number of OTUs with one individual found in A and less than or equal to 10 in B. fBRare = number of OTUs with one individual found in B and less than or equal to 10 in A. @@ -89,7 +89,7 @@ EstOutput SharedAce::getValues(vector shared) { if (isnan(Gamma3) || isinf(Gamma3)) { Gamma3 = 0; } if (isnan(part1) || isinf(part1)) { part1 = 0; } if (isnan(part2) || isinf(part2)) { part2 = 0; } - + part3 = fARare * Gamma1; part4 = fBRare * Gamma2; part5 = f11 * Gamma3; diff --git a/sharedanderbergs.cpp b/sharedanderbergs.cpp index 38eec63..09ec692 100644 --- a/sharedanderbergs.cpp +++ b/sharedanderbergs.cpp @@ -13,7 +13,7 @@ EstOutput Anderberg::getValues(vector shared) { try { - int S1, S2, S12, tempA, tempB; + double S1, S2, S12, tempA, tempB; S1 = 0; S2 = 0; S12 = 0; tempA = 0; tempB = 0; /*S1, S2 = number of OTUs observed or estimated in A and B diff --git a/sharedbraycurtis.cpp b/sharedbraycurtis.cpp index e8512d4..566df3b 100644 --- a/sharedbraycurtis.cpp +++ b/sharedbraycurtis.cpp @@ -15,7 +15,7 @@ EstOutput BrayCurtis::getValues(vector shared) { try { data.resize(1,0); - int sumSharedA, sumSharedB, sumSharedAB, tempA, tempB; + double sumSharedA, sumSharedB, sumSharedAB, tempA, tempB; sumSharedA = 0; sumSharedB = 0; sumSharedAB = 0; /*Xi, Yi = abundance of the ith shared OTU in A and B diff --git a/sharedjclass.cpp b/sharedjclass.cpp index d44f22f..ede4e31 100644 --- a/sharedjclass.cpp +++ b/sharedjclass.cpp @@ -13,7 +13,7 @@ EstOutput Jclass::getValues(vector shared) { try { - int S1, S2, S12, tempA, tempB; + double S1, S2, S12, tempA, tempB; S1 = 0; S2 = 0; S12 = 0; tempA = 0; tempB = 0; /*S1, S2 = number of OTUs observed or estimated in A and B diff --git a/sharedkulczynski.cpp b/sharedkulczynski.cpp index 2b4f723..c19df53 100644 --- a/sharedkulczynski.cpp +++ b/sharedkulczynski.cpp @@ -2,7 +2,7 @@ * sharedkulczynski.cpp * Mothur * - * Created by John Westcott on 3/24/09. + * Created by Sarah Westcott on 3/24/09. * Copyright 2009 Schloss Lab UMASS Amherst. All rights reserved. * */ @@ -13,7 +13,7 @@ EstOutput Kulczynski::getValues(vector shared) { try { - int S1, S2, S12, tempA, tempB; + double S1, S2, S12, tempA, tempB; S1 = 0; S2 = 0; S12 = 0; tempA = 0; tempB = 0; /*S1, S2 = number of OTUs observed or estimated in A and B diff --git a/sharedkulczynskicody.cpp b/sharedkulczynskicody.cpp index a2385d4..57fcb24 100644 --- a/sharedkulczynskicody.cpp +++ b/sharedkulczynskicody.cpp @@ -13,7 +13,7 @@ EstOutput KulczynskiCody::getValues(vector shared) { try { - int S1, S2, S12, tempA, tempB; + double S1, S2, S12, tempA, tempB; S1 = 0; S2 = 0; S12 = 0; tempA = 0; tempB = 0; /*S1, S2 = number of OTUs observed or estimated in A and B diff --git a/sharedlennon.cpp b/sharedlennon.cpp index a68da6b..0c0f833 100644 --- a/sharedlennon.cpp +++ b/sharedlennon.cpp @@ -13,7 +13,7 @@ EstOutput Lennon::getValues(vector shared) { try { - int S1, S2, S12, tempA, tempB, min; + double S1, S2, S12, tempA, tempB, min; S1 = 0; S2 = 0; S12 = 0; tempA = 0; tempB = 0; min = 0; /*S1, S2 = number of OTUs observed or estimated in A and B diff --git a/sharedmorisitahorn.cpp b/sharedmorisitahorn.cpp index 3439747..d52ead5 100644 --- a/sharedmorisitahorn.cpp +++ b/sharedmorisitahorn.cpp @@ -14,9 +14,9 @@ EstOutput MorHorn::getValues(vector shared) { try { data.resize(1,0); - float Atotal, Btotal, tempA, tempB; + double Atotal, Btotal, tempA, tempB; Atotal = 0; Btotal = 0; - float morhorn, sumSharedA, sumSharedB, a, b, d; + double morhorn, sumSharedA, sumSharedB, a, b, d; morhorn = 0.0; sumSharedA = 0.0; sumSharedB = 0.0; a = 0.0; b = 0.0; d = 0.0; //get the total values we need to calculate the theta denominator sums diff --git a/sharedochiai.cpp b/sharedochiai.cpp index a1991c9..004e535 100644 --- a/sharedochiai.cpp +++ b/sharedochiai.cpp @@ -13,7 +13,7 @@ EstOutput Ochiai::getValues(vector shared) { try { - int S1, S2, S12, tempA, tempB; + double S1, S2, S12, tempA, tempB; S1 = 0; S2 = 0; S12 = 0; tempA = 0; tempB = 0; /*S1, S2 = number of OTUs observed or estimated in A and B diff --git a/sharedordervector.h b/sharedordervector.h index d0a714c..3568450 100644 --- a/sharedordervector.h +++ b/sharedordervector.h @@ -26,6 +26,16 @@ struct individual { } }; +struct individualFloat { + string group; + int bin; + float abundance; + bool operator()(const individual& i1, const individual& i2) { + return (i1.abundance > i2.abundance); + } +}; + + #include "sabundvector.hpp" #include "rabundvector.hpp" #include "sharedrabundvector.h" diff --git a/sharedsobs.cpp b/sharedsobs.cpp index c8ecba7..36cacd6 100644 --- a/sharedsobs.cpp +++ b/sharedsobs.cpp @@ -16,7 +16,7 @@ EstOutput SharedSobs::getValues(vector shared){ try { data.resize(1,0); - int observed = 0; + double observed = 0; //loop through the species in each group for (int k = 0; k < shared[0]->size(); k++) { diff --git a/sharedsobscollectsummary.cpp b/sharedsobscollectsummary.cpp index c64c825..e2e169c 100644 --- a/sharedsobscollectsummary.cpp +++ b/sharedsobscollectsummary.cpp @@ -16,7 +16,7 @@ EstOutput SharedSobsCS::getValues(vector shared){ try { data.resize(1,0); - int observed = 0; + double observed = 0; int numGroups = shared.size(); for (int i = 0; i < shared[0]->size(); i++) { diff --git a/sharedsorclass.cpp b/sharedsorclass.cpp index 8f60650..2cebc7d 100644 --- a/sharedsorclass.cpp +++ b/sharedsorclass.cpp @@ -13,7 +13,7 @@ EstOutput SorClass::getValues(vector shared) { try { - int S1, S2, S12, tempA, tempB; + double S1, S2, S12, tempA, tempB; S1 = 0; S2 = 0; S12 = 0; tempA = 0; tempB = 0; /*S1, S2 = number of OTUs observed or estimated in A and B diff --git a/sharedthetan.cpp b/sharedthetan.cpp index d118be9..361c7b2 100644 --- a/sharedthetan.cpp +++ b/sharedthetan.cpp @@ -14,9 +14,9 @@ EstOutput ThetaN::getValues(vector shared) { try { data.resize(1,0); - int Atotal, Btotal, tempA, tempB; + double Atotal, Btotal, tempA, tempB; Atotal = 0; Btotal = 0; - float numerator, denominator, thetaN, sumSharedA, sumSharedB, a, b, d; + double numerator, denominator, thetaN, sumSharedA, sumSharedB, a, b, d; numerator = 0.0; denominator = 0.0; thetaN = 0.0; sumSharedA = 0.0; sumSharedB = 0.0; a = 0.0; b = 0.0; d = 0.0; //get the total values we need to calculate the theta denominator sums diff --git a/sharedthetayc.cpp b/sharedthetayc.cpp index ca066dd..4d8a5d9 100644 --- a/sharedthetayc.cpp +++ b/sharedthetayc.cpp @@ -14,19 +14,19 @@ EstOutput ThetaYC::getValues(vector shared) { try { data.resize(3,0.0000); - float Atotal = 0; - float Btotal = 0; - float thetaYC = 0; - float pi = 0; - float qi = 0; - float a = 0; - float b = 0; - float d = 0; + double Atotal = 0; + double Btotal = 0; + double thetaYC = 0; + double pi = 0; + double qi = 0; + double a = 0; + double b = 0; + double d = 0; - float sumPcubed = 0; - float sumQcubed = 0; - float sumPQsq = 0; - float sumPsqQ = 0; + double sumPcubed = 0; + double sumQcubed = 0; + double sumPQsq = 0; + double sumPsqQ = 0; //get the total values we need to calculate the theta denominator sums for (int i = 0; i < shared[0]->size(); i++) { @@ -55,16 +55,16 @@ EstOutput ThetaYC::getValues(vector shared) { if (isnan(thetaYC) || isinf(thetaYC)) { thetaYC = 0; } - float varA = 4 / Atotal * (sumPcubed - a * a); - float varB = 4 / Btotal * (sumQcubed - b * b); - float varD = sumPQsq / Atotal + sumPsqQ / Btotal - d * d * (1/Atotal + 1/Btotal); - float covAD = 2 / Atotal * (sumPsqQ - a * d); - float covBD = 2 / Btotal * (sumPQsq - b* d); + double varA = 4 / Atotal * (sumPcubed - a * a); + double varB = 4 / Btotal * (sumQcubed - b * b); + double varD = sumPQsq / Atotal + sumPsqQ / Btotal - d * d * (1/Atotal + 1/Btotal); + double covAD = 2 / Atotal * (sumPsqQ - a * d); + double covBD = 2 / Btotal * (sumPQsq - b* d); - float varT = d * d * (varA + varB) / pow(a + b - d, (float)4.0) + pow(a+b, (float)2.0) * varD / pow(a+b-d, (float)4.0) - - 2.0 * (a + b) * d / pow(a + b - d, (float)4.0) * (covAD + covBD); + double varT = d * d * (varA + varB) / pow(a + b - d, (double)4.0) + pow(a+b, (double)2.0) * varD / pow(a+b-d, (double)4.0) + - 2.0 * (a + b) * d / pow(a + b - d, (double)4.0) * (covAD + covBD); - float ci = 1.95 * sqrt(varT); + double ci = 1.95 * sqrt(varT); data[0] = thetaYC; data[1] = thetaYC - ci; diff --git a/splitabundcommand.cpp b/splitabundcommand.cpp index 891fb97..924ed69 100644 --- a/splitabundcommand.cpp +++ b/splitabundcommand.cpp @@ -184,15 +184,15 @@ int SplitAbundCommand::execute(){ remove((fileroot + "rare.list").c_str()); remove((fileroot + "abund.list").c_str()); - wroteListFile["rare"] = false; - wroteListFile["abund"] = false; + outputNames.push_back((fileroot + "rare.list")); + outputNames.push_back((fileroot + "abund.list")); }else{ for (int i=0; i::iterator itBool = wroteListFile.begin(); itBool != wroteListFile.end(); itBool++) { - string filename = fileroot + itBool->first; - if ((itBool->first == "rare") || (itBool->first == "abund")) { - filename = fileroot + itBool->first + ".list"; - } - if (itBool->second) { //we wrote to this file - outputNames.push_back(filename); - }else{ - remove(filename.c_str()); - } - } - if (m->control_pressed) { for (int i = 0; i < outputNames.size(); i++) { remove(outputNames[i].c_str()); } return 0; } - }else { //you are using the namefile to determine abundance if (outputDir == "") { outputDir = hasPath(namefile); } @@ -377,19 +364,16 @@ int SplitAbundCommand::writeList(ListVector* thisList) { ofstream aout; ofstream rout; - if (rareNames.size() != 0) { - string rare = outputDir + getRootName(getSimpleName(listfile)) + "rare.list"; - wroteListFile["rare"] = true; - openOutputFileAppend(rare, rout); - rout << thisList->getLabel() << '\t' << numRareBins << '\t'; - } + string rare = outputDir + getRootName(getSimpleName(listfile)) + "rare.list"; + openOutputFileAppend(rare, rout); + outputNames.push_back(rare); - if (abundNames.size() != 0) { - string abund = outputDir + getRootName(getSimpleName(listfile)) + "abund.list"; - wroteListFile["abund"] = true; - openOutputFileAppend(abund, aout); - aout << thisList->getLabel() << '\t' << numAbundBins << '\t'; - } + string abund = outputDir + getRootName(getSimpleName(listfile)) + "abund.list"; + openOutputFileAppend(abund, aout); + outputNames.push_back(abund); + + if (rareNames.size() != 0) { rout << thisList->getLabel() << '\t' << numRareBins << '\t'; } + if (abundNames.size() != 0) { aout << thisList->getLabel() << '\t' << numAbundBins << '\t'; } for (int i = 0; i < thisList->getNumBins(); i++) { if (m->control_pressed) { break; } @@ -402,14 +386,17 @@ int SplitAbundCommand::writeList(ListVector* thisList) { else { aout << bin << '\t'; } } - if (rareNames.size() != 0) { rout << endl; rout.close(); } - if (abundNames.size() != 0) { aout << endl; aout.close(); } - + if (rareNames.size() != 0) { rout << endl; } + if (abundNames.size() != 0) { aout << endl; } + + rout.close(); + aout.close(); + }else{ //parse names by abundance and group string fileroot = outputDir + getRootName(getSimpleName(listfile)); ofstream* temp; ofstream* temp2; - map wroteFile; + //map wroteFile; map filehandles; map::iterator it3; @@ -474,7 +461,6 @@ int SplitAbundCommand::writeList(ListVector* thisList) { //end list vector for (it3 = filehandles.begin(); it3 != filehandles.end(); it3++) { (*(filehandles[it3->first])) << thisList->getLabel() << '\t' << groupNumBins[it3->first] << '\t' << groupVector[it3->first] << endl; // label numBins listvector for that group - wroteListFile[it3->first] = true; (*(filehandles[it3->first])).close(); delete it3->second; } @@ -584,33 +570,32 @@ int SplitAbundCommand::writeNames() { //namefile ofstream aout; ofstream rout; + string rare = outputDir + getRootName(getSimpleName(namefile)) + "rare.names"; + openOutputFile(rare, rout); + outputNames.push_back(rare); + + string abund = outputDir + getRootName(getSimpleName(namefile)) + "abund.names"; + openOutputFile(abund, aout); + outputNames.push_back(abund); + if (rareNames.size() != 0) { - string rare = outputDir + getRootName(getSimpleName(namefile)) + "rare.names"; - openOutputFile(rare, rout); - outputNames.push_back(rare); - for (set::iterator itRare = rareNames.begin(); itRare != rareNames.end(); itRare++) { rout << (*itRare) << '\t' << nameMap[(*itRare)] << endl; } - rout.close(); } + rout.close(); if (abundNames.size() != 0) { - string abund = outputDir + getRootName(getSimpleName(namefile)) + "abund.names"; - openOutputFile(abund, aout); - outputNames.push_back(abund); - for (set::iterator itAbund = abundNames.begin(); itAbund != abundNames.end(); itAbund++) { aout << (*itAbund) << '\t' << nameMap[(*itAbund)] << endl; } - aout.close(); } - + aout.close(); + }else{ //parse names by abundance and group string fileroot = outputDir + getRootName(getSimpleName(namefile)); ofstream* temp; ofstream* temp2; - map wroteFile; map filehandles; map::iterator it3; @@ -622,9 +607,6 @@ int SplitAbundCommand::writeNames() { //namefile openOutputFile(fileroot + Groups[i] + ".rare.names", *(filehandles[Groups[i]+".rare"])); openOutputFile(fileroot + Groups[i] + ".abund.names", *(filehandles[Groups[i]+".abund"])); - - wroteFile[Groups[i] + ".rare"] = false; - wroteFile[Groups[i] + ".abund"] = false; } for (map::iterator itName = nameMap.begin(); itName != nameMap.end(); itName++) { @@ -654,17 +636,13 @@ int SplitAbundCommand::writeNames() { //namefile } } - for (itout = outputStrings.begin(); itout != outputStrings.end(); itout++) { - *(filehandles[itout->first]) << itout->second << endl; - wroteFile[itout->first] = true; - } + for (itout = outputStrings.begin(); itout != outputStrings.end(); itout++) { *(filehandles[itout->first]) << itout->second << endl; } } for (it3 = filehandles.begin(); it3 != filehandles.end(); it3++) { (*(filehandles[it3->first])).close(); - if (wroteFile[it3->first] == true) { outputNames.push_back(fileroot + it3->first + ".names"); } - else { remove((it3->first).c_str()); } + outputNames.push_back(fileroot + it3->first + ".names"); delete it3->second; } } @@ -688,32 +666,29 @@ int SplitAbundCommand::writeAccnos(string tag) { ofstream aout; ofstream rout; - if (rareNames.size() != 0) { - string rare = outputDir + getRootName(getSimpleName(inputFile)) + tag + "rare.accnos"; - openOutputFile(rare, rout); - outputNames.push_back(rare); - - for (set::iterator itRare = rareNames.begin(); itRare != rareNames.end(); itRare++) { - rout << (*itRare) << endl; - } - rout.close(); + + string rare = outputDir + getRootName(getSimpleName(inputFile)) + tag + "rare.accnos"; + openOutputFile(rare, rout); + outputNames.push_back(rare); + + for (set::iterator itRare = rareNames.begin(); itRare != rareNames.end(); itRare++) { + rout << (*itRare) << endl; } + rout.close(); + + string abund = outputDir + getRootName(getSimpleName(inputFile)) + tag + "abund.accnos"; + openOutputFile(abund, aout); + outputNames.push_back(abund); - if (abundNames.size() != 0) { - string abund = outputDir + getRootName(getSimpleName(inputFile)) + tag + "abund.accnos"; - openOutputFile(abund, aout); - outputNames.push_back(abund); - - for (set::iterator itAbund = abundNames.begin(); itAbund != abundNames.end(); itAbund++) { - aout << (*itAbund) << endl; - } - aout.close(); + for (set::iterator itAbund = abundNames.begin(); itAbund != abundNames.end(); itAbund++) { + aout << (*itAbund) << endl; } + aout.close(); + }else{ //parse names by abundance and group string fileroot = outputDir + getRootName(getSimpleName(inputFile)); ofstream* temp; ofstream* temp2; - map wroteFile; map filehandles; map::iterator it3; @@ -725,9 +700,6 @@ int SplitAbundCommand::writeAccnos(string tag) { openOutputFile(fileroot + tag + Groups[i] + ".rare.accnos", *(filehandles[Groups[i]+".rare"])); openOutputFile(fileroot + tag + Groups[i] + ".abund.accnos", *(filehandles[Groups[i]+".abund"])); - - wroteFile[Groups[i] + ".rare"] = false; - wroteFile[Groups[i] + ".abund"] = false; } //write rare @@ -736,7 +708,6 @@ int SplitAbundCommand::writeAccnos(string tag) { if (inUsersGroups(group, Groups)) { //only add if this is in a group we want *(filehandles[group+".rare"]) << *itRare << endl; - wroteFile[group+".rare"] = true; } } @@ -746,15 +717,13 @@ int SplitAbundCommand::writeAccnos(string tag) { if (inUsersGroups(group, Groups)) { //only add if this is in a group we want *(filehandles[group+".abund"]) << *itAbund << endl; - wroteFile[group+".abund"] = true; } } //close files for (it3 = filehandles.begin(); it3 != filehandles.end(); it3++) { (*(filehandles[it3->first])).close(); - if (wroteFile[it3->first] == true) { outputNames.push_back(fileroot + tag + it3->first + ".accnos"); } - else { remove((fileroot + tag + it3->first + ".accnos").c_str()); } + outputNames.push_back(fileroot + tag + it3->first + ".accnos"); delete it3->second; } } @@ -777,19 +746,14 @@ int SplitAbundCommand::parseGroup(string tag) { //namefile ofstream aout; ofstream rout; - if (rareNames.size() != 0) { - string rare = outputDir + getRootName(getSimpleName(groupfile)) + tag + "rare.groups"; - openOutputFile(rare, rout); - outputNames.push_back(rare); - } + string rare = outputDir + getRootName(getSimpleName(groupfile)) + tag + "rare.groups"; + openOutputFile(rare, rout); + outputNames.push_back(rare); + + string abund = outputDir + getRootName(getSimpleName(groupfile)) + tag + "abund.groups"; + openOutputFile(abund, aout); + outputNames.push_back(abund); - if (abundNames.size() != 0) { - string abund = outputDir + getRootName(getSimpleName(groupfile)) + tag + "abund.groups"; - openOutputFile(abund, aout); - outputNames.push_back(abund); - } - - for (map::iterator itName = nameMap.begin(); itName != nameMap.end(); itName++) { vector names; splitAtComma(itName->second, names); //parses bin into individual sequence names @@ -810,14 +774,13 @@ int SplitAbundCommand::parseGroup(string tag) { //namefile } } - if (rareNames.size() != 0) { rout.close(); } - if (abundNames.size() != 0) { aout.close(); } + rout.close(); + aout.close(); }else{ //parse names by abundance and group string fileroot = outputDir + getRootName(getSimpleName(groupfile)); ofstream* temp; ofstream* temp2; - map wroteFile; map filehandles; map::iterator it3; @@ -829,9 +792,6 @@ int SplitAbundCommand::parseGroup(string tag) { //namefile openOutputFile(fileroot + tag + Groups[i] + ".rare.groups", *(filehandles[Groups[i]+".rare"])); openOutputFile(fileroot + tag + Groups[i] + ".abund.groups", *(filehandles[Groups[i]+".abund"])); - - wroteFile[Groups[i] + ".rare"] = false; - wroteFile[Groups[i] + ".abund"] = false; } for (map::iterator itName = nameMap.begin(); itName != nameMap.end(); itName++) { @@ -851,15 +811,13 @@ int SplitAbundCommand::parseGroup(string tag) { //namefile if (inUsersGroups(group, Groups)) { //only add if this is in a group we want *(filehandles[group+rareAbund]) << names[i] << '\t' << group << endl; - wroteFile[group+rareAbund] = true; } } } for (it3 = filehandles.begin(); it3 != filehandles.end(); it3++) { (*(filehandles[it3->first])).close(); - if (wroteFile[it3->first] == true) { outputNames.push_back(fileroot + tag + it3->first + ".groups"); } - else { remove((fileroot + tag + it3->first + ".groups").c_str()); } + outputNames.push_back(fileroot + tag + it3->first + ".groups"); delete it3->second; } } @@ -882,19 +840,14 @@ int SplitAbundCommand::parseFasta(string tag) { //namefile ofstream aout; ofstream rout; - if (rareNames.size() != 0) { - string rare = outputDir + getRootName(getSimpleName(fastafile)) + tag + "rare.fasta"; - openOutputFile(rare, rout); - outputNames.push_back(rare); - } - - if (abundNames.size() != 0) { - string abund = outputDir + getRootName(getSimpleName(fastafile)) + tag + "abund.fasta"; - openOutputFile(abund, aout); - outputNames.push_back(abund); - } - - + string rare = outputDir + getRootName(getSimpleName(fastafile)) + tag + "rare.fasta"; + openOutputFile(rare, rout); + outputNames.push_back(rare); + + string abund = outputDir + getRootName(getSimpleName(fastafile)) + tag + "abund.fasta"; + openOutputFile(abund, aout); + outputNames.push_back(abund); + //open input file ifstream in; openInputFile(fastafile, in); @@ -922,14 +875,13 @@ int SplitAbundCommand::parseFasta(string tag) { //namefile } } in.close(); - if (rareNames.size() != 0) { rout.close(); } - if (abundNames.size() != 0) { aout.close(); } + rout.close(); + aout.close(); }else{ //parse names by abundance and group string fileroot = outputDir + getRootName(getSimpleName(fastafile)); ofstream* temp; ofstream* temp2; - map wroteFile; map filehandles; map::iterator it3; @@ -941,9 +893,6 @@ int SplitAbundCommand::parseFasta(string tag) { //namefile openOutputFile(fileroot + tag + Groups[i] + ".rare.fasta", *(filehandles[Groups[i]+".rare"])); openOutputFile(fileroot + tag + Groups[i] + ".abund.fasta", *(filehandles[Groups[i]+".abund"])); - - wroteFile[Groups[i] + ".rare"] = false; - wroteFile[Groups[i] + ".abund"] = false; } //open input file @@ -977,7 +926,6 @@ int SplitAbundCommand::parseFasta(string tag) { //namefile if (inUsersGroups(group, Groups)) { //only add if this is in a group we want seq.printSequence(*(filehandles[group+rareAbund])); - wroteFile[group+rareAbund] = true; }else if(group == "not found") { m->mothurOut(seq.getName() + " is not in your groupfile. Ignoring."); m->mothurOutEndLine(); } @@ -989,8 +937,7 @@ int SplitAbundCommand::parseFasta(string tag) { //namefile for (it3 = filehandles.begin(); it3 != filehandles.end(); it3++) { (*(filehandles[it3->first])).close(); - if (wroteFile[it3->first] == true) { outputNames.push_back(fileroot + tag + it3->first + ".fasta"); } - else { remove((fileroot + tag + it3->first + ".fasta").c_str()); } + outputNames.push_back(fileroot + tag + it3->first + ".fasta"); delete it3->second; } } diff --git a/splitabundcommand.h b/splitabundcommand.h index fb88380..df5bd55 100644 --- a/splitabundcommand.h +++ b/splitabundcommand.h @@ -54,7 +54,7 @@ private: vector Groups; bool abort, allLines, accnos; int cutoff; - map wroteListFile; + //map wroteListFile; map nameMap;