rarefactcommand.cpp

   1 /*
   2  *  rarefactcommand.cpp
   3  *  Dotur
   4  *
   5  *  Created by Sarah Westcott on 1/2/09.
   6  *  Copyright 2009 Schloss Lab UMASS Amherst. All rights reserved.
   7  *
   8  */
   9
  10 #include "rarefactcommand.h"
  11 #include "ace.h"
  12 #include "sobs.h"
  13 #include "nseqs.h"
  14 #include "chao1.h"
  15 #include "bootstrap.h"
  16 #include "simpson.h"
  17 #include "npshannon.h"
  18 #include "shannon.h"
  19 #include "jackknife.h"
  20 #include "coverage.h"
  21
  22 //**********************************************************************************************************************
  23
  24
  25 RareFactCommand::RareFactCommand(){
  26         try {
  27                 globaldata = GlobalData::getInstance();
  28                 string fileNameRoot;
  29                 fileNameRoot = getRootName(globaldata->inputFileName);
  30                 convert(globaldata->getFreq(), freq);
  31                 convert(globaldata->getIters(), nIters);
  32                 validCalculator = new ValidCalculators();
  33
  34                 int i;
  35                 for (i=0; i<globaldata->Estimators.size(); i++) {
  36                         if (validCalculator->isValidCalculator("rarefaction", globaldata->Estimators[i]) == true) {
  37                                 if (globaldata->Estimators[i] == "sobs") {
  38                                         rDisplays.push_back(new RareDisplay(new Sobs(), new ThreeColumnFile(fileNameRoot+"rarefaction")));
  39                                 }else if (globaldata->Estimators[i] == "chao") {
  40                                         rDisplays.push_back(new RareDisplay(new Chao1(), new ThreeColumnFile(fileNameRoot+"r_chao")));
  41                                 }else if (globaldata->Estimators[i] == "ace") {
  42                                         convert(globaldata->getAbund(), abund);
  43                                         if(abund < 5)
  44                                                 abund = 10;
  45                                         rDisplays.push_back(new RareDisplay(new Ace(abund), new ThreeColumnFile(fileNameRoot+"r_ace")));
  46                                 }else if (globaldata->Estimators[i] == "jack") {
  47                                         rDisplays.push_back(new RareDisplay(new Jackknife(), new ThreeColumnFile(fileNameRoot+"r_jack")));
  48                                 }else if (globaldata->Estimators[i] == "shannon") {
  49                                         rDisplays.push_back(new RareDisplay(new Shannon(), new ThreeColumnFile(fileNameRoot+"r_shannon")));
  50                                 }else if (globaldata->Estimators[i] == "npshannon") {
  51                                         rDisplays.push_back(new RareDisplay(new NPShannon(), new ThreeColumnFile(fileNameRoot+"r_npshannon")));
  52                                 }else if (globaldata->Estimators[i] == "simpson") {
  53                                         rDisplays.push_back(new RareDisplay(new Simpson(), new ThreeColumnFile(fileNameRoot+"r_simpson")));
  54                                 }else if (globaldata->Estimators[i] == "bootstrap") {
  55                                         rDisplays.push_back(new RareDisplay(new Bootstrap(), new ThreeColumnFile(fileNameRoot+"r_bootstrap")));
  56                                 }else if (globaldata->Estimators[i] == "coverage") {
  57                                         rDisplays.push_back(new RareDisplay(new Coverage(), new ThreeColumnFile(fileNameRoot+"r_coverage")));
  58                                 }else if (globaldata->Estimators[i] == "nseqs") {
  59                                         rDisplays.push_back(new RareDisplay(new NSeqs(), new ThreeColumnFile(fileNameRoot+"r_nseqs")));
  60                                 }
  61                         }
  62                 }
  63
  64                 //reset calc for next command
  65                 globaldata->setCalc("");
  66
  67         }
  68         catch(exception& e) {
  69                 cout << "Standard Error: " << e.what() << " has occurred in the RareFactCommand class Function RareFactCommand. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";
  70                 exit(1);
  71         }
  72         catch(...) {
  73                 cout << "An unknown error has occurred in the RareFactCommand class function RareFactCommand. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";
  74                 exit(1);
  75         }
  76
  77 }
  78
  79 //**********************************************************************************************************************
  80
  81 RareFactCommand::~RareFactCommand(){
  82         delete order;
  83         delete input;
  84         delete rCurve;
  85         delete read;
  86 }
  87
  88 //**********************************************************************************************************************
  89
  90 int RareFactCommand::execute(){
  91         try {
  92                 int count = 1;
  93
  94                 //if the users entered no valid calculators don't execute command
  95                 if (rDisplays.size() == 0) { return 0; }
  96
  97                 read = new ReadOTUFile(globaldata->inputFileName);
  98                 read->read(&*globaldata);
  99
 100                 order = globaldata->gorder;
 101                 OrderVector* lastOrder = order;
 102                 input = globaldata->ginput;
 103
 104                 //if the users enters label "0.06" and there is no "0.06" in their file use the next lowest label.
 105                 set<string> processedLabels;
 106                 set<string> userLabels = globaldata->labels;
 107
 108                 //as long as you are not at the end of the file or done wih the lines you want
 109                 while((order != NULL) && ((globaldata->allLines == 1) || (userLabels.size() != 0))) {
 110
 111                         if(globaldata->allLines == 1 || globaldata->lines.count(count) == 1 || globaldata->labels.count(order->getLabel()) == 1){
 112
 113                                 rCurve = new Rarefact(order, rDisplays);
 114                                 rCurve->getCurve(freq, nIters);
 115                                 delete rCurve;
 116
 117                                 cout << order->getLabel() << '\t' << count << endl;
 118                                 processedLabels.insert(order->getLabel());
 119                                 userLabels.erase(order->getLabel());
 120                         }
 121
 122                         if ((anyLabelsToProcess(order->getLabel(), userLabels, "") == true) && (processedLabels.count(lastOrder->getLabel()) != 1)) {
 123                                 rCurve = new Rarefact(lastOrder, rDisplays);
 124                                 rCurve->getCurve(freq, nIters);
 125                                 delete rCurve;
 126
 127                                 cout << lastOrder->getLabel() << '\t' << count << endl;
 128                                 processedLabels.insert(lastOrder->getLabel());
 129                                 userLabels.erase(lastOrder->getLabel());
 130                         }
 131
 132                         if (count != 1) { delete lastOrder; }
 133                         lastOrder = order;
 134
 135                         order = (input->getOrderVector());
 136                         count++;
 137                 }
 138
 139                 //output error messages about any remaining user labels
 140                 set<string>::iterator it;
 141                 bool needToRun = false;
 142                 for (it = userLabels.begin(); it != userLabels.end(); it++) {
 143                         cout << "Your file does not include the label "<< *it;
 144                         if (processedLabels.count(lastOrder->getLabel()) != 1) {
 145                                 cout << ". I will use " << lastOrder->getLabel() << "." << endl;
 146                                 needToRun = true;
 147                         }else {
 148                                 cout << ". Please refer to " << lastOrder->getLabel() << "." << endl;
 149                         }
 150                 }
 151
 152                 //run last line if you need to
 153                 if (needToRun == true)  {
 154                         rCurve = new Rarefact(lastOrder, rDisplays);
 155                         rCurve->getCurve(freq, nIters);
 156                         delete rCurve;
 157
 158                         cout << lastOrder->getLabel() << '\t' << count << endl;
 159                 }
 160
 161                 delete lastOrder;
 162
 163                 for(int i=0;i<rDisplays.size();i++){    delete rDisplays[i];    }
 164                 return 0;
 165         }
 166         catch(exception& e) {
 167                 cout << "Standard Error: " << e.what() << " has occurred in the RareFactCommand class Function execute. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";
 168                 exit(1);
 169         }
 170         catch(...) {
 171                 cout << "An unknown error has occurred in the RareFactCommand class function execute. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";
 172                 exit(1);
 173         }
 174 }
 175
 176 //**********************************************************************************************************************