r332: added output threshold

Otherwise there are far too many short hits
This commit is contained in:
Heng Li 2013-03-05 22:49:38 -05:00
parent 6476343a83
commit 5fbd454682
4 changed files with 8 additions and 3 deletions

View File

@ -44,6 +44,7 @@ mem_opt_t *mem_opt_init()
o = calloc(1, sizeof(mem_opt_t)); o = calloc(1, sizeof(mem_opt_t));
o->flag = 0; o->flag = 0;
o->a = 1; o->b = 4; o->q = 6; o->r = 1; o->w = 100; o->a = 1; o->b = 4; o->q = 6; o->r = 1; o->w = 100;
o->T = 30;
o->zdrop = 100; o->zdrop = 100;
o->pen_unpaired = 9; o->pen_unpaired = 9;
o->pen_clip = 5; o->pen_clip = 5;
@ -742,11 +743,12 @@ void mem_sam_se(const mem_opt_t *opt, const bntseq_t *bns, const uint8_t *pac, b
int k; int k;
kstring_t str; kstring_t str;
str.l = str.m = 0; str.s = 0; str.l = str.m = 0; str.s = 0;
if (a->n > 0) { if (a->n > 0 && a->a[0].score >= opt->T) {
int mapq0 = -1; int mapq0 = -1;
for (k = 0; k < a->n; ++k) { for (k = 0; k < a->n; ++k) {
bwahit_t h; bwahit_t h;
mem_alnreg_t *p = &a->a[k]; mem_alnreg_t *p = &a->a[k];
if (p->score < opt->T) continue;
if (p->secondary >= 0 && !(opt->flag&MEM_F_ALL)) continue; if (p->secondary >= 0 && !(opt->flag&MEM_F_ALL)) continue;
if (p->secondary >= 0 && p->score < a->a[p->secondary].score * .5) continue; if (p->secondary >= 0 && p->score < a->a[p->secondary].score * .5) continue;
mem_alnreg2hit(p, &h); mem_alnreg2hit(p, &h);

View File

@ -24,6 +24,7 @@ typedef struct {
int w; // band width int w; // band width
int zdrop; // Z-dropoff int zdrop; // Z-dropoff
int T; // output score threshold; only affecting output
int flag; // see MEM_F_* macros int flag; // see MEM_F_* macros
int min_seed_len; // minimum seed length int min_seed_len; // minimum seed length
float split_factor; // split into a seed if MEM is longer than min_seed_len*split_factor float split_factor; // split into a seed if MEM is longer than min_seed_len*split_factor

View File

@ -26,13 +26,14 @@ int main_mem(int argc, char *argv[])
void *ko = 0, *ko2 = 0; void *ko = 0, *ko2 = 0;
opt = mem_opt_init(); opt = mem_opt_init();
while ((c = getopt(argc, argv, "paMCPHk:c:v:s:r:t:R:A:B:O:E:U:w:L:d:")) >= 0) { while ((c = getopt(argc, argv, "paMCPHk:c:v:s:r:t:R:A:B:O:E:U:w:L:d:T:")) >= 0) {
if (c == 'k') opt->min_seed_len = atoi(optarg); if (c == 'k') opt->min_seed_len = atoi(optarg);
else if (c == 'w') opt->w = atoi(optarg); else if (c == 'w') opt->w = atoi(optarg);
else if (c == 'A') opt->a = atoi(optarg); else if (c == 'A') opt->a = atoi(optarg);
else if (c == 'B') opt->b = atoi(optarg); else if (c == 'B') opt->b = atoi(optarg);
else if (c == 'O') opt->q = atoi(optarg); else if (c == 'O') opt->q = atoi(optarg);
else if (c == 'E') opt->r = atoi(optarg); else if (c == 'E') opt->r = atoi(optarg);
else if (c == 'T') opt->T = atoi(optarg);
else if (c == 'L') opt->pen_clip = atoi(optarg); else if (c == 'L') opt->pen_clip = atoi(optarg);
else if (c == 'U') opt->pen_unpaired = atoi(optarg); else if (c == 'U') opt->pen_unpaired = atoi(optarg);
else if (c == 't') opt->n_threads = atoi(optarg), opt->n_threads = opt->n_threads > 1? opt->n_threads : 1; else if (c == 't') opt->n_threads = atoi(optarg), opt->n_threads = opt->n_threads > 1? opt->n_threads : 1;
@ -74,6 +75,7 @@ int main_mem(int argc, char *argv[])
fprintf(stderr, " -R STR read group header line such as '@RG\\tID:foo\\tSM:bar' [null]\n"); fprintf(stderr, " -R STR read group header line such as '@RG\\tID:foo\\tSM:bar' [null]\n");
fprintf(stderr, "\n"); fprintf(stderr, "\n");
fprintf(stderr, " -v INT verbose level: 1=error, 2=warning, 3=message, 4+=debugging [%d]\n", bwa_verbose); fprintf(stderr, " -v INT verbose level: 1=error, 2=warning, 3=message, 4+=debugging [%d]\n", bwa_verbose);
fprintf(stderr, " -T INT minimum score to output [%d]\n", opt->T);
fprintf(stderr, " -a output all alignments for SE or unpaired PE\n"); fprintf(stderr, " -a output all alignments for SE or unpaired PE\n");
fprintf(stderr, " -C append FASTA/FASTQ comment to SAM output\n"); fprintf(stderr, " -C append FASTA/FASTQ comment to SAM output\n");
fprintf(stderr, " -H hard clipping\n"); fprintf(stderr, " -H hard clipping\n");

2
main.c
View File

@ -4,7 +4,7 @@
#include "utils.h" #include "utils.h"
#ifndef PACKAGE_VERSION #ifndef PACKAGE_VERSION
#define PACKAGE_VERSION "0.7.0-r331-beta" #define PACKAGE_VERSION "0.7.0-r332-beta"
#endif #endif
static int usage() static int usage()