]> git.donarmstrong.com Git - mothur.git/blobdiff - getseqscommand.cpp
changed dups parameter to dereplicate in chimera.uchime.
[mothur.git] / getseqscommand.cpp
index 84fe0acab9abc0ec603df10896be1c0fe1e7f47c..6b16111fdcda83099c0e36396543d40ee345c865 100644 (file)
 #include "getseqscommand.h"
 #include "sequence.hpp"
 #include "listvector.hpp"
+#include "counttable.h"
 
 //**********************************************************************************************************************
 vector<string> GetSeqsCommand::setParameters(){        
        try {
                CommandParameter pfasta("fasta", "InputTypes", "", "", "none", "FNGLT", "none",false,false); parameters.push_back(pfasta);
-               CommandParameter pname("name", "InputTypes", "", "", "none", "FNGLT", "none",false,false); parameters.push_back(pname);
-               CommandParameter pgroup("group", "InputTypes", "", "", "none", "FNGLT", "none",false,false); parameters.push_back(pgroup);
+        CommandParameter pname("name", "InputTypes", "", "", "NameCount", "FNGLT", "none",false,false); parameters.push_back(pname);
+        CommandParameter pcount("count", "InputTypes", "", "", "NameCount-CountGroup", "FNGLT", "none",false,false); parameters.push_back(pcount);
+               CommandParameter pgroup("group", "InputTypes", "", "", "CountGroup", "FNGLT", "none",false,false); parameters.push_back(pgroup);
                CommandParameter plist("list", "InputTypes", "", "", "none", "FNGLT", "none",false,false); parameters.push_back(plist);
                CommandParameter ptaxonomy("taxonomy", "InputTypes", "", "", "none", "FNGLT", "none",false,false); parameters.push_back(ptaxonomy);
                CommandParameter palignreport("alignreport", "InputTypes", "", "", "none", "FNGLT", "none",false,false); parameters.push_back(palignreport);
@@ -40,7 +42,7 @@ vector<string> GetSeqsCommand::setParameters(){
 string GetSeqsCommand::getHelpString(){        
        try {
                string helpString = "";
-               helpString += "The get.seqs command reads an .accnos file and any of the following file types: fasta, name, group, list, taxonomy, quality or alignreport file.\n";
+               helpString += "The get.seqs command reads an .accnos file and any of the following file types: fasta, name, group, count, list, taxonomy, quality or alignreport file.\n";
                helpString += "It outputs a file containing only the sequences in the .accnos file.\n";
                helpString += "The get.seqs command parameters are accnos, fasta, name, group, list, taxonomy, qfile, alignreport and dups.  You must provide accnos unless you have a valid current accnos file, and at least one of the other parameters.\n";
                helpString += "The dups parameter allows you to add the entire line from a name file if you add any name from the line. default=false. \n";
@@ -68,6 +70,7 @@ GetSeqsCommand::GetSeqsCommand(){
                outputTypes["alignreport"] = tempOutNames;
                outputTypes["list"] = tempOutNames;
                outputTypes["qfile"] = tempOutNames;
+        outputTypes["count"] = tempOutNames;
                outputTypes["accnosreport"] = tempOutNames;
        }
        catch(exception& e) {
@@ -76,6 +79,34 @@ GetSeqsCommand::GetSeqsCommand(){
        }
 }
 //**********************************************************************************************************************
+string GetSeqsCommand::getOutputFileNameTag(string type, string inputName=""){ 
+       try {
+        string outputFileName = "";
+               map<string, vector<string> >::iterator it;
+        
+        //is this a type this command creates
+        it = outputTypes.find(type);
+        if (it == outputTypes.end()) {  m->mothurOut("[ERROR]: this command doesn't create a " + type + " output file.\n"); }
+        else {
+            if (type == "fasta")            {   outputFileName =  "pick" + m->getExtension(inputName);   }
+            else if (type == "taxonomy")    {   outputFileName =  "pick" + m->getExtension(inputName);   }
+            else if (type == "name")        {   outputFileName =  "pick" + m->getExtension(inputName);   }
+            else if (type == "count")       {   outputFileName =  "pick.count_table";   }
+            else if (type == "group")       {   outputFileName =  "pick" + m->getExtension(inputName);   }
+            else if (type == "list")        {   outputFileName =  "pick" + m->getExtension(inputName);   }
+            else if (type == "qfile")       {   outputFileName =  "pick" + m->getExtension(inputName);   }
+            else if (type == "accnosreport"){   outputFileName =  "accnos.report";                       }
+            else if (type == "alignreport") {   outputFileName =  "pick.align.report";                   }
+            else { m->mothurOut("[ERROR]: No definition for type " + type + " output file tag.\n"); m->control_pressed = true;  }
+        }
+        return outputFileName;
+       }
+       catch(exception& e) {
+               m->errorOut(e, "GetSeqsCommand", "getOutputFileNameTag");
+               exit(1);
+       }
+}
+//**********************************************************************************************************************
 GetSeqsCommand::GetSeqsCommand(string option)  {
        try {
                abort = false; calledHelp = false;   
@@ -108,6 +139,7 @@ GetSeqsCommand::GetSeqsCommand(string option)  {
                        outputTypes["list"] = tempOutNames;
                        outputTypes["qfile"] = tempOutNames;
                        outputTypes["accnosreport"] = tempOutNames;
+            outputTypes["count"] = tempOutNames;
                        
                        //if the user changes the output directory command factory will send this info to us in the output parameter 
                        outputDir = validParameter.validFile(parameters, "outputdir", false);           if (outputDir == "not found"){  outputDir = "";         }
@@ -188,6 +220,14 @@ GetSeqsCommand::GetSeqsCommand(string option)  {
                                        //if the user has not given a path then, add inputdir. else leave path alone.
                                        if (path == "") {       parameters["qfile"] = inputDir + it->second;            }
                                }
+                
+                it = parameters.find("count");
+                               //user has given a template file
+                               if(it != parameters.end()){ 
+                                       path = m->hasPath(it->second);
+                                       //if the user has not given a path then, add inputdir. else leave path alone.
+                                       if (path == "") {       parameters["count"] = inputDir + it->second;            }
+                               }
                        }
 
                        
@@ -243,17 +283,32 @@ GetSeqsCommand::GetSeqsCommand(string option)  {
                        if (accnosfile2 == "not open") { abort = true; }
                        else if (accnosfile2 == "not found") {  accnosfile2 = "";  }
                        
+            countfile = validParameter.validFile(parameters, "count", true);
+            if (countfile == "not open") { countfile = ""; abort = true; }
+            else if (countfile == "not found") { countfile = "";  }    
+            else { m->setCountTableFile(countfile); }
+            
+            if ((namefile != "") && (countfile != "")) {
+                m->mothurOut("[ERROR]: you may only use one of the following: name or count."); m->mothurOutEndLine(); abort = true;
+            }
+            
+            if ((groupfile != "") && (countfile != "")) {
+                m->mothurOut("[ERROR]: you may only use one of the following: group or count."); m->mothurOutEndLine(); abort=true;
+            }
+
                        
                        string usedDups = "true";
                        string temp = validParameter.validFile(parameters, "dups", false);      if (temp == "not found") { temp = "true"; usedDups = ""; }
                        dups = m->isTrue(temp);
                        
-                       if ((fastafile == "") && (namefile == "") && (groupfile == "") && (alignfile == "") && (listfile == "") && (taxfile == "") && (qualfile == "") && (accnosfile2 == ""))  { m->mothurOut("You must provide one of the following: fasta, name, group, alignreport, taxonomy, quality or listfile."); m->mothurOutEndLine(); abort = true; }
-               
-                       if ((namefile == "") && ((fastafile != "") || (taxfile != ""))){
-                               vector<string> files; files.push_back(fastafile); files.push_back(taxfile);
-                               parser.getNameFile(files);
-                       }
+                       if ((fastafile == "") && (namefile == "") && (groupfile == "") && (alignfile == "") && (listfile == "") && (taxfile == "") && (qualfile == "") && (accnosfile2 == "") && (countfile == ""))  { m->mothurOut("You must provide one of the following: fasta, name, group, count, alignreport, taxonomy, quality or listfile."); m->mothurOutEndLine(); abort = true; }
+            
+            if (countfile == "") {
+                if ((namefile == "") && ((fastafile != "") || (taxfile != ""))){
+                    vector<string> files; files.push_back(fastafile); files.push_back(taxfile);
+                    parser.getNameFile(files);
+                }
+            }
                }
 
        }
@@ -270,14 +325,21 @@ int GetSeqsCommand::execute(){
                if (abort == true) { if (calledHelp) { return 0; }  return 2;   }
                
                //get names you want to keep
-               readAccnos();
+               names = m->readAccnos(accnosfile);
                
                if (m->control_pressed) { return 0; }
+        
+        if (countfile != "") {
+            if ((fastafile != "") || (listfile != "") || (taxfile != "")) { 
+                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");
+            }
+        }
                
                //read through the correct file and output lines you want to keep
                if (namefile != "")                     {               readName();                     }
                if (fastafile != "")            {               readFasta();            }
                if (groupfile != "")            {               readGroup();            }
+        if (countfile != "")           {               readCount();            }
                if (alignfile != "")            {               readAlign();            }
                if (listfile != "")                     {               readList();                     }
                if (taxfile != "")                      {               readTax();                      }
@@ -327,6 +389,10 @@ int GetSeqsCommand::execute(){
                                if ((itTypes->second).size() != 0) { current = (itTypes->second)[0]; m->setQualFile(current); }
                        }
                        
+            itTypes = outputTypes.find("count");
+                       if (itTypes != outputTypes.end()) {
+                               if ((itTypes->second).size() != 0) { current = (itTypes->second)[0]; m->setCountTableFile(current); }
+                       }
                }
                
                return 0;               
@@ -343,7 +409,7 @@ int GetSeqsCommand::readFasta(){
        try {
                string thisOutputDir = outputDir;
                if (outputDir == "") {  thisOutputDir += m->hasPath(fastafile);  }
-               string outputFileName = thisOutputDir + m->getRootName(m->getSimpleName(fastafile)) + "pick" +  m->getExtension(fastafile);
+               string outputFileName = thisOutputDir + m->getRootName(m->getSimpleName(fastafile)) + getOutputFileNameTag("fasta", fastafile);
                ofstream out;
                m->openOutputFile(outputFileName, out);
                
@@ -399,7 +465,7 @@ int GetSeqsCommand::readQual(){
        try {
                string thisOutputDir = outputDir;
                if (outputDir == "") {  thisOutputDir += m->hasPath(qualfile);  }
-               string outputFileName = thisOutputDir + m->getRootName(m->getSimpleName(qualfile)) + "pick" +  m->getExtension(qualfile);
+               string outputFileName = thisOutputDir + m->getRootName(m->getSimpleName(qualfile)) + getOutputFileNameTag("qfile", qualfile);
                ofstream out;
                m->openOutputFile(outputFileName, out);
                
@@ -466,12 +532,70 @@ int GetSeqsCommand::readQual(){
                exit(1);
        }
 }
+//**********************************************************************************************************************
+int GetSeqsCommand::readCount(){
+       try {
+               string thisOutputDir = outputDir;
+               if (outputDir == "") {  thisOutputDir += m->hasPath(countfile);  }
+               string outputFileName = thisOutputDir + m->getRootName(m->getSimpleName(countfile)) + getOutputFileNameTag("count", countfile);
+               
+               ofstream out;
+               m->openOutputFile(outputFileName, out);
+               
+               ifstream in;
+               m->openInputFile(countfile, in);
+               
+               bool wroteSomething = false;
+               int selectedCount = 0;
+               
+        string headers = m->getline(in); m->gobble(in);
+        out << headers << endl;
+        
+        string name, rest; int thisTotal;
+        while (!in.eof()) {
+            
+            if (m->control_pressed) { in.close();  out.close();  m->mothurRemove(outputFileName);  return 0; }
+            
+            in >> name; m->gobble(in); 
+            in >> thisTotal; m->gobble(in);
+            rest = m->getline(in); m->gobble(in);
+            if (m->debug) { m->mothurOut("[DEBUG]: " + name + '\t' + rest + "\n"); }
+            
+            if (names.count(name) != 0) {
+                out << name << '\t' << thisTotal << '\t' << rest << endl;
+                wroteSomething = true;
+                selectedCount+= thisTotal;
+            }
+        }
+        in.close();
+               out.close();
+        
+        //check for groups that have been eliminated
+        CountTable ct;
+        if (ct.testGroups(outputFileName)) {
+            ct.readTable(outputFileName);
+            ct.printTable(outputFileName);
+        }
+               
+               if (wroteSomething == false) {  m->mothurOut("Your file does not contain any sequence from the .accnos file."); m->mothurOutEndLine();  }
+               outputTypes["count"].push_back(outputFileName); outputNames.push_back(outputFileName);
+               
+               m->mothurOut("Selected " + toString(selectedCount) + " sequences from your count file."); m->mothurOutEndLine();
+        
+               return 0;
+       }
+       catch(exception& e) {
+               m->errorOut(e, "GetSeqsCommand", "readCount");
+               exit(1);
+       }
+}
+
 //**********************************************************************************************************************
 int GetSeqsCommand::readList(){
        try {
                string thisOutputDir = outputDir;
                if (outputDir == "") {  thisOutputDir += m->hasPath(listfile);  }
-               string outputFileName = thisOutputDir + m->getRootName(m->getSimpleName(listfile)) + "pick" +  m->getExtension(listfile);
+               string outputFileName = thisOutputDir + m->getRootName(m->getSimpleName(listfile)) + getOutputFileNameTag("list", listfile);
                ofstream out;
                m->openOutputFile(outputFileName, out);
                
@@ -501,19 +625,16 @@ int GetSeqsCommand::readList(){
                        
                                //parse out names that are in accnos file
                                string binnames = list.get(i);
+                vector<string> bnames;
+                m->splitAtComma(binnames, bnames);
                                
                                string newNames = "";
-                               while (binnames.find_first_of(',') != -1) { 
-                                       string name = binnames.substr(0,binnames.find_first_of(','));
-                                       binnames = binnames.substr(binnames.find_first_of(',')+1, binnames.length());
-                                       
+                for (int i = 0; i < bnames.size(); i++) {
+                                       string name = bnames[i];
                                        //if that name is in the .accnos file, add it
                                        if (names.count(name) != 0) {  newNames += name + ",";  selectedCount++; if (m->debug) { sanity["list"].insert(name); } }
                                }
                        
-                               //get last name
-                               if (names.count(binnames) != 0) {  newNames += binnames + ",";  selectedCount++;  if (m->debug) { sanity["list"].insert(binnames); } }
-
                                //if there are names in this bin add to new list
                                if (newNames != "") { 
                                        newNames = newNames.substr(0, newNames.length()-1); //rip off extra comma
@@ -550,7 +671,7 @@ int GetSeqsCommand::readName(){
        try {
                string thisOutputDir = outputDir;
                if (outputDir == "") {  thisOutputDir += m->hasPath(namefile);  }
-               string outputFileName = thisOutputDir + m->getRootName(m->getSimpleName(namefile)) + "pick" +  m->getExtension(namefile);
+               string outputFileName = thisOutputDir + m->getRootName(m->getSimpleName(namefile)) + getOutputFileNameTag("name", namefile);
                ofstream out;
                m->openOutputFile(outputFileName, out);
                
@@ -650,7 +771,7 @@ int GetSeqsCommand::readGroup(){
        try {
                string thisOutputDir = outputDir;
                if (outputDir == "") {  thisOutputDir += m->hasPath(groupfile);  }
-               string outputFileName = thisOutputDir + m->getRootName(m->getSimpleName(groupfile)) + "pick" + m->getExtension(groupfile);
+               string outputFileName = thisOutputDir + m->getRootName(m->getSimpleName(groupfile)) + getOutputFileNameTag("group", groupfile);
                ofstream out;
                m->openOutputFile(outputFileName, out);
                
@@ -706,7 +827,7 @@ int GetSeqsCommand::readTax(){
        try {
                string thisOutputDir = outputDir;
                if (outputDir == "") {  thisOutputDir += m->hasPath(taxfile);  }
-               string outputFileName = thisOutputDir + m->getRootName(m->getSimpleName(taxfile)) + "pick" + m->getExtension(taxfile);
+               string outputFileName = thisOutputDir + m->getRootName(m->getSimpleName(taxfile)) + getOutputFileNameTag("taxonomy", taxfile);
                ofstream out;
                m->openOutputFile(outputFileName, out);
                
@@ -760,7 +881,7 @@ int GetSeqsCommand::readAlign(){
        try {
                string thisOutputDir = outputDir;
                if (outputDir == "") {  thisOutputDir += m->hasPath(alignfile);  }
-               string outputFileName = thisOutputDir + m->getRootName(m->getSimpleName(alignfile)) + "pick.align.report";
+               string outputFileName = thisOutputDir + m->getRootName(m->getSimpleName(alignfile)) + getOutputFileNameTag("alignreport");
                ofstream out;
                m->openOutputFile(outputFileName, out);
                
@@ -827,32 +948,6 @@ int GetSeqsCommand::readAlign(){
        }
 }
 //**********************************************************************************************************************
-
-int GetSeqsCommand::readAccnos(){
-       try {
-               
-               ifstream in;
-               m->openInputFile(accnosfile, in);
-               string name;
-               
-               while(!in.eof()){
-                       in >> name;
-                                               
-                       names.insert(name);
-                       
-                       m->gobble(in);
-               }
-               in.close();     
-               
-               return 0;
-
-       }
-       catch(exception& e) {
-               m->errorOut(e, "GetSeqsCommand", "readAccnos");
-               exit(1);
-       }
-}
-//**********************************************************************************************************************
 //just looking at common mistakes. 
 int GetSeqsCommand::runSanityCheck(){
        try {
@@ -952,7 +1047,7 @@ int GetSeqsCommand::compareAccnos(){
                
                string thisOutputDir = outputDir;
                if (outputDir == "") {  thisOutputDir += m->hasPath(accnosfile);  }
-               string outputFileName = thisOutputDir + m->getRootName(m->getSimpleName(accnosfile)) + "accnos.report";
+               string outputFileName = thisOutputDir + m->getRootName(m->getSimpleName(accnosfile)) + getOutputFileNameTag("accnosreport");
                ofstream out;
                m->openOutputFile(outputFileName, out);
                
@@ -1045,7 +1140,7 @@ int GetSeqsCommand::compareAccnos(){
                
        }
        catch(exception& e) {
-               m->errorOut(e, "GetSeqsCommand", "readAccnos");
+               m->errorOut(e, "GetSeqsCommand", "compareAccnos");
                exit(1);
        }
 }