01.#!/user/bin/perl -w

02.

03.$ARGV[0]= “H:\Pro\blast\bin\test1.2”; # blast文件路径赋值给ARGV[0]

04.open(FIN,"$ARGV[0]")||die; #打开blast文件

05.open(FIN1,"H:\Pro\blast\bin\data\ testdatabase.fa ")||die; #打开基因组文件

06.open(FOUT,"> “H:\Pro\perl\bin\out]")||die; # 输出文件

07.

08.

09.while(){ #逐行读取文件

10.chomp; #并加换行符

11.

12.$l=$_; #将每一行的文件赋值给$|

13.if($l=~/^>(\S+)/){ #匹配以>的开头

14.$id=$1; #提取序列名称

15.$sq{$id}=""; #构建哈希表%sq，以序列名称$id为key

16.}

17.else{

18.$sq{$id}.=$l; #逐行将序列写入字符串

19.}

20.}

21.

22.while(){ #逐行读取文件

23.chomp; #去掉换行符

24.next if/^#/; #去掉注释行 跳过 以#开头的行

25.@a=split("\t"); #以空格为边界打散成数据组，分割文件

26.next if(@a < 8);#跳过元素小于8的行 此步可以省略吗？

27.$st=$a[8]; $len=$a[9]-$a[8];#计算序列长度 s.end –s.start？

28.if($a[8] eq "-"){ #表示$a[8]”-”符号

29.#$st=$a[7]; $len=-$len; #基因不存在？

30.}

31.$seq=substr($sq{$a[0]},$st-3000,$len+3000);#将基因组中$st即s.start为起始的$len个碱基写入子串seq 上游加3000bp

32.if($a[3] eq "-"){

33.$seq=reverse($seq); #倒置序列，得到互补序列

34.$seq=lc($seq);

35.$seq=~s/a/T/g;

36.$seq=~s/t/A/g;

37.$seq=~s/c/G/g;

38.$seq=~s/g/C/g;

39.}

40.#$a[-1]=~s/intron Parent=//; #序列名称输出 对应关系不确定？

41.if($a[8]=~/^ID=([^;]+);/){

42.$name=$1;

43.}

44.else{

45.print "$a[8]\n";

46.}

47.print FOUT ">$name\n$seq\n";#输出

#!/user/bin/perl -w

$ARGV[0]= "H:\\Pro\\blast\\bin\\test1.2"; # blast文件路径赋值给ARGV[0]
open(FIN,"$ARGV[0]")||die; #打开blast文件
open(FIN1,"H:\\Pro\\blast\\bin\\data\\testdatabase.fa")||die; #打开基因组文件
open(FOUT, "> H:\\Pro\\perl\\bin\\out")||die; # 输出文件


while(){ #逐行读取文件
        chomp; #并加换行符

        $l=$_; #将每一行的文件赋值给$|
        if($l=~/^>(\S+)/){ #匹配以>的开头
                $id=$1; #提取序列名称
                $sq{$id}=""; #构建哈希表%sq，以序列名称$id为key
        }
        else{
                $sq{$id}.=$l; #逐行将序列写入字符串
        }
}
while(){ #逐行读取文件
        chomp; #去掉换行符
next if/^#/;  #去掉注释行 跳过 以#开头的行
        @a=split("\t"); #以空格为边界打散成数据组，分割文件
        next if(@a < 8);#跳过元素小于8的行 此步可以省略吗？
        $st=$a[8]; $len=$a[9]-$a[8];#计算序列长度 s.end –s.start？
        if($a[8] eq "-"){ #表示$a[8]”-”符号
                #$st=$a[7]; $len=-$len; #基因不存在？
        }
        $seq=substr($sq{$a[0]},$st-3000,$len+3000);#将基因组中$st即s.start为起始的$len个碱基写入子串seq 上游加3000bp
        if($a[3] eq "-"){
                $seq=reverse($seq); #倒置序列，得到互补序列
                $seq=lc($seq);
                $seq=~s/a/T/g;
                $seq=~s/t/A/g;
                $seq=~s/c/G/g;
                $seq=~s/g/C/g;
        }
        #$a[-1]=~s/intron Parent=//; #序列名称输出 对应关系不确定？
        if($a[8]=~/^ID=([^;]+);/){
                $name=$1;
        }
        else{
                print "$a[8]\n";
        }
        print FOUT ">$name\n$seq\n";#输出

}