利用perl脚本统计基因组中的GAP

perl find_gap.pl genome.fa > gap.info

### find_gap.pl
local $/ = ">";
<>;
while (<>) {
    s/>//;
    my @a = split/\n/,$_,2;
    $a[1] =~ s/\n//g;
    my @seq = split//,$a[1];
    foreach  (0..$#seq) {
        if ($seq[$_] =~ /N/) {
            push @{$hash{$a[0]}},$_;
        }
    }
}

foreach my $chr (sort{$a cmp $b}keys %hash) {
    my @gap = group(@{$hash{$chr}});
    foreach  (@gap) {
        print "$chr\t$_->[0]\t$_->[1]\n";
    }
}

sub group {
    my @tmp = ();
    my $tmp = $_[0];
    my $cnt = 0;
    my @return;
    foreach  (0..$#_-1) {
        if ($_[$_+1] - $_[$_] != 1) {
            push @{$tmp[$cnt]},$_[$_];
            $cnt++;
        }else{
            push @{$tmp[$cnt]},$_[$_];
        }
    }
    push @{$tmp[$cnt]},$_[-1];
    foreach  (@tmp) {
        push @return,[(sort{$a<=>$b}@{$_})[0],scalar(@{$_})];
    }
    return @return;
}