zzh
/
fast-bwa
1
0
Fork 0
Code Issues Pull Requests Packages Projects Releases Wiki Activity

extend more seeds (and thus slower...)

This commit is contained in:
Heng Li 2013-02-21 12:52:00 -05:00
parent f8829318cf
commit 54da54ffd4
3 changed files with 9 additions and 3 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

4
bwamem.c
View File

@ -62,6 +62,7 @@ mem_opt_t *mem_opt_init()
o->max_ins = 10000; o->max_ins = 10000;
o->mask_level = 0.50; o->mask_level = 0.50;
o->chain_drop_ratio = 0.50; o->chain_drop_ratio = 0.50;
o->split_factor = 1.5;
o->chunk_size = 10000000; o->chunk_size = 10000000;
o->n_threads = 1; o->n_threads = 1;
o->pe_dir = 0<<1|1; o->pe_dir = 0<<1|1;
@ -186,7 +187,8 @@ static int test_and_merge(const mem_opt_t *opt, mem_chain_t *c, const mem_seed_t
static void mem_insert_seed(const mem_opt_t *opt, kbtree_t(chn) *tree, smem_i *itr) static void mem_insert_seed(const mem_opt_t *opt, kbtree_t(chn) *tree, smem_i *itr)
{ {
const bwtintv_v *a; const bwtintv_v *a;
while ((a = smem_next(itr, opt->min_seed_len<<1, opt->split_width)) != 0) { // to find all SMEM and some internal MEM int split_len = (int)(opt->min_seed_len * opt->split_factor + .499);
while ((a = smem_next(itr, split_len, opt->split_width)) != 0) { // to find all SMEM and some internal MEM
int i; int i;
for (i = 0; i < a->n; ++i) { // go through each SMEM/MEM up to itr->start for (i = 0; i < a->n; ++i) { // go through each SMEM/MEM up to itr->start
bwtintv_t *p = &a->a[i]; bwtintv_t *p = &a->a[i];

4
bwamem.h
View File

@ -27,7 +27,9 @@ typedef struct {
int min_seed_len, max_occ, max_chain_gap; int min_seed_len, max_occ, max_chain_gap;
int n_threads, chunk_size; int n_threads, chunk_size;
int pe_dir; int pe_dir;
float mask_level, chain_drop_ratio; float mask_level;
float chain_drop_ratio;
float split_factor; // split into a seed if MEM is longer than min_seed_len*split_factor
int pen_unpaired; // phred-scaled penalty for unpaired reads int pen_unpaired; // phred-scaled penalty for unpaired reads
int max_ins; // maximum insert size int max_ins; // maximum insert size
int8_t mat[25]; // scoring matrix; mat[0] == 0 if unset int8_t mat[25]; // scoring matrix; mat[0] == 0 if unset

4
fastmap.c
View File

@ -24,12 +24,13 @@ int main_mem(int argc, char *argv[])
bseq1_t *seqs; bseq1_t *seqs;
opt = mem_opt_init(); opt = mem_opt_init();
while ((c = getopt(argc, argv, "PHk:c:v:s:")) >= 0) { while ((c = getopt(argc, argv, "PHk:c:v:s:r:")) >= 0) {
if (c == 'k') opt->min_seed_len = atoi(optarg); if (c == 'k') opt->min_seed_len = atoi(optarg);
else if (c == 'P') opt->flag |= MEM_F_NOPAIRING; else if (c == 'P') opt->flag |= MEM_F_NOPAIRING;
else if (c == 'H') opt->flag |= MEM_F_HARDCLIP; else if (c == 'H') opt->flag |= MEM_F_HARDCLIP;
else if (c == 'c') opt->max_occ = atoi(optarg); else if (c == 'c') opt->max_occ = atoi(optarg);
else if (c == 'v') mem_verbose = atoi(optarg); else if (c == 'v') mem_verbose = atoi(optarg);
else if (c == 'r') opt->split_factor = atof(optarg);
else if (c == 's') opt->split_width = atoi(optarg); else if (c == 's') opt->split_width = atoi(optarg);
} }
if (optind + 1 >= argc) { if (optind + 1 >= argc) {
@ -38,6 +39,7 @@ int main_mem(int argc, char *argv[])
fprintf(stderr, "Options: -k INT minimum seed length [%d]\n", opt->min_seed_len); fprintf(stderr, "Options: -k INT minimum seed length [%d]\n", opt->min_seed_len);
fprintf(stderr, " -c INT skip seeds with more than INT occurrences [%d]\n", opt->max_occ); fprintf(stderr, " -c INT skip seeds with more than INT occurrences [%d]\n", opt->max_occ);
fprintf(stderr, " -s INT look for internal seeds inside a seed with less than INT occ [%d]\n", opt->split_width); fprintf(stderr, " -s INT look for internal seeds inside a seed with less than INT occ [%d]\n", opt->split_width);
fprintf(stderr, " -r FLOAT look for internal seeds inside a seed longer than {-k} * FLOAT [%g]\n", opt->split_factor);
fprintf(stderr, " -v INT verbose level [%d]\n", mem_verbose); fprintf(stderr, " -v INT verbose level [%d]\n", mem_verbose);
fprintf(stderr, "\n"); fprintf(stderr, "\n");
free(opt); free(opt);