forked from khmiga/nanopore_kmer
-
Notifications
You must be signed in to change notification settings - Fork 0
/
Copy pathkmer.pl
57 lines (53 loc) · 1.01 KB
/
kmer.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
#! /usr/bin/perl -w
use strict;
my $t=0;
my $seq;
my $c;
my $h;
die "must provide input file\nusage: kmer.pl fastaFile k-merSize\n" unless ($ARGV[0]);
die "must provide kmer size\nusage: kmer.pl fastaFile k-merSize\n" unless($ARGV[1]);
my $outFile=$ARGV[0] . '.' . $ARGV[1] . 'mer.fa';
open(OUT,">$outFile");
open(FA,"$ARGV[0]");
while(<FA>){
chomp;
if($_=~ />(\S+)/){
my $tmpH=$1;
if($t==0){
$t=1;
$h=$tmpH;
}
else{
kmer($seq,$h);
$seq='';
$h=$tmpH;
}
}
else{
$seq.=$_;
}
}
kmer($seq,$h);
sub kmer{
my ($s,$header)=@_;
my $len=length($s);
my $start = 0;
my $slide = 1;
my $windowSize=12;
my $end = $start + $windowSize;
my $o=0;
while ($end <= $len) {
my $oligo = substr($s,$start,$windowSize);
$start += $slide;
$end = $start + $slide;
if($oligo!~ /N/){
my $oL=length($oligo);
if($oL ==$windowSize){
my $rc=reverse($oligo);
$rc=~ tr/ATCGN/TAGCN/;
print OUT "$header\t$start\t$oligo\t$rc\n";
}
}
}
}
close OUT;