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