X-Git-Url: https://git.donarmstrong.com/?a=blobdiff_plain;f=seqsummarycommand.cpp;h=8730997883ac5da6f8588e216f8311285303160d;hb=32abf5d58aed8087cd0bcb32d7aa9053f103f2c1;hp=f551eddacd00dade4aff87526f0b5aca59ad2b58;hpb=4a2d841cb97fb02351022efe9d7068b1dc212bf9;p=mothur.git diff --git a/seqsummarycommand.cpp b/seqsummarycommand.cpp index f551edd..8730997 100644 --- a/seqsummarycommand.cpp +++ b/seqsummarycommand.cpp @@ -10,24 +10,71 @@ #include "seqsummarycommand.h" #include "sequence.hpp" +//********************************************************************************************************************** +vector SeqSummaryCommand::getValidParameters(){ + try { + string Array[] = {"fasta","name","processors","outputdir","inputdir"}; + vector myArray (Array, Array+(sizeof(Array)/sizeof(string))); + return myArray; + } + catch(exception& e) { + m->errorOut(e, "SeqSummaryCommand", "getValidParameters"); + exit(1); + } +} +//********************************************************************************************************************** +SeqSummaryCommand::SeqSummaryCommand(){ + try { + abort = true; calledHelp = true; + vector tempOutNames; + outputTypes["summary"] = tempOutNames; + } + catch(exception& e) { + m->errorOut(e, "SeqSummaryCommand", "SeqSummaryCommand"); + exit(1); + } +} +//********************************************************************************************************************** +vector SeqSummaryCommand::getRequiredParameters(){ + try { + string Array[] = {"fasta"}; + vector myArray (Array, Array+(sizeof(Array)/sizeof(string))); + return myArray; + } + catch(exception& e) { + m->errorOut(e, "SeqSummaryCommand", "getRequiredParameters"); + exit(1); + } +} +//********************************************************************************************************************** +vector SeqSummaryCommand::getRequiredFiles(){ + try { + vector myArray; + return myArray; + } + catch(exception& e) { + m->errorOut(e, "SeqSummaryCommand", "getRequiredFiles"); + exit(1); + } +} //*************************************************************************************************************** SeqSummaryCommand::SeqSummaryCommand(string option) { try { - abort = false; + abort = false; calledHelp = false; //allow user to run help - if(option == "help") { help(); abort = true; } + if(option == "help") { help(); abort = true; calledHelp = true; } else { //valid paramters for this command - string Array[] = {"fasta","processors","outputdir","inputdir"}; + string Array[] = {"fasta","name","processors","outputdir","inputdir"}; vector myArray (Array, Array+(sizeof(Array)/sizeof(string))); OptionParser parser(option); map parameters = parser.getParameters(); - ValidParameters validParameter; + ValidParameters validParameter("summary.seqs"); map::iterator it; //check to make sure all parameters are valid for command @@ -43,21 +90,37 @@ SeqSummaryCommand::SeqSummaryCommand(string option) { it = parameters.find("fasta"); //user has given a template file if(it != parameters.end()){ - path = hasPath(it->second); + path = m->hasPath(it->second); //if the user has not given a path then, add inputdir. else leave path alone. if (path == "") { parameters["fasta"] = inputDir + it->second; } } + + it = parameters.find("name"); + //user has given a template file + if(it != parameters.end()){ + path = m->hasPath(it->second); + //if the user has not given a path then, add inputdir. else leave path alone. + if (path == "") { parameters["name"] = inputDir + it->second; } + } } + //initialize outputTypes + vector tempOutNames; + outputTypes["summary"] = tempOutNames; + //check for required parameters fastafile = validParameter.validFile(parameters, "fasta", true); if (fastafile == "not open") { abort = true; } else if (fastafile == "not found") { fastafile = ""; m->mothurOut("fasta is a required parameter for the summary.seqs command."); m->mothurOutEndLine(); abort = true; } + namefile = validParameter.validFile(parameters, "name", true); + if (namefile == "not open") { namefile = ""; abort = true; } + else if (namefile == "not found") { namefile = ""; } + //if the user changes the output directory command factory will send this info to us in the output parameter outputDir = validParameter.validFile(parameters, "outputdir", false); if (outputDir == "not found"){ outputDir = ""; - outputDir += hasPath(fastafile); //if user entered a file with a path then preserve it + outputDir += m->hasPath(fastafile); //if user entered a file with a path then preserve it } string temp = validParameter.validFile(parameters, "processors", false); if (temp == "not found"){ temp = "1"; } @@ -76,7 +139,8 @@ SeqSummaryCommand::SeqSummaryCommand(string option) { void SeqSummaryCommand::help(){ try { m->mothurOut("The summary.seqs command reads a fastafile and summarizes the sequences.\n"); - m->mothurOut("The summary.seqs command parameters are fasta and processors, fasta is required.\n"); + m->mothurOut("The summary.seqs command parameters are fasta, name and processors, fasta is required.\n"); + m->mothurOut("The name parameter allows you to enter a name file associated with your fasta file. \n"); m->mothurOut("The summary.seqs command should be in the following format: \n"); m->mothurOut("summary.seqs(fasta=yourFastaFile, processors=2) \n"); m->mothurOut("Note: No spaces between parameter labels (i.e. fasta), '=' and parameters (i.e.yourFastaFile).\n\n"); @@ -96,9 +160,9 @@ SeqSummaryCommand::~SeqSummaryCommand(){ /* do nothing */ } int SeqSummaryCommand::execute(){ try{ - if (abort == true) { return 0; } + if (abort == true) { if (calledHelp) { return 0; } return 2; } - string summaryFile = outputDir + getSimpleName(fastafile) + ".summary"; + string summaryFile = outputDir + m->getSimpleName(fastafile) + ".summary"; int numSeqs = 0; @@ -107,13 +171,17 @@ int SeqSummaryCommand::execute(){ vector seqLength; vector ambigBases; vector longHomoPolymer; + + if (namefile != "") { readNames(); } + + if (m->control_pressed) { return 0; } #ifdef USE_MPI int pid, numSeqsPerProcessor; int tag = 2001; int startTag = 1; int endTag = 2; int lengthTag = 3; int baseTag = 4; int lhomoTag = 5; int outMode=MPI_MODE_CREATE|MPI_MODE_WRONLY; - vector MPIPos; + vector MPIPos; MPI_Status status; MPI_Status statusOut; @@ -135,7 +203,7 @@ int SeqSummaryCommand::execute(){ if (pid == 0) { //you are the root process //print header - string outputString = "seqname\tstart\tend\tnbases\tambigs\tpolymer\n"; + string outputString = "seqname\tstart\tend\tnbases\tambigs\tpolymer\tnumSeqs\n"; int length = outputString.length(); char* buf2 = new char[length]; memcpy(buf2, outputString.c_str(), length); @@ -143,12 +211,13 @@ int SeqSummaryCommand::execute(){ MPI_File_write_shared(outMPI, buf2, length, MPI_CHAR, &statusOut); delete buf2; - MPIPos = setFilePosFasta(fastafile, numSeqs); //fills MPIPos, returns numSeqs + MPIPos = m->setFilePosFasta(fastafile, numSeqs); //fills MPIPos, returns numSeqs + + for(int i = 1; i < processors; i++) { + MPI_Send(&numSeqs, 1, MPI_INT, i, tag, MPI_COMM_WORLD); + MPI_Send(&MPIPos[0], (numSeqs+1), MPI_LONG, i, tag, MPI_COMM_WORLD); + } - //send file positions to all processes - MPI_Bcast(&numSeqs, 1, MPI_INT, 0, MPI_COMM_WORLD); //send numSeqs - MPI_Bcast(&MPIPos[0], (numSeqs+1), MPI_LONG, 0, MPI_COMM_WORLD); //send file pos - //figure out how many sequences you have to do numSeqsPerProcessor = numSeqs / processors; int startIndex = pid * numSeqsPerProcessor; @@ -156,54 +225,55 @@ int SeqSummaryCommand::execute(){ //do your part MPICreateSummary(startIndex, numSeqsPerProcessor, startPosition, endPosition, seqLength, ambigBases, longHomoPolymer, inMPI, outMPI, MPIPos); - - if (m->control_pressed) { MPI_File_close(&inMPI); MPI_File_close(&outMPI); return 0; } - - //get the info from the child processes - for(int i = 1; i < processors; i++) { - - int size; - MPI_Recv(&size, 1, MPI_INT, i, tag, MPI_COMM_WORLD, &status); - - vector temp; temp.resize(size+1); - - for(int j = 0; j < 5; j++) { - - MPI_Recv(&temp[0], (size+1), MPI_INT, i, 2001, MPI_COMM_WORLD, &status); - int receiveTag = temp[temp.size()-1]; //child process added a int to the end to indicate what count this is for - - if (receiveTag == startTag) { - for (int k = 0; k < size; k++) { startPosition.push_back(temp[k]); } - }else if (receiveTag == endTag) { - for (int k = 0; k < size; k++) { endPosition.push_back(temp[k]); } - }else if (receiveTag == lengthTag) { - for (int k = 0; k < size; k++) { seqLength.push_back(temp[k]); } - }else if (receiveTag == baseTag) { - for (int k = 0; k < size; k++) { ambigBases.push_back(temp[k]); } - }else if (receiveTag == lhomoTag) { - for (int k = 0; k < size; k++) { longHomoPolymer.push_back(temp[k]); } - } - } - } - - }else { //i am the child process - MPI_Bcast(&numSeqs, 1, MPI_INT, 0, MPI_COMM_WORLD); //get numSeqs + MPI_Recv(&numSeqs, 1, MPI_INT, 0, tag, MPI_COMM_WORLD, &status); MPIPos.resize(numSeqs+1); - MPI_Bcast(&MPIPos[0], (numSeqs+1), MPI_LONG, 0, MPI_COMM_WORLD); //get file positions - + MPI_Recv(&MPIPos[0], (numSeqs+1), MPI_LONG, 0, tag, MPI_COMM_WORLD, &status); + //figure out how many sequences you have to align numSeqsPerProcessor = numSeqs / processors; int startIndex = pid * numSeqsPerProcessor; if(pid == (processors - 1)){ numSeqsPerProcessor = numSeqs - pid * numSeqsPerProcessor; } - + //do your part MPICreateSummary(startIndex, numSeqsPerProcessor, startPosition, endPosition, seqLength, ambigBases, longHomoPolymer, inMPI, outMPI, MPIPos); - - if (m->control_pressed) { MPI_File_close(&inMPI); MPI_File_close(&outMPI); return 0; } - + } + + MPI_File_close(&inMPI); + MPI_File_close(&outMPI); + MPI_Barrier(MPI_COMM_WORLD); //make everyone wait - just in case + + if (pid == 0) { + //get the info from the child processes + for(int i = 1; i < processors; i++) { + int size; + MPI_Recv(&size, 1, MPI_INT, i, tag, MPI_COMM_WORLD, &status); + + vector temp; temp.resize(size+1); + + for(int j = 0; j < 5; j++) { + + MPI_Recv(&temp[0], (size+1), MPI_INT, i, 2001, MPI_COMM_WORLD, &status); + int receiveTag = temp[temp.size()-1]; //child process added a int to the end to indicate what count this is for + + if (receiveTag == startTag) { + for (int k = 0; k < size; k++) { startPosition.push_back(temp[k]); } + }else if (receiveTag == endTag) { + for (int k = 0; k < size; k++) { endPosition.push_back(temp[k]); } + }else if (receiveTag == lengthTag) { + for (int k = 0; k < size; k++) { seqLength.push_back(temp[k]); } + }else if (receiveTag == baseTag) { + for (int k = 0; k < size; k++) { ambigBases.push_back(temp[k]); } + }else if (receiveTag == lhomoTag) { + for (int k = 0; k < size; k++) { longHomoPolymer.push_back(temp[k]); } + } + } + } + + }else{ + //send my counts int size = startPosition.size(); MPI_Send(&size, 1, MPI_INT, 0, tag, MPI_COMM_WORLD); @@ -218,45 +288,33 @@ int SeqSummaryCommand::execute(){ ierr = MPI_Send(&(ambigBases[0]), (size+1), MPI_INT, 0, 2001, MPI_COMM_WORLD); longHomoPolymer.push_back(lhomoTag); ierr = MPI_Send(&(longHomoPolymer[0]), (size+1), MPI_INT, 0, 2001, MPI_COMM_WORLD); - } - MPI_File_close(&inMPI); - MPI_File_close(&outMPI); - + MPI_Barrier(MPI_COMM_WORLD); //make everyone wait - just in case #else + vector positions = m->divideFile(fastafile, processors); + + for (int i = 0; i < (positions.size()-1); i++) { + lines.push_back(new linePair(positions[i], positions[(i+1)])); + } + #if defined (__APPLE__) || (__MACH__) || (linux) || (__linux) if(processors == 1){ - ifstream inFASTA; - openInputFile(fastafile, inFASTA); - numSeqs=count(istreambuf_iterator(inFASTA),istreambuf_iterator(), '>'); - inFASTA.close(); - - lines.push_back(new linePair(0, numSeqs)); - - driverCreateSummary(startPosition, endPosition, seqLength, ambigBases, longHomoPolymer, fastafile, summaryFile, lines[0]); + numSeqs = driverCreateSummary(startPosition, endPosition, seqLength, ambigBases, longHomoPolymer, fastafile, summaryFile, lines[0]); }else{ - numSeqs = setLines(fastafile); - createProcessesCreateSummary(startPosition, endPosition, seqLength, ambigBases, longHomoPolymer, fastafile, summaryFile); + numSeqs = createProcessesCreateSummary(startPosition, endPosition, seqLength, ambigBases, longHomoPolymer, fastafile, summaryFile); rename((summaryFile + toString(processIDS[0]) + ".temp").c_str(), summaryFile.c_str()); //append files for(int i=1;iappendFiles((summaryFile + toString(processIDS[i]) + ".temp"), summaryFile); remove((summaryFile + toString(processIDS[i]) + ".temp").c_str()); } } if (m->control_pressed) { return 0; } #else - ifstream inFASTA; - openInputFile(fastafile, inFASTA); - numSeqs=count(istreambuf_iterator(inFASTA),istreambuf_iterator(), '>'); - inFASTA.close(); - - lines.push_back(new linePair(0, numSeqs)); - - driverCreateSummary(startPosition, endPosition, seqLength, ambigBases, longHomoPolymer, fastafile, summaryFile, lines[0]); + numSeqs = driverCreateSummary(startPosition, endPosition, seqLength, ambigBases, longHomoPolymer, fastafile, summaryFile, lines[0]); if (m->control_pressed) { return 0; } #endif #endif @@ -293,13 +351,14 @@ int SeqSummaryCommand::execute(){ m->mothurOut("75%-tile:\t" + toString(startPosition[ptile75]) + "\t" + toString(endPosition[ptile75]) + "\t" + toString(seqLength[ptile75]) + "\t" + toString(ambigBases[ptile75]) + "\t" + toString(longHomoPolymer[ptile75])); m->mothurOutEndLine(); m->mothurOut("97.5%-tile:\t" + toString(startPosition[ptile97_5]) + "\t" + toString(endPosition[ptile97_5]) + "\t" + toString(seqLength[ptile97_5]) + "\t" + toString(ambigBases[ptile97_5]) + "\t" + toString(longHomoPolymer[ptile97_5])); m->mothurOutEndLine(); m->mothurOut("Maximum:\t" + toString(startPosition[ptile100]) + "\t" + toString(endPosition[ptile100]) + "\t" + toString(seqLength[ptile100]) + "\t" + toString(ambigBases[ptile100]) + "\t" + toString(longHomoPolymer[ptile100])); m->mothurOutEndLine(); - m->mothurOut("# of Seqs:\t" + toString(numSeqs)); m->mothurOutEndLine(); + if (namefile == "") { m->mothurOut("# of Seqs:\t" + toString(numSeqs)); m->mothurOutEndLine(); } + else { m->mothurOut("# of unique seqs:\t" + toString(numSeqs)); m->mothurOutEndLine(); m->mothurOut("total # of seqs:\t" + toString(startPosition.size())); m->mothurOutEndLine(); } if (m->control_pressed) { remove(summaryFile.c_str()); return 0; } m->mothurOutEndLine(); m->mothurOut("Output File Name: "); m->mothurOutEndLine(); - m->mothurOut(summaryFile); m->mothurOutEndLine(); + m->mothurOut(summaryFile); m->mothurOutEndLine(); outputNames.push_back(summaryFile); outputTypes["summary"].push_back(summaryFile); m->mothurOutEndLine(); #ifdef USE_MPI @@ -314,44 +373,74 @@ int SeqSummaryCommand::execute(){ } } /**************************************************************************************/ -int SeqSummaryCommand::driverCreateSummary(vector& startPosition, vector& endPosition, vector& seqLength, vector& ambigBases, vector& longHomoPolymer, string filename, string sumFile, linePair* line) { +int SeqSummaryCommand::driverCreateSummary(vector& startPosition, vector& endPosition, vector& seqLength, vector& ambigBases, vector& longHomoPolymer, string filename, string sumFile, linePair* filePos) { try { ofstream outSummary; - openOutputFile(sumFile, outSummary); + m->openOutputFile(sumFile, outSummary); //print header if you are process 0 - if (line->start == 0) { - outSummary << "seqname\tstart\tend\tnbases\tambigs\tpolymer" << endl; + if (filePos->start == 0) { + outSummary << "seqname\tstart\tend\tnbases\tambigs\tpolymer\tnumSeqs" << endl; } ifstream in; - openInputFile(filename, in); + m->openInputFile(filename, in); - in.seekg(line->start); - - for(int i=0;inum;i++){ + in.seekg(filePos->start); + + bool done = false; + int count = 0; + + while (!done) { if (m->control_pressed) { in.close(); outSummary.close(); return 1; } - Sequence current(in); + Sequence current(in); m->gobble(in); + if (current.getName() != "") { - startPosition.push_back(current.getStartPos()); - endPosition.push_back(current.getEndPos()); - seqLength.push_back(current.getNumBases()); - ambigBases.push_back(current.getAmbigBases()); - longHomoPolymer.push_back(current.getLongHomoPolymer()); + int num = 1; + if (namefile != "") { + //make sure this sequence is in the namefile, else error + map::iterator it = nameMap.find(current.getName()); + + if (it == nameMap.end()) { m->mothurOut("[ERROR]: " + current.getName() + " is not in your namefile, please correct."); m->mothurOutEndLine(); m->control_pressed = true; } + else { num = it->second; } + } + + //for each sequence this sequence represents + for (int i = 0; i < num; i++) { + startPosition.push_back(current.getStartPos()); + endPosition.push_back(current.getEndPos()); + seqLength.push_back(current.getNumBases()); + ambigBases.push_back(current.getAmbigBases()); + longHomoPolymer.push_back(current.getLongHomoPolymer()); + } + + count++; outSummary << current.getName() << '\t'; outSummary << current.getStartPos() << '\t' << current.getEndPos() << '\t'; outSummary << current.getNumBases() << '\t' << current.getAmbigBases() << '\t'; - outSummary << current.getLongHomoPolymer() << endl; + outSummary << current.getLongHomoPolymer() << '\t' << num << endl; } - gobble(in); + + #if defined (__APPLE__) || (__MACH__) || (linux) || (__linux) + unsigned long int pos = in.tellg(); + if ((pos == -1) || (pos >= filePos->end)) { break; } + #else + if (in.eof()) { break; } + #endif + + //report progress + //if((count) % 100 == 0){ m->mothurOut(toString(count)); m->mothurOutEndLine(); } } + //report progress + //if((count) % 100 != 0){ m->mothurOut(toString(count)); m->mothurOutEndLine(); } + in.close(); - return 0; + return count; } catch(exception& e) { m->errorOut(e, "SeqSummaryCommand", "driverCreateSummary"); @@ -360,11 +449,13 @@ int SeqSummaryCommand::driverCreateSummary(vector& startPosition, vector& startPosition, vector& endPosition, vector& seqLength, vector& ambigBases, vector& longHomoPolymer, MPI_File& inMPI, MPI_File& outMPI, vector& MPIPos) { +int SeqSummaryCommand::MPICreateSummary(int start, int num, vector& startPosition, vector& endPosition, vector& seqLength, vector& ambigBases, vector& longHomoPolymer, MPI_File& inMPI, MPI_File& outMPI, vector& MPIPos) { try { + int pid; MPI_Status status; - + MPI_Comm_rank(MPI_COMM_WORLD, &pid); + for(int i=0;icontrol_pressed) { return 0; } @@ -383,14 +474,27 @@ int SeqSummaryCommand::MPICreateSummary(int start, int num, vector& startPo Sequence current(iss); if (current.getName() != "") { - startPosition.push_back(current.getStartPos()); - endPosition.push_back(current.getEndPos()); - seqLength.push_back(current.getNumBases()); - ambigBases.push_back(current.getAmbigBases()); - longHomoPolymer.push_back(current.getLongHomoPolymer()); + + int num = 1; + if (namefile != "") { + //make sure this sequence is in the namefile, else error + map::iterator it = nameMap.find(current.getName()); + + if (it == nameMap.end()) { m->mothurOut("[ERROR]: " + current.getName() + " is not in your namefile, please correct."); m->mothurOutEndLine(); m->control_pressed = true; } + else { num = it->second; } + } + + //for each sequence this sequence represents + for (int i = 0; i < num; i++) { + startPosition.push_back(current.getStartPos()); + endPosition.push_back(current.getEndPos()); + seqLength.push_back(current.getNumBases()); + ambigBases.push_back(current.getAmbigBases()); + longHomoPolymer.push_back(current.getLongHomoPolymer()); + } string outputString = current.getName() + "\t" + toString(current.getStartPos()) + "\t" + toString(current.getEndPos()) + "\t"; - outputString += toString(current.getNumBases()) + "\t" + toString(current.getAmbigBases()) + "\t" + toString(current.getLongHomoPolymer()) + "\n"; + outputString += toString(current.getNumBases()) + "\t" + toString(current.getAmbigBases()) + "\t" + toString(current.getLongHomoPolymer()) + "\t" + toString(num) + "\n"; //output to file length = outputString.length(); @@ -399,7 +503,7 @@ int SeqSummaryCommand::MPICreateSummary(int start, int num, vector& startPo MPI_File_write_shared(outMPI, buf3, length, MPI_CHAR, &status); delete buf3; - } + } } return 0; @@ -415,20 +519,40 @@ int SeqSummaryCommand::createProcessesCreateSummary(vector& startPosition, try { #if defined (__APPLE__) || (__MACH__) || (linux) || (__linux) int process = 0; - int exitCommand = 1; + int num = 0; processIDS.clear(); //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){ - driverCreateSummary(startPosition, endPosition, seqLength, ambigBases, longHomoPolymer, fastafile, sumFile + toString(getpid()) + ".temp", lines[process]); + num = driverCreateSummary(startPosition, endPosition, seqLength, ambigBases, longHomoPolymer, fastafile, sumFile + toString(getpid()) + ".temp", lines[process]); + + //pass numSeqs to parent + ofstream out; + string tempFile = fastafile + toString(getpid()) + ".num.temp"; + m->openOutputFile(tempFile, out); + + out << num << endl; + out << startPosition.size() << endl; + for (int k = 0; k < startPosition.size(); k++) { out << startPosition[k] << '\t'; } out << endl; + for (int k = 0; k < endPosition.size(); k++) { out << endPosition[k] << '\t'; } out << endl; + for (int k = 0; k < seqLength.size(); k++) { out << seqLength[k] << '\t'; } out << endl; + for (int k = 0; k < ambigBases.size(); k++) { out << ambigBases[k] << '\t'; } out << endl; + for (int k = 0; k < longHomoPolymer.size(); k++) { out << longHomoPolymer[k] << '\t'; } out << endl; + + out.close(); + + exit(0); + }else { + m->mothurOut("[ERROR]: unable to spawn the necessary processes."); m->mothurOutEndLine(); + for (int i = 0; i < processIDS.size(); i++) { kill (processIDS[i], SIGINT); } exit(0); - }else { m->mothurOut("unable to spawn the necessary processes."); m->mothurOutEndLine(); exit(0); } + } } //force parent to wait until all the processes are done @@ -437,7 +561,26 @@ int SeqSummaryCommand::createProcessesCreateSummary(vector& startPosition, wait(&temp); } - return exitCommand; + //parent reads in and combine Filter info + for (int i = 0; i < processIDS.size(); i++) { + string tempFilename = fastafile + toString(processIDS[i]) + ".num.temp"; + ifstream in; + m->openInputFile(tempFilename, in); + + int temp, tempNum; + in >> tempNum; m->gobble(in); num += tempNum; + in >> tempNum; m->gobble(in); + for (int k = 0; k < tempNum; k++) { in >> temp; startPosition.push_back(temp); } m->gobble(in); + for (int k = 0; k < tempNum; k++) { in >> temp; endPosition.push_back(temp); } m->gobble(in); + for (int k = 0; k < tempNum; k++) { in >> temp; seqLength.push_back(temp); } m->gobble(in); + for (int k = 0; k < tempNum; k++) { in >> temp; ambigBases.push_back(temp); } m->gobble(in); + for (int k = 0; k < tempNum; k++) { in >> temp; longHomoPolymer.push_back(temp); } m->gobble(in); + + in.close(); + remove(tempFilename.c_str()); + } + + return num; #endif } catch(exception& e) { @@ -445,60 +588,35 @@ int SeqSummaryCommand::createProcessesCreateSummary(vector& startPosition, exit(1); } } -/**************************************************************************************************/ - -int SeqSummaryCommand::setLines(string filename) { +/**********************************************************************************************************************/ +int SeqSummaryCommand::readNames() { try { + //open input file + ifstream in; + m->openInputFile(namefile, in); - vector positions; - - ifstream inFASTA; - openInputFile(filename, inFASTA); + while (!in.eof()) { + if (m->control_pressed) { break; } - string input; - while(!inFASTA.eof()){ - input = getline(inFASTA); - - if (input.length() != 0) { - if(input[0] == '>'){ long int pos = inFASTA.tellg(); positions.push_back(pos - input.length() - 1); } - } - } - inFASTA.close(); - - int numFastaSeqs = positions.size(); - - FILE * pFile; - long size; - - //get num bytes in file - pFile = fopen (filename.c_str(),"rb"); - if (pFile==NULL) perror ("Error opening file"); - else{ - fseek (pFile, 0, SEEK_END); - size=ftell (pFile); - fclose (pFile); + string firstCol, secondCol; + in >> firstCol >> secondCol; m->gobble(in); + + int num = m->getNumNames(secondCol); + + nameMap[firstCol] = num; } + in.close(); - int numSeqsPerProcessor = numFastaSeqs / processors; - - for (int i = 0; i < processors; i++) { - - long int startPos = positions[ i * numSeqsPerProcessor ]; - if(i == processors - 1){ - numSeqsPerProcessor = numFastaSeqs - i * numSeqsPerProcessor; - }else{ - long int myEnd = positions[ (i+1) * numSeqsPerProcessor ]; - } - lines.push_back(new linePair(startPos, numSeqsPerProcessor)); - } + return 0; - return numFastaSeqs; } catch(exception& e) { - m->errorOut(e, "SeqSummaryCommand", "setLines"); + m->errorOut(e, "SeqSummaryCommand", "readNames"); exit(1); } } -//*************************************************************************************************************** + +/**********************************************************************************************************************/ +