]> git.donarmstrong.com Git - mothur.git/blob - removeseqscommand.cpp
Merge remote-tracking branch 'mothur/master'
[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 #include "counttable.h"
14
15 //**********************************************************************************************************************
16 vector<string> RemoveSeqsCommand::setParameters(){      
17         try {
18                 CommandParameter pfasta("fasta", "InputTypes", "", "", "none", "FNGLT", "none",false,false); parameters.push_back(pfasta);
19         CommandParameter pname("name", "InputTypes", "", "", "NameCount", "FNGLT", "none",false,false); parameters.push_back(pname);
20         CommandParameter pcount("count", "InputTypes", "", "", "NameCount-CountGroup", "FNGLT", "none",false,false); parameters.push_back(pcount);
21                 CommandParameter pgroup("group", "InputTypes", "", "", "CountGroup", "FNGLT", "none",false,false); parameters.push_back(pgroup);
22                 CommandParameter plist("list", "InputTypes", "", "", "none", "FNGLT", "none",false,false); parameters.push_back(plist);
23                 CommandParameter ptaxonomy("taxonomy", "InputTypes", "", "", "none", "FNGLT", "none",false,false); parameters.push_back(ptaxonomy);
24                 CommandParameter palignreport("alignreport", "InputTypes", "", "", "none", "FNGLT", "none",false,false); parameters.push_back(palignreport);
25                 CommandParameter pqfile("qfile", "InputTypes", "", "", "none", "FNGLT", "none",false,false); parameters.push_back(pqfile);
26                 CommandParameter paccnos("accnos", "InputTypes", "", "", "none", "none", "none",false,true); parameters.push_back(paccnos);
27                 CommandParameter pdups("dups", "Boolean", "", "T", "", "", "",false,false); parameters.push_back(pdups);
28                 CommandParameter pinputdir("inputdir", "String", "", "", "", "", "",false,false); parameters.push_back(pinputdir);
29                 CommandParameter poutputdir("outputdir", "String", "", "", "", "", "",false,false); parameters.push_back(poutputdir);
30                 
31                 vector<string> myArray;
32                 for (int i = 0; i < parameters.size(); i++) {   myArray.push_back(parameters[i].name);          }
33                 return myArray;
34         }
35         catch(exception& e) {
36                 m->errorOut(e, "RemoveSeqsCommand", "setParameters");
37                 exit(1);
38         }
39 }
40 //**********************************************************************************************************************
41 string RemoveSeqsCommand::getHelpString(){      
42         try {
43                 string helpString = "";
44                 helpString += "The remove.seqs command reads an .accnos file and at least one of the following file types: fasta, name, group, count, list, taxonomy, quality or alignreport file.\n";
45                 helpString += "It outputs a file containing the sequences NOT in the .accnos file.\n";
46                 helpString += "The remove.seqs command parameters are accnos, fasta, name, group, count, list, taxonomy, qfile, alignreport and dups.  You must provide accnos and at least one of the file parameters.\n";
47                 helpString += "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";
48                 helpString += "The remove.seqs command should be in the following format: remove.seqs(accnos=yourAccnos, fasta=yourFasta).\n";
49                 helpString += "Example remove.seqs(accnos=amazon.accnos, fasta=amazon.fasta).\n";
50                 helpString += "Note: No spaces between parameter labels (i.e. fasta), '=' and parameters (i.e.yourFasta).\n";
51                 return helpString;
52         }
53         catch(exception& e) {
54                 m->errorOut(e, "RemoveSeqsCommand", "getHelpString");
55                 exit(1);
56         }
57 }
58 //**********************************************************************************************************************
59 string RemoveSeqsCommand::getOutputFileNameTag(string type, string inputName=""){       
60         try {
61         string outputFileName = "";
62                 map<string, vector<string> >::iterator it;
63         
64         //is this a type this command creates
65         it = outputTypes.find(type);
66         if (it == outputTypes.end()) {  m->mothurOut("[ERROR]: this command doesn't create a " + type + " output file.\n"); }
67         else {
68             if (type == "fasta")            {   outputFileName =  "pick" + m->getExtension(inputName);   }
69             else if (type == "taxonomy")    {   outputFileName =  "pick" + m->getExtension(inputName);   }
70             else if (type == "name")        {   outputFileName =  "pick" + m->getExtension(inputName);   }
71             else if (type == "group")       {   outputFileName =  "pick" + m->getExtension(inputName);   }
72             else if (type == "list")        {   outputFileName =  "pick" + m->getExtension(inputName);   }
73             else if (type == "qfile")       {   outputFileName =  "pick" + m->getExtension(inputName);   }
74             else if (type == "alignreport") {   outputFileName =  "pick.align.report";                   }
75             else if (type == "count")       {   outputFileName =  "pick.count_table";   }
76             else { m->mothurOut("[ERROR]: No definition for type " + type + " output file tag.\n"); m->control_pressed = true;  }
77         }
78         return outputFileName;
79         }
80         catch(exception& e) {
81                 m->errorOut(e, "RemoveSeqsCommand", "getOutputFileNameTag");
82                 exit(1);
83         }
84 }
85
86 //**********************************************************************************************************************
87 RemoveSeqsCommand::RemoveSeqsCommand(){ 
88         try {
89                 abort = true; calledHelp = true; 
90                 setParameters();
91                 vector<string> tempOutNames;
92                 outputTypes["fasta"] = tempOutNames;
93                 outputTypes["taxonomy"] = tempOutNames;
94                 outputTypes["name"] = tempOutNames;
95                 outputTypes["group"] = tempOutNames;
96                 outputTypes["alignreport"] = tempOutNames;
97                 outputTypes["list"] = tempOutNames;
98                 outputTypes["qfile"] = tempOutNames;
99         outputTypes["count"] = tempOutNames;
100         }
101         catch(exception& e) {
102                 m->errorOut(e, "RemoveSeqsCommand", "RemoveSeqsCommand");
103                 exit(1);
104         }
105 }
106 //**********************************************************************************************************************
107 RemoveSeqsCommand::RemoveSeqsCommand(string option)  {
108         try {
109                 abort = false; calledHelp = false;   
110                 
111                 //allow user to run help
112                 if(option == "help") { help(); abort = true; calledHelp = true; }
113                 else if(option == "citation") { citation(); abort = true; calledHelp = true;}
114                 
115                 else {
116                         vector<string> myArray = setParameters();
117                         
118                         OptionParser parser(option);
119                         map<string,string> parameters = parser.getParameters();
120                         
121                         ValidParameters validParameter;
122                         map<string,string>::iterator it;
123                         
124                         //check to make sure all parameters are valid for command
125                         for (it = parameters.begin(); it != parameters.end(); it++) { 
126                                 if (validParameter.isValidParameter(it->first, myArray, it->second) != true) {  abort = true;  }
127                         }
128                         
129                         //initialize outputTypes
130                         vector<string> tempOutNames;
131                         outputTypes["fasta"] = tempOutNames;
132                         outputTypes["taxonomy"] = tempOutNames;
133                         outputTypes["name"] = tempOutNames;
134                         outputTypes["group"] = tempOutNames;
135                         outputTypes["alignreport"] = tempOutNames;
136                         outputTypes["list"] = tempOutNames;
137                         outputTypes["qfile"] = tempOutNames;
138             outputTypes["count"] = tempOutNames;
139                         
140                         //if the user changes the output directory command factory will send this info to us in the output parameter 
141                         outputDir = validParameter.validFile(parameters, "outputdir", false);           if (outputDir == "not found"){  outputDir = "";         }
142                         
143                         //if the user changes the input directory command factory will send this info to us in the output parameter 
144                         string inputDir = validParameter.validFile(parameters, "inputdir", false);              
145                         if (inputDir == "not found"){   inputDir = "";          }
146                         else {
147                                 string path;
148                                 it = parameters.find("alignreport");
149                                 //user has given a template file
150                                 if(it != parameters.end()){ 
151                                         path = m->hasPath(it->second);
152                                         //if the user has not given a path then, add inputdir. else leave path alone.
153                                         if (path == "") {       parameters["alignreport"] = inputDir + it->second;              }
154                                 }
155                                 
156                                 it = parameters.find("fasta");
157                                 //user has given a template file
158                                 if(it != parameters.end()){ 
159                                         path = m->hasPath(it->second);
160                                         //if the user has not given a path then, add inputdir. else leave path alone.
161                                         if (path == "") {       parameters["fasta"] = inputDir + it->second;            }
162                                 }
163                                 
164                                 it = parameters.find("accnos");
165                                 //user has given a template file
166                                 if(it != parameters.end()){ 
167                                         path = m->hasPath(it->second);
168                                         //if the user has not given a path then, add inputdir. else leave path alone.
169                                         if (path == "") {       parameters["accnos"] = inputDir + it->second;           }
170                                 }
171                                 
172                                 it = parameters.find("list");
173                                 //user has given a template file
174                                 if(it != parameters.end()){ 
175                                         path = m->hasPath(it->second);
176                                         //if the user has not given a path then, add inputdir. else leave path alone.
177                                         if (path == "") {       parameters["list"] = inputDir + it->second;             }
178                                 }
179                                 
180                                 it = parameters.find("name");
181                                 //user has given a template file
182                                 if(it != parameters.end()){ 
183                                         path = m->hasPath(it->second);
184                                         //if the user has not given a path then, add inputdir. else leave path alone.
185                                         if (path == "") {       parameters["name"] = inputDir + it->second;             }
186                                 }
187                                 
188                                 it = parameters.find("group");
189                                 //user has given a template file
190                                 if(it != parameters.end()){ 
191                                         path = m->hasPath(it->second);
192                                         //if the user has not given a path then, add inputdir. else leave path alone.
193                                         if (path == "") {       parameters["group"] = inputDir + it->second;            }
194                                 }
195                                 
196                                 it = parameters.find("taxonomy");
197                                 //user has given a template file
198                                 if(it != parameters.end()){ 
199                                         path = m->hasPath(it->second);
200                                         //if the user has not given a path then, add inputdir. else leave path alone.
201                                         if (path == "") {       parameters["taxonomy"] = inputDir + it->second;         }
202                                 }
203                                 
204                                 it = parameters.find("qfile");
205                                 //user has given a template file
206                                 if(it != parameters.end()){ 
207                                         path = m->hasPath(it->second);
208                                         //if the user has not given a path then, add inputdir. else leave path alone.
209                                         if (path == "") {       parameters["qfile"] = inputDir + it->second;            }
210                                 }
211                 
212                 it = parameters.find("count");
213                                 //user has given a template file
214                                 if(it != parameters.end()){ 
215                                         path = m->hasPath(it->second);
216                                         //if the user has not given a path then, add inputdir. else leave path alone.
217                                         if (path == "") {       parameters["count"] = inputDir + it->second;            }
218                                 }
219                         }
220
221                         
222                         //check for required parameters
223                         accnosfile = validParameter.validFile(parameters, "accnos", true);
224                         if (accnosfile == "not open") { abort = true; }
225                         else if (accnosfile == "not found") {  
226                                 accnosfile = m->getAccnosFile(); 
227                                 if (accnosfile != "") {  m->mothurOut("Using " + accnosfile + " as input file for the accnos parameter."); m->mothurOutEndLine(); }
228                                 else { 
229                                         m->mothurOut("You have no valid accnos file and accnos is required."); m->mothurOutEndLine(); 
230                                         abort = true;
231                                 }  
232                         }else { m->setAccnosFile(accnosfile); } 
233                         
234                         fastafile = validParameter.validFile(parameters, "fasta", true);
235                         if (fastafile == "not open") { fastafile = ""; abort = true; }
236                         else if (fastafile == "not found") {  fastafile = "";  }        
237                         else { m->setFastaFile(fastafile); }
238                                                                    
239                         namefile = validParameter.validFile(parameters, "name", true);
240                         if (namefile == "not open") { namefile = ""; abort = true; }
241                         else if (namefile == "not found") {  namefile = "";  }  
242                         else { m->setNameFile(namefile); } 
243                                                                    
244                         groupfile = validParameter.validFile(parameters, "group", true);
245                         if (groupfile == "not open") { abort = true; }
246                         else if (groupfile == "not found") {  groupfile = "";  }
247                         else { m->setGroupFile(groupfile); }
248                         
249                         alignfile = validParameter.validFile(parameters, "alignreport", true);
250                         if (alignfile == "not open") { abort = true; }
251                         else if (alignfile == "not found") {  alignfile = "";  }
252                         
253                         listfile = validParameter.validFile(parameters, "list", true);
254                         if (listfile == "not open") { abort = true; }
255                         else if (listfile == "not found") {  listfile = "";  }
256                         else { m->setListFile(listfile); }
257                         
258                         taxfile = validParameter.validFile(parameters, "taxonomy", true);
259                         if (taxfile == "not open") { abort = true; }
260                         else if (taxfile == "not found") {  taxfile = "";  }
261                         else { m->setTaxonomyFile(taxfile); }
262                         
263                         qualfile = validParameter.validFile(parameters, "qfile", true);
264                         if (qualfile == "not open") { abort = true; }
265                         else if (qualfile == "not found") {  qualfile = "";  }                  
266                         else { m->setQualFile(qualfile); }
267                         
268                         string usedDups = "true";
269                         string temp = validParameter.validFile(parameters, "dups", false);      
270                         if (temp == "not found") { 
271                                 if (namefile != "") {  temp = "true";                                   }
272                                 else                            {  temp = "false"; usedDups = "";       }
273                         }
274                         dups = m->isTrue(temp);
275             
276             countfile = validParameter.validFile(parameters, "count", true);
277             if (countfile == "not open") { countfile = ""; abort = true; }
278             else if (countfile == "not found") { countfile = "";  }     
279             else { m->setCountTableFile(countfile); }
280             
281             if ((namefile != "") && (countfile != "")) {
282                 m->mothurOut("[ERROR]: you may only use one of the following: name or count."); m->mothurOutEndLine(); abort = true;
283             }
284             
285             if ((groupfile != "") && (countfile != "")) {
286                 m->mothurOut("[ERROR]: you may only use one of the following: group or count."); m->mothurOutEndLine(); abort=true;
287             }
288                         
289                         if ((countfile == "") && (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; }
290                         
291             if (countfile == "") {
292                 if ((fastafile != "") && (namefile == "")) {
293                     vector<string> files; files.push_back(fastafile);
294                     parser.getNameFile(files);
295                 }
296             }
297                 }
298
299         }
300         catch(exception& e) {
301                 m->errorOut(e, "RemoveSeqsCommand", "RemoveSeqsCommand");
302                 exit(1);
303         }
304 }
305 //**********************************************************************************************************************
306
307 int RemoveSeqsCommand::execute(){
308         try {
309                 
310                 if (abort == true) { if (calledHelp) { return 0; }  return 2;   }
311                 
312                 //get names you want to keep
313                 names = m->readAccnos(accnosfile);
314                 
315                 if (m->control_pressed) { return 0; }
316         
317         if (countfile != "") {
318             if ((fastafile != "") || (listfile != "") || (taxfile != "")) { 
319                 m->mothurOut("\n[NOTE]: The count file should contain only unique names, so mothur assumes your fasta, list and taxonomy files also contain only uniques.\n\n");
320             }
321         }
322                 
323                 //read through the correct file and output lines you want to keep
324                 if (namefile != "")                     {               readName();             }
325                 if (fastafile != "")            {               readFasta();    }
326                 if (groupfile != "")            {               readGroup();    }
327                 if (alignfile != "")            {               readAlign();    }
328                 if (listfile != "")                     {               readList();             }
329                 if (taxfile != "")                      {               readTax();              }
330                 if (qualfile != "")                     {               readQual();             }
331         if (countfile != "")            {               readCount();            }
332                 
333                 if (m->control_pressed) { for (int i = 0; i < outputNames.size(); i++) {        m->mothurRemove(outputNames[i]); } return 0; }
334         
335                 if (outputNames.size() != 0) {
336                         m->mothurOutEndLine();
337                         m->mothurOut("Output File Names: "); m->mothurOutEndLine();
338                         for (int i = 0; i < outputNames.size(); i++) {  m->mothurOut(outputNames[i]); m->mothurOutEndLine();    }
339                         m->mothurOutEndLine();
340                         
341                         //set fasta file as new current fastafile
342                         string current = "";
343                         itTypes = outputTypes.find("fasta");
344                         if (itTypes != outputTypes.end()) {
345                                 if ((itTypes->second).size() != 0) { current = (itTypes->second)[0]; m->setFastaFile(current); }
346                         }
347                         
348                         itTypes = outputTypes.find("name");
349                         if (itTypes != outputTypes.end()) {
350                                 if ((itTypes->second).size() != 0) { current = (itTypes->second)[0]; m->setNameFile(current); }
351                         }
352                         
353                         itTypes = outputTypes.find("group");
354                         if (itTypes != outputTypes.end()) {
355                                 if ((itTypes->second).size() != 0) { current = (itTypes->second)[0]; m->setGroupFile(current); }
356                         }
357                         
358                         itTypes = outputTypes.find("list");
359                         if (itTypes != outputTypes.end()) {
360                                 if ((itTypes->second).size() != 0) { current = (itTypes->second)[0]; m->setListFile(current); }
361                         }
362                         
363                         itTypes = outputTypes.find("taxonomy");
364                         if (itTypes != outputTypes.end()) {
365                                 if ((itTypes->second).size() != 0) { current = (itTypes->second)[0]; m->setTaxonomyFile(current); }
366                         }
367                         
368                         itTypes = outputTypes.find("qfile");
369                         if (itTypes != outputTypes.end()) {
370                                 if ((itTypes->second).size() != 0) { current = (itTypes->second)[0]; m->setQualFile(current); }
371                         }       
372             
373             itTypes = outputTypes.find("count");
374                         if (itTypes != outputTypes.end()) {
375                                 if ((itTypes->second).size() != 0) { current = (itTypes->second)[0]; m->setCountTableFile(current); }
376                         }
377                 }
378                 
379                 return 0;               
380         }
381
382         catch(exception& e) {
383                 m->errorOut(e, "RemoveSeqsCommand", "execute");
384                 exit(1);
385         }
386 }
387
388 //**********************************************************************************************************************
389 int RemoveSeqsCommand::readFasta(){
390         try {
391                 string thisOutputDir = outputDir;
392                 if (outputDir == "") {  thisOutputDir += m->hasPath(fastafile);  }
393                 string outputFileName = thisOutputDir + m->getRootName(m->getSimpleName(fastafile)) + getOutputFileNameTag("fasta", fastafile);
394                 
395                 ofstream out;
396                 m->openOutputFile(outputFileName, out);
397                 
398                 ifstream in;
399                 m->openInputFile(fastafile, in);
400                 string name;
401                 
402                 bool wroteSomething = false;
403                 int removedCount = 0;
404                 
405                 while(!in.eof()){
406                         if (m->control_pressed) { in.close();  out.close();  m->mothurRemove(outputFileName);  return 0; }
407                         
408                         Sequence currSeq(in);
409             
410             if (!dups) {//adjust name if needed
411                 map<string, string>::iterator it = uniqueMap.find(currSeq.getName());
412                 if (it != uniqueMap.end()) { currSeq.setName(it->second); }
413             }
414
415                         name = currSeq.getName();
416                         
417                         if (name != "") {
418                                 //if this name is in the accnos file
419                                 if (names.count(name) == 0) {
420                                         wroteSomething = true;
421                                         
422                     currSeq.printSequence(out);
423                                 }else {  removedCount++;  }
424                         }
425                         m->gobble(in);
426                 }
427                 in.close();     
428                 out.close();
429                 
430                 if (wroteSomething == false) {  m->mothurOut("Your file contains only sequences from the .accnos file."); m->mothurOutEndLine();  }
431                 outputTypes["fasta"].push_back(outputFileName);  outputNames.push_back(outputFileName);
432                 
433                 m->mothurOut("Removed " + toString(removedCount) + " sequences from your fasta file."); m->mothurOutEndLine();
434                 
435                 return 0;
436                 
437         }
438         catch(exception& e) {
439                 m->errorOut(e, "RemoveSeqsCommand", "readFasta");
440                 exit(1);
441         }
442 }
443 //**********************************************************************************************************************
444 int RemoveSeqsCommand::readQual(){
445         try {
446                 string thisOutputDir = outputDir;
447                 if (outputDir == "") {  thisOutputDir += m->hasPath(qualfile);  }
448                 string outputFileName = thisOutputDir + m->getRootName(m->getSimpleName(qualfile)) + getOutputFileNameTag("qfile", qualfile);
449                 ofstream out;
450                 m->openOutputFile(outputFileName, out);
451                 
452                 
453                 ifstream in;
454                 m->openInputFile(qualfile, in);
455                 string name;
456                 
457                 bool wroteSomething = false;
458                 int removedCount = 0;
459                 
460                 
461                 while(!in.eof()){       
462                         string saveName = "";
463                         string name = "";
464                         string scores = "";
465                         
466                         in >> name; 
467                         
468                         if (name.length() != 0) { 
469                                 saveName = name.substr(1);
470                                 while (!in.eof())       {       
471                                         char c = in.get(); 
472                                         if (c == 10 || c == 13){        break;  }
473                                         else { name += c; }     
474                                 } 
475                                 m->gobble(in);
476                         }
477                         
478                         while(in){
479                                 char letter= in.get();
480                                 if(letter == '>'){      in.putback(letter);     break;  }
481                                 else{ scores += letter; }
482                         }
483                         
484                         m->gobble(in);
485                         
486             if (!dups) {//adjust name if needed
487                 map<string, string>::iterator it = uniqueMap.find(saveName);
488                 if (it != uniqueMap.end()) { name = ">" + it->second; saveName = it->second; }
489             }
490             
491                         if (names.count(saveName) == 0) {
492                                 wroteSomething = true;
493                                 
494                                 out << name << endl << scores;
495                         }else {  removedCount++;  }
496                         
497                         m->gobble(in);
498                 }
499                 in.close();
500                 out.close();
501                 
502                 
503                 if (wroteSomething == false) { m->mothurOut("Your file contains only sequences from the .accnos file."); m->mothurOutEndLine();  }
504                 outputNames.push_back(outputFileName);  outputTypes["qfile"].push_back(outputFileName); 
505                 
506                 m->mothurOut("Removed " + toString(removedCount) + " sequences from your quality file."); m->mothurOutEndLine();
507                 
508                 return 0;
509                 
510         }
511         catch(exception& e) {
512                 m->errorOut(e, "RemoveSeqsCommand", "readQual");
513                 exit(1);
514         }
515 }
516 //**********************************************************************************************************************
517 int RemoveSeqsCommand::readCount(){
518         try {
519         
520                 string thisOutputDir = outputDir;
521                 if (outputDir == "") {  thisOutputDir += m->hasPath(countfile);  }
522                 string outputFileName = thisOutputDir + m->getRootName(m->getSimpleName(countfile)) + getOutputFileNameTag("count", countfile);
523                 
524                 ofstream out;
525                 m->openOutputFile(outputFileName, out);
526                 
527                 ifstream in;
528                 m->openInputFile(countfile, in);
529                 
530                 bool wroteSomething = false;
531                 int removedCount = 0;
532                 
533         string headers = m->getline(in); m->gobble(in);
534         out << headers << endl;
535         
536         string name, rest; int thisTotal;
537         while (!in.eof()) {
538             
539             if (m->control_pressed) { in.close();  out.close();  m->mothurRemove(outputFileName);  return 0; }
540             
541             in >> name; m->gobble(in); 
542             in >> thisTotal; m->gobble(in);
543             rest = m->getline(in); m->gobble(in);
544             if (m->debug) { m->mothurOut("[DEBUG]: " + name + '\t' + rest + "\n"); }
545             
546             if (names.count(name) == 0) {
547                 out << name << '\t' << thisTotal << '\t' << rest << endl;
548                 wroteSomething = true;
549             }else { removedCount += thisTotal; }
550         }
551         in.close();
552                 out.close();
553         
554         //check for groups that have been eliminated
555         CountTable ct;
556         if (ct.testGroups(outputFileName)) {
557             ct.readTable(outputFileName);
558             ct.printTable(outputFileName);
559         }
560
561                 
562                 if (wroteSomething == false) {  m->mothurOut("Your file contains only sequences from the .accnos file."); m->mothurOutEndLine();  }
563                 outputTypes["count"].push_back(outputFileName); outputNames.push_back(outputFileName);
564                 
565                 m->mothurOut("Removed " + toString(removedCount) + " sequences from your count file."); m->mothurOutEndLine();
566         
567                 return 0;
568         }
569         catch(exception& e) {
570                 m->errorOut(e, "RemoveSeqsCommand", "readCount");
571                 exit(1);
572         }
573 }
574 //**********************************************************************************************************************
575 int RemoveSeqsCommand::readList(){
576         try {
577                 string thisOutputDir = outputDir;
578                 if (outputDir == "") {  thisOutputDir += m->hasPath(listfile);  }
579                 string outputFileName = thisOutputDir + m->getRootName(m->getSimpleName(listfile)) + getOutputFileNameTag("list", listfile);            
580                 ofstream out;
581                 m->openOutputFile(outputFileName, out);
582                 
583                 ifstream in;
584                 m->openInputFile(listfile, in);
585                 
586                 bool wroteSomething = false;
587                 int removedCount = 0;
588                 
589                 while(!in.eof()){
590                         
591                         removedCount = 0;
592                         
593                         //read in list vector
594                         ListVector list(in);
595                         
596                         //make a new list vector
597                         ListVector newList;
598                         newList.setLabel(list.getLabel());
599                         
600                         //for each bin
601                         for (int i = 0; i < list.getNumBins(); i++) {
602                                 if (m->control_pressed) { in.close();  out.close();  m->mothurRemove(outputFileName);  return 0; }
603                         
604                                 //parse out names that are in accnos file
605                                 string binnames = list.get(i);
606                                 
607                                 string newNames = "";
608                                 while (binnames.find_first_of(',') != -1) { 
609                                         string name = binnames.substr(0,binnames.find_first_of(','));
610                                         binnames = binnames.substr(binnames.find_first_of(',')+1, binnames.length());
611                                         
612                                         //if that name is in the .accnos file, add it
613                                         if (names.count(name) == 0) {  newNames += name + ",";  }
614                                         else {  removedCount++;  }
615                                 }
616                         
617                                 //get last name
618                                 if (names.count(binnames) == 0) {  newNames += binnames + ",";  }
619                                 else {  removedCount++;  }
620
621                                 //if there are names in this bin add to new list
622                                 if (newNames != "") {  
623                                         newNames = newNames.substr(0, newNames.length()-1); //rip off extra comma
624                                         newList.push_back(newNames);    
625                                 }
626                         }
627                                 
628                         //print new listvector
629                         if (newList.getNumBins() != 0) {
630                                 wroteSomething = true;
631                                 newList.print(out);
632                         }
633                         
634                         m->gobble(in);
635                 }
636                 in.close();     
637                 out.close();
638                 
639                 if (wroteSomething == false) {  m->mothurOut("Your file contains only sequences from the .accnos file."); m->mothurOutEndLine();  }
640                 outputTypes["list"].push_back(outputFileName); outputNames.push_back(outputFileName);
641                 
642                 m->mothurOut("Removed " + toString(removedCount) + " sequences from your list file."); m->mothurOutEndLine();
643                 
644                 return 0;
645
646         }
647         catch(exception& e) {
648                 m->errorOut(e, "RemoveSeqsCommand", "readList");
649                 exit(1);
650         }
651 }
652 //**********************************************************************************************************************
653 int RemoveSeqsCommand::readName(){
654         try {
655                 string thisOutputDir = outputDir;
656                 if (outputDir == "") {  thisOutputDir += m->hasPath(namefile);  }
657                 string outputFileName = thisOutputDir + m->getRootName(m->getSimpleName(namefile)) + getOutputFileNameTag("name", namefile);
658                 ofstream out;
659                 m->openOutputFile(outputFileName, out);
660
661                 ifstream in;
662                 m->openInputFile(namefile, in);
663                 string name, firstCol, secondCol;
664                 
665                 bool wroteSomething = false;
666                 int removedCount = 0;
667                 
668                 while(!in.eof()){
669                         if (m->control_pressed) { in.close();  out.close();  m->mothurRemove(outputFileName);  return 0; }
670                         
671                         in >> firstCol;         m->gobble(in);          
672                         in >> secondCol;                        
673                         
674                         vector<string> parsedNames;
675                         m->splitAtComma(secondCol, parsedNames);
676                         
677                         vector<string> validSecond;  validSecond.clear();
678                         for (int i = 0; i < parsedNames.size(); i++) {
679                                 if (names.count(parsedNames[i]) == 0) {
680                                         validSecond.push_back(parsedNames[i]);
681                                 }
682                         }
683                         
684                         if ((dups) && (validSecond.size() != parsedNames.size())) {  //if dups is true and we want to get rid of anyone, get rid of everyone
685                                 for (int i = 0; i < parsedNames.size(); i++) {  names.insert(parsedNames[i]);  }
686                                 removedCount += parsedNames.size();
687                         }else {
688                                 removedCount += parsedNames.size()-validSecond.size();
689                                 //if the name in the first column is in the set then print it and any other names in second column also in set
690                                 if (names.count(firstCol) == 0) {
691                                         
692                                         wroteSomething = true;
693                                         
694                                         out << firstCol << '\t';
695                                         
696                                         //you know you have at least one valid second since first column is valid
697                                         for (int i = 0; i < validSecond.size()-1; i++) {  out << validSecond[i] << ',';  }
698                                         out << validSecond[validSecond.size()-1] << endl;
699                                         
700                                         //make first name in set you come to first column and then add the remaining names to second column
701                                 }else {
702                                         
703                                         //you want part of this row
704                                         if (validSecond.size() != 0) {
705                                                 
706                                                 wroteSomething = true;
707                                                 
708                                                 out << validSecond[0] << '\t';
709                         //we are changing the unique name in the fasta file
710                         uniqueMap[firstCol] = validSecond[0];
711                                                 
712                                                 //you know you have at least one valid second since first column is valid
713                                                 for (int i = 0; i < validSecond.size()-1; i++) {  out << validSecond[i] << ',';  }
714                                                 out << validSecond[validSecond.size()-1] << endl;
715                                         }
716                                 }
717                         }
718                         m->gobble(in);
719                 }
720                 in.close();
721                 out.close();
722                 
723                 if (wroteSomething == false) {  m->mothurOut("Your file contains only sequences from the .accnos file."); m->mothurOutEndLine();  }
724                 outputTypes["name"].push_back(outputFileName); outputNames.push_back(outputFileName);
725                 
726                 m->mothurOut("Removed " + toString(removedCount) + " sequences from your name file."); m->mothurOutEndLine();
727                 
728                 return 0;
729         }
730         catch(exception& e) {
731                 m->errorOut(e, "RemoveSeqsCommand", "readName");
732                 exit(1);
733         }
734 }
735
736 //**********************************************************************************************************************
737 int RemoveSeqsCommand::readGroup(){
738         try {
739                 string thisOutputDir = outputDir;
740                 if (outputDir == "") {  thisOutputDir += m->hasPath(groupfile);  }
741                 string outputFileName = thisOutputDir + m->getRootName(m->getSimpleName(groupfile)) + getOutputFileNameTag("group", groupfile);         
742                 ofstream out;
743                 m->openOutputFile(outputFileName, out);
744
745                 ifstream in;
746                 m->openInputFile(groupfile, in);
747                 string name, group;
748                 
749                 bool wroteSomething = false;
750                 int removedCount = 0;
751                 
752                 while(!in.eof()){
753                         if (m->control_pressed) { in.close();  out.close();  m->mothurRemove(outputFileName);  return 0; }
754                         
755                         in >> name;                             //read from first column
756                         in >> group;                    //read from second column
757                         
758                         //if this name is in the accnos file
759                         if (names.count(name) == 0) {
760                                 wroteSomething = true;
761                                 out << name << '\t' << group << endl;
762                         }else {  removedCount++;  }
763                                         
764                         m->gobble(in);
765                 }
766                 in.close();
767                 out.close();
768                 
769                 if (wroteSomething == false) {  m->mothurOut("Your file contains only sequences from the .accnos file."); m->mothurOutEndLine();  }
770                 outputTypes["group"].push_back(outputFileName); outputNames.push_back(outputFileName);
771                 
772                 m->mothurOut("Removed " + toString(removedCount) + " sequences from your group file."); m->mothurOutEndLine();
773
774                 
775                 return 0;
776         }
777         catch(exception& e) {
778                 m->errorOut(e, "RemoveSeqsCommand", "readGroup");
779                 exit(1);
780         }
781 }
782 //**********************************************************************************************************************
783 int RemoveSeqsCommand::readTax(){
784         try {
785                 string thisOutputDir = outputDir;
786                 if (outputDir == "") {  thisOutputDir += m->hasPath(taxfile);  }
787                 string outputFileName = thisOutputDir + m->getRootName(m->getSimpleName(taxfile)) + getOutputFileNameTag("taxonomy", taxfile);
788                 ofstream out;
789                 m->openOutputFile(outputFileName, out);
790
791                 ifstream in;
792                 m->openInputFile(taxfile, in);
793                 string name, tax;
794                 
795                 bool wroteSomething = false;
796                 int removedCount = 0;
797                 
798                 while(!in.eof()){
799                         if (m->control_pressed) { in.close();  out.close();  m->mothurRemove(outputFileName);  return 0; }
800                         
801                         in >> name;                             //read from first column
802                         in >> tax;                      //read from second column
803                         
804             if (!dups) {//adjust name if needed
805                 map<string, string>::iterator it = uniqueMap.find(name);
806                 if (it != uniqueMap.end()) { name = it->second; }
807             }
808             
809                         //if this name is in the accnos file
810                         if (names.count(name) == 0) {
811                                 wroteSomething = true;
812             
813                                 out << name << '\t' << tax << endl;
814                         }else {  removedCount++;  }
815                                         
816                         m->gobble(in);
817                 }
818                 in.close();
819                 out.close();
820                 
821                 if (wroteSomething == false) {  m->mothurOut("Your file contains only sequences from the .accnos file."); m->mothurOutEndLine();  }
822                 outputTypes["taxonomy"].push_back(outputFileName); outputNames.push_back(outputFileName);
823                 
824                 m->mothurOut("Removed " + toString(removedCount) + " sequences from your taxonomy file."); m->mothurOutEndLine();
825                 
826                 return 0;
827         }
828         catch(exception& e) {
829                 m->errorOut(e, "RemoveSeqsCommand", "readTax");
830                 exit(1);
831         }
832 }
833 //**********************************************************************************************************************
834 //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
835 int RemoveSeqsCommand::readAlign(){
836         try {
837                 string thisOutputDir = outputDir;
838                 if (outputDir == "") {  thisOutputDir += m->hasPath(alignfile);  }
839                 string outputFileName = thisOutputDir + m->getRootName(m->getSimpleName(alignfile)) + getOutputFileNameTag("alignreport");
840                 
841                 ofstream out;
842                 m->openOutputFile(outputFileName, out);
843
844                 ifstream in;
845                 m->openInputFile(alignfile, in);
846                 string name, junk;
847                 
848                 bool wroteSomething = false;
849                 int removedCount = 0;
850                 
851                 //read column headers
852                 for (int i = 0; i < 16; i++) {  
853                         if (!in.eof())  {       in >> junk;      out << junk << '\t';   }
854                         else                    {       break;                  }
855                 }
856                 out << endl;
857                 
858                 while(!in.eof()){
859                         if (m->control_pressed) { in.close();  out.close();  m->mothurRemove(outputFileName);  return 0; }
860                         
861                         in >> name;                             //read from first column
862             
863             if (!dups) {//adjust name if needed
864                 map<string, string>::iterator it = uniqueMap.find(name);
865                 if (it != uniqueMap.end()) { name = it->second; }
866             }
867                         
868                         //if this name is in the accnos file
869                         if (names.count(name) == 0) {
870                                 wroteSomething = true;
871                                 
872                                 out << name << '\t';
873                                 
874                                 //read rest
875                                 for (int i = 0; i < 15; i++) {  
876                                         if (!in.eof())  {       in >> junk;      out << junk << '\t';   }
877                                         else                    {       break;                  }
878                                 }
879                                 out << endl;
880                                 
881                         }else {//still read just don't do anything with it
882                                 removedCount++;  
883                                 
884                                 //read rest
885                                 for (int i = 0; i < 15; i++) {  
886                                         if (!in.eof())  {       in >> junk;             }
887                                         else                    {       break;                  }
888                                 }
889                         }
890                         
891                         m->gobble(in);
892                 }
893                 in.close();
894                 out.close();
895                 
896                 if (wroteSomething == false) {  m->mothurOut("Your file contains only sequences from the .accnos file."); m->mothurOutEndLine();  }
897                 outputTypes["alignreport"].push_back(outputFileName); outputNames.push_back(outputFileName);
898                 
899                 m->mothurOut("Removed " + toString(removedCount) + " sequences from your alignreport file."); m->mothurOutEndLine();
900
901                 
902                 return 0;
903                 
904         }
905         catch(exception& e) {
906                 m->errorOut(e, "RemoveSeqsCommand", "readAlign");
907                 exit(1);
908         }
909 }
910 //**********************************************************************************************************************
911
912