X-Git-Url: https://git.donarmstrong.com/?a=blobdiff_plain;f=getlineagecommand.cpp;h=8a6f47c0998202943e398fce09a9ae52edfa5eda;hb=b83fd3896a0575b61dabfb5ea82af11eae83ddce;hp=a1d44857531ddea86b025ac1b6d0d888f21df702;hpb=8bc3e5b38c2317a1715f53be22fa96455868c281;p=mothur.git

diff --git a/getlineagecommand.cpp b/getlineagecommand.cpp
index a1d4485..8a6f47c 100644
--- a/getlineagecommand.cpp
+++ b/getlineagecommand.cpp
@@ -27,6 +27,7 @@ vector<string> GetLineageCommand::getValidParameters(){
 //**********************************************************************************************************************
 GetLineageCommand::GetLineageCommand(){	
 	try {
+		abort = true;
 		//initialize outputTypes
 		vector<string> tempOutNames;
 		outputTypes["fasta"] = tempOutNames;
@@ -44,7 +45,7 @@ GetLineageCommand::GetLineageCommand(){
 //**********************************************************************************************************************
 vector<string> GetLineageCommand::getRequiredParameters(){	
 	try {
-		string Array[] =  {"taxonomy"};
+		string Array[] =  {"taxonomy","taxon"};
 		vector<string> myArray (Array, Array+(sizeof(Array)/sizeof(string)));
 		return myArray;
 	}
@@ -181,7 +182,11 @@ GetLineageCommand::GetLineageCommand(string option)  {
 			else if (taxfile == "not found") {  taxfile = ""; m->mothurOut("The taxonomy parameter is required for the get.lineage command."); m->mothurOutEndLine();  abort = true; }
 			
 			string usedDups = "true";
-			string temp = validParameter.validFile(parameters, "dups", false);	if (temp == "not found") { temp = "false"; usedDups = ""; }
+			string temp = validParameter.validFile(parameters, "dups", false);	
+			if (temp == "not found") { 
+				if (namefile != "") {  temp = "true";					}
+				else				{  temp = "false"; usedDups = "";	}
+			}
 			dups = m->isTrue(temp);
 			
 			taxons = validParameter.validFile(parameters, "taxon", false);	
@@ -214,6 +219,8 @@ void GetLineageCommand::help(){
 		m->mothurOut("The get.lineage command parameters are taxon, fasta, name, group, list, taxonomy, alignreport and dups.  You must provide taxonomy and taxon.\n");
 		m->mothurOut("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");
 		m->mothurOut("The taxon parameter allows you to select the taxons you would like to get.\n");
+		m->mothurOut("You may enter your taxons with confidence scores, doing so will get only those sequences that belong to the taxonomy and whose cofidence scores is above the scores you give.\n");
+		m->mothurOut("If they belong to the taxonomy and have confidences below those you provide the sequence will not be selected.\n");
 		m->mothurOut("The get.lineage command should be in the following format: get.lineage(taxonomy=yourTaxonomyFile, taxon=yourTaxons).\n");
 		m->mothurOut("Example get.lineage(taxonomy=amazon.silva.taxonomy, taxon=Bacteria;Firmicutes;Bacilli;Lactobacillales;).\n");
 		m->mothurOut("Note: If you are running mothur in script mode you must wrap the taxon in ' characters so mothur will ignore the ; in the taxon.\n");
@@ -535,7 +542,18 @@ int GetLineageCommand::readTax(){
 		m->openInputFile(taxfile, in);
 		string name, tax;
 		
-		bool wroteSomething = false;
+		//bool wroteSomething = false;
+		
+		bool taxonsHasConfidence = false;
+		vector< map<string, int> > searchTaxons;
+		string noConfidenceTaxons = taxons;
+		int hasConPos = taxons.find_first_of('(');
+		if (hasConPos != string::npos) {  
+			taxonsHasConfidence = true; 
+			searchTaxons = getTaxons(taxons); 
+			noConfidenceTaxons = removeConfidences(taxons);
+		}
+		
 		
 		while(!in.eof()){
 
@@ -546,18 +564,89 @@ int GetLineageCommand::readTax(){
 			
 			string newtax = tax;
 			
-			//if the users file contains confidence scores we want to ignore them when searching for the taxons
-			int hasConfidences = tax.find_first_of('(');
-			if (hasConfidences != string::npos) { 
-				newtax = removeConfidences(tax);
+			
+			//if the users file contains confidence scores we want to ignore them when searching for the taxons, unless the taxon has them
+			if (!taxonsHasConfidence) {
+				int hasConfidences = tax.find_first_of('(');
+				if (hasConfidences != string::npos) { 
+					newtax = removeConfidences(tax);
+				}
+				
+				int pos = newtax.find(taxons);
+				
+				if (pos != string::npos) { //this sequence contains the taxon the user wants
+					names.insert(name);
+					out << name << '\t' << tax << endl;
+				}
+				
+			}else{//if taxons has them and you don't them remove taxons
+				int hasConfidences = tax.find_first_of('(');
+				if (hasConfidences == string::npos) { 
+					
+					int pos = newtax.find(noConfidenceTaxons);
+					
+					if (pos != string::npos) { //this sequence contains the taxon the user wants
+						names.insert(name);
+						out << name << '\t' << tax << endl;
+					}
+				}else { //both have confidences so we want to make sure the users confidences are greater then or equal to the taxons
+					//first remove confidences from both and see if the taxonomy exists
+					
+					string noNewTax = tax;
+					int hasConfidences = tax.find_first_of('(');
+					if (hasConfidences != string::npos) { 
+						noNewTax = removeConfidences(tax);
+					}
+					
+					int pos = noNewTax.find(noConfidenceTaxons);
+					
+					if (pos != string::npos) { //if yes, then are the confidences okay
+						
+						bool good = true;
+						vector< map<string, int> > usersTaxon = getTaxons(newtax);
+						
+						//the usersTaxon is most likely longer than the searchTaxons, and searchTaxon[0] may relate to userTaxon[4]
+						//we want to "line them up", so we will find the the index where the searchstring starts
+						int index = 0;
+						for (int i = 0; i < usersTaxon.size(); i++) {
+							
+							if (usersTaxon[i].begin()->first == searchTaxons[0].begin()->first) { 
+								index = i;  
+								int spot = 0;
+								bool goodspot = true;
+								//is this really the start, or are we dealing with a taxon of the same name?
+								while ((spot < searchTaxons.size()) && ((i+spot) < usersTaxon.size())) {
+									if (usersTaxon[i+spot].begin()->first != searchTaxons[spot].begin()->first) { goodspot = false; break; }
+									else { spot++; }
+								}
+								
+								if (goodspot) { break; }
+							}
+						}
+						
+						for (int i = 0; i < searchTaxons.size(); i++) {
+							
+							if ((i+index) < usersTaxon.size()) { //just in case, should never be false
+								if (usersTaxon[i+index].begin()->second < searchTaxons[i].begin()->second) { //is the users cutoff less than the search taxons
+									good = false;
+									break;
+								}
+							}else {
+								good = false;
+								break;
+							}
+						}
+						
+						//passed the test so add you
+						if (good) {
+							names.insert(name);
+							out << name << '\t' << tax << endl;
+						}
+					}
+				}
 			}
 			
-			int pos = newtax.find(taxons);
 			
-			if (pos != string::npos) { //this sequence contains the taxon the user wants
-				names.insert(name);
-				out << name << '\t' << tax << endl;
-			}
 			
 			m->gobble(in);
 		}
@@ -576,6 +665,48 @@ int GetLineageCommand::readTax(){
 	}
 }
 /**************************************************************************************************/
+vector< map<string, int> > GetLineageCommand::getTaxons(string tax) {
+	try {
+		
+		vector< map<string, int> > t;
+		string taxon = "";
+		int taxLength = tax.length();
+		for(int i=0;i<taxLength;i++){
+			if(tax[i] == ';'){
+		
+				int openParen = taxon.find_first_of('(');
+				int closeParen = taxon.find_last_of(')');
+				
+				string newtaxon, confidence;
+				if ((openParen != string::npos) && (closeParen != string::npos)) {
+					newtaxon = taxon.substr(0, openParen); //rip off confidence
+					confidence = taxon.substr((openParen+1), (closeParen-openParen-1));  
+				}else{
+					newtaxon = taxon;
+					confidence = "0";
+				} 
+				int con = 0;
+				convert(confidence, con);
+				
+				map<string, int> temp;
+				temp[newtaxon] = con;
+				t.push_back(temp);
+				
+				taxon = "";
+			}
+			else{
+				taxon += tax[i];
+			}
+		}
+		
+		return t;
+	}
+	catch(exception& e) {
+		m->errorOut(e, "GetLineageCommand", "getTaxons");
+		exit(1);
+	}
+}
+/**************************************************************************************************/
 string GetLineageCommand::removeConfidences(string tax) {
 	try {