From 89d4d219cd93c778c9053c322babfe1f6e6c1329 Mon Sep 17 00:00:00 2001 From: Heng Li Date: Wed, 20 Mar 2024 15:29:54 -0400 Subject: [PATCH] r1206: enabled RMQ for lr:hqae Also fixed a bug in determining inner_dist for RMQ. It should have no effect on previous presets. --- align.c | 4 ++-- lchain.c | 5 +++-- minimap.h | 2 +- options.c | 1 + 4 files changed, 7 insertions(+), 5 deletions(-) diff --git a/align.c b/align.c index a1d90ea6..9586dd80 100644 --- a/align.c +++ b/align.c @@ -933,14 +933,14 @@ double mm_event_identity(const mm_reg1_t *r) static int32_t mm_recal_max_dp(const mm_reg1_t *r, double b2, int32_t match_sc) { uint32_t i; - int32_t n_gap = 0, n_gapo = 0, n_mis; + int32_t n_gap = 0, n_mis; double gap_cost = 0.0; if (r->p == 0) return -1; for (i = 0; i < r->p->n_cigar; ++i) { int32_t op = r->p->cigar[i] & 0xf, len = r->p->cigar[i] >> 4; if (op == MM_CIGAR_INS || op == MM_CIGAR_DEL) { gap_cost += b2 + (double)mg_log2(1.0 + len); - ++n_gapo, n_gap += len; + n_gap += len; } } n_mis = r->blen + r->p->n_ambi - r->mlen - n_gap; diff --git a/lchain.c b/lchain.c index 9c76517b..d20aa3ab 100644 --- a/lchain.c +++ b/lchain.c @@ -203,7 +203,7 @@ mm128_t *mg_lchain_dp(int max_dist_x, int max_dist_y, int bw, int max_skip, int if (max_ii < 0 || (a[i].x - a[max_ii].x <= (int64_t)max_dist_x && f[max_ii] < f[i])) max_ii = i; if (mmax_f < max_f) mmax_f = max_f; - //fprintf(stderr, "X1\t%ld\t%ld:%d\t%ld\t%ld:%d\t%ld\t%ld\t%ld\n", (long)i, (long)(a[i].x>>32), (int32_t)a[i].x, (long)max_j, max_j<0?-1L:(long)(a[max_j].x>>32), max_j<0?-1:(int32_t)a[max_j].x, (long)max_f, (long)v[i], (long)mmax_f); + //fprintf(stderr, "X1\t%ld\t%ld:%d\t%ld\t%ld:%d\t%ld\t%ld\n", (long)i, (long)(a[i].x>>32), (int32_t)a[i].x, (long)max_j, max_j<0?-1L:(long)(a[max_j].x>>32), max_j<0?-1:(int32_t)a[max_j].x, (long)max_f, (long)v[i]); } u = mg_chain_backtrack(km, n, f, p, v, t, min_cnt, min_sc, max_drop, &n_u, &n_v); @@ -263,7 +263,8 @@ mm128_t *mg_lchain_rmq(int max_dist, int max_dist_inner, int bw, int max_chn_ski return 0; } if (max_dist < bw) max_dist = bw; - if (max_dist_inner <= 0 || max_dist_inner >= max_dist) max_dist_inner = 0; + if (max_dist_inner < 0) max_dist_inner = 0; + if (max_dist_inner > max_dist) max_dist_inner = max_dist; p = Kmalloc(km, int64_t, n); f = Kmalloc(km, int32_t, n); t = Kcalloc(km, int32_t, n); diff --git a/minimap.h b/minimap.h index 54f85965..8930cd5d 100644 --- a/minimap.h +++ b/minimap.h @@ -5,7 +5,7 @@ #include #include -#define MM_VERSION "2.27-r1205-dirty" +#define MM_VERSION "2.27-r1206-dirty" #define MM_F_NO_DIAG (0x001LL) // no exact diagonal hit #define MM_F_NO_DUAL (0x002LL) // skip pairs where query name is lexicographically larger than target name diff --git a/options.c b/options.c index cd2990d1..5dd6af95 100644 --- a/options.c +++ b/options.c @@ -116,6 +116,7 @@ int mm_set_opt(const char *preset, mm_idxopt_t *io, mm_mapopt_t *mo) } } else if (strcmp(preset, "lr:hqae") == 0) { // high-quality assembly evaluation io->flag = 0, io->k = 25, io->w = 51; + mo->flag |= MM_F_RMQ; mo->min_mid_occ = 50, mo->max_mid_occ = 500; mo->rmq_inner_dist = 5000; mo->occ_dist = 200;