1 #ifndef SUMMARYQUALCOMMAND_H
2 #define SUMMARYQUALCOMMAND_H
8 * Created by westcott on 11/28/11.
9 * Copyright 2011 Schloss Lab. All rights reserved.
14 #include "command.hpp"
15 #include "qualityscores.h"
17 /**************************************************************************************************/
19 class SummaryQualCommand : public Command {
21 SummaryQualCommand(string);
23 ~SummaryQualCommand(){}
25 vector<string> setParameters();
26 string getCommandName() { return "summary.qual"; }
27 string getCommandCategory() { return "Sequence Processing"; }
29 string getHelpString();
30 string getOutputPattern(string);
31 string getCitation() { return "http://www.mothur.org/wiki/Summary.qual"; }
32 string getDescription() { return "summarize the quality of a set of sequences"; }
35 void help() { m->mothurOut(getHelpString()); }
39 string qualfile, outputDir, namefile, countfile;
40 vector<string> outputNames;
41 map<string, int> nameMap;
45 unsigned long long start;
46 unsigned long long end;
47 linePair(unsigned long long i, unsigned long long j) : start(i), end(j) {}
50 vector<linePair> lines;
51 vector<int> processIDS;
53 int createProcessesCreateSummary(vector<int>&, vector<int>&, vector< vector<int> >&, string);
54 int driverCreateSummary(vector<int>&, vector<int>&, vector< vector<int> >&, string, linePair);
55 int printQual(string, vector<int>&, vector<int>&, vector< vector<int> >&);
58 /**************************************************************************************************/
59 //custom data structure for threads to use.
60 // This is passed by void pointer so it can be any data type
61 // that can be passed using a single void pointer (LPVOID).
62 struct seqSumQualData {
65 vector< vector<int> > scores;
67 unsigned long long start;
68 unsigned long long end;
72 map<string, int> nameMap;
75 seqSumQualData(string f, MothurOut* mout, unsigned long long st, unsigned long long en, bool n, map<string, int> nam) {
86 /**************************************************************************************************/
87 #if defined (__APPLE__) || (__MACH__) || (linux) || (__linux) || (__linux__) || (__unix__) || (__unix)
89 static DWORD WINAPI MySeqSumQualThreadFunction(LPVOID lpParam){
90 seqSumQualData* pDataArray;
91 pDataArray = (seqSumQualData*)lpParam;
95 pDataArray->m->openInputFile(pDataArray->filename, in);
97 //print header if you are process 0
98 if ((pDataArray->start == 0) || (pDataArray->start == 1)) {
100 }else { //this accounts for the difference in line endings.
101 in.seekg(pDataArray->start-1); pDataArray->m->gobble(in);
104 pDataArray->count = 0;
105 pDataArray->numSeqs = 0;
106 for(int i = 0; i < pDataArray->end; i++){ //end is the number of sequences to process
108 if (pDataArray->m->control_pressed) { in.close(); pDataArray->count = 1; return 1; }
110 QualityScores current(in); pDataArray->m->gobble(in);
112 if (current.getName() != "") {
115 if (pDataArray->hasNameMap) {
116 //make sure this sequence is in the namefile, else error
117 map<string, int>::iterator it = pDataArray->nameMap.find(current.getName());
119 if (it == pDataArray->nameMap.end()) { pDataArray->m->mothurOut("[ERROR]: " + current.getName() + " is not in your namefile, please correct."); pDataArray->m->mothurOutEndLine(); pDataArray->m->control_pressed = true; }
120 else { num = it->second; }
123 vector<int> thisScores = current.getQualityScores();
125 //resize to num of positions setting number of seqs with that size to 1
126 if (pDataArray->position.size() < thisScores.size()) { pDataArray->position.resize(thisScores.size(), 0); }
127 if (pDataArray->averageQ.size() < thisScores.size()) { pDataArray->averageQ.resize(thisScores.size(), 0); }
128 if (pDataArray->scores.size() < thisScores.size()) {
129 pDataArray->scores.resize(thisScores.size());
130 for (int i = 0; i < pDataArray->scores.size(); i++) { pDataArray->scores.at(i).resize(41, 0); }
133 //increase counts of number of seqs with this position
134 //average is really the total, we will average in execute
135 for (int i = 0; i < thisScores.size(); i++) {
136 pDataArray->position.at(i) += num;
137 pDataArray->averageQ.at(i) += (thisScores[i] * num); //weighting for namesfile
138 if (thisScores[i] > 40) { pDataArray->m->mothurOut("[ERROR]: " + current.getName() + " has a quality scores of " + toString(thisScores[i]) + ", expecting values to be less than 40."); pDataArray->m->mothurOutEndLine(); pDataArray->m->control_pressed = true; }
139 else { pDataArray->scores.at(i)[thisScores[i]] += num; }
142 pDataArray->numSeqs += num;
152 catch(exception& e) {
153 pDataArray->m->errorOut(e, "SummaryQualCommand", "MySeqSumQualThreadFunction");
160 /**************************************************************************************************/