mergesfffilecommand.cpp

   1 //
   2 //  mergesfffilecommand.cpp
   3 //  Mothur
   4 //
   5 //  Created by Sarah Westcott on 1/31/14.
   6 //  Copyright (c) 2014 Schloss Lab. All rights reserved.
   7 //
   8
   9 #include "mergesfffilecommand.h"
  10 #include "endiannessmacros.h"
  11
  12 //**********************************************************************************************************************
  13 vector<string> MergeSfffilesCommand::setParameters(){
  14         try {
  15                 CommandParameter psff("sff", "InputTypes", "", "", "sffFile", "sffFile", "none","sff",false,false); parameters.push_back(psff);
  16         CommandParameter pfile("file", "InputTypes", "", "", "sffFile", "sffFile", "none","sff",false,false); parameters.push_back(pfile);
  17                 CommandParameter pinputdir("inputdir", "String", "", "", "", "", "","",false,false); parameters.push_back(pinputdir);
  18                 CommandParameter poutputdir("outputdir", "String", "", "", "", "", "","",false,false); parameters.push_back(poutputdir);
  19
  20                 vector<string> myArray;
  21                 for (int i = 0; i < parameters.size(); i++) {   myArray.push_back(parameters[i].name);          }
  22                 return myArray;
  23         }
  24         catch(exception& e) {
  25                 m->errorOut(e, "MergeSfffilesCommand", "setParameters");
  26                 exit(1);
  27         }
  28 }
  29 //**********************************************************************************************************************
  30 string MergeSfffilesCommand::getHelpString(){
  31         try {
  32                 string helpString = "";
  33                 helpString += "The merge.sfffiles command reads a sff set of files or a file containing a list of sff files and merges the individual files into a single sff file.\n";
  34                 helpString += "The merge.sfffiles command parameters are sff and file. sff or file is required. \n";
  35                 helpString += "The sff parameter allows you to enter the sff list of sff files separated by -'s.\n";
  36                 helpString += "The file parameter allows you to provide a file containing a list of sff files to merge.  \n";
  37                 helpString += "Example sffinfo(sff=mySffFile.sff-mySecond.sff).\n";
  38                 helpString += "Note: No spaces between parameter labels (i.e. sff), '=' and parameters (i.e.yourSffFileName).\n";
  39                 return helpString;
  40         }
  41         catch(exception& e) {
  42                 m->errorOut(e, "MergeSfffilesCommand", "getHelpString");
  43                 exit(1);
  44         }
  45 }
  46
  47 //**********************************************************************************************************************
  48 string MergeSfffilesCommand::getOutputPattern(string type) {
  49     try {
  50         string pattern = "";
  51
  52         if (type == "sff")            {   pattern =  "[filename],merge.sff";   }
  53         else { m->mothurOut("[ERROR]: No definition for type " + type + " output pattern.\n"); m->control_pressed = true;  }
  54
  55         return pattern;
  56     }
  57     catch(exception& e) {
  58         m->errorOut(e, "MergeSfffilesCommand", "getOutputPattern");
  59         exit(1);
  60     }
  61 }
  62 //**********************************************************************************************************************
  63 MergeSfffilesCommand::MergeSfffilesCommand(){
  64         try {
  65                 abort = true; calledHelp = true;
  66                 setParameters();
  67                 vector<string> tempOutNames;
  68         outputTypes["sff"] = tempOutNames;
  69         }
  70         catch(exception& e) {
  71                 m->errorOut(e, "MergeSfffilesCommand", "MergeSfffilesCommand");
  72                 exit(1);
  73         }
  74 }
  75 //**********************************************************************************************************************
  76
  77 MergeSfffilesCommand::MergeSfffilesCommand(string option)  {
  78         try {
  79                 abort = false; calledHelp = false;
  80
  81                 //allow user to run help
  82                 if(option == "help") { help(); abort = true; calledHelp = true; }
  83                 else if(option == "citation") { citation(); abort = true; calledHelp = true;}
  84
  85                 else {
  86                         //valid paramters for this command
  87                         vector<string> myArray = setParameters();
  88
  89                         OptionParser parser(option);
  90                         map<string, string> parameters = parser.getParameters();
  91             map<string,string>::iterator it;
  92
  93                         ValidParameters validParameter;
  94                         //check to make sure all parameters are valid for command
  95                         for (map<string,string>::iterator it = parameters.begin(); it != parameters.end(); it++) {
  96                                 if (validParameter.isValidParameter(it->first, myArray, it->second) != true) {  abort = true;  }
  97                         }
  98
  99                         //initialize outputTypes
 100                         vector<string> tempOutNames;
 101             outputTypes["sff"] = tempOutNames;
 102
 103                         //if the user changes the output directory command factory will send this info to us in the output parameter
 104                         outputDir = validParameter.validFile(parameters, "outputdir", false);           if (outputDir == "not found"){  outputDir = "";         }
 105
 106                         //if the user changes the input directory command factory will send this info to us in the output parameter
 107                         string inputDir = validParameter.validFile(parameters, "inputdir", false);        if (inputDir == "not found"){ inputDir = "";          }
 108             else {
 109                 it = parameters.find("file");
 110                                 //user has given a template file
 111                                 if(it != parameters.end()){
 112                                         string path = m->hasPath(it->second);
 113                                         //if the user has not given a path then, add inputdir. else leave path alone.
 114                                         if (path == "") {       parameters["file"] = inputDir + it->second;             }
 115                                 }
 116             }
 117
 118                         sffFilename = validParameter.validFile(parameters, "sff", false);
 119                         if (sffFilename == "not found") { sffFilename = "";  }
 120                         else {
 121                                 m->splitAtDash(sffFilename, filenames);
 122
 123                                 //go through files and make sure they are good, if not, then disregard them
 124                                 for (int i = 0; i < filenames.size(); i++) {
 125                                         bool ignore = false;
 126                                         if (filenames[i] == "current") {
 127                                                 filenames[i] = m->getSFFFile();
 128                                                 if (filenames[i] != "") {  m->mothurOut("Using " + filenames[i] + " as input file for the sff parameter where you had given current."); m->mothurOutEndLine(); }
 129                                                 else {
 130                                                         m->mothurOut("You have no current sfffile, ignoring current."); m->mothurOutEndLine(); ignore=true;
 131                                                         //erase from file list
 132                                                         filenames.erase(filenames.begin()+i);
 133                                                         i--;
 134                                                 }
 135                                         }
 136
 137                                         if (!ignore) {
 138                                                 if (inputDir != "") {
 139                                                         string path = m->hasPath(filenames[i]);
 140                                                         //if the user has not given a path then, add inputdir. else leave path alone.
 141                                                         if (path == "") {       filenames[i] = inputDir + filenames[i];         }
 142                                                 }
 143
 144                                                 ifstream in;
 145                                                 int ableToOpen = m->openInputFile(filenames[i], in, "noerror");
 146
 147                                                 //if you can't open it, try default location
 148                                                 if (ableToOpen == 1) {
 149                                                         if (m->getDefaultPath() != "") { //default path is set
 150                                                                 string tryPath = m->getDefaultPath() + m->getSimpleName(filenames[i]);
 151                                                                 m->mothurOut("Unable to open " + filenames[i] + ". Trying default " + tryPath); m->mothurOutEndLine();
 152                                                                 ifstream in2;
 153                                                                 ableToOpen = m->openInputFile(tryPath, in2, "noerror");
 154                                                                 in2.close();
 155                                                                 filenames[i] = tryPath;
 156                                                         }
 157                                                 }
 158
 159                                                 //if you can't open it, try default location
 160                                                 if (ableToOpen == 1) {
 161                                                         if (m->getOutputDir() != "") { //default path is set
 162                                                                 string tryPath = m->getOutputDir() + m->getSimpleName(filenames[i]);
 163                                                                 m->mothurOut("Unable to open " + filenames[i] + ". Trying output directory " + tryPath); m->mothurOutEndLine();
 164                                                                 ifstream in2;
 165                                                                 ableToOpen = m->openInputFile(tryPath, in2, "noerror");
 166                                                                 in2.close();
 167                                                                 filenames[i] = tryPath;
 168                                                         }
 169                                                 }
 170
 171                                                 in.close();
 172
 173                                                 if (ableToOpen == 1) {
 174                                                         m->mothurOut("Unable to open " + filenames[i] + ". It will be disregarded."); m->mothurOutEndLine();
 175                                                         //erase from file list
 176                                                         filenames.erase(filenames.begin()+i);
 177                                                         i--;
 178                                                 }else { m->setSFFFile(filenames[i]); }
 179                                         }
 180                                 }
 181                         }
 182
 183                         file = validParameter.validFile(parameters, "file", true);
 184                         if (file == "not open") {  abort = true; }
 185                         else if (file == "not found") { file = "";  }
 186
 187             if ((file == "") && (filenames.size() == 0)) {
 188                 m->mothurOut("[ERROR]: no valid files."); m->mothurOutEndLine(); abort = true;
 189             }
 190
 191             if ((file != "") && (filenames.size() != 0)) { //both are given
 192                 m->mothurOut("[ERROR]: cannot use file option and sff option at the same time, choose one."); m->mothurOutEndLine(); abort = true;
 193             }
 194
 195                 }
 196         }
 197         catch(exception& e) {
 198                 m->errorOut(e, "MergeSfffilesCommand", "MergeSfffilesCommand");
 199                 exit(1);
 200         }
 201 }
 202 //**********************************************************************************************************************
 203 int MergeSfffilesCommand::execute(){
 204         try {
 205                 if (abort == true) { if (calledHelp) { return 0; }  return 2;   }
 206
 207         map<string, string> variables;
 208         if (file != "") {
 209             string thisOutputDir = outputDir;
 210             if (outputDir == "") {  thisOutputDir += m->hasPath(file);  }
 211             variables["[filename]"] = thisOutputDir + m->getRootName(m->getSimpleName(file));
 212             readFile();
 213         }else {
 214             string thisOutputDir = outputDir;
 215             if (outputDir == "") {  thisOutputDir += m->hasPath(filenames[0]);  }
 216             variables["[filename]"] = thisOutputDir + "mergedSff";
 217         }
 218         ofstream out;
 219                 outputFile = getOutputFileName("sff",variables);
 220         m->openOutputFile(outputFile, out);
 221         outputNames.push_back(outputFile); outputTypes["sff"].push_back(outputFile);
 222         outputFileHeader = outputFile + ".headers";
 223         numTotalReads = 0;
 224
 225                 for (int s = 0; s < filenames.size(); s++) {
 226
 227                         if (m->control_pressed) {  for (int i = 0; i < outputNames.size(); i++) {       m->mothurRemove(outputNames[i]);        } return 0; }
 228
 229                         int start = time(NULL);
 230
 231             filenames[s] = m->getFullPathName(filenames[s]);
 232                         m->mothurOut("\nMerging info from " + filenames[s] + " ..." ); m->mothurOutEndLine();
 233
 234                         int numReads = mergeSffInfo(filenames[s], out);
 235
 236                         m->mothurOut("It took " + toString(time(NULL) - start) + " secs to merge " + toString(numReads) + ".\n");
 237                 }
 238         out.close();
 239
 240         //create new common header and add to merged file
 241         adjustCommonHeader();
 242
 243                 if (m->control_pressed) {  for (int i = 0; i < outputNames.size(); i++) {       m->mothurRemove(outputNames[i]);        } return 0; }
 244
 245                 //set sff file as new current sff file
 246                 string current = "";
 247                 itTypes = outputTypes.find("sff");
 248                 if (itTypes != outputTypes.end()) {
 249                         if ((itTypes->second).size() != 0) { current = (itTypes->second)[0]; m->setSFFFile(current); }
 250                 }
 251
 252                 //report output filenames
 253                 m->mothurOutEndLine();
 254                 m->mothurOut("Output File Names: "); m->mothurOutEndLine();
 255                 for (int i = 0; i < outputNames.size(); i++) {  m->mothurOut(outputNames[i]); m->mothurOutEndLine();    }
 256                 m->mothurOutEndLine();
 257
 258                 return 0;
 259         }
 260         catch(exception& e) {
 261                 m->errorOut(e, "MergeSfffilesCommand", "execute");
 262                 exit(1);
 263         }
 264 }
 265 //**********************************************************************************************************************
 266 int MergeSfffilesCommand::mergeSffInfo(string input, ofstream& out){
 267         try {
 268                 currentFileName = input;
 269
 270                 ifstream in;
 271                 m->openInputFileBinary(input, in);
 272
 273                 CommonHeader header;
 274                 readCommonHeader(in, header);
 275
 276                 int count = 0;
 277
 278                 //check magic number and version
 279                 if (header.magicNumber != 779314790) { m->mothurOut("Magic Number is not correct, not a valid .sff file"); m->mothurOutEndLine(); return count; }
 280                 if (header.version != "0001") { m->mothurOut("Version is not supported, only support version 0001."); m->mothurOutEndLine(); return count; }
 281
 282         //save for adjustHeader sanity check
 283         commonHeaders.push_back(header);
 284
 285                 //read through the sff file
 286                 while (!in.eof()) {
 287
 288                         //read data
 289                         seqRead read;  Header readheader;
 290             readSeqData(in, read, header.numFlowsPerRead, readheader, out);
 291
 292             bool okay = sanityCheck(readheader, read);
 293             if (!okay) { break; }
 294
 295                         count++;
 296
 297                         //report progress
 298                         if((count+1) % 10000 == 0){     m->mothurOut(toString(count+1)); m->mothurOutEndLine();         }
 299
 300                         if (m->control_pressed) { count = 0; break;   }
 301
 302                         if (count >= header.numReads) { break; }
 303                 }
 304
 305                 //report progress
 306                 if (!m->control_pressed) {   if((count) % 10000 != 0){  m->mothurOut(toString(count)); m->mothurOutEndLine();           }  }
 307
 308                 in.close();
 309
 310                 return count;
 311         }
 312         catch(exception& e) {
 313                 m->errorOut(e, "MergeSfffilesCommand", "mergeSffInfo");
 314                 exit(1);
 315         }
 316 }
 317 //**********************************************************************************************************************
 318 int MergeSfffilesCommand::readCommonHeader(ifstream& in, CommonHeader& header){
 319         try {
 320
 321                 if (!in.eof()) {
 322
 323                         //read magic number
 324                         char buffer[4];
 325                         in.read(buffer, 4);
 326                         header.magicNumber = be_int4(*(unsigned int *)(&buffer));
 327
 328                         //read version
 329                         char buffer9[4];
 330                         in.read(buffer9, 4);
 331                         header.version = "";
 332                         for (int i = 0; i < 4; i++) {  header.version += toString((int)(buffer9[i]));  }
 333
 334                         //read offset
 335                         char buffer2 [8];
 336                         in.read(buffer2, 8);
 337                         header.indexOffset =  be_int8(*(unsigned long long *)(&buffer2));
 338
 339                         //read index length
 340                         char buffer3 [4];
 341                         in.read(buffer3, 4);
 342                         header.indexLength =  be_int4(*(unsigned int *)(&buffer3));
 343
 344                         //read num reads
 345                         char buffer4 [4];
 346                         in.read(buffer4, 4);
 347                         header.numReads =  be_int4(*(unsigned int *)(&buffer4));
 348
 349             if (m->debug) { m->mothurOut("[DEBUG]: numReads = " + toString(header.numReads) + "\n"); }
 350
 351                         //read header length
 352                         char buffer5 [2];
 353                         in.read(buffer5, 2);
 354                         header.headerLength =  be_int2(*(unsigned short *)(&buffer5));
 355
 356                         //read key length
 357                         char buffer6 [2];
 358                         in.read(buffer6, 2);
 359                         header.keyLength = be_int2(*(unsigned short *)(&buffer6));
 360
 361                         //read number of flow reads
 362                         char buffer7 [2];
 363                         in.read(buffer7, 2);
 364                         header.numFlowsPerRead =  be_int2(*(unsigned short *)(&buffer7));
 365
 366                         //read format code
 367                         char buffer8 [1];
 368                         in.read(buffer8, 1);
 369                         header.flogramFormatCode = (int)(buffer8[0]);
 370
 371                         //read flow chars
 372                         char* tempBuffer = new char[header.numFlowsPerRead];
 373                         in.read(&(*tempBuffer), header.numFlowsPerRead);
 374                         header.flowChars = tempBuffer;
 375                         if (header.flowChars.length() > header.numFlowsPerRead) { header.flowChars = header.flowChars.substr(0, header.numFlowsPerRead);  }
 376                         delete[] tempBuffer;
 377
 378                         //read key
 379                         char* tempBuffer2 = new char[header.keyLength];
 380                         in.read(&(*tempBuffer2), header.keyLength);
 381                         header.keySequence = tempBuffer2;
 382                         if (header.keySequence.length() > header.keyLength) { header.keySequence = header.keySequence.substr(0, header.keyLength);  }
 383                         delete[] tempBuffer2;
 384
 385                         /* Pad to 8 chars */
 386                         unsigned long long spotInFile = in.tellg();
 387                         unsigned long long spot = (spotInFile + 7)& ~7;  // ~ inverts
 388                         in.seekg(spot);
 389
 390         }else{
 391                         m->mothurOut("Error reading sff common header."); m->mothurOutEndLine();
 392                 }
 393
 394                 return 0;
 395
 396         }
 397         catch(exception& e) {
 398                 m->errorOut(e, "MergeSfffilesCommand", "readCommonHeader");
 399                 exit(1);
 400         }
 401 }
 402 //**********************************************************************************************************************
 403 int MergeSfffilesCommand::adjustCommonHeader(){
 404         try {
 405         //sanity check
 406         bool okay = true;
 407         if (commonHeaders.size() != 0) {
 408             unsigned int magicN = commonHeaders[0].magicNumber;
 409             string version = commonHeaders[0].version;
 410             unsigned short headerLength = commonHeaders[0].headerLength;
 411             unsigned short keyLength = commonHeaders[0].keyLength;
 412             unsigned short numFlows = commonHeaders[0].numFlowsPerRead;
 413             int flowCode = commonHeaders[0].flogramFormatCode;
 414             string flowChars = commonHeaders[0].flowChars;
 415             string keySeq = commonHeaders[0].keySequence;
 416
 417             for (int i = 1; i < commonHeaders.size(); i++) {
 418                 if (commonHeaders[i].magicNumber != magicN)             { okay = false;  m->mothurOut("[ERROR]: merge issue with common headers. Magic numbers do not match. " + filenames[0] + " magic number is " + toString(commonHeaders[0].magicNumber) + ", but " + filenames[i] + " magic number is " + toString(commonHeaders[i].magicNumber) + ".\n");  }
 419                 if (commonHeaders[i].version != version)                { okay = false;   m->mothurOut("[ERROR]: merge issue with common headers. Versions do not match. " + filenames[0] + " version is " + commonHeaders[0].version + ", but " + filenames[i] + " version is " + commonHeaders[i].version + ".\n");     }
 420                 if (commonHeaders[i].headerLength != headerLength)      { okay = false;    m->mothurOut("[ERROR]: merge issue with common headers. Header lengths do not match. " + filenames[0] + " header length is " + toString(commonHeaders[0].headerLength) + ", but " + filenames[i] + " header length is " + toString(commonHeaders[i].headerLength) + ".\n");    }
 421                 if (commonHeaders[i].keyLength != keyLength)            { okay = false;  m->mothurOut("[ERROR]: merge issue with common headers. Key Lengths do not match. " + filenames[0] + " Key length is " + toString(commonHeaders[0].keyLength) + ", but " + filenames[i] + " key length is " + toString(commonHeaders[i].keyLength) + ".\n");    }
 422                 if (commonHeaders[i].numFlowsPerRead != numFlows)       { okay = false;   m->mothurOut("[ERROR]: merge issue with common headers. Number of flows per read do not match. " + filenames[0] + " number of flows is " + toString(commonHeaders[0].numFlowsPerRead) + ", but " + filenames[i] + " number of flows is " + toString(commonHeaders[i].numFlowsPerRead) + ".\n");     }
 423                 if (commonHeaders[i].flogramFormatCode != flowCode)     { okay = false;    m->mothurOut("[ERROR]: merge issue with common headers. Flow format codes do not match. " + filenames[0] + " Flow format code is " + toString(commonHeaders[0].flogramFormatCode) + ", but " + filenames[i] + " flow format code is " + toString(commonHeaders[i].flogramFormatCode) + ".\n");    }
 424                 if (commonHeaders[i].flowChars != flowChars)            { okay = false;   m->mothurOut("[ERROR]: merge issue with common headers. Flow characters do not match. " + filenames[0] + " Flow characters are " + commonHeaders[0].flowChars + ", but " + filenames[i] + " flow characters are " + commonHeaders[i].flowChars + ".\n");    }
 425                 if (commonHeaders[i].keySequence != keySeq)             { okay = false;    m->mothurOut("[ERROR]: merge issue with common headers. Key sequences do not match. " + filenames[0] + " Key sequence is " + commonHeaders[0].keySequence + ", but " + filenames[i] + " key sequence is " + commonHeaders[i].keySequence + ".\n");     }
 426             }
 427         }else { m->control_pressed = true; return 0; } //should never get here
 428
 429         if (!okay) { m->control_pressed = true; return 0; }
 430
 431         string endian = m->findEdianness();
 432         char* mybuffer = new char[4];
 433         ifstream in;
 434         m->openInputFileBinary(currentFileName, in);
 435
 436         //magic number
 437         in.read(mybuffer,4);
 438         ofstream out;
 439         m->openOutputFileBinaryAppend(outputFileHeader, out);
 440         out.write(mybuffer, in.gcount());
 441         delete[] mybuffer;
 442
 443         //version
 444         mybuffer = new char[4];
 445         in.read(mybuffer,4);
 446         out.write(mybuffer, in.gcount());
 447         delete[] mybuffer;
 448
 449         //offset
 450         mybuffer = new char[8];
 451         in.read(mybuffer,8);
 452         unsigned long long offset = 0;
 453         char* thisbuffer = new char[8];
 454         thisbuffer[0] = (offset >> 56) & 0xFF;
 455         thisbuffer[1] = (offset >> 48) & 0xFF;
 456         thisbuffer[2] = (offset >> 40) & 0xFF;
 457         thisbuffer[3] = (offset >> 32) & 0xFF;
 458         thisbuffer[4] = (offset >> 24) & 0xFF;
 459         thisbuffer[5] = (offset >> 16) & 0xFF;
 460         thisbuffer[6] = (offset >> 8) & 0xFF;
 461         thisbuffer[7] = offset & 0xFF;
 462         out.write(thisbuffer, 8);
 463         delete[] thisbuffer;
 464         delete[] mybuffer;
 465
 466         //read index length
 467                 mybuffer = new char[4];
 468         in.read(mybuffer,4);
 469         offset = 0;
 470         char* thisbuffer2 = new char[4];
 471         thisbuffer2[0] = (offset >> 24) & 0xFF;
 472         thisbuffer2[1] = (offset >> 16) & 0xFF;
 473         thisbuffer2[2] = (offset >> 8) & 0xFF;
 474         thisbuffer2[3] = offset & 0xFF;
 475         out.write(thisbuffer2, 4);
 476         delete[] thisbuffer2;
 477         delete[] mybuffer;
 478
 479         //change num reads
 480         mybuffer = new char[4];
 481         in.read(mybuffer,4);
 482         delete[] mybuffer;
 483         thisbuffer2 = new char[4];
 484         if (endian == "BIG_ENDIAN") {
 485             thisbuffer2[0] = (numTotalReads >> 24) & 0xFF;
 486             thisbuffer2[1] = (numTotalReads >> 16) & 0xFF;
 487             thisbuffer2[2] = (numTotalReads >> 8) & 0xFF;
 488             thisbuffer2[3] = numTotalReads & 0xFF;
 489         }else {
 490             thisbuffer2[0] = numTotalReads & 0xFF;
 491             thisbuffer2[1] = (numTotalReads >> 8) & 0xFF;
 492             thisbuffer2[2] = (numTotalReads >> 16) & 0xFF;
 493             thisbuffer2[3] = (numTotalReads >> 24) & 0xFF;
 494         }
 495         out.write(thisbuffer2, 4);
 496         delete[] thisbuffer2;
 497
 498
 499         //read header length
 500         mybuffer = new char[2];
 501         in.read(mybuffer,2);
 502         out.write(mybuffer, in.gcount());
 503         delete[] mybuffer;
 504
 505         //read key length
 506         mybuffer = new char[2];
 507         in.read(mybuffer,2);
 508         out.write(mybuffer, in.gcount());
 509         delete[] mybuffer;
 510
 511         //read number of flow reads
 512         mybuffer = new char[2];
 513         in.read(mybuffer,2);
 514         out.write(mybuffer, in.gcount());
 515         delete[] mybuffer;
 516
 517         //read format code
 518         mybuffer = new char[1];
 519         in.read(mybuffer,1);
 520         out.write(mybuffer, in.gcount());
 521         delete[] mybuffer;
 522
 523         //read flow chars
 524         mybuffer = new char[commonHeaders[0].numFlowsPerRead];
 525         in.read(mybuffer,commonHeaders[0].numFlowsPerRead);
 526         out.write(mybuffer, in.gcount());
 527         delete[] mybuffer;
 528
 529         //read key
 530         mybuffer = new char[commonHeaders[0].keyLength];
 531         in.read(mybuffer,commonHeaders[0].keyLength);
 532         out.write(mybuffer, in.gcount());
 533         delete[] mybuffer;
 534
 535         /* Pad to 8 chars */
 536         unsigned long long spotInFile = in.tellg();
 537         unsigned long long spot = (spotInFile + 7)& ~7;  // ~ inverts
 538         in.seekg(spot);
 539
 540         mybuffer = new char[spot-spotInFile];
 541         out.write(mybuffer, spot-spotInFile);
 542         delete[] mybuffer;
 543         in.close();
 544         out.close();
 545
 546         m->appendBinaryFiles(outputFile, outputFileHeader);
 547         m->renameFile(outputFileHeader, outputFile);
 548         m->mothurRemove(outputFileHeader);
 549
 550                 return 0;
 551
 552         }
 553         catch(exception& e) {
 554                 m->errorOut(e, "MergeSfffilesCommand", "adjustCommonHeader");
 555                 exit(1);
 556         }
 557 }
 558 //**********************************************************************************************************************
 559 bool MergeSfffilesCommand::readSeqData(ifstream& in, seqRead& read, int numFlowReads, Header& header, ofstream& out){
 560         try {
 561         unsigned long long startSpotInFile = in.tellg();
 562                 if (!in.eof()) {
 563
 564             /*****************************************/
 565             //read header
 566
 567             //read header length
 568                         char buffer [2];
 569                         in.read(buffer, 2);
 570                         header.headerLength = be_int2(*(unsigned short *)(&buffer));
 571
 572                         //read name length
 573                         char buffer2 [2];
 574                         in.read(buffer2, 2);
 575                         header.nameLength = be_int2(*(unsigned short *)(&buffer2));
 576
 577                         //read num bases
 578                         char buffer3 [4];
 579                         in.read(buffer3, 4);
 580                         header.numBases =  be_int4(*(unsigned int *)(&buffer3));
 581
 582
 583                         //read clip qual left
 584                         char buffer4 [2];
 585                         in.read(buffer4, 2);
 586                         header.clipQualLeft =  be_int2(*(unsigned short *)(&buffer4));
 587                         header.clipQualLeft = 5;
 588
 589
 590                         //read clip qual right
 591                         char buffer5 [2];
 592                         in.read(buffer5, 2);
 593                         header.clipQualRight =  be_int2(*(unsigned short *)(&buffer5));
 594
 595
 596                         //read clipAdapterLeft
 597                         char buffer6 [2];
 598                         in.read(buffer6, 2);
 599                         header.clipAdapterLeft = be_int2(*(unsigned short *)(&buffer6));
 600
 601
 602                         //read clipAdapterRight
 603                         char buffer7 [2];
 604                         in.read(buffer7, 2);
 605                         header.clipAdapterRight = be_int2(*(unsigned short *)(&buffer7));
 606
 607
 608                         //read name
 609                         char* tempBuffer = new char[header.nameLength];
 610                         in.read(&(*tempBuffer), header.nameLength);
 611                         header.name = tempBuffer;
 612                         if (header.name.length() > header.nameLength) { header.name = header.name.substr(0, header.nameLength);  }
 613
 614                         delete[] tempBuffer;
 615
 616                         /* Pad to 8 chars */
 617                         unsigned long long spotInFile = in.tellg();
 618                         unsigned long long spot = (spotInFile + 7)& ~7;
 619                         in.seekg(spot);
 620
 621             /*****************************************/
 622             //sequence read
 623
 624                         //read flowgram
 625                         read.flowgram.resize(numFlowReads);
 626                         for (int i = 0; i < numFlowReads; i++) {
 627                                 char buffer [2];
 628                                 in.read(buffer, 2);
 629                                 read.flowgram[i] = be_int2(*(unsigned short *)(&buffer));
 630                         }
 631
 632                         //read flowIndex
 633                         read.flowIndex.resize(header.numBases);
 634                         for (int i = 0; i < header.numBases; i++) {
 635                                 char temp[1];
 636                                 in.read(temp, 1);
 637                                 read.flowIndex[i] = be_int1(*(unsigned char *)(&temp));
 638                         }
 639
 640                         //read bases
 641                         char* tempBuffer6 = new char[header.numBases];
 642                         in.read(&(*tempBuffer6), header.numBases);
 643                         read.bases = tempBuffer6;
 644                         if (read.bases.length() > header.numBases) { read.bases = read.bases.substr(0, header.numBases);  }
 645                         delete[] tempBuffer6;
 646
 647                         //read qual scores
 648                         read.qualScores.resize(header.numBases);
 649                         for (int i = 0; i < header.numBases; i++) {
 650                                 char temp[1];
 651                                 in.read(temp, 1);
 652                                 read.qualScores[i] = be_int1(*(unsigned char *)(&temp));
 653                         }
 654
 655                         /* Pad to 8 chars */
 656                         spotInFile = in.tellg();
 657                         spot = (spotInFile + 7)& ~7;
 658                         in.seekg(spot);
 659
 660             char * mybuffer;
 661             mybuffer = new char [spot-startSpotInFile];
 662
 663             ifstream in2;
 664             m->openInputFileBinary(currentFileName, in2);
 665             in2.seekg(startSpotInFile);
 666             in2.read(mybuffer,spot-startSpotInFile);
 667
 668             out.write(mybuffer, in2.gcount());
 669             numTotalReads++;
 670
 671             delete[] mybuffer;
 672             in2.close();
 673
 674                 }else{
 675                         m->mothurOut("Error reading."); m->mothurOutEndLine();
 676                 }
 677
 678         if (in.eof()) {  return true; }
 679
 680                 return false;
 681         }
 682         catch(exception& e) {
 683                 m->errorOut(e, "MergeSfffilesCommand", "readSeqData");
 684                 exit(1);
 685         }
 686 }
 687 //**********************************************************************************************************************
 688 bool MergeSfffilesCommand::sanityCheck(Header& header, seqRead& read) {
 689         try {
 690         bool okay = true;
 691         string message = "[WARNING]: Your sff file may be corrupted! Sequence: " + header.name + "\n";
 692
 693         if (header.clipQualLeft > read.bases.length()) {
 694             okay = false; message += "Clip Qual Left = " + toString(header.clipQualLeft) + ", but we only read " + toString(read.bases.length()) + " bases.\n";
 695         }
 696         if (header.clipQualRight > read.bases.length()) {
 697             okay = false; message += "Clip Qual Right = " + toString(header.clipQualRight) + ", but we only read " + toString(read.bases.length()) + " bases.\n";
 698         }
 699         if (header.clipQualLeft > read.qualScores.size()) {
 700             okay = false; message += "Clip Qual Left = " + toString(header.clipQualLeft) + ", but we only read " + toString(read.qualScores.size()) + " quality scores.\n";
 701         }
 702         if (header.clipQualRight > read.qualScores.size()) {
 703             okay = false; message += "Clip Qual Right = " + toString(header.clipQualRight) + ", but we only read " + toString(read.qualScores.size()) + " quality scores.\n";
 704         }
 705
 706         if (okay == false) {
 707             m->mothurOut(message); m->mothurOutEndLine();
 708         }
 709
 710                 return okay;
 711         }
 712         catch(exception& e) {
 713                 m->errorOut(e, "MergeSfffilesCommand", "sanityCheck");
 714                 exit(1);
 715         }
 716 }
 717 //**********************************************************************************************************************
 718 int MergeSfffilesCommand::readFile(){
 719         try {
 720
 721         string filename;
 722
 723         ifstream in;
 724         m->openInputFile(file, in);
 725
 726         while(!in.eof()) {
 727
 728             if (m->control_pressed) { return 0; }
 729
 730             in >> filename; m->gobble(in);
 731
 732             if (m->debug) { m->mothurOut("[DEBUG]: filename = " + filename + ".\n"); }
 733
 734             //check to make sure both are able to be opened
 735             ifstream in2;
 736             int openForward = m->openInputFile(filename, in2, "noerror");
 737
 738             //if you can't open it, try default location
 739             if (openForward == 1) {
 740                 if (m->getDefaultPath() != "") { //default path is set
 741                     string tryPath = m->getDefaultPath() + m->getSimpleName(filename);
 742                     m->mothurOut("Unable to open " + filename + ". Trying default " + tryPath); m->mothurOutEndLine();
 743                     ifstream in3;
 744                     openForward = m->openInputFile(tryPath, in3, "noerror");
 745                     in3.close();
 746                     filename = tryPath;
 747                 }
 748             }
 749
 750             //if you can't open it, try output location
 751             if (openForward == 1) {
 752                 if (m->getOutputDir() != "") { //default path is set
 753                     string tryPath = m->getOutputDir() + m->getSimpleName(filename);
 754                     m->mothurOut("Unable to open " + filename + ". Trying output directory " + tryPath); m->mothurOutEndLine();
 755                     ifstream in4;
 756                     openForward = m->openInputFile(tryPath, in4, "noerror");
 757                     filename = tryPath;
 758                     in4.close();
 759                 }
 760             }
 761
 762             if (openForward == 1) { //can't find it
 763                 m->mothurOut("[WARNING]: can't find " + filename + ", ignoring.\n");
 764             }else{  filenames.push_back(filename); }
 765
 766         }
 767         in.close();
 768
 769         return 0;
 770     }
 771     catch(exception& e) {
 772         m->errorOut(e, "MergeSfffilesCommand", "readFileNames");
 773         exit(1);
 774     }
 775 }
 776 //**********************************************************************************************************************
 777
 778
 779
 780