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;
94 Sequence seq(in); gobble(in);
95 if (seq.getName() != "") {
96 statData data = getStats(seq.getAligned());
98 out << seq.getName() << '\t' << data.pound << '\t' << data.dash << '\t' << data.plus << '\t' << data.equal << '\t';
99 out << data.loop << '\t' << data.tilde << '\t' << data.total << endl;
109 catch(exception& e) {
110 errorOut(e, "AlignCheckCommand", "execute");
114 //**********************************************************************************************************************
115 void AlignCheckCommand::readMap(){
118 structMap.resize(1, 0);
121 openInputFile(mapfile, in);
126 structMap.push_back(position);
131 seqLength = structMap.size();
134 //check you make sure is structMap[10] = 380 then structMap[380] = 10.
135 for(int i=0;i<seqLength;i++){
136 if(structMap[i] != 0){
137 if(structMap[structMap[i]] != i){
138 mothurOut("Your map file contains an error: line " + toString(i) + " does not match line " + toString(structMap[i]) + "."); mothurOutEndLine();
145 catch(exception& e) {
146 errorOut(e, "AlignCheckCommand", "readMap");
150 /**************************************************************************************************/
152 statData AlignCheckCommand::getStats(string sequence){
156 sequence = "*" + sequence; // need to pad the sequence so we can index it by 1
158 int seqLength = sequence.length();
159 for(int i=1;i<seqLength;i++){
160 if(structMap[i] != 0){
161 if(sequence[i] == 'A'){
162 if(sequence[structMap[i]] == 'T') { data.tilde++; }
163 else if(sequence[structMap[i]] == 'A') { data.pound++; }
164 else if(sequence[structMap[i]] == 'G') { data.equal++; }
165 else if(sequence[structMap[i]] == 'C') { data.pound++; }
166 else if(sequence[structMap[i]] == '-') { data.pound++; }
169 else if(sequence[i] == 'T'){
170 if(sequence[structMap[i]] == 'T') { data.plus++; }
171 else if(sequence[structMap[i]] == 'A') { data.tilde++; }
172 else if(sequence[structMap[i]] == 'G') { data.dash++; }
173 else if(sequence[structMap[i]] == 'C') { data.pound++; }
174 else if(sequence[structMap[i]] == '-') { data.pound++; }
177 else if(sequence[i] == 'G'){
178 if(sequence[structMap[i]] == 'T') { data.dash++; }
179 else if(sequence[structMap[i]] == 'A') { data.equal++; }
180 else if(sequence[structMap[i]] == 'G') { data.pound++; }
181 else if(sequence[structMap[i]] == 'C') { data.tilde++; }
182 else if(sequence[structMap[i]] == '-') { data.pound++; }
185 else if(sequence[i] == 'C'){
186 if(sequence[structMap[i]] == 'T') { data.pound++; }
187 else if(sequence[structMap[i]] == 'A') { data.pound++; }
188 else if(sequence[structMap[i]] == 'G') { data.tilde++; }
189 else if(sequence[structMap[i]] == 'C') { data.pound++; }
190 else if(sequence[structMap[i]] == '-') { data.pound++; }
193 else if(sequence[i] == '-'){
194 if(sequence[structMap[i]] == 'T') { data.pound++; data.total++; }
195 else if(sequence[structMap[i]] == 'A') { data.pound++; data.total++; }
196 else if(sequence[structMap[i]] == 'G') { data.pound++; data.total++; }
197 else if(sequence[structMap[i]] == 'C') { data.pound++; data.total++; }
198 else if(sequence[structMap[i]] == '-') { /*donothing*/ }
201 else if(isalnum(sequence[i])){
209 catch(exception& e) {
210 errorOut(e, "AlignCheckCommand", "getStats");
216 //**********************************************************************************************************************