Mercurial > repos > big-tiandm > mirplant2
comparison matching.pl @ 11:861eb4fbfbd0 draft
Uploaded
author | big-tiandm |
---|---|
date | Fri, 25 Jul 2014 05:20:44 -0400 |
parents | |
children |
comparison
equal
deleted
inserted
replaced
10:5ac42e277c76 | 11:861eb4fbfbd0 |
---|---|
1 #!/usr/bin/perl -w | |
2 #Filename: | |
3 #Author: Tian Dongmei | |
4 #Email: tiandm@big.ac.cn | |
5 #Date: 2013/7/19 | |
6 #Modified: | |
7 #Description: | |
8 my $version=1.00; | |
9 | |
10 use strict; | |
11 use Getopt::Long; | |
12 | |
13 my %opts; | |
14 GetOptions(\%opts,"i=s","g=s","index:s","v:i","p:i","r:s","o=s","time:s","h"); | |
15 if (!(defined $opts{i} and defined $opts{o} ) || defined $opts{h}) { #necessary arguments | |
16 &usage; | |
17 } | |
18 | |
19 my $filein=$opts{'i'}; | |
20 my $fileout=$opts{'o'}; | |
21 unless ($fileout=~/\/$/) {$fileout.="/";} | |
22 my $genome=$opts{'g'}; | |
23 my $mis=defined $opts{'v'}? $opts{'v'} : 0; | |
24 my $hits=defined $opts{'r'}? $opts{'r'} : 25; | |
25 my $index=defined $opts{'index'} ? $opts{'index'} : ""; | |
26 my $threads=defined $opts{'p'} ? $opts{'p'} : 1; | |
27 | |
28 | |
29 my $time=&Time(); | |
30 if (defined $opts{'time'}) { | |
31 $time=$opts{'time'}; | |
32 } | |
33 | |
34 my $mapdir=$fileout."/genome_match_".$time; | |
35 if(not -d $mapdir){ | |
36 mkdir $mapdir; | |
37 } | |
38 chdir $mapdir; | |
39 ###check genome index | |
40 if (-s $index.".1.ebwt") { | |
41 }else{ | |
42 `bowtie-build $genome genome`; | |
43 $index="genome"; | |
44 } | |
45 | |
46 ### genome mapping | |
47 `bowtie -v $mis -f -p $threads -m $hits -a --best --strata $index $filein --al genome_mapped.fa --un genome_not_mapped.fa --max genome_mapped_Mlimit.fa > genome_mapped.bwt`; | |
48 | |
49 #`convert_bowtie_to_blast.pl genome_mapped.bwt genome_mapped.fa $genome > genome_mapped.bst`; | |
50 | |
51 sub Time{ | |
52 my $time=time(); | |
53 my ($sec,$min,$hour,$day,$month,$year) = (localtime($time))[0,1,2,3,4,5,6]; | |
54 $month++; | |
55 $year+=1900; | |
56 if (length($sec) == 1) {$sec = "0"."$sec";} | |
57 if (length($min) == 1) {$min = "0"."$min";} | |
58 if (length($hour) == 1) {$hour = "0"."$hour";} | |
59 if (length($day) == 1) {$day = "0"."$day";} | |
60 if (length($month) == 1) {$month = "0"."$month";} | |
61 #print "$year-$month-$day $hour:$min:$sec\n"; | |
62 return("$year-$month-$day-$hour-$min-$sec"); | |
63 } | |
64 | |
65 sub usage{ | |
66 print <<"USAGE"; | |
67 Version $version | |
68 Usage: | |
69 $0 -i -o | |
70 options: | |
71 -i input file# input reads fasta/fastq file | |
72 -g input file# genome file | |
73 -index file-prefix #(must be indexed by bowtie-build) The parameter | |
74 string must be the prefix of the bowtie index. For instance, if | |
75 the first indexed file is called 'h_sapiens_37_asm.1.ebwt' then | |
76 the prefix is 'h_sapiens_37_asm'.##can be null | |
77 -v <int> report end-to-end hits w/ <=v mismatches; ignore qualities,default 0; | |
78 | |
79 -p/--threads <int> number of alignment threads to launch (default: 1) | |
80 | |
81 -r int a read is allowed to map up to this number of positions in the genome | |
82 default is 25 | |
83 | |
84 -o output directory | |
85 | |
86 -h help | |
87 USAGE | |
88 exit(1); | |
89 } | |
90 |