]> git.donarmstrong.com Git - mothur.git/blobdiff - getoturepcommand.cpp
added get.rabund and get.sabund command and fixed bug introduced by line by line...
[mothur.git] / getoturepcommand.cpp
index 064455f49e28538c95858ba574dc9ef3b0cc9289..eb88fbf5ab807b451b00025d55a1ce719a6c6ccd 100644 (file)
@@ -20,15 +20,31 @@ GetOTURepCommand::GetOTURepCommand(){
                if(globaldata->gListVector != NULL)             {       
                        listOfNames = new ListVector(*globaldata->gListVector); 
                        
+                       vector<string> names;
+                       string binnames;
                        //map names to rows in sparsematrix
                        for (int i = 0; i < listOfNames->size(); i++) {
-                               nameToIndex[listOfNames->get(i)] = i;
+                               names.clear();
+                               binnames = listOfNames->get(i);
+                               splitAtComma(binnames, names);
+                               
+                               for (int j = 0; j < names.size(); j++) {
+                                       nameToIndex[names[j]] = i;
+                               }
                        }
-               }else { cout << "error" << endl; }
+               }else { cout << "error, no listvector." << endl; }
 
                
                fastafile = globaldata->getFastaFile();
                namesfile = globaldata->getNameFile();
+               groupfile = globaldata->getGroupFile();
+               
+               if (groupfile != "") {
+                       //read in group map info.
+                       groupMap = new GroupMap(groupfile);
+                       groupMap->readMap();
+               }
+
                openInputFile(fastafile, in);
                
                fasta = new FastaMap();
