mothur.h

   1 #ifndef MOTHUR_H
   2 #define MOTHUR_H
   3
   4
   5
   6 /*
   7  *  mothur.h
   8  *  Mothur
   9  *
  10  *  Created by Sarah Westcott on 2/19/09.
  11  *  Copyright 2009 Schloss Lab UMASS Amherst. All rights reserved.
  12  *
  13  */
  14
  15 /* This file contains all the standard incudes we use in the project as well as some common utilities. */
  16
  17 //#include <cstddef>
  18
  19 //io libraries
  20 #include <iostream>
  21 #include <iomanip>
  22 #include <fstream>
  23 #include <sstream>
  24
  25 //exception
  26 #include <stdexcept>
  27 #include <exception>
  28 #include <cstdlib>
  29
  30
  31 //containers
  32 #include <vector>
  33 #include <set>
  34 #include <map>
  35 #include <string>
  36 #include <list>
  37
  38 //math
  39 #include <cmath>
  40 #include <math.h>
  41 #include <algorithm>
  42
  43 //misc
  44 #include <cerrno>
  45 #include <ctime>
  46 #include <limits>
  47 #include <unistd.h>
  48
  49 #if defined (__APPLE__) || (__MACH__) || (linux) || (__linux)
  50         #include <sys/wait.h>
  51 #endif
  52
  53 using namespace std;
  54
  55 #define exp(x) (exp((double) x))
  56 #define sqrt(x) (sqrt((double) x))
  57 #define log10(x) (log10((double) x))
  58 #define log2(x) (log10(x)/log10(2))
  59 #define isnan(x) ((x) != (x))
  60 #define isinf(x) (fabs(x) == std::numeric_limits<double>::infinity())
  61
  62
  63 typedef unsigned long ull;
  64
  65 struct IntNode {
  66         int lvalue;
  67         int rvalue;
  68         int lcoef;
  69         int rcoef;
  70         IntNode* left;
  71         IntNode* right;
  72 };
  73
  74 struct ThreadNode {
  75         int* pid;
  76         IntNode* left;
  77         IntNode* right;
  78 };
  79
  80 /***********************************************************************/
  81
  82 // snagged from http://www.parashift.com/c++-faq-lite/misc-technical-issues.html#faq-39.2
  83 // works for now, but there should be a way to do it without killing the whole program
  84
  85 class BadConversion : public runtime_error {
  86 public:
  87         BadConversion(const string& s) : runtime_error(s){ }
  88 };
  89
  90 //**********************************************************************************************************************
  91
  92 template<typename T>
  93 inline void convert(const string& s, T& x, bool failIfLeftoverChars = true){
  94         istringstream i(s);
  95         char c;
  96         if (!(i >> x) || (failIfLeftoverChars && i.get(c)))
  97                 throw BadConversion(s);
  98 }
  99
 100 //**********************************************************************************************************************
 101
 102 template<typename T>
 103 inline bool convertTestFloat(const string& s, T& x, bool failIfLeftoverChars = true){
 104         istringstream i(s);
 105         char c;
 106         if (!(i >> x) || (failIfLeftoverChars && i.get(c)))
 107         {
 108                 return false;
 109         }
 110         return true;
 111 }
 112
 113 //**********************************************************************************************************************
 114
 115 template<typename T>
 116 inline bool convertTest(const string& s, T& x, bool failIfLeftoverChars = true){
 117         istringstream i(s);
 118         char c;
 119         if (!(i >> x) || (failIfLeftoverChars && i.get(c)))
 120         {
 121                 cout << "unable to be converted into an integer.\n" << endl;
 122                 return false;
 123         }
 124         return true;
 125 }
 126
 127 //**********************************************************************************************************************
 128
 129 template<typename T>
 130 string toString(const T&x){
 131     stringstream output;
 132     output << x;
 133     return output.str();
 134 }
 135
 136 //**********************************************************************************************************************
 137
 138 template<typename T>
 139 string toHex(const T&x){
 140         stringstream output;
 141
 142         output << hex << x;
 143
 144     return output.str();
 145 }
 146 //**********************************************************************************************************************
 147
 148 template<typename T>
 149 string toString(const T&x, int i){
 150         stringstream output;
 151
 152         output.precision(i);
 153     output << fixed << x;
 154
 155     return output.str();
 156 }
 157 /***********************************************************************/
 158
 159 inline int openOutputFileAppend(string fileName, ofstream& fileHandle){
 160
 161         fileHandle.open(fileName.c_str(), ios::app);
 162         if(!fileHandle) {
 163                 cout << "Error: Could not open " <<  fileName << endl;
 164                 return 1;
 165         }
 166         else {
 167                 return 0;
 168         }
 169
 170 }
 171 /***********************************************************************/
 172
 173 inline string getline(ifstream& fileHandle) {
 174         try {
 175
 176                 string line = "";
 177
 178                 while (!fileHandle.eof())       {
 179                         //get next character
 180                         char c = fileHandle.get();
 181
 182                         //are you at the end of the line
 183                         if (c == 10 || c == 13){   break;       }
 184                         else {          line += c;              }
 185                 }
 186
 187                 return line;
 188
 189         }
 190         catch(exception& e) {
 191                 cout << "Error in mothur function getline" << endl;
 192                 exit(1);
 193         }
 194 }
 195
 196 /**************************************************************************************************/
 197
 198 inline void mothurOut(string message) {
 199         try{
 200                 ofstream out;
 201                 string logFileName = "mothur.logFile";
 202                 openOutputFileAppend(logFileName, out);
 203
 204                 cout << message;
 205                 out << message;
 206
 207                 out.close();
 208         }
 209         catch(exception& e) {
 210                 cout << "Error in mothur class mothurOut" << endl;
 211                 exit(1);
 212         }
 213 }
 214 /**************************************************************************************************/
 215
 216 inline void mothurOut(string message, string precision) {
 217         try{
 218                 ofstream out;
 219                 string logFileName = "mothur.logFile";
 220                 openOutputFileAppend(logFileName, out);
 221
 222                 cout << precision << message;
 223                 out << precision << message;
 224
 225                 out.close();
 226         }
 227         catch(exception& e) {
 228                 cout << "Error in mothur class mothurOut" << endl;
 229                 exit(1);
 230         }
 231 }
 232
 233 /**************************************************************************************************/
 234
 235 inline void mothurOutEndLine() {
 236         try {
 237                 ofstream out;
 238                 string logFileName = "mothur.logFile";
 239                 openOutputFileAppend(logFileName, out);
 240
 241                 cout << endl;
 242                 out << endl;
 243
 244                 out.close();
 245         }
 246         catch(exception& e) {
 247                 cout << "error in mothur mothurOutEndLine" << endl;
 248                 exit(1);
 249         }
 250 }
 251
 252
 253 /**************************************************************************************************/
 254
 255 inline void errorOut(exception& e, string object, string function) {
 256
 257                 mothurOut("Error: ");
 258                 mothurOut(toString(e.what()));
 259                 mothurOut(" has occurred in the " + object + " class function " + function + ". Please contact Pat Schloss at pschloss@microbio.umass.edu, and be sure to include the mothur.logFile with your inquiry.");
 260                 mothurOutEndLine();
 261
 262 }
 263
 264
 265
 266
 267 /***********************************************************************/
 268
 269 inline void gobble(istream& f){
 270
 271         char d;
 272     while(isspace(d=f.get()))           {;}
 273         f.putback(d);
 274
 275 }
 276 /***********************************************************************/
 277
 278 inline bool isTrue(string f){
 279
 280         if ((f == "TRUE") || (f == "T") || (f == "true") || (f == "t")) {       return true;    }
 281         else {  return false;  }
 282 }
 283
 284 /***********************************************************************/
 285
 286 inline float roundDist(float dist, int precision){
 287
 288         return int(dist * precision + 0.5)/float(precision);
 289
 290 }
 291
 292 /***********************************************************************/
 293
 294 inline int getNumNames(string names){
 295
 296         int count = 0;
 297
 298         if(names != ""){
 299                 count = 1;
 300                 for(int i=0;i<names.size();i++){
 301                         if(names[i] == ','){
 302                                 count++;
 303                         }
 304                 }
 305         }
 306
 307         return count;
 308
 309 }
 310
 311 /**************************************************************************************************/
 312
 313 inline vector<vector<double> > binomial(int maxOrder){
 314
 315         vector<vector<double> > binomial(maxOrder+1);
 316
 317     for(int i=0;i<=maxOrder;i++){
 318                 binomial[i].resize(maxOrder+1);
 319                 binomial[i][0]=1;
 320                 binomial[0][i]=0;
 321     }
 322     binomial[0][0]=1;
 323
 324     binomial[1][0]=1;
 325     binomial[1][1]=1;
 326
 327     for(int i=2;i<=maxOrder;i++){
 328                 binomial[1][i]=0;
 329     }
 330
 331     for(int i=2;i<=maxOrder;i++){
 332                 for(int j=1;j<=maxOrder;j++){
 333                         if(i==j){       binomial[i][j]=1;                                                                       }
 334                         if(j>i) {       binomial[i][j]=0;                                                                       }
 335                         else    {       binomial[i][j]=binomial[i-1][j-1]+binomial[i-1][j];     }
 336                 }
 337     }
 338
 339         return binomial;
 340 }
 341
 342 /***********************************************************************/
 343
 344 inline string getRootName(string longName){
 345
 346         string rootName = longName;
 347
 348         if(longName.find_last_of(".") != longName.npos){
 349                 int pos = longName.find_last_of('.')+1;
 350                 rootName = longName.substr(0, pos);
 351         }
 352
 353         return rootName;
 354 }
 355 /***********************************************************************/
 356
 357 inline string getSimpleName(string longName){
 358
 359         string simpleName = longName;
 360
 361         if(longName.find_last_of("/") != longName.npos){
 362                 int pos = longName.find_last_of('/')+1;
 363                 simpleName = longName.substr(pos, longName.length());
 364         }
 365
 366         return simpleName;
 367 }
 368
 369 /***********************************************************************/
 370
 371 inline int factorial(int num){
 372         int total = 1;
 373
 374         for (int i = 1; i <= num; i++) {
 375                 total *= i;
 376         }
 377
 378         return total;
 379 }
 380 /**************************************************************************************************
 381
 382 double min(double x, double y)
 383 {
 384     if(x<y){    return x;    }
 385     else   {    return y;    }
 386 }
 387
 388 /***********************************************************************/
 389
 390 inline string getPathName(string longName){
 391
 392         string rootPathName = longName;
 393
 394         if(longName.find_last_of('/') != longName.npos){
 395                 int pos = longName.find_last_of('/')+1;
 396                 rootPathName = longName.substr(0, pos);
 397         }
 398
 399         return rootPathName;
 400 }
 401
 402 /***********************************************************************/
 403
 404 inline string getExtension(string longName){
 405
 406         string extension = longName;
 407
 408         if(longName.find_last_of('.') != longName.npos){
 409                 int pos = longName.find_last_of('.');
 410                 extension = longName.substr(pos, longName.length());
 411         }
 412
 413         return extension;
 414 }
 415
 416 /***********************************************************************/
 417
 418 inline int openInputFile(string fileName, ifstream& fileHandle){
 419
 420         fileHandle.open(fileName.c_str());
 421         if(!fileHandle) {
 422                 mothurOut("Error: Could not open " + fileName);  mothurOutEndLine();
 423                 return 1;
 424         }
 425         else {
 426                 return 0;
 427         }
 428
 429 }
 430
 431 /***********************************************************************/
 432
 433 inline int openOutputFile(string fileName, ofstream& fileHandle){
 434
 435         fileHandle.open(fileName.c_str(), ios::trunc);
 436         if(!fileHandle) {
 437                 mothurOut("Error: Could not open " + fileName);  mothurOutEndLine();
 438                 return 1;
 439         }
 440         else {
 441                 return 0;
 442         }
 443
 444 }
 445
 446 /***********************************************************************/
 447
 448 inline int getNumSeqs(ifstream& file){
 449
 450         int numSeqs = count(istreambuf_iterator<char>(file),istreambuf_iterator<char>(), '>');
 451         file.seekg(0);
 452         return numSeqs;
 453
 454 }
 455
 456 /***********************************************************************/
 457
 458 //This function parses the estimator options and puts them in a vector
 459 inline void splitAtDash(string& estim, vector<string>& container) {
 460         try {
 461                 string individual;
 462
 463                 while (estim.find_first_of('-') != -1) {
 464                         individual = estim.substr(0,estim.find_first_of('-'));
 465                         if ((estim.find_first_of('-')+1) <= estim.length()) { //checks to make sure you don't have dash at end of string
 466                                 estim = estim.substr(estim.find_first_of('-')+1, estim.length());
 467                                 container.push_back(individual);
 468                         }
 469                 }
 470                 //get last one
 471                 container.push_back(estim);
 472         }
 473         catch(exception& e) {
 474                 errorOut(e, "mothur", "splitAtDash");
 475                 exit(1);
 476         }
 477 }
 478
 479 /***********************************************************************/
 480 //This function parses the label options and puts them in a set
 481 inline void splitAtDash(string& estim, set<string>& container) {
 482         try {
 483                 string individual;
 484
 485                 while (estim.find_first_of('-') != -1) {
 486                         individual = estim.substr(0,estim.find_first_of('-'));
 487                         if ((estim.find_first_of('-')+1) <= estim.length()) { //checks to make sure you don't have dash at end of string
 488                                 estim = estim.substr(estim.find_first_of('-')+1, estim.length());
 489                                 container.insert(individual);
 490                         }
 491                 }
 492                 //get last one
 493                 container.insert(estim);
 494         }
 495         catch(exception& e) {
 496                 errorOut(e, "mothur", "splitAtDash");
 497                 exit(1);
 498         }
 499 }
 500 /***********************************************************************/
 501 //This function parses the line options and puts them in a set
 502 inline void splitAtDash(string& estim, set<int>& container) {
 503         try {
 504                 string individual;
 505                 int lineNum;
 506
 507                 while (estim.find_first_of('-') != -1) {
 508                         individual = estim.substr(0,estim.find_first_of('-'));
 509                         if ((estim.find_first_of('-')+1) <= estim.length()) { //checks to make sure you don't have dash at end of string
 510                                 estim = estim.substr(estim.find_first_of('-')+1, estim.length());
 511                                 convert(individual, lineNum); //convert the string to int
 512                                 container.insert(lineNum);
 513                         }
 514                 }
 515                 //get last one
 516                 convert(estim, lineNum); //convert the string to int
 517                 container.insert(lineNum);
 518         }
 519         catch(exception& e) {
 520                 errorOut(e, "mothur", "splitAtDash");
 521                 exit(1);
 522         }
 523 }
 524 /***********************************************************************/
 525 //This function parses the a string and puts peices in a vector
 526 inline void splitAtComma(string& estim, vector<string>& container) {
 527         try {
 528                 string individual;
 529
 530                 while (estim.find_first_of(',') != -1) {
 531                         individual = estim.substr(0,estim.find_first_of(','));
 532                         if ((estim.find_first_of(',')+1) <= estim.length()) { //checks to make sure you don't have comma at end of string
 533                                 estim = estim.substr(estim.find_first_of(',')+1, estim.length());
 534                                 container.push_back(individual);
 535                         }
 536                 }
 537                 //get last one
 538                 container.push_back(estim);
 539         }
 540         catch(exception& e) {
 541                 errorOut(e, "mothur", "splitAtComma");
 542                 exit(1);
 543         }
 544 }
 545 /***********************************************************************/
 546
 547 //This function splits up the various option parameters
 548 inline void splitAtComma(string& prefix, string& suffix){
 549         try {
 550                 prefix = suffix.substr(0,suffix.find_first_of(','));
 551                 if ((suffix.find_first_of(',')+2) <= suffix.length()) {  //checks to make sure you don't have comma at end of string
 552                         suffix = suffix.substr(suffix.find_first_of(',')+1, suffix.length());
 553                         string space = " ";
 554                         while(suffix.at(0) == ' ')
 555                                 suffix = suffix.substr(1, suffix.length());
 556                 }
 557
 558         }
 559         catch(exception& e) {
 560                 errorOut(e, "mothur", "splitAtComma");
 561                 exit(1);
 562         }
 563 }
 564 /***********************************************************************/
 565
 566 //This function separates the key value from the option value i.e. dist=96_...
 567 inline void splitAtEquals(string& key, string& value){
 568         try {
 569                 if(value.find_first_of('=') != -1){
 570                         key = value.substr(0,value.find_first_of('='));
 571                         if ((value.find_first_of('=')+1) <= value.length()) {
 572                                 value = value.substr(value.find_first_of('=')+1, value.length());
 573                         }
 574                 }else{
 575                         key = value;
 576                         value = 1;
 577                 }
 578         }
 579         catch(exception& e) {
 580                 errorOut(e, "mothur", "splitAtEquals");
 581                 exit(1);
 582         }
 583 }
 584 /**************************************************************************************************/
 585
 586 inline bool inUsersGroups(string groupname, vector<string> Groups) {
 587         try {
 588                 for (int i = 0; i < Groups.size(); i++) {
 589                         if (groupname == Groups[i]) { return true; }
 590                 }
 591                 return false;
 592         }
 593         catch(exception& e) {
 594                 errorOut(e, "mothur", "inUsersGroups");
 595                 exit(1);
 596         }
 597 }
 598
 599 /**************************************************************************************************/
 600
 601 inline void mothurOutJustToLog(string message) {
 602         try {
 603                 ofstream out;
 604                 string logFileName = "mothur.logFile";
 605                 openOutputFileAppend(logFileName, out);
 606
 607                 out << message;
 608
 609                 out.close();
 610         }
 611         catch(exception& e) {
 612                 errorOut(e, "mothur", "mothurOutJustToLog");
 613                 exit(1);
 614         }
 615 }
 616
 617
 618 /**************************************************************************************************/
 619
 620 inline void mothurOut(float num) {
 621         try {
 622                 ofstream out;
 623                 string logFileName = "mothur.logFile";
 624                 openOutputFileAppend(logFileName, out);
 625
 626                 cout << num;
 627                 out << num;
 628
 629                 out.close();
 630         }
 631         catch(exception& e) {
 632                 cout << "Error in mothur class mothurOut float" << endl;
 633                 exit(1);
 634         }
 635 }
 636 /***********************************************************************/
 637 inline void mothurOut(double value) {
 638         try {
 639                 ofstream out;
 640                 string logFileName = "mothur.logFile";
 641                 openOutputFileAppend(logFileName, out);
 642
 643                 cout << value;
 644                 out << value;
 645
 646                 out.close();
 647         }
 648         catch(exception& e) {
 649                 cout << "Error in mothur class mothurOut double" << endl;
 650                 exit(1);
 651         }
 652 }
 653
 654 /***********************************************************************/
 655 //this function determines if the user has given us labels that are smaller than the given label.
 656 //if so then it returns true so that the calling function can run the previous valid distance.
 657 //it's a "smart" distance function.  It also checks for invalid labels.
 658 inline bool anyLabelsToProcess(string label, set<string>& userLabels, string errorOff) {
 659         try {
 660                 set<string>::iterator it;
 661                 vector<float> orderFloat;
 662                 map<string, float> userMap;  //the conversion process removes trailing 0's which we need to put back
 663                 map<string, float>::iterator it2;
 664                 float labelFloat;
 665                 bool smaller = false;
 666
 667                 //unique is the smallest line
 668                 if (label == "unique") {  return false;  }
 669                 else { convert(label, labelFloat); }
 670
 671                 //go through users set and make them floats
 672                 for(it = userLabels.begin(); it != userLabels.end(); ++it) {
 673
 674                         float temp;
 675                         if ((*it != "unique") && (convertTestFloat(*it, temp) == true)){
 676                                 convert(*it, temp);
 677                                 orderFloat.push_back(temp);
 678                                 userMap[*it] = temp;
 679                         }else if (*it == "unique") {
 680                                 orderFloat.push_back(-1.0);
 681                                 userMap["unique"] = -1.0;
 682                         }else {
 683                                 if (errorOff == "") {  mothurOut(*it + " is not a valid label."); mothurOutEndLine();  }
 684                                 userLabels.erase(*it);
 685                                 it--;
 686                         }
 687                 }
 688
 689                 //sort order
 690                 sort(orderFloat.begin(), orderFloat.end());
 691
 692                 /*************************************************/
 693                 //is this label bigger than any of the users labels
 694                 /*************************************************/
 695
 696                 //loop through order until you find a label greater than label
 697                 for (int i = 0; i < orderFloat.size(); i++) {
 698                         if (orderFloat[i] < labelFloat) {
 699                                 smaller = true;
 700                                 if (orderFloat[i] == -1) {
 701                                         if (errorOff == "") { mothurOut("Your file does not include the label unique."); mothurOutEndLine(); }
 702                                         userLabels.erase("unique");
 703                                 }
 704                                 else {
 705                                         if (errorOff == "") { mothurOut("Your file does not include the label "); mothurOutEndLine(); }
 706                                         string s = "";
 707                                         for (it2 = userMap.begin(); it2!= userMap.end(); it2++) {
 708                                                 if (it2->second == orderFloat[i]) {
 709                                                         s = it2->first;
 710                                                         //remove small labels
 711                                                         userLabels.erase(s);
 712                                                         break;
 713                                                 }
 714                                         }
 715                                         if (errorOff == "") { mothurOut(s + ". I will use the next smallest distance. "); mothurOutEndLine(); }
 716                                 }
 717                         //since they are sorted once you find a bigger one stop looking
 718                         }else { break; }
 719                 }
 720
 721                 return smaller;
 722
 723         }
 724         catch(exception& e) {
 725                 errorOut(e, "mothur", "anyLabelsToProcess");
 726                 exit(1);
 727         }
 728 }
 729
 730 /**************************************************************************************************/
 731 #endif
 732