forked from sanahmed/PhaME
-
Notifications
You must be signed in to change notification settings - Fork 0
/
checkNUCmer.pl
61 lines (54 loc) · 1001 Bytes
/
checkNUCmer.pl
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
#!/usr/bin/perl -w
use strict;
use Getopt::Long;
use FileHandle;
my $infile;
my $reffile;
my $reference;
my $length=0;
GetOptions(
'i=s' => \$infile,
'r=s' => \$reffile,
'h|help' => sub{usage()},
);
read_reference($reffile);
read_gaps($infile);
sub read_reference
{
my ($header,@seq,$sequence);
my $fh=FileHandle->new($reffile)|| die "$!";
if ($fh->open("< $reffile")){
$/=">";
while(<$fh>){
$_=~ s/\>//g;
unless ($_){next;};
($header,@seq)=split /\n/,$_;
$sequence=join "",@seq;
$reference=length $sequence;
}
}
}
sub read_gaps
{
my $line;
my $fh=FileHandle->new($infile)|| die "$!";
if ($fh->open("< $infile")){
$/="\n";
while (<$fh>){
unless ($_){next;};
my ($ref,$start,$end,$len,$query)=split /\t/,$_;
$length+=$len;
}
}
my $difference =$length/$reference;
if ($difference <0.75){print "0";}
if ($difference >=0.75){print "1";}
}
sub usage
{
print STDERR"
Usage:
$0 -i <gaps file> -r <reference file>
";
exit;
}