Package: samtools
Version: 0.1.18-1
Severity: important
Tag: patch

If a read is aligned outside of the reference (possibly due to
reference mismatch or truncation), samtools mpileup will segfault. It
should instead warn, and continue on. The attached patch fixes this.


Don Armstrong

-- 
"Ban cryptography! Yes. Let's also ban pencils, pens and paper, since
criminals can use them to draw plans of the joint they are casing or
even, god forbid, create one time pads to pass uncrackable codes to
each other. Ban open spaces since criminals could use them to converse
with each other out of earshot of the police. Let's ban flags since
they could be used to pass secret messages in semaphore. In fact let's
just ban all forms of verbal and non-verbal communication -- let's see
those criminals make plans now!"

http://www.donarmstrong.com              http://rzlab.ucr.edu
Description: Fix segfault if position is outside of the reference sequence length
Origin: Don Armstrong <d...@donarmstrong.com>

Index: samtools-0.1.18/bam_plcmd.c
===================================================================
--- samtools-0.1.18.orig/bam_plcmd.c	2011-07-06 20:41:26.000000000 -0700
+++ samtools-0.1.18/bam_plcmd.c	2011-12-22 14:44:05.000000000 -0800
@@ -93,6 +93,7 @@
 	bam_iter_t iter;
 	bam_header_t *h;
 	int ref_id;
+        int ref_len;
 	char *ref;
 	const mplp_conf_t *conf;
 } mplp_aux_t;
@@ -134,6 +135,11 @@
 				qual[i] = qual[i] > 31? qual[i] - 31 : 0;
 		}
 		has_ref = (ma->ref && ma->ref_id == b->core.tid)? 1 : 0;
+		if (has_ref && (ma->ref_len <= b->core.pos)) { // exclude reads outside of the reference sequence
+		  fprintf(stderr,"[%s] Skipping because %d is outside of %d [ref:%s]\n",__func__,b->core.pos,ma->ref_len,ma->ref_id);
+		  skip = 1;
+		  continue;
+		}
 		skip = 0;
 		if (has_ref && (ma->conf->flag&MPLP_REALN)) bam_prob_realn_core(b, ma->ref, (ma->conf->flag & MPLP_EXT_BAQ)? 3 : 1);
 		if (has_ref && ma->conf->capQ_thres > 10) {
@@ -277,7 +283,11 @@
 	if (tid0 >= 0 && conf->fai) { // region is set
 		ref = faidx_fetch_seq(conf->fai, h->target_name[tid0], 0, 0x7fffffff, &ref_len);
 		ref_tid = tid0;
-		for (i = 0; i < n; ++i) data[i]->ref = ref, data[i]->ref_id = tid0;
+		for (i = 0; i < n; ++i) {
+		  data[i]->ref = ref;
+                  data[i]->ref_id = tid0;
+		  data[i]->ref_len = ref_len;
+                }
 	} else ref_tid = -1, ref = 0;
 	iter = bam_mplp_init(n, mplp_func, (void**)data);
 	max_depth = conf->max_depth;
@@ -295,7 +305,11 @@
 		if (tid != ref_tid) {
 			free(ref); ref = 0;
 			if (conf->fai) ref = faidx_fetch_seq(conf->fai, h->target_name[tid], 0, 0x7fffffff, &ref_len);
-			for (i = 0; i < n; ++i) data[i]->ref = ref, data[i]->ref_id = tid;
+			for (i = 0; i < n; ++i) {
+			  data[i]->ref = ref;
+			  data[i]->ref_id = tid;
+			  data[i]->ref_len = ref_len;
+			}
 			ref_tid = tid;
 		}
 		if (conf->flag & MPLP_GLF) {

Reply via email to