X-Git-Url: http://woldlab.caltech.edu/gitweb/?p=samtools.git;a=blobdiff_plain;f=bcftools%2Fvcfutils.pl;h=6cc168fb363f59541511429b6e6538e0263acbe1;hp=bbc479bfffa47aa2b1089b5b6e56414b696be413;hb=6a0c6f060a60789b48e10a72b1381f6e54599302;hpb=e3b3a0177339fb8c099346986e965e3bd5b85999

diff --git a/bcftools/vcfutils.pl b/bcftools/vcfutils.pl
index bbc479b..6cc168f 100755
--- a/bcftools/vcfutils.pl
+++ b/bcftools/vcfutils.pl
@@ -231,13 +231,17 @@ Note: Some of the filters rely on annotations generated by SAMtools\/BCFtools.
 	}
 	next if ($t[4] eq '.'); # skip non-var sites
 	# check if the site is a SNP
-	my $is_snp = 1;
+	my $type = 1; # SNP
 	if (length($t[3]) > 1) {
-	  $is_snp = 0;
+	  $type = 2; # MNP
+	  my @s = split(',', $t[4]);
+	  for (@s) {
+		$type = 3 if (length != length($t[3]));
+	  }
 	} else {
 	  my @s = split(',', $t[4]);
 	  for (@s) {
-		$is_snp = 0 if (length > 1);
+		$type = 3 if (length > 1);
 	  }
 	}
 	# clear the out-of-range elements
@@ -270,35 +274,42 @@ Note: Some of the filters rely on annotations generated by SAMtools\/BCFtools.
 	$flt = 7 if ($flt == 0 && /PV4=([^,]+),([^,]+),([^,]+),([^,;\t]+)/
 				 && ($1<$opts{1} || $2<$opts{2} || $3<$opts{3} || $4<$opts{4}));
 
-	# site dependent filters
-	my ($rlen, $indel_score) = (0, -1); # $indel_score<0 for SNPs
+	my $score = $t[5] * 100 + $dp_alt;
+	my $rlen = length($t[3]) - 1; # $indel_score<0 for SNPs
 	if ($flt == 0) {
-	  if (!$is_snp) { # an indel
-		$rlen = length($t[3]) - 1;
-		$indel_score = $t[5] * 100 + $dp_alt;
-		# filtering SNPs
+	  if ($type == 3) { # an indel
+		# filtering SNPs and MNPs
 		for my $x (@staging) {
-		  next if ($x->[0] >= 0 || $x->[1] || $x->[4] + $x->[2] + $ow < $t[1]);
+		  next if (($x->[0]&3) == 3 || $x->[1] || $x->[4] + $x->[2] + $ow < $t[1]);
 		  $x->[1] = 5;
 		}
 		# check the staging list for indel filtering
 		for my $x (@staging) {
-		  next if ($x->[0] < 0 || $x->[1] || $x->[4] + $x->[2] + $ol < $t[1]);
-		  if ($x->[0] < $indel_score) {
+		  next if (($x->[0]&3) != 3 || $x->[1] || $x->[4] + $x->[2] + $ol < $t[1]);
+		  if ($x->[0]>>2 < $score) {
 			$x->[1] = 6;
 		  } else {
 			$flt = 6; last;
 		  }
 		}
-	  } else { # a SNP
+	  } else { # SNP or MNP
 		for my $x (@staging) {
-		  next if ($x->[0] < 0 || $x->[1] || $x->[4] + $x->[2] + $ow < $t[1]);
+		  next if (($x->[0]&3) != 3 || $x->[4] + $x->[2] + $ow < $t[1]);
 		  $flt = 5;
 		  last;
 		}
+		# check MNP
+		for my $x (@staging) {
+		  next if (($x->[0]&3) == 3 || $x->[4] + $x->[2] < $t[1]);
+		  if ($x->[0]>>2 < $score) {
+			$x->[1] = 8;
+		  } else {
+			$flt = 8; last;
+		  }
+		}
 	  }
 	}
-	push(@staging, [$indel_score, $flt, $rlen, @t]);
+	push(@staging, [$score<<2|$type, $flt, $rlen, @t]);
   }
   # output the last few elements in the staging list
   while (@staging) {
@@ -311,7 +322,7 @@ sub varFilter_aux {
   if ($first->[1] == 0) {
 	print join("\t", @$first[3 .. @$first-1]), "\n";
   } elsif ($is_print) {
-	print STDERR join("\t", substr("UQdDaGgP", $first->[1], 1), @$first[3 .. @$first-1]), "\n";
+	print STDERR join("\t", substr("UQdDaGgPM", $first->[1], 1), @$first[3 .. @$first-1]), "\n";
   }
 }