-void ScreenSeqsCommand::help(){
- try {
- cout << "The screen.seqs command reads a fastafile and creates ....." << "\n";
- cout << "The screen.seqs command parameters are fasta, start, end, maxambig, maxhomop, minlength, maxlength, name, and group." << "\n";
- cout << "The fasta parameter is required." << "\n";
- cout << "The start parameter .... The default is -1." << "\n";
- cout << "The end parameter .... The default is -1." << "\n";
- cout << "The maxambig parameter .... The default is -1." << "\n";
- cout << "The maxhomop parameter .... The default is -1." << "\n";
- cout << "The minlength parameter .... The default is -1." << "\n";
- cout << "The maxlength parameter .... The default is -1." << "\n";
- cout << "The name parameter allows you to provide a namesfile, and the group parameter allows you to provide a groupfile." << "\n";
- cout << "The screen.seqs command should be in the following format: " << "\n";
- cout << "screen.seqs(fasta=yourFastaFile, name=youNameFile, group=yourGroupFIle, start=yourStart, end=yourEnd, maxambig=yourMaxambig, " << "\n";
- cout << "maxhomop=yourMaxhomop, minlength=youMinlength, maxlength=yourMaxlength) " << "\n";
- cout << "Example screen.seqs(fasta=abrecovery.fasta, name=abrecovery.names, group=abrecovery.groups, start=..., end=..., maxambig=..., maxhomop=..., minlength=..., maxlength=...)." << "\n";
- cout << "Note: No spaces between parameter labels (i.e. fasta), '=' and parameters (i.e.yourFasta)." << "\n" << "\n";
+int ScreenSeqsCommand::execute(){
+ try{
+
+ if (abort == true) { if (calledHelp) { return 0; } return 2; }
+
+ //if the user want to optimize we need to know the 90% mark
+ vector<unsigned long long> positions;
+ if (optimize.size() != 0) { //get summary is paralellized so we need to divideFile, no need to do this step twice so I moved it here
+ //use the namefile to optimize correctly
+ if (namefile != "") { nameMap = m->readNames(namefile); }
+ getSummary(positions);
+ }
+ else {
+ #if defined (__APPLE__) || (__MACH__) || (linux) || (__linux)
+ positions = m->divideFile(fastafile, processors);
+ for (int i = 0; i < (positions.size()-1); i++) { lines.push_back(linePair(positions[i], positions[(i+1)])); }
+ #else
+ if(processors == 1){ lines.push_back(linePair(0, 1000)); }
+ else {
+ int numFastaSeqs = 0;
+ positions = m->setFilePosFasta(fastafile, numFastaSeqs);
+
+ //figure out how many sequences you have to process
+ int numSeqsPerProcessor = numFastaSeqs / processors;
+ for (int i = 0; i < processors; i++) {
+ int startIndex = i * numSeqsPerProcessor;
+ if(i == (processors - 1)){ numSeqsPerProcessor = numFastaSeqs - i * numSeqsPerProcessor; }
+ lines.push_back(linePair(positions[startIndex], numSeqsPerProcessor));
+ }
+ }
+ #endif
+ }
+
+ string goodSeqFile = outputDir + m->getRootName(m->getSimpleName(fastafile)) + "good" + m->getExtension(fastafile);
+ string badAccnosFile = outputDir + m->getRootName(m->getSimpleName(fastafile)) + "bad.accnos";
+
+ int numFastaSeqs = 0;
+ set<string> badSeqNames;
+ int start = time(NULL);
+
+#ifdef USE_MPI
+ int pid, numSeqsPerProcessor;
+ int tag = 2001;
+ vector<unsigned long long> MPIPos;
+
+ MPI_Status status;
+ MPI_Comm_rank(MPI_COMM_WORLD, &pid); //find out who we are
+ MPI_Comm_size(MPI_COMM_WORLD, &processors);
+
+ MPI_File inMPI;
+ MPI_File outMPIGood;
+ MPI_File outMPIBadAccnos;
+
+ int outMode=MPI_MODE_CREATE|MPI_MODE_WRONLY;
+ int inMode=MPI_MODE_RDONLY;
+
+ char outGoodFilename[1024];
+ strcpy(outGoodFilename, goodSeqFile.c_str());
+
+ char outBadAccnosFilename[1024];
+ strcpy(outBadAccnosFilename, badAccnosFile.c_str());
+
+ char inFileName[1024];
+ strcpy(inFileName, fastafile.c_str());
+
+ MPI_File_open(MPI_COMM_WORLD, inFileName, inMode, MPI_INFO_NULL, &inMPI); //comm, filename, mode, info, filepointer
+ MPI_File_open(MPI_COMM_WORLD, outGoodFilename, outMode, MPI_INFO_NULL, &outMPIGood);
+ MPI_File_open(MPI_COMM_WORLD, outBadAccnosFilename, outMode, MPI_INFO_NULL, &outMPIBadAccnos);
+
+ if (m->control_pressed) { MPI_File_close(&inMPI); MPI_File_close(&outMPIGood); MPI_File_close(&outMPIBadAccnos); return 0; }
+
+ if (pid == 0) { //you are the root process
+
+ MPIPos = m->setFilePosFasta(fastafile, numFastaSeqs); //fills MPIPos, returns numSeqs
+
+ //send file positions to all processes
+ for(int i = 1; i < processors; i++) {
+ MPI_Send(&numFastaSeqs, 1, MPI_INT, i, tag, MPI_COMM_WORLD);
+ MPI_Send(&MPIPos[0], (numFastaSeqs+1), MPI_LONG, i, tag, MPI_COMM_WORLD);
+ }
+
+ //figure out how many sequences you have to align
+ numSeqsPerProcessor = numFastaSeqs / processors;
+ int startIndex = pid * numSeqsPerProcessor;
+ if(pid == (processors - 1)){ numSeqsPerProcessor = numFastaSeqs - pid * numSeqsPerProcessor; }
+
+ //align your part
+ driverMPI(startIndex, numSeqsPerProcessor, inMPI, outMPIGood, outMPIBadAccnos, MPIPos, badSeqNames);
+
+ if (m->control_pressed) { MPI_File_close(&inMPI); MPI_File_close(&outMPIGood); MPI_File_close(&outMPIBadAccnos); return 0; }
+
+ for (int i = 1; i < processors; i++) {
+ //get bad lists
+ int badSize;
+ MPI_Recv(&badSize, 1, MPI_INT, i, tag, MPI_COMM_WORLD, &status);
+ }
+ }else{ //you are a child process
+ MPI_Recv(&numFastaSeqs, 1, MPI_INT, 0, tag, MPI_COMM_WORLD, &status);
+ MPIPos.resize(numFastaSeqs+1);
+ MPI_Recv(&MPIPos[0], (numFastaSeqs+1), MPI_LONG, 0, tag, MPI_COMM_WORLD, &status);
+
+ //figure out how many sequences you have to align
+ numSeqsPerProcessor = numFastaSeqs / processors;
+ int startIndex = pid * numSeqsPerProcessor;
+ if(pid == (processors - 1)){ numSeqsPerProcessor = numFastaSeqs - pid * numSeqsPerProcessor; }
+
+ //align your part
+ driverMPI(startIndex, numSeqsPerProcessor, inMPI, outMPIGood, outMPIBadAccnos, MPIPos, badSeqNames);
+
+ if (m->control_pressed) { MPI_File_close(&inMPI); MPI_File_close(&outMPIGood); MPI_File_close(&outMPIBadAccnos); return 0; }
+
+ //send bad list
+ int badSize = badSeqNames.size();
+ MPI_Send(&badSize, 1, MPI_INT, 0, tag, MPI_COMM_WORLD);
+ }
+
+ //close files
+ MPI_File_close(&inMPI);
+ MPI_File_close(&outMPIGood);
+ MPI_File_close(&outMPIBadAccnos);
+ MPI_Barrier(MPI_COMM_WORLD); //make everyone wait - just in case
+
+#else
+
+ //#if defined (__APPLE__) || (__MACH__) || (linux) || (__linux)
+ if(processors == 1){
+ numFastaSeqs = driver(lines[0], goodSeqFile, badAccnosFile, fastafile, badSeqNames);
+ }else{
+ numFastaSeqs = createProcesses(goodSeqFile, badAccnosFile, fastafile, badSeqNames);
+ }
+ //#else
+ // numFastaSeqs = driver(lines[0], goodSeqFile, badAccnosFile, fastafile, badSeqNames);
+ //#endif
+ if (m->control_pressed) { m->mothurRemove(goodSeqFile); return 0; }
+#endif
+
+ #ifdef USE_MPI
+ MPI_Comm_rank(MPI_COMM_WORLD, &pid);
+
+ if (pid == 0) { //only one process should fix files
+
+ //read accnos file with all names in it, process 0 just has its names
+ MPI_File inMPIAccnos;
+ MPI_Offset size;
+
+ char inFileName[1024];
+ strcpy(inFileName, badAccnosFile.c_str());
+
+ MPI_File_open(MPI_COMM_SELF, inFileName, inMode, MPI_INFO_NULL, &inMPIAccnos); //comm, filename, mode, info, filepointer
+ MPI_File_get_size(inMPIAccnos, &size);
+
+ char* buffer = new char[size];
+ MPI_File_read(inMPIAccnos, buffer, size, MPI_CHAR, &status);
+
+ string tempBuf = buffer;
+ if (tempBuf.length() > size) { tempBuf = tempBuf.substr(0, size); }
+ istringstream iss (tempBuf,istringstream::in);
+
+ delete buffer;
+ MPI_File_close(&inMPIAccnos);
+
+ badSeqNames.clear();
+ string tempName;
+ while (!iss.eof()) {
+ iss >> tempName; m->gobble(iss);
+ badSeqNames.insert(tempName);
+ }
+ #endif
+
+ if(namefile != "" && groupfile != "") {
+ screenNameGroupFile(badSeqNames);
+ if (m->control_pressed) { m->mothurRemove(goodSeqFile); return 0; }
+ }else if(namefile != "") {
+ screenNameGroupFile(badSeqNames);
+ if (m->control_pressed) { m->mothurRemove(goodSeqFile); return 0; }
+ }else if(groupfile != "") { screenGroupFile(badSeqNames); } // this screens just the group
+
+ if (m->control_pressed) { m->mothurRemove(goodSeqFile); return 0; }