From 5fbd4546829d552eec52b1e45ea9ea86b7d5419c Mon Sep 17 00:00:00 2001
From: Heng Li <lh3@me.com>
Date: Tue, 5 Mar 2013 22:49:38 -0500
Subject: [PATCH] r332: added output threshold

Otherwise there are far too many short hits
---
 bwamem.c  | 4 +++-
 bwamem.h  | 1 +
 fastmap.c | 4 +++-
 main.c    | 2 +-
 4 files changed, 8 insertions(+), 3 deletions(-)

diff --git a/bwamem.c b/bwamem.c
index cd3a1f5..7a6d175 100644
--- a/bwamem.c
+++ b/bwamem.c
@@ -44,6 +44,7 @@ mem_opt_t *mem_opt_init()
 	o = calloc(1, sizeof(mem_opt_t));
 	o->flag = 0;
 	o->a = 1; o->b = 4; o->q = 6; o->r = 1; o->w = 100;
+	o->T = 30;
 	o->zdrop = 100;
 	o->pen_unpaired = 9;
 	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;
 	kstring_t str;
 	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;
 		for (k = 0; k < a->n; ++k) {
 			bwahit_t h;
 			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 && p->score < a->a[p->secondary].score * .5) continue;
 			mem_alnreg2hit(p, &h);
diff --git a/bwamem.h b/bwamem.h
index 96a3308..14f04d5 100644
--- a/bwamem.h
+++ b/bwamem.h
@@ -24,6 +24,7 @@ typedef struct {
 	int w;                  // band width
 	int zdrop;              // Z-dropoff
 
+	int T;                  // output score threshold; only affecting output
 	int flag;               // see MEM_F_* macros
 	int min_seed_len;       // minimum seed length
 	float split_factor;     // split into a seed if MEM is longer than min_seed_len*split_factor
diff --git a/fastmap.c b/fastmap.c
index d4e5626..e6c2b1e 100644
--- a/fastmap.c
+++ b/fastmap.c
@@ -26,13 +26,14 @@ int main_mem(int argc, char *argv[])
 	void *ko = 0, *ko2 = 0;
 
 	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);
 		else if (c == 'w') opt->w = atoi(optarg);
 		else if (c == 'A') opt->a = atoi(optarg);
 		else if (c == 'B') opt->b = atoi(optarg);
 		else if (c == 'O') opt->q = 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 == '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;
@@ -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, "\n");
 		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, "       -C         append FASTA/FASTQ comment to SAM output\n");
 		fprintf(stderr, "       -H         hard clipping\n");
diff --git a/main.c b/main.c
index 0954e01..c342fd4 100644
--- a/main.c
+++ b/main.c
@@ -4,7 +4,7 @@
 #include "utils.h"
 
 #ifndef PACKAGE_VERSION
-#define PACKAGE_VERSION "0.7.0-r331-beta"
+#define PACKAGE_VERSION "0.7.0-r332-beta"
 #endif
 
 static int usage()