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

View of /FigKernelScripts/find_special_proteins.pl

Parent Directory Parent Directory | Revision Log Revision Log


Revision 1.4 - (download) (as text) (annotate)
Fri Jul 18 16:53:59 2008 UTC (11 years, 4 months ago) by gdpusch
Branch: MAIN
CVS Tags: mgrast_dev_08112011, rast_rel_2009_05_18, mgrast_dev_08022011, rast_rel_2014_0912, myrast_rel40, mgrast_dev_05262011, rast_rel_2008_12_18, mgrast_dev_04082011, rast_rel_2010_0928, rast_2008_0924, mgrast_version_3_2, mgrast_dev_12152011, mgrast_dev_06072011, rast_rel_2008_09_30, rast_rel_2009_0925, rast_rel_2010_0526, rast_rel_2014_0729, mgrast_dev_02212011, rast_rel_2010_1206, mgrast_release_3_0, mgrast_dev_03252011, rast_rel_2010_0118, mgrast_rel_2008_0924, mgrast_rel_2008_1110_v2, rast_rel_2009_02_05, rast_rel_2011_0119, mgrast_rel_2008_0625, mgrast_release_3_0_4, mgrast_release_3_0_2, mgrast_release_3_0_3, mgrast_release_3_0_1, mgrast_dev_03312011, mgrast_release_3_1_2, mgrast_release_3_1_1, mgrast_release_3_1_0, mgrast_dev_04132011, rast_rel_2008_10_09, mgrast_dev_04012011, rast_release_2008_09_29, mgrast_rel_2008_0806, mgrast_rel_2008_0923, mgrast_rel_2008_0919, rast_rel_2009_07_09, rast_rel_2010_0827, mgrast_rel_2008_1110, myrast_33, rast_rel_2011_0928, rast_rel_2008_09_29, mgrast_rel_2008_0917, rast_rel_2008_10_29, mgrast_dev_04052011, mgrast_dev_02222011, rast_rel_2009_03_26, mgrast_dev_10262011, rast_rel_2008_11_24, rast_rel_2008_08_07, HEAD
Changes since 1.3: +3 -3 lines
Suppress writing "called_by" and "special_pegs" when "-keep" is selected. -- /gdp

# -*- perl -*-
########################################################################
# Copyright (c) 2003-2006 University of Chicago and Fellowship
# for Interpretations of Genomes. All Rights Reserved.
#
# This file is part of the SEED Toolkit.
# 
# The SEED Toolkit is free software. You can redistribute
# it and/or modify it under the terms of the SEED Toolkit
# Public License. 
#
# You should have received a copy of the SEED Toolkit Public License
# along with this program; if not write to the University of Chicago
# at info@ci.uchicago.edu or the Fellowship for Interpretation of
# Genomes at veronika@thefig.info or download a copy from
# http://www.theseed.org/LICENSE.TXT.
########################################################################

use strict;
use warnings;

$0 =~ m/([^\/]+)$/o;
my $self  = $1;
my $usage = "usage: find_special_proteins [-h(elp)] To_Call_Dir [old]";

use FIG;
# use FIGV;
my $fig = new FIG;

use NewGenome;
use ToCall;

if (@ARGV && ($ARGV[0] =~ m/-h(elp)?/)) {
    print STDERR "\n   $usage\n\n";
    exit (0);
}

my $to_call_dir;
($to_call_dir = shift @ARGV) || die "\n   $usage\n\n";

### This is a bit complex.  Normally, one uses this program to call genes
# in a newly-sequenced genome.  In this case the functionality of NewGenome.pm
# is needed.  To call genes in an existing genome (to, for example, see how well
# we are doing) you would use ToCall.pm.

my ($to_call, $keep);
if ((@ARGV > 0) && ($ARGV[0] =~ /^old/i)) {
    $keep    = 1;
    $to_call = new ToCall($to_call_dir);
}
else {
    $to_call = NewGenome->new($to_call_dir);
}

my $special_pegs_file = $keep ? qq(/dev/null) : qq($to_call_dir/special_pegs);
open(SPECIAL_PEGS, ">>$special_pegs_file")
    || die "Could not append-open $special_pegs_file";

my $called_by_file = $keep ? qq(/dev/null) : qq($to_call_dir/called_by);
open(CALLED_BY,    ">>$called_by_file")
    || die "Could not append-open $called_by_file";

my @results = ();
my $contigs_file = $to_call->get_genome_dir() . qq(/contigs);

push @results, `find_selenoproteins    < $contigs_file`;
push @results, `find_selenoproteins -p < $contigs_file`;

chomp @results;
foreach my $result (@results) {
    if ($keep) {
	#...Do nothing.
    }
    else {
	
	my ($loc, $seq, $type, $func) = split /\t/o, $result;
	my $fid = $to_call->add_feature(
					-type  => qq(peg),
					-loc   => $to_call->from_seed_loc($loc),
					-seq   => $seq,
					-func  => ($func || $type),
					-annot => [qq(RAST), qq($func\nCalled by "$self".)]
					);
	
	print SPECIAL_PEGS "$fid\t$type\n";
	print CALLED_BY    "$fid\t$self\n";
    }
}

close(CALLED_BY)    || die "Could not close $called_by_file";

close(SPECIAL_PEGS) || die "Could not close $special_pegs_file";
if (!-s $special_pegs_file) {
    system("rm -f $special_pegs_file");
}

if (not $keep) {
    $to_call->export_features();
}

exit(0);
	       

MCS Webmaster
ViewVC Help
Powered by ViewVC 1.0.3