seqerrorcommand.h

   1 #ifndef SEQERRORCOMMAND
   2 #define SEQERRORCOMMAND
   3
   4 /*
   5  *  seqerrorcommand.h
   6  *  Mothur
   7  *
   8  *  Created by Pat Schloss on 7/15/10.
   9  *  Copyright 2010 Schloss Lab. All rights reserved.
  10  *
  11  */
  12
  13 #include "mothur.h"
  14 #include "command.hpp"
  15 #include "sequence.hpp"
  16 #include "referencedb.h"
  17
  18 struct Compare {
  19         int AA, AT, AG, AC,     TA, TT, TG, TC, GA, GT, GG, GC, CA, CT, CG, CC, NA, NT, NG, NC, Ai, Ti, Gi, Ci, Ni, dA, dT, dG, dC;
  20         string refName, queryName, sequence;
  21         double errorRate;
  22         int weight, matches, mismatches, total;
  23
  24         Compare(){
  25                 AA=0; AT=0; AG=0; AC=0;
  26                 TA=0; TT=0; TG=0; TC=0;
  27                 GA=0; GT=0; GG=0; GC=0;
  28                 CA=0; CT=0; CG=0; CC=0;
  29                 NA=0; NT=0; NG=0; NC=0;
  30                 Ai=0; Ti=0; Gi=0; Ci=0; Ni=0;
  31                 dA=0; dT=0; dG=0; dC=0;
  32                 refName = "";
  33                 queryName = "";
  34                 weight = 1;
  35                 matches = 0;
  36                 mismatches = 0;
  37                 total = 0;
  38                 errorRate = 1.0000;
  39                 sequence = "";
  40         }
  41 };
  42
  43 class SeqErrorCommand : public Command {
  44 public:
  45         SeqErrorCommand(string);
  46         SeqErrorCommand();
  47         ~SeqErrorCommand() {}
  48
  49         vector<string> setParameters();
  50         string getCommandName()                 { return "seq.error";                           }
  51         string getCommandCategory()             { return "Sequence Processing";         }
  52         string getHelpString();
  53         string getCitation() { return "http://www.mothur.org/wiki/Seq.error"; }
  54         string getDescription()         { return "seq.error"; }
  55
  56
  57         int execute();
  58         void help() { m->mothurOut(getHelpString()); }
  59
  60 private:
  61         bool abort;
  62         ReferenceDB* rdb;
  63
  64         struct linePair {
  65                 unsigned long int start;
  66                 unsigned long int end;
  67                 linePair(unsigned long int i, unsigned long int j) : start(i), end(j) {}
  68         };
  69
  70         vector<int> processIDS;   //processid
  71         vector<linePair> lines;
  72         vector<linePair> qLines;
  73         vector<linePair> rLines;
  74
  75         void getReferences();
  76         map<string,int> getWeights();
  77         Compare getErrors(Sequence, Sequence);
  78         void printErrorHeader(ofstream&);
  79         void printErrorData(Compare, int, ofstream&, ofstream&);
  80         void printSubMatrix();
  81         void printErrorFRFile(map<char, vector<int> >, map<char, vector<int> >);
  82         void printErrorQuality(map<char, vector<int> >);
  83         void printQualityFR(vector<vector<int> >, vector<vector<int> >);
  84
  85         int setLines(string, string, string, vector<unsigned long int>&, vector<unsigned long int>&, vector<unsigned long int>&);
  86         int driver(string, string, string, string, string, string, linePair, linePair, linePair);
  87         int createProcesses(string, string, string, string, string, string);
  88
  89         string queryFileName, referenceFileName, qualFileName, reportFileName, namesFileName, outputDir;
  90         double threshold;
  91         bool ignoreChimeras, save;
  92         int numRefs, processors;
  93         int maxLength, totalBases, totalMatches;
  94         //ofstream errorSummaryFile, errorSeqFile;
  95         vector<string> outputNames;
  96
  97         vector<Sequence> referenceSeqs;
  98         vector<vector<int> > substitutionMatrix;
  99         vector<vector<int> > qualForwardMap;
 100         vector<vector<int> > qualReverseMap;
 101         vector<int> misMatchCounts;
 102         map<char, vector<int> > qScoreErrorMap;
 103         map<char, vector<int> > errorForward;
 104         map<char, vector<int> > errorReverse;
 105         map<string, int> weights;
 106         vector<string> megaAlignVector;
 107
 108 };
 109
 110 #endif