]> git.donarmstrong.com Git - mothur.git/blob - removeseqscommand.cpp
modified trim.seqs to speed up the append with multiple processors
[mothur.git] / removeseqscommand.cpp
1 /*
2  *  removeseqscommand.cpp
3  *  Mothur
4  *
5  *  Created by Sarah Westcott on 7/8/09.
6  *  Copyright 2009 Schloss Lab UMASS Amherst. All rights reserved.
7  *
8  */
9
10 #include "removeseqscommand.h"
11 #include "sequence.hpp"
12 #include "listvector.hpp"
13
14 //**********************************************************************************************************************
15 vector<string> RemoveSeqsCommand::getValidParameters(){ 
16         try {
17                 string Array[] =  {"fasta","name", "group", "alignreport", "accnos", "qfile","list","taxonomy","outputdir","inputdir", "dups" };
18                 vector<string> myArray (Array, Array+(sizeof(Array)/sizeof(string)));
19                 return myArray;
20         }
21         catch(exception& e) {
22                 m->errorOut(e, "RemoveSeqsCommand", "getValidParameters");
23                 exit(1);
24         }
25 }
26 //**********************************************************************************************************************
27 RemoveSeqsCommand::RemoveSeqsCommand(){ 
28         try {
29                 abort = true;
30                 //initialize outputTypes
31                 vector<string> tempOutNames;
32                 outputTypes["fasta"] = tempOutNames;
33                 outputTypes["taxonomy"] = tempOutNames;
34                 outputTypes["name"] = tempOutNames;
35                 outputTypes["group"] = tempOutNames;
36                 outputTypes["alignreport"] = tempOutNames;
37                 outputTypes["list"] = tempOutNames;
38                 outputTypes["qfile"] = tempOutNames;
39         }
40         catch(exception& e) {
41                 m->errorOut(e, "RemoveSeqsCommand", "RemoveSeqsCommand");
42                 exit(1);
43         }
44 }
45 //**********************************************************************************************************************
46 vector<string> RemoveSeqsCommand::getRequiredParameters(){      
47         try {
48                 string Array[] =  {"accnos"};
49                 vector<string> myArray (Array, Array+(sizeof(Array)/sizeof(string)));
50                 return myArray;
51         }
52         catch(exception& e) {
53                 m->errorOut(e, "RemoveSeqsCommand", "getRequiredParameters");
54                 exit(1);
55         }
56 }
57 //**********************************************************************************************************************
58 vector<string> RemoveSeqsCommand::getRequiredFiles(){   
59         try {
60                 vector<string> myArray;
61                 return myArray;
62         }
63         catch(exception& e) {
64                 m->errorOut(e, "RemoveSeqsCommand", "getRequiredFiles");
65                 exit(1);
66         }
67 }
68 //**********************************************************************************************************************
69 RemoveSeqsCommand::RemoveSeqsCommand(string option)  {
70         try {
71                 abort = false;
72                 
73                 //allow user to run help
74                 if(option == "help") { help(); abort = true; }
75                 
76                 else {
77                         //valid paramters for this command
78                         string Array[] =  {"fasta","name", "group", "alignreport", "accnos", "qfile", "list","taxonomy","outputdir","inputdir", "dups" };
79                         vector<string> myArray (Array, Array+(sizeof(Array)/sizeof(string)));
80                         
81                         OptionParser parser(option);
82                         map<string,string> parameters = parser.getParameters();
83                         
84                         ValidParameters validParameter;
85                         map<string,string>::iterator it;
86                         
87                         //check to make sure all parameters are valid for command
88                         for (it = parameters.begin(); it != parameters.end(); it++) { 
89                                 if (validParameter.isValidParameter(it->first, myArray, it->second) != true) {  abort = true;  }
90                         }
91                         
92                         //initialize outputTypes
93                         vector<string> tempOutNames;
94                         outputTypes["fasta"] = tempOutNames;
95                         outputTypes["taxonomy"] = tempOutNames;
96                         outputTypes["name"] = tempOutNames;
97                         outputTypes["group"] = tempOutNames;
98                         outputTypes["alignreport"] = tempOutNames;
99                         outputTypes["list"] = tempOutNames;
100                         outputTypes["qfile"] = tempOutNames;
101                         
102                         //if the user changes the output directory command factory will send this info to us in the output parameter 
103                         outputDir = validParameter.validFile(parameters, "outputdir", false);           if (outputDir == "not found"){  outputDir = "";         }
104                         
105                         //if the user changes the input directory command factory will send this info to us in the output parameter 
106                         string inputDir = validParameter.validFile(parameters, "inputdir", false);              
107                         if (inputDir == "not found"){   inputDir = "";          }
108                         else {
109                                 string path;
110                                 it = parameters.find("alignreport");
111                                 //user has given a template file
112                                 if(it != parameters.end()){ 
113                                         path = m->hasPath(it->second);
114                                         //if the user has not given a path then, add inputdir. else leave path alone.
115                                         if (path == "") {       parameters["alignreport"] = inputDir + it->second;              }
116                                 }
117                                 
118                                 it = parameters.find("fasta");
119                                 //user has given a template file
120                                 if(it != parameters.end()){ 
121                                         path = m->hasPath(it->second);
122                                         //if the user has not given a path then, add inputdir. else leave path alone.
123                                         if (path == "") {       parameters["fasta"] = inputDir + it->second;            }
124                                 }
125                                 
126                                 it = parameters.find("accnos");
127                                 //user has given a template file
128                                 if(it != parameters.end()){ 
129                                         path = m->hasPath(it->second);
130                                         //if the user has not given a path then, add inputdir. else leave path alone.
131                                         if (path == "") {       parameters["accnos"] = inputDir + it->second;           }
132                                 }
133                                 
134                                 it = parameters.find("list");
135                                 //user has given a template file
136                                 if(it != parameters.end()){ 
137                                         path = m->hasPath(it->second);
138                                         //if the user has not given a path then, add inputdir. else leave path alone.
139                                         if (path == "") {       parameters["list"] = inputDir + it->second;             }
140                                 }
141                                 
142                                 it = parameters.find("name");
143                                 //user has given a template file
144                                 if(it != parameters.end()){ 
145                                         path = m->hasPath(it->second);
146                                         //if the user has not given a path then, add inputdir. else leave path alone.
147                                         if (path == "") {       parameters["name"] = inputDir + it->second;             }
148                                 }
149                                 
150                                 it = parameters.find("group");
151                                 //user has given a template file
152                                 if(it != parameters.end()){ 
153                                         path = m->hasPath(it->second);
154                                         //if the user has not given a path then, add inputdir. else leave path alone.
155                                         if (path == "") {       parameters["group"] = inputDir + it->second;            }
156                                 }
157                                 
158                                 it = parameters.find("taxonomy");
159                                 //user has given a template file
160                                 if(it != parameters.end()){ 
161                                         path = m->hasPath(it->second);
162                                         //if the user has not given a path then, add inputdir. else leave path alone.
163                                         if (path == "") {       parameters["taxonomy"] = inputDir + it->second;         }
164                                 }
165                                 
166                                 it = parameters.find("qfile");
167                                 //user has given a template file
168                                 if(it != parameters.end()){ 
169                                         path = m->hasPath(it->second);
170                                         //if the user has not given a path then, add inputdir. else leave path alone.
171                                         if (path == "") {       parameters["qfile"] = inputDir + it->second;            }
172                                 }
173                         }
174
175                         
176                         //check for required parameters
177                         accnosfile = validParameter.validFile(parameters, "accnos", true);
178                         if (accnosfile == "not open") { abort = true; }
179                         else if (accnosfile == "not found") {  accnosfile = "";  m->mothurOut("You must provide an accnos file."); m->mothurOutEndLine(); abort = true; }       
180                         
181                         fastafile = validParameter.validFile(parameters, "fasta", true);
182                         if (fastafile == "not open") { abort = true; }
183                         else if (fastafile == "not found") {  fastafile = "";  }        
184                         
185                         namefile = validParameter.validFile(parameters, "name", true);
186                         if (namefile == "not open") { abort = true; }
187                         else if (namefile == "not found") {  namefile = "";  }  
188                         
189                         groupfile = validParameter.validFile(parameters, "group", true);
190                         if (groupfile == "not open") { abort = true; }
191                         else if (groupfile == "not found") {  groupfile = "";  }        
192                         
193                         alignfile = validParameter.validFile(parameters, "alignreport", true);
194                         if (alignfile == "not open") { abort = true; }
195                         else if (alignfile == "not found") {  alignfile = "";  }
196                         
197                         listfile = validParameter.validFile(parameters, "list", true);
198                         if (listfile == "not open") { abort = true; }
199                         else if (listfile == "not found") {  listfile = "";  }
200                         
201                         taxfile = validParameter.validFile(parameters, "taxonomy", true);
202                         if (taxfile == "not open") { abort = true; }
203                         else if (taxfile == "not found") {  taxfile = "";  }
204                         
205                         qualfile = validParameter.validFile(parameters, "qfile", true);
206                         if (qualfile == "not open") { abort = true; }
207                         else if (qualfile == "not found") {  qualfile = "";  }                  
208
209                         
210                         string usedDups = "true";
211                         string temp = validParameter.validFile(parameters, "dups", false);      
212                         if (temp == "not found") { 
213                                 if (namefile != "") {  temp = "true";                                   }
214                                 else                            {  temp = "false"; usedDups = "";       }
215                         }
216                         dups = m->isTrue(temp);
217                         
218                         if ((fastafile == "") && (namefile == "") && (groupfile == "") && (alignfile == "") && (listfile == "") && (taxfile == "") && (qualfile == ""))  { m->mothurOut("You must provide at least one of the following: fasta, name, group, taxonomy, quality, alignreport or list."); m->mothurOutEndLine(); abort = true; }
219                         
220                         if ((usedDups != "") && (namefile == "")) {  m->mothurOut("You may only use dups with the name option."); m->mothurOutEndLine();  abort = true; }                       
221                 }
222
223         }
224         catch(exception& e) {
225                 m->errorOut(e, "RemoveSeqsCommand", "RemoveSeqsCommand");
226                 exit(1);
227         }
228 }
229 //**********************************************************************************************************************
230
231 void RemoveSeqsCommand::help(){
232         try {
233                 m->mothurOut("The remove.seqs command reads an .accnos file and at least one of the following file types: fasta, name, group, list, taxonomy, quality or alignreport file.\n");
234                 m->mothurOut("It outputs a file containing the sequences NOT in the .accnos file.\n");
235                 m->mothurOut("The remove.seqs command parameters are accnos, fasta, name, group, list, taxonomy, qfile, alignreport and dups.  You must provide accnos and at least one of the file parameters.\n");
236                 m->mothurOut("The dups parameter allows you to remove the entire line from a name file if you remove any name from the line. default=true. \n");
237                 m->mothurOut("The remove.seqs command should be in the following format: remove.seqs(accnos=yourAccnos, fasta=yourFasta).\n");
238                 m->mothurOut("Example remove.seqs(accnos=amazon.accnos, fasta=amazon.fasta).\n");
239                 m->mothurOut("Note: No spaces between parameter labels (i.e. fasta), '=' and parameters (i.e.yourFasta).\n\n");
240         }
241         catch(exception& e) {
242                 m->errorOut(e, "RemoveSeqsCommand", "help");
243                 exit(1);
244         }
245 }
246
247 //**********************************************************************************************************************
248
249 int RemoveSeqsCommand::execute(){
250         try {
251                 
252                 if (abort == true) { return 0; }
253                 
254                 //get names you want to keep
255                 readAccnos();
256                 
257                 if (m->control_pressed) { return 0; }
258                 
259                 //read through the correct file and output lines you want to keep
260                 if (namefile != "")                     {               readName();             }
261                 if (fastafile != "")            {               readFasta();    }
262                 if (groupfile != "")            {               readGroup();    }
263                 if (alignfile != "")            {               readAlign();    }
264                 if (listfile != "")                     {               readList();             }
265                 if (taxfile != "")                      {               readTax();              }
266                 if (qualfile != "")                     {               readQual();             }
267                 
268                 if (m->control_pressed) { for (int i = 0; i < outputNames.size(); i++) {        remove(outputNames[i].c_str()); } return 0; }
269                 
270                 m->mothurOut("Removed " + toString(names.size()) + " sequences."); m->mothurOutEndLine();
271                 
272                 if (outputNames.size() != 0) {
273                         m->mothurOutEndLine();
274                         m->mothurOut("Output File Names: "); m->mothurOutEndLine();
275                         for (int i = 0; i < outputNames.size(); i++) {  m->mothurOut(outputNames[i]); m->mothurOutEndLine();    }
276                         m->mothurOutEndLine();
277                 }
278                 
279                 return 0;               
280         }
281
282         catch(exception& e) {
283                 m->errorOut(e, "RemoveSeqsCommand", "execute");
284                 exit(1);
285         }
286 }
287
288 //**********************************************************************************************************************
289 int RemoveSeqsCommand::readFasta(){
290         try {
291                 string thisOutputDir = outputDir;
292                 if (outputDir == "") {  thisOutputDir += m->hasPath(fastafile);  }
293                 string outputFileName = thisOutputDir + m->getRootName(m->getSimpleName(fastafile)) + "pick" + m->getExtension(fastafile);
294                 
295                 ofstream out;
296                 m->openOutputFile(outputFileName, out);
297                 
298                 ifstream in;
299                 m->openInputFile(fastafile, in);
300                 string name;
301                 
302                 bool wroteSomething = false;
303                 
304                 while(!in.eof()){
305                         if (m->control_pressed) { in.close();  out.close();  remove(outputFileName.c_str());  return 0; }
306                         
307                         Sequence currSeq(in);
308                         name = currSeq.getName();
309                         
310                         if (name != "") {
311                                 //if this name is in the accnos file
312                                 if (names.count(name) == 0) {
313                                         wroteSomething = true;
314                                         
315                                         currSeq.printSequence(out);
316                                 }
317                         }
318                         m->gobble(in);
319                 }
320                 in.close();     
321                 out.close();
322                 
323                 if (wroteSomething == false) {  m->mothurOut("Your file contains only sequences from the .accnos file."); m->mothurOutEndLine();  }
324                 outputTypes["fasta"].push_back(outputFileName);  outputNames.push_back(outputFileName);
325                 
326                 return 0;
327                 
328         }
329         catch(exception& e) {
330                 m->errorOut(e, "RemoveSeqsCommand", "readFasta");
331                 exit(1);
332         }
333 }
334 //**********************************************************************************************************************
335 int RemoveSeqsCommand::readQual(){
336         try {
337                 string thisOutputDir = outputDir;
338                 if (outputDir == "") {  thisOutputDir += m->hasPath(qualfile);  }
339                 string outputFileName = thisOutputDir + m->getRootName(m->getSimpleName(qualfile)) + "pick" +  m->getExtension(qualfile);
340                 ofstream out;
341                 m->openOutputFile(outputFileName, out);
342                 
343                 
344                 ifstream in;
345                 m->openInputFile(qualfile, in);
346                 string name;
347                 
348                 bool wroteSomething = false;
349                 
350                 
351                 while(!in.eof()){       
352                         string saveName = "";
353                         string name = "";
354                         string scores = "";
355                         
356                         in >> name; 
357                         
358                         if (name.length() != 0) { 
359                                 saveName = name.substr(1);
360                                 while (!in.eof())       {       
361                                         char c = in.get(); 
362                                         if (c == 10 || c == 13){        break;  }
363                                         else { name += c; }     
364                                 } 
365                                 m->gobble(in);
366                         }
367                         
368                         while(in){
369                                 char letter= in.get();
370                                 if(letter == '>'){      in.putback(letter);     break;  }
371                                 else{ scores += letter; }
372                         }
373                         
374                         m->gobble(in);
375                         
376                         if (names.count(saveName) == 0) {
377                                 wroteSomething = true;
378                                 
379                                 out << name << endl << scores;
380                         }
381                         
382                         m->gobble(in);
383                 }
384                 in.close();
385                 out.close();
386                 
387                 
388                 if (wroteSomething == false) { m->mothurOut("Your file contains only sequences from the .accnos file."); m->mothurOutEndLine();  }
389                 outputNames.push_back(outputFileName);  outputTypes["qfile"].push_back(outputFileName); 
390                 
391                 return 0;
392                 
393         }
394         catch(exception& e) {
395                 m->errorOut(e, "RemoveSeqsCommand", "readQual");
396                 exit(1);
397         }
398 }
399 //**********************************************************************************************************************
400 int RemoveSeqsCommand::readList(){
401         try {
402                 string thisOutputDir = outputDir;
403                 if (outputDir == "") {  thisOutputDir += m->hasPath(listfile);  }
404                 string outputFileName = thisOutputDir + m->getRootName(m->getSimpleName(listfile)) + "pick" +  m->getExtension(listfile);
405                 
406                 ofstream out;
407                 m->openOutputFile(outputFileName, out);
408                 
409                 ifstream in;
410                 m->openInputFile(listfile, in);
411                 
412                 bool wroteSomething = false;
413                 
414                 while(!in.eof()){
415                         //read in list vector
416                         ListVector list(in);
417                         
418                         //make a new list vector
419                         ListVector newList;
420                         newList.setLabel(list.getLabel());
421                         
422                         //for each bin
423                         for (int i = 0; i < list.getNumBins(); i++) {
424                                 if (m->control_pressed) { in.close();  out.close();  remove(outputFileName.c_str());  return 0; }
425                         
426                                 //parse out names that are in accnos file
427                                 string binnames = list.get(i);
428                                 
429                                 string newNames = "";
430                                 while (binnames.find_first_of(',') != -1) { 
431                                         string name = binnames.substr(0,binnames.find_first_of(','));
432                                         binnames = binnames.substr(binnames.find_first_of(',')+1, binnames.length());
433                                         
434                                         //if that name is in the .accnos file, add it
435                                         if (names.count(name) == 0) {  newNames += name + ",";  }
436                                 }
437                         
438                                 //get last name
439                                 if (names.count(binnames) == 0) {  newNames += binnames + ",";  }
440
441                                 //if there are names in this bin add to new list
442                                 if (newNames != "") {  
443                                         newNames = newNames.substr(0, newNames.length()-1); //rip off extra comma
444                                         newList.push_back(newNames);    
445                                 }
446                         }
447                                 
448                         //print new listvector
449                         if (newList.getNumBins() != 0) {
450                                 wroteSomething = true;
451                                 newList.print(out);
452                         }
453                         
454                         m->gobble(in);
455                 }
456                 in.close();     
457                 out.close();
458                 
459                 if (wroteSomething == false) {  m->mothurOut("Your file contains only sequences from the .accnos file."); m->mothurOutEndLine();  }
460                 outputTypes["list"].push_back(outputFileName); outputNames.push_back(outputFileName);
461                                 
462                 return 0;
463
464         }
465         catch(exception& e) {
466                 m->errorOut(e, "RemoveSeqsCommand", "readList");
467                 exit(1);
468         }
469 }
470 //**********************************************************************************************************************
471 int RemoveSeqsCommand::readName(){
472         try {
473                 string thisOutputDir = outputDir;
474                 if (outputDir == "") {  thisOutputDir += m->hasPath(namefile);  }
475                 string outputFileName = thisOutputDir + m->getRootName(m->getSimpleName(namefile)) + "pick" + m->getExtension(namefile);
476
477                 ofstream out;
478                 m->openOutputFile(outputFileName, out);
479
480                 ifstream in;
481                 m->openInputFile(namefile, in);
482                 string name, firstCol, secondCol;
483                 
484                 bool wroteSomething = false;
485                 
486                 while(!in.eof()){
487                         if (m->control_pressed) { in.close();  out.close();  remove(outputFileName.c_str());  return 0; }
488
489                         in >> firstCol;         m->gobble(in);          
490                         in >> secondCol;                        
491
492                         vector<string> parsedNames;
493                         //parse second column saving each name
494                         while (secondCol.find_first_of(',') != -1) { 
495                                 name = secondCol.substr(0,secondCol.find_first_of(','));
496                                 secondCol = secondCol.substr(secondCol.find_first_of(',')+1, secondCol.length());
497                                 parsedNames.push_back(name);
498                         }
499                         
500                         //get name after last ,
501                         parsedNames.push_back(secondCol);
502
503                         vector<string> validSecond;  validSecond.clear();
504                         for (int i = 0; i < parsedNames.size(); i++) {
505                                 if (names.count(parsedNames[i]) == 0) {
506                                         validSecond.push_back(parsedNames[i]);
507                                 }
508                         }
509                         
510                         if ((dups) && (validSecond.size() != parsedNames.size())) {  //if dups is true and we want to get rid of anyone, get rid of everyone
511                                 for (int i = 0; i < parsedNames.size(); i++) {  names.insert(parsedNames[i]);  }
512                         }else {
513                                         //if the name in the first column is in the set then print it and any other names in second column also in set
514                                 if (names.count(firstCol) == 0) {
515                                         
516                                         wroteSomething = true;
517                                         
518                                         out << firstCol << '\t';
519                                         
520                                         //you know you have at least one valid second since first column is valid
521                                         for (int i = 0; i < validSecond.size()-1; i++) {  out << validSecond[i] << ',';  }
522                                         out << validSecond[validSecond.size()-1] << endl;
523                                         
524                                         //make first name in set you come to first column and then add the remaining names to second column
525                                 }else {
526                                         
527                                         //you want part of this row
528                                         if (validSecond.size() != 0) {
529                                                 
530                                                 wroteSomething = true;
531                                                 
532                                                 out << validSecond[0] << '\t';
533                                                 
534                                                 //you know you have at least one valid second since first column is valid
535                                                 for (int i = 0; i < validSecond.size()-1; i++) {  out << validSecond[i] << ',';  }
536                                                 out << validSecond[validSecond.size()-1] << endl;
537                                         }
538                                 }
539                         }
540                         m->gobble(in);
541                 }
542                 in.close();
543                 out.close();
544
545                 if (wroteSomething == false) {  m->mothurOut("Your file contains only sequences from the .accnos file."); m->mothurOutEndLine();  }
546                 outputTypes["name"].push_back(outputFileName); outputNames.push_back(outputFileName);
547                                 
548                 return 0;
549         }
550         catch(exception& e) {
551                 m->errorOut(e, "RemoveSeqsCommand", "readName");
552                 exit(1);
553         }
554 }
555
556 //**********************************************************************************************************************
557 int RemoveSeqsCommand::readGroup(){
558         try {
559                 string thisOutputDir = outputDir;
560                 if (outputDir == "") {  thisOutputDir += m->hasPath(groupfile);  }
561                 string outputFileName = thisOutputDir + m->getRootName(m->getSimpleName(groupfile)) + "pick" + m->getExtension(groupfile);
562                 
563                 ofstream out;
564                 m->openOutputFile(outputFileName, out);
565
566                 ifstream in;
567                 m->openInputFile(groupfile, in);
568                 string name, group;
569                 
570                 bool wroteSomething = false;
571                 
572                 while(!in.eof()){
573                         if (m->control_pressed) { in.close();  out.close();  remove(outputFileName.c_str());  return 0; }
574                         
575                         in >> name;                             //read from first column
576                         in >> group;                    //read from second column
577                         
578                         //if this name is in the accnos file
579                         if (names.count(name) == 0) {
580                                 wroteSomething = true;
581                                 out << name << '\t' << group << endl;
582                         }
583                                         
584                         m->gobble(in);
585                 }
586                 in.close();
587                 out.close();
588                 
589                 if (wroteSomething == false) {  m->mothurOut("Your file contains only sequences from the .accnos file."); m->mothurOutEndLine();  }
590                 outputTypes["group"].push_back(outputFileName); outputNames.push_back(outputFileName);
591                 
592                 return 0;
593         }
594         catch(exception& e) {
595                 m->errorOut(e, "RemoveSeqsCommand", "readGroup");
596                 exit(1);
597         }
598 }
599 //**********************************************************************************************************************
600 int RemoveSeqsCommand::readTax(){
601         try {
602                 string thisOutputDir = outputDir;
603                 if (outputDir == "") {  thisOutputDir += m->hasPath(taxfile);  }
604                 string outputFileName = thisOutputDir + m->getRootName(m->getSimpleName(taxfile)) + "pick" + m->getExtension(taxfile);
605                 ofstream out;
606                 m->openOutputFile(outputFileName, out);
607
608                 ifstream in;
609                 m->openInputFile(taxfile, in);
610                 string name, tax;
611                 
612                 bool wroteSomething = false;
613                 
614                 while(!in.eof()){
615                         if (m->control_pressed) { in.close();  out.close();  remove(outputFileName.c_str());  return 0; }
616                         
617                         in >> name;                             //read from first column
618                         in >> tax;                      //read from second column
619                         
620                         //if this name is in the accnos file
621                         if (names.count(name) == 0) {
622                                 wroteSomething = true;
623                                 out << name << '\t' << tax << endl;
624                         }
625                                         
626                         m->gobble(in);
627                 }
628                 in.close();
629                 out.close();
630                 
631                 if (wroteSomething == false) {  m->mothurOut("Your file contains only sequences from the .accnos file."); m->mothurOutEndLine();  }
632                 outputTypes["taxonomy"].push_back(outputFileName); outputNames.push_back(outputFileName);
633                 
634                 return 0;
635         }
636         catch(exception& e) {
637                 m->errorOut(e, "RemoveSeqsCommand", "readTax");
638                 exit(1);
639         }
640 }
641 //**********************************************************************************************************************
642 //alignreport file has a column header line then all other lines contain 16 columns.  we just want the first column since that contains the name
643 int RemoveSeqsCommand::readAlign(){
644         try {
645                 string thisOutputDir = outputDir;
646                 if (outputDir == "") {  thisOutputDir += m->hasPath(alignfile);  }
647                 string outputFileName = thisOutputDir + m->getRootName(m->getSimpleName(alignfile)) + "pick.align.report";
648                 
649                 ofstream out;
650                 m->openOutputFile(outputFileName, out);
651
652                 ifstream in;
653                 m->openInputFile(alignfile, in);
654                 string name, junk;
655                 
656                 bool wroteSomething = false;
657                 
658                 //read column headers
659                 for (int i = 0; i < 16; i++) {  
660                         if (!in.eof())  {       in >> junk;      out << junk << '\t';   }
661                         else                    {       break;                  }
662                 }
663                 out << endl;
664                 
665                 while(!in.eof()){
666                         if (m->control_pressed) { in.close();  out.close();  remove(outputFileName.c_str());  return 0; }
667                         
668                         in >> name;                             //read from first column
669                         
670                         //if this name is in the accnos file
671                         if (names.count(name) == 0) {
672                                 wroteSomething = true;
673                                 
674                                 out << name << '\t';
675                                 
676                                 //read rest
677                                 for (int i = 0; i < 15; i++) {  
678                                         if (!in.eof())  {       in >> junk;      out << junk << '\t';   }
679                                         else                    {       break;                  }
680                                 }
681                                 out << endl;
682                                 
683                         }else {//still read just don't do anything with it
684                                 
685                                 //read rest
686                                 for (int i = 0; i < 15; i++) {  
687                                         if (!in.eof())  {       in >> junk;             }
688                                         else                    {       break;                  }
689                                 }
690                         }
691                         
692                         m->gobble(in);
693                 }
694                 in.close();
695                 out.close();
696                 
697                 if (wroteSomething == false) {  m->mothurOut("Your file contains only sequences from the .accnos file."); m->mothurOutEndLine();  }
698                 outputTypes["alignreport"].push_back(outputFileName); outputNames.push_back(outputFileName);
699                 
700                 return 0;
701                 
702         }
703         catch(exception& e) {
704                 m->errorOut(e, "RemoveSeqsCommand", "readAlign");
705                 exit(1);
706         }
707 }
708 //**********************************************************************************************************************
709 void RemoveSeqsCommand::readAccnos(){
710         try {
711                 
712                 ifstream in;
713                 m->openInputFile(accnosfile, in);
714                 string name;
715                 
716                 while(!in.eof()){
717                         in >> name;
718                                                 
719                         names.insert(name);
720                         
721                         m->gobble(in);
722                 }
723                 in.close();             
724
725         }
726         catch(exception& e) {
727                 m->errorOut(e, "RemoveSeqsCommand", "readAccnos");
728                 exit(1);
729         }
730 }
731
732 //**********************************************************************************************************************
733
734