]> git.donarmstrong.com Git - mothur.git/blobdiff - secondarystructurecommand.cpp
fixed mempearson and structpearson
[mothur.git] / secondarystructurecommand.cpp
index 48ea1549aa78537966745dbe11fbc9d1afbf8f51..26bd35444d866132fd37c34fc3e0d8fd4cf84c40 100644 (file)
 #include "secondarystructurecommand.h"
 #include "sequence.hpp"
 
+//**********************************************************************************************************************
+vector<string> AlignCheckCommand::getValidParameters(){        
+       try {
+               string Array[] =  {"fasta","map", "outputdir","inputdir"};
+               vector<string> myArray (Array, Array+(sizeof(Array)/sizeof(string)));
+               return myArray;
+       }
+       catch(exception& e) {
+               m->errorOut(e, "AlignCheckCommand", "getValidParameters");
+               exit(1);
+       }
+}
+//**********************************************************************************************************************
+AlignCheckCommand::AlignCheckCommand(){        
+       try {
+               abort = true;
+               //initialize outputTypes
+               vector<string> tempOutNames;
+               outputTypes["aligncheck"] = tempOutNames;
+       }
+       catch(exception& e) {
+               m->errorOut(e, "AlignCheckCommand", "AlignCheckCommand");
+               exit(1);
+       }
+}
+//**********************************************************************************************************************
+vector<string> AlignCheckCommand::getRequiredParameters(){     
+       try {
+               string Array[] =  {"fasta","map"};
+               vector<string> myArray (Array, Array+(sizeof(Array)/sizeof(string)));
+               return myArray;
+       }
+       catch(exception& e) {
+               m->errorOut(e, "AlignCheckCommand", "getRequiredParameters");
+               exit(1);
+       }
+}
+//**********************************************************************************************************************
+vector<string> AlignCheckCommand::getRequiredFiles(){  
+       try {
+               vector<string> myArray;
+               return myArray;
+       }
+       catch(exception& e) {
+               m->errorOut(e, "AlignCheckCommand", "getRequiredFiles");
+               exit(1);
+       }
+}
 //**********************************************************************************************************************
 
-AlignCheckCommand::AlignCheckCommand(string option){
+AlignCheckCommand::AlignCheckCommand(string option)  {
        try {
                abort = false;
-               
+               haderror = 0;
+                       
                //allow user to run help
                if(option == "help") { help(); abort = true; }
                
                else {
                        //valid paramters for this command
-                       string Array[] =  {"fasta","map"};
+                       string Array[] =  {"fasta","map", "outputdir","inputdir"};
                        vector<string> myArray (Array, Array+(sizeof(Array)/sizeof(string)));
                        
                        OptionParser parser(option);
                        map<string,string> parameters = parser.getParameters();
                        
                        ValidParameters validParameter;
+                       map<string,string>::iterator it;
                        
                        //check to make sure all parameters are valid for command
-                       for (map<string,string>::iterator it = parameters.begin(); it != parameters.end(); it++) { 
+                       for (it = parameters.begin(); it != parameters.end(); it++) { 
                                if (validParameter.isValidParameter(it->first, myArray, it->second) != true) {  abort = true;  }
                        }
                        
+                       //initialize outputTypes
+                       vector<string> tempOutNames;
+                       outputTypes["aligncheck"] = tempOutNames;
+                       
+                       //if the user changes the input directory command factory will send this info to us in the output parameter 
+                       string inputDir = validParameter.validFile(parameters, "inputdir", false);              
+                       if (inputDir == "not found"){   inputDir = "";          }
+                       else {
+                               string path;
+                               it = parameters.find("fasta");
+                               //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["fasta"] = inputDir + it->second;            }
+                               }
+                               
+                               it = parameters.find("map");
+                               //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["map"] = inputDir + it->second;              }
+                               }
+                       }
+
                        //check for required parameters
                        mapfile = validParameter.validFile(parameters, "map", true);
                        if (mapfile == "not open") { abort = true; }
-                       else if (mapfile == "not found") {  mapfile = "";  mothurOut("You must provide an map file."); mothurOutEndLine(); abort = true; }      
+                       else if (mapfile == "not found") {  mapfile = "";  m->mothurOut("You must provide an map file."); m->mothurOutEndLine(); abort = true; }        
                        
                        fastafile = validParameter.validFile(parameters, "fasta", true);
                        if (fastafile == "not open") { abort = true; }
-                       else if (fastafile == "not found") {  fastafile = "";  mothurOut("You must provide an fasta file."); mothurOutEndLine(); abort = true;  }       
+                       else if (fastafile == "not found") {  fastafile = "";  m->mothurOut("You must provide an fasta file."); m->mothurOutEndLine(); abort = true;  } 
                        
+                       //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 = ""; 
+                               outputDir += m->hasPath(fastafile); //if user entered a file with a path then preserve it       
+                       }
+
                }
 
        }
        catch(exception& e) {
-               errorOut(e, "AlignCheckCommand", "RemoveSeqsCommand");
+               m->errorOut(e, "AlignCheckCommand", "RemoveSeqsCommand");
                exit(1);
        }
 }
