zzh
/
minimap2
1
0
Fork 0
Code Issues Pull Requests Packages Projects Releases Wiki Activity

r205: when computing mapq, consider min_chain_sc 

Not doing this was a mistake.
This commit is contained in:
Heng Li 2017-07-26 11:34:14 -04:00
parent a8ad53ee81
commit e9dc1ce2b6
4 changed files with 8 additions and 7 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

9
hit.c
View File

@ -288,17 +288,18 @@ void mm_join_long(void *km, const mm_mapopt_t *opt, int qlen, int *n_regs_, mm_r
}
}
void mm_set_mapq(int n_regs, mm_reg1_t *regs)
void mm_set_mapq(int n_regs, mm_reg1_t *regs, int min_chain_sc)
{
int i;
for (i = 0; i < n_regs; ++i) {
mm_reg1_t *r = &regs[i];
if (r->parent == r->id) {
int mapq;
int mapq, subsc;
subsc = r->subsc > min_chain_sc? r->subsc : min_chain_sc;
if (r->p && r->p->dp_max2 > 0 && r->p->dp_max > 0) {
float identity = (float)(r->p->blen - r->p->n_diff - r->p->n_ambi) / (r->p->blen - r->p->n_ambi);
mapq = (int)(identity * 30.0 * (1. - (float)r->p->dp_max2 * r->subsc / r->p->dp_max / r->score) * logf(r->score));
} else mapq = (int)(30.0 * (1. - (float)r->subsc / r->score) * logf(r->score));
mapq = (int)(identity * 30.0 * (1. - (float)r->p->dp_max2 * subsc / r->p->dp_max / r->score) * logf(r->score));
} else mapq = (int)(30.0 * (1. - (float)subsc / r->score) * logf(r->score));
mapq = mapq > 0? mapq : 0;
r->mapq = mapq < 60? mapq : 60;
} else r->mapq = 0;

2
main.c
View File

@ -8,7 +8,7 @@
#include "minimap.h"
#include "mmpriv.h"
#define MM_VERSION "2.0-r203-dirty"
#define MM_VERSION "2.0-r205-dirty"
void liftrlimit()
{

2
map.c
View File

@ -255,7 +255,7 @@ mm_reg1_t *mm_map_frag(const mm_mapopt_t *opt, const mm_idx_t *mi, mm_tbuf_t *b,
mm_set_sam_pri(*n_regs, regs);
}
}
mm_set_mapq(*n_regs, regs);
mm_set_mapq(*n_regs, regs, opt->min_chain_score);
// free
kfree(b->km, a);

2
mmpriv.h
View File

@ -49,7 +49,7 @@ void mm_select_sub(void *km, float mask_level, float pri_ratio, int min_diff, in
void mm_filter_regs(void *km, const mm_mapopt_t *opt, int *n_regs, mm_reg1_t *regs);
void mm_join_long(void *km, const mm_mapopt_t *opt, int qlen, int *n_regs, mm_reg1_t *regs, mm128_t *a);
void mm_hit_sort_by_dp(void *km, int *n_regs, mm_reg1_t *r);
void mm_set_mapq(int n_regs, mm_reg1_t *regs);
void mm_set_mapq(int n_regs, mm_reg1_t *regs, int min_chain_sc);
#ifdef __cplusplus
}