removerarecommand.cpp

   1 /*
   2  *  removerarecommand.cpp
   3  *  mothur
   4  *
   5  *  Created by westcott on 1/21/11.
   6  *  Copyright 2011 Schloss Lab. All rights reserved.
   7  *
   8  */
   9
  10 #include "removerarecommand.h"
  11 #include "sequence.hpp"
  12 #include "groupmap.h"
  13 #include "sharedutilities.h"
  14 #include "inputdata.h"
  15
  16 //**********************************************************************************************************************
  17 vector<string> RemoveRareCommand::getValidParameters(){
  18         try {
  19                 string Array[] =  {"rabund","sabund", "group", "list", "shared","bygroup","nseqs","groups","label","outputdir","inputdir"};
  20                 vector<string> myArray (Array, Array+(sizeof(Array)/sizeof(string)));
  21                 return myArray;
  22         }
  23         catch(exception& e) {
  24                 m->errorOut(e, "RemoveRareCommand", "getValidParameters");
  25                 exit(1);
  26         }
  27 }
  28 //**********************************************************************************************************************
  29 RemoveRareCommand::RemoveRareCommand(){
  30         try {
  31                 abort = true; calledHelp = true;
  32                 vector<string> tempOutNames;
  33                 outputTypes["rabund"] = tempOutNames;
  34                 outputTypes["sabund"] = tempOutNames;
  35                 outputTypes["list"] = tempOutNames;
  36                 outputTypes["group"] = tempOutNames;
  37                 outputTypes["shared"] = tempOutNames;
  38         }
  39         catch(exception& e) {
  40                 m->errorOut(e, "RemoveRareCommand", "RemoveRareCommand");
  41                 exit(1);
  42         }
  43 }
  44 //**********************************************************************************************************************
  45 vector<string> RemoveRareCommand::getRequiredParameters(){
  46         try {
  47                 string Array[] =  {"nseqs"};
  48                 vector<string> myArray (Array, Array+(sizeof(Array)/sizeof(string)));
  49                 return myArray;
  50         }
  51         catch(exception& e) {
  52                 m->errorOut(e, "RemoveRareCommand", "getRequiredParameters");
  53                 exit(1);
  54         }
  55 }
  56 //**********************************************************************************************************************
  57 vector<string> RemoveRareCommand::getRequiredFiles(){
  58         try {
  59                 vector<string> myArray;
  60                 return myArray;
  61         }
  62         catch(exception& e) {
  63                 m->errorOut(e, "RemoveRareCommand", "getRequiredFiles");
  64                 exit(1);
  65         }
  66 }
  67 //**********************************************************************************************************************
  68 RemoveRareCommand::RemoveRareCommand(string option)  {
  69         try {
  70                 globaldata = GlobalData::getInstance();
  71                 abort = false; calledHelp = false;
  72                 allLines = 1;
  73
  74                 //allow user to run help
  75                 if(option == "help") { help(); abort = true; calledHelp = true; }
  76
  77                 else {
  78                         //valid paramters for this command
  79                         string Array[] =  {"rabund","sabund", "group", "list", "shared", "bygroup", "nseqs","groups","label","outputdir","inputdir"};
  80                         vector<string> myArray (Array, Array+(sizeof(Array)/sizeof(string)));
  81
  82                         OptionParser parser(option);
  83                         map<string,string> parameters = parser.getParameters();
  84
  85                         ValidParameters validParameter;
  86                         map<string,string>::iterator it;
  87
  88                         //check to make sure all parameters are valid for command
  89                         for (it = parameters.begin(); it != parameters.end(); it++) {
  90                                 if (validParameter.isValidParameter(it->first, myArray, it->second) != true) {  abort = true;  }
  91                         }
  92
  93                         //initialize outputTypes
  94                         vector<string> tempOutNames;
  95                         outputTypes["rabund"] = tempOutNames;
  96                         outputTypes["sabund"] = tempOutNames;
  97                         outputTypes["list"] = tempOutNames;
  98                         outputTypes["group"] = tempOutNames;
  99                         outputTypes["shared"] = tempOutNames;
 100
 101                         //if the user changes the output directory command factory will send this info to us in the output parameter
 102                         outputDir = validParameter.validFile(parameters, "outputdir", false);           if (outputDir == "not found"){  outputDir = "";         }
 103
 104                         //if the user changes the input directory command factory will send this info to us in the output parameter
 105                         string inputDir = validParameter.validFile(parameters, "inputdir", false);
 106                         if (inputDir == "not found"){   inputDir = "";          }
 107                         else {
 108                                 string path;
 109                                 it = parameters.find("list");
 110                                 //user has given a template file
 111                                 if(it != parameters.end()){
 112                                         path = m->hasPath(it->second);
 113                                         //if the user has not given a path then, add inputdir. else leave path alone.
 114                                         if (path == "") {       parameters["list"] = inputDir + it->second;             }
 115                                 }
 116
 117                                 it = parameters.find("group");
 118                                 //user has given a template file
 119                                 if(it != parameters.end()){
 120                                         path = m->hasPath(it->second);
 121                                         //if the user has not given a path then, add inputdir. else leave path alone.
 122                                         if (path == "") {       parameters["group"] = inputDir + it->second;            }
 123                                 }
 124
 125                                 it = parameters.find("sabund");
 126                                 //user has given a template file
 127                                 if(it != parameters.end()){
 128                                         path = m->hasPath(it->second);
 129                                         //if the user has not given a path then, add inputdir. else leave path alone.
 130                                         if (path == "") {       parameters["sabund"] = inputDir + it->second;           }
 131                                 }
 132
 133                                 it = parameters.find("rabund");
 134                                 //user has given a template file
 135                                 if(it != parameters.end()){
 136                                         path = m->hasPath(it->second);
 137                                         //if the user has not given a path then, add inputdir. else leave path alone.
 138                                         if (path == "") {       parameters["rabund"] = inputDir + it->second;           }
 139                                 }
 140
 141                                 it = parameters.find("shared");
 142                                 //user has given a template file
 143                                 if(it != parameters.end()){
 144                                         path = m->hasPath(it->second);
 145                                         //if the user has not given a path then, add inputdir. else leave path alone.
 146                                         if (path == "") {       parameters["shared"] = inputDir + it->second;           }
 147                                 }
 148                         }
 149
 150
 151                         //check for file parameters
 152                         listfile = validParameter.validFile(parameters, "list", true);
 153                         if (listfile == "not open") { abort = true; }
 154                         else if (listfile == "not found") {  listfile = "";  }
 155
 156                         sabundfile = validParameter.validFile(parameters, "sabund", true);
 157                         if (sabundfile == "not open") { abort = true; }
 158                         else if (sabundfile == "not found") {  sabundfile = "";  }
 159
 160                         rabundfile = validParameter.validFile(parameters, "rabund", true);
 161                         if (rabundfile == "not open") { abort = true; }
 162                         else if (rabundfile == "not found") {  rabundfile = "";  }
 163
 164                         groupfile = validParameter.validFile(parameters, "group", true);
 165                         if (groupfile == "not open") { groupfile = ""; abort = true; }
 166                         else if (groupfile == "not found") {  groupfile = "";  }
 167
 168                         sharedfile = validParameter.validFile(parameters, "shared", true);
 169                         if (sharedfile == "not open") { sharedfile = "";  abort = true; }
 170                         else if (sharedfile == "not found") {  sharedfile = "";  }
 171
 172                         if ((sabundfile == "") && (rabundfile == "")  && (sharedfile == "") && (listfile == ""))  { m->mothurOut("You must provide at least one of the following: rabund, sabund, shared or list."); m->mothurOutEndLine(); abort = true; }
 173
 174                         groups = validParameter.validFile(parameters, "groups", false);
 175                         if (groups == "not found") { groups = "all"; }
 176                         m->splitAtDash(groups, Groups);
 177
 178                         label = validParameter.validFile(parameters, "label", false);
 179                         if (label == "not found") { label = ""; }
 180                         else {
 181                                 if(label != "all") {  m->splitAtDash(label, labels);  allLines = 0;  }
 182                                 else { allLines = 1;  }
 183                         }
 184
 185                         string temp = validParameter.validFile(parameters, "nseqs", false);
 186                         if (temp == "not found") { m->mothurOut("nseqs is a required parameter."); m->mothurOutEndLine(); abort = true; }
 187                         else { convert(temp, nseqs); }
 188
 189                         temp = validParameter.validFile(parameters, "bygroup", false);   if (temp == "not found") { temp = "f"; }
 190                         byGroup = m->isTrue(temp);
 191
 192                         if (byGroup && (sharedfile == "")) { m->mothurOut("The byGroup parameter is only valid with a shared file."); m->mothurOutEndLine(); }
 193
 194                         if ((groupfile != "") && (listfile == "")) { m->mothurOut("A groupfile is only valid with a list file."); m->mothurOutEndLine(); groupfile = ""; }
 195                 }
 196
 197         }
 198         catch(exception& e) {
 199                 m->errorOut(e, "RemoveRareCommand", "RemoveRareCommand");
 200                 exit(1);
 201         }
 202 }
 203 //**********************************************************************************************************************
 204
 205 void RemoveRareCommand::help(){
 206         try {
 207                 m->mothurOut("The remove.rare command parameters are list, rabund, sabund, shared, group, label, groups, bygroup and nseqs.\n");
 208                 m->mothurOut("The remove.rare command reads one of the following file types: list, rabund, sabund or shared file. It outputs a new file after removing the rare otus.\n");
 209                 m->mothurOut("The groups parameter allows you to specify which of the groups you would like analyzed.  Default=all. You may separate group names with dashes.\n");
 210                 m->mothurOut("The label parameter is used to analyze specific labels in your input. default=all. You may separate label names with dashes.\n");
 211                 m->mothurOut("The bygroup parameter is only valid with the shared file. default=f, meaning remove any OTU that has nseqs or fewer sequences across all groups.\n");
 212                 m->mothurOut("bygroups=T means remove any OTU that has nseqs or fewer sequences in each group (if groupA has 1 sequence and group B has 100 sequences in OTUZ and nseqs=1, then set the groupA count for OTUZ to 0 and keep groupB's count at 100.) \n");
 213                 m->mothurOut("The nseqs parameter allows you to set the cutoff for an otu to be deemed rare. It is required.\n");
 214                 m->mothurOut("The remove.rare command should be in the following format: remove.rare(shared=yourSharedFile, nseqs=yourRareCutoff).\n");
 215                 m->mothurOut("Example remove.rare(shared=amazon.fn.shared, nseqs=2).\n");
 216                 m->mothurOut("Note: No spaces between parameter labels (i.e. shared), '=' and parameters (i.e.yourSharedFile).\n\n");
 217         }
 218         catch(exception& e) {
 219                 m->errorOut(e, "RemoveRareCommand", "help");
 220                 exit(1);
 221         }
 222 }
 223
 224 //**********************************************************************************************************************
 225
 226 int RemoveRareCommand::execute(){
 227         try {
 228
 229                 if (abort == true) { if (calledHelp) { return 0; }  return 2;   }
 230
 231                 if (m->control_pressed) { return 0; }
 232
 233                 //read through the correct file and output lines you want to keep
 234                 if (sabundfile != "")           {               processSabund();        }
 235                 if (rabundfile != "")           {               processRabund();        }
 236                 if (listfile != "")                     {               processList();          }
 237                 if (sharedfile != "")           {               processShared();        }
 238
 239                 if (m->control_pressed) { for (int i = 0; i < outputNames.size(); i++) {        remove(outputNames[i].c_str()); } return 0; }
 240
 241                 if (outputNames.size() != 0) {
 242                         m->mothurOutEndLine();
 243                         m->mothurOut("Output File Names: "); m->mothurOutEndLine();
 244                         for (int i = 0; i < outputNames.size(); i++) {  m->mothurOut(outputNames[i]); m->mothurOutEndLine();    }
 245                         m->mothurOutEndLine();
 246
 247                         //set rabund file as new current rabundfile
 248                         string current = "";
 249                         itTypes = outputTypes.find("rabund");
 250                         if (itTypes != outputTypes.end()) {
 251                                 if ((itTypes->second).size() != 0) { current = (itTypes->second)[0]; m->setRabundFile(current); }
 252                         }
 253
 254                         itTypes = outputTypes.find("sabund");
 255                         if (itTypes != outputTypes.end()) {
 256                                 if ((itTypes->second).size() != 0) { current = (itTypes->second)[0]; m->setSabundFile(current); }
 257                         }
 258
 259                         itTypes = outputTypes.find("group");
 260                         if (itTypes != outputTypes.end()) {
 261                                 if ((itTypes->second).size() != 0) { current = (itTypes->second)[0]; m->setGroupFile(current); }
 262                         }
 263
 264                         itTypes = outputTypes.find("list");
 265                         if (itTypes != outputTypes.end()) {
 266                                 if ((itTypes->second).size() != 0) { current = (itTypes->second)[0]; m->setListFile(current); }
 267                         }
 268
 269                         itTypes = outputTypes.find("shared");
 270                         if (itTypes != outputTypes.end()) {
 271                                 if ((itTypes->second).size() != 0) { current = (itTypes->second)[0]; m->setSharedFile(current); }
 272                         }
 273                 }
 274
 275                 return 0;
 276         }
 277
 278         catch(exception& e) {
 279                 m->errorOut(e, "RemoveRareCommand", "execute");
 280                 exit(1);
 281         }
 282 }
 283
 284 //**********************************************************************************************************************
 285 int RemoveRareCommand::processList(){
 286         try {
 287                 string thisOutputDir = outputDir;
 288                 if (outputDir == "") {  thisOutputDir += m->hasPath(listfile);  }
 289                 string outputFileName = thisOutputDir + m->getRootName(m->getSimpleName(listfile)) + "pick" +  m->getExtension(listfile);
 290                 string outputGroupFileName = thisOutputDir + m->getRootName(m->getSimpleName(groupfile)) + "pick" +  m->getExtension(groupfile);
 291
 292                 ofstream out, outGroup;
 293                 m->openOutputFile(outputFileName, out);
 294
 295                 bool wroteSomething = false;
 296
 297                 //you must provide a label because the names in the listfile need to be consistent
 298                 string thisLabel = "";
 299                 if (allLines) { m->mothurOut("For the listfile you must select one label, using first label in your listfile."); m->mothurOutEndLine(); }
 300                 else if (labels.size() > 1) { m->mothurOut("For the listfile you must select one label, using " + (*labels.begin()) + "."); m->mothurOutEndLine(); thisLabel = *labels.begin(); }
 301                 else { thisLabel = *labels.begin(); }
 302
 303                 InputData input(listfile, "list");
 304                 ListVector* list = input.getListVector();
 305
 306                 //get first one or the one we want
 307                 if (thisLabel != "") {
 308                         //use smart distancing
 309                         set<string> userLabels; userLabels.insert(thisLabel);
 310                         set<string> processedLabels;
 311                         string lastLabel = list->getLabel();
 312                         while((list != NULL) && (userLabels.size() != 0)) {
 313                                 if(userLabels.count(list->getLabel()) == 1){
 314                                         processedLabels.insert(list->getLabel());
 315                                         userLabels.erase(list->getLabel());
 316                                         break;
 317                                 }
 318
 319                                 if ((m->anyLabelsToProcess(list->getLabel(), userLabels, "") == true) && (processedLabels.count(lastLabel) != 1)) {
 320                                         processedLabels.insert(list->getLabel());
 321                                         userLabels.erase(list->getLabel());
 322                                         delete list;
 323                                         list = input.getListVector(lastLabel);
 324                                         break;
 325                                 }
 326                                 lastLabel = list->getLabel();
 327                                 delete list;
 328                                 list = input.getListVector();
 329                         }
 330                         if (userLabels.size() != 0) {
 331                                 m->mothurOut("Your file does not include the label " + thisLabel + ". I will use " + lastLabel + ".");  m->mothurOutEndLine();
 332                                 list = input.getListVector(lastLabel);
 333                         }
 334                 }
 335
 336                 //if groupfile is given then use it
 337                 GroupMap* groupMap;
 338                 if (groupfile != "") {
 339                         groupMap = new GroupMap(groupfile); groupMap->readMap();
 340                         SharedUtil util;
 341                         util.setGroups(Groups, groupMap->namesOfGroups);
 342                         m->openOutputFile(outputGroupFileName, outGroup);
 343                 }
 344
 345
 346                 if (list != NULL) {
 347                         //make a new list vector
 348                         ListVector newList;
 349                         newList.setLabel(list->getLabel());
 350
 351                         //for each bin
 352                         for (int i = 0; i < list->getNumBins(); i++) {
 353                                 if (m->control_pressed) {  if (groupfile != "") { delete groupMap; outGroup.close(); remove(outputGroupFileName.c_str()); } out.close();  remove(outputFileName.c_str());  return 0; }
 354
 355                                 //parse out names that are in accnos file
 356                                 string binnames = list->get(i);
 357                                 vector<string> names;
 358                                 string saveBinNames = binnames;
 359                                 m->splitAtComma(binnames, names);
 360
 361                                 vector<string> newGroupFile;
 362                                 if (groupfile != "") {
 363                                         vector<string> newNames;
 364                                         saveBinNames = "";
 365                                         for(int k = 0; k < names.size(); k++) {
 366                                                 string group = groupMap->getGroup(names[k]);
 367
 368                                                 if (m->inUsersGroups(group, Groups)) {
 369                                                         newGroupFile.push_back(names[k] + "\t" + group);
 370
 371                                                         newNames.push_back(names[k]);
 372                                                         saveBinNames += names[k] + ",";
 373                                                 }
 374                                         }
 375                                         names = newNames;
 376                                         saveBinNames = saveBinNames.substr(0, saveBinNames.length()-1);
 377                                 }
 378
 379                                 if (names.size() > nseqs) { //keep bin
 380                                         newList.push_back(saveBinNames);
 381                                         for(int k = 0; k < newGroupFile.size(); k++) { outGroup << newGroupFile[k] << endl; }
 382                                 }
 383                         }
 384
 385                         //print new listvector
 386                         if (newList.getNumBins() != 0) {
 387                                 wroteSomething = true;
 388                                 newList.print(out);
 389                         }
 390                 }
 391
 392                 out.close();
 393                 if (groupfile != "") { outGroup.close(); outputTypes["group"].push_back(outputGroupFileName); outputNames.push_back(outputGroupFileName); }
 394
 395                 if (wroteSomething == false) {  m->mothurOut("Your file contains only rare sequences."); m->mothurOutEndLine();  }
 396                 outputTypes["list"].push_back(outputFileName); outputNames.push_back(outputFileName);
 397
 398                 return 0;
 399         }
 400         catch(exception& e) {
 401                 m->errorOut(e, "RemoveRareCommand", "processList");
 402                 exit(1);
 403         }
 404 }
 405 //**********************************************************************************************************************
 406 int RemoveRareCommand::processSabund(){
 407         try {
 408                 string thisOutputDir = outputDir;
 409                 if (outputDir == "") {  thisOutputDir += m->hasPath(sabundfile);  }
 410                 string outputFileName = thisOutputDir + m->getRootName(m->getSimpleName(sabundfile)) + "pick" +  m->getExtension(sabundfile);
 411                 outputTypes["sabund"].push_back(outputFileName); outputNames.push_back(outputFileName);
 412
 413                 ofstream out;
 414                 m->openOutputFile(outputFileName, out);
 415
 416                 //if the users enters label "0.06" and there is no "0.06" in their file use the next lowest label.
 417                 InputData input(sabundfile, "sabund");
 418                 SAbundVector* sabund = input.getSAbundVector();
 419                 string lastLabel = sabund->getLabel();
 420                 set<string> processedLabels;
 421                 set<string> userLabels = labels;
 422
 423                 while((sabund != NULL) && ((allLines == 1) || (userLabels.size() != 0))) {
 424
 425                         if (m->control_pressed) { delete sabund; out.close(); return 0; }
 426
 427                         if(allLines == 1 || labels.count(sabund->getLabel()) == 1){
 428
 429                                 m->mothurOut(sabund->getLabel()); m->mothurOutEndLine();
 430                                 processedLabels.insert(sabund->getLabel());
 431                                 userLabels.erase(sabund->getLabel());
 432
 433                                 if (sabund->getMaxRank() > nseqs) {
 434                                         for(int i = 1; i <=nseqs; i++) {  sabund->set(i, 0); }
 435                                 }else { sabund->clear(); }
 436
 437                                 if (sabund->getNumBins() > 0) { sabund->print(out); }
 438                         }
 439
 440                         if ((m->anyLabelsToProcess(sabund->getLabel(), userLabels, "") == true) && (processedLabels.count(lastLabel) != 1)) {
 441                                 string saveLabel = sabund->getLabel();
 442
 443                                 delete sabund;
 444                                 sabund = input.getSAbundVector(lastLabel);
 445
 446                                 m->mothurOut(sabund->getLabel()); m->mothurOutEndLine();
 447                                 processedLabels.insert(sabund->getLabel());
 448                                 userLabels.erase(sabund->getLabel());
 449
 450                                 if (sabund->getMaxRank() > nseqs) {
 451                                         for(int i = 1; i <=nseqs; i++) {  sabund->set(i, 0); }
 452                                 }else { sabund->clear(); }
 453
 454                                 if (sabund->getNumBins() > 0) { sabund->print(out); }
 455
 456                                 //restore real lastlabel to save below
 457                                 sabund->setLabel(saveLabel);
 458                         }
 459
 460                         lastLabel = sabund->getLabel();
 461
 462                         delete sabund;
 463                         sabund = input.getSAbundVector();
 464                 }
 465
 466                 if (m->control_pressed) {  out.close(); return 0; }
 467
 468                 //output error messages about any remaining user labels
 469                 set<string>::iterator it;
 470                 bool needToRun = false;
 471                 for (it = userLabels.begin(); it != userLabels.end(); it++) {
 472                         m->mothurOut("Your file does not include the label " + *it);
 473                         if (processedLabels.count(lastLabel) != 1) {
 474                                 m->mothurOut(". I will use " + lastLabel + "."); m->mothurOutEndLine();
 475                                 needToRun = true;
 476                         }else {
 477                                 m->mothurOut(". Please refer to " + lastLabel + "."); m->mothurOutEndLine();
 478                         }
 479                 }
 480
 481                 //run last label if you need to
 482                 if (needToRun == true)  {
 483                         if (sabund != NULL) {   delete sabund;  }
 484                         sabund = input.getSAbundVector(lastLabel);
 485
 486                         m->mothurOut(sabund->getLabel()); m->mothurOutEndLine();
 487
 488                         if (sabund->getMaxRank() > nseqs) {
 489                                 for(int i = 1; i <=nseqs; i++) {  sabund->set(i, 0); }
 490                         }else { sabund->clear(); }
 491
 492                         if (sabund->getNumBins() > 0) { sabund->print(out); }
 493
 494                         delete sabund;
 495                 }
 496
 497                 return 0;
 498         }
 499         catch(exception& e) {
 500                 m->errorOut(e, "RemoveRareCommand", "processSabund");
 501                 exit(1);
 502         }
 503 }
 504 //**********************************************************************************************************************
 505 int RemoveRareCommand::processRabund(){
 506         try {
 507                 string thisOutputDir = outputDir;
 508                 if (outputDir == "") {  thisOutputDir += m->hasPath(rabundfile);  }
 509                 string outputFileName = thisOutputDir + m->getRootName(m->getSimpleName(rabundfile)) + "pick" +  m->getExtension(rabundfile);
 510                 outputTypes["rabund"].push_back(outputFileName); outputNames.push_back(outputFileName);
 511
 512                 ofstream out;
 513                 m->openOutputFile(outputFileName, out);
 514
 515                 //if the users enters label "0.06" and there is no "0.06" in their file use the next lowest label.
 516                 InputData input(rabundfile, "rabund");
 517                 RAbundVector* rabund = input.getRAbundVector();
 518                 string lastLabel = rabund->getLabel();
 519                 set<string> processedLabels;
 520                 set<string> userLabels = labels;
 521
 522                 while((rabund != NULL) && ((allLines == 1) || (userLabels.size() != 0))) {
 523
 524                         if (m->control_pressed) { delete rabund; out.close(); return 0; }
 525
 526                         if(allLines == 1 || labels.count(rabund->getLabel()) == 1){
 527
 528                                 m->mothurOut(rabund->getLabel()); m->mothurOutEndLine();
 529                                 processedLabels.insert(rabund->getLabel());
 530                                 userLabels.erase(rabund->getLabel());
 531
 532                                 RAbundVector newRabund; newRabund.setLabel(rabund->getLabel());
 533                                 for (int i = 0; i < rabund->getNumBins(); i++) {
 534                                         if (rabund->get(i) > nseqs) {
 535                                                 newRabund.push_back(rabund->get(i));
 536                                         }
 537                                 }
 538                                 if (newRabund.getNumBins() > 0) { newRabund.print(out); }
 539                         }
 540
 541                         if ((m->anyLabelsToProcess(rabund->getLabel(), userLabels, "") == true) && (processedLabels.count(lastLabel) != 1)) {
 542                                 string saveLabel = rabund->getLabel();
 543
 544                                 delete rabund;
 545                                 rabund = input.getRAbundVector(lastLabel);
 546
 547                                 m->mothurOut(rabund->getLabel()); m->mothurOutEndLine();
 548                                 processedLabels.insert(rabund->getLabel());
 549                                 userLabels.erase(rabund->getLabel());
 550
 551                                 RAbundVector newRabund; newRabund.setLabel(rabund->getLabel());
 552                                 for (int i = 0; i < rabund->getNumBins(); i++) {
 553                                         if (rabund->get(i) > nseqs) {
 554                                                 newRabund.push_back(rabund->get(i));
 555                                         }
 556                                 }
 557                                 if (newRabund.getNumBins() > 0) { newRabund.print(out); }
 558
 559                                 //restore real lastlabel to save below
 560                                 rabund->setLabel(saveLabel);
 561                         }
 562
 563                         lastLabel = rabund->getLabel();
 564
 565                         delete rabund;
 566                         rabund = input.getRAbundVector();
 567                 }
 568
 569                 if (m->control_pressed) {  out.close(); return 0; }
 570
 571                 //output error messages about any remaining user labels
 572                 set<string>::iterator it;
 573                 bool needToRun = false;
 574                 for (it = userLabels.begin(); it != userLabels.end(); it++) {
 575                         m->mothurOut("Your file does not include the label " + *it);
 576                         if (processedLabels.count(lastLabel) != 1) {
 577                                 m->mothurOut(". I will use " + lastLabel + "."); m->mothurOutEndLine();
 578                                 needToRun = true;
 579                         }else {
 580                                 m->mothurOut(". Please refer to " + lastLabel + "."); m->mothurOutEndLine();
 581                         }
 582                 }
 583
 584                 //run last label if you need to
 585                 if (needToRun == true)  {
 586                         if (rabund != NULL) {   delete rabund;  }
 587                         rabund = input.getRAbundVector(lastLabel);
 588
 589                         m->mothurOut(rabund->getLabel()); m->mothurOutEndLine();
 590
 591                         RAbundVector newRabund; newRabund.setLabel(rabund->getLabel());
 592                         for (int i = 0; i < rabund->getNumBins(); i++) {
 593                                 if (rabund->get(i) > nseqs) {
 594                                         newRabund.push_back(rabund->get(i));
 595                                 }
 596                         }
 597                         if (newRabund.getNumBins() > 0) { newRabund.print(out); }
 598
 599                         delete rabund;
 600                 }
 601
 602                 return 0;
 603         }
 604         catch(exception& e) {
 605                 m->errorOut(e, "RemoveRareCommand", "processRabund");
 606                 exit(1);
 607         }
 608 }
 609 //**********************************************************************************************************************
 610 int RemoveRareCommand::processShared(){
 611         try {
 612                 globaldata->Groups = Groups;
 613
 614                 string thisOutputDir = outputDir;
 615                 if (outputDir == "") {  thisOutputDir += m->hasPath(sharedfile);  }
 616                 string outputFileName = thisOutputDir + m->getRootName(m->getSimpleName(sharedfile)) + "pick" +  m->getExtension(sharedfile);
 617                 outputTypes["shared"].push_back(outputFileName); outputNames.push_back(outputFileName);
 618
 619                 ofstream out;
 620                 m->openOutputFile(outputFileName, out);
 621
 622                 //if the users enters label "0.06" and there is no "0.06" in their file use the next lowest label.
 623                 InputData input(sharedfile, "sharedfile");
 624                 vector<SharedRAbundVector*> lookup = input.getSharedRAbundVectors();
 625                 string lastLabel = lookup[0]->getLabel();
 626                 set<string> processedLabels;
 627                 set<string> userLabels = labels;
 628
 629                 while((lookup[0] != NULL) && ((allLines == 1) || (userLabels.size() != 0))) {
 630
 631                         if (m->control_pressed) { for (int i = 0; i < lookup.size(); i++) {  delete lookup[i];  }  out.close(); return 0; }
 632
 633                         if(allLines == 1 || labels.count(lookup[0]->getLabel()) == 1){
 634
 635                                 m->mothurOut(lookup[0]->getLabel()); m->mothurOutEndLine();
 636                                 processedLabels.insert(lookup[0]->getLabel());
 637                                 userLabels.erase(lookup[0]->getLabel());
 638
 639                                 processLookup(lookup, out);
 640                         }
 641
 642                         if ((m->anyLabelsToProcess(lookup[0]->getLabel(), userLabels, "") == true) && (processedLabels.count(lastLabel) != 1)) {
 643                                 string saveLabel = lookup[0]->getLabel();
 644
 645                                 for (int i = 0; i < lookup.size(); i++) {  delete lookup[i];  }
 646                                 lookup = input.getSharedRAbundVectors(lastLabel);
 647
 648                                 m->mothurOut(lookup[0]->getLabel()); m->mothurOutEndLine();
 649                                 processedLabels.insert(lookup[0]->getLabel());
 650                                 userLabels.erase(lookup[0]->getLabel());
 651
 652                                 processLookup(lookup, out);
 653
 654                                 //restore real lastlabel to save below
 655                                 lookup[0]->setLabel(saveLabel);
 656                         }
 657
 658                         lastLabel = lookup[0]->getLabel();
 659
 660                         for (int i = 0; i < lookup.size(); i++) {  delete lookup[i];  }
 661                         lookup = input.getSharedRAbundVectors();
 662                 }
 663
 664                 if (m->control_pressed) {  out.close(); return 0; }
 665
 666                 //output error messages about any remaining user labels
 667                 set<string>::iterator it;
 668                 bool needToRun = false;
 669                 for (it = userLabels.begin(); it != userLabels.end(); it++) {
 670                         m->mothurOut("Your file does not include the label " + *it);
 671                         if (processedLabels.count(lastLabel) != 1) {
 672                                 m->mothurOut(". I will use " + lastLabel + "."); m->mothurOutEndLine();
 673                                 needToRun = true;
 674                         }else {
 675                                 m->mothurOut(". Please refer to " + lastLabel + "."); m->mothurOutEndLine();
 676                         }
 677                 }
 678
 679                 //run last label if you need to
 680                 if (needToRun == true)  {
 681                         for (int i = 0; i < lookup.size(); i++) {  if (lookup[i] != NULL) {     delete lookup[i];       }  }
 682                         lookup = input.getSharedRAbundVectors(lastLabel);
 683
 684                         m->mothurOut(lookup[0]->getLabel()); m->mothurOutEndLine();
 685
 686                         processLookup(lookup, out);
 687
 688                         for (int i = 0; i < lookup.size(); i++) {  delete lookup[i];  }
 689                 }
 690
 691                 return 0;
 692         }
 693         catch(exception& e) {
 694                 m->errorOut(e, "RemoveRareCommand", "processSabund");
 695                 exit(1);
 696         }
 697 }
 698 //**********************************************************************************************************************
 699 int RemoveRareCommand::processLookup(vector<SharedRAbundVector*>& lookup, ofstream& out){
 700         try {
 701
 702                 vector<SharedRAbundVector> newRabunds;  newRabunds.resize(lookup.size());
 703                 for (int i = 0; i < lookup.size(); i++) {
 704                         newRabunds[i].setGroup(lookup[i]->getGroup());
 705                         newRabunds[i].setLabel(lookup[i]->getLabel());
 706                 }
 707
 708                 if (byGroup) {
 709
 710                         //for each otu
 711                         for (int i = 0; i < lookup[0]->getNumBins(); i++) {
 712                                 bool allZero = true;
 713
 714                                 if (m->control_pressed) { return 0; }
 715
 716                                 //for each group
 717                                 for (int j = 0; j < lookup.size(); j++) {
 718
 719                                         //are you rare?
 720                                         if (lookup[j]->getAbundance(i) > nseqs) {
 721                                                 newRabunds[j].push_back(lookup[j]->getAbundance(i), newRabunds[j].getGroup());
 722                                                 allZero = false;
 723                                         }else {
 724                                                 newRabunds[j].push_back(0, newRabunds[j].getGroup());
 725                                         }
 726                                 }
 727
 728                                 //eliminates zero otus
 729                                 if (allZero) { for (int j = 0; j < newRabunds.size(); j++) {  newRabunds[j].pop_back(); } }
 730                         }
 731                 }else {
 732                         //for each otu
 733                         for (int i = 0; i < lookup[0]->getNumBins(); i++) {
 734
 735                                 if (m->control_pressed) { return 0; }
 736
 737                                 int totalAbund = 0;
 738                                 //get total otu abundance
 739                                 for (int j = 0; j < lookup.size(); j++) {
 740                                         newRabunds[j].push_back(lookup[j]->getAbundance(i), newRabunds[j].getGroup());
 741                                         totalAbund += lookup[j]->getAbundance(i);
 742                                 }
 743
 744                                 //eliminates otus below rare cutoff
 745                                 if (totalAbund <= nseqs) { for (int j = 0; j < newRabunds.size(); j++) {  newRabunds[j].pop_back(); } }
 746                         }
 747                 }
 748
 749                 //do we have any otus above the rare cutoff
 750                 if (newRabunds[0].getNumBins() != 0) {
 751                         for (int j = 0; j < newRabunds.size(); j++) {
 752                                 out << newRabunds[j].getLabel() << '\t' << newRabunds[j].getGroup() << '\t';
 753                                 newRabunds[j].print(out);
 754                         }
 755                 }
 756
 757                 return 0;
 758         }
 759         catch(exception& e) {
 760                 m->errorOut(e, "RemoveRareCommand", "processLookup");
 761                 exit(1);
 762         }
 763 }
 764 //**********************************************************************************************************************
 765
 766
 767
 768