X-Git-Url: http://woldlab.caltech.edu/gitweb/?p=samtools.git;a=blobdiff_plain;f=misc%2Fsamtools.pl;h=d03c1c7f148dcca835946e7ccd9ee52f3b6dc00b;hp=ce8449de622884847003ba538a2cf60c57118c99;hb=HEAD;hpb=8d2494d1fb7cd0fa7c63be5ffba8dd1a11457522

diff --git a/misc/samtools.pl b/misc/samtools.pl
index ce8449d..d03c1c7 100755
--- a/misc/samtools.pl
+++ b/misc/samtools.pl
@@ -10,7 +10,7 @@ my $version = '0.3.3';
 &usage if (@ARGV < 1);
 
 my $command = shift(@ARGV);
-my %func = (showALEN=>\&showALEN, pileup2fq=>\&pileup2fq, varFilter=>\&varFilter,
+my %func = (showALEN=>\&showALEN, pileup2fq=>\&pileup2fq, varFilter=>\&varFilter, plp2vcf=>\&plp2vcf,
 			unique=>\&unique, uniqcmp=>\&uniqcmp, sra2hdr=>\&sra2hdr, sam2fq=>\&sam2fq);
 
 die("Unknown command \"$command\".\n") if (!defined($func{$command}));
@@ -473,6 +473,44 @@ sub uniqcmp_aux {
   close($fh);
 }
 
+sub plp2vcf {
+  while (<>) {
+	my @t = split;
+	next if ($t[3] eq '*/*');
+	if ($t[2] eq '*') { # indel
+	  my @s = split("/", $t[3]);
+	  my (@a, @b);
+	  my ($ref, $alt);
+	  for (@s) {
+		next if ($_ eq '*');
+		if (/^-/) {
+		  push(@a, 'N'.substr($_, 1));
+		  push(@b, 'N');
+		} elsif (/^\+/) {
+		  push(@a, 'N');
+		  push(@b, 'N'.substr($_, 1));
+		}
+	  }
+	  if ($a[0] && $a[1]) {
+		if (length($a[0]) < length($a[1])) {
+		  $ref = $a[1];
+		  $alt = ($b[0] . ('N' x (length($a[1]) - length($a[0])))) . ",$b[1]";
+		} elsif (length($a[0]) > length($a[1])) {
+		  $ref = $a[0];
+		  $alt = ($b[1] . ('N' x (length($a[0]) - length($a[1])))) . ",$b[0]";
+		} else {
+		  $ref = $a[0];
+		  $alt = ($b[0] eq $b[1])? $b[0] : "$b[0],$b[1]";
+		}
+	  } else {
+		$ref = $a[0]; $alt = $b[0];
+	  }
+	  print join("\t", @t[0,1], '.', $ref, $alt, $t[5], '.', '.'), "\n";
+	} else { # SNP
+	}
+  }
+}
+
 #
 # Usage
 #