@@ -55,15 +137,15 @@ AlignCheckCommand::AlignCheckCommand(string option){
 
 void AlignCheckCommand::help(){
        try {
-               //mothurOut("The remove.seqs command reads an .accnos file and one of the following file types: fasta, name, group or alignreport file.\n");
-               //mothurOut("It outputs a file containing the sequences NOT in the .accnos file.\n");
-               //mothurOut("The remove.seqs command parameters are accnos, fasta, name, group and alignreport.  You must provide accnos and one of the other parameters.\n");
-               //mothurOut("The remove.seqs command should be in the following format: remove.seqs(accnos=yourAccnos, fasta=yourFasta).\n");
-               //mothurOut("Example remove.seqs(accnos=amazon.accnos, fasta=amazon.fasta).\n");
-               //mothurOut("Note: No spaces between parameter labels (i.e. fasta), '=' and parameters (i.e.yourFasta).\n\n");
+               m->mothurOut("The align.check command reads a fasta file and map file.\n");
+               m->mothurOut("It outputs a file containing the secondary structure matches in the .align.check file.\n");
+               m->mothurOut("The align.check command parameters are fasta and map, both are required.\n");
+               m->mothurOut("The align.check command should be in the following format: align.check(fasta=yourFasta, map=yourMap).\n");
+               m->mothurOut("Example align.check(map=silva.ss.map, fasta=amazon.fasta).\n");
+               m->mothurOut("Note: No spaces between parameter labels (i.e. fasta), '=' and parameters (i.e.yourFasta).\n\n");
        }
        catch(exception& e) {
-               errorOut(e, "AlignCheckCommand", "help");
+               m->errorOut(e, "AlignCheckCommand", "help");
                exit(1);
        }
 }
@@ -78,13 +160,46 @@ int AlignCheckCommand::execute(){
                //get secondary structure info.
                readMap();
                
-       
+               ifstream in;
+               m->openInputFile(fastafile, in);
+               
+               ofstream out;
+               string outfile = outputDir + m->getRootName(m->getSimpleName(fastafile)) + "align.check";
+               m->openOutputFile(outfile, out);
+               
+               out << "name" << '\t' << "pound" << '\t' << "dash" << '\t' << "plus" << '\t' << "equal" << '\t';
+               out << "loop" << '\t' << "tilde" << '\t' << "total" << endl;
+
+               
+               while(!in.eof()){
+                       if (m->control_pressed) { in.close(); out.close(); remove(outfile.c_str()); return 0; }
+                       
+                       Sequence seq(in);  m->gobble(in);
+                       if (seq.getName() != "") {
+                               statData data = getStats(seq.getAligned());
+                               
+                               if (haderror == 1) { break; }
+                               
+                               out << seq.getName() << '\t' << data.pound << '\t' << data.dash << '\t' << data.plus << '\t' << data.equal << '\t';
+                               out << data.loop << '\t' << data.tilde << '\t' << data.total << endl;
+                       }
+               }
+
+               in.close();
+               out.close();
+               
+               if (m->control_pressed) {  remove(outfile.c_str()); return 0; }
+               
+               m->mothurOutEndLine();
+               m->mothurOut("Output File Name: "); m->mothurOutEndLine();
+               m->mothurOut(outfile); m->mothurOutEndLine();   outputNames.push_back(outfile); outputTypes["aligncheck"].push_back(outfile);
+               m->mothurOutEndLine();
                
                return 0;               
        }
 
        catch(exception& e) {
-               errorOut(e, "AlignCheckCommand", "execute");
+               m->errorOut(e, "AlignCheckCommand", "execute");
                exit(1);
        }
 }
