ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/gclib/scripts/gffsummary.pl
Revision: 24
Committed: Tue Jul 26 21:46:39 2011 UTC (8 years, 1 month ago) by gpertea
File size: 1907 byte(s)
Log Message:
Line File contents
1 #!/usr/bin/perl
2 use strict;
3 use Getopt::Std;
4
5 my $usage = q/Usage:
6 gffstat [-a] [-t <ftype>] [-f <subfeature_name>] input.gffs..
7
8 Returns statistics about the number, length and composition of
9 'mRNA' features and their 'CDS' segments.
10 Options:
11 -q omit header line
12 -a treat all input files as a single input stream,
13 so only the overall counts (totals) are shown
14 -t only consider entries with <ftype> in the
15 second column (e.g. "jigsaw")
16 -f use <subfeature_name> instead of default 'CDS'
17 (needed for exon statistics)
18 /;
19 getopts('qat:f:') || die($usage."\n");
20
21 my ($readall, $tflt)=($Getopt::Std::opt_a, $Getopt::Std::opt_t);
22 my $subfname= $Getopt::Std::opt_f || 'CDS';
23 $subfname='CDS' unless $subfname;
24 print join("\t",' ','G_count', 'G_cov', 'Ex_count', 'Ex_cov')."\n"
25 unless ($Getopt::Std::opt_q);
26 if ($readall || @ARGV==0) {
27 processInput();
28 }
29 else {
30 foreach my $f (@ARGV) {
31 processInput($f);
32 }
33 }
34
35 sub processInput {
36 my $inf=$_[0];
37 my ($gcount, $gextent, $exoncount, $exoncov);
38 my $rlinefunc;
39 if ($inf) {
40 open(INFILE, $inf)
41 || die("Error opening file $inf!\n");
42 $rlinefunc = \&readfile;
43 } else {
44 $rlinefunc=\&readstdin;
45 }
46
47 while (&$rlinefunc()) {
48 next if m/^\s*#/;
49 my ($gobj, $ftype, $fname, $fstart, $fend, $fscore, $fstrand, $frame, $fdesrc)=split(/\t/);
50 next if ($tflt && $tflt ne $ftype);
51 ($fstart, $fend)=($fend, $fstart) if $fstart>$fend;
52 if ($fname eq 'mRNA') {
53 $gcount++;
54 $gextent+= ($fend-$fstart+1);
55 next;
56 }
57 if ($fname eq $subfname) {
58 $exoncount++;
59 $exoncov+=($fend-$fstart+1);
60 }
61 }#while parsing
62 my $label='Total';
63 if ($inf) {
64 close(INFILE);
65 ($label)=($inf=~m/([^\/]+)$/);
66 }
67 print join("\t", $label.':', $gcount, $gextent, $exoncount, $exoncov)."\n";
68 }
69
70 sub readfile {
71 $_=<INFILE>;
72 return $_;
73 }
74
75 sub readstdin {
76 $_=<>;
77 return $_;
78 }

Properties

Name Value
svn:executable *