diff --git a/src/util/convertalignments.cpp b/src/util/convertalignments.cpp index a2b9a4e06..637e1790c 100644 --- a/src/util/convertalignments.cpp +++ b/src/util/convertalignments.cpp @@ -349,6 +349,9 @@ int convertalignments(int argc, const char **argv, const Command &command) { std::string newBacktrace; newBacktrace.reserve(1024); + std::string complementBuffer; + complementBuffer.reserve(1024); + const TaxonNode * taxonNode = NULL; #pragma omp for schedule(dynamic, 10) @@ -714,32 +717,45 @@ int convertalignments(int argc, const char **argv, const Command &command) { break; } case Parameters::FORMAT_ALIGNMENT_SAM: { - bool strand = res.qEndPos > res.qStartPos; + bool forward = res.qEndPos > res.qStartPos; int rawScore = static_cast(evaluer->computeRawScoreFromBitScore(res.score) + 0.5); uint32_t mapq = -4.343 * log(exp(static_cast(-rawScore))); mapq = (uint32_t) (mapq + 4.99); mapq = mapq < 254 ? mapq : 254; - int count = snprintf(buffer, sizeof(buffer), "%s\t%d\t%s\t%d\t%d\t", queryId.c_str(), (strand) ? 16: 0, targetId.c_str(), res.dbStartPos + 1, mapq); + int count = snprintf(buffer, sizeof(buffer), "%s\t%d\t%s\t%d\t%d\t", queryId.c_str(), (forward) ? 0 : 16, targetId.c_str(), res.dbStartPos + 1, mapq); if (count < 0 || static_cast(count) >= sizeof(buffer)) { Debug(Debug::WARNING) << "Truncated line in entry" << i << "!\n"; continue; } result.append(buffer, count); - if (isTranslatedSearch == true && targetNucs == true && queryNucs == true) { + if (isTranslatedSearch == true && queryNucs == true) { Matcher::result_t::protein2nucl(res.backtrace, newBacktrace); result.append(newBacktrace); newBacktrace.clear(); - } else { result.append(res.backtrace); } result.append("\t*\t0\t0\t"); int start = std::min(res.qStartPos, res.qEndPos); int end = std::max(res.qStartPos, res.qEndPos); - if (queryProfile) { - result.append(queryProfData.c_str() + start, (end + 1) - start); + if (forward == false && queryNucs == true) { + if (queryProfile) { + complementBuffer.assign(queryProfData.c_str() + start, (end + 1) - start); + } else { + complementBuffer.assign(querySeqData + start, (end + 1) - start); + } + std::reverse(complementBuffer.begin(), complementBuffer.end()); + for (size_t i = 0; i < complementBuffer.size(); i++) { + complementBuffer[i] = Orf::complement(complementBuffer[i]); + } + result.append(complementBuffer); + complementBuffer.clear(); } else { - result.append(querySeqData + start, (end + 1) - start); + if (queryProfile) { + result.append(queryProfData.c_str() + start, (end + 1) - start); + } else { + result.append(querySeqData + start, (end + 1) - start); + } } count = snprintf(buffer, sizeof(buffer), "\t*\tAS:i:%d\tNM:i:%d\n", rawScore, missMatchCount); if (count < 0 || static_cast(count) >= sizeof(buffer)) {