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

View of /FigKernelScripts/make_PHOBs_for_subsystems.pl

Parent Directory Parent Directory | Revision Log Revision Log


Revision 1.4 - (download) (as text) (annotate)
Thu May 12 20:38:37 2005 UTC (14 years, 7 months ago) by overbeek
Branch: MAIN
CVS Tags: merge-trunktag-bobdev_news-2, Root-bobdev_news, merge-bobdev_news-1, caBIG-dataload-0, merge-trunktag-bobdev_news-1, merge-bodev_news-3, caBIG-00-00-00, merge-bobdev_news-2, merge-trunktag-bodev_news-3
Branch point for: Branch-bobdev_news
Changes since 1.3: +11 -4 lines
add tool to generate gene blocks for variation analysis

# -*- perl -*-

use FIG;
my $fig = new FIG;

$usage = "usage: make_PHOBs_for_subsystems Dir [SimilarityCutoff FracCoverage] < FileOf_SS-Role-PEG-Func";

(
 ($dir = shift @ARGV) && (! -d $dir)
)
    || die $usage;

if (@ARGV == 2)
{
    $sim_cutoff = shift @ARGV;
    $frac_cov   = shift @ARGV;
}

mkdir($dir,0777) || die "could not make $dir";
while (defined($_ = <STDIN>))
{
    chop;
    ($ss,$role,$peg) = split(/\t/,$_);
    push(@{$subH{$ss}->{$role}},$peg);
}

$tmpF = "/tmp/pegs$$.fasta";

open(INDEXSS,">$dir/index") || die "could not open $dir/index";
$ssN = 0;
foreach $ss (sort keys(%subH))
{
    ++$ssN;
    mkdir("$dir/$ssN",0777) || die "could not make $dir/$ssN";
    print INDEXSS "$ssN\t$ss\n";

    open(INDEXR,">$dir/$ssN/index") || die "could not open $dir/$ssN/index";
    $roleN = 0;
    foreach $role (sort keys(%{$subH{$ss}}))
    {
	++$roleN;
	
	$pegs = $subH{$ss}->{$role};
	if (@$pegs < 8)
	{
	    print STDERR "Skipping $dir/$ssN/$roleN, role $role --- only ", (scalar @$pegs), " pegs\n";
	    next;
	}
	
	print INDEXR "$roleN\t$role\n";
	
	open(TMP,">$tmpF") || die "could not open $tmpF";
	foreach $peg (@$pegs)
	{
	    if ($pseq = $fig->get_translation($peg))
	    {
		print TMP ">$peg\n$pseq\n";
	    }
	    else
	    {
		print STDERR "failed to get translation for $peg\n";
	    }
	}
	close(TMP);
	
	mkdir("$dir/$ssN/$roleN",0777) || die "could not make $dir/$ssN/$roleN";
	&FIG::run("split_and_trim_sequences $dir/$ssN/$roleN/split_info $sim_cutoff $frac_cov < $tmpF");
	
	if (-s "$dir/$ssN/$roleN/split_info/set.sizes")
	{
	    open(SZ,"<$dir/$ssN/$roleN/split_info/set.sizes") || die " could not open $dir/$ssN/$roleN/split_info/set.sizes";
	    while (defined($_ = <SZ>))
	    {
		if (($_ =~ /^(\d+)\t(\d+)/) && ($2 > 3))
		{
		    $n = $1;
		    &FIG::run("make_phob_from_seqs $dir/$ssN/$roleN/$n < $dir/$ssN/$roleN/split_info/$n");
		}
	    }
	    close(SZ);
	}
	else
	{
	    system("rmdir $dir/$ssN/$roleN") && warn "Could not rmdir $dir/$ssN/$roleN";
	    print STDERR "Skipping role $role in $dir/$ssN/$roleN\n";
	}
    }
    close(INDEXR);
}
close(INDEXSS);
unlink($tmpF);

MCS Webmaster
ViewVC Help
Powered by ViewVC 1.0.3