m = MothurOut::getInstance();
seqName = "";
seqLength = -1;
+
}
catch(exception& e) {
m->errorOut(e, "QualityScores", "QualityScores");
/**************************************************************************************************/
-QualityScores::QualityScores(ifstream& qFile, int l){
+QualityScores::QualityScores(ifstream& qFile){
try {
m = MothurOut::getInstance();
-
+
seqName = "";
- seqLength = l;
int score;
- string line;
- getline(qFile, line); gobble(qFile);
- istringstream nameStream(line);
-
- nameStream >> seqName;
- seqName = seqName.substr(1);
-
- //getline(qFile, line);
- //istringstream qualStream(line);
-
- //while(qualStream){
- // qualStream >> score;
- // qScores.push_back(score);
- //}
- //qScores.pop_back();
+ qFile >> seqName;
+ m->getline(qFile);
+
+ if (seqName == "") {
+ m->mothurOut("Error reading quality file, name blank at position, " + toString(qFile.tellg()));
+ m->mothurOutEndLine();
+ }
+ else{
+ seqName = seqName.substr(1);
+ }
- //seqLength = qScores.size();
+ string qScoreString = m->getline(qFile);
- for(int i=0;i<seqLength;i++){
- qFile >> score;
+ istringstream qScoreStringStream(qScoreString);
+ while(!qScoreStringStream.eof()){
+ qScoreStringStream >> score;
qScores.push_back(score);
}
- gobble(qFile);
-
+ qScores.pop_back();
+ seqLength = qScores.size();
}
catch(exception& e) {
m->errorOut(e, "QualityScores", "QualityScores");
exit(1);
}
-
+
}
/**************************************************************************************************/
qScores = hold;
}
if(start == -1){
- hold = vector<int>(qScores.begin(), qScores.begin()+end); //not sure if indexing is correct
- qScores = hold;
+ if(qScores.size() > end){
+ hold = vector<int>(qScores.begin(), qScores.begin()+end);
+ qScores = hold;
+ }
}
seqLength = qScores.size();
}
}
+ //every score passed
+ if (end == (seqLength-1)) { end = seqLength; }
+
sequence.setUnaligned(rawSequence.substr(0,end));
trimQScores(-1, end);
-
+
return 1;
}
catch(exception& e) {
if(end == -1){ end = seqLength; }
+
sequence.setUnaligned(rawSequence.substr(0,end));
trimQScores(-1, end);
+
return 1;
}
catch(exception& e) {
int seqLength = sequence.getNumBases();
if(seqName != sequence.getName()){
- m->mothurOut("sequence name mismatch btwn fasta: " + sequence.getName() + " and qual file: " + seqName);
- m->mothurOutEndLine();
+ m->mothurOut("sequence name mismatch between fasta: " + sequence.getName() + " and qual file: " + seqName);
+ m->mothurOutEndLine();
}
int end = windowSize;
int start = 0;
-
+ if(seqLength < windowSize) { return 0; }
+
while(start < seqLength){
double windowSum = 0.0000;
for(int i=start;i<end;i++){
- windowSum += qScores[i];
+ windowSum += qScores[i];
}
double windowAverage = windowSum / (double)(end-start);
if(windowAverage < qThreshold){
- end = start;
+ end = end - stepSize;
break;
}
start += stepSize;
if(end >= seqLength){ end = seqLength - 1; }
}
-
if(end == -1){ end = seqLength; }
sequence.setUnaligned(rawSequence.substr(0,end));
return success;
}
catch(exception& e) {
- m->errorOut(e, "TrimSeqsCommand", "cullQualAverage");
+ m->errorOut(e, "QualityScores", "cullQualAverage");
+ exit(1);
+ }
+}
+
+/**************************************************************************************************/
+
+void QualityScores::updateQScoreErrorMap(map<char, vector<int> >& qualErrorMap, string errorSeq, int start, int stop, int weight){
+ try {
+
+ int seqLength = errorSeq.size();
+ int qIndex = start - 1;
+ for(int i=0;i<seqLength;i++){
+
+ if(errorSeq[i] == 'm') { qualErrorMap['m'][qScores[qIndex]] += weight; }
+ else if(errorSeq[i] == 's') { qualErrorMap['s'][qScores[qIndex]] += weight; }
+ else if(errorSeq[i] == 'i') { qualErrorMap['i'][qScores[qIndex]] += weight; }
+ else if(errorSeq[i] == 'a') { qualErrorMap['a'][qScores[qIndex]] += weight; }
+ else if(errorSeq[i] == 'd') { /* there are no qScores for deletions */ }
+
+ if(errorSeq[i] != 'd') { qIndex++; }
+
+ }
+ }
+ catch(exception& e) {
+ m->errorOut(e, "QualityScores", "updateQScoreErrorMap");
+ exit(1);
+ }
+}
+
+/**************************************************************************************************/
+
+void QualityScores::updateForwardMap(vector<vector<int> >& forwardMap, int start, int stop, int weight){
+ try {
+
+ int index = 0;
+ for(int i=start-1;i<stop;i++){
+ forwardMap[index++][qScores[i]] += weight;
+ }
+
+ }
+ catch(exception& e) {
+ m->errorOut(e, "QualityScores", "updateForwardMap");
+ exit(1);
+ }
+}
+
+/**************************************************************************************************/
+
+void QualityScores::updateReverseMap(vector<vector<int> >& reverseMap, int start, int stop, int weight){
+ try {
+
+ int index = 0;
+ for(int i=stop-1;i>=start;i--){
+ reverseMap[index++][qScores[i]] += weight;
+ }
+
+ }
+ catch(exception& e) {
+ m->errorOut(e, "QualityScores", "updateForwardMap");
exit(1);
}
}