@@ -46,11 +62,12 @@ GetOTURepCommand::GetOTURepCommand(){
 //**********************************************************************************************************************
 
 GetOTURepCommand::~GetOTURepCommand(){
-       delete matrix;
-       delete list;
        delete input;
        delete read;
        delete fasta;
+       if (groupfile != "") {
+               delete groupMap;
+       }
 }
 
 //**********************************************************************************************************************
@@ -58,7 +75,7 @@ GetOTURepCommand::~GetOTURepCommand(){
 int GetOTURepCommand::execute(){
        try {
                int count = 1;
-               string nameRep, name, sequence;
+               int error;
                
                //read fastafile
                fasta->readFastaFile(in);
@@ -72,48 +89,73 @@ int GetOTURepCommand::execute(){
                }
                
                //read list file
-               read = new ReadPhilFile(globaldata->getListFile());     
+               read = new ReadOTUFile(globaldata->getListFile());      
                read->read(&*globaldata); 
                
                input = globaldata->ginput;
                list = globaldata->gListVector;
+               ListVector* lastList = list;
+               
+               //if the users enters label "0.06" and there is no "0.06" in their file use the next lowest label.
+               set<string> processedLabels;
+               set<string> userLabels = globaldata->labels;
+               set<int> userLines = globaldata->lines;
+
                
-               while(list != NULL){
+               while((list != NULL) && ((globaldata->allLines == 1) || (userLabels.size() != 0) || (userLines.size() != 0))) {
                        
                        if(globaldata->allLines == 1 || globaldata->lines.count(count) == 1 || globaldata->labels.count(list->getLabel()) == 1){
-                               
-                               //create output file
-                               string outputFileName = getRootName(globaldata->getListFile()) + list->getLabel() + ".rep.fasta";
-                               openOutputFile(outputFileName, out);
-
-                               cout << list->getLabel() << '\t' << count << endl;
-                               
-                               //for each bin in the list vector
-                               for (int i = 0; i < list->size(); i++) {
-                                       nameRep = FindRep(i);
+                                       cout << list->getLabel() << '\t' << count << endl;
+                                       error = process(list);
+                                       if (error == 1) { return 0; } //there is an error in hte input files, abort command
                                        
-                                       //print out name and sequence for that bin
-                                       sequence = fasta->getSequence(nameRep);
-
-                                       if (sequence != "not found") {
-                                               nameRep = nameRep + "|" + toString(i+1);
-                                               out << ">" << nameRep << endl;
-                                               out << sequence << endl;
-                                       }else { 
-                                               cout << nameRep << " is missing from your fasta or name file. Please correct. " << endl; 
-                                               remove(outputFileName.c_str());
-                                               return 0;
-                                       }
-                               }
-                               
-                               out.close();
+                                       processedLabels.insert(list->getLabel());
+                                       userLabels.erase(list->getLabel());
+                                       userLines.erase(count);
                        }
                        
+                       if ((anyLabelsToProcess(list->getLabel(), userLabels, "") == true) && (processedLabels.count(lastList->getLabel()) != 1)) {
+                                       cout << lastList->getLabel() << '\t' << count << endl;
+                                       error = process(lastList);
+                                       if (error == 1) { return 0; } //there is an error in hte input files, abort command
+                                       
+                                       processedLabels.insert(lastList->getLabel());
+                                       userLabels.erase(lastList->getLabel());
+                       }
+                       
+                       if (count != 1) { delete lastList; }
+                       lastList = list;                        
+                       
                        list = input->getListVector();
                        count++;
                }
-
                
+               //output error messages about any remaining user labels
+               set<string>::iterator it;
+               bool needToRun = false;
+               for (it = userLabels.begin(); it != userLabels.end(); it++) {  
+                       cout << "Your file does not include the label "<< *it; 
+                       if (processedLabels.count(lastList->getLabel()) != 1) {
+                               cout << ". I will use " << lastList->getLabel() << "." << endl;
+                               needToRun = true;
+                       }else {
+                               cout << ". Please refer to " << lastList->getLabel() << "." << endl;
+                       }
+               }
+               
+               //run last line if you need to
+               if (needToRun == true)  {
+                       cout << lastList->getLabel() << '\t' << count << endl;
+                       error = process(lastList);
+                       if (error == 1) { return 0; } //there is an error in hte input files, abort command
+               }
+               delete lastList;
+               
+               delete matrix;
+               globaldata->gSparseMatrix = NULL;
+               delete list;
+               globaldata->gListVector = NULL;
+
                return 0;
        }
        catch(exception& e) {
@@ -165,7 +207,7 @@ void GetOTURepCommand::readNamesFile() {
        }       
 }
 //**********************************************************************************************************************
-string GetOTURepCommand::FindRep(int bin) {
+string GetOTURepCommand::FindRep(int bin, string& group, ListVector* thisList) {
        try{
                vector<string> names;
                map<string, float> sums;
@@ -174,18 +216,41 @@ string GetOTURepCommand::FindRep(int bin) {
                string binnames;
                float min = 10000;
                string minName;
+               map<string, string> groups;
+               map<string, string>::iterator groupIt;
                
-               binnames = list->get(bin);
-               
+               binnames = thisList->get(bin);
+       
                //parse names into vector
                splitAtComma(binnames, names);
                
+               //if you have a groupfile
+               if(groupfile != "") {
+                       //find the groups that are in this bin
+                       for (int i = 0; i < names.size(); i++) {
+                               string groupName = groupMap->getGroup(names[i]);
+                               if (groupName == "not found") {  
+                                       cout << names[i] << " is missing from your group file. Please correct. " << endl;
+                                       groupError = true;
+                               }else{
+                                       groups[groupName] = groupName;
+                               }
+                       }
+                       
+                       //turn the groups into a string
+                       for(groupIt = groups.begin(); groupIt != groups.end(); groupIt++) { group += groupIt->first + "-"; }
+                       
+                       //rip off last dash
+                       group = group.substr(0, group.length()-1);
+               }
+               
                //if only 1 sequence in bin then that's the rep
                if (names.size() == 1) { return names[0]; }
                else {
                        //fill binMap
                        for (int i = 0; i < names.size(); i++) {
                                for (it3 = nameToIndex.begin(); it3 != nameToIndex.end(); it3++) {
+
                                        if (it3->first == names[i]) {  
                                                binMap[it3->second] = it3->first;
 
@@ -233,6 +298,56 @@ string GetOTURepCommand::FindRep(int bin) {
        }       
 }
 
+//**********************************************************************************************************************
+int GetOTURepCommand::process(ListVector* processList) {
+       try{
+                               string nameRep, name, sequence;
+
+                               //create output file
+                               string outputFileName = getRootName(globaldata->getListFile()) + processList->getLabel() + ".rep.fasta";
+                               openOutputFile(outputFileName, out);
+                               
+                               //for each bin in the list vector
+                               for (int i = 0; i < processList->size(); i++) {
+                                       string groups;
+                                       nameRep = FindRep(i, groups, processList);
+                                       
+                                       //print out name and sequence for that bin
+                                       sequence = fasta->getSequence(nameRep);
+
+                                       if (sequence != "not found") {
+                                               if (groupfile == "") {
+                                                       nameRep = nameRep + "|" + toString(i+1);
+                                                       out << ">" << nameRep << endl;
+                                                       out << sequence << endl;
+                                               }else {
+                                                       nameRep = nameRep + "|" + groups + "|" + toString(i+1);
+                                                       out << ">" << nameRep << endl;
+                                                       out << sequence << endl;
+                                               }
+                                       }else { 
+                                               cout << nameRep << " is missing from your fasta or name file. Please correct. " << endl; 
+                                               remove(outputFileName.c_str());
+                                               return 1;
+                                       }
+                               }
+                               
+                               out.close();
+                               return 0;
+       
+       }
+       catch(exception& e) {
+               cout << "Standard Error: " << e.what() << " has occurred in the GetOTURepCommand class Function process. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";
+               exit(1);
+       }
+       catch(...) {
+               cout << "An unknown error has occurred in the GetOTURepCommand class function process. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";
+               exit(1);
+       }       
+}
+
+//**********************************************************************************************************************
+