forked from vcflib/vcflib
-
Notifications
You must be signed in to change notification settings - Fork 3
/
Copy pathvcfindelproximity
executable file
·82 lines (67 loc) · 1.75 KB
/
vcfindelproximity
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
#!/usr/bin/env perl
# Show SNPs around an INDEL
# for line in the vcf
# stuff the line into a queue
# when you reach an indel
# record the position
# pop lines from the back of the queue until we are at the current position
#
my @lines;
my $prox = $ARGV[0];
my $lastchrom = "";
my $indelpos = 0;
while (<STDIN>) {
if ($_ =~ /^#/) {
print $_;
next;
}
$_ =~ /^(.+?)\t(.+?)\t(.+?)\t(.+?)\t(.+?)\t/;
my $chrom = $1;
my $pos = $2;
my $tag = $3;
my $ref = $4;
my $alt = $5;
#print "chrom: $chrom, pos: $pos, ref: $ref, alt: $alt\n";
# if new chrom, print out everything from last one
if ($lastchrom == "") {
$lastchrom = $chrom;
}
if ($chrom != $lastchrom) {
while ($lines) {
print pop(@lines);
}
}
unshift(@lines, $_);
my $diff = length($ref) - length($alt);
if ($diff != 0) {
# insertion
if ($indelpos == 0) {
$indelpos = $pos;
}
$nextindelpos = $pos;
#print "last $indelpos next $nextindelpos\n";
while (@lines) {
my $line = pop(@lines);
$line =~ /^(.+?)\t(.+?)\t(.+?)\t(.+?)\t(.+?)\t/;
my $c = $1;
my $p = $2;
my $t = $3;
my $r = $4;
my $a = $5;
# print indels
if (length($r) - length($a) != 0) {
print $line;
} else {
# print other events which are more than prox away from indels
if (abs($indelpos - $p) >= $prox and abs($nextindelpos - $p) >= $prox) {
print $line;
}
}
}
$indelpos = $pos;
}
}
# flush lines end of file
while ($lines) {
print pop(@lines);
}