Merge branch 'master' into dev-lh3

This commit is contained in:
Heng Li
2020-04-12 18:57:42 -04:00
2 changed files with 3 additions and 263 deletions
-262
View File
@@ -7654,257 +7654,6 @@ ma_hit_t_alloc* reverse_sources, long long miniedgeLen, R_to_U* ruIndex)
}
int asg_arc_del_tri_link(asg_t *g, int max_dist)
{
double startTime = Get_T();
kvec_t(uint64_t) b;
memset(&b, 0, sizeof(b));
kvec_t(uint32_t) b_f;
memset(&b_f, 0, sizeof(b_f));
kvec_t(uint32_t) b_r;
memset(&b_r, 0, sizeof(b_r));
uint32_t v, w, n_vtx = g->n_seq * 2, n_cut = 0;
buf_t bub;
if (!g->is_symm) asg_symm(g);
memset(&bub, 0, sizeof(buf_t));
bub.a = (binfo_t*)calloc(n_vtx, sizeof(binfo_t));
uint32_t Ns_first[3];
uint32_t Ns_second[3];
long long NodeLen_first[3];
long long NodeLen_second[3];
for (v = 0; v < n_vtx; ++v)
{
if(g->seq_vis[v] == 0)
{
asg_arc_t *av = asg_arc_a(g, v);
uint32_t nv = asg_arc_n(g, v);
if(nv != 2) continue;
if(av[0].v == av[1].v) continue;
/**********************test first node************************/
NodeLen_first[0] = NodeLen_first[1] = NodeLen_first[2] = -1;
if(asg_is_single_edge(g, av[0].v, v>>1) <= 2 && asg_is_single_edge(g, av[1].v, v>>1) <= 2)
{
NodeLen_first[asg_is_single_edge(g, av[0].v, v>>1)] = 0;
NodeLen_first[asg_is_single_edge(g, av[1].v, v>>1)] = 1;
}
///one node has one out-edge, another node has two out-edges
if(NodeLen_first[1] == -1 || NodeLen_first[2] == -1)
{
continue;
}
/**********************test first node************************/
/**********************test second node************************/
w = av[NodeLen_first[2]].v^1;
asg_arc_t *aw = asg_arc_a(g, w);
uint32_t nw = asg_arc_n(g, w);
if(nw != 2)
{
fprintf(stderr, "error\n");
}
NodeLen_second[0] = NodeLen_second[1] = NodeLen_second[2] = -1;
if(asg_is_single_edge(g, aw[0].v, w>>1) <= 2 && asg_is_single_edge(g, aw[1].v, w>>1) <= 2)
{
NodeLen_second[asg_is_single_edge(g, aw[0].v, w>>1)] = 0;
NodeLen_second[asg_is_single_edge(g, aw[1].v, w>>1)] = 1;
}
///one node has one out-edge, another node has two out-edges
if(NodeLen_second[1] == -1 || NodeLen_second[2] == -1)
{
continue;
}
/**********************test second node************************/
uint64_t t_ol = av[NodeLen_first[2]].ol;
kv_push(uint64_t, b, (uint64_t)(t_ol << 32 | v));
}
}
radix_sort_arch64(b.a, b.a + b.n);
uint64_t k;
for (k = 0; k < b.n; k++)
{
///v is the node
v = (uint32_t)b.a[k];
if (g->seq[v>>1].del) continue;
uint32_t nv = asg_arc_n(g, v), nw, to_del;
uint32_t kv = get_real_length(g, v, NULL), kw;
///at the begining, the nv of all nodes must be == 2;
///here kv == 2, that means all edges are kept
///so we can use normal method to delete edges
if (nv != 2) continue;
if (kv != 2) continue;
asg_arc_t *av = asg_arc_a(g, v), *aw;
if(av[0].v == av[1].v)
{
continue;
}
/**********************test first node************************/
NodeLen_first[0] = NodeLen_first[1] = NodeLen_first[2] = -1;
if(asg_is_single_edge(g, av[0].v, v>>1) <= 2 && asg_is_single_edge(g, av[1].v, v>>1) <= 2)
{
NodeLen_first[asg_is_single_edge(g, av[0].v, v>>1)] = 0;
NodeLen_first[asg_is_single_edge(g, av[1].v, v>>1)] = 1;
}
///one node has one out-edge, another node has two out-edges
if(NodeLen_first[1] == -1 || NodeLen_first[2] == -1)
{
continue;
}
/**********************test first node************************/
/**********************test second node************************/
w = av[NodeLen_first[2]].v^1;
aw = asg_arc_a(g, w);
nw = asg_arc_n(g, w);
kw = get_real_length(g, w, NULL);
///at the begining, the nw of all nodes must be == 2;
///here kw == 2, that means all edges are kept
///so we can use normal method to delete edges
if(nw != 2) continue;
if(kw != 2) continue;
NodeLen_second[0] = NodeLen_second[1] = NodeLen_second[2] = -1;
if(asg_is_single_edge(g, aw[0].v, w>>1) <= 2 && asg_is_single_edge(g, aw[1].v, w>>1) <= 2)
{
NodeLen_second[asg_is_single_edge(g, aw[0].v, w>>1)] = 0;
NodeLen_second[asg_is_single_edge(g, aw[1].v, w>>1)] = 1;
}
///one node has one out-edge, another node has two out-edges
if(NodeLen_second[1] == -1 || NodeLen_second[2] == -1)
{
continue;
}
/**********************test second node************************/
uint32_t convex1, convex2, f1, f2;
long long l1, l2;
to_del = 0;
f1 = detect_bubble_end_with_bubbles(g, av[0].v, av[1].v, &convex1, &l1, NULL);
f2 = detect_bubble_end_with_bubbles(g, aw[0].v, aw[1].v, &convex2, &l2, NULL);
if(f1 && f2)
{
if((l1 <= min_thres) || (l2 <= min_thres))
{
continue;
}
to_del = 1;
}
else if(f1)
{
if(l1 <= min_thres)
{
continue;
}
to_del = 1;
}
else if(f2)
{
if(l2 <= min_thres)
{
continue;
}
to_del = 1;
}
///here both f1 and f2 must be false
if(to_del == 0)
{
if(!f1)
{
Ns_first[0] = av[0].v; Ns_first[1] = av[1].v;
f1 = asg_bub_end_finder_with_del_advance(g, Ns_first, 2, max_dist,
&bub, 0, (u_int32_t)-1, &convex1);
l1 = min_thres + 10;
}
if(!f2)
{
Ns_second[0] = aw[0].v; Ns_second[1] = aw[1].v;
f2 = asg_bub_end_finder_with_del_advance(g, Ns_second, 2, max_dist,
&bub, 0, (u_int32_t)-1, &convex2);
l2 = min_thres + 10;
}
///since both f1 and f2 must be false before 'if(to_del == 0)'
///here l1 = min_thres + 10, l2 = min_thres + 10
if(f1 && f2)
{
// if(l1 <= min_thres || l2 <= min_thres) // TODO: check this block: it is always false
// {
// continue;
// }
to_del = 1;
}
else if(f1)
{
if(l1 <= min_thres)
{
continue;
}
to_del = 1;
}
else if(f2)
{
if(l2 <= min_thres)
{
continue;
}
to_del = 1;
}
}
if (to_del)
{
++n_cut;
av[NodeLen_first[2]].del = 1;
asg_arc_del(g, av[NodeLen_first[2]].v^1, av[NodeLen_first[2]].ul>>32^1, 1);
}
}
free(b.a); free(b_f.a); free(b_r.a);
free(bub.a); free(bub.S.a); free(bub.T.a); free(bub.b.a); free(bub.e.a);
if (n_cut)
{
asg_cleanup(g);
asg_symm(g);
}
fprintf(stderr, "[M::%s] removed %d false overlaps\n", __func__, n_cut);
fprintf(stderr, "[M::%s] takes %0.2f s\n\n", __func__, Get_T()-startTime);
return n_cut;
}
int asg_arc_del_complex_false_link(asg_t *g, float drop_ratio, float o_drop_ratio, int max_dist,
ma_hit_t_alloc* reverse_sources, long long miniedgeLen)
{
@@ -13454,11 +13203,6 @@ kvec_asg_arc_t_warp* new_rtg_edges)
delete_useless_nodes(ug);
// purge_dups(*ug, read_g, coverage_cut, reverse_sources, ruIndex, new_rtg_edges, 0.75, 50, 50, 0.5, max_hang,
// min_ovlp, bubble_dist, drop_ratio, 1);
// deduplicate(*ug, read_g, reverse_sources, 20, 100, 0.05, 0.2, ruIndex, 1);
// delete_useless_nodes(ug);
update_unitig_graph((*ug), read_g, reverse_sources, ruIndex, flag, drop_rate);
renew_utg(ug, read_g, new_rtg_edges);
@@ -22264,19 +22008,13 @@ kvec_asg_arc_t_warp* new_rtg_edges)
renew_utg(ug, read_g, new_rtg_edges);
///enable_debug_mode(1);
purge_dups(*ug, read_g, coverage_cut, reverse_sources, ruIndex, new_rtg_edges, 0.75, 50, 50, 0.5, max_hang,
min_ovlp, bubble_dist, drop_ratio, 0);
delete_useless_nodes(ug);
renew_utg(ug, read_g, new_rtg_edges);
rescue_missing_overlaps_aggressive(*ug, read_g, sources, coverage_cut, ruIndex, max_hang,
min_ovlp, 0, 0, 1, NULL);
renew_utg(ug, read_g, new_rtg_edges);
+3 -1
View File
@@ -1599,6 +1599,7 @@ hap_overlaps_list* all_ovlp)
asg_arc_t t;
asg_arc_t_offset t_offset;
hap_candidates hap_can;
memset(&hap_can, 0, sizeof(hap_candidates));
hap_overlaps hap_align;
xUid = Input_uId;
if(nsg->seq[xUid].del || nsg->seq[xUid].c == ALTER_LABLE) return;
@@ -2161,7 +2162,8 @@ void hap_alignment_worker(void *_data, long eid, int tid)
int32_t r;
asg_arc_t t;
asg_arc_t_offset t_offset;
hap_candidates hap_can;
hap_candidates hap_can;
memset(&hap_can, 0, sizeof(hap_candidates));
hap_overlaps hap_align;
xUid = Input_uId;
if(nsg->seq[xUid].del || nsg->seq[xUid].c == ALTER_LABLE) return;