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