]> git.donarmstrong.com Git - mothur.git/blobdiff - secondarystructurecommand.cpp
added summary outputs for align.check
[mothur.git] / secondarystructurecommand.cpp
index c9002cab4a458c1ef9fbb3a2acb85706e01513dc..ece2de0fee12364effae1b23065ef30903c22485 100644 (file)
@@ -13,7 +13,7 @@
 //**********************************************************************************************************************
 vector<string> AlignCheckCommand::getValidParameters(){        
        try {
-               string Array[] =  {"fasta","map", "outputdir","inputdir"};
+               string Array[] =  {"fasta", "name","map", "outputdir","inputdir"};
                vector<string> myArray (Array, Array+(sizeof(Array)/sizeof(string)));
                return myArray;
        }
@@ -69,7 +69,7 @@ AlignCheckCommand::AlignCheckCommand(string option)  {
                
                else {
                        //valid paramters for this command
-                       string Array[] =  {"fasta","map", "outputdir","inputdir"};
+                       string Array[] =  {"fasta","name","map", "outputdir","inputdir"};
                        vector<string> myArray (Array, Array+(sizeof(Array)/sizeof(string)));
                        
                        OptionParser parser(option);
@@ -107,6 +107,14 @@ AlignCheckCommand::AlignCheckCommand(string option)  {
                                        //if the user has not given a path then, add inputdir. else leave path alone.
                                        if (path == "") {       parameters["map"] = inputDir + it->second;              }
                                }
+                               
+                               it = parameters.find("name");
+                               //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["name"] = inputDir + it->second;             }
+                               }
                        }
 
                        //check for required parameters
@@ -118,6 +126,10 @@ AlignCheckCommand::AlignCheckCommand(string option)  {
                        if (fastafile == "not open") { abort = true; }
                        else if (fastafile == "not found") {  fastafile = "";  m->mothurOut("You must provide an fasta file."); m->mothurOutEndLine(); abort = true;  } 
                        
+                       namefile = validParameter.validFile(parameters, "name", true);
+                       if (namefile == "not open") { namefile = ""; abort = true; }
+                       else if (namefile == "not found") { namefile = "";  }   
+                       
                        //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 = ""; 
@@ -159,6 +171,10 @@ int AlignCheckCommand::execute(){
                //get secondary structure info.
                readMap();
                
+               if (namefile != "") { nameMap = m->readNames(namefile); }
+               
+               if (m->control_pressed) { return 0; }
+               
                ifstream in;
                m->openInputFile(fastafile, in);
                
@@ -166,10 +182,19 @@ int AlignCheckCommand::execute(){
                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;
+               out << "loop" << '\t' << "tilde" << '\t' << "total"  << '\t' << "numseqs" << endl;
 
+               vector<int> pound;
+               vector<int> dash;
+               vector<int> plus;
+               vector<int> equal;
+               vector<int> loop;
+               vector<int> tilde;
+               vector<int> total;
                
+               int count = 0;
                while(!in.eof()){
                        if (m->control_pressed) { in.close(); out.close(); remove(outfile.c_str()); return 0; }
                        
@@ -177,10 +202,31 @@ int AlignCheckCommand::execute(){
                        if (seq.getName() != "") {
                                statData data = getStats(seq.getAligned());
                                
-                               if (haderror == 1) { break; }
+                               if (haderror == 1) { m->control_pressed = true; break; }
+                               
+                               int num = 1;
+                               if (namefile != "") {
+                                       //make sure this sequence is in the namefile, else error 
+                                       map<string, int>::iterator it = nameMap.find(seq.getName());
+                                       
+                                       if (it == nameMap.end()) { cout << "[ERROR]: " << seq.getName() << " is not in your namefile, please correct." << endl; m->control_pressed = true; }
+                                       else { num = it->second; }
+                               }
+                               
+                               //for each sequence this sequence represents
+                               for (int i = 0; i < num; i++) {
+                                       pound.push_back(data.pound);
+                                       dash.push_back(data.dash);
+                                       plus.push_back(data.plus);
+                                       equal.push_back(data.equal);
+                                       loop.push_back(data.loop);
+                                       tilde.push_back(data.tilde);
+                                       total.push_back(data.total);
+                               }       
+                               count++;
                                
                                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;
+                               out << data.loop << '\t' << data.tilde << '\t' << data.total << '\t' << num << endl;
                        }
                }
 
@@ -189,6 +235,37 @@ int AlignCheckCommand::execute(){
                
                if (m->control_pressed) {  remove(outfile.c_str()); return 0; }
                
+               sort(pound.begin(), pound.end());
+               sort(dash.begin(), dash.end());
+               sort(plus.begin(), plus.end());
+               sort(equal.begin(), equal.end());
+               sort(loop.begin(), loop.end());
+               sort(tilde.begin(), tilde.end());
+               sort(total.begin(), total.end());
+               int size = pound.size();
+               
+               int ptile0_25   = int(size * 0.025);
+               int ptile25             = int(size * 0.250);
+               int ptile50             = int(size * 0.500);
+               int ptile75             = int(size * 0.750);
+               int ptile97_5   = int(size * 0.975);
+               int ptile100    = size - 1;
+               
+               if (m->control_pressed) {  remove(outfile.c_str()); return 0; }
+               
+               m->mothurOutEndLine();
+               m->mothurOut("\t\tPound\tDash\tPlus\tEqual\tLoop\tTilde\tTotal"); m->mothurOutEndLine();
+               m->mothurOut("Minimum:\t" + toString(pound[0]) + "\t" + toString(dash[0]) + "\t" + toString(plus[0]) + "\t" + toString(equal[0]) + "\t" + toString(loop[0]) + "\t" + toString(tilde[0]) + "\t" + toString(total[0])); m->mothurOutEndLine();
+               m->mothurOut("2.5%-tile:\t" + toString(pound[ptile0_25]) + "\t" + toString(dash[ptile0_25]) + "\t" + toString(plus[ptile0_25]) + "\t" + toString(equal[ptile0_25]) + "\t"+ toString(loop[ptile0_25]) + "\t"+ toString(tilde[ptile0_25]) + "\t"+ toString(total[ptile0_25])); m->mothurOutEndLine();
+               m->mothurOut("25%-tile:\t" + toString(pound[ptile25]) + "\t" + toString(dash[ptile25]) + "\t" + toString(plus[ptile25]) + "\t" + toString(equal[ptile25]) + "\t" + toString(loop[ptile25]) + "\t" + toString(tilde[ptile25]) + "\t" + toString(total[ptile25])); m->mothurOutEndLine();
+               m->mothurOut("Median: \t" + toString(pound[ptile50]) + "\t" + toString(dash[ptile50]) + "\t" + toString(plus[ptile50]) + "\t" + toString(equal[ptile50]) + "\t" + toString(loop[ptile50]) + "\t" + toString(tilde[ptile50]) + "\t" + toString(total[ptile50])); m->mothurOutEndLine();
+               m->mothurOut("75%-tile:\t" + toString(pound[ptile75]) + "\t" + toString(dash[ptile75]) + "\t" + toString(plus[ptile75]) + "\t" + toString(equal[ptile75]) + "\t" + toString(loop[ptile75]) + "\t" + toString(tilde[ptile75]) + "\t" + toString(total[ptile75])); m->mothurOutEndLine();
+               m->mothurOut("97.5%-tile:\t" + toString(pound[ptile97_5]) + "\t" + toString(dash[ptile97_5]) + "\t" + toString(plus[ptile97_5]) + "\t" + toString(equal[ptile97_5]) + "\t" + toString(loop[ptile97_5]) + "\t" + toString(tilde[ptile97_5]) + "\t" + toString(total[ptile97_5])); m->mothurOutEndLine();
+               m->mothurOut("Maximum:\t" + toString(pound[ptile100]) + "\t" + toString(dash[ptile100]) + "\t" + toString(plus[ptile100]) + "\t" + toString(equal[ptile100]) + "\t" + toString(loop[ptile100]) + "\t" + toString(tilde[ptile100]) + "\t" + toString(total[ptile100])); m->mothurOutEndLine();
+               if (namefile == "") {  m->mothurOut("# of Seqs:\t" + toString(count)); m->mothurOutEndLine(); }
+               else { m->mothurOut("# of unique seqs:\t" + toString(count)); m->mothurOutEndLine(); m->mothurOut("total # of seqs:\t" + toString(size)); m->mothurOutEndLine(); }
+               
+               
                m->mothurOutEndLine();
                m->mothurOut("Output File Name: "); m->mothurOutEndLine();
                m->mothurOut(outfile); m->mothurOutEndLine();   outputNames.push_back(outfile); outputTypes["aligncheck"].push_back(outfile);
@@ -305,6 +382,4 @@ statData AlignCheckCommand::getStats(string sequence){
                exit(1);
        }
 }
-
-
 //**********************************************************************************************************************