r759: bugfix - frac_rep not working
Also added commented code for a 3rd round seeding. Not used.
This commit is contained in:
parent
c9b33502f3
commit
43b498a37e
22
bwamem.c
22
bwamem.c
|
|
@ -131,7 +131,27 @@ static void mem_collect_intv(const mem_opt_t *opt, const bwt_t *bwt, int len, co
|
||||||
if (end - start < split_len || p->x[2] > opt->split_width) continue;
|
if (end - start < split_len || p->x[2] > opt->split_width) continue;
|
||||||
bwt_smem1(bwt, len, seq, (start + end)>>1, p->x[2]+1, &a->mem1, a->tmpv);
|
bwt_smem1(bwt, len, seq, (start + end)>>1, p->x[2]+1, &a->mem1, a->tmpv);
|
||||||
for (i = 0; i < a->mem1.n; ++i)
|
for (i = 0; i < a->mem1.n; ++i)
|
||||||
kv_push(bwtintv_t, a->mem, a->mem1.a[i]);
|
if ((a->mem1.a[i].info>>32) - (uint32_t)a->mem1.a[i].info >= opt->min_seed_len)
|
||||||
|
kv_push(bwtintv_t, a->mem, a->mem1.a[i]);
|
||||||
|
/*
|
||||||
|
{ // third pass:
|
||||||
|
int max_len = 0, max_i = -1;
|
||||||
|
for (i = 0; i < a->mem1.n; ++i) {
|
||||||
|
bwtintv_t *q = &a->mem1.a[i];
|
||||||
|
int start = q->info>>32, end = (int32_t)q->info;
|
||||||
|
if (end - start > max_len)
|
||||||
|
max_len = end - start, max_i = i;
|
||||||
|
}
|
||||||
|
if (max_i >= 0 && max_len > .9 * (end - start) && max_len >= split_len) {
|
||||||
|
bwtintv_t t = a->mem1.a[max_i];
|
||||||
|
int start = t.info>>32, end = (int32_t)t.info;
|
||||||
|
bwt_smem1(bwt, len, seq, (start + end)>>1, t.x[2]+1, &a->mem1, a->tmpv);
|
||||||
|
for (i = 0; i < a->mem1.n; ++i)
|
||||||
|
if ((a->mem1.a[i].info>>32) - (uint32_t)a->mem1.a[i].info >= opt->min_seed_len)
|
||||||
|
kv_push(bwtintv_t, a->mem, a->mem1.a[i]);
|
||||||
|
}
|
||||||
|
}
|
||||||
|
*/
|
||||||
}
|
}
|
||||||
// sort
|
// sort
|
||||||
ks_introsort(mem_intv, a->mem.n, a->mem.a);
|
ks_introsort(mem_intv, a->mem.n, a->mem.a);
|
||||||
|
|
|
||||||
2
bwamem.h
2
bwamem.h
|
|
@ -65,7 +65,7 @@ typedef struct {
|
||||||
int secondary; // index of the parent hit shadowing the current hit; <0 if primary
|
int secondary; // index of the parent hit shadowing the current hit; <0 if primary
|
||||||
int seedlen0; // length of the starting seed
|
int seedlen0; // length of the starting seed
|
||||||
int n_comp; // number of sub-alignments chained together
|
int n_comp; // number of sub-alignments chained together
|
||||||
int frac_rep;
|
float frac_rep;
|
||||||
uint64_t hash;
|
uint64_t hash;
|
||||||
} mem_alnreg_t;
|
} mem_alnreg_t;
|
||||||
|
|
||||||
|
|
|
||||||
|
|
@ -289,6 +289,7 @@ int mem_sam_pe(const mem_opt_t *opt, const bntseq_t *bns, const uint8_t *pac, co
|
||||||
if (n_sub > 0) q_pe -= (int)(4.343 * log(n_sub+1) + .499);
|
if (n_sub > 0) q_pe -= (int)(4.343 * log(n_sub+1) + .499);
|
||||||
if (q_pe < 0) q_pe = 0;
|
if (q_pe < 0) q_pe = 0;
|
||||||
if (q_pe > 60) q_pe = 60;
|
if (q_pe > 60) q_pe = 60;
|
||||||
|
q_pe = (int)(q_pe * (1. - .5 * (a[0].a[0].frac_rep + a[1].a[0].frac_rep)) + .499);
|
||||||
// the following assumes no split hits
|
// the following assumes no split hits
|
||||||
if (o > score_un) { // paired alignment is preferred
|
if (o > score_un) { // paired alignment is preferred
|
||||||
mem_alnreg_t *c[2];
|
mem_alnreg_t *c[2];
|
||||||
|
|
|
||||||
Loading…
Reference in New Issue