2 * secondarystructurecommand.cpp
5 * Created by westcott on 9/18/09.
6 * Copyright 2009 Schloss Lab. All rights reserved.
10 #include "secondarystructurecommand.h"
11 #include "sequence.hpp"
13 //**********************************************************************************************************************
15 AlignCheckCommand::AlignCheckCommand(string option){
19 //allow user to run help
20 if(option == "help") { help(); abort = true; }
23 //valid paramters for this command
24 string Array[] = {"fasta","map"};
25 vector<string> myArray (Array, Array+(sizeof(Array)/sizeof(string)));
27 OptionParser parser(option);
28 map<string,string> parameters = parser.getParameters();
30 ValidParameters validParameter;
32 //check to make sure all parameters are valid for command
33 for (map<string,string>::iterator it = parameters.begin(); it != parameters.end(); it++) {
34 if (validParameter.isValidParameter(it->first, myArray, it->second) != true) { abort = true; }
37 //check for required parameters
38 mapfile = validParameter.validFile(parameters, "map", true);
39 if (mapfile == "not open") { abort = true; }
40 else if (mapfile == "not found") { mapfile = ""; mothurOut("You must provide an map file."); mothurOutEndLine(); abort = true; }
42 fastafile = validParameter.validFile(parameters, "fasta", true);
43 if (fastafile == "not open") { abort = true; }
44 else if (fastafile == "not found") { fastafile = ""; mothurOut("You must provide an fasta file."); mothurOutEndLine(); abort = true; }
50 errorOut(e, "AlignCheckCommand", "RemoveSeqsCommand");
54 //**********************************************************************************************************************
56 void AlignCheckCommand::help(){
58 mothurOut("The align.check command reads a fasta file and map file.\n");
59 mothurOut("It outputs a file containing the secondary structure matches in the .align.check file.\n");
60 mothurOut("The align.check command parameters are fasta and map, both are required.\n");
61 mothurOut("The align.check command should be in the following format: align.check(fasta=yourFasta, map=yourMap).\n");
62 mothurOut("Example align.check(map=silva.ss.map, fasta=amazon.fasta).\n");
63 mothurOut("Note: No spaces between parameter labels (i.e. fasta), '=' and parameters (i.e.yourFasta).\n\n");
66 errorOut(e, "AlignCheckCommand", "help");
71 //**********************************************************************************************************************
73 int AlignCheckCommand::execute(){
76 if (abort == true) { return 0; }
78 //get secondary structure info.
82 openInputFile(fastafile, in);
85 string outfile = getRootName(fastafile) + "align.check";
86 openOutputFile(outfile, out);
88 out << "name" << '\t' << "pound" << '\t' << "dash" << '\t' << "plus" << '\t' << "equal" << '\t';
89 out << "loop" << '\t' << "tilde" << '\t' << "total" << endl;
95 statData data = getStats(seq.getAligned());
97 out << seq.getName() << '\t' << data.pound << '\t' << data.dash << '\t' << data.plus << '\t' << data.equal << '\t';
98 out << data.loop << '\t' << data.tilde << '\t' << data.total << endl;
107 catch(exception& e) {
108 errorOut(e, "AlignCheckCommand", "execute");
112 //**********************************************************************************************************************
113 void AlignCheckCommand::readMap(){
116 structMap.resize(1, 0);
119 openInputFile(mapfile, in);
124 structMap.push_back(position);
129 seqLength = structMap.size();
132 //check you make sure is structMap[10] = 380 then structMap[380] = 10.
133 for(int i=0;i<seqLength;i++){
134 if(structMap[i] != 0){
135 if(structMap[structMap[i]] != i){
136 mothurOut("Your map file contains an error: line " + toString(i) + " does not match line " + toString(structMap[i]) + "."); mothurOutEndLine();
143 catch(exception& e) {
144 errorOut(e, "AlignCheckCommand", "readMap");
148 /**************************************************************************************************/
150 statData AlignCheckCommand::getStats(string sequence){
154 sequence = "*" + sequence; // need to pad the sequence so we can index it by 1
156 int seqLength = sequence.length();
157 for(int i=1;i<seqLength;i++){
158 if(structMap[i] != 0){
159 if(sequence[i] == 'A'){
160 if(sequence[structMap[i]] == 'T') { data.tilde++; }
161 else if(sequence[structMap[i]] == 'A') { data.pound++; }
162 else if(sequence[structMap[i]] == 'G') { data.equal++; }
163 else if(sequence[structMap[i]] == 'C') { data.pound++; }
164 else if(sequence[structMap[i]] == '-') { data.pound++; }
167 else if(sequence[i] == 'T'){
168 if(sequence[structMap[i]] == 'T') { data.plus++; }
169 else if(sequence[structMap[i]] == 'A') { data.tilde++; }
170 else if(sequence[structMap[i]] == 'G') { data.dash++; }
171 else if(sequence[structMap[i]] == 'C') { data.pound++; }
172 else if(sequence[structMap[i]] == '-') { data.pound++; }
175 else if(sequence[i] == 'G'){
176 if(sequence[structMap[i]] == 'T') { data.dash++; }
177 else if(sequence[structMap[i]] == 'A') { data.equal++; }
178 else if(sequence[structMap[i]] == 'G') { data.pound++; }
179 else if(sequence[structMap[i]] == 'C') { data.tilde++; }
180 else if(sequence[structMap[i]] == '-') { data.pound++; }
183 else if(sequence[i] == 'C'){
184 if(sequence[structMap[i]] == 'T') { data.pound++; }
185 else if(sequence[structMap[i]] == 'A') { data.pound++; }
186 else if(sequence[structMap[i]] == 'G') { data.tilde++; }
187 else if(sequence[structMap[i]] == 'C') { data.pound++; }
188 else if(sequence[structMap[i]] == '-') { data.pound++; }
191 else if(sequence[i] == '-'){
192 if(sequence[structMap[i]] == 'T') { data.pound++; data.total++; }
193 else if(sequence[structMap[i]] == 'A') { data.pound++; data.total++; }
194 else if(sequence[structMap[i]] == 'G') { data.pound++; data.total++; }
195 else if(sequence[structMap[i]] == 'C') { data.pound++; data.total++; }
196 else if(sequence[structMap[i]] == '-') { /*donothing*/ }
199 else if(isalnum(sequence[i])){
207 catch(exception& e) {
208 errorOut(e, "AlignCheckCommand", "getStats");
214 //**********************************************************************************************************************