r719: more stringent overlap merge
Will consider to make it the default
This commit is contained in:
parent
4e22270eba
commit
bdb7b000cd
18
bwamem.c
18
bwamem.c
|
|
@ -371,28 +371,32 @@ int mem_sort_and_dedup2(int n, mem_alnreg_t *a, float mask_level_redun, int merg
|
||||||
int64_t or, oq, mr, mq;
|
int64_t or, oq, mr, mq;
|
||||||
int is_merged = 0;
|
int is_merged = 0;
|
||||||
if (q->qe == q->qb) continue; // a[j] has been excluded
|
if (q->qe == q->qb) continue; // a[j] has been excluded
|
||||||
if (merge_bw > 0 && q->qb < p->qb && q->qe < p->qe && q->re < p->re) {
|
or = q->re - p->rb; // overlap length on the reference
|
||||||
int flag = 0;
|
oq = q->qb < p->qb? q->qe - p->qb : p->qe - q->qb; // overlap length on the query
|
||||||
|
mr = q->re - q->rb < p->re - p->rb? q->re - q->rb : p->re - p->rb; // min ref len in alignment
|
||||||
|
mq = q->qe - q->qb < p->qe - p->qb? q->qe - q->qb : p->qe - p->qb; // min qry len in alignment
|
||||||
|
if (merge_bw > 0 && q->qb < p->qb && q->qe < p->qe && q->re < p->re && fabs((float)oq / (p->qe - q->qb) - (float)or / (p->re - q->rb)) < .05) {
|
||||||
|
int flag = 0, q_s, r_s;
|
||||||
int64_t l1, l2;
|
int64_t l1, l2;
|
||||||
l1 = p->qb - q->qb, l2 = p->rb - q->rb;
|
l1 = p->qb - q->qb, l2 = p->rb - q->rb;
|
||||||
if (l1 - l2 < merge_bw && l2 - l1 < merge_bw) ++flag;
|
if (l1 - l2 < merge_bw && l2 - l1 < merge_bw) ++flag;
|
||||||
l1 = p->qe - q->qe, l2 = p->re - q->re;
|
l1 = p->qe - q->qe, l2 = p->re - q->re;
|
||||||
if (l1 - l2 < merge_bw && l2 - l1 < merge_bw) ++flag;
|
if (l1 - l2 < merge_bw && l2 - l1 < merge_bw) ++flag;
|
||||||
|
l1 = p->qe - q->qb, l2 = p->re - q->rb;
|
||||||
|
q_s = (int)((double)(p->qe - q->qb) / ((p->qe - p->qb) + (q->qe - q->qb)) * (p->score + q->score) + .499);
|
||||||
|
r_s = (int)((double)(p->re - q->rb) / ((p->re - p->rb) + (q->re - q->rb)) * (p->score + q->score) + .499);
|
||||||
if (flag == 2) { // merge q into p
|
if (flag == 2) { // merge q into p
|
||||||
mem_alnreg_t t = *p;
|
mem_alnreg_t t = *p;
|
||||||
*p = *q;
|
*p = *q;
|
||||||
p->qe = t.qe, p->re = t.re;
|
p->qe = t.qe, p->re = t.re;
|
||||||
p->score = p->score > t.score? p->score : t.score;
|
p->truesc = p->score = (q_s + r_s) >> 1;
|
||||||
p->w = p->w > t.w? p->w : t.w;
|
p->w = p->w > t.w? p->w : t.w;
|
||||||
|
p->w = p->w > abs(l1 - l2) + 50? p->w : abs(l1 - l2) + 50;
|
||||||
q->qb = q->qe;
|
q->qb = q->qe;
|
||||||
is_merged = 1;
|
is_merged = 1;
|
||||||
}
|
}
|
||||||
}
|
}
|
||||||
if (is_merged == 0) {
|
if (is_merged == 0) {
|
||||||
or = q->re - p->rb; // overlap length on the reference
|
|
||||||
oq = q->qb < p->qb? q->qe - p->qb : p->qe - q->qb; // overlap length on the query
|
|
||||||
mr = q->re - q->rb < p->re - p->rb? q->re - q->rb : p->re - p->rb; // min ref len in alignment
|
|
||||||
mq = q->qe - q->qb < p->qe - p->qb? q->qe - q->qb : p->qe - p->qb; // min qry len in alignment
|
|
||||||
if (or > mask_level_redun * mr && oq > mask_level_redun * mq) { // one of the hits is redundant
|
if (or > mask_level_redun * mr && oq > mask_level_redun * mq) { // one of the hits is redundant
|
||||||
if (p->score < q->score) {
|
if (p->score < q->score) {
|
||||||
p->qe = p->qb;
|
p->qe = p->qb;
|
||||||
|
|
|
||||||
Loading…
Reference in New Issue