[Bio] / FigKernelPackages / FIG.pm Repository:
ViewVC logotype

Diff of /FigKernelPackages/FIG.pm

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

revision 1.620, Mon Sep 17 19:52:10 2007 UTC revision 1.621, Wed Sep 19 21:08:12 2007 UTC
# Line 931  Line 931 
931          }          }
932      }      }
933    
934        #
935        # If this is a RAST replacement genome, perform subsystem salvage.
936        #
937        my $replaces = &FIG::file_head("$genome_dir/REPLACES", 1);
938        chomp $replaces;
939        if (-f "$genome_dir/RAST" and $replaces ne '')
940        {
941            if (open(MAP, "$genome_dir/peg_maps"))
942            {
943                my %map;
944                while (<MAP>)
945                {
946                    chomp;
947                    my($f, $t) = split(/\t/);
948                    $map{$f} = $t;
949                }
950                close(MAP);
951    
952                $self->perform_subsystem_salvage([[$replaces, $genome]], \%map);
953            }
954            else
955            {
956                warn "Genome $genome that replaces $replaces is missing a peg_maps file: $!";
957            }
958        }
959    
960      return $rc;      return $rc;
961  }  }
962    
# Line 15342  Line 15368 
15368   );   );
15369  }  }
15370    
15371    =head3 perform_subsystem_salvage
15372    
15373    C<< my $glist = [['273035.1', '273035.4']]; >>
15374    C<< my $pmap = { 'fig|273035.1.peg.1' => 'fig|273035.4.peg.4', ... }; >>
15375    C<< $fig->perform_subsystem_salvage($glist, $pmap); >>
15376    
15377    For each subsystem in this SEED, perform a subsystem salvage operation for each old-genome / new-genome pair in $glist.
15378    This operation will determine if the old genome exists in the subsystem. If it does, the new genome is
15379    added to the subsystem, and we attempt to map the pegs from the cells in the old subsystem's row to the new
15380    subsystem. If all pegs map, we copy the variant code for the genome. If all cells did not map, we prepend
15381    a * to the variant code before copying.
15382    
15383    =cut
15384    
15385    sub perform_subsystem_salvage
15386    {
15387        my($fig, $genome_pairs, $map) = @_;
15388    
15389        for my $ssname ($fig->all_subsystems())
15390        {
15391            my $ss = new Subsystem($ssname, $fig);
15392    
15393            for my $gpair (@$genome_pairs)
15394            {
15395                my($g, $ng) = @$gpair;
15396    
15397                my $idx = $ss->get_genome_index($g);
15398                if (defined($idx))
15399                {
15400                    print "Salvaging $ssname for old_genome=$g new_genome=$ng\n";
15401                    my $row = $ss->get_row($idx);
15402    
15403                    my $new_idx = $ss->get_genome_index($ng);
15404    
15405                    if (defined($new_idx))
15406                    {
15407                        warn "$ng was already present as $new_idx\n";
15408                        next;
15409                    }
15410    
15411                    $new_idx = $ss->add_genome($ng);
15412    
15413                    if (!defined($new_idx))
15414                    {
15415                        die "Subsystem $ss add_genome($ng) failed";
15416                    }
15417    
15418                    my @nrow;
15419                    my $pegcount = 0;
15420                    my $mapcount = 0;
15421                    for (my $ci = 0; $ci < @$row; $ci++)
15422                    {
15423                        my $c = $row->[$ci];
15424    
15425                        my @nc;
15426    
15427                        for my $p (@$c)
15428                        {
15429                            $pegcount++;
15430                            my $new = $map->{$p};
15431                            if ($new)
15432                            {
15433                                $mapcount++;
15434                                push(@nc, $new);
15435                            }
15436                        }
15437    
15438                        $ss->set_pegs_in_cell($new_idx, $ci, \@nc);
15439                    }
15440    
15441                    my $vc = $ss->get_variant_code($idx);
15442                    if ($pegcount < $mapcount)
15443                    {
15444                        $vc = "*$vc";
15445                    }
15446                    $ss->set_variant_code($new_idx, $vc);
15447                    $ss->write_subsystem(1);
15448                }
15449            }
15450        }
15451    }
15452    
15453    
15454    
15455  =head3 all_constructs  =head3 all_constructs
15456    

Legend:
Removed from v.1.620  
changed lines
  Added in v.1.621

MCS Webmaster
ViewVC Help
Powered by ViewVC 1.0.3