trimseqscommand.h

   1 #ifndef TRIMSEQSCOMMAND_H
   2 #define TRIMSEQSCOMMAND_H
   3
   4 /*
   5  *  trimseqscommand.h
   6  *  Mothur
   7  *
   8  *  Created by Pat Schloss on 6/6/09.
   9  *  Copyright 2009 Patrick D. Schloss. All rights reserved.
  10  *
  11  */
  12
  13 #include "mothur.h"
  14 #include "command.hpp"
  15 #include "sequence.hpp"
  16 #include "qualityscores.h"
  17 #include "trimoligos.h"
  18 #include "counttable.h"
  19 #include "oligos.h"
  20
  21
  22 class TrimSeqsCommand : public Command {
  23 public:
  24         TrimSeqsCommand(string);
  25         TrimSeqsCommand();
  26         ~TrimSeqsCommand(){}
  27
  28         vector<string> setParameters();
  29         string getCommandName()                 { return "trim.seqs";   }
  30         string getCommandCategory()             { return "Sequence Processing";         }
  31
  32         string getHelpString();
  33     string getOutputPattern(string);
  34         string getCitation() { return "http://www.mothur.org/wiki/Trim.seqs"; }
  35         string getDescription()         { return "provides the preprocessing features needed to screen and sort pyrosequences"; }
  36
  37         int execute();
  38         void help() { m->mothurOut(getHelpString()); }
  39
  40 private:
  41     struct linePair {
  42         unsigned long long start;
  43         unsigned long long end;
  44         linePair(unsigned long long i, unsigned long long j) : start(i), end(j) {}
  45         linePair() {}
  46     };
  47
  48     Oligos oligos;
  49         bool getOligos(vector<vector<string> >&, vector<vector<string> >&, vector<vector<string> >&, map<string, string>&);
  50         bool keepFirstTrim(Sequence&, QualityScores&);
  51         bool removeLastTrim(Sequence&, QualityScores&);
  52         bool cullLength(Sequence&);
  53         bool cullHomoP(Sequence&);
  54         bool cullAmbigs(Sequence&);
  55         bool abort, createGroup;
  56         string fastaFile, oligoFile, qFileName, groupfile, nameFile, countfile, outputDir;
  57
  58         bool flip, allFiles, qtrim, keepforward, pairedOligos, reorient, logtransform;
  59         int numBarcodes, numFPrimers, numRPrimers, numLinkers, numSpacers, maxAmbig, maxHomoP, minLength, maxLength, processors, tdiffs, bdiffs, pdiffs, ldiffs, sdiffs, comboStarts;
  60         int qWindowSize, qWindowStep, keepFirst, removeLast;
  61         double qRollAverage, qThreshold, qWindowAverage, qAverage;
  62         set<string> filesToRemove;
  63     vector<string> groupVector,outputNames;
  64         map<string, int> groupCounts;
  65         map<string, string> nameMap;
  66     map<string, int> nameCount; //for countfile name -> repCount
  67     map<string, string> groupMap; //for countfile name -> group
  68
  69         vector<int> processIDS;   //processid
  70         vector<linePair> lines;
  71         vector<linePair> qLines;
  72
  73         int driverCreateTrim(string, string, string, string, string, string, string, string, string, string, string, vector<vector<string> >, vector<vector<string> >, vector<vector<string> >, linePair, linePair);
  74         int createProcessesCreateTrim(string, string, string, string, string, string, string, string, string, string, string, vector<vector<string> >, vector<vector<string> >, vector<vector<string> >);
  75         int setLines(string, string);
  76 };
  77
  78 /**************************************************************************************************/
  79 //custom data structure for threads to use.
  80 // This is passed by void pointer so it can be any data type
  81 // that can be passed using a single void pointer (LPVOID).
  82 struct trimData {
  83     unsigned long long start, end;
  84     MothurOut* m;
  85     string filename, qFileName, oligosfile, trimFileName, scrapFileName, trimQFileName, scrapQFileName, trimNFileName, scrapNFileName, trimCFileName, scrapCFileName, groupFileName, nameFile, countfile;
  86         vector<vector<string> > fastaFileNames;
  87     vector<vector<string> > qualFileNames;
  88     vector<vector<string> > nameFileNames;
  89     unsigned long long lineStart, lineEnd, qlineStart, qlineEnd;
  90     bool flip, allFiles, qtrim, keepforward, createGroup, reorient, logtransform;
  91         int maxAmbig, maxHomoP, minLength, maxLength, tdiffs, bdiffs, pdiffs, ldiffs, sdiffs;
  92         int qWindowSize, qWindowStep, keepFirst, removeLast, count;
  93         double qRollAverage, qThreshold, qWindowAverage, qAverage;
  94     map<string, int> nameCount;
  95         map<string, int> groupCounts;
  96         map<string, string> nameMap;
  97     map<string, string> groupMap;
  98
  99         trimData(){}
 100         trimData(string fn, string qn, string nf, string cf, string tn, string sn, string tqn, string sqn, string tnn, string snn, string tcn, string scn,string gn, vector<vector<string> > ffn, vector<vector<string> > qfn, vector<vector<string> > nfn, unsigned long long lstart, unsigned long long lend, unsigned long long qstart, unsigned long long qend,  MothurOut* mout,
 101                       int pd, int bd, int ld, int sd, int td, string o, bool cGroup, bool aFiles, bool keepF, int keepfi, int removeL,
 102                       int WindowStep, int WindowSize, int WindowAverage, bool trim, double Threshold, double Average, double RollAverage, bool lt,
 103                       int minL, int maxA, int maxH, int maxL, bool fli, bool reo, map<string, string> nm, map<string, int> ncount) {
 104         filename = fn;
 105         qFileName = qn;
 106         nameFile = nf;
 107         countfile = cf;
 108         trimFileName = tn;
 109         scrapFileName = sn;
 110         trimQFileName = tqn;
 111         scrapQFileName = sqn;
 112         trimNFileName = tnn;
 113         scrapNFileName = snn;
 114         trimCFileName = tcn;
 115         scrapCFileName = scn;
 116         groupFileName = gn;
 117         fastaFileNames = ffn;
 118         qualFileNames = qfn;
 119         nameFileNames = nfn;
 120         lineStart = lstart;
 121         lineEnd = lend;
 122         qlineStart = qstart;
 123         qlineEnd = qend;
 124                 m = mout;
 125         nameCount = ncount;
 126         oligosfile = o;
 127
 128         pdiffs = pd;
 129         bdiffs = bd;
 130         ldiffs = ld;
 131         sdiffs = sd;
 132         tdiffs = td;
 133
 134         createGroup = cGroup;
 135         allFiles = aFiles;
 136         keepforward = keepF;
 137         keepFirst = keepfi;
 138         removeLast = removeL;
 139         qWindowStep = WindowStep;
 140         qWindowSize = WindowSize;
 141         qWindowAverage = WindowAverage;
 142         qtrim = trim;
 143         qThreshold = Threshold;
 144         qAverage = Average;
 145         qRollAverage = RollAverage;
 146         logtransform = lt;
 147         minLength = minL;
 148         maxAmbig = maxA;
 149         maxHomoP = maxH;
 150         maxLength = maxL;
 151         flip = fli;
 152         reorient = reo;
 153         nameMap = nm;
 154         count = 0;
 155         }
 156 };
 157 /**************************************************************************************************/
 158 #if defined (__APPLE__) || (__MACH__) || (linux) || (__linux) || (__linux__) || (__unix__) || (__unix)
 159 #else
 160 static DWORD WINAPI MyTrimThreadFunction(LPVOID lpParam){
 161         trimData* pDataArray;
 162         pDataArray = (trimData*)lpParam;
 163
 164         try {
 165         ofstream trimFASTAFile;
 166                 pDataArray->m->openOutputFile(pDataArray->trimFileName, trimFASTAFile);
 167
 168                 ofstream scrapFASTAFile;
 169                 pDataArray->m->openOutputFile(pDataArray->scrapFileName, scrapFASTAFile);
 170
 171                 ofstream trimQualFile;
 172                 ofstream scrapQualFile;
 173                 if(pDataArray->qFileName != ""){
 174                         pDataArray->m->openOutputFile(pDataArray->trimQFileName, trimQualFile);
 175                         pDataArray->m->openOutputFile(pDataArray->scrapQFileName, scrapQualFile);
 176                 }
 177
 178                 ofstream trimNameFile;
 179                 ofstream scrapNameFile;
 180                 if(pDataArray->nameFile != ""){
 181                         pDataArray->m->openOutputFile(pDataArray->trimNFileName, trimNameFile);
 182                         pDataArray->m->openOutputFile(pDataArray->scrapNFileName, scrapNameFile);
 183                 }
 184
 185
 186                 ofstream outGroupsFile;
 187                 if ((pDataArray->createGroup) && (pDataArray->countfile == "")){        pDataArray->m->openOutputFile(pDataArray->groupFileName, outGroupsFile);   }
 188                 if(pDataArray->allFiles){
 189                         for (int i = 0; i < pDataArray->fastaFileNames.size(); i++) { //clears old file
 190                                 for (int j = 0; j < pDataArray->fastaFileNames[i].size(); j++) { //clears old file
 191                                         if (pDataArray->fastaFileNames[i][j] != "") {
 192                                                 ofstream temp;
 193                                                 pDataArray->m->openOutputFile(pDataArray->fastaFileNames[i][j], temp);                  temp.close();
 194                                                 if(pDataArray->qFileName != ""){
 195                                                         pDataArray->m->openOutputFile(pDataArray->qualFileNames[i][j], temp);                   temp.close();
 196                                                 }
 197
 198                                                 if(pDataArray->nameFile != ""){
 199                                                         pDataArray->m->openOutputFile(pDataArray->nameFileNames[i][j], temp);                   temp.close();
 200                                                 }
 201                                         }
 202                                 }
 203                         }
 204                 }
 205
 206         ofstream trimCountFile;
 207                 ofstream scrapCountFile;
 208                 if(pDataArray->countfile != ""){
 209                         pDataArray->m->openOutputFile(pDataArray->trimCFileName, trimCountFile);
 210                         pDataArray->m->openOutputFile(pDataArray->scrapCFileName, scrapCountFile);
 211             if ((pDataArray->lineStart == 0) || (pDataArray->lineStart == 1)) { trimCountFile << "Representative_Sequence\ttotal" << endl; scrapCountFile << "Representative_Sequence\ttotal" << endl; }
 212                 }
 213
 214                 ifstream inFASTA;
 215                 pDataArray->m->openInputFile(pDataArray->filename, inFASTA);
 216                 if ((pDataArray->lineStart == 0) || (pDataArray->lineStart == 1)) {
 217                         inFASTA.seekg(0);
 218                 }else { //this accounts for the difference in line endings.
 219                         inFASTA.seekg(pDataArray->lineStart-1); pDataArray->m->gobble(inFASTA);
 220                 }
 221
 222                 ifstream qFile;
 223                 if(pDataArray->qFileName != "") {
 224                         pDataArray->m->openInputFile(pDataArray->qFileName, qFile);
 225                         if ((pDataArray->qlineStart == 0) || (pDataArray->qlineStart == 1)) {
 226                 qFile.seekg(0);
 227             }else { //this accounts for the difference in line endings.
 228                 qFile.seekg(pDataArray->qlineStart-1); pDataArray->m->gobble(qFile);
 229             }
 230                 }
 231
 232         int numBarcodes, numLinkers, numSpacers, numRPrimers, numFPrimers;
 233         bool pairedOligos;
 234         Oligos oligos(pDataArray->oligosfile);
 235         if (oligos.hasPairedBarcodes()) {
 236             pairedOligos = true;
 237             numFPrimers = oligos.getPairedPrimers().size();
 238             numBarcodes = oligos.getPairedBarcodes().size();
 239         }else {
 240             pairedOligos = false;
 241             numFPrimers = oligos.getPrimers().size();
 242             numBarcodes = oligos.getBarcodes().size();
 243         }
 244
 245         numLinkers = oligos.getLinkers().size();
 246         numSpacers = oligos.getSpacers().size();
 247         numRPrimers = oligos.getReversePrimers().size();
 248
 249         TrimOligos* trimOligos = NULL;
 250         if (pairedOligos)   {   trimOligos = new TrimOligos(pDataArray->pdiffs, pDataArray->bdiffs, 0, 0, oligos.getPairedPrimers(), oligos.getPairedBarcodes());    }
 251         else                {   trimOligos = new TrimOligos(pDataArray->pdiffs, pDataArray->bdiffs, pDataArray->ldiffs, pDataArray->sdiffs, oligos.getPrimers(), oligos.getBarcodes(), oligos.getReversePrimers(), oligos.getLinkers(), oligos.getSpacers());  }
 252
 253         TrimOligos* rtrimOligos = NULL;
 254         if (pDataArray->reorient) {
 255             rtrimOligos = new TrimOligos(pDataArray->pdiffs, pDataArray->bdiffs, 0, 0,  oligos.getReorientedPairedPrimers(), oligos.getReorientedPairedBarcodes()); numBarcodes = oligos.getReorientedPairedBarcodes().size();
 256         }
 257
 258                 pDataArray->count = 0;
 259                 for(int i = 0; i < pDataArray->lineEnd; i++){ //end is the number of sequences to process
 260
 261                         if (pDataArray->m->control_pressed) {
 262                 delete trimOligos; if (pDataArray->reorient) { delete rtrimOligos; }
 263                                 inFASTA.close(); trimFASTAFile.close(); scrapFASTAFile.close();
 264                                 if ((pDataArray->createGroup) && (pDataArray->countfile == "")) {        outGroupsFile.close();   }
 265                 if(pDataArray->qFileName != "") {       qFile.close();  scrapQualFile.close(); trimQualFile.close();    }
 266                 if(pDataArray->nameFile != "")  {       scrapNameFile.close(); trimNameFile.close();    }
 267                 if(pDataArray->countfile != "") {       scrapCountFile.close(); trimCountFile.close();  }
 268
 269                                 if(pDataArray->qFileName != ""){ qFile.close(); }
 270                                 return 0;
 271                         }
 272
 273                         int success = 1;
 274                         string trashCode = "";
 275                         int currentSeqsDiffs = 0;
 276
 277                         Sequence currSeq(inFASTA); pDataArray->m->gobble(inFASTA);
 278             Sequence savedSeq(currSeq.getName(), currSeq.getAligned());
 279
 280                         QualityScores currQual; QualityScores savedQual;
 281                         if(pDataArray->qFileName != ""){
 282                                 currQual = QualityScores(qFile);  pDataArray->m->gobble(qFile);
 283                 savedQual.setName(currQual.getName()); savedQual.setScores(currQual.getScores());
 284                         }
 285
 286
 287                         string origSeq = currSeq.getUnaligned();
 288                         if (origSeq != "") {
 289                 pDataArray->count++;
 290
 291                                 int barcodeIndex = 0;
 292                                 int primerIndex = 0;
 293
 294                 if(numLinkers != 0){
 295                                         success = trimOligos->stripLinker(currSeq, currQual);
 296                                         if(success > pDataArray->ldiffs)                {       trashCode += 'k';       }
 297                                         else{ currentSeqsDiffs += success;  }
 298                                 }
 299
 300                                 if(numBarcodes != 0){
 301                                         success = trimOligos->stripBarcode(currSeq, currQual, barcodeIndex);
 302                                         if(success > pDataArray->bdiffs)                {       trashCode += 'b';       }
 303                                         else{ currentSeqsDiffs += success;  }
 304                                 }
 305
 306                 if(numSpacers != 0){
 307                                         success = trimOligos->stripSpacer(currSeq, currQual);
 308                                         if(success > pDataArray->sdiffs)                {       trashCode += 's';       }
 309                                         else{ currentSeqsDiffs += success;  }
 310
 311                                 }
 312
 313                                 if(numFPrimers != 0){
 314                                         success = trimOligos->stripForward(currSeq, currQual, primerIndex, pDataArray->keepforward);
 315                                         if(success > pDataArray->pdiffs)                {       trashCode += 'f';       }
 316                                         else{ currentSeqsDiffs += success;  }
 317                                 }
 318
 319                                 if (currentSeqsDiffs > pDataArray->tdiffs)      {       trashCode += 't';   }
 320
 321                                 if(numRPrimers != 0){
 322                                         success = trimOligos->stripReverse(currSeq, currQual);
 323                                         if(!success)                            {       trashCode += 'r';       }
 324                                 }
 325
 326                 if (pDataArray->reorient && (trashCode != "")) { //if you failed and want to check the reverse
 327                     int thisSuccess = 0;
 328                     string thisTrashCode = "";
 329                     int thisCurrentSeqsDiffs = 0;
 330
 331                     int thisBarcodeIndex = 0;
 332                     int thisPrimerIndex = 0;
 333
 334                     if(numBarcodes != 0){
 335                         thisSuccess = rtrimOligos->stripBarcode(savedSeq, savedQual, thisBarcodeIndex);
 336                         if(thisSuccess > pDataArray->bdiffs)            {       thisTrashCode += 'b';   }
 337                         else{ thisCurrentSeqsDiffs += thisSuccess;  }
 338                     }
 339
 340                     if(numFPrimers != 0){
 341                         thisSuccess = rtrimOligos->stripForward(savedSeq, savedQual, thisPrimerIndex, pDataArray->keepforward);
 342                         if(thisSuccess > pDataArray->pdiffs)            {       thisTrashCode += 'f';   }
 343                         else{ thisCurrentSeqsDiffs += thisSuccess;  }
 344                     }
 345
 346                     if (thisCurrentSeqsDiffs > pDataArray->tdiffs)      {       thisTrashCode += 't';   }
 347
 348                     if (thisTrashCode == "") {
 349                         trashCode = thisTrashCode;
 350                         success = thisSuccess;
 351                         currentSeqsDiffs = thisCurrentSeqsDiffs;
 352                         barcodeIndex = thisBarcodeIndex;
 353                         primerIndex = thisPrimerIndex;
 354                         savedSeq.reverseComplement();
 355                         currSeq.setAligned(savedSeq.getAligned());
 356                         if(pDataArray->qFileName != ""){
 357                             savedQual.flipQScores();
 358                             currQual.setScores(savedQual.getScores());
 359                         }
 360                     }else { trashCode += "(" + thisTrashCode + ")";  }
 361                 }
 362
 363
 364                                 if(pDataArray->keepFirst != 0){
 365                                         //success = keepFirstTrim(currSeq, currQual);
 366                     success = 1;
 367                     if(currQual.getName() != ""){
 368                         currQual.trimQScores(-1, pDataArray->keepFirst);
 369                     }
 370                     currSeq.trim(pDataArray->keepFirst);
 371                                 }
 372
 373                                 if(pDataArray->removeLast != 0){
 374                                         //success = removeLastTrim(currSeq, currQual);
 375                     success = 0;
 376                     int length = currSeq.getNumBases() - pDataArray->removeLast;
 377
 378                     if(length > 0){
 379                         if(currQual.getName() != ""){
 380                             currQual.trimQScores(-1, length);
 381                         }
 382                         currSeq.trim(length);
 383                         success = 1;
 384                     }
 385                     else{ success = 0; }
 386
 387                                         if(!success)                            {       trashCode += 'l';       }
 388                                 }
 389
 390
 391                                 if(pDataArray->qFileName != ""){
 392                                         int origLength = currSeq.getNumBases();
 393
 394                                         if(pDataArray->qThreshold != 0)                 {       success = currQual.stripQualThreshold(currSeq, pDataArray->qThreshold);                 }
 395                                         else if(pDataArray->qAverage != 0)              {       success = currQual.cullQualAverage(currSeq, pDataArray->qAverage, pDataArray->logtransform);                            }
 396                                         else if(pDataArray->qRollAverage != 0)  {       success = currQual.stripQualRollingAverage(currSeq, pDataArray->qRollAverage, pDataArray->logtransform);        }
 397                                         else if(pDataArray->qWindowAverage != 0){       success = currQual.stripQualWindowAverage(currSeq, pDataArray->qWindowStep, pDataArray->qWindowSize, pDataArray->qWindowAverage, pDataArray->logtransform);     }
 398                                         else                                            {       success = 1;                            }
 399
 400                                         //you don't want to trim, if it fails above then scrap it
 401                                         if ((!pDataArray->qtrim) && (origLength != currSeq.getNumBases())) {  success = 0; }
 402
 403                                         if(!success)                            {       trashCode += 'q';       }
 404                                 }
 405
 406                                 if(pDataArray->minLength > 0 || pDataArray->maxLength > 0){
 407                                         //success = cullLength(currSeq);
 408                     int length = currSeq.getNumBases();
 409                     success = 0;        //guilty until proven innocent
 410                     if(length >= pDataArray->minLength && pDataArray->maxLength == 0)                   {       success = 1;    }
 411                     else if(length >= pDataArray->minLength && length <= pDataArray->maxLength) {       success = 1;    }
 412                     else                                                                                                {       success = 0;    }
 413
 414                                         if(!success)                            {       trashCode += 'l';       }
 415                                 }
 416                                 if(pDataArray->maxHomoP > 0){
 417                                         //success = cullHomoP(currSeq);
 418                     int longHomoP = currSeq.getLongHomoPolymer();
 419                     success = 0;        //guilty until proven innocent
 420                     if(longHomoP <= pDataArray->maxHomoP){      success = 1;    }
 421                     else                                        {       success = 0;    }
 422
 423                                         if(!success)                            {       trashCode += 'h';       }
 424                                 }
 425                                 if(pDataArray->maxAmbig != -1){
 426                                         //success = cullAmbigs(currSeq);
 427                     int numNs = currSeq.getAmbigBases();
 428                     success = 0;        //guilty until proven innocent
 429                     if(numNs <= pDataArray->maxAmbig)   {       success = 1;    }
 430                     else                                        {       success = 0;    }
 431                                         if(!success)                            {       trashCode += 'n';       }
 432                                 }
 433
 434                                 if(pDataArray->flip){           // should go last
 435                                         currSeq.reverseComplement();
 436                                         if(pDataArray->qFileName != ""){
 437                                                 currQual.flipQScores();
 438                                         }
 439                                 }
 440
 441                                 if(trashCode.length() == 0){
 442                     string thisGroup = "";
 443                     if (pDataArray->createGroup) {    thisGroup = oligos.getGroupName(barcodeIndex, primerIndex);                     }
 444
 445                     int pos = thisGroup.find("ignore");
 446                     if (pos == string::npos) {
 447
 448                         currSeq.setAligned(currSeq.getUnaligned());
 449                         currSeq.printSequence(trimFASTAFile);
 450
 451                         if(pDataArray->qFileName != ""){
 452                             currQual.printQScores(trimQualFile);
 453                         }
 454
 455                         if(pDataArray->nameFile != ""){
 456                             map<string, string>::iterator itName = pDataArray->nameMap.find(currSeq.getName());
 457                             if (itName != pDataArray->nameMap.end()) {  trimNameFile << itName->first << '\t' << itName->second << endl; }
 458                             else { pDataArray->m->mothurOut("[ERROR]: " + currSeq.getName() + " is not in your namefile, please correct."); pDataArray->m->mothurOutEndLine(); }
 459                         }
 460
 461                         int numRedundants = 0;
 462                         if (pDataArray->countfile != "") {
 463                             map<string, int>::iterator itCount = pDataArray->nameCount.find(currSeq.getName());
 464                             if (itCount != pDataArray->nameCount.end()) {
 465                                 trimCountFile << itCount->first << '\t' << itCount->second << endl;
 466                                 numRedundants = itCount->second-1;
 467                             }else { pDataArray->m->mothurOut("[ERROR]: " + currSeq.getName() + " is not in your count file, please correct."); pDataArray->m->mothurOutEndLine(); }
 468                         }
 469
 470                         if (pDataArray->createGroup) {
 471                             if(numBarcodes != 0){
 472
 473                                 if (pDataArray->countfile == "") { outGroupsFile << currSeq.getName() << '\t' << thisGroup << endl; }
 474                                 else {   pDataArray->groupMap[currSeq.getName()] = thisGroup; }
 475
 476                                 if (pDataArray->nameFile != "") {
 477                                     map<string, string>::iterator itName = pDataArray->nameMap.find(currSeq.getName());
 478                                     if (itName != pDataArray->nameMap.end()) {
 479                                         vector<string> thisSeqsNames;
 480                                         pDataArray->m->splitAtChar(itName->second, thisSeqsNames, ',');
 481                                         numRedundants = thisSeqsNames.size()-1; //we already include ourselves below
 482                                         for (int k = 1; k < thisSeqsNames.size(); k++) { //start at 1 to skip self
 483                                             outGroupsFile << thisSeqsNames[k] << '\t' << thisGroup << endl;
 484                                         }
 485                                     }else { pDataArray->m->mothurOut("[ERROR]: " + currSeq.getName() + " is not in your namefile, please correct."); pDataArray->m->mothurOutEndLine(); }
 486                                 }
 487
 488                                 map<string, int>::iterator it = pDataArray->groupCounts.find(thisGroup);
 489                                 if (it == pDataArray->groupCounts.end()) {      pDataArray->groupCounts[thisGroup] = 1 + numRedundants; }
 490                                 else { pDataArray->groupCounts[it->first] += (1 + numRedundants); }
 491
 492                             }
 493                         }
 494
 495                         if(pDataArray->allFiles){
 496                             ofstream output;
 497                             pDataArray->m->openOutputFileAppend(pDataArray->fastaFileNames[barcodeIndex][primerIndex], output);
 498                             currSeq.printSequence(output);
 499                             output.close();
 500
 501                             if(pDataArray->qFileName != ""){
 502                                 pDataArray->m->openOutputFileAppend(pDataArray->qualFileNames[barcodeIndex][primerIndex], output);
 503                                 currQual.printQScores(output);
 504                                 output.close();
 505                             }
 506
 507                             if(pDataArray->nameFile != ""){
 508                                 map<string, string>::iterator itName = pDataArray->nameMap.find(currSeq.getName());
 509                                 if (itName != pDataArray->nameMap.end()) {
 510                                     pDataArray->m->openOutputFileAppend(pDataArray->nameFileNames[barcodeIndex][primerIndex], output);
 511                                     output << itName->first << '\t' << itName->second << endl;
 512                                     output.close();
 513                                 }else { pDataArray->m->mothurOut("[ERROR]: " + currSeq.getName() + " is not in your namefile, please correct."); pDataArray->m->mothurOutEndLine(); }
 514                             }
 515                         }
 516                     }
 517                                 }
 518                                 else{
 519                                         if(pDataArray->nameFile != ""){ //needs to be before the currSeq name is changed
 520                                                 map<string, string>::iterator itName = pDataArray->nameMap.find(currSeq.getName());
 521                                                 if (itName != pDataArray->nameMap.end()) {  scrapNameFile << itName->first << '\t' << itName->second << endl; }
 522                                                 else { pDataArray->m->mothurOut("[ERROR]: " + currSeq.getName() + " is not in your namefile, please correct."); pDataArray->m->mothurOutEndLine(); }
 523                                         }
 524                     if (pDataArray->countfile != "") {
 525                         map<string, int>::iterator itCount = pDataArray->nameCount.find(currSeq.getName());
 526                         if (itCount != pDataArray->nameCount.end()) {
 527                             trimCountFile << itCount->first << '\t' << itCount->second << endl;
 528                         }else { pDataArray->m->mothurOut("[ERROR]: " + currSeq.getName() + " is not in your count file, please correct."); pDataArray->m->mothurOutEndLine(); }
 529                     }
 530                                         currSeq.setName(currSeq.getName() + '|' + trashCode);
 531                                         currSeq.setUnaligned(origSeq);
 532                                         currSeq.setAligned(origSeq);
 533                                         currSeq.printSequence(scrapFASTAFile);
 534                                         if(pDataArray->qFileName != ""){
 535                                                 currQual.printQScores(scrapQualFile);
 536                                         }
 537                                 }
 538
 539                         }
 540
 541                         //report progress
 542                         if((pDataArray->count) % 1000 == 0){    pDataArray->m->mothurOut(toString(pDataArray->count)); pDataArray->m->mothurOutEndLine();               }
 543
 544                 }
 545                 //report progress
 546                 if((pDataArray->count) % 1000 != 0){    pDataArray->m->mothurOut(toString(pDataArray->count)); pDataArray->m->mothurOutEndLine();               }
 547
 548         if (pDataArray->reorient) { delete rtrimOligos; }
 549                 delete trimOligos;
 550                 inFASTA.close();
 551                 trimFASTAFile.close();
 552                 scrapFASTAFile.close();
 553                 if (pDataArray->createGroup) {   outGroupsFile.close();   }
 554                 if(pDataArray->qFileName != "") {       qFile.close();  scrapQualFile.close(); trimQualFile.close();    }
 555                 if(pDataArray->nameFile != "")  {       scrapNameFile.close(); trimNameFile.close();    }
 556
 557         return 0;
 558
 559         }
 560         catch(exception& e) {
 561             pDataArray->m->errorOut(e, "TrimSeqsCommand", "MyTrimThreadFunction");
 562             exit(1);
 563         }
 564     }
 565 #endif
 566
 567
 568 /**************************************************************************************************/
 569
 570 #endif