perl/etc/WoodyMiaoLin/PbePopGen/addRG.pl

   1 #!/usr/bin/perl
   2 use strict;
   3 use warnings;
   4
   5 die "This program add HEADER and TAG of READ GROUP to bam file.\nAuther: Woody\nUsage: $0 <in.bam> <out.bam>\n" if @ARGV<2;
   6
   7 my $in = shift;
   8 my $out = shift;
   9 $in =~ /^.+\/(\w+)\.rmdup\.bam$/;
  10 my $ip = $1;
  11 $out =~ /^(.+)\.bam$/;
  12 my $op = $1;
  13
  14 open I, "-|", "samtools-0.1.7 view -h $in";
  15 open HEAD, ">", "$op.header.sam";
  16 open READ, "|-", "samtools-0.1.7 view -bS - >$op.read.bam";
  17
  18 warn "$ip adding TAG...\n";
  19 my %RG;
  20 while (<I>) {
  21         if (/^@/) {
  22                 print HEAD;
  23                 print READ;
  24         } else {
  25                 chomp;
  26                 my @a = split /\t/;
  27                 my @b = split /:/, $a[0];
  28                 my $rg = "$b[0]_$b[1]_$b[2]_$b[3]";
  29                 $RG{$rg} = 1;
  30                 print READ "$_\tRG:Z:$rg\n";
  31         }
  32 }
  33
  34 close I;
  35 close READ;
  36
  37 foreach (sort keys %RG) {
  38         print HEAD "\@RG\tID:$_\tPL:ILLUMINA\tSM:$ip\n";
  39 }
  40 close HEAD;
  41
  42 warn "$ip reheadering...\n";
  43 system "samtools reheader -P $op.header.sam $op.read.bam >$op.bam";
  44 system "rm $op.header.sam $op.read.bam";
  45 warn "Done!\n";