temp = temp.substr(0,alnLength);
int numDiff = countDiffs(oligo, temp);
+ if (m->debug) { m->mothurOut("[DEBUG]: " + seq.getName() + " aligned fragment =" + temp + ", barcode =" + oligo + ", numDiffs = " + toString(numDiff) + "\n"); }
+
if(numDiff < minDiff){
minDiff = numDiff;
minCount = 1;
break;
}
- if((compareDNASeq(foligo, rawFSequence.substr(0,foligo.length()))) && (compareDNASeq(roligo, rawRSequence.substr((rawRSequence.length()-roligo.length()),roligo.length())))) {
+ if((compareDNASeq(foligo, rawFSequence.substr(0,foligo.length()))) && (compareDNASeq(roligo, rawRSequence.substr(0,roligo.length())))) {
group = it->first;
forwardSeq.setUnaligned(rawFSequence.substr(foligo.length()));
- reverseSeq.setUnaligned(rawRSequence.substr(0,(rawRSequence.length()-roligo.length())));
+ reverseSeq.setUnaligned(rawRSequence.substr(roligo.length()));
success = 0;
break;
}
temp = temp.substr(0,alnLength);
int numDiff = countDiffs(oligo, temp);
+ if (m->debug) { m->mothurOut("[DEBUG]: forward " + forwardSeq.getName() + " aligned fragment=" + temp + ", barcode=" + oligo + ", numDiffs=" + toString(numDiff) + ".\n"); }
+
if (alnLength == 0) { numDiff = bdiffs + 100; }
//cout << "after = " << oligo << '\t' << temp << '\t' << numDiff << endl;
oligo = oligo.substr(0,alnLength);
temp = temp.substr(0,alnLength);
int numDiff = countDiffs(oligo, temp);
+
+ if (m->debug) { m->mothurOut("[DEBUG]: reverse " + forwardSeq.getName() + " aligned fragment=" + temp + ", barcode=" + oligo + ", numDiffs=" + toString(numDiff) + ".\n"); }
+
if (alnLength == 0) { numDiff = bdiffs + 100; }
//cout << "after = " << oligo << '\t' << temp << '\t' << numDiff << endl;
}
- /*cout << minDiff << '\t' << minCount << '\t' << endl;
- for (int i = 0; i < minFGroup.size(); i++) {
- cout << i << '\t';
- for (int j = 0; j < minFGroup[i].size(); j++) { cout << minFGroup[i][j] << " "; }
- cout << endl;
- }
- cout << endl;
- for (int i = 0; i < minRGroup.size(); i++) {
- cout << i << '\t';
- for (int j = 0; j < minRGroup[i].size(); j++) { cout << minRGroup[i][j] << " "; }
- cout << endl;
- }
- cout << endl;*/
if(minDiff > bdiffs) { success = minDiff; } //no good matches
else {
bool foundMatch = false;
break;
}
- if((compareDNASeq(foligo, rawFSequence.substr(0,foligo.length()))) && (compareDNASeq(roligo, rawRSequence.substr((rawRSequence.length()-roligo.length()),roligo.length())))) {
+ if((compareDNASeq(foligo, rawFSequence.substr(0,foligo.length()))) && (compareDNASeq(roligo, rawRSequence.substr(0,roligo.length())))) {
group = it->first;
forwardSeq.setUnaligned(rawFSequence.substr(foligo.length()));
- reverseSeq.setUnaligned(rawRSequence.substr(0,(rawRSequence.length()-roligo.length())));
+ reverseSeq.setUnaligned(rawRSequence.substr(roligo.length()));
forwardQual.trimQScores(foligo.length(), -1);
- reverseQual.trimQScores(-1, rawRSequence.length()-roligo.length());
+ reverseQual.trimQScores(roligo.length(), -1);
success = 0;
break;
}
int numDiff = countDiffs(oligo, temp);
if (alnLength == 0) { numDiff = bdiffs + 100; }
- //cout << "after = " << oligo << '\t' << temp << '\t' << numDiff << endl;
+ if (m->debug) { m->mothurOut("[DEBUG]: forward " + forwardSeq.getName() + " aligned fragment=" + temp + ", barcode=" + oligo + ", numDiffs=" + toString(numDiff) + ".\n"); }
if(numDiff < minDiff){
minDiff = numDiff;
int numDiff = countDiffs(oligo, temp);
if (alnLength == 0) { numDiff = bdiffs + 100; }
+ if (m->debug) { m->mothurOut("[DEBUG]: reverse " + reverseSeq.getName() + " aligned fragment=" + temp + ", barcode=" + oligo + ", numDiffs=" + toString(numDiff) + ".\n"); }
+
//cout << "after = " << oligo << '\t' << temp << '\t' << numDiff << endl;
if(numDiff < minDiff){
minDiff = numDiff;
}
- /*cout << minDiff << '\t' << minCount << '\t' << endl;
- for (int i = 0; i < minFGroup.size(); i++) {
- cout << i << '\t';
- for (int j = 0; j < minFGroup[i].size(); j++) { cout << minFGroup[i][j] << " "; }
- cout << endl;
- }
- cout << endl;
- for (int i = 0; i < minRGroup.size(); i++) {
- cout << i << '\t';
- for (int j = 0; j < minRGroup[i].size(); j++) { cout << minRGroup[i][j] << " "; }
- cout << endl;
- }
- cout << endl;*/
if(minDiff > bdiffs) { success = minDiff; } //no good matches
else {
bool foundMatch = false;
break;
}
+ if (rawSeq.length() < (foligo.length() + roligo.length())) { success = pdiffs + 10; break; }
+
if((compareDNASeq(foligo, rawSeq.substr(0,foligo.length()))) && (compareDNASeq(roligo, rawSeq.substr(rawSeq.length()-roligo.length(),roligo.length())))) {
group = it->first;
string trimmedSeq = rawSeq.substr(foligo.length()); //trim forward barcode
}
- /*cout << minDiff << '\t' << minCount << '\t' << endl;
- for (int i = 0; i < minFGroup.size(); i++) {
- cout << i << '\t';
- for (int j = 0; j < minFGroup[i].size(); j++) { cout << minFGroup[i][j] << " "; }
- cout << endl;
- }
- cout << endl;
- for (int i = 0; i < minRGroup.size(); i++) {
- cout << i << '\t';
- for (int j = 0; j < minRGroup[i].size(); j++) { cout << minRGroup[i][j] << " "; }
- cout << endl;
- }
- cout << endl;*/
if(minDiff > bdiffs) { success = minDiff; } //no good matches
else {
bool foundMatch = false;
group = matches[0];
fStart = minFPos[0];
rStart = rawSeq.length() - minRPos[0];
+ if (fStart > rStart) { foundMatch = false; } //only barcodes not a good sequence
}
//we have an acceptable match for the forward and reverse, but do they match?
break;
}
+ if (rawSeq.length() < (foligo.length() + roligo.length())) { success = pdiffs + 10; break; }
+
if((compareDNASeq(foligo, rawSeq.substr(0,foligo.length()))) && (compareDNASeq(roligo, rawSeq.substr(rawSeq.length()-roligo.length(),roligo.length())))) {
group = it->first;
if (!keepForward) {
success = pdiffs + 10;
break;
}
- //cout << "before = " << oligo << '\t' << rawSeq.substr(0,oligo.length()+bdiffs) << endl;
+ //cout << "before = " << oligo << '\t' << rawSeq.substr(0,oligo.length()+pdiffs) << endl;
//use needleman to align first barcode.length()+numdiffs of sequence to each barcode
alignment->alignPrimer(oligo, rawSeq.substr(0,oligo.length()+pdiffs));
oligo = alignment->getSeqAAln();
}
- /*cout << minDiff << '\t' << minCount << '\t' << endl;
- for (int i = 0; i < minFGroup.size(); i++) {
- cout << i << '\t';
- for (int j = 0; j < minFGroup[i].size(); j++) { cout << minFGroup[i][j] << " "; }
- cout << endl;
- }
- cout << endl;
- for (int i = 0; i < minRGroup.size(); i++) {
- cout << i << '\t';
- for (int j = 0; j < minRGroup[i].size(); j++) { cout << minRGroup[i][j] << " "; }
- cout << endl;
- }
- cout << endl;*/
if(minDiff > pdiffs) { success = minDiff; } //no good matches
else {
bool foundMatch = false;
group = matches[0];
fStart = minFPos[0];
rStart = rawSeq.length() - minRPos[0];
+ if (fStart > rStart) { foundMatch = false; } //only primers not a good sequence
}
//we have an acceptable match for the forward and reverse, but do they match?
break;
}
- if((compareDNASeq(foligo, rawFSequence.substr(0,foligo.length()))) && (compareDNASeq(roligo, rawRSequence.substr((rawRSequence.length()-roligo.length()),roligo.length())))) {
+ if((compareDNASeq(foligo, rawFSequence.substr(0,foligo.length()))) && (compareDNASeq(roligo, rawRSequence.substr(0,roligo.length())))) {
group = it->first;
forwardSeq.setUnaligned(rawFSequence.substr(foligo.length()));
- reverseSeq.setUnaligned(rawRSequence.substr(0,(rawRSequence.length()-roligo.length())));
+ reverseSeq.setUnaligned(rawRSequence.substr(roligo.length()));
forwardQual.trimQScores(foligo.length(), -1);
- reverseQual.trimQScores(-1, rawRSequence.length()-roligo.length());
+ reverseQual.trimQScores(roligo.length(), -1);
success = 0;
break;
}
if (alnLength == 0) { numDiff = pdiffs + 100; }
//cout << "after = " << oligo << '\t' << temp << '\t' << numDiff << endl;
+ if (m->debug) { m->mothurOut("[DEBUG]: forward " + forwardSeq.getName() + " aligned fragment=" + temp + ", primer=" + oligo + ", numDiffs=" + toString(numDiff) + ".\n"); }
+
if(numDiff < minDiff){
minDiff = numDiff;
minCount = 1;
int numDiff = countDiffs(oligo, temp);
if (alnLength == 0) { numDiff = pdiffs + 100; }
+ if (m->debug) { m->mothurOut("[DEBUG]: reverse " + forwardSeq.getName() + " aligned fragment=" + temp + ", primer=" + oligo + ", numDiffs=" + toString(numDiff) + ".\n"); }
+
//cout << "after = " << oligo << '\t' << temp << '\t' << numDiff << endl;
if(numDiff < minDiff){
minDiff = numDiff;
}
- /*cout << minDiff << '\t' << minCount << '\t' << endl;
- for (int i = 0; i < minFGroup.size(); i++) {
- cout << i << '\t';
- for (int j = 0; j < minFGroup[i].size(); j++) { cout << minFGroup[i][j] << " "; }
- cout << endl;
- }
- cout << endl;
- for (int i = 0; i < minRGroup.size(); i++) {
- cout << i << '\t';
- for (int j = 0; j < minRGroup[i].size(); j++) { cout << minRGroup[i][j] << " "; }
- cout << endl;
- }
- cout << endl;*/
if(minDiff > pdiffs) { success = minDiff; } //no good matches
else {
bool foundMatch = false;
if((compareDNASeq(foligo, rawFSequence.substr(0,foligo.length()))) && (compareDNASeq(roligo, rawRSequence.substr((rawRSequence.length()-roligo.length()),roligo.length())))) {
group = it->first;
forwardSeq.setUnaligned(rawFSequence.substr(foligo.length()));
- reverseSeq.setUnaligned(rawRSequence.substr(0,(rawRSequence.length()-roligo.length())));
+ reverseSeq.setUnaligned(rawRSequence.substr(roligo.length()));
success = 0;
break;
}
temp = temp.substr(0,alnLength);
int numDiff = countDiffs(oligo, temp);
+ if (m->debug) { m->mothurOut("[DEBUG]: forward " + forwardSeq.getName() + " aligned fragment=" + temp + ", primer=" + oligo + ", numDiffs=" + toString(numDiff) + ".\n"); }
+
if (alnLength == 0) { numDiff = pdiffs + 100; }
//cout << "after = " << oligo << '\t' << temp << '\t' << numDiff << endl;
oligo = oligo.substr(0,alnLength);
temp = temp.substr(0,alnLength);
int numDiff = countDiffs(oligo, temp);
+
+ if (m->debug) { m->mothurOut("[DEBUG]: reverse " + forwardSeq.getName() + " aligned fragment=" + temp + ", primer=" + oligo + ", numDiffs=" + toString(numDiff) + ".\n"); }
+
if (alnLength == 0) { numDiff = pdiffs + 100; }
//cout << "after = " << oligo << '\t' << temp << '\t' << numDiff << endl;
}
- /*cout << minDiff << '\t' << minCount << '\t' << endl;
- for (int i = 0; i < minFGroup.size(); i++) {
- cout << i << '\t';
- for (int j = 0; j < minFGroup[i].size(); j++) { cout << minFGroup[i][j] << " "; }
- cout << endl;
- }
- cout << endl;
- for (int i = 0; i < minRGroup.size(); i++) {
- cout << i << '\t';
- for (int j = 0; j < minRGroup[i].size(); j++) { cout << minRGroup[i][j] << " "; }
- cout << endl;
- }
- cout << endl;*/
if(minDiff > pdiffs) { success = minDiff; } //no good matches
else {
bool foundMatch = false;
int numDiff = countDiffs(oligo, temp);
+ if (m->debug) { m->mothurOut("[DEBUG]: " + seq.getName() + " aligned fragment =" + temp + ", barcode =" + oligo + ", numDiffs = " + toString(numDiff) + "\n"); }
+
if(numDiff < minDiff){
minDiff = numDiff;
minCount = 1;