pipelinepdscommand.cpp

   1 /*
   2  *  pipelinepdscommand.cpp
   3  *  Mothur
   4  *
   5  *  Created by westcott on 10/5/10.
   6  *  Copyright 2010 Schloss Lab. All rights reserved.
   7  *
   8  */
   9
  10 #include "pipelinepdscommand.h"
  11 #include "sffinfocommand.h"
  12 #include "commandoptionparser.hpp"
  13
  14 //**********************************************************************************************************************
  15 vector<string> PipelineCommand::getValidParameters(){
  16         try {
  17                 string Array[] =  {"sff","oligos","align","chimera","classify","taxonomy","pipeline","processors","outputdir","inputdir"};
  18                 vector<string> myArray (Array, Array+(sizeof(Array)/sizeof(string)));
  19                 return myArray;
  20         }
  21         catch(exception& e) {
  22                 m->errorOut(e, "PipelineCommand", "getValidParameters");
  23                 exit(1);
  24         }
  25 }
  26 //**********************************************************************************************************************
  27 vector<string> PipelineCommand::getRequiredParameters(){
  28         try {
  29                 vector<string> myArray;
  30                 return myArray;
  31         }
  32         catch(exception& e) {
  33                 m->errorOut(e, "PipelineCommand", "getRequiredParameters");
  34                 exit(1);
  35         }
  36 }
  37 //**********************************************************************************************************************
  38 vector<string> PipelineCommand::getRequiredFiles(){
  39         try {
  40                 vector<string> myArray;
  41                 return myArray;
  42         }
  43         catch(exception& e) {
  44                 m->errorOut(e, "PipelineCommand", "getRequiredFiles");
  45                 exit(1);
  46         }
  47 }
  48 //**********************************************************************************************************************
  49 PipelineCommand::PipelineCommand(string option) {
  50         try {
  51                 cFactory = CommandFactory::getInstance();
  52                 abort = false; calledHelp = false;
  53
  54                 //allow user to run help
  55                 if(option == "help") { help(); abort = true; calledHelp = true; }
  56
  57                 else {
  58
  59                         //valid paramters for this command
  60                         string AlignArray[] =  {"sff","oligos","align","chimera","classify","taxonomy","pipeline","processors","outputdir","inputdir"};
  61                         vector<string> myArray (AlignArray, AlignArray+(sizeof(AlignArray)/sizeof(string)));
  62
  63                         OptionParser parser(option);
  64                         map<string, string> parameters = parser.getParameters();
  65
  66                         ValidParameters validParameter;
  67                         map<string, string>::iterator it;
  68
  69                         //check to make sure all parameters are valid for command
  70                         for (it = parameters.begin(); it != parameters.end(); it++) {
  71                                 if (validParameter.isValidParameter(it->first, myArray, it->second) != true) {  abort = true;  }
  72                         }
  73
  74                         //if the user changes the input directory command factory will send this info to us in the output parameter
  75                         string inputDir = validParameter.validFile(parameters, "inputdir", false);
  76                         if (inputDir == "not found"){   inputDir = "";          }
  77                         else {
  78                                 string path;
  79                                 it = parameters.find("sff");
  80                                 //user has given a template file
  81                                 if(it != parameters.end()){
  82                                         path = m->hasPath(it->second);
  83                                         //if the user has not given a path then, add inputdir. else leave path alone.
  84                                         if (path == "") {       parameters["sff"] = inputDir + it->second;              }
  85                                 }
  86
  87                                 it = parameters.find("oligos");
  88                                 //user has given a template file
  89                                 if(it != parameters.end()){
  90                                         path = m->hasPath(it->second);
  91                                         //if the user has not given a path then, add inputdir. else leave path alone.
  92                                         if (path == "") {       parameters["oligos"] = inputDir + it->second;           }
  93                                 }
  94
  95                                 it = parameters.find("align");
  96                                 //user has given a template file
  97                                 if(it != parameters.end()){
  98                                         path = m->hasPath(it->second);
  99                                         //if the user has not given a path then, add inputdir. else leave path alone.
 100                                         if (path == "") {       parameters["align"] = inputDir + it->second;            }
 101                                 }
 102
 103                                 it = parameters.find("chimera");
 104                                 //user has given a template file
 105                                 if(it != parameters.end()){
 106                                         path = m->hasPath(it->second);
 107                                         //if the user has not given a path then, add inputdir. else leave path alone.
 108                                         if (path == "") {       parameters["chimera"] = inputDir + it->second;          }
 109                                 }
 110
 111                                 it = parameters.find("classify");
 112                                 //user has given a template file
 113                                 if(it != parameters.end()){
 114                                         path = m->hasPath(it->second);
 115                                         //if the user has not given a path then, add inputdir. else leave path alone.
 116                                         if (path == "") {       parameters["classify"] = inputDir + it->second;         }
 117                                 }
 118
 119                                 it = parameters.find("taxonomy");
 120                                 //user has given a template file
 121                                 if(it != parameters.end()){
 122                                         path = m->hasPath(it->second);
 123                                         //if the user has not given a path then, add inputdir. else leave path alone.
 124                                         if (path == "") {       parameters["taxonomy"] = inputDir + it->second;         }
 125                                 }
 126
 127                                 it = parameters.find("pipeline");
 128                                 //user has given a template file
 129                                 if(it != parameters.end()){
 130                                         path = m->hasPath(it->second);
 131                                         //if the user has not given a path then, add inputdir. else leave path alone.
 132                                         if (path == "") {       parameters["pipeline"] = inputDir + it->second;         }
 133                                 }
 134                         }
 135
 136                         outputDir = validParameter.validFile(parameters, "outputdir", false);           if (outputDir == "not found"){  outputDir = ""; }
 137
 138                         pipeFilename = validParameter.validFile(parameters, "pipeline", true);
 139                         if (pipeFilename == "not found") { pipeFilename = "";  }
 140                         else if (pipeFilename == "not open") { pipeFilename = ""; abort = true; }
 141
 142                         string temp = validParameter.validFile(parameters, "processors", false);        if (temp == "not found"){       temp = "1";                             }
 143                         convert(temp, processors);
 144
 145                         if (pipeFilename != "") {
 146                                 abort = readUsersPipeline();
 147                         }else{
 148                                 sffFile = validParameter.validFile(parameters, "sff", true);
 149                                 if (sffFile == "not found") { m->mothurOut("sff is a required parameter for the pipeline command."); m->mothurOutEndLine(); abort = true;  }
 150                                 else if (sffFile == "not open") { sffFile = ""; abort = true; }
 151
 152                                 oligosFile = validParameter.validFile(parameters, "oligos", true);
 153                                 if (oligosFile == "not found") { m->mothurOut("oligos is a required parameter for the pipeline command."); m->mothurOutEndLine(); abort = true;  }
 154                                 else if (oligosFile == "not open") { oligosFile = ""; abort = true; }
 155
 156                                 alignFile = validParameter.validFile(parameters, "align", true);
 157                                 if (alignFile == "not found") { m->mothurOut("align is a required parameter for the pipeline command. Please provide the template to align with."); m->mothurOutEndLine(); abort = true;  }
 158                                 else if (alignFile == "not open") { alignFile = ""; abort = true; }
 159
 160                                 chimeraFile = validParameter.validFile(parameters, "chimera", true);
 161                                 if (chimeraFile == "not found") { m->mothurOut("chimera is a required parameter for the pipeline command. Please provide the template to check for chimeras with."); m->mothurOutEndLine(); abort = true;  }
 162                                 else if (chimeraFile == "not open") { chimeraFile = ""; abort = true; }
 163
 164                                 classifyFile = validParameter.validFile(parameters, "classify", true);
 165                                 if (classifyFile == "not found") { m->mothurOut("classify is a required parameter for the pipeline command. Please provide the template to use with the classifier."); m->mothurOutEndLine(); abort = true;  }
 166                                 else if (classifyFile == "not open") { classifyFile = ""; abort = true; }
 167
 168                                 taxonomyFile = validParameter.validFile(parameters, "taxonomy", true);
 169                                 if (taxonomyFile == "not found") { m->mothurOut("taxonomy is a required parameter for the pipeline command."); m->mothurOutEndLine(); abort = true;  }
 170                                 else if (taxonomyFile == "not open") { taxonomyFile = ""; abort = true; }
 171                         }
 172                 }
 173
 174         }
 175         catch(exception& e) {
 176                 m->errorOut(e, "PipelineCommand", "PipelineCommand");
 177                 exit(1);
 178         }
 179 }
 180 //**********************************************************************************************************************
 181
 182 void PipelineCommand::help(){
 183         try {
 184                  m->mothurOut("The pipeline command is designed to guide you through your analysis using mothur.\n");
 185                  m->mothurOut("The pipeline command parameters are pipeline, sff, oligos, align, chimera, classify, taxonomy and processors.\n");
 186                  m->mothurOut("The sff parameter allows you to enter your sff file. It is required, if not using pipeline parameter.\n");
 187                  m->mothurOut("The oligos parameter allows you to enter your oligos file. It is required, if not using pipeline parameter.\n");
 188                  m->mothurOut("The align parameter allows you to enter a template to use with the aligner. It is required, if not using pipeline parameter.\n");
 189                  m->mothurOut("The chimera parameter allows you to enter a template to use for chimera detection. It is required, if not using pipeline parameter.\n");
 190                  m->mothurOut("The classify parameter allows you to enter a template to use for classification. It is required, if not using pipeline parameter.\n");
 191                  m->mothurOut("The taxonomy parameter allows you to enter a taxonomy file for the classify template to use for classification. It is required, if not using pipeline parameter.\n");
 192                  m->mothurOut("The processors parameter allows you to specify the number of processors to use. The default is 1.\n");
 193                  m->mothurOut("The pipeline parameter allows you to enter your own pipeline file. This file should look like a mothur batchfile, but where you would be using a mothur generated file, you can use mothurmade instead.\n");
 194                  m->mothurOut("Example: trim.seqs(processors=8, allfiles=T, maxambig=0, maxhomop=8, flip=T, bdiffs=1, pdiffs=2, qwindowaverage=35, qwindowsize=50, fasta=may1.v13.fasta, oligos=may1.v13.oligos, qfile=may1.v13.qual)\n");
 195                  m->mothurOut("then, you could enter unique.seqs(fasta=mothurmade), and mothur would use the .trim.fasta file from the trim.seqs command. \n");
 196                  m->mothurOut("then you could enter align.seqs(candidate=mothurmade, template=silva.v13.align, processors=8). , and mothur would use the .trim.unique.fasta file from the unique.seqs command. \n");
 197                  m->mothurOut("If no pipeline file is given then mothur will use Pat's pipeline. \n\n");
 198                  m->mothurOut("Here is a list of the commands used in Pat's pipeline.\n");
 199                  m->mothurOut("All paralellized commands will use the processors you entered.\n");
 200                  m->mothurOut("The sffinfo command takes your sff file and extracts the fasta and quality files.\n");
 201                  m->mothurOut("The trim.seqs command uses your oligos file and the quality and fasta files generated by sffinfo.\n");
 202                  m->mothurOut("The trim.seqs command sets the following parameters: allfiles=T, maxambig=0, maxhomop=8, flip=T, bdiffs=1, pdiffs=2, qwindowaverage=35, qwindowsize=50.\n");
 203                  m->mothurOut("The unique.seqs command uses the trimmed fasta file and removes redundant sequences, don't worry the names file generated by unique.seqs will be used in the pipeline to make sure they are included.\n");
 204                  m->mothurOut("The align.seqs command aligns the unique sequences using the aligners default options. \n");
 205                  m->mothurOut("The screen.seqs command screens the sequences using optimize=end-minlength. \n");
 206                  m->mothurOut("The pipeline uses chimera.slayer to detect chimeras using the default options. \n");
 207                  m->mothurOut("The pipeline removes all sequences determined to be chimeric by chimera.slayer. \n");
 208                  m->mothurOut("The filter.seqs command filters the sequences using vertical=T, trump=. \n");
 209                  m->mothurOut("The unique.seqs command uses the filtered fasta file and name file to remove sequences that have become redundant after filtering.\n");
 210                  m->mothurOut("The pre.cluster command clusters sequences that have no more than 2 differences.\n");
 211                  m->mothurOut("The dist.seqs command is used to generate a column and phylip formatted distance matrix using cutoff=0.20 for column.\n");
 212                  m->mothurOut("The pipeline uses cluster with method=average, hard=T. \n");
 213                  m->mothurOut("The classify.seqs command is used to classify the sequences using the bayesian method with a cutoff of 80.\n");
 214                  m->mothurOut("The phylotype command is used to cluster the sequences based on their classification.\n");
 215                  m->mothurOut("The clearcut command is used to generate a tree using neighbor=T. \n");
 216                  m->mothurOut("The summary.single and summary.shared commands are run on the otu files from cluster and phylotype commands. \n");
 217                  m->mothurOut("The summary.shared command uses calc=sharednseqs-sharedsobs-sharedchao-sharedace-anderberg-jclass-jest-kulczynski-kulczynskicody-lennon-ochiai-sorclass-sorest-whittaker-braycurtis-jabund-morisitahorn-sorabund-thetan-thetayc. \n");
 218                  m->mothurOut("The summary.single command uses calc=nseqs-sobs-coverage-bergerparker-chao-ace-jack-bootstrap-boneh-efron-shen-solow-shannon-npshannon-invsimpson-qstat-simpsoneven-shannoneven-heip-smithwilson. \n");
 219                  m->mothurOut("The classify.otu command is used to get the concensus taxonomy for otu files from cluster and phylotype commands. \n");
 220                  m->mothurOut("The phylo.diversity command run on the tree generated by clearcut with rarefy=T, iters=100. \n");
 221                  m->mothurOut("The unifrac commands are also run on the tree generated by clearcut with random=F, distance=T. \n");
 222                  m->mothurOut("\n\n");
 223         }
 224         catch(exception& e) {
 225                 m->errorOut(e, "PipelineCommand", "help");
 226                 exit(1);
 227         }
 228 }
 229
 230 //**********************************************************************************************************************
 231
 232 PipelineCommand::~PipelineCommand(){}
 233
 234 //**********************************************************************************************************************
 235
 236 int PipelineCommand::execute(){
 237         try {
 238                 if (abort == true) { if (calledHelp) { return 0; }  return 2;   }
 239
 240                 int start = time(NULL);
 241
 242                 if (pipeFilename == "") {
 243                         createPatsPipeline();
 244
 245                         //run Pats pipeline
 246                         for (int i = 0; i < commands.size(); i++) {
 247                                 m->mothurOutEndLine(); m->mothurOut("mothur > " + commands[i]); m->mothurOutEndLine();
 248
 249                                 if (m->control_pressed) { return 0; }
 250
 251                                 CommandOptionParser parser(commands[i]);
 252                                 string commandName = parser.getCommandString();
 253                                 string options = parser.getOptionString();
 254
 255                                 #ifdef USE_MPI
 256                                         int pid;
 257                                         MPI_Comm_rank(MPI_COMM_WORLD, &pid);
 258
 259                                         if ((cFactory->MPIEnabled(commandName)) || (pid == 0)) {
 260                                 #endif
 261
 262                                 //executes valid command
 263                                 Command* command = cFactory->getCommand(commandName, options, "pipe");
 264                                 command->execute();
 265
 266                                 //add output files to list
 267                                 map<string, vector<string> > thisCommandsFile = command->getOutputFiles();
 268                                 map<string, vector<string> >::iterator itMade;
 269                                 for (itMade = thisCommandsFile.begin(); itMade != thisCommandsFile.end(); itMade++) {
 270                                         vector<string> temp = itMade->second;
 271                                         for (int j = 0; j < temp.size(); j++) { outputNames.push_back(temp[j]); }
 272                                 }
 273
 274                                 #ifdef USE_MPI
 275                                         }
 276                                 #endif
 277                         }
 278
 279                 }else {  runUsersPipeline(); }
 280
 281                 if (m->control_pressed) { return 0; }
 282
 283                 m->mothurOut("It took " + toString(time(NULL) - start) + " secs to run the pipeline analysis."); m->mothurOutEndLine(); m->mothurOutEndLine();
 284
 285                 m->mothurOutEndLine();
 286                 m->mothurOut("Output File Names: "); m->mothurOutEndLine();
 287                 for (int i = 0; i < outputNames.size(); i++) {  m->mothurOut(outputNames[i]); m->mothurOutEndLine();    }
 288                 m->mothurOutEndLine();
 289
 290                 return 0;
 291         }
 292         catch(exception& e) {
 293                 m->errorOut(e, "PipelineCommand", "execute");
 294                 exit(1);
 295         }
 296 }
 297 //**********************************************************************************************************************
 298
 299 bool PipelineCommand::readUsersPipeline(){
 300         try {
 301
 302                 ifstream in;
 303                 m->openInputFile(pipeFilename, in);
 304
 305                 string nextCommand = "";
 306
 307                 map<string, vector<string> > mothurMadeFiles;
 308
 309                 while(!in.eof()) {
 310                         nextCommand = m->getline(in); m->gobble(in);
 311
 312                         if (nextCommand[0] != '#') {
 313                                 bool error = false;
 314
 315                                 string commandName, options;
 316                                 error = parseCommand(nextCommand, commandName, options);
 317
 318                                 if (error) { in.close(); return error; }
 319                                 if (commandName == "pipeline.pds") { m->mothurOut("Cannot run the pipeline.pds command from inside the pipeline.pds command."); m->mothurOutEndLine(); in.close(); return true; }
 320
 321                                 error = checkForValidAndRequiredParameters(commandName, options, mothurMadeFiles);
 322
 323                                 if (error) { in.close(); return error; }
 324                         }
 325                 }
 326
 327                 in.close();
 328
 329                 return false;
 330         }
 331         catch(exception& e) {
 332                 m->errorOut(e, "PipelineCommand", "readUsersPipeline");
 333                 exit(1);
 334         }
 335 }
 336 //**********************************************************************************************************************
 337
 338 bool PipelineCommand::parseCommand(string nextCommand, string& name, string& options){
 339         try {
 340                 CommandOptionParser parser(nextCommand);
 341                 name = parser.getCommandString();
 342                 options = parser.getOptionString();
 343
 344                 if (name == "") { return true; } //name == "" if () are not right
 345
 346                 return false;
 347         }
 348         catch(exception& e) {
 349                 m->errorOut(e, "PipelineCommand", "parseCommand");
 350                 exit(1);
 351         }
 352 }
 353 //**********************************************************************************************************************
 354
 355 bool PipelineCommand::checkForValidAndRequiredParameters(string name, string options, map<string, vector<string> >& mothurMadeFiles){
 356         try {
 357
 358                 if (name == "system") { return false; }
 359
 360                 if (name == "system") { return false; }
 361
 362                 //get shell of the command so we can check to make sure its valid without running it
 363                 Command* command = cFactory->getCommand(name);
 364
 365                 //check to make sure all parameters are valid for command
 366                 vector<string> validParameters = command->getValidParameters();
 367
 368                 OptionParser parser(options);
 369                 map<string, string> parameters = parser.getParameters();
 370
 371                 ValidParameters validParameter;
 372                 map<string, string>::iterator it;
 373                 map<string, vector<string> >::iterator itMade;
 374
 375                 for (it = parameters.begin(); it != parameters.end(); it++) {
 376
 377                         if (validParameter.isValidParameter(it->first, validParameters, it->second) != true) {  return true;  } // not valid
 378                         if (it->second == "mothurmade") {
 379                                 itMade = mothurMadeFiles.find(it->first);
 380
 381                                 if (itMade == mothurMadeFiles.end()) {
 382                                         if ((name == "align.seqs") && (it->first == "candidate")) {} //do nothing about candidate
 383                                         else {
 384                                                 m->mothurOut("You have the " + it->first + " listed as a mothurmade file for the " + name + " command, but it seems mothur will not make that file in your current pipeline, please correct."); m->mothurOutEndLine();
 385                                                 return true;
 386                                         }
 387                                 }
 388                         }
 389                 }
 390
 391                 //is the command missing any required
 392                 vector<string> requiredParameters = command->getRequiredParameters();
 393
 394                 //check for or
 395                 bool hasOr = false;
 396                 int numFound = 0;
 397                 if (requiredParameters.size() > 2) {
 398                         if (requiredParameters[(requiredParameters.size()-1)] == "or") { hasOr = true; }
 399                 }
 400
 401                 for (int i = 0; i < requiredParameters.size(); i++) {
 402                         it = parameters.find(requiredParameters[i]);
 403
 404                         if (it != parameters.end()) { numFound++; }
 405                         else {
 406                                 if (!hasOr) { m->mothurOut(name + " requires the " + requiredParameters[i] + " parameter, please correct."); m->mothurOutEndLine(); }
 407                         }
 408                 }
 409
 410                 // if all are needed and not all are found
 411                 if ((!hasOr) && (numFound != requiredParameters.size())) { return true; }
 412                 //if one is needed and none are found
 413                 else if ((hasOr) && (numFound == 0)) { return true; }
 414
 415                 //update MothurMade
 416                 map<string, vector<string> > thisCommandsFile = command->getOutputFiles();
 417                 for (itMade = thisCommandsFile.begin(); itMade != thisCommandsFile.end(); itMade++) {
 418                         mothurMadeFiles[itMade->first] = itMade->second; //adds any new types
 419                 }
 420
 421                 return false;
 422         }
 423         catch(exception& e) {
 424                 m->errorOut(e, "PipelineCommand", "checkForValidAndRequiredParameters");
 425                 exit(1);
 426         }
 427 }
 428 //**********************************************************************************************************************
 429 int PipelineCommand::runUsersPipeline(){
 430         try {
 431                 ifstream in;
 432                 m->openInputFile(pipeFilename, in);
 433
 434                 string nextCommand = "";
 435
 436                 map<string, vector<string> > mothurMadeFiles;
 437
 438                 while(!in.eof()) {
 439                         nextCommand = m->getline(in);  m->gobble(in);
 440
 441                         if (nextCommand[0] != '#') {
 442                                 CommandOptionParser parser(nextCommand);
 443                                 string commandName = parser.getCommandString();
 444                                 string options = parser.getOptionString();
 445
 446                                 if ((options != "") && (commandName != "system")) {
 447                                         bool error = fillInMothurMade(options, mothurMadeFiles);
 448                                         if (error) { in.close(); return 0; }
 449                                 }
 450
 451                                 m->mothurOutEndLine(); m->mothurOut("mothur > " + commandName + "(" + options + ")"); m->mothurOutEndLine();
 452
 453                                 if (m->control_pressed) { return 0; }
 454
 455                                 #ifdef USE_MPI
 456                                         int pid;
 457                                         MPI_Comm_rank(MPI_COMM_WORLD, &pid);
 458
 459                                         if ((cFactory->MPIEnabled(commandName)) || (pid == 0)) {
 460                                 #endif
 461
 462                                 //executes valid command
 463                                 Command* command = cFactory->getCommand(commandName, options, "pipe");
 464                                 command->execute();
 465
 466                                 //add output files to list
 467                                 map<string, vector<string> > thisCommandsFile = command->getOutputFiles();
 468                                 map<string, vector<string> >::iterator itMade;
 469                                 map<string, vector<string> >::iterator it;
 470                                 for (itMade = thisCommandsFile.begin(); itMade != thisCommandsFile.end(); itMade++) {
 471
 472                                         vector<string> temp = itMade->second;
 473                                         for (int k = 0; k < temp.size(); k++) { outputNames.push_back(temp[k]); }  //
 474
 475                                         //update Mothur Made for each file
 476                                         it = mothurMadeFiles.find(itMade->first);
 477
 478                                         if (it == mothurMadeFiles.end()) { //new type
 479
 480                                                 mothurMadeFiles[itMade->first] = temp;
 481
 482                                         }else{ //update existing type
 483                                                 vector<string> oldFileNames = it->second;
 484                                                 //look at new files, see if an old version of the file exists, if so update, else just add.
 485                                                 //for example you may have abrecovery.fasta and amazon.fasta as old files and you created a new amazon.trim.fasta.
 486
 487                                                 for (int k = 0; k < temp.size(); k++) {
 488
 489                                                         //get base name
 490                                                         string root = m->getSimpleName(temp[k]);
 491                                                         string individual = "";
 492                                                         for(int i=0;i<root.length();i++){
 493                                                                 if(root[i] == '.'){
 494                                                                         root = individual;
 495                                                                         break;
 496                                                                 }else{
 497                                                                         individual += root[i];
 498                                                                 }
 499                                                         }
 500
 501                                                         //look for that base name in oldfiles
 502                                                         int spot = -1;
 503                                                         for (int l = 0; l < oldFileNames.size(); l++) {
 504                                                                 int pos = oldFileNames[l].find(root);
 505                                                                 if (pos != string::npos) {
 506                                                                         spot = l;
 507                                                                         break;
 508                                                                 }
 509                                                         }
 510
 511                                                         //if you found it update it, else add it
 512                                                         if (spot != -1) {
 513                                                                 mothurMadeFiles[it->first][spot] = temp[k];
 514                                                         }else{
 515                                                                 mothurMadeFiles[it->first].push_back(temp[k]);
 516                                                         }
 517                                                 }
 518                                         }
 519                                 }
 520
 521                                 #ifdef USE_MPI
 522                                         }
 523                                 #endif
 524                         }
 525                 }
 526
 527                 in.close();
 528
 529                 return 0;
 530         }
 531         catch(exception& e) {
 532                 m->errorOut(e, "PipelineCommand", "runUsersPipeline");
 533                 exit(1);
 534         }
 535 }
 536 //**********************************************************************************************************************
 537 bool PipelineCommand::fillInMothurMade(string& options, map<string, vector<string> >& mothurMadeFiles){
 538         try {
 539
 540                 OptionParser parser(options);
 541                 map<string, string> parameters = parser.getParameters();
 542                 map<string, string>::iterator it;
 543                 map<string, vector<string> >::iterator itMade;
 544
 545                 options = "";
 546
 547                 //fill in mothurmade filenames
 548                 for (it = parameters.begin(); it != parameters.end(); it++) {
 549                         string paraType = it->first;
 550                         string tempOption = it->second;
 551
 552                         if (tempOption == "mothurmade") {
 553
 554                                 if (it->first == "candidate") { paraType = "fasta"; }
 555
 556                                 itMade = mothurMadeFiles.find(paraType);
 557
 558                                 if (itMade == mothurMadeFiles.end()) {
 559                                         m->mothurOut("Looking for a mothurmade " + paraType + " file, but it seems mothur has not made that file type in your current pipeline, please correct."); m->mothurOutEndLine();
 560                                         return true;
 561                                 }else{
 562                                         vector<string> temp = itMade->second;
 563
 564                                         if (temp.size() > 1) {
 565                                                 //ask user which file to use
 566                                                 m->mothurOut("More than one file has been created for the " + paraType + " parameter. "); m->mothurOutEndLine();
 567                                                 for (int i = 0; i < temp.size(); i++) {
 568                                                         m->mothurOut(toString(i) + " - " + temp[i]); m->mothurOutEndLine();
 569                                                 }
 570
 571                                                 m->mothurOut("Please select the number of the file you would like to use: ");
 572                                                 int num = 0;
 573                                                 cin >> num;
 574                                                 m->mothurOutJustToLog(toString(num)); m->mothurOutEndLine();
 575
 576                                                 if ((num < 0) || (num > (temp.size()-1))) { m->mothurOut("Not a valid response, quitting."); m->mothurOutEndLine(); return true; }
 577                                                 else {
 578                                                         tempOption = temp[num];
 579                                                 }
 580
 581                                                 //clears buffer so next command doesn't have error
 582                                                 string s;
 583                                                 getline(cin, s);
 584
 585                                                 vector<string> newTemp;
 586                                                 for (int i = 0; i < temp.size(); i++) {
 587                                                         if (i == num) { newTemp.push_back(temp[i]); }
 588                                                         else {
 589                                                                 m->mothurOut("Would you like to remove " + temp[i] + " as an option for " + paraType + ", (y/n): "); m->mothurOutEndLine();
 590                                                                 string response;
 591                                                                 cin >> response;
 592                                                                 m->mothurOutJustToLog(response); m->mothurOutEndLine();
 593
 594                                                                 if (response == "n") {  newTemp.push_back(temp[i]); }
 595
 596                                                                 //clears buffer so next command doesn't have error
 597                                                                 string s;
 598                                                                 getline(cin, s);
 599                                                         }
 600                                                 }
 601
 602                                                 mothurMadeFiles[paraType] = newTemp;
 603
 604
 605                                         }else if (temp.size() == 0){
 606                                                 m->mothurOut("Sorry, we seem to think you created a " + paraType + " file, but it seems mothur doesn't have a filename."); m->mothurOutEndLine();
 607                                                 return true;
 608                                         }else{
 609                                                 tempOption = temp[0];
 610                                         }
 611                                 }
 612                         }
 613
 614                         options += it->first + "=" + tempOption + ", ";
 615                 }
 616
 617                 //rip off extra comma
 618                 options = options.substr(0, (options.length()-2));
 619
 620                 return false;
 621         }
 622         catch(exception& e) {
 623                 m->errorOut(e, "PipelineCommand", "fillInMothurMade");
 624                 exit(1);
 625         }
 626 }
 627
 628 //**********************************************************************************************************************
 629 void PipelineCommand::createPatsPipeline(){
 630         try {
 631
 632                 //sff.info command
 633                 string thisCommand = "sffinfo(sff=" + sffFile + ")";
 634                 commands.push_back(thisCommand);
 635
 636                 //trim.seqs command
 637                 string fastaFile = m->getRootName(m->getSimpleName(sffFile)) + "fasta";
 638                 string qualFile = m->getRootName(m->getSimpleName(sffFile)) + "qual";
 639                 thisCommand = "trim.seqs(processors=" + toString(processors) + ", fasta=" + fastaFile + ", allfiles=T, maxambig=0, maxhomop=8, flip=T, bdiffs=1, pdiffs=2, qwindowaverage=35, qwindowsize=50, oligos=" + oligosFile + ", qfile=" + qualFile + ")";
 640                 commands.push_back(thisCommand);
 641
 642                 //unique.seqs
 643                 string groupFile = m->getRootName(m->getSimpleName(fastaFile)) + "groups";
 644                 qualFile =  m->getRootName(m->getSimpleName(fastaFile)) + "trim.qual";
 645                 fastaFile =  m->getRootName(m->getSimpleName(fastaFile)) + "trim.fasta";
 646                 thisCommand = "unique.seqs(fasta=" + fastaFile + ")";
 647                 commands.push_back(thisCommand);
 648
 649                 //align.seqs
 650                 string nameFile = m->getRootName(m->getSimpleName(fastaFile)) + "names";
 651                 fastaFile = m->getRootName(m->getSimpleName(fastaFile)) + "unique" + m->getExtension(fastaFile);
 652                 thisCommand = "align.seqs(processors=" + toString(processors) + ", candidate=" + fastaFile + ", template=" + alignFile + ")";
 653                 commands.push_back(thisCommand);
 654
 655                 //screen.seqs
 656                 fastaFile = m->getRootName(m->getSimpleName(fastaFile)) + "align";
 657                 thisCommand = "screen.seqs(processors=" + toString(processors) + ", fasta=" + fastaFile + ", name=" + nameFile + ", group=" + groupFile + ", optimize=end-minlength)";
 658                 commands.push_back(thisCommand);
 659
 660                 //chimera.slayer
 661                 fastaFile = m->getRootName(m->getSimpleName(fastaFile)) + "good" + m->getExtension(fastaFile);
 662                 nameFile = m->getRootName(m->getSimpleName(nameFile)) + "good" + m->getExtension(nameFile);
 663                 groupFile = m->getRootName(m->getSimpleName(groupFile)) + "good" + m->getExtension(groupFile);
 664                 thisCommand = "chimera.slayer(processors=" + toString(processors) + ", fasta=" + fastaFile + ", template=" + chimeraFile + ")";
 665                 commands.push_back(thisCommand);
 666
 667                 //remove.seqs
 668                 string accnosFile = m->getRootName(m->getSimpleName(fastaFile))  + "slayer.accnos";
 669                 thisCommand = "remove.seqs(fasta=" + fastaFile + ", name=" + nameFile + ", group=" + groupFile + ", accnos=" + accnosFile + ", dups=T)";
 670                 commands.push_back(thisCommand);
 671
 672                 //filter.seqs
 673                 nameFile = m->getRootName(m->getSimpleName(nameFile)) + "pick" + m->getExtension(nameFile);
 674                 groupFile = m->getRootName(m->getSimpleName(groupFile)) + "pick" + m->getExtension(groupFile);
 675                 fastaFile = m->getRootName(m->getSimpleName(fastaFile)) + "pick" + m->getExtension(fastaFile);
 676                 thisCommand = "filter.seqs(processors=" + toString(processors) + ", fasta=" + fastaFile + ", vertical=T, trump=.)";
 677                 commands.push_back(thisCommand);
 678
 679                 //unique.seqs
 680                 fastaFile =  m->getRootName(m->getSimpleName(fastaFile)) + "filter.fasta";
 681                 thisCommand = "unique.seqs(fasta=" + fastaFile + ", name=" + nameFile + ")";
 682                 commands.push_back(thisCommand);
 683
 684                 //pre.cluster
 685                 nameFile = m->getRootName(m->getSimpleName(fastaFile)) + "names";
 686                 fastaFile = m->getRootName(m->getSimpleName(fastaFile)) + "unique" + m->getExtension(fastaFile);
 687                 thisCommand = "pre.cluster(fasta=" + fastaFile + ", name=" + nameFile + ", diffs=2)";
 688                 commands.push_back(thisCommand);
 689
 690                 //dist.seqs
 691                 nameFile = m->getRootName(m->getSimpleName(fastaFile)) + "precluster.names";
 692                 fastaFile = m->getRootName(m->getSimpleName(fastaFile)) + "precluster" + m->getExtension(fastaFile);
 693                 thisCommand = "dist.seqs(processors=" + toString(processors) + ", fasta=" + fastaFile + ", cutoff=0.20)";
 694                 commands.push_back(thisCommand);
 695
 696                 //dist.seqs
 697                 string columnFile = m->getRootName(m->getSimpleName(fastaFile)) + "dist";
 698                 thisCommand = "dist.seqs(processors=" + toString(processors) + ", fasta=" + fastaFile + ", output=lt)";
 699                 commands.push_back(thisCommand);
 700
 701                 //read.dist
 702                 string phylipFile = m->getRootName(m->getSimpleName(fastaFile)) + "phylip.dist";
 703                 thisCommand = "read.dist(column=" + columnFile + ", name=" + nameFile + ")";
 704                 commands.push_back(thisCommand);
 705
 706                 //cluster
 707                 thisCommand = "cluster(method=average, hard=T)";
 708                 commands.push_back(thisCommand);
 709
 710                 string listFile = m->getRootName(m->getSimpleName(columnFile)) + "an.list";
 711                 string rabundFile = m->getRootName(m->getSimpleName(columnFile)) + "an.rabund";
 712
 713                 //degap.seqs
 714                 thisCommand = "degap.seqs(fasta=" + fastaFile + ")";
 715                 commands.push_back(thisCommand);
 716
 717                 //classify.seqs
 718                 fastaFile = m->getRootName(m->getSimpleName(fastaFile)) + "ng.fasta";
 719                 thisCommand = "classify.seqs(processors=" + toString(processors) + ", fasta=" + fastaFile + ", name=" + nameFile + ", template=" + classifyFile + ", taxonomy=" + taxonomyFile + ", cutoff=80)";
 720                 commands.push_back(thisCommand);
 721
 722                 string RippedTaxName = m->getRootName(m->getSimpleName(taxonomyFile));
 723                 RippedTaxName = m->getExtension(RippedTaxName.substr(0, RippedTaxName.length()-1));
 724                 if (RippedTaxName[0] == '.') { RippedTaxName = RippedTaxName.substr(1, RippedTaxName.length()); }
 725                 RippedTaxName +=  ".";
 726
 727                 string fastaTaxFile = m->getRootName(m->getSimpleName(fastaFile)) + RippedTaxName + "taxonomy";
 728                 string taxSummaryFile = m->getRootName(m->getSimpleName(fastaFile)) + RippedTaxName + "tax.summary";
 729
 730                 //phylotype
 731                 thisCommand = "phylotype(taxonomy=" + fastaTaxFile + ", name=" + nameFile + ")";
 732                 commands.push_back(thisCommand);
 733
 734                 string phyloListFile = m->getRootName(m->getSimpleName(fastaTaxFile)) + "tx.list";
 735                 string phyloRabundFile = m->getRootName(m->getSimpleName(fastaTaxFile)) + "tx.rabund";
 736
 737                 //clearcut
 738                 thisCommand = "clearcut(phylip=" + phylipFile + ", neighbor=T)";
 739                 commands.push_back(thisCommand);
 740
 741                 string treeFile = m->getRootName(m->getSimpleName(phylipFile)) + "tre";
 742
 743                 //read.otu
 744                 thisCommand = "read.otu(list=" + listFile + ", group=" + groupFile + ", label=0.03)";
 745                 commands.push_back(thisCommand);
 746
 747                 string sharedFile = m->getRootName(m->getSimpleName(listFile)) + "shared";
 748
 749                 //read.otu
 750                 thisCommand = "read.otu(list=" + phyloListFile + ", group=" + groupFile + ", label=1)";
 751                 commands.push_back(thisCommand);
 752
 753                 string phyloSharedFile = m->getRootName(m->getSimpleName(phyloListFile)) + "shared";
 754
 755                 //read.otu
 756                 thisCommand = "read.otu(shared=" + sharedFile + ")";
 757                 commands.push_back(thisCommand);
 758
 759                 //summary.single
 760                 thisCommand = "summary.single(calc=nseqs-sobs-coverage-bergerparker-chao-ace-jack-bootstrap-boneh-efron-shen-solow-shannon-npshannon-invsimpson-qstat-simpsoneven-shannoneven-heip-smithwilson, size=5000)";
 761                 commands.push_back(thisCommand);
 762
 763                 //summary.shared
 764                 thisCommand = "summary.shared(calc=sharednseqs-sharedsobs-sharedchao-sharedace-anderberg-jclass-jest-kulczynski-kulczynskicody-lennon-ochiai-sorclass-sorest-whittaker-braycurtis-jabund-morisitahorn-sorabund-thetan-thetayc)";
 765                 commands.push_back(thisCommand);
 766
 767                 //read.otu
 768                 thisCommand = "read.otu(rabund=" + rabundFile + ", label=0.03)";
 769                 commands.push_back(thisCommand);
 770
 771                 //summary.single
 772                 thisCommand = "summary.single(calc=nseqs-sobs-coverage-bergerparker-chao-ace-jack-bootstrap-boneh-efron-shen-solow-shannon-npshannon-invsimpson-qstat-simpsoneven-shannoneven-heip-smithwilson, size=5000)";
 773                 commands.push_back(thisCommand);
 774
 775                 //read.otu
 776                 thisCommand = "read.otu(shared=" + phyloSharedFile + ")";
 777                 commands.push_back(thisCommand);
 778
 779                 //summary.single
 780                 thisCommand = "summary.single(calc=nseqs-sobs-coverage-bergerparker-chao-ace-jack-bootstrap-boneh-efron-shen-solow-shannon-npshannon-invsimpson-qstat-simpsoneven-shannoneven-heip-smithwilson, size=5000)";
 781                 commands.push_back(thisCommand);
 782
 783                 //summary.shared
 784                 thisCommand = "summary.shared(calc=sharednseqs-sharedsobs-sharedchao-sharedace-anderberg-jclass-jest-kulczynski-kulczynskicody-lennon-ochiai-sorclass-sorest-whittaker-braycurtis-jabund-morisitahorn-sorabund-thetan-thetayc)";
 785                 commands.push_back(thisCommand);
 786
 787                 //read.otu
 788                 thisCommand = "read.otu(rabund=" + phyloRabundFile + ", label=1)";
 789                 commands.push_back(thisCommand);
 790
 791                 //summary.single
 792                 thisCommand = "summary.single(calc=nseqs-sobs-coverage-bergerparker-chao-ace-jack-bootstrap-boneh-efron-shen-solow-shannon-npshannon-invsimpson-qstat-simpsoneven-shannoneven-heip-smithwilson, size=5000)";
 793                 commands.push_back(thisCommand);
 794
 795                 //classify.otu
 796                 thisCommand = "classify.otu(taxonomy=" + fastaTaxFile + ", name=" + nameFile + ", list=" + listFile + ", cutoff=51, label=0.03)";
 797                 commands.push_back(thisCommand);
 798
 799                 //classify.otu
 800                 thisCommand = "classify.otu(taxonomy=" + fastaTaxFile + ", name=" + nameFile + ", list=" + phyloListFile + ", cutoff=51, label=1)";
 801                 commands.push_back(thisCommand);
 802
 803                 //read.tree
 804                 thisCommand = "read.tree(tree=" + treeFile + ", name=" + nameFile + ", group=" + groupFile + ")";
 805                 commands.push_back(thisCommand);
 806
 807                 //phylo.diversity
 808                 thisCommand = "phylo.diversity(iters=100,rarefy=T)";
 809                 commands.push_back(thisCommand);
 810
 811                 //unifrac.weighted
 812                 thisCommand = "unifrac.weighted(random=false, distance=true, groups=all, processors=" + toString(processors) + ")";
 813                 commands.push_back(thisCommand);
 814
 815                 //unifrac.unweighted
 816                 thisCommand = "unifrac.unweighted(random=false, distance=true, processors=" + toString(processors) + ")";
 817                 commands.push_back(thisCommand);
 818
 819
 820         }
 821         catch(exception& e) {
 822                 m->errorOut(e, "PipelineCommand", "createPatsPipeline");
 823                 exit(1);
 824         }
 825 }
 826
 827 //**********************************************************************************************************************