diff --git a/src/common/wflign/src/wflign_alignment.cpp b/src/common/wflign/src/wflign_alignment.cpp index 3bb42731..98471ef8 100644 --- a/src/common/wflign/src/wflign_alignment.cpp +++ b/src/common/wflign/src/wflign_alignment.cpp @@ -14,7 +14,7 @@ // Default constructor alignment_t::alignment_t() - : j(0), i(0), query_length(0), target_length(0), ok(false), keep(false) { + : j(0), i(0), query_length(0), target_length(0), score(std::numeric_limits::max()), ok(false), keep(false) { edit_cigar = {nullptr, 0, 0}; } diff --git a/src/common/wflign/src/wflign_alignment.hpp b/src/common/wflign/src/wflign_alignment.hpp index 23f00445..8a2127c5 100644 --- a/src/common/wflign/src/wflign_alignment.hpp +++ b/src/common/wflign/src/wflign_alignment.hpp @@ -35,6 +35,7 @@ class alignment_t { int i; int query_length; int target_length; + int score; bool ok; bool keep; bool is_rev; diff --git a/src/common/wflign/src/wflign_patch.cpp b/src/common/wflign/src/wflign_patch.cpp index ff960406..fa050908 100644 --- a/src/common/wflign/src/wflign_patch.cpp +++ b/src/common/wflign/src/wflign_patch.cpp @@ -166,8 +166,8 @@ void do_wfa_patch_alignment( wf_aligner.setMaxAlignmentSteps(max_score); - int fwd_score = std::numeric_limits::max(); - int rev_score = std::numeric_limits::max(); + //int fwd_score = std::numeric_limits::max(); + //int rev_score = std::numeric_limits::max(); const int status = wf_aligner.alignEnd2End(target + i, target_length, query + j, query_length); aln.ok = (status == WF_STATUS_ALG_COMPLETED); @@ -191,13 +191,13 @@ void do_wfa_patch_alignment( #endif wflign_edit_cigar_copy(wf_aligner, &aln.edit_cigar); - fwd_score = calculate_alignment_score(aln.edit_cigar, convex_penalties); + aln.score = calculate_alignment_score(aln.edit_cigar, convex_penalties); //std::cerr << "forward score is " << fwd_score << std::endl; } if (query_length >= min_inversion_length && target_length >= min_inversion_length) { if (aln.ok) { - wf_aligner.setMaxAlignmentSteps(fwd_score * 2 + 64); + wf_aligner.setMaxAlignmentSteps(aln.score); } // Try reverse complement alignment std::string rev_comp_query = reverse_complement(std::string(query + j, query_length)); @@ -225,7 +225,7 @@ void do_wfa_patch_alignment( } #endif - rev_score = calculate_alignment_score(rev_aln.edit_cigar, convex_penalties); + rev_aln.score = calculate_alignment_score(rev_aln.edit_cigar, convex_penalties); //std::cerr << "reverse score is " << rev_score << std::endl; rev_aln.j = j; rev_aln.i = i; @@ -236,7 +236,7 @@ void do_wfa_patch_alignment( rev_aln.ok = false; } - if (rev_aln.ok && rev_score < fwd_score) { + if (rev_aln.ok && rev_aln.score < aln.score) { rev_aln.ok = true; aln.ok = false; #ifdef WFLIGN_DEBUG @@ -246,8 +246,8 @@ void do_wfa_patch_alignment( << " chain_gap " << chain_gap << " max_patching_score " << max_patching_score << " max_score " << max_score - << " fwd_score " << fwd_score - << " rev_score " << rev_score + << " fwd_score " << aln.score + << " rev_score " << rev_aln.score << std::endl << "query " << std::string(query + j, query_length) << " target " << std::string(target + i, target_length) @@ -470,7 +470,7 @@ std::vector do_progressive_wfa_patch_alignment( min_inversion_length); if (aln.ok || rev_aln.ok) { - alignment_t& chosen_aln = (rev_aln.ok && (!aln.ok || calculate_alignment_score(rev_aln.edit_cigar, convex_penalties) < calculate_alignment_score(aln.edit_cigar, convex_penalties))) ? rev_aln : aln; + alignment_t& chosen_aln = (rev_aln.ok && (!aln.ok || rev_aln.score < aln.score)) ? rev_aln : aln; // Erode short matches throughout the alignment //erode_alignment(chosen_aln, erode_k); if (chosen_aln.query_length > 0 && chosen_aln.target_length > 0) { @@ -1677,13 +1677,13 @@ void write_merged_alignment( //std::cerr << "FIRST PATCH ROUND" << std::endl; // +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++\n"; - patching(erodev, pre_tracev, 4096, 32, 512, false); // In the 1st round, we patch more aggressively and don't save inversions + patching(erodev, tracev, 4096, 32, 512, true); // one patchin gpass #ifdef VALIDATE_WFA_WFLIGN - if (!validate_trace(pre_tracev, query, + if (!validate_trace(tracev, query, target - target_pointer_shift, query_length, target_length_mut, query_start, target_start)) { - std::cerr << "cigar failure in pre_tracev " + std::cerr << "cigar failure in tracev " << "\t" << query_name << "\t" << query_total_length << "\t" << query_offset + (query_is_rev @@ -1702,12 +1702,13 @@ void write_merged_alignment( #endif // normalize: sort so that I