xxay/hifiasm
1
0
Fork 0
mirror of https://github.com/chhylp123/hifiasm.git synced 2026-09-15 12:47:57 +08:00
Code Issues Packages Projects Releases Wiki Activity

adjusted chain sc

Browse Source
This commit is contained in:
chhylp123
2022-10-09 18:14:07 -04:00
parent 1d0562ebd1
commit 3d65bef9de
3 changed files with 31 additions and 21 deletions
Download Patch File Download Diff File Expand all files Collapse all files

38
Correct.cpp
View File

@@ -14970,13 +14970,15 @@ uint64_t query_gen_gov_idx(asg64_v *ovidx, uint64_t v, uint64_t w)
return 0;
}
// #define id_mm ((uint64_t)0x7fffffffffffffff)
#define id_set ((uint64_t)0x8000000000000000)
#define id_get(a) ((uint32_t)(a))
#define err_get(a) (((a)&((uint64_t)0x7fffffffffffffff))>>32)
///[s, e)
uint64_t gen_region_phase(overlap_region* ol, uint64_t *id_a, uint64_t id_n, uint64_t s, uint64_t e, uint64_t dp, ul_ov_t *c_idx, uint64_t *buf, asg64_v *ovidx)
{
#define id_mm ((uint64_t)0x7fffffffffffffff)
#define id_set ((uint64_t)0x8000000000000000)
if(!id_n) return id_n;
uint64_t k, m, mn, q[2], buf_n, rm_n, i; int64_t err, msc, msc_k, msc_n;
uint64_t k, m, mn, q[2], buf_n, rm_n, i, oid; int64_t err, msc, msc_k, msc_n;
overlap_region *z; ul_ov_t *p;
for (k = buf_n = rm_n = 0; k < id_n; k++) {
p = &(c_idx[id_a[k]]);
@@ -15030,25 +15032,26 @@ uint64_t gen_region_phase(overlap_region* ol, uint64_t *id_a, uint64_t id_n, uin
}
// fprintf(stderr, "[M::%s] buf_n::%ld, msc_n::%ld, mn::%lu\n", __func__, buf_n, msc_n, mn);
for (k = 0; k < buf_n; k++) {
buf[k] = ovlp_id((c_idx[(uint32_t)buf[k]]));
oid = ovlp_id((c_idx[(uint32_t)buf[k]]));
buf[k] >>= 32; buf[k] <<= 32; buf[k] |= oid;
// if(k < mn) fprintf(stderr, "d::bst::[M::%s::utg%.6dl]\n", __func__, (int32_t)ol[buf[k]].y_id+1);
}
for (k = 0; k < mn; k++) {
if(buf[k]&id_set) continue;
for (i = 0; i < k; i++) {
if(query_gen_gov_idx(ovidx, buf[k]&id_mm, buf[i]&id_mm)) break;
if(query_gen_gov_idx(ovidx, id_get(buf[k]), id_get(buf[i]))) break;
}
if(i < k) {
if(!(buf[k]&id_set)) {
ol[buf[k]&id_mm].overlapLen -= e - s;
ol[buf[k]&id_mm].align_length -= e - s;
ol[id_get(buf[k])].overlapLen -= e - s;
ol[id_get(buf[k])].align_length -= e - s;
}
buf[k] |= id_set;
if(!(buf[i]&id_set)) {
ol[buf[i]&id_mm].overlapLen -= e - s;
ol[buf[i]&id_mm].align_length -= e - s;
ol[id_get(buf[i])].overlapLen -= e - s;
ol[id_get(buf[i])].align_length -= e - s;
}
buf[i] |= id_set;
}
@@ -15057,24 +15060,27 @@ uint64_t gen_region_phase(overlap_region* ol, uint64_t *id_a, uint64_t id_n, uin
for (k = mn; k < buf_n; k++) {
z = &(ol[buf[k]&id_mm]);
z = &(ol[id_get(buf[k])]);
z->align_length -= e - s;
for (i = 0; i < mn; i++) {
if(query_gen_gov_idx(ovidx, buf[k]&id_mm, buf[i]&id_mm)) break;
if(query_gen_gov_idx(ovidx, id_get(buf[k]), id_get(buf[i]))) break;
}
if(i < mn) {
z->align_length += e - s;
if(!(buf[k]&id_set)) {
ol[buf[k]&id_mm].overlapLen -= e - s;
ol[buf[k]&id_mm].align_length -= e - s;
ol[id_get(buf[k])].overlapLen -= e - s;
ol[id_get(buf[k])].align_length -= e - s;
}
buf[k] |= id_set;
if(!(buf[i]&id_set)) {
ol[buf[i]&id_mm].overlapLen -= e - s;
ol[buf[i]&id_mm].align_length -= e - s;
ol[id_get(buf[i])].overlapLen -= e - s;
ol[id_get(buf[i])].align_length -= e - s;
}
buf[i] |= id_set;
// fprintf(stderr, "i::bst::[M::%s::utg%.6dl]\n", __func__, (int32_t)ol[buf[k]].y_id+1);
} else {
assert(err_get(buf[i]) > (uint64_t)msc);
z->non_homopolymer_errors += err_get(buf[i])-msc;
}
}
}
@@ -15287,7 +15293,7 @@ void region_phase(overlap_region_alloc* ol, const ul_idx_t *uref, const ug_opt_t
kv_resize(ul_ov_t, *c_idx, ol->length);
for (k = idx->n = c_idx->n = 0; k < on; k++) {
z = &(ol->list[k]); zwn = z->w_list.n;
z->align_length = z->overlapLen = z->x_pos_e+1-z->x_pos_s;
z->align_length = z->overlapLen = z->x_pos_e+1-z->x_pos_s; z->non_homopolymer_errors = 0;
if(!zwn) continue;
q[0] = q[1] = t[0] = t[1] = w[0] = w[1] = INT32_MIN;
for (i = 0; i < zwn; i++) {

13
inter.cpp
View File

@@ -5073,17 +5073,17 @@ int64_t debug_i, int64_t tid, void *km)
}
// #define aln_sc(a, w) (((int64_t)((a).sec))-((int64_t)(((a).qe-(a).qs-(a).sec)*(w))))
#define aln_sc(a, w) (((int64_t)((a).align_length))-((int64_t)(((a).overlapLen-(a).align_length)*(w))))
#define aln_sc(a, w) (((int64_t)((a).align_length))-((int64_t)(((a).overlapLen-(a).align_length)*(w)))-(((int64_t)((a).non_homopolymer_errors))*UG_TRANS_ERR_W))
void gen_gl_aln(overlap_region_alloc* olist, const ul_idx_t *uref, kv_ul_ov_t *res)
{
uint64_t k; ul_ov_t *p = NULL;
res->n = 0; kv_resize(ul_ov_t, *res, olist->length);
for (k = 0; k < olist->length; k++) {
// fprintf(stderr, "+++[M::%s::utg%.6dl] q[%d, %d), t[%d, %d), tot::%u, cis::%u\n", __func__,
// fprintf(stderr, "+++[M::%s::utg%.6dl] q[%d, %d), t[%d, %d), tot::%u, cis::%u, err::%u\n", __func__,
// (int32_t)olist->list[k].y_id+1, olist->list[k].x_pos_s, olist->list[k].x_pos_e+1,
// olist->list[k].y_pos_s, olist->list[k].y_pos_e+1,
// olist->list[k].overlapLen, olist->list[k].align_length);
// olist->list[k].overlapLen, olist->list[k].align_length, olist->list[k].non_homopolymer_errors);
p = &(res->a[res->n++]); assert(olist->list[k].overlapLen >= olist->list[k].align_length);
p->qn = k; p->qs = olist->list[k].x_pos_s; p->qe = olist->list[k].x_pos_e+1;
p->tn = olist->list[k].y_id; p->sec = olist->list[k].align_length;
@@ -6001,8 +6001,10 @@ uint32_t ck_w_err(overlap_region *z, uint64_t *w_idx, int64_t wl, int64_t ql)
// fprintf(stderr, "[M::%s::utg%.6dl] x::[%u, %u), ol::%ld, e[0]::%ld, e[1]::%ld\n",
// __func__, (int32_t)z->y_id+1, z->x_pos_s, z->x_pos_e+1, ol, e[0], e[1]);
if(e[1] > (e[0]+64)) {
if(e[1] > (e[0]+(ol*0.01))) return 0;
if(e[1] > (e[0]+(e[0]*0.03))) return 0;
if((e[1] > (e[0]+(ol*0.01)))||(e[1] > (e[0]+(e[0]*0.03)))) {
z->non_homopolymer_errors = e[1] - e[0];
return 0;
}
}
// if((e[1] > (e[0]+16)) && (e[1] > (e[0]+(ol*0.01)))) return 0;
return 1;
@@ -6667,6 +6669,7 @@ static void worker_for_ul_rescall_alignment(void *data, long i, int tid) // call
for (k = 0; k < b->olist.length; k++) {
b->olist.list[k].align_length = b->olist.list[k].overlapLen =
b->olist.list[k].x_pos_e+1-b->olist.list[k].x_pos_s;
b->olist.list[k].non_homopolymer_errors = 0;
}
}

1
inter.h
View File

@@ -18,6 +18,7 @@
#define UG_ITER_N 5000
#define UG_DIS_N 50000
#define UG_TRANS_W 2
#define UG_TRANS_ERR_W 512
#define G_CHAIN_TRANS_RATE 0.25
#define G_CHAIN_TRANS_WEIGHT -1
#define G_CHAIN_INDEL 128