]> git.donarmstrong.com Git - samtools.git/commitdiff
* move indel_filter.pl to samtools.pl
authorHeng Li <lh3@live.co.uk>
Tue, 17 Feb 2009 10:49:55 +0000 (10:49 +0000)
committerHeng Li <lh3@live.co.uk>
Tue, 17 Feb 2009 10:49:55 +0000 (10:49 +0000)
misc/indel_filter.pl [deleted file]
misc/samtools.pl [new file with mode: 0755]

diff --git a/misc/indel_filter.pl b/misc/indel_filter.pl
deleted file mode 100755 (executable)
index d357fcc..0000000
+++ /dev/null
@@ -1,42 +0,0 @@
-#!/usr/bin/perl -w
-
-# Contact: lh3
-
-use strict;
-use warnings;
-use Getopt::Std;
-
-my %opts = (D=>100, m=>10, r=>undef);
-getopts('D:m:r', \%opts);
-
-die(qq/
-Usage:   indel_filter.pl [options] <in.indel>\n
-Options: -D INT    maximum read depth [$opts{D}]
-         -m INT    minimum distance between two adjacent indels [$opts{m}]
-\n/) if (@ARGV == 0 && -t STDIN);
-
-my (@arr1, @arr2);
-my ($curr, $last) = (\@arr1, \@arr2);
-my $is_ref = defined($opts{r})? 1 : 0;
-while (<>) {
-  my @t = split;
-  next if ($t[2] ne '*');
-  if (!$is_ref) {
-       next if ($t[3] eq '*/*');
-       next if ($t[5] == 0);
-  }
-  next if ($t[7] > $opts{D});
-  @$curr = ($t[0], $t[1], $t[5], $_);
-  my $do_swap = 1;
-  if (defined $last->[0]) {
-       if ($curr->[0] eq $last->[0] && $last->[1] + $opts{m} > $curr->[1]) {
-         $do_swap = 0 if ($last->[2] > $curr->[2]);
-       } else { # then print
-         print $last->[3];
-       }
-  }
-  if ($do_swap) {
-       my $tmp = $curr; $curr = $last; $last = $tmp;
-  }
-}
-print $last->[3] if (defined $last->[0]);
diff --git a/misc/samtools.pl b/misc/samtools.pl
new file mode 100755 (executable)
index 0000000..0b77782
--- /dev/null
@@ -0,0 +1,81 @@
+#!/usr/bin/perl -w
+
+# Author: lh3
+
+use strict;
+use warnings;
+use Getopt::Std;
+
+my $version = '0.1.0';
+&usage if (@ARGV < 1);
+
+my $command = shift(@ARGV);
+my %func = (indelFilter=>\&indelFilter, showALEN=>\&showALEN);
+
+die("Unknown command \"$command\".\n") if (!defined($func{$command}));
+&{$func{$command}};
+exit(0);
+
+#
+# showALEN
+#
+
+sub showALEN {
+  die(qq/Usage: samtools.pl showALEN <in.sam>\n/) if (@ARGV == 0 && -t STDIN);
+  while (<>) {
+       my @t = split;
+       my $l = 0;
+       $_ = $t[5];
+       s/(\d+)[SMI]/$l+=$1/eg;
+       print join("\t", @t[0..5]), "\t$l\t", join("\t", @t[6..$#t]), "\n";
+  }
+}
+
+#
+# indelFilter
+#
+
+sub indelFilter {
+  my %opts = (D=>100, m=>10, r=>undef);
+  getopts('D:m:r', \%opts);
+
+  die(qq/
+Usage:   samtools.pl indelFilter [options] <in.indel>\n
+Options: -D INT    maximum read depth [$opts{D}]
+         -m INT    minimum distance between two adjacent indels [$opts{m}]
+\n/) if (@ARGV == 0 && -t STDIN);
+
+  my (@arr1, @arr2);
+  my ($curr, $last) = (\@arr1, \@arr2);
+  my $is_ref = defined($opts{r})? 1 : 0;
+  while (<>) {
+       my @t = split;
+       next if ($t[2] ne '*');
+       if (!$is_ref) {
+         next if ($t[3] eq '*/*');
+         next if ($t[5] == 0);
+       }
+       next if ($t[7] > $opts{D});
+       @$curr = ($t[0], $t[1], $t[5], $_);
+       my $do_swap = 1;
+       if (defined $last->[0]) {
+         if ($curr->[0] eq $last->[0] && $last->[1] + $opts{m} > $curr->[1]) {
+               $do_swap = 0 if ($last->[2] > $curr->[2]);
+         } else { # then print
+               print $last->[3];
+         }
+       }
+       if ($do_swap) {
+         my $tmp = $curr; $curr = $last; $last = $tmp;
+       }
+  }
+  print $last->[3] if (defined $last->[0]);
+}
+
+sub usage {
+  die(qq/
+Usage:   samtools.pl <command> [<arguments>]\n
+Command: indelFilter   filter indels generated by `pileup -c'
+         showALEN      print alignment length (ALEN) following CIGAR
+\n/);
+}