-
Notifications
You must be signed in to change notification settings - Fork 3
/
Copy pathoptions.c
154 lines (147 loc) · 4.59 KB
/
options.c
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
#include <string.h>
#include "mgpriv.h"
#include "sys.h"
void mg_idxopt_init(mg_idxopt_t *io)
{
memset(io, 0, sizeof(mg_idxopt_t));
io->k = 17;
io->w = 11;
io->bucket_bits = 14;
}
void mg_mapopt_init(mg_mapopt_t *mo)
{
memset(mo, 0, sizeof(mg_mapopt_t));
mo->seed = 11;
mo->occ_max1 = 50, mo->occ_max1_cap = 250;
mo->occ_max1_frac = 2e-4f;
mo->max_gap = 5000;
mo->max_gap_ref = -1;
mo->max_gap_pre = 1000;
mo->max_lc_skip = 25, mo->max_gc_skip = 25;
mo->max_lc_iter = 5000;
mo->bw = 500, mo->bw_long = 20000;
mo->rmq_size_cap = 100000;
mo->rmq_rescue_size = 1000;
mo->rmq_rescue_ratio = 0.1f;
mo->mini_batch_size = 500000000;
mo->div = 0.1f;
mo->l_chn_pen_gap = 0.0f, mo->l_chn_pen_skip = 0.0f; // Added for linear chaining
mo->chn_pen_gap = 1.0f, mo->chn_pen_skip = 0.05f;
mo->min_lc_cnt = 1, mo->min_lc_score = 0;
mo->min_gc_cnt = 5, mo->min_gc_score = 50;
mo->gdp_max_ed = 10000;
mo->lc_max_trim = 50;
mo->lc_max_occ = 2;
mo->mask_level = 0.5f;
mo->sub_diff = 6;
mo->best_n = 5;
mo->pri_ratio = 0.8f;
mo->ref_bonus = 0;
mo->pe_ori = 0; // FF
mo->min_cov_mapq = 20;
mo->min_cov_blen = 1000;
mo->cap_kalloc = 1000000000;
// tau_1 : threshold to compute "intra cid" disjoint set of chains
mo->tau_1 = 0.99f;
// tau_2 : threshold to pick "inter cid" set of chains
mo->tau_2 = 0.95f;
mo->is_ggen = 0;
mo->G = 5000; // for long read mapping
mo->max_itr = 100;
}
void mg_ggopt_init(mg_ggopt_t *go)
{
memset(go, 0, sizeof(mg_ggopt_t));
go->algo = MG_G_NONE;
go->flag |= MG_G_NO_QOVLP;
go->min_map_len = 100000;
go->min_depth_len = 20000;
go->min_mapq = 5;
go->min_var_len = 50;
go->match_pen = 10;
// for ggs
go->ggs_shrink_pen = 9;
go->ggs_min_end_cnt = 10;
go->ggs_min_end_frac = 0.1f;
go->ggs_max_iden = 0.80f;
go->ggs_min_inv_iden = 0.95f;
}
int mg_opt_set(const char *preset, mg_idxopt_t *io, mg_mapopt_t *mo, mg_ggopt_t *go)
{
if (preset == 0) {
mg_idxopt_init(io);
mg_mapopt_init(mo);
mg_ggopt_init(go);
} else if (strcmp(preset, "lr") == 0) { // this is the default
} else if (strcmp(preset, "asm") == 0 || strcmp(preset, "ggs") == 0) {
io->k = 19, io->w = 10;
mo->flag |= MG_M_RMQ;
mo->occ_max1 = 10, mo->occ_max1_cap = 100;
mo->bw = 1000, mo->bw_long = 150000;
mo->max_gap = 150000, mo->max_gap_pre = 1000;
mo->min_lc_cnt = 1, mo->min_lc_score = 0;
mo->min_gc_cnt = 5, mo->min_gc_score = 1000;
mo->min_cov_mapq = 5;
mo->min_cov_blen = 100000;
mo->max_lc_skip = mo->max_gc_skip = 50;
mo->div = 0.01f;
mo->mini_batch_size = 4000000000LL;
// tau_1 : threshold to compute "intra cid" disjoint set of chains
mo->tau_1 = 0.0f;
// tau_2 : threshold to pick "inter cid" set of chains
mo->tau_2 = 0.0f;
mo->is_ggen = 1;
mo->G = 150000; // for graph generation
mo->max_itr = 500; // works reasonably well
if (strcmp(preset, "ggs") == 0)
go->algo = MG_G_GGSIMPLE, mo->best_n = 0;
} else if (strcmp(preset, "se") == 0 || strcmp(preset, "sr") == 0) {
io->k = 21, io->w = 10;
mo->flag |= MG_M_SR | MG_M_HEAP_SORT | MG_M_2_IO_THREADS;
mo->occ_max1 = 1000;
mo->occ_max1_cap = 2500;
mo->max_gap = 100;
mo->bw = mo->bw_long = 100;
mo->max_frag_len = 800;
mo->pri_ratio = 0.5f;
mo->min_lc_cnt = 2, mo->min_lc_score = 0;
mo->min_gc_cnt = 3, mo->min_gc_score = 40;
mo->mini_batch_size = 50000000;
mo->min_cov_blen = 50;
mo->chn_pen_gap = 0.2f;
mo->ref_bonus = 1;
// tau_1 : threshold to compute "intra cid" disjoint set of chains
mo->tau_1 = 0.99;
// tau_2 : threshold to pick "inter cid" set of chains
mo->tau_2 = 0.95;
if (strcmp(preset, "sr") == 0) {
mo->flag |= MG_M_FRAG_MODE | MG_M_FRAG_MERGE;
mo->pe_ori = 0<<1|1; // FR
}
} else return -1;
return 0;
}
int mg_opt_check(const mg_idxopt_t *io, const mg_mapopt_t *mo, const mg_ggopt_t *go)
{
if ((mo->flag & MG_M_FRAG_MODE) && !(mo->flag & MG_M_FRAG_MERGE)) {
if (mg_verbose >= 1)
fprintf(stderr, "[ERROR]\033[1;31m the fragment-without-merge mode is not implemented\033[0m\n");
return -1;
}
return 0;
}
void mg_opt_update(const mg_idx_t *gi, mg_mapopt_t *mo, mg_ggopt_t *go)
{
float f[2];
int32_t q[2];
f[0] = 0.1f, f[1] = mo->occ_max1_frac;
mg_idx_cal_quantile(gi, 2, f, q);
if (q[0] > mo->lc_max_occ) mo->lc_max_occ = q[0];
if (mo->lc_max_occ > mo->occ_max1_cap) mo->lc_max_occ = mo->occ_max1_cap;
if (q[1] > mo->occ_max1) mo->occ_max1 = q[1];
if (mo->occ_max1 > mo->occ_max1_cap) mo->occ_max1 = mo->occ_max1_cap;
if (mo->bw_long < mo->bw) mo->bw_long = mo->bw;
if (mg_verbose >= 3)
fprintf(stderr, "[M::%s::%.3f*%.2f] occ_max1=%d; lc_max_occ=%d\n", __func__,
realtime() - mg_realtime0, cputime() / (realtime() - mg_realtime0), mo->occ_max1, mo->lc_max_occ);
}