[Bio] / FigKernelScripts / svr_pegs_in_subsystems.pl Repository:
ViewVC logotype

Diff of /FigKernelScripts/svr_pegs_in_subsystems.pl

Parent Directory Parent Directory | Revision Log Revision Log | View Patch Patch

revision 1.6, Thu Jul 29 14:48:54 2010 UTC revision 1.7, Mon Apr 4 19:25:34 2011 UTC
# Line 1  Line 1 
1  #!/usr/bin/perl -w  #!/usr/bin/perl -w
2    
3    use strict;
4  use SAPserver;  use SAPserver;
5  use Getopt::Long;  use Getopt::Long;
6    use ScriptThing;
7    
8  # This is a SAS Component  # This is a SAS Component
9    
# Line 41  Line 43 
43    
44  The URL for the Sapling server, if it is to be different from the default.  The URL for the Sapling server, if it is to be different from the default.
45    
46    =item c
47    
48    Number (1-based) of the column in the input file containing the subsystem name. If omitted,
49    the last column is used.
50    
51  =back  =back
52    
53  =cut  =cut
# Line 48  Line 55 
55  my $noroles = 0;  my $noroles = 0;
56  my $group = 0;  my $group = 0;
57  my $show_owner = 0;  my $show_owner = 0;
58    my $column;
59  my $oldid = "";  my $oldid = "";
60  my $url = "";  my $url = "";
61    my $inputFile = "-";
62    
63  $0 =~ m/([^\/]+)$/;  $0 =~ m/([^\/]+)$/;
64  my $self = $1;  my $self = $1;
65  my $usage = "$self [--noroles --group --url=http://...] GenomeF < SubsystemIDs";  my $usage = "$self [--noroles --group --url=http://...] GenomeF < SubsystemIDs";
66    
67  my $rc = GetOptions("noroles" => \$noroles, "group" => \$group, "url=s" => \$url);  my $rc = GetOptions("noroles" => \$noroles, "group" => \$group, "url=s" => \$url, "i=s" => \$inputFile, "c=i" => \$column);
68    
69  if (!$rc) {  if (!$rc) {
70      die "\n   usage: $usage\n\n";      die "\n   usage: $usage\n\n";
# Line 64  Line 73 
73  my $roles = $noroles ? 0 : 1;  my $roles = $noroles ? 0 : 1;
74  my $ss = SAPserver->new(url => $url);  my $ss = SAPserver->new(url => $url);
75    
76  open GENOMES, "<" . $ARGV[$#ARGV] || die "Genome file error: $!";  open my $gh, "<" . $ARGV[$#ARGV] || die "Genome file error: $!";
   
 my @genomes;  
 my @subs;  
 while (<GENOMES>) {  
     chomp;  
     push (@genomes, $_);  
 }  
 while (<STDIN>) {  
     chomp;  
     push (@subs, $_);  
 }  
77    
78  my $pegs_inss = $ss->pegs_in_subsystems(\@genomes, \@subs);  my @genomes = ScriptThing::GetList($gh);
79    close $gh;
80    open(my $ih, "<$inputFile") || die "Error opening input: $!";
81    while (my @tuples = ScriptThing::GetBatch($ih, 10, $column)) {
82        my @subs = map { $_->[0] } @tuples;
83        my $pegs_inss = $ss->pegs_in_subsystems(-genomes => \@genomes,
84                                                -subsystems => \@subs);
85        # Loop through the incoming lines, and pair the results with the inputs.
86        for my $tuple (@tuples) {
87            # Get the current line and its subsystem ID.
88            my ($sub, $line) = @$tuple;
89            # Get the role hash for this subsystem.
90            my $roleHash = $pegs_inss->{$sub};
91            # Only proceed if we found something.
92            if ($roleHash) {
93                # Are we including roles in the output?
94  if ($roles) {  if ($roles) {
95      foreach my $ss_role (@{$pegs_inss}) { #foreach subsystem/role                  # Yes. Loop through the roles.
96          #(ss, role, (peg))                  for my $role (sort keys %$roleHash) {
97                        # Get the features for this role.
98                        my $fids = $roleHash->{$role};
99                        # Are we in group mode?
100          if ($group) {          if ($group) {
101              print $ss_role->[0], "\t", $ss_role->[1]->[0], "\t", join (",", @{$ss_role->[1]->[1]}), "\n";                          # Yes. Put all the pegs on a single line.
102                            print "$line\t$role\t" . join(", ", @$fids) . "\n";
103          } else {          } else {
104              foreach my $peg (@{$ss_role->[1]->[1]}) { #foreach peg in this peg list                          # No. Put each peg on a line by itself.
105                  print join("\t", ($ss_role->[0], $ss_role->[1]->[0], $peg)), "\n";                          for my $fid (@$fids) {
106                                print "$line\t$role\t$fid\n";
107              }              }
108          }          }
109      }      }
110                } else {
111  } else { # no roles                  # Roles are not included in the output. We need to create
112      my %ss_pegs;                  # a list of the features for all the roles.
113      foreach my $ss_role (@{$pegs_inss}) {                  my %fids;
114          foreach my $peg (@{$ss_role->[1]->[1]}) {                  for my $role (keys %$roleHash) {
115              $ss_pegs{$ss_role->[0]}{$peg} = 1;                      my $fidList = $roleHash->{$role};
116                        for my $fid (@$fidList) {
117                            $fids{$fid} = 1;
118          }          }
119      }      }
120      for my $ss (keys %ss_pegs) {                  my @fids = sort keys %fids;
121                    # Are we grouping the features?
122          if ($group) {          if ($group) {
123              print $ss, "\t", join(",", keys %{$ss_pegs{$ss}}), "\n";                      # Yes. Put all the pegs on one line.
124                        print "$line\t" . join(", ", @fids) . "\n";
125          } else {          } else {
126              foreach my $peg (keys %{$ss_pegs{$ss}}) {                      # No. Put each peg on its own line.
127                  print join("\t", $ss, $peg), "\n";                      for my $fid (@fids) {
128                            print "$line\t$fid\n";
129                        }
130                    }
131              }              }
132          }          }
133      }      }

Legend:
Removed from v.1.6  
changed lines
  Added in v.1.7

MCS Webmaster
ViewVC Help
Powered by ViewVC 1.0.3