chimerapintailcommand.cpp

   1 /*
   2  *  chimerapintailcommand.cpp
   3  *  Mothur
   4  *
   5  *  Created by westcott on 4/1/10.
   6  *  Copyright 2010 Schloss Lab. All rights reserved.
   7  *
   8  */
   9
  10 #include "chimerapintailcommand.h"
  11 #include "pintail.h"
  12
  13 //***************************************************************************************************************
  14
  15 ChimeraPintailCommand::ChimeraPintailCommand(string option)  {
  16         try {
  17                 abort = false;
  18
  19                 //allow user to run help
  20                 if(option == "help") { help(); abort = true; }
  21
  22                 else {
  23                         //valid paramters for this command
  24                         string Array[] =  {"fasta","filter","processors","window" ,"increment","template","conservation","quantile","mask","outputdir","inputdir"};
  25                         vector<string> myArray (Array, Array+(sizeof(Array)/sizeof(string)));
  26
  27                         OptionParser parser(option);
  28                         map<string,string> parameters = parser.getParameters();
  29
  30                         ValidParameters validParameter("chimera.pintail");
  31                         map<string,string>::iterator it;
  32
  33                         //check to make sure all parameters are valid for command
  34                         for (it = parameters.begin(); it != parameters.end(); it++) {
  35                                 if (validParameter.isValidParameter(it->first, myArray, it->second) != true) {  abort = true;  }
  36                         }
  37
  38                         //if the user changes the input directory command factory will send this info to us in the output parameter
  39                         inputDir = validParameter.validFile(parameters, "inputdir", false);
  40                         if (inputDir == "not found"){   inputDir = "";          }
  41                         else {
  42                                 string path;
  43                                 it = parameters.find("template");
  44                                 //user has given a template file
  45                                 if(it != parameters.end()){
  46                                         path = hasPath(it->second);
  47                                         //if the user has not given a path then, add inputdir. else leave path alone.
  48                                         if (path == "") {       parameters["template"] = inputDir + it->second;         }
  49                                 }
  50
  51                                 it = parameters.find("conservation");
  52                                 //user has given a template file
  53                                 if(it != parameters.end()){
  54                                         path = hasPath(it->second);
  55                                         //if the user has not given a path then, add inputdir. else leave path alone.
  56                                         if (path == "") {       parameters["conservation"] = inputDir + it->second;             }
  57                                 }
  58
  59                                 it = parameters.find("quantile");
  60                                 //user has given a template file
  61                                 if(it != parameters.end()){
  62                                         path = hasPath(it->second);
  63                                         //if the user has not given a path then, add inputdir. else leave path alone.
  64                                         if (path == "") {       parameters["quantile"] = inputDir + it->second;         }
  65                                 }
  66                         }
  67
  68
  69                         //check for required parameters
  70                         fastafile = validParameter.validFile(parameters, "fasta", false);
  71                         if (fastafile == "not found") { fastafile = ""; m->mothurOut("fasta is a required parameter for the chimera.pintail command."); m->mothurOutEndLine(); abort = true;  }
  72                         else {
  73                                 splitAtDash(fastafile, fastaFileNames);
  74
  75                                 //go through files and make sure they are good, if not, then disregard them
  76                                 for (int i = 0; i < fastaFileNames.size(); i++) {
  77                                         if (inputDir != "") {
  78                                                 string path = hasPath(fastaFileNames[i]);
  79                                                 //if the user has not given a path then, add inputdir. else leave path alone.
  80                                                 if (path == "") {       fastaFileNames[i] = inputDir + fastaFileNames[i];               }
  81                                         }
  82
  83                                         int ableToOpen;
  84                                         ifstream in;
  85
  86                                         ableToOpen = openInputFile(fastaFileNames[i], in, "noerror");
  87
  88                                         //if you can't open it, try default location
  89                                         if (ableToOpen == 1) {
  90                                                 if (m->getDefaultPath() != "") { //default path is set
  91                                                         string tryPath = m->getDefaultPath() + getSimpleName(fastaFileNames[i]);
  92                                                         m->mothurOut("Unable to open " + fastaFileNames[i] + ". Trying default " + tryPath); m->mothurOutEndLine();
  93                                                         ableToOpen = openInputFile(tryPath, in, "noerror");
  94                                                         fastaFileNames[i] = tryPath;
  95                                                 }
  96                                         }
  97                                         in.close();
  98
  99                                         if (ableToOpen == 1) {
 100                                                 m->mothurOut("Unable to open " + fastaFileNames[i] + ". It will be disregarded."); m->mothurOutEndLine();
 101                                                 //erase from file list
 102                                                 fastaFileNames.erase(fastaFileNames.begin()+i);
 103                                                 i--;
 104                                         }
 105                                 }
 106
 107                                 //make sure there is at least one valid file left
 108                                 if (fastaFileNames.size() == 0) { m->mothurOut("no valid files."); m->mothurOutEndLine(); abort = true; }
 109                         }
 110
 111                         string temp;
 112                         temp = validParameter.validFile(parameters, "filter", false);                   if (temp == "not found") { temp = "F"; }
 113                         filter = isTrue(temp);
 114
 115                         temp = validParameter.validFile(parameters, "processors", false);               if (temp == "not found") { temp = "1"; }
 116                         convert(temp, processors);
 117
 118                         temp = validParameter.validFile(parameters, "window", false);                   if (temp == "not found") { temp = "0"; }
 119                         convert(temp, window);
 120
 121                         temp = validParameter.validFile(parameters, "increment", false);                if (temp == "not found") { temp = "25"; }
 122                         convert(temp, increment);
 123
 124                         maskfile = validParameter.validFile(parameters, "mask", false);
 125                         if (maskfile == "not found") { maskfile = "";  }
 126                         else if (maskfile != "default")  {
 127                                 if (inputDir != "") {
 128                                         string path = hasPath(maskfile);
 129                                         //if the user has not given a path then, add inputdir. else leave path alone.
 130                                         if (path == "") {       maskfile = inputDir + maskfile;         }
 131                                 }
 132
 133                                 ifstream in;
 134                                 int     ableToOpen = openInputFile(maskfile, in);
 135                                 if (ableToOpen == 1) { abort = true; }
 136                                 in.close();
 137                         }
 138
 139
 140                         //if the user changes the output directory command factory will send this info to us in the output parameter
 141                         outputDir = validParameter.validFile(parameters, "outputdir", false);           if (outputDir == "not found"){
 142                                 outputDir = "";
 143                                 outputDir += hasPath(fastafile); //if user entered a file with a path then preserve it
 144                         }
 145
 146                         templatefile = validParameter.validFile(parameters, "template", true);
 147                         if (templatefile == "not open") { abort = true; }
 148                         else if (templatefile == "not found") { templatefile = "";  m->mothurOut("template is a required parameter for the chimera.pintail command."); m->mothurOutEndLine(); abort = true;  }
 149
 150                         consfile = validParameter.validFile(parameters, "conservation", true);
 151                         if (consfile == "not open") { abort = true; }
 152                         else if (consfile == "not found") {
 153                                 consfile = "";
 154                                 //check for consfile
 155                                 string tempConsFile = getRootName(inputDir + getSimpleName(templatefile)) + "freq";
 156                                 ifstream FileTest(tempConsFile.c_str());
 157                                 if(FileTest){
 158                                         bool GoodFile = checkReleaseVersion(FileTest, m->getVersion());
 159                                         if (GoodFile) {
 160                                                 m->mothurOut("I found " + tempConsFile + " in your input file directory. I will use it to save time."); m->mothurOutEndLine();  consfile = tempConsFile;  FileTest.close();
 161                                         }
 162                                 }
 163                         }
 164
 165                         quanfile = validParameter.validFile(parameters, "quantile", true);
 166                         if (quanfile == "not open") { abort = true; }
 167                         else if (quanfile == "not found") { quanfile = ""; }
 168                 }
 169         }
 170         catch(exception& e) {
 171                 m->errorOut(e, "ChimeraPintailCommand", "ChimeraPintailCommand");
 172                 exit(1);
 173         }
 174 }
 175 //**********************************************************************************************************************
 176
 177 void ChimeraPintailCommand::help(){
 178         try {
 179
 180                 m->mothurOut("The chimera.pintail command reads a fastafile and templatefile and outputs potentially chimeric sequences.\n");
 181                 m->mothurOut("This command was created using the algorythms described in the 'At Least 1 in 20 16S rRNA Sequence Records Currently Held in the Public Repositories is Estimated To Contain Substantial Anomalies' paper by Kevin E. Ashelford 1, Nadia A. Chuzhanova 3, John C. Fry 1, Antonia J. Jones 2 and Andrew J. Weightman 1.\n");
 182                 m->mothurOut("The chimera.pintail command parameters are fasta, template, filter, mask, processors, window, increment, conservation and quantile.\n");
 183                 m->mothurOut("The fasta parameter allows you to enter the fasta file containing your potentially chimeric sequences, and is required. \n");
 184                 m->mothurOut("You may enter multiple fasta files by separating their names with dashes. ie. fasta=abrecovery.fasta-amzon.fasta \n");
 185                 m->mothurOut("The template parameter allows you to enter a template file containing known non-chimeric sequences, and is required. \n");
 186                 m->mothurOut("The filter parameter allows you to specify if you would like to apply a vertical and 50% soft filter. \n");
 187                 m->mothurOut("The mask parameter allows you to specify a file containing one sequence you wish to use as a mask for the your sequences, by default no mask is applied.  You can apply an ecoli mask by typing, mask=default. \n");
 188                 m->mothurOut("The processors parameter allows you to specify how many processors you would like to use.  The default is 1. \n");
 189                 #ifdef USE_MPI
 190                 m->mothurOut("When using MPI, the processors parameter is set to the number of MPI processes running. \n");
 191                 #endif
 192                 m->mothurOut("The window parameter allows you to specify the window size for searching for chimeras, default=300. \n");
 193                 m->mothurOut("The increment parameter allows you to specify how far you move each window while finding chimeric sequences, default=25.\n");
 194                 m->mothurOut("The conservation parameter allows you to enter a frequency file containing the highest bases frequency at each place in the alignment.\n");
 195                 m->mothurOut("The quantile parameter allows you to enter a file containing quantiles for a template files sequences, if you use the filter the quantile file generated becomes unique to the fasta file you used.\n");
 196                 m->mothurOut("The chimera.pintail command should be in the following format: \n");
 197                 m->mothurOut("chimera.pintail(fasta=yourFastaFile, template=yourTemplate) \n");
 198                 m->mothurOut("Example: chimera.pintail(fasta=AD.align, template=silva.bacteria.fasta) \n");
 199                 m->mothurOut("Note: No spaces between parameter labels (i.e. fasta), '=' and parameters (i.e.yourFastaFile).\n\n");
 200         }
 201         catch(exception& e) {
 202                 m->errorOut(e, "ChimeraPintailCommand", "help");
 203                 exit(1);
 204         }
 205 }
 206
 207 //***************************************************************************************************************
 208
 209 ChimeraPintailCommand::~ChimeraPintailCommand(){        /*      do nothing      */      }
 210
 211 //***************************************************************************************************************
 212
 213 int ChimeraPintailCommand::execute(){
 214         try{
 215
 216                 if (abort == true) { return 0; }
 217
 218                 for (int s = 0; s < fastaFileNames.size(); s++) {
 219
 220                         m->mothurOut("Checking sequences from " + fastaFileNames[s] + " ..." ); m->mothurOutEndLine();
 221
 222                         int start = time(NULL);
 223
 224                         //set user options
 225                         if (maskfile == "default") { m->mothurOut("I am using the default 236627 EU009184.1 Shigella dysenteriae str. FBD013."); m->mothurOutEndLine();  }
 226
 227                         //check for quantile to save the time
 228                         string tempQuan = "";
 229                         if ((!filter) && (maskfile == "")) {
 230                                 tempQuan = inputDir + getRootName(getSimpleName(templatefile)) + "pintail.quan";
 231                         }else if ((!filter) && (maskfile != "")) {
 232                                 tempQuan = inputDir + getRootName(getSimpleName(templatefile)) + "pintail.masked.quan";
 233                         }else if ((filter) && (maskfile != "")) {
 234                                 tempQuan = inputDir + getRootName(getSimpleName(templatefile)) + "pintail.filtered." + getSimpleName(getRootName(fastaFileNames[s])) + "masked.quan";
 235                         }else if ((filter) && (maskfile == "")) {
 236                                 tempQuan = inputDir + getRootName(getSimpleName(templatefile)) + "pintail.filtered." + getSimpleName(getRootName(fastaFileNames[s])) + "quan";
 237                         }
 238
 239                         ifstream FileTest(tempQuan.c_str());
 240                         if(FileTest){
 241                                 bool GoodFile = checkReleaseVersion(FileTest, m->getVersion());
 242                                 if (GoodFile) {
 243                                         m->mothurOut("I found " + tempQuan + " in your input file directory. I will use it to save time."); m->mothurOutEndLine();  quanfile = tempQuan;  FileTest.close();
 244                                 }
 245                         }
 246
 247                         chimera = new Pintail(fastaFileNames[s], templatefile, filter, processors, maskfile, consfile, quanfile, window, increment, outputDir);
 248
 249                         string outputFileName, accnosFileName;
 250                         if (maskfile != "") {
 251                                 outputFileName = outputDir + getRootName(getSimpleName(fastaFileNames[s])) + maskfile + ".pintail.chimeras";
 252                                 accnosFileName = outputDir + getRootName(getSimpleName(fastaFileNames[s])) + maskfile + ".pintail.accnos";
 253                         }else {
 254                                 outputFileName = outputDir + getRootName(getSimpleName(fastaFileNames[s]))  + "pintail.chimeras";
 255                                 accnosFileName = outputDir + getRootName(getSimpleName(fastaFileNames[s]))  + "pintail.accnos";
 256                         }
 257
 258                         if (m->control_pressed) { delete chimera; for (int j = 0; j < outputNames.size(); j++) {        remove(outputNames[j].c_str()); }  return 0;    }
 259
 260                         if (chimera->getUnaligned()) {
 261                                 m->mothurOut("Your template sequences are different lengths, please correct."); m->mothurOutEndLine();
 262                                 delete chimera;
 263                                 return 0;
 264                         }
 265                         templateSeqsLength = chimera->getLength();
 266
 267                 #ifdef USE_MPI
 268                         int pid, end, numSeqsPerProcessor;
 269                                 int tag = 2001;
 270                                 vector<unsigned long int> MPIPos;
 271
 272                                 MPI_Status status;
 273                                 MPI_Comm_rank(MPI_COMM_WORLD, &pid); //find out who we are
 274                                 MPI_Comm_size(MPI_COMM_WORLD, &processors);
 275
 276                                 MPI_File inMPI;
 277                                 MPI_File outMPI;
 278                                 MPI_File outMPIAccnos;
 279
 280                                 int outMode=MPI_MODE_CREATE|MPI_MODE_WRONLY;
 281                                 int inMode=MPI_MODE_RDONLY;
 282
 283                                 char outFilename[1024];
 284                                 strcpy(outFilename, outputFileName.c_str());
 285
 286                                 char outAccnosFilename[1024];
 287                                 strcpy(outAccnosFilename, accnosFileName.c_str());
 288
 289                                 char inFileName[1024];
 290                                 strcpy(inFileName, fastaFileNames[s].c_str());
 291
 292                                 MPI_File_open(MPI_COMM_WORLD, inFileName, inMode, MPI_INFO_NULL, &inMPI);  //comm, filename, mode, info, filepointer
 293                                 MPI_File_open(MPI_COMM_WORLD, outFilename, outMode, MPI_INFO_NULL, &outMPI);
 294                                 MPI_File_open(MPI_COMM_WORLD, outAccnosFilename, outMode, MPI_INFO_NULL, &outMPIAccnos);
 295
 296                                 if (m->control_pressed) {  MPI_File_close(&inMPI);  MPI_File_close(&outMPI);   MPI_File_close(&outMPIAccnos);  for (int j = 0; j < outputNames.size(); j++) {   remove(outputNames[j].c_str()); }  delete chimera; return 0;  }
 297
 298                                 if (pid == 0) { //you are the root process
 299
 300                                         MPIPos = setFilePosFasta(fastaFileNames[s], numSeqs); //fills MPIPos, returns numSeqs
 301
 302                                         //send file positions to all processes
 303                                         for(int i = 1; i < processors; i++) {
 304                                                 MPI_Send(&numSeqs, 1, MPI_INT, i, tag, MPI_COMM_WORLD);
 305                                                 MPI_Send(&MPIPos[0], (numSeqs+1), MPI_LONG, i, tag, MPI_COMM_WORLD);
 306                                         }
 307
 308                                         //figure out how many sequences you have to align
 309                                         numSeqsPerProcessor = numSeqs / processors;
 310                                         int startIndex =  pid * numSeqsPerProcessor;
 311                                         if(pid == (processors - 1)){    numSeqsPerProcessor = numSeqs - pid * numSeqsPerProcessor;      }
 312
 313                                         //do your part
 314                                         driverMPI(startIndex, numSeqsPerProcessor, inMPI, outMPI, outMPIAccnos, MPIPos);
 315
 316                                         if (m->control_pressed) {  MPI_File_close(&inMPI);  MPI_File_close(&outMPI);   MPI_File_close(&outMPIAccnos);  remove(outputFileName.c_str());  remove(accnosFileName.c_str());  for (int j = 0; j < outputNames.size(); j++) { remove(outputNames[j].c_str()); }  delete chimera; return 0;  }
 317
 318                                 }else{ //you are a child process
 319                                         MPI_Recv(&numSeqs, 1, MPI_INT, 0, tag, MPI_COMM_WORLD, &status);
 320                                         MPIPos.resize(numSeqs+1);
 321                                         MPI_Recv(&MPIPos[0], (numSeqs+1), MPI_LONG, 0, tag, MPI_COMM_WORLD, &status);
 322
 323                                         //figure out how many sequences you have to align
 324                                         numSeqsPerProcessor = numSeqs / processors;
 325                                         int startIndex =  pid * numSeqsPerProcessor;
 326                                         if(pid == (processors - 1)){    numSeqsPerProcessor = numSeqs - pid * numSeqsPerProcessor;      }
 327
 328                                         //do your part
 329                                         driverMPI(startIndex, numSeqsPerProcessor, inMPI, outMPI, outMPIAccnos, MPIPos);
 330
 331                                         if (m->control_pressed) {  MPI_File_close(&inMPI);  MPI_File_close(&outMPI);   MPI_File_close(&outMPIAccnos);  for (int j = 0; j < outputNames.size(); j++) {   remove(outputNames[j].c_str()); }  delete chimera; return 0;  }
 332                                 }
 333
 334                                 //close files
 335                                 MPI_File_close(&inMPI);
 336                                 MPI_File_close(&outMPI);
 337                                 MPI_File_close(&outMPIAccnos);
 338                                 MPI_Barrier(MPI_COMM_WORLD); //make everyone wait - just in case
 339                 #else
 340
 341                         //break up file
 342                         #if defined (__APPLE__) || (__MACH__) || (linux) || (__linux)
 343                                 if(processors == 1){
 344                                         ifstream inFASTA;
 345                                         openInputFile(fastaFileNames[s], inFASTA);
 346                                         getNumSeqs(inFASTA, numSeqs);
 347                                         inFASTA.close();
 348
 349                                         lines.push_back(new linePair(0, numSeqs));
 350
 351                                         driver(lines[0], outputFileName, fastaFileNames[s], accnosFileName);
 352
 353                                         if (m->control_pressed) {
 354                                                 remove(outputFileName.c_str());
 355                                                 remove(accnosFileName.c_str());
 356                                                 for (int j = 0; j < outputNames.size(); j++) {  remove(outputNames[j].c_str()); }
 357                                                 for (int i = 0; i < lines.size(); i++) {  delete lines[i];  }  lines.clear();
 358                                                 delete chimera;
 359                                                 return 0;
 360                                         }
 361
 362                                 }else{
 363                                         vector<unsigned long int> positions;
 364                                         processIDS.resize(0);
 365
 366                                         ifstream inFASTA;
 367                                         openInputFile(fastaFileNames[s], inFASTA);
 368
 369                                         string input;
 370                                         while(!inFASTA.eof()){
 371                                                 input = getline(inFASTA);
 372                                                 if (input.length() != 0) {
 373                                                         if(input[0] == '>'){    unsigned long int pos = inFASTA.tellg(); positions.push_back(pos - input.length() - 1); }
 374                                                 }
 375                                         }
 376                                         inFASTA.close();
 377
 378                                         numSeqs = positions.size();
 379
 380                                         int numSeqsPerProcessor = numSeqs / processors;
 381
 382                                         for (int i = 0; i < processors; i++) {
 383                                                 unsigned long int startPos = positions[ i * numSeqsPerProcessor ];
 384                                                 if(i == processors - 1){
 385                                                         numSeqsPerProcessor = numSeqs - i * numSeqsPerProcessor;
 386                                                 }
 387                                                 lines.push_back(new linePair(startPos, numSeqsPerProcessor));
 388                                         }
 389
 390                                         createProcesses(outputFileName, fastaFileNames[s], accnosFileName);
 391
 392                                         rename((outputFileName + toString(processIDS[0]) + ".temp").c_str(), outputFileName.c_str());
 393                                         rename((accnosFileName + toString(processIDS[0]) + ".temp").c_str(), accnosFileName.c_str());
 394
 395                                         //append output files
 396                                         for(int i=1;i<processors;i++){
 397                                                 appendFiles((outputFileName + toString(processIDS[i]) + ".temp"), outputFileName);
 398                                                 remove((outputFileName + toString(processIDS[i]) + ".temp").c_str());
 399                                         }
 400
 401                                         //append output files
 402                                         for(int i=1;i<processors;i++){
 403                                                 appendFiles((accnosFileName + toString(processIDS[i]) + ".temp"), accnosFileName);
 404                                                 remove((accnosFileName + toString(processIDS[i]) + ".temp").c_str());
 405                                         }
 406
 407                                         if (m->control_pressed) {
 408                                                 remove(outputFileName.c_str());
 409                                                 remove(accnosFileName.c_str());
 410                                                 for (int j = 0; j < outputNames.size(); j++) {  remove(outputNames[j].c_str()); }
 411                                                 for (int i = 0; i < lines.size(); i++) {  delete lines[i];  }  lines.clear();
 412                                                 delete chimera;
 413                                                 return 0;
 414                                         }
 415                                 }
 416
 417                         #else
 418                                 ifstream inFASTA;
 419                                 openInputFile(fastaFileNames[s], inFASTA);
 420                                 getNumSeqs(inFASTA, numSeqs);
 421                                 inFASTA.close();
 422                                 lines.push_back(new linePair(0, numSeqs));
 423
 424                                 driver(lines[0], outputFileName, fastaFileNames[s], accnosFileName);
 425
 426                                 if (m->control_pressed) {
 427                                                 remove(outputFileName.c_str());
 428                                                 remove(accnosFileName.c_str());
 429                                                 for (int j = 0; j < outputNames.size(); j++) {  remove(outputNames[j].c_str()); }
 430                                                 for (int i = 0; i < lines.size(); i++) {  delete lines[i];  }  lines.clear();
 431                                                 delete chimera;
 432                                                 return 0;
 433                                 }
 434                         #endif
 435
 436                 #endif
 437
 438                         delete chimera;
 439                         for (int i = 0; i < lines.size(); i++) {  delete lines[i];  }  lines.clear();
 440
 441                         outputNames.push_back(outputFileName);
 442                         outputNames.push_back(accnosFileName);
 443
 444                         m->mothurOutEndLine();
 445                         m->mothurOutEndLine(); m->mothurOut("It took " + toString(time(NULL) - start) + " secs to check " + toString(numSeqs) + " sequences."); m->mothurOutEndLine();
 446                 }
 447
 448                 m->mothurOutEndLine();
 449                 m->mothurOut("Output File Names: "); m->mothurOutEndLine();
 450                 for (int i = 0; i < outputNames.size(); i++) {  m->mothurOut(outputNames[i]); m->mothurOutEndLine();    }
 451                 m->mothurOutEndLine();
 452
 453                 return 0;
 454
 455         }
 456         catch(exception& e) {
 457                 m->errorOut(e, "ChimeraPintailCommand", "execute");
 458                 exit(1);
 459         }
 460 }
 461 //**********************************************************************************************************************
 462
 463 int ChimeraPintailCommand::driver(linePair* line, string outputFName, string filename, string accnos){
 464         try {
 465                 ofstream out;
 466                 openOutputFile(outputFName, out);
 467
 468                 ofstream out2;
 469                 openOutputFile(accnos, out2);
 470
 471                 ifstream inFASTA;
 472                 openInputFile(filename, inFASTA);
 473
 474                 inFASTA.seekg(line->start);
 475
 476                 for(int i=0;i<line->numSeqs;i++){
 477
 478                         if (m->control_pressed) {       return 1;       }
 479
 480                         Sequence* candidateSeq = new Sequence(inFASTA);  gobble(inFASTA);
 481
 482                         if (candidateSeq->getName() != "") { //incase there is a commented sequence at the end of a file
 483
 484                                 if (candidateSeq->getAligned().length() != templateSeqsLength)  {  //chimeracheck does not require seqs to be aligned
 485                                         m->mothurOut(candidateSeq->getName() + " is not the same length as the template sequences. Skipping."); m->mothurOutEndLine();
 486                                 }else{
 487                                         //find chimeras
 488                                         chimera->getChimeras(candidateSeq);
 489
 490                                         if (m->control_pressed) {       delete candidateSeq; return 1;  }
 491
 492                                         //print results
 493                                         chimera->print(out, out2);
 494                                 }
 495                         }
 496                         delete candidateSeq;
 497
 498                         //report progress
 499                         if((i+1) % 100 == 0){   m->mothurOut("Processing sequence: " + toString(i+1)); m->mothurOutEndLine();           }
 500                 }
 501                 //report progress
 502                 if((line->numSeqs) % 100 != 0){ m->mothurOut("Processing sequence: " + toString(line->numSeqs)); m->mothurOutEndLine();         }
 503
 504                 out.close();
 505                 out2.close();
 506                 inFASTA.close();
 507
 508                 return 0;
 509         }
 510         catch(exception& e) {
 511                 m->errorOut(e, "ChimeraPintailCommand", "driver");
 512                 exit(1);
 513         }
 514 }
 515 //**********************************************************************************************************************
 516 #ifdef USE_MPI
 517 int ChimeraPintailCommand::driverMPI(int start, int num, MPI_File& inMPI, MPI_File& outMPI, MPI_File& outAccMPI, vector<unsigned long int>& MPIPos){
 518         try {
 519
 520                 MPI_Status status;
 521                 int pid;
 522                 MPI_Comm_rank(MPI_COMM_WORLD, &pid); //find out who we are
 523
 524                 for(int i=0;i<num;i++){
 525
 526                         if (m->control_pressed) {       return 1;       }
 527
 528                         //read next sequence
 529                         int length = MPIPos[start+i+1] - MPIPos[start+i];
 530
 531                         char* buf4 = new char[length];
 532                         MPI_File_read_at(inMPI, MPIPos[start+i], buf4, length, MPI_CHAR, &status);
 533
 534                         string tempBuf = buf4;
 535                         if (tempBuf.length() > length) { tempBuf = tempBuf.substr(0, length);  }
 536                         istringstream iss (tempBuf,istringstream::in);
 537                         delete buf4;
 538
 539                         Sequence* candidateSeq = new Sequence(iss);  gobble(iss);
 540
 541                         if (candidateSeq->getName() != "") { //incase there is a commented sequence at the end of a file
 542
 543                                 if      (candidateSeq->getAligned().length() != templateSeqsLength) {  //chimeracheck does not require seqs to be aligned
 544                                         m->mothurOut(candidateSeq->getName() + " is not the same length as the template sequences. Skipping."); m->mothurOutEndLine();
 545                                 }else{
 546                                         //find chimeras
 547                                         chimera->getChimeras(candidateSeq);
 548
 549                                         if (m->control_pressed) {       delete candidateSeq; return 1;  }
 550
 551                                         //print results
 552                                         bool isChimeric = chimera->print(outMPI, outAccMPI);
 553                                 }
 554                         }
 555                         delete candidateSeq;
 556
 557                         //report progress
 558                         if((i+1) % 100 == 0){  cout << "Processing sequence: " << (i+1) << endl;        m->mothurOutJustToLog("Processing sequence: " + toString(i+1) + "\n");          }
 559                 }
 560                 //report progress
 561                 if(num % 100 != 0){             cout << "Processing sequence: " << num << endl; m->mothurOutJustToLog("Processing sequence: " + toString(num) + "\n");  }
 562
 563
 564                 return 0;
 565         }
 566         catch(exception& e) {
 567                 m->errorOut(e, "ChimeraPintailCommand", "driverMPI");
 568                 exit(1);
 569         }
 570 }
 571 #endif
 572
 573 /**************************************************************************************************/
 574
 575 int ChimeraPintailCommand::createProcesses(string outputFileName, string filename, string accnos) {
 576         try {
 577 #if defined (__APPLE__) || (__MACH__) || (linux) || (__linux)
 578                 int process = 0;
 579                 //              processIDS.resize(0);
 580
 581                 //loop through and create all the processes you want
 582                 while (process != processors) {
 583                         int pid = fork();
 584
 585                         if (pid > 0) {
 586                                 processIDS.push_back(pid);  //create map from line number to pid so you can append files in correct order later
 587                                 process++;
 588                         }else if (pid == 0){
 589                                 driver(lines[process], outputFileName + toString(getpid()) + ".temp", filename, accnos + toString(getpid()) + ".temp");
 590                                 exit(0);
 591                         }else { m->mothurOut("unable to spawn the necessary processes."); m->mothurOutEndLine(); exit(0); }
 592                 }
 593
 594                 //force parent to wait until all the processes are done
 595                 for (int i=0;i<processors;i++) {
 596                         int temp = processIDS[i];
 597                         wait(&temp);
 598                 }
 599
 600                 return 0;
 601 #endif
 602         }
 603         catch(exception& e) {
 604                 m->errorOut(e, "ChimeraPintailCommand", "createProcesses");
 605                 exit(1);
 606         }
 607 }
 608
 609 /**************************************************************************************************/
 610
 611