comparison scripts/pickUniqIntervalPos.pl @ 0:28d1a6f8143f draft

planemo upload for repository https://github.com/portiahollyoak/Tools commit 132bb96bba8e7aed66a102ed93b7744f36d10d37-dirty
author portiahollyoak
date Mon, 25 Apr 2016 13:08:56 -0400
parents
children 9672fe07a232
comparison
equal deleted inserted replaced
-1:000000000000 0:28d1a6f8143f
1 #!/share/bin/perl
2 use Bio::Seq;
3 use List::Util qw(sum);
4
5 die "perl $0 <sam> <fragment_size>\n" if @ARGV<1;
6 open in,$ARGV[0];
7 my %pe;
8 while(<in>)
9 {
10 chomp;
11 my @f=split/\t/,$_,12;
12 ## read number 1 or 2
13 my ($rnum)=$f[1]=~/(\d)$/;
14
15 ## XT:A:*
16 my ($xt)=$f[11]=~/XT:A:(.)/;
17
18 my $strand="+";
19
20 ## parse CIGAR
21 if(($f[1]=~/R/)&&($f[8] > $ARGV[1])&&($f[8] <= 10000))
22 {
23 # CIGAR
24 my (@cigar_m)=$f[5]=~/(\d+)M/g;
25 my (@cigar_d)=$f[5]=~/(\d+)D/g;
26 my (@cigar_s)=$f[5]=~/(\d+)S/g;
27 my (@cigar_i)=$f[5]=~/(\d+)I/g;
28 my $aln_ln=sum(@cigar_m,@cigar_d);
29
30 # print $f[2],"\t",$f[3]-1+$aln_ln,"\t",$f[3]+$f[8],"\t$f[0]/$rnum\t","\n";
31 if ($f[2] =~ /^\d{1,2}$/) {$f[2]="chr$f[2]";}
32 print $f[2],"\t",$f[3]-6+$aln_ln,"\t",$f[7]+5,"\t$f[0]/$rnum\t","\n";
33 }
34 }
35 close in;
36