@@ -95,16 +210,16 @@ void AlignCheckCommand::readMap(){
                structMap.resize(1, 0);
                ifstream in;
                
-               openInputFile(mapfile, in);
+               m->openInputFile(mapfile, in);
                
                while(!in.eof()){
                        int position;
                        in >> position;
                        structMap.push_back(position);  
-                       gobble(in);
+                       m->gobble(in);
                }
                in.close();
-               
+
                seqLength = structMap.size();
                
                
@@ -112,7 +227,7 @@ void AlignCheckCommand::readMap(){
                for(int i=0;i<seqLength;i++){
                        if(structMap[i] != 0){
                                if(structMap[structMap[i]] != i){
-                                       mothurOut("Your map file contains an error:  line " + toString(i) + " does not match line " + toString(structMap[i]) + "."); mothurOutEndLine();
+                                       m->mothurOut("Your map file contains an error:  line " + toString(i) + " does not match line " + toString(structMap[i]) + "."); m->mothurOutEndLine();
                                }
                        }
                }
@@ -120,9 +235,77 @@ void AlignCheckCommand::readMap(){
                
        }
        catch(exception& e) {
-               errorOut(e, "AlignCheckCommand", "readFasta");
+               m->errorOut(e, "AlignCheckCommand", "readMap");
                exit(1);
        }
 }
+/**************************************************************************************************/
+
+statData AlignCheckCommand::getStats(string sequence){
+       try {
+       
+               statData data;
+               sequence = "*" + sequence; // need to pad the sequence so we can index it by 1
+               
+               int length = sequence.length();
+               
+               if (length != seqLength) { m->mothurOut("your sequences are " + toString(length) + " long, but your map file only contains " + toString(seqLength) + " entries. please correct."); m->mothurOutEndLine(); haderror = 1; return data;  }
+               
+               for(int i=1;i<length;i++){
+                       if(structMap[i] != 0){
+                               if(sequence[i] == 'A'){
+                                       if(sequence[structMap[i]] == 'T')               {       data.tilde++;   }
+                                       else if(sequence[structMap[i]] == 'A')  {       data.pound++;   }
+                                       else if(sequence[structMap[i]] == 'G')  {       data.equal++;   }
+                                       else if(sequence[structMap[i]] == 'C')  {       data.pound++;   }
+                                       else if(sequence[structMap[i]] == '-')  {       data.pound++;   }
+                                       data.total++;
+                               }
+                               else if(sequence[i] == 'T'){
+                                       if(sequence[structMap[i]] == 'T')               {       data.plus++;    }
+                                       else if(sequence[structMap[i]] == 'A')  {       data.tilde++;   }
+                                       else if(sequence[structMap[i]] == 'G')  {       data.dash++;    }
+                                       else if(sequence[structMap[i]] == 'C')  {       data.pound++;   }
+                                       else if(sequence[structMap[i]] == '-')  {       data.pound++;   }
+                                       data.total++;
+                               }
+                               else if(sequence[i] == 'G'){
+                                       if(sequence[structMap[i]] == 'T')               {       data.dash++;    }
+                                       else if(sequence[structMap[i]] == 'A')  {       data.equal++;   }
+                                       else if(sequence[structMap[i]] == 'G')  {       data.pound++;   }
+                                       else if(sequence[structMap[i]] == 'C')  {       data.tilde++;   }
+                                       else if(sequence[structMap[i]] == '-')  {       data.pound++;   }
+                                       data.total++;
+                               }
+                               else if(sequence[i] == 'C'){
+                                       if(sequence[structMap[i]] == 'T')               {       data.pound++;   }
+                                       else if(sequence[structMap[i]] == 'A')  {       data.pound++;   }
+                                       else if(sequence[structMap[i]] == 'G')  {       data.tilde++;   }
+                                       else if(sequence[structMap[i]] == 'C')  {       data.pound++;   }
+                                       else if(sequence[structMap[i]] == '-')  {       data.pound++;   }
+                                       data.total++;
+                               }
+                               else if(sequence[i] == '-'){
+                                       if(sequence[structMap[i]] == 'T')               {       data.pound++;   data.total++;   }
+                                       else if(sequence[structMap[i]] == 'A')  {       data.pound++;   data.total++;   }
+                                       else if(sequence[structMap[i]] == 'G')  {       data.pound++;   data.total++;   }
+                                       else if(sequence[structMap[i]] == 'C')  {       data.pound++;   data.total++;   }
+                                       else if(sequence[structMap[i]] == '-')  {               /*donothing*/                           }
+                               }                       
+                       }
+                       else if(isalnum(sequence[i])){
+                               data.loop++;
+                               data.total++;
+                       }
+               }
+               return data;
+               
+       }
+       catch(exception& e) {
+               m->errorOut(e, "AlignCheckCommand", "getStats");
+               exit(1);
+       }
+}
+
 
 //**********************************************************************************************************************