]> git.donarmstrong.com Git - rsem.git/blob - rsem_perl_utils.pm
Modified rsem-tbam2gbam so that the original alignment quality MAPQ will be preserved...
[rsem.git] / rsem_perl_utils.pm
1 #!/usr/bin/perl
2
3 package rsem_perl_utils;
4
5 use strict;
6
7 require Exporter;
8 our @ISA = qw(Exporter);
9 our @EXPORT = qw(runCommand);
10 our @EXPORT_OK = qw(runCommand collectResults showVersionInfo);
11
12 # command, {err_msg}
13 sub runCommand {
14     print $_[0]."\n";
15     my $status = system($_[0]);
16
17     if ($? == -1) {
18         my @arr = split(/[ \t]+/, $_[0]);
19         print "$arr[0] : $!!\n";
20         print "Please check if you have compiled the associated codes by typing related \"make\" commands and/or made related executables ready to use.\n";
21         exit(-1);
22     }
23
24     if ($status != 0) {
25         my $errmsg = "";
26         if (scalar(@_) > 1) { $errmsg .= $_[1]."\n"; }
27         $errmsg .= "\"$_[0]\" failed! Plase check if you provide correct parameters/options for the pipeline!\n";
28         print $errmsg;
29         exit(-1);
30     }
31     print "\n";
32 }
33
34
35 my @transcript_title = ("transcript_id", "gene_id", "length", "effective_length", "expected_count", "TPM", "FPKM", "IsoPct", "pme_expected_count", "pme_TPM", "pme_FPKM", "IsoPct_from_pme_TPM", "TPM_ci_lower_bound", "TPM_ci_upper_bound", "FPKM_ci_lower_bound", "FPKM_ci_upper_bound");
36
37 my @gene_title = ("gene_id", "transcript_id(s)", "length", "effective_length", "expected_count", "TPM", "FPKM", "pme_expected_count", "pme_TPM", "pme_FPKM", "TPM_ci_lower_bound", "TPM_ci_upper_bound", "FPKM_ci_lower_bound", "FPKM_ci_upper_bound");
38
39 # inpF, outF
40 sub collectResults {
41     my $local_status;
42     my ($inpF, $outF);
43     my @results = ();
44     my $line;
45
46     $inpF = $_[1];
47     $outF = $_[2];
48
49     $local_status = open(INPUT, $inpF);
50     if ($local_status == 0) { print "Fail to open file $inpF!\n"; exit(-1); }
51     
52     @results = ();
53     
54     while ($line = <INPUT>) {
55         chomp($line);
56         my @local_arr = split(/\t/, $line);
57         push(@results, \@local_arr); 
58     }
59
60     close(INPUT);
61
62     $local_status = open(OUTPUT, ">$outF");
63     if ($local_status == 0) { print "Fail to create file $outF!\n"; exit(-1); }
64
65     my $n = scalar(@results);
66     my $m = scalar(@{$results[0]});
67
68     $" = "\t";
69
70     my @out_arr = ();
71     for (my $i = 0; $i < $n; $i++) {
72         if ($_[0] eq "isoform") { push(@out_arr, $transcript_title[$i]); }
73         elsif ($_[0] eq "gene") { push(@out_arr, $gene_title[$i]); }
74         else { print "A bug on 'collectResults' is detected!\n"; exit(-1); }
75     }
76     print OUTPUT "@out_arr\n";
77
78     for (my $i = 0; $i < $m; $i++) {
79         @out_arr = ();
80         for (my $j = 0; $j < $n; $j++) { push(@out_arr, $results[$j][$i]); }
81         print OUTPUT "@out_arr\n"; 
82     }
83
84     close(OUTPUT);
85 }
86
87 # dir
88 sub showVersionInfo {
89     open(INPUT, "$_[0]\WHAT_IS_NEW");
90     my $line = <INPUT>;
91     chomp($line);
92     close(INPUT);
93     print "Current version is $line\n";
94     exit(0);
95 }
96
97 1;