Mercurial > repos > portiahollyoak > temp
comparison scripts/pickOverlapPair.ex_MEM.pl @ 21:9672fe07a232 draft default tip
planemo upload for repository https://github.com/portiahollyoak/Tools commit 0fea84d05f8976b8360a8b4943ecb01b87e3ade0-dirty
author | mvdbeek |
---|---|
date | Mon, 05 Dec 2016 09:58:47 -0500 |
parents | ca36262102d8 |
children |
comparison
equal
deleted
inserted
replaced
20:6e02b9179a24 | 21:9672fe07a232 |
---|---|
54 my $rightlower=$right-500; | 54 my $rightlower=$right-500; |
55 my $rightupper=$right+500; | 55 my $rightupper=$right+500; |
56 my $chr_num=$a[0]; | 56 my $chr_num=$a[0]; |
57 $chr_num =~ s/chr//; | 57 $chr_num =~ s/chr//; |
58 if (($chrs{$a[0]} == 1) && (! defined $chrs{$chr_num})) {$chr_num=$a[0];} | 58 if (($chrs{$a[0]} == 1) && (! defined $chrs{$chr_num})) {$chr_num=$a[0];} |
59 system("samtools view -Xf 0x2 $title $chr_num\:$leftlower\-$leftupper $chr_num\:$rightlower\-$rightupper > temp.sam"); | 59 system("samtools view -f 0x2 $title $chr_num\:$leftlower\-$leftupper $chr_num\:$rightlower\-$rightupper > temp.sam"); |
60 | 60 |
61 open in,"temp.sam"; | 61 open in,"temp.sam"; |
62 my %ps=(); | 62 my %ps=(); |
63 my %me=(); | 63 my %me=(); |
64 my %uniqp=(); | 64 my %uniqp=(); |
68 while(<in>) | 68 while(<in>) |
69 { | 69 { |
70 chomp; | 70 chomp; |
71 my @f=split/\t/,$_,12; | 71 my @f=split/\t/,$_,12; |
72 ## read number 1 or 2 | 72 ## read number 1 or 2 |
73 my ($rnum)=$f[1]=~/(\d)$/; | 73 #my ($rnum)=$f[1]=~/(\d)$/; |
74 my $rnum=1; | |
75 if (($f[1] & 128) == 128) {$rnum=2;} | |
74 | 76 |
75 ## XT:A:* | 77 ## XT:A:* |
76 my $xt=""; | 78 my $xt=""; |
77 my @a=split(/\s+/, $_); | 79 my @a=split(/\s+/, $_); |
78 my $as=0; | 80 my $as=0; |
90 elsif ($as > 0) {$xt="U";} | 92 elsif ($as > 0) {$xt="U";} |
91 } | 93 } |
92 | 94 |
93 ## Coordinate | 95 ## Coordinate |
94 my $coor=$f[3]; | 96 my $coor=$f[3]; |
95 if ($f[1]=~/r/) | 97 if (($f[1] & 16) == 16) |
96 { | 98 { |
97 if ($xt eq "U") {$uniqm{$f[0]}=1;} | 99 if ($xt eq "U") {$uniqm{$f[0]}=1;} |
98 my (@cigar_m)=$f[5]=~/(\d+)M/g; | 100 my (@cigar_m)=$f[5]=~/(\d+)M/g; |
99 my (@cigar_d)=$f[5]=~/(\d+)D/g; | 101 my (@cigar_d)=$f[5]=~/(\d+)D/g; |
100 my (@cigar_s)=$f[5]=~/(\d+)S/g; | 102 my (@cigar_s)=$f[5]=~/(\d+)S/g; |
101 my (@cigar_i)=$f[5]=~/(\d+)I/g; | 103 my (@cigar_i)=$f[5]=~/(\d+)I/g; |
102 my $aln_ln=sum(@cigar_m,@cigar_d); | 104 my $aln_ln=sum(@cigar_m,@cigar_d); |
103 $me{$f[0]}=$f[3]+$aln_ln-1; | 105 $me{$f[0]}=$f[3]+$aln_ln-1; |
104 } | 106 } |
105 elsif ($f[1]=~/R/) { | 107 elsif (($f[1] & 32) == 32) { |
106 $ps{$f[0]}=$f[3]; | 108 $ps{$f[0]}=$f[3]; |
107 if ($xt eq "U") {$uniqp{$f[0]}=1;} | 109 if ($xt eq "U") {$uniqp{$f[0]}=1;} |
108 } | 110 } |
109 | 111 |
110 # ${$pe{$f[0]}}[$rnum-1]=[$xt,$coor]; | 112 # ${$pe{$f[0]}}[$rnum-1]=[$xt,$coor]; |