mothurout.cpp

   1 /*
   2  *  mothurOut.cpp
   3  *  Mothur
   4  *
   5  *  Created by westcott on 2/25/10.
   6  *  Copyright 2010 Schloss Lab. All rights reserved.
   7  *
   8  */
   9
  10 #include "mothurout.h"
  11
  12
  13 /******************************************************/
  14 MothurOut* MothurOut::getInstance() {
  15         if( _uniqueInstance == 0) {
  16                 _uniqueInstance = new MothurOut();
  17         }
  18         return _uniqueInstance;
  19 }
  20 /*********************************************************************************************/
  21 set<string> MothurOut::getCurrentTypes()  {
  22         try {
  23
  24         set<string> types;
  25         types.insert("fasta");
  26         types.insert("summary");
  27         types.insert("accnos");
  28         types.insert("column");
  29         types.insert("design");
  30         types.insert("group");
  31         types.insert("list");
  32         types.insert("name");
  33         types.insert("oligos");
  34         types.insert("order");
  35         types.insert("ordergroup");
  36         types.insert("phylip");
  37         types.insert("qfile");
  38         types.insert("relabund");
  39         types.insert("sabund");
  40         types.insert("rabund");
  41         types.insert("sff");
  42         types.insert("shared");
  43         types.insert("taxonomy");
  44         types.insert("tree");
  45         types.insert("flow");
  46         types.insert("biom");
  47         types.insert("count");
  48         types.insert("processors");
  49
  50                 return types;
  51         }
  52         catch(exception& e) {
  53                 errorOut(e, "MothurOut", "getCurrentTypes");
  54                 exit(1);
  55         }
  56 }
  57 /*********************************************************************************************/
  58 void MothurOut::printCurrentFiles()  {
  59         try {
  60
  61
  62                 if (accnosfile != "")           {  mothurOut("accnos=" + accnosfile); mothurOutEndLine();                       }
  63                 if (columnfile != "")           {  mothurOut("column=" + columnfile); mothurOutEndLine();                       }
  64                 if (designfile != "")           {  mothurOut("design=" + designfile); mothurOutEndLine();                       }
  65                 if (fastafile != "")            {  mothurOut("fasta=" + fastafile); mothurOutEndLine();                         }
  66                 if (groupfile != "")            {  mothurOut("group=" + groupfile); mothurOutEndLine();                         }
  67                 if (listfile != "")                     {  mothurOut("list=" + listfile); mothurOutEndLine();                           }
  68                 if (namefile != "")                     {  mothurOut("name=" + namefile); mothurOutEndLine();                           }
  69                 if (oligosfile != "")           {  mothurOut("oligos=" + oligosfile); mothurOutEndLine();                       }
  70                 if (orderfile != "")            {  mothurOut("order=" + orderfile); mothurOutEndLine();                         }
  71                 if (ordergroupfile != "")       {  mothurOut("ordergroup=" + ordergroupfile); mothurOutEndLine();       }
  72                 if (phylipfile != "")           {  mothurOut("phylip=" + phylipfile); mothurOutEndLine();                       }
  73                 if (qualfile != "")                     {  mothurOut("qfile=" + qualfile); mothurOutEndLine();                          }
  74                 if (rabundfile != "")           {  mothurOut("rabund=" + rabundfile); mothurOutEndLine();                       }
  75                 if (relabundfile != "")         {  mothurOut("relabund=" + relabundfile); mothurOutEndLine();           }
  76                 if (sabundfile != "")           {  mothurOut("sabund=" + sabundfile); mothurOutEndLine();                       }
  77                 if (sfffile != "")                      {  mothurOut("sff=" + sfffile); mothurOutEndLine();                                     }
  78                 if (sharedfile != "")           {  mothurOut("shared=" + sharedfile); mothurOutEndLine();                       }
  79                 if (taxonomyfile != "")         {  mothurOut("taxonomy=" + taxonomyfile); mothurOutEndLine();           }
  80                 if (treefile != "")                     {  mothurOut("tree=" + treefile); mothurOutEndLine();                           }
  81                 if (flowfile != "")                     {  mothurOut("flow=" + flowfile); mothurOutEndLine();                           }
  82         if (biomfile != "")                     {  mothurOut("biom=" + biomfile); mothurOutEndLine();                           }
  83         if (counttablefile != "")       {  mothurOut("count=" + counttablefile); mothurOutEndLine();    }
  84                 if (processors != "1")          {  mothurOut("processors=" + processors); mothurOutEndLine();           }
  85         if (summaryfile != "")          {  mothurOut("summary=" + summaryfile); mothurOutEndLine();             }
  86
  87         }
  88         catch(exception& e) {
  89                 errorOut(e, "MothurOut", "printCurrentFiles");
  90                 exit(1);
  91         }
  92 }
  93 /*********************************************************************************************/
  94 bool MothurOut::hasCurrentFiles()  {
  95         try {
  96                 bool hasCurrent = false;
  97
  98                 if (accnosfile != "")           {  return true;                 }
  99                 if (columnfile != "")           {  return true;                 }
 100                 if (designfile != "")           {  return true;                 }
 101                 if (fastafile != "")            {  return true;                 }
 102                 if (groupfile != "")            {  return true;                 }
 103                 if (listfile != "")                     {  return true;                 }
 104                 if (namefile != "")                     {  return true;                 }
 105                 if (oligosfile != "")           {  return true;                 }
 106                 if (orderfile != "")            {  return true;                 }
 107                 if (ordergroupfile != "")       {  return true;                 }
 108                 if (phylipfile != "")           {  return true;                 }
 109                 if (qualfile != "")                     {  return true;                 }
 110                 if (rabundfile != "")           {  return true;                 }
 111                 if (relabundfile != "")         {  return true;                 }
 112                 if (sabundfile != "")           {  return true;                 }
 113                 if (sfffile != "")                      {  return true;                 }
 114                 if (sharedfile != "")           {  return true;                 }
 115                 if (taxonomyfile != "")         {  return true;                 }
 116                 if (treefile != "")                     {  return true;                 }
 117                 if (flowfile != "")                     {  return true;                 }
 118         if (biomfile != "")                     {  return true;                 }
 119         if (counttablefile != "")       {  return true;                 }
 120         if (summaryfile != "")  {  return true;                 }
 121                 if (processors != "1")          {  return true;                 }
 122
 123                 return hasCurrent;
 124
 125         }
 126         catch(exception& e) {
 127                 errorOut(e, "MothurOut", "hasCurrentFiles");
 128                 exit(1);
 129         }
 130 }
 131
 132 /*********************************************************************************************/
 133 void MothurOut::clearCurrentFiles()  {
 134         try {
 135                 phylipfile = "";
 136                 columnfile = "";
 137                 listfile = "";
 138                 rabundfile = "";
 139                 sabundfile = "";
 140                 namefile = "";
 141                 groupfile = "";
 142                 designfile = "";
 143                 orderfile = "";
 144                 treefile = "";
 145                 sharedfile = "";
 146                 ordergroupfile = "";
 147                 relabundfile = "";
 148                 fastafile = "";
 149                 qualfile = "";
 150                 sfffile = "";
 151                 oligosfile = "";
 152                 accnosfile = "";
 153                 taxonomyfile = "";
 154                 flowfile = "";
 155         biomfile = "";
 156         counttablefile = "";
 157         summaryfile = "";
 158                 processors = "1";
 159         }
 160         catch(exception& e) {
 161                 errorOut(e, "MothurOut", "clearCurrentFiles");
 162                 exit(1);
 163         }
 164 }
 165 /***********************************************************************/
 166 string MothurOut::findProgramPath(string programName){
 167         try {
 168
 169                 string envPath = getenv("PATH");
 170                 string pPath = "";
 171
 172                 //delimiting path char
 173                 char delim;
 174 #if defined (__APPLE__) || (__MACH__) || (linux) || (__linux) || (__linux__) || (__unix__) || (__unix)
 175         delim = ':';
 176 #else
 177         delim = ';';
 178 #endif
 179
 180                 //break apart path variable by ':'
 181                 vector<string> dirs;
 182                 splitAtChar(envPath, dirs, delim);
 183
 184         if (debug) { mothurOut("[DEBUG]: dir's in path: \n"); }
 185
 186                 //get path related to mothur
 187                 for (int i = 0; i < dirs.size(); i++) {
 188
 189             if (debug) { mothurOut("[DEBUG]: " + dirs[i] + "\n"); }
 190
 191                         //to lower so we can find it
 192                         string tempLower = "";
 193                         for (int j = 0; j < dirs[i].length(); j++) {  tempLower += tolower(dirs[i][j]);  }
 194
 195                         //is this mothurs path?
 196                         if (tempLower.find(programName) != -1) {  pPath = dirs[i]; break;  }
 197                 }
 198
 199                 if (debug) { mothurOut("[DEBUG]: programPath = " + pPath + "\n"); }
 200
 201                 if (pPath != "") {
 202                         //add programName so it looks like what argv would look like
 203 #if defined (__APPLE__) || (__MACH__) || (linux) || (__linux) || (__linux__) || (__unix__) || (__unix)
 204             pPath += "/" + programName;
 205 #else
 206             pPath += "\\" + programName;
 207 #endif
 208                 }else {
 209                         //okay programName is not in the path, so the folder programName is in must be in the path
 210                         //lets find out which one
 211
 212                         //get path related to the program
 213                         for (int i = 0; i < dirs.size(); i++) {
 214
 215                 if (debug) { mothurOut("[DEBUG]: looking in " + dirs[i] + " for " + programName + " \n"); }
 216
 217                                 //is this the programs path?
 218                                 ifstream in;
 219                                 string tempIn = dirs[i];
 220 #if defined (__APPLE__) || (__MACH__) || (linux) || (__linux) || (__linux__) || (__unix__) || (__unix)
 221                 tempIn += "/" + programName;
 222 #else
 223                 tempIn += "\\" + programName;
 224 #endif
 225                                 openInputFile(tempIn, in, "");
 226
 227                                 //if this file exists
 228                                 if (in) { in.close(); pPath = tempIn; if (debug) { mothurOut("[DEBUG]: found it, programPath = " + pPath + "\n"); } break;   }
 229                         }
 230                 }
 231
 232                 return pPath;
 233
 234         }
 235         catch(exception& e) {
 236                 errorOut(e, "MothurOut", "findProgramPath");
 237                 exit(1);
 238         }
 239 }
 240 /*********************************************************************************************/
 241 void MothurOut::setFileName(string filename)  {
 242         try {
 243                 logFileName = filename;
 244
 245                 #ifdef USE_MPI
 246                         int pid;
 247                         MPI_Comm_rank(MPI_COMM_WORLD, &pid);
 248
 249                         if (pid == 0) { //only one process should output to screen
 250                 #endif
 251
 252                 openOutputFile(filename, out);
 253
 254                 #ifdef USE_MPI
 255                         }
 256                 #endif
 257         }
 258         catch(exception& e) {
 259                 errorOut(e, "MothurOut", "setFileName");
 260                 exit(1);
 261         }
 262 }
 263 /*********************************************************************************************/
 264 void MothurOut::setDefaultPath(string pathname)  {
 265         try {
 266
 267                 //add / to name if needed
 268                 string lastChar = pathname.substr(pathname.length()-1);
 269                 #if defined (__APPLE__) || (__MACH__) || (linux) || (__linux) || (__linux__) || (__unix__) || (__unix)
 270                         if (lastChar != "/") { pathname += "/"; }
 271                 #else
 272                         if (lastChar != "\\") { pathname += "\\"; }
 273                 #endif
 274
 275                 defaultPath = pathname;
 276
 277         }
 278         catch(exception& e) {
 279                 errorOut(e, "MothurOut", "setDefaultPath");
 280                 exit(1);
 281         }
 282 }
 283 /*********************************************************************************************/
 284 void MothurOut::setOutputDir(string pathname)  {
 285         try {
 286                 outputDir = pathname;
 287         }
 288         catch(exception& e) {
 289                 errorOut(e, "MothurOut", "setOutputDir");
 290                 exit(1);
 291         }
 292 }
 293 /*********************************************************************************************/
 294 void MothurOut::closeLog()  {
 295         try {
 296
 297                 #ifdef USE_MPI
 298                         int pid;
 299                         MPI_Comm_rank(MPI_COMM_WORLD, &pid);
 300
 301                         if (pid == 0) { //only one process should output to screen
 302                 #endif
 303
 304                 out.close();
 305
 306                 #ifdef USE_MPI
 307                         }
 308                 #endif
 309         }
 310         catch(exception& e) {
 311                 errorOut(e, "MothurOut", "closeLog");
 312                 exit(1);
 313         }
 314 }
 315
 316 /*********************************************************************************************/
 317 MothurOut::~MothurOut() {
 318         try {
 319                 _uniqueInstance = 0;
 320
 321         }
 322         catch(exception& e) {
 323                 errorOut(e, "MothurOut", "MothurOut");
 324                 exit(1);
 325         }
 326 }
 327 /*********************************************************************************************/
 328 void MothurOut::mothurOut(string output) {
 329         try {
 330
 331                 #ifdef USE_MPI
 332                         int pid;
 333                         MPI_Comm_rank(MPI_COMM_WORLD, &pid);
 334
 335                         if (pid == 0) { //only one process should output to screen
 336                 #endif
 337
 338                 out << output;
 339         logger() << output;
 340
 341                 #ifdef USE_MPI
 342                         }
 343                 #endif
 344         }
 345         catch(exception& e) {
 346                 errorOut(e, "MothurOut", "MothurOut");
 347                 exit(1);
 348         }
 349 }
 350 /*********************************************************************************************/
 351 void MothurOut::mothurOutJustToScreen(string output) {
 352         try {
 353
 354 #ifdef USE_MPI
 355         int pid;
 356         MPI_Comm_rank(MPI_COMM_WORLD, &pid);
 357
 358         if (pid == 0) { //only one process should output to screen
 359 #endif
 360             logger() << output;
 361
 362 #ifdef USE_MPI
 363         }
 364 #endif
 365         }
 366         catch(exception& e) {
 367                 errorOut(e, "MothurOut", "MothurOut");
 368                 exit(1);
 369         }
 370 }
 371 /*********************************************************************************************/
 372 void MothurOut::mothurOutEndLine() {
 373         try {
 374                 #ifdef USE_MPI
 375                         int pid;
 376                         MPI_Comm_rank(MPI_COMM_WORLD, &pid);
 377
 378                         if (pid == 0) { //only one process should output to screen
 379                 #endif
 380
 381                 out << endl;
 382         logger() << endl;
 383
 384                 #ifdef USE_MPI
 385                         }
 386                 #endif
 387         }
 388         catch(exception& e) {
 389                 errorOut(e, "MothurOut", "MothurOutEndLine");
 390                 exit(1);
 391         }
 392 }
 393 /*********************************************************************************************/
 394 void MothurOut::mothurOut(string output, ofstream& outputFile) {
 395         try {
 396
 397 #ifdef USE_MPI
 398                 int pid;
 399                 MPI_Comm_rank(MPI_COMM_WORLD, &pid);
 400
 401                 if (pid == 0) { //only one process should output to screen
 402 #endif
 403
 404
 405                         out << output;
 406                         outputFile << output;
 407             logger() << output;
 408
 409 #ifdef USE_MPI
 410                 }
 411 #endif
 412
 413         }
 414         catch(exception& e) {
 415                 errorOut(e, "MothurOut", "MothurOut");
 416                 exit(1);
 417         }
 418 }
 419 /*********************************************************************************************/
 420 void MothurOut::mothurOutEndLine(ofstream& outputFile) {
 421         try {
 422 #ifdef USE_MPI
 423                 int pid;
 424                 MPI_Comm_rank(MPI_COMM_WORLD, &pid);
 425
 426                 if (pid == 0) { //only one process should output to screen
 427 #endif
 428
 429                         out << endl;
 430                         outputFile << endl;
 431             logger() << endl;
 432
 433 #ifdef USE_MPI
 434                 }
 435 #endif
 436         }
 437         catch(exception& e) {
 438                 errorOut(e, "MothurOut", "MothurOutEndLine");
 439                 exit(1);
 440         }
 441 }
 442 /*********************************************************************************************/
 443 void MothurOut::mothurOutJustToLog(string output) {
 444         try {
 445                 #ifdef USE_MPI
 446                         int pid;
 447                         MPI_Comm_rank(MPI_COMM_WORLD, &pid);
 448
 449                         if (pid == 0) { //only one process should output to screen
 450                 #endif
 451
 452                 out << output;
 453
 454                 #ifdef USE_MPI
 455                         }
 456                 #endif
 457         }
 458         catch(exception& e) {
 459                 errorOut(e, "MothurOut", "MothurOutJustToLog");
 460                 exit(1);
 461         }
 462 }
 463 /*********************************************************************************************/
 464 void MothurOut::errorOut(exception& e, string object, string function) {
 465         //double vm, rss;
 466         //mem_usage(vm, rss);
 467
 468     string errorType = toString(e.what());
 469
 470     int pos = errorType.find("bad_alloc");
 471     mothurOut("[ERROR]: ");
 472     mothurOut(errorType);
 473
 474     if (pos == string::npos) { //not bad_alloc
 475         mothurOut(" has occurred in the " + object + " class function " + function + ". Please contact Pat Schloss at mothur.bugs@gmail.com, and be sure to include the mothur.logFile with your inquiry.");
 476         mothurOutEndLine();
 477     }else { //bad alloc
 478         if (object == "cluster"){
 479             mothurOut(" has occurred in the " + object + " class function " + function + ". This error indicates your computer is running out of memory.  There are two common causes for this, file size and format.\n\nFile Size:\nThe cluster command loads your distance matrix into RAM, and your distance file is most likely too large to fit in RAM. There are two options to help with this. The first is to use a cutoff. By using a cutoff mothur will only load distances that are below the cutoff. If that is still not enough, there is a command called cluster.split, http://www.mothur.org/wiki/cluster.split which divides the distance matrix, and clusters the smaller pieces separately. You may also be able to reduce the size of the original distance matrix by using the commands outlined in the Schloss SOP, http://www.mothur.org/wiki/Schloss_SOP. \n\nWrong Format:\nThis error can be caused by trying to read a column formatted distance matrix using the phylip parameter. By default, the dist.seqs command generates a column formatted distance matrix. To make a phylip formatted matrix set the dist.seqs command parameter output to lt.  \n\nIf you are uable to resolve the issue, please contact Pat Schloss at mothur.bugs@gmail.com, and be sure to include the mothur.logFile with your inquiry.");
 480         }else if (object == "shhh.flows"){
 481                 mothurOut(" has occurred in the " + object + " class function " + function + ". This error indicates your computer is running out of memory. The shhh.flows command is very memory intensive. This error is most commonly caused by trying to process a dataset too large, using multiple processors, or failing to run trim.flows before shhh.flows. If you are running our 32bit version, your memory usage is limited to 4G.  If you have more than 4G of RAM and are running a 64bit OS, using our 64bit version may resolve your issue.  If you are using multiple processors, try running the command with processors=1, the more processors you use the more memory is required. Running trim.flows with an oligos file, and then shhh.flows with the file option may also resolve the issue. If for some reason you are unable to run shhh.flows with your data, a good alternative is to use the trim.seqs command using a 50-bp sliding window and to trim the sequence when the average quality score over that window drops below 35. Our results suggest that the sequencing error rates by this method are very good, but not quite as good as by shhh.flows and that the resulting sequences tend to be a bit shorter. If you are uable to resolve the issue, please contact Pat Schloss at mothur.bugs@gmail.com, and be sure to include the mothur.logFile with your inquiry. ");
 482         }else {
 483             mothurOut(" has occurred in the " + object + " class function " + function + ". This error indicates your computer is running out of memory.  This is most commonly caused by trying to process a dataset too large, using multiple processors, or a file format issue. If you are running our 32bit version, your memory usage is limited to 4G.  If you have more than 4G of RAM and are running a 64bit OS, using our 64bit version may resolve your issue.  If you are using multiple processors, try running the command with processors=1, the more processors you use the more memory is required. Also, you may be able to reduce the size of your dataset by using the commands outlined in the Schloss SOP, http://www.mothur.org/wiki/Schloss_SOP. If you are uable to resolve the issue, please contact Pat Schloss at mothur.bugs@gmail.com, and be sure to include the mothur.logFile with your inquiry.");
 484         }
 485     }
 486 }
 487 /*********************************************************************************************/
 488 //The following was originally from http://stackoverflow.com/questions/669438/how-to-get-memory-usage-at-run-time-in-c
 489 // process_mem_usage(double &, double &) - takes two doubles by reference,
 490 // attempts to read the system-dependent data for a process' virtual memory
 491 // size and resident set size, and return the results in KB.
 492 //
 493 // On failure, returns 0.0, 0.0
 494 int MothurOut::mem_usage(double& vm_usage, double& resident_set) {
 495   #if defined (__APPLE__) || (__MACH__) || (linux) || (__linux) || (__linux__) || (__unix__) || (__unix)
 496
 497            vm_usage     = 0.0;
 498            resident_set = 0.0;
 499
 500            // 'file' stat seems to give the most reliable results
 501            //
 502            ifstream stat_stream("/proc/self/stat",ios_base::in);
 503
 504            // dummy vars for leading entries in stat that we don't care about
 505            //
 506            string pid, comm, state, ppid, pgrp, session, tty_nr;
 507            string tpgid, flags, minflt, cminflt, majflt, cmajflt;
 508            string utime, stime, cutime, cstime, priority, nice;
 509            string O, itrealvalue, starttime;
 510
 511            // the two fields we want
 512            //
 513            unsigned long vsize;
 514            long rss;
 515
 516            stat_stream >> pid >> comm >> state >> ppid >> pgrp >> session >> tty_nr
 517                                    >> tpgid >> flags >> minflt >> cminflt >> majflt >> cmajflt
 518                                    >> utime >> stime >> cutime >> cstime >> priority >> nice
 519                                    >> O >> itrealvalue >> starttime >> vsize >> rss; // don't care about the rest
 520
 521            long page_size_kb = sysconf(_SC_PAGE_SIZE) / 1024; // in case x86-64 is configured to use 2MB pages
 522            vm_usage     = vsize / 1024.0;
 523            resident_set = rss * page_size_kb;
 524
 525            mothurOut("Memory Usage: vm = " + toString(vm_usage) + " rss = " + toString(resident_set) + "\n");
 526                 return 0;
 527
 528         #else
 529 /*              //windows memory usage
 530                 // Get the list of process identifiers.
 531                 DWORD aProcesses[1024], cbNeeded, cProcesses;
 532
 533                 if ( !EnumProcesses( aProcesses, sizeof(aProcesses), &cbNeeded ) ){ return 1; }
 534
 535                 // Calculate how many process identifiers were returned.
 536                 cProcesses = cbNeeded / sizeof(DWORD);
 537
 538                 // Print the memory usage for each process
 539                 for (int i = 0; i < cProcesses; i++ ) {
 540                         DWORD processID = aProcesses[i];
 541
 542                         PROCESS_MEMORY_COUNTERS pmc;
 543
 544                         HANDLE hProcess = OpenProcess((PROCESS_QUERY_INFORMATION | PROCESS_VM_READ), FALSE, processID);
 545
 546                         // Print the process identifier.
 547                         printf( "\nProcess ID: %u\n", processID);
 548
 549                         if (NULL != hProcess) {
 550
 551                                 if ( GetProcessMemoryInfo( hProcess, &pmc, sizeof(pmc)) ) {
 552                                         printf( "\tPageFaultCount: 0x%08X\n", pmc.PageFaultCount );
 553                                         printf( "\tPeakWorkingSetSize: 0x%08X\n", pmc.PeakWorkingSetSize );
 554                                         printf( "\tWorkingSetSize: 0x%08X\n", pmc.WorkingSetSize );
 555                                         printf( "\tQuotaPeakPagedPoolUsage: 0x%08X\n", pmc.QuotaPeakPagedPoolUsage );
 556                                         printf( "\tQuotaPagedPoolUsage: 0x%08X\n", pmc.QuotaPagedPoolUsage );
 557                                         printf( "\tQuotaPeakNonPagedPoolUsage: 0x%08X\n", pmc.QuotaPeakNonPagedPoolUsage );
 558                                         printf( "\tQuotaNonPagedPoolUsage: 0x%08X\n", pmc.QuotaNonPagedPoolUsage );
 559                                         printf( "\tPagefileUsage: 0x%08X\n", pmc.PagefileUsage );
 560                                         printf( "\tPeakPagefileUsage: 0x%08X\n", pmc.PeakPagefileUsage );
 561                                 }
 562                                 CloseHandle(hProcess);
 563                         }
 564                 }
 565 */
 566                         return 0;
 567
 568         #endif
 569 }
 570
 571
 572 /***********************************************************************/
 573 int MothurOut::openOutputFileAppend(string fileName, ofstream& fileHandle){
 574         try {
 575                 fileName = getFullPathName(fileName);
 576
 577                 fileHandle.open(fileName.c_str(), ios::app);
 578                 if(!fileHandle) {
 579                         mothurOut("[ERROR]: Could not open " + fileName); mothurOutEndLine();
 580                         return 1;
 581                 }
 582                 else {
 583                         return 0;
 584                 }
 585         }
 586         catch(exception& e) {
 587                 errorOut(e, "MothurOut", "openOutputFileAppend");
 588                 exit(1);
 589         }
 590 }
 591 /***********************************************************************/
 592 void MothurOut::gobble(istream& f){
 593         try {
 594
 595                 char d;
 596                 while(isspace(d=f.get()))               { ;}
 597                 if(!f.eof()) { f.putback(d); }
 598         }
 599         catch(exception& e) {
 600                 errorOut(e, "MothurOut", "gobble");
 601                 exit(1);
 602         }
 603 }
 604 /***********************************************************************/
 605 void MothurOut::gobble(istringstream& f){
 606         try {
 607                 char d;
 608                 while(isspace(d=f.get()))               {;}
 609                 if(!f.eof()) { f.putback(d); }
 610         }
 611         catch(exception& e) {
 612                 errorOut(e, "MothurOut", "gobble");
 613                 exit(1);
 614         }
 615 }
 616
 617 /***********************************************************************/
 618
 619 string MothurOut::getline(istringstream& fileHandle) {
 620         try {
 621
 622                 string line = "";
 623
 624                 while (!fileHandle.eof())       {
 625                         //get next character
 626                         char c = fileHandle.get();
 627
 628                         //are you at the end of the line
 629                         if ((c == '\n') || (c == '\r') || (c == '\f')){  break; }
 630                         else {          line += c;              }
 631                 }
 632
 633                 return line;
 634
 635         }
 636         catch(exception& e) {
 637                 errorOut(e, "MothurOut", "getline");
 638                 exit(1);
 639         }
 640 }
 641 /***********************************************************************/
 642
 643 string MothurOut::getline(ifstream& fileHandle) {
 644         try {
 645
 646                 string line = "";
 647
 648                 while (fileHandle)      {
 649                         //get next character
 650                         char c = fileHandle.get();
 651
 652                         //are you at the end of the line
 653                         if ((c == '\n') || (c == '\r') || (c == '\f') || (c == EOF)){  break;   }
 654                         else {          line += c;              }
 655                 }
 656
 657                 return line;
 658
 659         }
 660         catch(exception& e) {
 661                 errorOut(e, "MothurOut", "getline");
 662                 exit(1);
 663         }
 664 }
 665 /***********************************************************************/
 666
 667 #if defined (__APPLE__) || (__MACH__) || (linux) || (__linux) || (__linux__) || (__unix__) || (__unix)
 668 #ifdef USE_COMPRESSION
 669 inline bool endsWith(string s, const char * suffix){
 670   size_t suffixLength = strlen(suffix);
 671   return s.size() >= suffixLength && s.substr(s.size() - suffixLength, suffixLength).compare(suffix) == 0;
 672 }
 673 #endif
 674 #endif
 675
 676 string MothurOut::getRootName(string longName){
 677         try {
 678
 679                 string rootName = longName;
 680
 681 #if defined (__APPLE__) || (__MACH__) || (linux) || (__linux) || (__linux__) || (__unix__) || (__unix)
 682 #ifdef USE_COMPRESSION
 683     if (endsWith(rootName, ".gz") || endsWith(rootName, ".bz2")) {
 684       int pos = rootName.find_last_of('.');
 685       rootName = rootName.substr(0, pos);
 686       cerr << "shortening " << longName << " to " << rootName << "\n";
 687     }
 688 #endif
 689 #endif
 690                 if(rootName.find_last_of(".") != rootName.npos){
 691                         int pos = rootName.find_last_of('.')+1;
 692                         rootName = rootName.substr(0, pos);
 693                 }
 694
 695                 return rootName;
 696         }
 697         catch(exception& e) {
 698                 errorOut(e, "MothurOut", "getRootName");
 699                 exit(1);
 700         }
 701 }
 702 /***********************************************************************/
 703
 704 string MothurOut::getSimpleName(string longName){
 705         try {
 706                 string simpleName = longName;
 707
 708                 size_t found;
 709                 found=longName.find_last_of("/\\");
 710
 711                 if(found != longName.npos){
 712                         simpleName = longName.substr(found+1);
 713                 }
 714
 715                 return simpleName;
 716         }
 717         catch(exception& e) {
 718                 errorOut(e, "MothurOut", "getSimpleName");
 719                 exit(1);
 720         }
 721 }
 722
 723 /***********************************************************************/
 724
 725 int MothurOut::getRandomIndex(int highest){
 726         try {
 727
 728                 int random = (int) ((float)(highest+1) * (float)(rand()) / ((float)RAND_MAX+1.0));
 729
 730                 return random;
 731         }
 732         catch(exception& e) {
 733                 errorOut(e, "MothurOut", "getRandomIndex");
 734                 exit(1);
 735         }
 736
 737 }
 738 /**********************************************************************/
 739
 740 string MothurOut::getPathName(string longName){
 741         try {
 742                 string rootPathName = longName;
 743
 744                 if(longName.find_last_of("/\\") != longName.npos){
 745                         int pos = longName.find_last_of("/\\")+1;
 746                         rootPathName = longName.substr(0, pos);
 747                 }
 748
 749                 return rootPathName;
 750         }
 751         catch(exception& e) {
 752                 errorOut(e, "MothurOut", "getPathName");
 753                 exit(1);
 754         }
 755
 756 }
 757 /***********************************************************************/
 758
 759 bool MothurOut::dirCheck(string& dirName){
 760         try {
 761
 762         string tag = "";
 763         #ifdef USE_MPI
 764             int pid;
 765             MPI_Comm_rank(MPI_COMM_WORLD, &pid); //find out who we are
 766
 767             tag = toString(pid);
 768         #endif
 769
 770         //add / to name if needed
 771         string lastChar = dirName.substr(dirName.length()-1);
 772         #if defined (__APPLE__) || (__MACH__) || (linux) || (__linux) || (__linux__) || (__unix__) || (__unix)
 773         if (lastChar != "/") { dirName += "/"; }
 774         #else
 775         if (lastChar != "\\") { dirName += "\\"; }
 776         #endif
 777
 778         //test to make sure directory exists
 779         dirName = getFullPathName(dirName);
 780         string outTemp = dirName + tag + "temp";
 781         ofstream out;
 782         out.open(outTemp.c_str(), ios::trunc);
 783         if(!out) {
 784             mothurOut(dirName + " directory does not exist or is not writable."); mothurOutEndLine();
 785         }else{
 786             out.close();
 787             mothurRemove(outTemp);
 788             return true;
 789         }
 790
 791         return false;
 792     }
 793         catch(exception& e) {
 794                 errorOut(e, "MothurOut", "dirCheck");
 795                 exit(1);
 796         }
 797
 798 }
 799 /***********************************************************************/
 800
 801 string MothurOut::hasPath(string longName){
 802         try {
 803                 string path = "";
 804
 805                 size_t found;
 806                 found=longName.find_last_of("~/\\");
 807
 808                 if(found != longName.npos){
 809                         path = longName.substr(0, found+1);
 810                 }
 811
 812                 return path;
 813         }
 814         catch(exception& e) {
 815                 errorOut(e, "MothurOut", "hasPath");
 816                 exit(1);
 817         }
 818 }
 819
 820 /***********************************************************************/
 821
 822 string MothurOut::getExtension(string longName){
 823         try {
 824                 string extension = "";
 825
 826                 if(longName.find_last_of('.') != longName.npos){
 827                         int pos = longName.find_last_of('.');
 828                         extension = longName.substr(pos, longName.length());
 829                 }
 830
 831                 return extension;
 832         }
 833         catch(exception& e) {
 834                 errorOut(e, "MothurOut", "getExtension");
 835                 exit(1);
 836         }
 837 }
 838 /***********************************************************************/
 839 bool MothurOut::isBlank(string fileName){
 840         try {
 841
 842                 fileName = getFullPathName(fileName);
 843
 844                 ifstream fileHandle;
 845                 fileHandle.open(fileName.c_str());
 846                 if(!fileHandle) {
 847                         mothurOut("[ERROR]: Could not open " + fileName); mothurOutEndLine();
 848                         return false;
 849                 }else {
 850                         //check for blank file
 851                         gobble(fileHandle);
 852                         if (fileHandle.eof()) { fileHandle.close(); return true;  }
 853                         fileHandle.close();
 854                 }
 855                 return false;
 856         }
 857         catch(exception& e) {
 858                 errorOut(e, "MothurOut", "isBlank");
 859                 exit(1);
 860         }
 861 }
 862 /***********************************************************************/
 863
 864 string MothurOut::getFullPathName(string fileName){
 865         try{
 866
 867         string path = hasPath(fileName);
 868         string newFileName;
 869         int pos;
 870
 871         if (path == "") { return fileName; } //its a simple name
 872         else { //we need to complete the pathname
 873                 // ex. ../../../filename
 874                 // cwd = /user/work/desktop
 875
 876                 string cwd;
 877                 //get current working directory
 878                 #if defined (__APPLE__) || (__MACH__) || (linux) || (__linux) || (__linux__) || (__unix__) || (__unix)
 879
 880                         if (path.find("~") != -1) { //go to home directory
 881                                 string homeDir;
 882
 883                                 char *homepath = NULL;
 884                                 homepath = getenv ("HOME");
 885                                 if ( homepath != NULL) { homeDir = homepath; }
 886                                 else { homeDir = "";  }
 887
 888                                 newFileName = homeDir + fileName.substr(fileName.find("~")+1);
 889                                 return newFileName;
 890                         }else { //find path
 891                                 if (path.rfind("./") == string::npos) { return fileName; } //already complete name
 892                                 else { newFileName = fileName.substr(fileName.rfind("./")+2); } //save the complete part of the name
 893
 894                                 //char* cwdpath = new char[1024];
 895                                 //size_t size;
 896                                 //cwdpath=getcwd(cwdpath,size);
 897                                 //cwd = cwdpath;
 898
 899                                 char *cwdpath = NULL;
 900                                 cwdpath = getcwd(NULL, 0); // or _getcwd
 901                                 if ( cwdpath != NULL) { cwd = cwdpath; }
 902                                 else { cwd = "";  }
 903
 904
 905                                 //rip off first '/'
 906                                 string simpleCWD;
 907                                 if (cwd.length() > 0) { simpleCWD = cwd.substr(1); }
 908
 909                                 //break apart the current working directory
 910                                 vector<string> dirs;
 911                                 while (simpleCWD.find_first_of('/') != string::npos) {
 912                                         string dir = simpleCWD.substr(0,simpleCWD.find_first_of('/'));
 913                                         simpleCWD = simpleCWD.substr(simpleCWD.find_first_of('/')+1, simpleCWD.length());
 914                                         dirs.push_back(dir);
 915                                 }
 916                                 //get last one              // ex. ../../../filename = /user/work/desktop/filename
 917                                 dirs.push_back(simpleCWD);  //ex. dirs[0] = user, dirs[1] = work, dirs[2] = desktop
 918
 919
 920                                 int index = dirs.size()-1;
 921
 922                                 while((pos = path.rfind("./")) != string::npos) { //while you don't have a complete path
 923                                         if (pos == 0) { break;  //you are at the end
 924                                         }else if (path[(pos-1)] == '.') { //you want your parent directory ../
 925                                                 path = path.substr(0, pos-1);
 926                                                 index--;
 927                                                 if (index == 0) {  break; }
 928                                         }else if (path[(pos-1)] == '/') { //you want the current working dir ./
 929                                                 path = path.substr(0, pos);
 930                                         }else if (pos == 1) { break;  //you are at the end
 931                                         }else { mothurOut("cannot resolve path for " +  fileName + "\n"); return fileName; }
 932                                 }
 933
 934                                 for (int i = index; i >= 0; i--) {
 935                                         newFileName = dirs[i] +  "/" + newFileName;
 936                                 }
 937
 938                                 newFileName =  "/" +  newFileName;
 939                                 return newFileName;
 940                         }
 941                 #else
 942                         if (path.find("~") != string::npos) { //go to home directory
 943                                 string homeDir = getenv ("HOMEPATH");
 944                                 newFileName = homeDir + fileName.substr(fileName.find("~")+1);
 945                                 return newFileName;
 946                         }else { //find path
 947                                 if (path.rfind(".\\") == string::npos) { return fileName; } //already complete name
 948                                 else { newFileName = fileName.substr(fileName.rfind(".\\")+2); } //save the complete part of the name
 949
 950                                 char *cwdpath = NULL;
 951                                 cwdpath = getcwd(NULL, 0); // or _getcwd
 952                                 if ( cwdpath != NULL) { cwd = cwdpath; }
 953                                 else { cwd = "";  }
 954
 955                                 //break apart the current working directory
 956                                 vector<string> dirs;
 957                                 while (cwd.find_first_of('\\') != -1) {
 958                                         string dir = cwd.substr(0,cwd.find_first_of('\\'));
 959                                         cwd = cwd.substr(cwd.find_first_of('\\')+1, cwd.length());
 960                                         dirs.push_back(dir);
 961
 962                                 }
 963                                 //get last one
 964                                 dirs.push_back(cwd);  //ex. dirs[0] = user, dirs[1] = work, dirs[2] = desktop
 965
 966                                 int index = dirs.size()-1;
 967
 968                                 while((pos = path.rfind(".\\")) != string::npos) { //while you don't have a complete path
 969                                         if (pos == 0) { break;  //you are at the end
 970                                         }else if (path[(pos-1)] == '.') { //you want your parent directory ../
 971                                                 path = path.substr(0, pos-1);
 972                                                 index--;
 973                                                 if (index == 0) {  break; }
 974                                         }else if (path[(pos-1)] == '\\') { //you want the current working dir ./
 975                                                 path = path.substr(0, pos);
 976                                         }else if (pos == 1) { break;  //you are at the end
 977                                         }else { mothurOut("cannot resolve path for " +  fileName + "\n"); return fileName; }
 978                                 }
 979
 980                                 for (int i = index; i >= 0; i--) {
 981                                         newFileName = dirs[i] +  "\\" + newFileName;
 982                                 }
 983
 984                                 return newFileName;
 985                         }
 986
 987                 #endif
 988         }
 989         }
 990         catch(exception& e) {
 991                 errorOut(e, "MothurOut", "getFullPathName");
 992                 exit(1);
 993         }
 994 }
 995 /***********************************************************************/
 996
 997 int MothurOut::openInputFile(string fileName, ifstream& fileHandle, string m){
 998         try {
 999                         //get full path name
1000                         string completeFileName = getFullPathName(fileName);
1001 #if defined (__APPLE__) || (__MACH__) || (linux) || (__linux) || (__linux__) || (__unix__) || (__unix)
1002 #ifdef USE_COMPRESSION
1003       // check for gzipped or bzipped file
1004       if (endsWith(completeFileName, ".gz") || endsWith(completeFileName, ".bz2")) {
1005         string tempName = string(tmpnam(0));
1006         mkfifo(tempName.c_str(), 0666);
1007         int fork_result = fork();
1008         if (fork_result < 0) {
1009           cerr << "Error forking.\n";
1010           exit(1);
1011         } else if (fork_result == 0) {
1012           string command = (endsWith(completeFileName, ".gz") ? "zcat " : "bzcat ") + completeFileName + string(" > ") + tempName;
1013           cerr << "Decompressing " << completeFileName << " via temporary named pipe " << tempName << "\n";
1014           system(command.c_str());
1015           cerr << "Done decompressing " << completeFileName << "\n";
1016           mothurRemove(tempName);
1017           exit(EXIT_SUCCESS);
1018         } else {
1019           cerr << "waiting on child process " << fork_result << "\n";
1020           completeFileName = tempName;
1021         }
1022       }
1023 #endif
1024 #endif
1025                         fileHandle.open(completeFileName.c_str());
1026                         if(!fileHandle) {
1027                                 //mothurOut("[ERROR]: Could not open " + completeFileName); mothurOutEndLine();
1028                                 return 1;
1029                         }else {
1030                                 //check for blank file
1031                                 gobble(fileHandle);
1032                                 return 0;
1033                         }
1034         }
1035         catch(exception& e) {
1036                 errorOut(e, "MothurOut", "openInputFile - no Error");
1037                 exit(1);
1038         }
1039 }
1040 /***********************************************************************/
1041
1042 int MothurOut::openInputFile(string fileName, ifstream& fileHandle){
1043         try {
1044
1045                 //get full path name
1046                 string completeFileName = getFullPathName(fileName);
1047 #if defined (__APPLE__) || (__MACH__) || (linux) || (__linux) || (__linux__) || (__unix__) || (__unix)
1048 #ifdef USE_COMPRESSION
1049   // check for gzipped or bzipped file
1050   if (endsWith(completeFileName, ".gz") || endsWith(completeFileName, ".bz2")) {
1051     string tempName = string(tmpnam(0));
1052     mkfifo(tempName.c_str(), 0666);
1053     int fork_result = fork();
1054     if (fork_result < 0) {
1055       cerr << "Error forking.\n";
1056       exit(1);
1057     } else if (fork_result == 0) {
1058       string command = (endsWith(completeFileName, ".gz") ? "zcat " : "bzcat ") + completeFileName + string(" > ") + tempName;
1059       cerr << "Decompressing " << completeFileName << " via temporary named pipe " << tempName << "\n";
1060       system(command.c_str());
1061       cerr << "Done decompressing " << completeFileName << "\n";
1062       mothurRemove(tempName);
1063       exit(EXIT_SUCCESS);
1064     } else {
1065       cerr << "waiting on child process " << fork_result << "\n";
1066       completeFileName = tempName;
1067     }
1068   }
1069 #endif
1070 #endif
1071
1072                 fileHandle.open(completeFileName.c_str());
1073                 if(!fileHandle) {
1074                         mothurOut("[ERROR]: Could not open " + completeFileName); mothurOutEndLine();
1075                         return 1;
1076                 }
1077                 else {
1078                         //check for blank file
1079                         gobble(fileHandle);
1080                         if (fileHandle.eof()) { mothurOut("[ERROR]: " + completeFileName + " is blank. Please correct."); mothurOutEndLine();  }
1081
1082                         return 0;
1083                 }
1084         }
1085         catch(exception& e) {
1086                 errorOut(e, "MothurOut", "openInputFile");
1087                 exit(1);
1088         }
1089 }
1090 /***********************************************************************/
1091
1092 int MothurOut::renameFile(string oldName, string newName){
1093         try {
1094
1095         if (oldName == newName) { return 0; }
1096
1097                 ifstream inTest;
1098                 int exist = openInputFile(newName, inTest, "");
1099                 inTest.close();
1100
1101         #if defined (__APPLE__) || (__MACH__) || (linux) || (__linux) || (__linux__) || (__unix__) || (__unix)
1102                 if (exist == 0) { //you could open it so you want to delete it
1103                         string command = "rm " + newName;
1104                         system(command.c_str());
1105                 }
1106
1107                 string command = "mv " + oldName + " " + newName;
1108                 system(command.c_str());
1109         #else
1110                 mothurRemove(newName);
1111                 int renameOk = rename(oldName.c_str(), newName.c_str());
1112         #endif
1113                 return 0;
1114
1115         }
1116         catch(exception& e) {
1117                 errorOut(e, "MothurOut", "renameFile");
1118                 exit(1);
1119         }
1120 }
1121
1122 /***********************************************************************/
1123
1124 int MothurOut::openOutputFile(string fileName, ofstream& fileHandle){
1125         try {
1126
1127                 string completeFileName = getFullPathName(fileName);
1128 #if defined (__APPLE__) || (__MACH__) || (linux) || (__linux) || (__linux__) || (__unix__) || (__unix)
1129 #ifdef USE_COMPRESSION
1130     // check for gzipped file
1131     if (endsWith(completeFileName, ".gz") || endsWith(completeFileName, ".bz2")) {
1132       string tempName = string(tmpnam(0));
1133       mkfifo(tempName.c_str(), 0666);
1134       cerr << "Compressing " << completeFileName << " via temporary named pipe " << tempName << "\n";
1135       int fork_result = fork();
1136       if (fork_result < 0) {
1137         cerr << "Error forking.\n";
1138         exit(1);
1139       } else if (fork_result == 0) {
1140         string command = string(endsWith(completeFileName, ".gz") ?  "gzip" : "bzip2") + " -v > " + completeFileName + string(" < ") + tempName;
1141         system(command.c_str());
1142         exit(0);
1143       } else {
1144         completeFileName = tempName;
1145       }
1146     }
1147 #endif
1148 #endif
1149                 fileHandle.open(completeFileName.c_str(), ios::trunc);
1150                 if(!fileHandle) {
1151                         mothurOut("[ERROR]: Could not open " + completeFileName); mothurOutEndLine();
1152                         return 1;
1153                 }
1154                 else {
1155                         return 0;
1156                 }
1157         }
1158         catch(exception& e) {
1159                 errorOut(e, "MothurOut", "openOutputFile");
1160                 exit(1);
1161         }
1162
1163 }
1164
1165 /**************************************************************************************************/
1166 int MothurOut::appendFiles(string temp, string filename) {
1167         try{
1168                 ofstream output;
1169                 ifstream input;
1170
1171                 //open output file in append mode
1172                 openOutputFileAppend(filename, output);
1173                 int ableToOpen = openInputFile(temp, input, "no error");
1174                 //int ableToOpen = openInputFile(temp, input);
1175
1176                 int numLines = 0;
1177                 if (ableToOpen == 0) { //you opened it
1178
1179             char buffer[4096];
1180             while (!input.eof()) {
1181                 input.read(buffer, 4096);
1182                 output.write(buffer, input.gcount());
1183                 //count number of lines
1184                 for (int i = 0; i < input.gcount(); i++) {  if (buffer[i] == '\n') {numLines++;} }
1185             }
1186                         input.close();
1187                 }
1188
1189                 output.close();
1190
1191                 return numLines;
1192         }
1193         catch(exception& e) {
1194                 errorOut(e, "MothurOut", "appendFiles");
1195                 exit(1);
1196         }
1197 }
1198 /**************************************************************************************************/
1199 int MothurOut::appendFilesWithoutHeaders(string temp, string filename) {
1200         try{
1201                 ofstream output;
1202                 ifstream input;
1203
1204                 //open output file in append mode
1205                 openOutputFileAppend(filename, output);
1206                 int ableToOpen = openInputFile(temp, input, "no error");
1207                 //int ableToOpen = openInputFile(temp, input);
1208
1209                 int numLines = 0;
1210                 if (ableToOpen == 0) { //you opened it
1211
1212             string headers = getline(input); gobble(input);
1213             if (debug) { mothurOut("[DEBUG]: skipping headers " + headers +'\n'); }
1214
1215             char buffer[4096];
1216             while (!input.eof()) {
1217                 input.read(buffer, 4096);
1218                 output.write(buffer, input.gcount());
1219                 //count number of lines
1220                 for (int i = 0; i < input.gcount(); i++) {  if (buffer[i] == '\n') {numLines++;} }
1221             }
1222                         input.close();
1223                 }
1224
1225                 output.close();
1226
1227                 return numLines;
1228         }
1229         catch(exception& e) {
1230                 errorOut(e, "MothurOut", "appendFiles");
1231                 exit(1);
1232         }
1233 }
1234 /**************************************************************************************************/
1235 string MothurOut::sortFile(string distFile, string outputDir){
1236         try {
1237
1238                 //if (outputDir == "") {  outputDir += hasPath(distFile);  }
1239                 string outfile = getRootName(distFile) + "sorted.dist";
1240
1241
1242                 //if you can, use the unix sort since its been optimized for years
1243                 #if defined (__APPLE__) || (__MACH__) || (linux) || (__linux) || (__linux__) || (__unix__) || (__unix)
1244                         string command = "sort -n -k +3 " + distFile + " -o " + outfile;
1245                         system(command.c_str());
1246                 #else //you are stuck with my best attempt...
1247                         //windows sort does not have a way to specify a column, only a character in the line
1248                         //since we cannot assume that the distance will always be at the the same character location on each line
1249                         //due to variable sequence name lengths, I chose to force the distance into first position, then sort and then put it back.
1250
1251                         //read in file line by file and put distance first
1252                         string tempDistFile = distFile + ".temp";
1253                         ifstream input;
1254                         ofstream output;
1255                         openInputFile(distFile, input);
1256                         openOutputFile(tempDistFile, output);
1257
1258                         string firstName, secondName;
1259                         float dist;
1260                         while (!input.eof()) {
1261                                 input >> firstName >> secondName >> dist;
1262                                 output << dist << '\t' << firstName << '\t' << secondName << endl;
1263                                 gobble(input);
1264                         }
1265                         input.close();
1266                         output.close();
1267
1268
1269                         //sort using windows sort
1270                         string tempOutfile = outfile + ".temp";
1271                         string command = "sort " + tempDistFile + " /O " + tempOutfile;
1272                         system(command.c_str());
1273
1274                         //read in sorted file and put distance at end again
1275                         ifstream input2;
1276             ofstream output2;
1277                         openInputFile(tempOutfile, input2);
1278                         openOutputFile(outfile, output2);
1279
1280             while (!input2.eof()) {
1281                                 input2 >> dist >> firstName >> secondName;
1282                                 output2 << firstName << '\t' << secondName << '\t' << dist << endl;
1283                                 gobble(input2);
1284                         }
1285                         input2.close();
1286                         output2.close();
1287
1288                         //remove temp files
1289                         mothurRemove(tempDistFile);
1290                         mothurRemove(tempOutfile);
1291                 #endif
1292
1293                 return outfile;
1294         }
1295         catch(exception& e) {
1296                 errorOut(e, "MothurOut", "sortFile");
1297                 exit(1);
1298         }
1299 }
1300 /**************************************************************************************************/
1301 vector<unsigned long long> MothurOut::setFilePosFasta(string filename, int& num) {
1302         try {
1303                         vector<unsigned long long> positions;
1304                         ifstream inFASTA;
1305                         //openInputFile(filename, inFASTA);
1306                         inFASTA.open(filename.c_str(), ios::binary);
1307
1308                         string input;
1309                         unsigned long long count = 0;
1310                         while(!inFASTA.eof()){
1311                                 //input = getline(inFASTA);
1312                                 //cout << input << '\t' << inFASTA.tellg() << endl;
1313                                 //if (input.length() != 0) {
1314                                 //      if(input[0] == '>'){    unsigned long int pos = inFASTA.tellg(); positions.push_back(pos - input.length() - 1);  cout << (pos - input.length() - 1) << endl; }
1315                                 //}
1316                                 //gobble(inFASTA); //has to be here since windows line endings are 2 characters and mess up the positions
1317                                 char c = inFASTA.get(); count++;
1318                                 if (c == '>') {
1319                                         positions.push_back(count-1);
1320                                         if (debug) { mothurOut("[DEBUG]: numSeqs = " + toString(positions.size()) +  " count = " + toString(count) + ".\n"); }
1321                                 }
1322                         }
1323                         inFASTA.close();
1324
1325                         num = positions.size();
1326             if (debug) { mothurOut("[DEBUG]: num = " + toString(num) + ".\n"); }
1327                         FILE * pFile;
1328                         unsigned long long size;
1329
1330                         //get num bytes in file
1331                         pFile = fopen (filename.c_str(),"rb");
1332                         if (pFile==NULL) perror ("Error opening file");
1333                         else{
1334                                 fseek (pFile, 0, SEEK_END);
1335                                 size=ftell (pFile);
1336                                 fclose (pFile);
1337                         }
1338
1339                         /*unsigned long long size = positions[(positions.size()-1)];
1340                         ifstream in;
1341                         openInputFile(filename, in);
1342
1343                         in.seekg(size);
1344
1345                         while(in.get()){
1346                                 if(in.eof())            {       break;  }
1347                                 else                            {       size++; }
1348                         }
1349                         in.close();*/
1350
1351             if (debug) { mothurOut("[DEBUG]: size = " + toString(size) + ".\n"); }
1352
1353                         positions.push_back(size);
1354                         positions[0] = 0;
1355
1356                         return positions;
1357         }
1358         catch(exception& e) {
1359                 errorOut(e, "MothurOut", "setFilePosFasta");
1360                 exit(1);
1361         }
1362 }
1363 //**********************************************************************************************************************
1364 vector<consTax> MothurOut::readConsTax(string inputfile){
1365         try {
1366
1367         vector<consTax> taxes;
1368
1369         ifstream in;
1370         openInputFile(inputfile, in);
1371
1372         //read headers
1373         getline(in);
1374
1375         while (!in.eof()) {
1376
1377             if (control_pressed) { break; }
1378
1379             string otu = ""; string tax = "unknown";
1380             int size = 0;
1381
1382             in >> otu >> size >> tax; gobble(in);
1383             consTax temp(otu, tax, size);
1384             taxes.push_back(temp);
1385         }
1386         in.close();
1387
1388         return taxes;
1389     }
1390         catch(exception& e) {
1391                 errorOut(e, "MothurOut", "readConsTax");
1392                 exit(1);
1393         }
1394 }
1395 //**********************************************************************************************************************
1396 int MothurOut::readConsTax(string inputfile, map<string, consTax2>& taxes){
1397         try {
1398         ifstream in;
1399         openInputFile(inputfile, in);
1400
1401         //read headers
1402         getline(in);
1403
1404         while (!in.eof()) {
1405
1406             if (control_pressed) { break; }
1407
1408             string otu = ""; string tax = "unknown";
1409             int size = 0;
1410
1411             in >> otu >> size >> tax; gobble(in);
1412             consTax2 temp(tax, size);
1413             taxes[otu] = temp;
1414         }
1415         in.close();
1416
1417         return 0;
1418     }
1419         catch(exception& e) {
1420                 errorOut(e, "MothurOut", "readConsTax");
1421                 exit(1);
1422         }
1423 }
1424 /**************************************************************************************************/
1425 vector<unsigned long long> MothurOut::setFilePosEachLine(string filename, int& num) {
1426         try {
1427                         filename = getFullPathName(filename);
1428
1429                         vector<unsigned long long> positions;
1430                         ifstream in;
1431                         //openInputFile(filename, in);
1432                         in.open(filename.c_str(), ios::binary);
1433
1434                         string input;
1435                         unsigned long long count = 0;
1436                         positions.push_back(0);
1437
1438                         while(!in.eof()){
1439                                 //getline counting reads
1440                                 char d = in.get(); count++;
1441                                 while ((d != '\n') && (d != '\r') && (d != '\f') && (d != in.eof()))    {
1442                                         //get next character
1443                                         d = in.get();
1444                                         count++;
1445                                 }
1446
1447                                 if (!in.eof()) {
1448                                         d=in.get(); count++;
1449                                         while(isspace(d) && (d != in.eof()))            { d=in.get(); count++;}
1450                                 }
1451                                 positions.push_back(count-1);
1452                                 //cout << count-1 << endl;
1453                         }
1454                         in.close();
1455
1456                         num = positions.size()-1;
1457
1458                         FILE * pFile;
1459                         unsigned long long size;
1460
1461                         //get num bytes in file
1462                         pFile = fopen (filename.c_str(),"rb");
1463                         if (pFile==NULL) perror ("Error opening file");
1464                         else{
1465                                 fseek (pFile, 0, SEEK_END);
1466                                 size=ftell (pFile);
1467                                 fclose (pFile);
1468                         }
1469
1470                         positions[(positions.size()-1)] = size;
1471
1472                         return positions;
1473         }
1474         catch(exception& e) {
1475                 errorOut(e, "MothurOut", "setFilePosEachLine");
1476                 exit(1);
1477         }
1478 }
1479 /**************************************************************************************************/
1480
1481 vector<unsigned long long> MothurOut::divideFile(string filename, int& proc) {
1482         try{
1483                 vector<unsigned long long> filePos;
1484                 filePos.push_back(0);
1485
1486                 FILE * pFile;
1487                 unsigned long long size;
1488
1489                 filename = getFullPathName(filename);
1490
1491                 //get num bytes in file
1492                 pFile = fopen (filename.c_str(),"rb");
1493                 if (pFile==NULL) perror ("Error opening file");
1494                 else{
1495                         fseek (pFile, 0, SEEK_END);
1496                         size=ftell (pFile);
1497                         fclose (pFile);
1498                 }
1499
1500         #if defined (__APPLE__) || (__MACH__) || (linux) || (__linux) || (__linux__) || (__unix__) || (__unix)
1501
1502                 //estimate file breaks
1503                 unsigned long long chunkSize = 0;
1504                 chunkSize = size / proc;
1505
1506                 //file to small to divide by processors
1507                 if (chunkSize == 0)  {  proc = 1;       filePos.push_back(size); return filePos;        }
1508
1509                 //for each process seekg to closest file break and search for next '>' char. make that the filebreak
1510                 for (int i = 0; i < proc; i++) {
1511                         unsigned long long spot = (i+1) * chunkSize;
1512
1513                         ifstream in;
1514                         openInputFile(filename, in);
1515                         in.seekg(spot);
1516
1517                         //look for next '>'
1518                         unsigned long long newSpot = spot;
1519                         while (!in.eof()) {
1520                            char c = in.get();
1521
1522                            if (c == '>') {   in.putback(c); newSpot = in.tellg(); break;  }
1523                            else if (int(c) == -1) { break; }
1524
1525                         }
1526
1527                         //there was not another sequence before the end of the file
1528                         unsigned long long sanityPos = in.tellg();
1529
1530                         if (sanityPos == -1) {  break;  }
1531                         else {  filePos.push_back(newSpot);  }
1532
1533                         in.close();
1534                 }
1535
1536                 //save end pos
1537                 filePos.push_back(size);
1538
1539                 //sanity check filePos
1540                 for (int i = 0; i < (filePos.size()-1); i++) {
1541                         if (filePos[(i+1)] <= filePos[i]) {  filePos.erase(filePos.begin()+(i+1)); i--; }
1542                 }
1543
1544                 proc = (filePos.size() - 1);
1545 #else
1546                 mothurOut("[ERROR]: Windows version should not be calling the divideFile function."); mothurOutEndLine();
1547                 proc=1;
1548                 filePos.push_back(size);
1549 #endif
1550                 return filePos;
1551         }
1552         catch(exception& e) {
1553                 errorOut(e, "MothurOut", "divideFile");
1554                 exit(1);
1555         }
1556 }
1557 /**************************************************************************************************/
1558
1559 vector<unsigned long long> MothurOut::divideFilePerLine(string filename, int& proc) {
1560         try{
1561                 vector<unsigned long long> filePos;
1562                 filePos.push_back(0);
1563
1564                 FILE * pFile;
1565                 unsigned long long size;
1566
1567                 filename = getFullPathName(filename);
1568
1569                 //get num bytes in file
1570                 pFile = fopen (filename.c_str(),"rb");
1571                 if (pFile==NULL) perror ("Error opening file");
1572                 else{
1573                         fseek (pFile, 0, SEEK_END);
1574                         size=ftell (pFile);
1575                         fclose (pFile);
1576                 }
1577
1578 #if defined (__APPLE__) || (__MACH__) || (linux) || (__linux) || (__linux__) || (__unix__) || (__unix)
1579
1580                 //estimate file breaks
1581                 unsigned long long chunkSize = 0;
1582                 chunkSize = size / proc;
1583
1584                 //file to small to divide by processors
1585                 if (chunkSize == 0)  {  proc = 1;       filePos.push_back(size); return filePos;        }
1586
1587                 //for each process seekg to closest file break and search for next '>' char. make that the filebreak
1588                 for (int i = 0; i < proc; i++) {
1589                         unsigned long long spot = (i+1) * chunkSize;
1590
1591                         ifstream in;
1592                         openInputFile(filename, in);
1593                         in.seekg(spot);
1594
1595                         //look for next line break
1596                         unsigned long long newSpot = spot;
1597                         while (!in.eof()) {
1598                 char c = in.get();
1599
1600                                 if ((c == '\n') || (c == '\r') || (c == '\f'))  { gobble(in); newSpot = in.tellg(); break; }
1601                 else if (int(c) == -1) { break; }
1602             }
1603
1604                         //there was not another line before the end of the file
1605                         unsigned long long sanityPos = in.tellg();
1606
1607                         if (sanityPos == -1) {  break;  }
1608                         else {  filePos.push_back(newSpot);  }
1609
1610                         in.close();
1611                 }
1612
1613                 //save end pos
1614                 filePos.push_back(size);
1615
1616                 //sanity check filePos
1617                 for (int i = 0; i < (filePos.size()-1); i++) {
1618                         if (filePos[(i+1)] <= filePos[i]) {  filePos.erase(filePos.begin()+(i+1)); i--; }
1619                 }
1620
1621                 proc = (filePos.size() - 1);
1622 #else
1623                 mothurOut("[ERROR]: Windows version should not be calling the divideFile function."); mothurOutEndLine();
1624                 proc=1;
1625                 filePos.push_back(size);
1626 #endif
1627                 return filePos;
1628         }
1629         catch(exception& e) {
1630                 errorOut(e, "MothurOut", "divideFile");
1631                 exit(1);
1632         }
1633 }
1634 /**************************************************************************************************/
1635 int MothurOut::divideFile(string filename, int& proc, vector<string>& files) {
1636         try{
1637
1638                 vector<unsigned long long> filePos = divideFile(filename, proc);
1639
1640                 for (int i = 0; i < (filePos.size()-1); i++) {
1641
1642                         //read file chunk
1643                         ifstream in;
1644                         openInputFile(filename, in);
1645                         in.seekg(filePos[i]);
1646                         unsigned long long size = filePos[(i+1)] - filePos[i];
1647                         char* chunk = new char[size];
1648                         in.read(chunk, size);
1649                         in.close();
1650
1651                         //open new file
1652                         string fileChunkName = filename + "." + toString(i) + ".tmp";
1653                         ofstream out;
1654                         openOutputFile(fileChunkName, out);
1655
1656                         out << chunk << endl;
1657                         out.close();
1658                         delete[] chunk;
1659
1660                         //save name
1661                         files.push_back(fileChunkName);
1662                 }
1663
1664                 return 0;
1665         }
1666         catch(exception& e) {
1667                 errorOut(e, "MothurOut", "divideFile");
1668                 exit(1);
1669         }
1670 }
1671 /***********************************************************************/
1672
1673 bool MothurOut::isTrue(string f){
1674         try {
1675
1676                 for (int i = 0; i < f.length(); i++) { f[i] = toupper(f[i]); }
1677
1678                 if ((f == "TRUE") || (f == "T")) {      return true;    }
1679                 else {  return false;  }
1680         }
1681         catch(exception& e) {
1682                 errorOut(e, "MothurOut", "isTrue");
1683                 exit(1);
1684         }
1685 }
1686
1687 /***********************************************************************/
1688
1689 float MothurOut::roundDist(float dist, int precision){
1690         try {
1691                 return int(dist * precision + 0.5)/float(precision);
1692         }
1693         catch(exception& e) {
1694                 errorOut(e, "MothurOut", "roundDist");
1695                 exit(1);
1696         }
1697 }
1698 /***********************************************************************/
1699
1700 float MothurOut::ceilDist(float dist, int precision){
1701         try {
1702                 return int(ceil(dist * precision))/float(precision);
1703         }
1704         catch(exception& e) {
1705                 errorOut(e, "MothurOut", "ceilDist");
1706                 exit(1);
1707         }
1708 }
1709 /***********************************************************************/
1710
1711 vector<string> MothurOut::splitWhiteSpace(string& rest, char buffer[], int size){
1712         try {
1713         vector<string> pieces;
1714
1715         for (int i = 0; i < size; i++) {
1716             if (!isspace(buffer[i]))  { rest += buffer[i];  }
1717             else {
1718                 if (rest != "") { pieces.push_back(rest);  rest = ""; }
1719                 while (i < size) {  //gobble white space
1720                     if (isspace(buffer[i])) { i++; }
1721                     else { rest = buffer[i];  break; } //cout << "next piece buffer = " << nextPiece << endl;
1722                 }
1723             }
1724         }
1725
1726         return pieces;
1727         }
1728         catch(exception& e) {
1729                 errorOut(e, "MothurOut", "splitWhiteSpace");
1730                 exit(1);
1731         }
1732 }
1733 /***********************************************************************/
1734 vector<string> MothurOut::splitWhiteSpace(string input){
1735         try {
1736         vector<string> pieces;
1737         string rest = "";
1738
1739         for (int i = 0; i < input.length(); i++) {
1740             if (!isspace(input[i]))  { rest += input[i];  }
1741             else {
1742                 if (rest != "") { pieces.push_back(rest);  rest = ""; }
1743                 while (i < input.length()) {  //gobble white space
1744                     if (isspace(input[i])) { i++; }
1745                     else { rest = input[i];  break; } //cout << "next piece buffer = " << nextPiece << endl;
1746                 }
1747             }
1748         }
1749
1750         if (rest != "") { pieces.push_back(rest); }
1751
1752         return pieces;
1753         }
1754         catch(exception& e) {
1755                 errorOut(e, "MothurOut", "splitWhiteSpace");
1756                 exit(1);
1757         }
1758 }
1759 /***********************************************************************/
1760 vector<string> MothurOut::splitWhiteSpaceWithQuotes(string input){
1761         try {
1762         vector<string> pieces;
1763         string rest = "";
1764
1765         int pos = input.find('\'');
1766         int pos2 = input.find('\"');
1767
1768         if ((pos == string::npos) && (pos2 == string::npos)) { return splitWhiteSpace(input); } //no quotes to worry about
1769         else {
1770             for (int i = 0; i < input.length(); i++) {
1771                 if ((input[i] == '\'') || (input[i] == '\"') || (rest == "\'") || (rest == "\"")) { //grab everything til end or next ' or "
1772                     rest += input[i];
1773                     for (int j = i+1; j < input.length(); j++) {
1774                         if ((input[j] == '\'') || (input[j] == '\"')) {  //then quit
1775                             rest += input[j];
1776                             i = j+1;
1777                             j+=input.length();
1778                         }else { rest += input[j]; }
1779                     }
1780                 }else if (!isspace(input[i]))  { rest += input[i];  }
1781                 else {
1782                     if (rest != "") { pieces.push_back(rest);  rest = ""; }
1783                     while (i < input.length()) {  //gobble white space
1784                         if (isspace(input[i])) { i++; }
1785                         else { rest = input[i];  break; } //cout << "next piece buffer = " << nextPiece << endl;
1786                     }
1787                 }
1788             }
1789
1790             if (rest != "") { pieces.push_back(rest); }
1791         }
1792         return pieces;
1793         }
1794         catch(exception& e) {
1795                 errorOut(e, "MothurOut", "splitWhiteSpace");
1796                 exit(1);
1797         }
1798 }
1799 //**********************************************************************************************************************
1800 int MothurOut::readTax(string namefile, map<string, string>& taxMap) {
1801         try {
1802         //open input file
1803                 ifstream in;
1804                 openInputFile(namefile, in);
1805
1806         string rest = "";
1807         char buffer[4096];
1808         bool pairDone = false;
1809         bool columnOne = true;
1810         string firstCol, secondCol;
1811
1812                 while (!in.eof()) {
1813                         if (control_pressed) { break; }
1814
1815             in.read(buffer, 4096);
1816             vector<string> pieces = splitWhiteSpace(rest, buffer, in.gcount());
1817
1818             for (int i = 0; i < pieces.size(); i++) {
1819                 if (columnOne) {  firstCol = pieces[i]; columnOne=false; }
1820                 else  { secondCol = pieces[i]; pairDone = true; columnOne=true; }
1821
1822                 if (pairDone) {
1823                     checkName(firstCol);
1824                     //are there confidence scores, if so remove them
1825                     if (secondCol.find_first_of('(') != -1) {  removeConfidences(secondCol);    }
1826                     map<string, string>::iterator itTax = taxMap.find(firstCol);
1827
1828                     if(itTax == taxMap.end()) {
1829                         bool ignore = false;
1830                         if (secondCol != "") { if (secondCol[secondCol.length()-1] != ';') { mothurOut("[ERROR]: " + firstCol + " is missing the final ';', ignoring.\n"); ignore=true; }
1831                         }
1832                         if (!ignore) { taxMap[firstCol] = secondCol; }
1833                         if (debug) {  mothurOut("[DEBUG]: name = '" + firstCol + "' tax = '" + secondCol + "'\n");  }
1834                     }else {
1835                         mothurOut("[ERROR]: " + firstCol + " is already in your taxonomy file, names must be unique./n"); control_pressed = true;
1836                     }
1837                     pairDone = false;
1838                 }
1839             }
1840                 }
1841                 in.close();
1842
1843         if (rest != "") {
1844             vector<string> pieces = splitWhiteSpace(rest);
1845
1846             for (int i = 0; i < pieces.size(); i++) {
1847                 if (columnOne) {  firstCol = pieces[i]; columnOne=false; }
1848                 else  { secondCol = pieces[i]; pairDone = true; columnOne=true; }
1849
1850                 if (pairDone) {
1851                     checkName(firstCol);
1852                     //are there confidence scores, if so remove them
1853                     if (secondCol.find_first_of('(') != -1) {  removeConfidences(secondCol);    }
1854                     map<string, string>::iterator itTax = taxMap.find(firstCol);
1855
1856                     if(itTax == taxMap.end()) {
1857                         bool ignore = false;
1858                         if (secondCol != "") { if (secondCol[secondCol.length()-1] != ';') { mothurOut("[ERROR]: " + firstCol + " is missing the final ';', ignoring.\n"); ignore=true; }
1859                         }
1860                         if (!ignore) { taxMap[firstCol] = secondCol; }
1861                         if (debug) {  mothurOut("[DEBUG]: name = '" + firstCol + "' tax = '" + secondCol + "'\n");  }
1862                     }else {
1863                         mothurOut("[ERROR]: " + firstCol + " is already in your taxonomy file, names must be unique./n"); control_pressed = true;
1864                     }
1865
1866                     pairDone = false;
1867                 }
1868             }
1869         }
1870
1871                 return taxMap.size();
1872
1873         }
1874         catch(exception& e) {
1875                 errorOut(e, "MothurOut", "readTax");
1876                 exit(1);
1877         }
1878 }
1879 /**********************************************************************************************************************/
1880 int MothurOut::readNames(string namefile, map<string, string>& nameMap, bool redund) {
1881         try {
1882                 //open input file
1883                 ifstream in;
1884                 openInputFile(namefile, in);
1885
1886         string rest = "";
1887         char buffer[4096];
1888         bool pairDone = false;
1889         bool columnOne = true;
1890         string firstCol, secondCol;
1891
1892                 while (!in.eof()) {
1893                         if (control_pressed) { break; }
1894
1895             in.read(buffer, 4096);
1896             vector<string> pieces = splitWhiteSpace(rest, buffer, in.gcount());
1897
1898             for (int i = 0; i < pieces.size(); i++) {
1899                 if (columnOne) {  firstCol = pieces[i]; columnOne=false; }
1900                 else  { secondCol = pieces[i]; pairDone = true; columnOne=true; }
1901
1902                 if (pairDone) {
1903                     checkName(firstCol);
1904                     checkName(secondCol);
1905
1906                     //parse names into vector
1907                     vector<string> theseNames;
1908                     splitAtComma(secondCol, theseNames);
1909                     for (int i = 0; i < theseNames.size(); i++) {  nameMap[theseNames[i]] = firstCol;  }
1910                     pairDone = false;
1911                 }
1912             }
1913                 }
1914                 in.close();
1915
1916         if (rest != "") {
1917             vector<string> pieces = splitWhiteSpace(rest);
1918
1919             for (int i = 0; i < pieces.size(); i++) {
1920                 if (columnOne) {  firstCol = pieces[i]; columnOne=false; }
1921                 else  { secondCol = pieces[i]; pairDone = true; columnOne=true; }
1922
1923                 if (pairDone) {
1924                     checkName(firstCol);
1925                     checkName(secondCol);
1926
1927                     //parse names into vector
1928                     vector<string> theseNames;
1929                     splitAtComma(secondCol, theseNames);
1930                     for (int i = 0; i < theseNames.size(); i++) {   nameMap[theseNames[i]] = firstCol;  }
1931                     pairDone = false;
1932                 }
1933             }
1934         }
1935
1936                 return nameMap.size();
1937
1938         }
1939         catch(exception& e) {
1940                 errorOut(e, "MothurOut", "readNames");
1941                 exit(1);
1942         }
1943 }
1944 /**********************************************************************************************************************/
1945 int MothurOut::readNames(string namefile, map<string, string>& nameMap, int flip) {
1946         try {
1947                 //open input file
1948                 ifstream in;
1949                 openInputFile(namefile, in);
1950
1951         string rest = "";
1952         char buffer[4096];
1953         bool pairDone = false;
1954         bool columnOne = true;
1955         string firstCol, secondCol;
1956
1957                 while (!in.eof()) {
1958                         if (control_pressed) { break; }
1959
1960             in.read(buffer, 4096);
1961             vector<string> pieces = splitWhiteSpace(rest, buffer, in.gcount());
1962
1963             for (int i = 0; i < pieces.size(); i++) {
1964                 if (columnOne) {  firstCol = pieces[i]; columnOne=false; }
1965                 else  { secondCol = pieces[i]; pairDone = true; columnOne=true; }
1966
1967                 if (pairDone) {
1968                     checkName(firstCol);
1969                     checkName(secondCol);
1970                     nameMap[secondCol] = firstCol;
1971                     pairDone = false;
1972                 }
1973             }
1974                 }
1975                 in.close();
1976
1977         if (rest != "") {
1978             vector<string> pieces = splitWhiteSpace(rest);
1979
1980             for (int i = 0; i < pieces.size(); i++) {
1981                 if (columnOne) {  firstCol = pieces[i]; columnOne=false; }
1982                 else  { secondCol = pieces[i]; pairDone = true; columnOne=true; }
1983
1984                 if (pairDone) {
1985                     checkName(firstCol);
1986                     checkName(secondCol);
1987                     nameMap[secondCol] = firstCol;
1988                     pairDone = false;
1989                 }
1990             }
1991         }
1992
1993                 return nameMap.size();
1994
1995         }
1996         catch(exception& e) {
1997                 errorOut(e, "MothurOut", "readNames");
1998                 exit(1);
1999         }
2000 }
2001 /**********************************************************************************************************************/
2002 int MothurOut::readNames(string namefile, map<string, string>& nameMap, map<string, int>& nameCount) {
2003         try {
2004                 nameMap.clear(); nameCount.clear();
2005                 //open input file
2006                 ifstream in;
2007                 openInputFile(namefile, in);
2008
2009         string rest = "";
2010         char buffer[4096];
2011         bool pairDone = false;
2012         bool columnOne = true;
2013         string firstCol, secondCol;
2014
2015                 while (!in.eof()) {
2016                         if (control_pressed) { break; }
2017
2018             in.read(buffer, 4096);
2019             vector<string> pieces = splitWhiteSpace(rest, buffer, in.gcount());
2020
2021             for (int i = 0; i < pieces.size(); i++) {
2022                 if (columnOne) {  firstCol = pieces[i]; columnOne=false; }
2023                 else  { secondCol = pieces[i]; pairDone = true; columnOne=true; }
2024
2025                 if (pairDone) {
2026                     checkName(firstCol);
2027                     checkName(secondCol);
2028                     //parse names into vector
2029                     vector<string> theseNames;
2030                     splitAtComma(secondCol, theseNames);
2031                     for (int i = 0; i < theseNames.size(); i++) {  nameMap[theseNames[i]] = firstCol;  }
2032                     nameCount[firstCol] = theseNames.size();
2033                     pairDone = false;
2034                 }
2035             }
2036                 }
2037                 in.close();
2038
2039         if (rest != "") {
2040             vector<string> pieces = splitWhiteSpace(rest);
2041
2042             for (int i = 0; i < pieces.size(); i++) {
2043                 if (columnOne) {  firstCol = pieces[i]; columnOne=false; }
2044                 else  { secondCol = pieces[i]; pairDone = true; columnOne=true; }
2045
2046                 if (pairDone) {
2047                     checkName(firstCol);
2048                     checkName(secondCol);
2049                     //parse names into vector
2050                     vector<string> theseNames;
2051                     splitAtComma(secondCol, theseNames);
2052                     for (int i = 0; i < theseNames.size(); i++) {  nameMap[theseNames[i]] = firstCol;  }
2053                     nameCount[firstCol] = theseNames.size();
2054                     pairDone = false;
2055                 }
2056             }
2057
2058         }
2059                 return nameMap.size();
2060
2061         }
2062         catch(exception& e) {
2063                 errorOut(e, "MothurOut", "readNames");
2064                 exit(1);
2065         }
2066 }
2067 /**********************************************************************************************************************/
2068 int MothurOut::readNames(string namefile, map<string, string>& nameMap) {
2069         try {
2070                 //open input file
2071                 ifstream in;
2072                 openInputFile(namefile, in);
2073
2074         string rest = "";
2075         char buffer[4096];
2076         bool pairDone = false;
2077         bool columnOne = true;
2078         string firstCol, secondCol;
2079
2080                 while (!in.eof()) {
2081                         if (control_pressed) { break; }
2082
2083             in.read(buffer, 4096);
2084             vector<string> pieces = splitWhiteSpace(rest, buffer, in.gcount());
2085
2086             for (int i = 0; i < pieces.size(); i++) {
2087                 if (columnOne) {  firstCol = pieces[i]; columnOne=false; }
2088                 else  { secondCol = pieces[i]; pairDone = true; columnOne=true; }
2089
2090                 if (pairDone) {
2091                     checkName(firstCol);
2092                     checkName(secondCol);
2093                     nameMap[firstCol] = secondCol; pairDone = false; }
2094             }
2095                 }
2096                 in.close();
2097
2098         if (rest != "") {
2099             vector<string> pieces = splitWhiteSpace(rest);
2100
2101             for (int i = 0; i < pieces.size(); i++) {
2102                 if (columnOne) {  firstCol = pieces[i]; columnOne=false; }
2103                 else  { secondCol = pieces[i]; pairDone = true; columnOne=true; }
2104
2105                 if (pairDone) {
2106                     checkName(firstCol);
2107                     checkName(secondCol);
2108                     nameMap[firstCol] = secondCol; pairDone = false; }
2109             }
2110         }
2111
2112                 return nameMap.size();
2113
2114         }
2115         catch(exception& e) {
2116                 errorOut(e, "MothurOut", "readNames");
2117                 exit(1);
2118         }
2119 }
2120 /**********************************************************************************************************************/
2121 int MothurOut::readNames(string namefile, map<string, vector<string> >& nameMap) {
2122         try {
2123                 //open input file
2124                 ifstream in;
2125                 openInputFile(namefile, in);
2126
2127         string rest = "";
2128         char buffer[4096];
2129         bool pairDone = false;
2130         bool columnOne = true;
2131         string firstCol, secondCol;
2132
2133                 while (!in.eof()) {
2134                         if (control_pressed) { break; }
2135
2136             in.read(buffer, 4096);
2137             vector<string> pieces = splitWhiteSpace(rest, buffer, in.gcount());
2138
2139             for (int i = 0; i < pieces.size(); i++) {
2140                 if (columnOne) {  firstCol = pieces[i]; columnOne=false; }
2141                 else  { secondCol = pieces[i]; pairDone = true; columnOne=true; }
2142
2143                 if (pairDone) {
2144                     checkName(firstCol);
2145                     checkName(secondCol);
2146                     vector<string> temp;
2147                     splitAtComma(secondCol, temp);
2148                     nameMap[firstCol] = temp;
2149                     pairDone = false;
2150                 }
2151             }
2152                 }
2153                 in.close();
2154
2155         if (rest != "") {
2156             vector<string> pieces = splitWhiteSpace(rest);
2157
2158             for (int i = 0; i < pieces.size(); i++) {
2159                 if (columnOne) {  firstCol = pieces[i]; columnOne=false; }
2160                 else  { secondCol = pieces[i]; pairDone = true; columnOne=true; }
2161
2162                 if (pairDone) {
2163                     checkName(firstCol);
2164                     checkName(secondCol);
2165                     vector<string> temp;
2166                     splitAtComma(secondCol, temp);
2167                     nameMap[firstCol] = temp;
2168                     pairDone = false;
2169                 }
2170             }
2171         }
2172
2173                 return nameMap.size();
2174         }
2175         catch(exception& e) {
2176                 errorOut(e, "MothurOut", "readNames");
2177                 exit(1);
2178         }
2179 }
2180 /**********************************************************************************************************************/
2181 map<string, int> MothurOut::readNames(string namefile) {
2182         try {
2183                 map<string, int> nameMap;
2184
2185                 //open input file
2186                 ifstream in;
2187                 openInputFile(namefile, in);
2188
2189         string rest = "";
2190         char buffer[4096];
2191         bool pairDone = false;
2192         bool columnOne = true;
2193         string firstCol, secondCol;
2194
2195                 while (!in.eof()) {
2196                         if (control_pressed) { break; }
2197
2198             in.read(buffer, 4096);
2199             vector<string> pieces = splitWhiteSpace(rest, buffer, in.gcount());
2200
2201             for (int i = 0; i < pieces.size(); i++) {
2202                 if (columnOne) {  firstCol = pieces[i]; columnOne=false; }
2203                 else  { secondCol = pieces[i]; pairDone = true; columnOne=true; }
2204
2205                 if (pairDone) {
2206                     checkName(firstCol);
2207                     checkName(secondCol);
2208                     int num = getNumNames(secondCol);
2209                     nameMap[firstCol] = num;
2210                     pairDone = false;
2211                 }
2212             }
2213                 }
2214         in.close();
2215
2216         if (rest != "") {
2217             vector<string> pieces = splitWhiteSpace(rest);
2218             for (int i = 0; i < pieces.size(); i++) {
2219                 if (columnOne) {  firstCol = pieces[i]; columnOne=false; }
2220                 else  { secondCol = pieces[i]; pairDone = true; columnOne=true; }
2221
2222                 if (pairDone) {
2223                     checkName(firstCol);
2224                     checkName(secondCol);
2225                     int num = getNumNames(secondCol);
2226                     nameMap[firstCol] = num;
2227                     pairDone = false;
2228                 }
2229             }
2230         }
2231
2232                 return nameMap;
2233
2234         }
2235         catch(exception& e) {
2236                 errorOut(e, "MothurOut", "readNames");
2237                 exit(1);
2238         }
2239 }
2240 /**********************************************************************************************************************/
2241 map<string, int> MothurOut::readNames(string namefile, unsigned long int& numSeqs) {
2242         try {
2243                 map<string, int> nameMap;
2244         numSeqs = 0;
2245
2246                 //open input file
2247                 ifstream in;
2248                 openInputFile(namefile, in);
2249
2250         string rest = "";
2251         char buffer[4096];
2252         bool pairDone = false;
2253         bool columnOne = true;
2254         string firstCol, secondCol;
2255
2256                 while (!in.eof()) {
2257                         if (control_pressed) { break; }
2258
2259             in.read(buffer, 4096);
2260             vector<string> pieces = splitWhiteSpace(rest, buffer, in.gcount());
2261
2262             for (int i = 0; i < pieces.size(); i++) {
2263                 if (columnOne) {  firstCol = pieces[i]; columnOne=false; }
2264                 else  { secondCol = pieces[i]; pairDone = true; columnOne=true; }
2265
2266                 if (pairDone) {
2267                     checkName(firstCol);
2268                     checkName(secondCol);
2269                     int num = getNumNames(secondCol);
2270                     nameMap[firstCol] = num;
2271                     pairDone = false;
2272                     numSeqs += num;
2273                 }
2274             }
2275                 }
2276         in.close();
2277
2278         if (rest != "") {
2279             vector<string> pieces = splitWhiteSpace(rest);
2280             for (int i = 0; i < pieces.size(); i++) {
2281                 if (columnOne) {  firstCol = pieces[i]; columnOne=false; }
2282                 else  { secondCol = pieces[i]; pairDone = true; columnOne=true; }
2283
2284                 if (pairDone) {
2285                     checkName(firstCol);
2286                     checkName(secondCol);
2287                     int num = getNumNames(secondCol);
2288                     nameMap[firstCol] = num;
2289                     pairDone = false;
2290                     numSeqs += num;
2291                 }
2292             }
2293         }
2294
2295                 return nameMap;
2296
2297         }
2298         catch(exception& e) {
2299                 errorOut(e, "MothurOut", "readNames");
2300                 exit(1);
2301         }
2302 }
2303 /************************************************************/
2304 int MothurOut::checkName(string& name) {
2305     try {
2306         if (modifyNames) {
2307             for (int i = 0; i < name.length(); i++) {
2308                 if (name[i] == ':') { name[i] = '_'; changedSeqNames = true; }
2309             }
2310         }
2311         return 0;
2312     }
2313         catch(exception& e) {
2314                 errorOut(e, "MothurOut", "checkName");
2315                 exit(1);
2316         }
2317 }
2318 /**********************************************************************************************************************/
2319 int MothurOut::readNames(string namefile, vector<seqPriorityNode>& nameVector, map<string, string>& fastamap) {
2320         try {
2321                 int error = 0;
2322
2323                 //open input file
2324                 ifstream in;
2325                 openInputFile(namefile, in);
2326
2327         string rest = "";
2328         char buffer[4096];
2329         bool pairDone = false;
2330         bool columnOne = true;
2331         string firstCol, secondCol;
2332
2333                 while (!in.eof()) {
2334                         if (control_pressed) { break; }
2335
2336             in.read(buffer, 4096);
2337             vector<string> pieces = splitWhiteSpace(rest, buffer, in.gcount());
2338
2339             for (int i = 0; i < pieces.size(); i++) {
2340                 if (columnOne) {  firstCol = pieces[i]; columnOne=false; }
2341                 else  { secondCol = pieces[i]; pairDone = true; columnOne=true; }
2342
2343                 if (pairDone) {
2344                     checkName(firstCol);
2345                     checkName(secondCol);
2346                     int num = getNumNames(secondCol);
2347
2348                     map<string, string>::iterator it = fastamap.find(firstCol);
2349                     if (it == fastamap.end()) {
2350                         error = 1;
2351                         mothurOut("[ERROR]: " + firstCol + " is not in your fastafile, but is in your namesfile, please correct."); mothurOutEndLine();
2352                     }else {
2353                         seqPriorityNode temp(num, it->second, firstCol);
2354                         nameVector.push_back(temp);
2355                     }
2356
2357                     pairDone = false;
2358                 }
2359             }
2360                 }
2361         in.close();
2362
2363         if (rest != "") {
2364             vector<string> pieces = splitWhiteSpace(rest);
2365
2366             for (int i = 0; i < pieces.size(); i++) {
2367                 if (columnOne) {  firstCol = pieces[i]; columnOne=false; }
2368                 else  { secondCol = pieces[i]; pairDone = true; columnOne=true; }
2369
2370                 if (pairDone) {
2371                     checkName(firstCol);
2372                     checkName(secondCol);
2373                     int num = getNumNames(secondCol);
2374
2375                     map<string, string>::iterator it = fastamap.find(firstCol);
2376                     if (it == fastamap.end()) {
2377                         error = 1;
2378                         mothurOut("[ERROR]: " + firstCol + " is not in your fastafile, but is in your namesfile, please correct."); mothurOutEndLine();
2379                     }else {
2380                         seqPriorityNode temp(num, it->second, firstCol);
2381                         nameVector.push_back(temp);
2382                     }
2383
2384                     pairDone = false;
2385                 }
2386             }
2387         }
2388                 return error;
2389         }
2390         catch(exception& e) {
2391                 errorOut(e, "MothurOut", "readNames");
2392                 exit(1);
2393         }
2394 }
2395 //**********************************************************************************************************************
2396 set<string> MothurOut::readAccnos(string accnosfile){
2397         try {
2398                 set<string> names;
2399                 ifstream in;
2400                 openInputFile(accnosfile, in);
2401                 string name;
2402
2403         string rest = "";
2404         char buffer[4096];
2405
2406                 while (!in.eof()) {
2407                         if (control_pressed) { break; }
2408
2409             in.read(buffer, 4096);
2410             vector<string> pieces = splitWhiteSpace(rest, buffer, in.gcount());
2411
2412             for (int i = 0; i < pieces.size(); i++) {  checkName(pieces[i]);
2413                 names.insert(pieces[i]);
2414             }
2415         }
2416                 in.close();
2417
2418         if (rest != "") {
2419             vector<string> pieces = splitWhiteSpace(rest);
2420             for (int i = 0; i < pieces.size(); i++) {  checkName(pieces[i]); names.insert(pieces[i]);  }
2421         }
2422                 return names;
2423         }
2424         catch(exception& e) {
2425                 errorOut(e, "MothurOut", "readAccnos");
2426                 exit(1);
2427         }
2428 }
2429 //**********************************************************************************************************************
2430 int MothurOut::readAccnos(string accnosfile, vector<string>& names){
2431         try {
2432         names.clear();
2433                 ifstream in;
2434                 openInputFile(accnosfile, in);
2435                 string name;
2436
2437         string rest = "";
2438         char buffer[4096];
2439
2440                 while (!in.eof()) {
2441                         if (control_pressed) { break; }
2442
2443             in.read(buffer, 4096);
2444             vector<string> pieces = splitWhiteSpace(rest, buffer, in.gcount());
2445
2446             for (int i = 0; i < pieces.size(); i++) {  checkName(pieces[i]); names.push_back(pieces[i]);  }
2447         }
2448                 in.close();
2449
2450         if (rest != "") {
2451             vector<string> pieces = splitWhiteSpace(rest);
2452             for (int i = 0; i < pieces.size(); i++) {  checkName(pieces[i]); names.push_back(pieces[i]);  }
2453         }
2454
2455                 return 0;
2456         }
2457         catch(exception& e) {
2458                 errorOut(e, "MothurOut", "readAccnos");
2459                 exit(1);
2460         }
2461 }
2462 /***********************************************************************/
2463
2464 int MothurOut::getNumNames(string names){
2465         try {
2466                 int count = 0;
2467
2468                 if(names != ""){
2469                         count = 1;
2470                         for(int i=0;i<names.size();i++){
2471                                 if(names[i] == ','){
2472                                         count++;
2473                                 }
2474                         }
2475                 }
2476
2477                 return count;
2478         }
2479         catch(exception& e) {
2480                 errorOut(e, "MothurOut", "getNumNames");
2481                 exit(1);
2482         }
2483 }
2484 /***********************************************************************/
2485
2486 int MothurOut::getNumChar(string line, char c){
2487         try {
2488                 int count = 0;
2489
2490                 if(line != ""){
2491                         for(int i=0;i<line.size();i++){
2492                                 if(line[i] == c){
2493                                         count++;
2494                                 }
2495                         }
2496                 }
2497
2498                 return count;
2499         }
2500         catch(exception& e) {
2501                 errorOut(e, "MothurOut", "getNumChar");
2502                 exit(1);
2503         }
2504 }
2505 //**********************************************************************************************************************
2506 bool MothurOut::isSubset(vector<string> bigset, vector<string> subset) {
2507         try {
2508
2509
2510                 if (subset.size() > bigset.size()) { return false;  }
2511
2512                 //check if each guy in suset is also in bigset
2513                 for (int i = 0; i < subset.size(); i++) {
2514                         bool match = false;
2515                         for (int j = 0; j < bigset.size(); j++) {
2516                                 if (subset[i] == bigset[j]) { match = true; break; }
2517                         }
2518
2519                         //you have a guy in subset that had no match in bigset
2520                         if (match == false) { return false; }
2521                 }
2522
2523                 return true;
2524
2525         }
2526         catch(exception& e) {
2527                 errorOut(e, "MothurOut", "isSubset");
2528                 exit(1);
2529         }
2530 }
2531 /***********************************************************************/
2532 int MothurOut::mothurRemove(string filename){
2533         try {
2534                 filename = getFullPathName(filename);
2535                 int error = remove(filename.c_str());
2536                 //if (error != 0) {
2537                 //      if (errno != ENOENT) { //ENOENT == file does not exist
2538                 //              string message = "Error deleting file " + filename;
2539                 //              perror(message.c_str());
2540                 //      }
2541                 //}
2542                 return error;
2543         }
2544         catch(exception& e) {
2545                 errorOut(e, "MothurOut", "mothurRemove");
2546                 exit(1);
2547         }
2548 }
2549 /***********************************************************************/
2550 bool MothurOut::mothurConvert(string item, int& num){
2551         try {
2552                 bool error = false;
2553
2554                 if (isNumeric1(item)) {
2555                         convert(item, num);
2556                 }else {
2557                         num = 0;
2558                         error = true;
2559                         mothurOut("[ERROR]: cannot convert " + item + " to an integer."); mothurOutEndLine();
2560                         commandInputsConvertError = true;
2561                 }
2562
2563                 return error;
2564         }
2565         catch(exception& e) {
2566                 errorOut(e, "MothurOut", "mothurConvert");
2567                 exit(1);
2568         }
2569 }
2570 /***********************************************************************/
2571 bool MothurOut::mothurConvert(string item, intDist& num){
2572         try {
2573                 bool error = false;
2574
2575                 if (isNumeric1(item)) {
2576                         convert(item, num);
2577                 }else {
2578                         num = 0;
2579                         error = true;
2580                         mothurOut("[ERROR]: cannot convert " + item + " to an integer."); mothurOutEndLine();
2581                         commandInputsConvertError = true;
2582                 }
2583
2584                 return error;
2585         }
2586         catch(exception& e) {
2587                 errorOut(e, "MothurOut", "mothurConvert");
2588                 exit(1);
2589         }
2590 }
2591
2592 /***********************************************************************/
2593 bool MothurOut::isNumeric1(string stringToCheck){
2594         try {
2595                 bool numeric = false;
2596
2597                 if(stringToCheck.find_first_not_of("0123456789.-") == string::npos) { numeric = true; }
2598
2599                 return numeric;
2600         }
2601         catch(exception& e) {
2602                 errorOut(e, "MothurOut", "isNumeric1");
2603                 exit(1);
2604         }
2605
2606 }
2607 /***********************************************************************/
2608 bool MothurOut::mothurConvert(string item, float& num){
2609         try {
2610                 bool error = false;
2611
2612                 if (isNumeric1(item)) {
2613                         convert(item, num);
2614                 }else {
2615                         num = 0;
2616                         error = true;
2617                         mothurOut("[ERROR]: cannot convert " + item + " to a float."); mothurOutEndLine();
2618                         commandInputsConvertError = true;
2619                 }
2620
2621                 return error;
2622         }
2623         catch(exception& e) {
2624                 errorOut(e, "MothurOut", "mothurConvert");
2625                 exit(1);
2626         }
2627 }
2628 /***********************************************************************/
2629 bool MothurOut::mothurConvert(string item, double& num){
2630         try {
2631                 bool error = false;
2632
2633                 if (isNumeric1(item)) {
2634                         convert(item, num);
2635                 }else {
2636                         num = 0;
2637                         error = true;
2638                         mothurOut("[ERROR]: cannot convert " + item + " to a double."); mothurOutEndLine();
2639                         commandInputsConvertError = true;
2640                 }
2641
2642                 return error;
2643         }
2644         catch(exception& e) {
2645                 errorOut(e, "MothurOut", "mothurConvert");
2646                 exit(1);
2647         }
2648 }
2649 /**************************************************************************************************/
2650
2651 vector<vector<double> > MothurOut::binomial(int maxOrder){
2652         try {
2653         vector<vector<double> > binomial(maxOrder+1);
2654
2655     for(int i=0;i<=maxOrder;i++){
2656                 binomial[i].resize(maxOrder+1);
2657                 binomial[i][0]=1;
2658                 binomial[0][i]=0;
2659     }
2660     binomial[0][0]=1;
2661
2662     binomial[1][0]=1;
2663     binomial[1][1]=1;
2664
2665     for(int i=2;i<=maxOrder;i++){
2666                 binomial[1][i]=0;
2667     }
2668
2669     for(int i=2;i<=maxOrder;i++){
2670                 for(int j=1;j<=maxOrder;j++){
2671                         if(i==j){       binomial[i][j]=1;                                                                       }
2672                         if(j>i) {       binomial[i][j]=0;                                                                       }
2673                         else    {       binomial[i][j]=binomial[i-1][j-1]+binomial[i-1][j];     }
2674                 }
2675     }
2676
2677         return binomial;
2678
2679         }
2680         catch(exception& e) {
2681                 errorOut(e, "MothurOut", "binomial");
2682                 exit(1);
2683         }
2684 }
2685 /**************************************************************************************************/
2686 unsigned int MothurOut::fromBase36(string base36){
2687         try {
2688                 unsigned int num = 0;
2689
2690                 map<char, int> converts;
2691                 converts['A'] = 0;
2692                 converts['a'] = 0;
2693                 converts['B'] = 1;
2694                 converts['b'] = 1;
2695                 converts['C'] = 2;
2696                 converts['c'] = 2;
2697                 converts['D'] = 3;
2698                 converts['d'] = 3;
2699                 converts['E'] = 4;
2700                 converts['e'] = 4;
2701                 converts['F'] = 5;
2702                 converts['f'] = 5;
2703                 converts['G'] = 6;
2704                 converts['g'] = 6;
2705                 converts['H'] = 7;
2706                 converts['h'] = 7;
2707                 converts['I'] = 8;
2708                 converts['i'] = 8;
2709                 converts['J'] = 9;
2710                 converts['j'] = 9;
2711                 converts['K'] = 10;
2712                 converts['k'] = 10;
2713                 converts['L'] = 11;
2714                 converts['l'] = 11;
2715                 converts['M'] = 12;
2716                 converts['m'] = 12;
2717                 converts['N'] = 13;
2718                 converts['n'] = 13;
2719                 converts['O'] = 14;
2720                 converts['o'] = 14;
2721                 converts['P'] = 15;
2722                 converts['p'] = 15;
2723                 converts['Q'] = 16;
2724                 converts['q'] = 16;
2725                 converts['R'] = 17;
2726                 converts['r'] = 17;
2727                 converts['S'] = 18;
2728                 converts['s'] = 18;
2729                 converts['T'] = 19;
2730                 converts['t'] = 19;
2731                 converts['U'] = 20;
2732                 converts['u'] = 20;
2733                 converts['V'] = 21;
2734                 converts['v'] = 21;
2735                 converts['W'] = 22;
2736                 converts['w'] = 22;
2737                 converts['X'] = 23;
2738                 converts['x'] = 23;
2739                 converts['Y'] = 24;
2740                 converts['y'] = 24;
2741                 converts['Z'] = 25;
2742                 converts['z'] = 25;
2743                 converts['0'] = 26;
2744                 converts['1'] = 27;
2745                 converts['2'] = 28;
2746                 converts['3'] = 29;
2747                 converts['4'] = 30;
2748                 converts['5'] = 31;
2749                 converts['6'] = 32;
2750                 converts['7'] = 33;
2751                 converts['8'] = 34;
2752                 converts['9'] = 35;
2753
2754                 int i = 0;
2755                 while (i < base36.length()) {
2756                         char c = base36[i];
2757                         num = 36 * num + converts[c];
2758                         i++;
2759                 }
2760
2761                 return num;
2762
2763         }
2764         catch(exception& e) {
2765                 errorOut(e, "MothurOut", "fromBase36");
2766                 exit(1);
2767         }
2768 }
2769 /***********************************************************************/
2770
2771 int MothurOut::factorial(int num){
2772         try {
2773                 int total = 1;
2774
2775                 for (int i = 1; i <= num; i++) {
2776                         total *= i;
2777                 }
2778
2779                 return total;
2780         }
2781         catch(exception& e) {
2782                 errorOut(e, "MothurOut", "factorial");
2783                 exit(1);
2784         }
2785 }
2786 /***********************************************************************/
2787
2788 int MothurOut::getNumSeqs(ifstream& file){
2789         try {
2790                 int numSeqs = count(istreambuf_iterator<char>(file),istreambuf_iterator<char>(), '>');
2791                 file.seekg(0);
2792                 return numSeqs;
2793         }
2794         catch(exception& e) {
2795                 errorOut(e, "MothurOut", "getNumSeqs");
2796                 exit(1);
2797         }
2798 }
2799 /***********************************************************************/
2800 void MothurOut::getNumSeqs(ifstream& file, int& numSeqs){
2801         try {
2802                 string input;
2803                 numSeqs = 0;
2804                 while(!file.eof()){
2805                         input = getline(file);
2806                         if (input.length() != 0) {
2807                                 if(input[0] == '>'){ numSeqs++; }
2808                         }
2809                 }
2810         }
2811         catch(exception& e) {
2812                 errorOut(e, "MothurOut", "getNumSeqs");
2813                 exit(1);
2814         }
2815 }
2816 /***********************************************************************/
2817
2818 //This function parses the estimator options and puts them in a vector
2819 void MothurOut::splitAtChar(string& estim, vector<string>& container, char symbol) {
2820         try {
2821
2822         if (symbol == '-') { splitAtDash(estim, container); return; }
2823
2824                 string individual = "";
2825                 int estimLength = estim.size();
2826                 for(int i=0;i<estimLength;i++){
2827                         if(estim[i] == symbol){
2828                                 container.push_back(individual);
2829                                 individual = "";
2830                         }
2831                         else{
2832                                 individual += estim[i];
2833                         }
2834                 }
2835                 container.push_back(individual);
2836
2837         }
2838         catch(exception& e) {
2839                 errorOut(e, "MothurOut", "splitAtChar");
2840                 exit(1);
2841         }
2842 }
2843
2844 /***********************************************************************/
2845
2846 //This function parses the estimator options and puts them in a vector
2847 void MothurOut::splitAtDash(string& estim, vector<string>& container) {
2848         try {
2849                 string individual = "";
2850                 int estimLength = estim.size();
2851                 bool prevEscape = false;
2852                 /*for(int i=0;i<estimLength;i++){
2853                         if(prevEscape){
2854                                 individual += estim[i];
2855                                 prevEscape = false;
2856                         }
2857                         else{
2858                                 if(estim[i] == '\\'){
2859                                         prevEscape = true;
2860                                 }
2861                                 else if(estim[i] == '-'){
2862                                         container.push_back(individual);
2863                                         individual = "";
2864                                         prevEscape = false;
2865                                 }
2866                                 else{
2867                                         individual += estim[i];
2868                                         prevEscape = false;
2869                                 }
2870                         }
2871                 }*/
2872
2873
2874         for(int i=0;i<estimLength;i++){
2875             if(estim[i] == '-'){
2876                 if (prevEscape) {  individual += estim[i]; prevEscape = false;  } //add in dash because it was escaped.
2877                 else {
2878                     container.push_back(individual);
2879                     individual = "";
2880                 }
2881             }else if(estim[i] == '\\'){
2882                 if (i < estimLength-1) {
2883                     if (estim[i+1] == '-') { prevEscape=true; }  //are you a backslash before a dash, if yes ignore
2884                     else { individual += estim[i]; prevEscape = false;  } //if no, add in
2885                 }else { individual += estim[i]; }
2886             }else {
2887                 individual += estim[i];
2888             }
2889         }
2890
2891
2892
2893                 container.push_back(individual);
2894         }
2895         catch(exception& e) {
2896                 errorOut(e, "MothurOut", "splitAtDash");
2897                 exit(1);
2898         }
2899 }
2900
2901 /***********************************************************************/
2902 //This function parses the label options and puts them in a set
2903 void MothurOut::splitAtDash(string& estim, set<string>& container) {
2904         try {
2905                 string individual = "";
2906                 int estimLength = estim.size();
2907                 bool prevEscape = false;
2908         /*
2909                 for(int i=0;i<estimLength;i++){
2910                         if(prevEscape){
2911                                 individual += estim[i];
2912                                 prevEscape = false;
2913                         }
2914                         else{
2915                                 if(estim[i] == '\\'){
2916                                         prevEscape = true;
2917                                 }
2918                                 else if(estim[i] == '-'){
2919                                         container.insert(individual);
2920                                         individual = "";
2921                                         prevEscape = false;
2922                                 }
2923                                 else{
2924                                         individual += estim[i];
2925                                         prevEscape = false;
2926                                 }
2927                         }
2928                 }
2929                 */
2930
2931         for(int i=0;i<estimLength;i++){
2932             if(estim[i] == '-'){
2933                 if (prevEscape) {  individual += estim[i]; prevEscape = false;  } //add in dash because it was escaped.
2934                 else {
2935                     container.insert(individual);
2936                     individual = "";
2937                 }
2938             }else if(estim[i] == '\\'){
2939                 if (i < estimLength-1) {
2940                     if (estim[i+1] == '-') { prevEscape=true; }  //are you a backslash before a dash, if yes ignore
2941                     else { individual += estim[i]; prevEscape = false;  } //if no, add in
2942                 }else { individual += estim[i]; }
2943             }else {
2944                 individual += estim[i];
2945             }
2946         }
2947         container.insert(individual);
2948
2949         }
2950         catch(exception& e) {
2951                 errorOut(e, "MothurOut", "splitAtDash");
2952                 exit(1);
2953         }
2954 }
2955 /***********************************************************************/
2956 //This function parses the line options and puts them in a set
2957 void MothurOut::splitAtDash(string& estim, set<int>& container) {
2958         try {
2959                 string individual = "";
2960                 int lineNum;
2961                 int estimLength = estim.size();
2962                 bool prevEscape = false;
2963         /*
2964                 for(int i=0;i<estimLength;i++){
2965                         if(prevEscape){
2966                                 individual += estim[i];
2967                                 prevEscape = false;
2968                         }
2969                         else{
2970                                 if(estim[i] == '\\'){
2971                                         prevEscape = true;
2972                                 }
2973                                 else if(estim[i] == '-'){
2974                                         convert(individual, lineNum); //convert the string to int
2975                                         container.insert(lineNum);
2976                                         individual = "";
2977                                         prevEscape = false;
2978                                 }
2979                                 else{
2980                                         individual += estim[i];
2981                                         prevEscape = false;
2982                                 }
2983                         }
2984                 }*/
2985
2986         for(int i=0;i<estimLength;i++){
2987             if(estim[i] == '-'){
2988                 if (prevEscape) {  individual += estim[i]; prevEscape = false;  } //add in dash because it was escaped.
2989                 else {
2990                     convert(individual, lineNum); //convert the string to int
2991                     container.insert(lineNum);
2992                     individual = "";
2993                 }
2994             }else if(estim[i] == '\\'){
2995                 if (i < estimLength-1) {
2996                     if (estim[i+1] == '-') { prevEscape=true; }  //are you a backslash before a dash, if yes ignore
2997                     else { individual += estim[i]; prevEscape = false;  } //if no, add in
2998                 }else { individual += estim[i]; }
2999             }else {
3000                 individual += estim[i];
3001             }
3002         }
3003
3004                 convert(individual, lineNum); //convert the string to int
3005                 container.insert(lineNum);
3006         }
3007         catch(exception& e) {
3008                 errorOut(e, "MothurOut", "splitAtDash");
3009                 exit(1);
3010         }
3011 }
3012
3013 /***********************************************************************/
3014 string MothurOut::makeList(vector<string>& names) {
3015         try {
3016                 string list = "";
3017
3018         if (names.size() == 0) { return list; }
3019
3020         for (int i = 0; i < names.size()-1; i++) { list += names[i] + ",";  }
3021
3022         //get last name
3023         list += names[names.size()-1];
3024
3025         return list;
3026     }
3027         catch(exception& e) {
3028                 errorOut(e, "MothurOut", "makeList");
3029                 exit(1);
3030         }
3031 }
3032
3033 /***********************************************************************/
3034 //This function parses the a string and puts peices in a vector
3035 void MothurOut::splitAtComma(string& estim, vector<string>& container) {
3036         try {
3037                 string individual = "";
3038                 int estimLength = estim.size();
3039                 for(int i=0;i<estimLength;i++){
3040                         if(estim[i] == ','){
3041                                 container.push_back(individual);
3042                                 individual = "";
3043                         }
3044                         else{
3045                                 individual += estim[i];
3046                         }
3047                 }
3048                 container.push_back(individual);
3049
3050
3051
3052
3053 //              string individual;
3054 //
3055 //              while (estim.find_first_of(',') != -1) {
3056 //                      individual = estim.substr(0,estim.find_first_of(','));
3057 //                      if ((estim.find_first_of(',')+1) <= estim.length()) { //checks to make sure you don't have comma at end of string
3058 //                              estim = estim.substr(estim.find_first_of(',')+1, estim.length());
3059 //                              container.push_back(individual);
3060 //                      }
3061 //              }
3062 //              //get last one
3063 //              container.push_back(estim);
3064         }
3065         catch(exception& e) {
3066                 errorOut(e, "MothurOut", "splitAtComma");
3067                 exit(1);
3068         }
3069 }
3070 /***********************************************************************/
3071 //This function splits up the various option parameters
3072 void MothurOut::splitAtChar(string& prefix, string& suffix, char c){
3073         try {
3074                 prefix = suffix.substr(0,suffix.find_first_of(c));
3075                 if ((suffix.find_first_of(c)+2) <= suffix.length()) {  //checks to make sure you don't have comma at end of string
3076                         suffix = suffix.substr(suffix.find_first_of(c)+1, suffix.length());
3077                         string space = " ";
3078                         while(suffix.at(0) == ' ')
3079                                 suffix = suffix.substr(1, suffix.length());
3080                 }else {  suffix = "";  }
3081
3082     }
3083         catch(exception& e) {
3084                 errorOut(e, "MothurOut", "splitAtChar");
3085                 exit(1);
3086         }
3087 }
3088
3089 /***********************************************************************/
3090
3091 //This function splits up the various option parameters
3092 void MothurOut::splitAtComma(string& prefix, string& suffix){
3093         try {
3094                 prefix = suffix.substr(0,suffix.find_first_of(','));
3095                 if ((suffix.find_first_of(',')+2) <= suffix.length()) {  //checks to make sure you don't have comma at end of string
3096                         suffix = suffix.substr(suffix.find_first_of(',')+1, suffix.length());
3097                         string space = " ";
3098                         while(suffix.at(0) == ' ')
3099                                 suffix = suffix.substr(1, suffix.length());
3100                 }else {  suffix = "";  }
3101
3102         }
3103         catch(exception& e) {
3104                 errorOut(e, "MothurOut", "splitAtComma");
3105                 exit(1);
3106         }
3107 }
3108 /***********************************************************************/
3109
3110 //This function separates the key value from the option value i.e. dist=96_...
3111 void MothurOut::splitAtEquals(string& key, string& value){
3112         try {
3113                 if(value.find_first_of('=') != -1){
3114                         key = value.substr(0,value.find_first_of('='));
3115                         if ((value.find_first_of('=')+1) <= value.length()) {
3116                                 value = value.substr(value.find_first_of('=')+1, value.length());
3117                         }
3118                 }else{
3119                         key = value;
3120                         value = 1;
3121                 }
3122         }
3123         catch(exception& e) {
3124                 errorOut(e, "MothurOut", "splitAtEquals");
3125                 exit(1);
3126         }
3127 }
3128
3129 /**************************************************************************************************/
3130
3131 bool MothurOut::inUsersGroups(string groupname, vector<string> Groups) {
3132         try {
3133                 for (int i = 0; i < Groups.size(); i++) {
3134                         if (groupname == Groups[i]) { return true; }
3135                 }
3136                 return false;
3137         }
3138         catch(exception& e) {
3139                 errorOut(e, "MothurOut", "inUsersGroups");
3140                 exit(1);
3141         }
3142 }
3143 /**************************************************************************************************/
3144
3145 bool MothurOut::inUsersGroups(vector<int> set, vector< vector<int> > sets) {
3146         try {
3147                 for (int i = 0; i < sets.size(); i++) {
3148                         if (set == sets[i]) { return true; }
3149                 }
3150                 return false;
3151         }
3152         catch(exception& e) {
3153                 errorOut(e, "MothurOut", "inUsersGroups");
3154                 exit(1);
3155         }
3156 }
3157 /**************************************************************************************************/
3158
3159 bool MothurOut::inUsersGroups(int groupname, vector<int> Groups) {
3160         try {
3161                 for (int i = 0; i < Groups.size(); i++) {
3162                         if (groupname == Groups[i]) { return true; }
3163                 }
3164                 return false;
3165         }
3166         catch(exception& e) {
3167                 errorOut(e, "MothurOut", "inUsersGroups");
3168                 exit(1);
3169         }
3170 }
3171
3172 /**************************************************************************************************/
3173 //returns true if any of the strings in first vector are in second vector
3174 bool MothurOut::inUsersGroups(vector<string> groupnames, vector<string> Groups) {
3175         try {
3176
3177                 for (int i = 0; i < groupnames.size(); i++) {
3178                         if (inUsersGroups(groupnames[i], Groups)) { return true; }
3179                 }
3180                 return false;
3181         }
3182         catch(exception& e) {
3183                 errorOut(e, "MothurOut", "inUsersGroups");
3184                 exit(1);
3185         }
3186 }
3187 /***********************************************************************/
3188 //this function determines if the user has given us labels that are smaller than the given label.
3189 //if so then it returns true so that the calling function can run the previous valid distance.
3190 //it's a "smart" distance function.  It also checks for invalid labels.
3191 bool MothurOut::anyLabelsToProcess(string label, set<string>& userLabels, string errorOff) {
3192         try {
3193
3194                 set<string>::iterator it;
3195                 vector<float> orderFloat;
3196                 map<string, float> userMap;  //the conversion process removes trailing 0's which we need to put back
3197                 map<string, float>::iterator it2;
3198                 float labelFloat;
3199                 bool smaller = false;
3200
3201                 //unique is the smallest line
3202                 if (label == "unique") {  return false;  }
3203                 else {
3204                         if (convertTestFloat(label, labelFloat)) {
3205                                 convert(label, labelFloat);
3206                         }else { //cant convert
3207                                 return false;
3208                         }
3209                 }
3210
3211                 //go through users set and make them floats
3212                 for(it = userLabels.begin(); it != userLabels.end();) {
3213
3214                         float temp;
3215                         if ((*it != "unique") && (convertTestFloat(*it, temp) == true)){
3216                                 convert(*it, temp);
3217                                 orderFloat.push_back(temp);
3218                                 userMap[*it] = temp;
3219                                 it++;
3220                         }else if (*it == "unique") {
3221                                 orderFloat.push_back(-1.0);
3222                                 userMap["unique"] = -1.0;
3223                                 it++;
3224                         }else {
3225                                 if (errorOff == "") {  mothurOut(*it + " is not a valid label."); mothurOutEndLine();  }
3226                                 userLabels.erase(it++);
3227                         }
3228                 }
3229
3230                 //sort order
3231                 sort(orderFloat.begin(), orderFloat.end());
3232
3233                 /*************************************************/
3234                 //is this label bigger than any of the users labels
3235                 /*************************************************/
3236
3237                 //loop through order until you find a label greater than label
3238                 for (int i = 0; i < orderFloat.size(); i++) {
3239                         if (orderFloat[i] < labelFloat) {
3240                                 smaller = true;
3241                                 if (orderFloat[i] == -1) {
3242                                         if (errorOff == "") { mothurOut("Your file does not include the label unique."); mothurOutEndLine(); }
3243                                         userLabels.erase("unique");
3244                                 }
3245                                 else {
3246                                         if (errorOff == "") { mothurOut("Your file does not include the label "); mothurOutEndLine(); }
3247                                         string s = "";
3248                                         for (it2 = userMap.begin(); it2!= userMap.end(); it2++) {
3249                                                 if (it2->second == orderFloat[i]) {
3250                                                         s = it2->first;
3251                                                         //remove small labels
3252                                                         userLabels.erase(s);
3253                                                         break;
3254                                                 }
3255                                         }
3256                                         if (errorOff == "") {mothurOut( s +  ". I will use the next smallest distance. "); mothurOutEndLine(); }
3257                                 }
3258                         //since they are sorted once you find a bigger one stop looking
3259                         }else { break; }
3260                 }
3261
3262                 return smaller;
3263
3264         }
3265         catch(exception& e) {
3266                 errorOut(e, "MothurOut", "anyLabelsToProcess");
3267                 exit(1);
3268         }
3269 }
3270
3271 /**************************************************************************************************/
3272 bool MothurOut::checkReleaseVersion(ifstream& file, string version) {
3273         try {
3274
3275                 bool good = true;
3276
3277                 string line = getline(file);
3278
3279                 //before we added this check
3280                 if (line[0] != '#') {  good = false;  }
3281                 else {
3282                         //rip off #
3283                         line = line.substr(1);
3284
3285                         vector<string> versionVector;
3286                         splitAtChar(version, versionVector, '.');
3287
3288                         //check file version
3289                         vector<string> linesVector;
3290                         splitAtChar(line, linesVector, '.');
3291
3292                         if (versionVector.size() != linesVector.size()) { good = false; }
3293                         else {
3294                                 for (int j = 0; j < versionVector.size(); j++) {
3295                                         int num1, num2;
3296                                         convert(versionVector[j], num1);
3297                                         convert(linesVector[j], num2);
3298
3299                                         //if mothurs version is newer than this files version, then we want to remake it
3300                                         if (num1 > num2) {  good = false; break;  }
3301                                 }
3302                         }
3303
3304                 }
3305
3306                 if (!good) {  file.close();  }
3307                 else { file.seekg(0);  }
3308
3309                 return good;
3310         }
3311         catch(exception& e) {
3312                 errorOut(e, "MothurOut", "checkReleaseVersion");
3313                 exit(1);
3314         }
3315 }
3316 /**************************************************************************************************/
3317 vector<double> MothurOut::getAverages(vector< vector<double> >& dists) {
3318         try{
3319         vector<double> averages; //averages.resize(numComp, 0.0);
3320         for (int i = 0; i < dists[0].size(); i++) { averages.push_back(0.0); }
3321
3322         for (int thisIter = 0; thisIter < dists.size(); thisIter++) {
3323             for (int i = 0; i < dists[thisIter].size(); i++) {
3324                 averages[i] += dists[thisIter][i];
3325             }
3326         }
3327
3328         //finds average.
3329         for (int i = 0; i < averages.size(); i++) {  averages[i] /= (double) dists.size(); }
3330
3331         return averages;
3332     }
3333         catch(exception& e) {
3334                 errorOut(e, "MothurOut", "getAverages");
3335                 exit(1);
3336         }
3337 }
3338 /**************************************************************************************************/
3339 double MothurOut::getAverage(vector<double> dists) {
3340         try{
3341         double average = 0;
3342
3343         for (int i = 0; i < dists.size(); i++) {
3344             average += dists[i];
3345         }
3346
3347         //finds average.
3348         average /= (double) dists.size();
3349
3350         return average;
3351     }
3352         catch(exception& e) {
3353                 errorOut(e, "MothurOut", "getAverage");
3354                 exit(1);
3355         }
3356 }
3357
3358 /**************************************************************************************************/
3359 vector<double> MothurOut::getStandardDeviation(vector< vector<double> >& dists) {
3360         try{
3361
3362         vector<double> averages = getAverages(dists);
3363
3364         //find standard deviation
3365         vector<double> stdDev; //stdDev.resize(numComp, 0.0);
3366         for (int i = 0; i < dists[0].size(); i++) { stdDev.push_back(0.0); }
3367
3368         for (int thisIter = 0; thisIter < dists.size(); thisIter++) { //compute the difference of each dist from the mean, and square the result of each
3369             for (int j = 0; j < dists[thisIter].size(); j++) {
3370                 stdDev[j] += ((dists[thisIter][j] - averages[j]) * (dists[thisIter][j] - averages[j]));
3371             }
3372         }
3373         for (int i = 0; i < stdDev.size(); i++) {
3374             stdDev[i] /= (double) dists.size();
3375             stdDev[i] = sqrt(stdDev[i]);
3376         }
3377
3378         return stdDev;
3379     }
3380         catch(exception& e) {
3381                 errorOut(e, "MothurOut", "getAverages");
3382                 exit(1);
3383         }
3384 }
3385 /**************************************************************************************************/
3386 vector<double> MothurOut::getStandardDeviation(vector< vector<double> >& dists, vector<double>& averages) {
3387         try{
3388         //find standard deviation
3389         vector<double> stdDev; //stdDev.resize(numComp, 0.0);
3390         for (int i = 0; i < dists[0].size(); i++) { stdDev.push_back(0.0); }
3391
3392         for (int thisIter = 0; thisIter < dists.size(); thisIter++) { //compute the difference of each dist from the mean, and square the result of each
3393             for (int j = 0; j < dists[thisIter].size(); j++) {
3394                 stdDev[j] += ((dists[thisIter][j] - averages[j]) * (dists[thisIter][j] - averages[j]));
3395             }
3396         }
3397         for (int i = 0; i < stdDev.size(); i++) {
3398             stdDev[i] /= (double) dists.size();
3399             stdDev[i] = sqrt(stdDev[i]);
3400         }
3401
3402         return stdDev;
3403     }
3404         catch(exception& e) {
3405                 errorOut(e, "MothurOut", "getAverages");
3406                 exit(1);
3407         }
3408 }
3409 /**************************************************************************************************/
3410 vector< vector<seqDist> > MothurOut::getAverages(vector< vector< vector<seqDist> > >& calcDistsTotals, string mode) {
3411         try{
3412
3413         vector< vector<seqDist>  > calcAverages; //calcAverages.resize(calcDistsTotals[0].size());
3414         for (int i = 0; i < calcDistsTotals[0].size(); i++) {  //initialize sums to zero.
3415             //calcAverages[i].resize(calcDistsTotals[0][i].size());
3416             vector<seqDist> temp;
3417             for (int j = 0; j < calcDistsTotals[0][i].size(); j++) {
3418                 seqDist tempDist;
3419                 tempDist.seq1 = calcDistsTotals[0][i][j].seq1;
3420                 tempDist.seq2 = calcDistsTotals[0][i][j].seq2;
3421                 tempDist.dist = 0.0;
3422                 temp.push_back(tempDist);
3423             }
3424             calcAverages.push_back(temp);
3425         }
3426
3427         if (mode == "average") {
3428             for (int thisIter = 0; thisIter < calcDistsTotals.size(); thisIter++) { //sum all groups dists for each calculator
3429                 for (int i = 0; i < calcAverages.size(); i++) {  //initialize sums to zero.
3430                     for (int j = 0; j < calcAverages[i].size(); j++) {
3431                         calcAverages[i][j].dist += calcDistsTotals[thisIter][i][j].dist;
3432                     }
3433                 }
3434             }
3435
3436             for (int i = 0; i < calcAverages.size(); i++) {  //finds average.
3437                 for (int j = 0; j < calcAverages[i].size(); j++) {
3438                     calcAverages[i][j].dist /= (float) calcDistsTotals.size();
3439                 }
3440             }
3441         }else { //find median
3442             for (int i = 0; i < calcAverages.size(); i++) { //for each calc
3443                 for (int j = 0; j < calcAverages[i].size(); j++) {  //for each comparison
3444                     vector<double> dists;
3445                     for (int thisIter = 0; thisIter < calcDistsTotals.size(); thisIter++) { //for each subsample
3446                         dists.push_back(calcDistsTotals[thisIter][i][j].dist);
3447                     }
3448                     sort(dists.begin(), dists.end());
3449                     calcAverages[i][j].dist = dists[(calcDistsTotals.size()/2)];
3450                 }
3451             }
3452         }
3453
3454         return calcAverages;
3455     }
3456         catch(exception& e) {
3457                 errorOut(e, "MothurOut", "getAverages");
3458                 exit(1);
3459         }
3460 }
3461 /**************************************************************************************************/
3462 vector< vector<seqDist> > MothurOut::getAverages(vector< vector< vector<seqDist> > >& calcDistsTotals) {
3463         try{
3464
3465         vector< vector<seqDist>  > calcAverages; //calcAverages.resize(calcDistsTotals[0].size());
3466         for (int i = 0; i < calcDistsTotals[0].size(); i++) {  //initialize sums to zero.
3467             //calcAverages[i].resize(calcDistsTotals[0][i].size());
3468             vector<seqDist> temp;
3469             for (int j = 0; j < calcDistsTotals[0][i].size(); j++) {
3470                 seqDist tempDist;
3471                 tempDist.seq1 = calcDistsTotals[0][i][j].seq1;
3472                 tempDist.seq2 = calcDistsTotals[0][i][j].seq2;
3473                 tempDist.dist = 0.0;
3474                 temp.push_back(tempDist);
3475             }
3476             calcAverages.push_back(temp);
3477         }
3478
3479
3480         for (int thisIter = 0; thisIter < calcDistsTotals.size(); thisIter++) { //sum all groups dists for each calculator
3481                 for (int i = 0; i < calcAverages.size(); i++) {  //initialize sums to zero.
3482                     for (int j = 0; j < calcAverages[i].size(); j++) {
3483                         calcAverages[i][j].dist += calcDistsTotals[thisIter][i][j].dist;
3484                     }
3485                 }
3486         }
3487
3488         for (int i = 0; i < calcAverages.size(); i++) {  //finds average.
3489                 for (int j = 0; j < calcAverages[i].size(); j++) {
3490                     calcAverages[i][j].dist /= (float) calcDistsTotals.size();
3491                 }
3492         }
3493
3494         return calcAverages;
3495     }
3496         catch(exception& e) {
3497                 errorOut(e, "MothurOut", "getAverages");
3498                 exit(1);
3499         }
3500 }
3501 /**************************************************************************************************/
3502 vector< vector<seqDist> > MothurOut::getStandardDeviation(vector< vector< vector<seqDist> > >& calcDistsTotals) {
3503         try{
3504
3505         vector< vector<seqDist> > calcAverages = getAverages(calcDistsTotals);
3506
3507         //find standard deviation
3508         vector< vector<seqDist>  > stdDev;
3509         for (int i = 0; i < calcDistsTotals[0].size(); i++) {  //initialize sums to zero.
3510             vector<seqDist> temp;
3511             for (int j = 0; j < calcDistsTotals[0][i].size(); j++) {
3512                 seqDist tempDist;
3513                 tempDist.seq1 = calcDistsTotals[0][i][j].seq1;
3514                 tempDist.seq2 = calcDistsTotals[0][i][j].seq2;
3515                 tempDist.dist = 0.0;
3516                 temp.push_back(tempDist);
3517             }
3518             stdDev.push_back(temp);
3519         }
3520
3521         for (int thisIter = 0; thisIter < calcDistsTotals.size(); thisIter++) { //compute the difference of each dist from the mean, and square the result of each
3522             for (int i = 0; i < stdDev.size(); i++) {
3523                 for (int j = 0; j < stdDev[i].size(); j++) {
3524                     stdDev[i][j].dist += ((calcDistsTotals[thisIter][i][j].dist - calcAverages[i][j].dist) * (calcDistsTotals[thisIter][i][j].dist - calcAverages[i][j].dist));
3525                 }
3526             }
3527         }
3528
3529         for (int i = 0; i < stdDev.size(); i++) {  //finds average.
3530             for (int j = 0; j < stdDev[i].size(); j++) {
3531                 stdDev[i][j].dist /= (float) calcDistsTotals.size();
3532                 stdDev[i][j].dist = sqrt(stdDev[i][j].dist);
3533             }
3534         }
3535
3536         return stdDev;
3537     }
3538         catch(exception& e) {
3539                 errorOut(e, "MothurOut", "getAverages");
3540                 exit(1);
3541         }
3542 }
3543 /**************************************************************************************************/
3544 vector< vector<seqDist> > MothurOut::getStandardDeviation(vector< vector< vector<seqDist> > >& calcDistsTotals, vector< vector<seqDist> >& calcAverages) {
3545         try{
3546         //find standard deviation
3547         vector< vector<seqDist>  > stdDev;
3548         for (int i = 0; i < calcDistsTotals[0].size(); i++) {  //initialize sums to zero.
3549             vector<seqDist> temp;
3550             for (int j = 0; j < calcDistsTotals[0][i].size(); j++) {
3551                 seqDist tempDist;
3552                 tempDist.seq1 = calcDistsTotals[0][i][j].seq1;
3553                 tempDist.seq2 = calcDistsTotals[0][i][j].seq2;
3554                 tempDist.dist = 0.0;
3555                 temp.push_back(tempDist);
3556             }
3557             stdDev.push_back(temp);
3558         }
3559
3560         for (int thisIter = 0; thisIter < calcDistsTotals.size(); thisIter++) { //compute the difference of each dist from the mean, and square the result of each
3561             for (int i = 0; i < stdDev.size(); i++) {
3562                 for (int j = 0; j < stdDev[i].size(); j++) {
3563                     stdDev[i][j].dist += ((calcDistsTotals[thisIter][i][j].dist - calcAverages[i][j].dist) * (calcDistsTotals[thisIter][i][j].dist - calcAverages[i][j].dist));
3564                 }
3565             }
3566         }
3567
3568         for (int i = 0; i < stdDev.size(); i++) {  //finds average.
3569             for (int j = 0; j < stdDev[i].size(); j++) {
3570                 stdDev[i][j].dist /= (float) calcDistsTotals.size();
3571                 stdDev[i][j].dist = sqrt(stdDev[i][j].dist);
3572             }
3573         }
3574
3575         return stdDev;
3576     }
3577         catch(exception& e) {
3578                 errorOut(e, "MothurOut", "getAverages");
3579                 exit(1);
3580         }
3581 }
3582
3583 /**************************************************************************************************/
3584 bool MothurOut::isContainingOnlyDigits(string input) {
3585         try{
3586
3587                 //are you a digit in ascii code
3588                 for (int i = 0;i < input.length(); i++){
3589                         if( input[i]>47 && input[i]<58){}
3590                         else { return false; }
3591                 }
3592
3593                 return true;
3594         }
3595         catch(exception& e) {
3596                 errorOut(e, "MothurOut", "isContainingOnlyDigits");
3597                 exit(1);
3598         }
3599 }
3600 /**************************************************************************************************/
3601 int MothurOut::removeConfidences(string& tax) {
3602         try {
3603
3604                 string taxon;
3605                 string newTax = "";
3606
3607                 while (tax.find_first_of(';') != -1) {
3608
3609                         if (control_pressed) { return 0; }
3610
3611                         //get taxon
3612                         taxon = tax.substr(0,tax.find_first_of(';'));
3613
3614                         int pos = taxon.find_last_of('(');
3615                         if (pos != -1) {
3616                                 //is it a number?
3617                                 int pos2 = taxon.find_last_of(')');
3618                                 if (pos2 != -1) {
3619                                         string confidenceScore = taxon.substr(pos+1, (pos2-(pos+1)));
3620                                         if (isNumeric1(confidenceScore)) {
3621                                                 taxon = taxon.substr(0, pos); //rip off confidence
3622                                         }
3623                                 }
3624                         }
3625                         taxon += ";";
3626
3627                         tax = tax.substr(tax.find_first_of(';')+1, tax.length());
3628                         newTax += taxon;
3629                 }
3630
3631                 tax = newTax;
3632
3633                 return 0;
3634         }
3635         catch(exception& e) {
3636                 errorOut(e, "MothurOut", "removeConfidences");
3637                 exit(1);
3638         }
3639 }
3640 /**************************************************************************************************/
3641 string MothurOut::removeQuotes(string tax) {
3642         try {
3643
3644                 string taxon;
3645                 string newTax = "";
3646
3647                 for (int i = 0; i < tax.length(); i++) {
3648
3649                         if (control_pressed) { return newTax; }
3650
3651             if ((tax[i] != '\'') && (tax[i] != '\"')) { newTax += tax[i]; }
3652
3653         }
3654
3655                 return newTax;
3656         }
3657         catch(exception& e) {
3658                 errorOut(e, "MothurOut", "removeQuotes");
3659                 exit(1);
3660         }
3661 }
3662 /**************************************************************************************************/
3663 // function for calculating standard deviation
3664 double MothurOut::getStandardDeviation(vector<int>& featureVector){
3665     try {
3666         //finds sum
3667         double average = 0;
3668         for (int i = 0; i < featureVector.size(); i++) { average += featureVector[i]; }
3669         average /= (double) featureVector.size();
3670
3671         //find standard deviation
3672         double stdDev = 0;
3673         for (int i = 0; i < featureVector.size(); i++) { //compute the difference of each dist from the mean, and square the result of each
3674             stdDev += ((featureVector[i] - average) * (featureVector[i] - average));
3675         }
3676
3677         stdDev /= (double) featureVector.size();
3678         stdDev = sqrt(stdDev);
3679
3680         return stdDev;
3681     }
3682         catch(exception& e) {
3683                 errorOut(e, "MothurOut", "getStandardDeviation");
3684                 exit(1);
3685         }
3686 }
3687 /**************************************************************************************************/
3688
3689
3690