X-Git-Url: http://source.jalview.org/gitweb/?a=blobdiff_plain;f=binaries%2Fsrc%2FViennaRNA%2Finterfaces%2FPerl%2FRNAfold.pl;fp=binaries%2Fsrc%2FViennaRNA%2Finterfaces%2FPerl%2FRNAfold.pl;h=58d60108e9e0ea75743457c2abcc738cac546c65;hb=7522ace91fc0804a9719dbac9f68bc8154da3132;hp=0000000000000000000000000000000000000000;hpb=8116c0444fe98e8eb21bcdd8ded06e1429085823;p=jabaws.git diff --git a/binaries/src/ViennaRNA/interfaces/Perl/RNAfold.pl b/binaries/src/ViennaRNA/interfaces/Perl/RNAfold.pl new file mode 100644 index 0000000..58d6010 --- /dev/null +++ b/binaries/src/ViennaRNA/interfaces/Perl/RNAfold.pl @@ -0,0 +1,144 @@ +#!/usr/bin/perl +# -*-Perl-*- +# Last changed Time-stamp: <2006-04-04 14:51:06 ivo> + +# Perl implementation of RNAfold +# This is an example script demonstrating how to use the RNA Perl module + +use RNA; +use Getopt::Long; +use strict; +use warnings; + + Getopt::Long::config("no_ignore_case"); + +use vars qw/$opt_debug $opt_v $ParamFile $pf $ns_bases/; +my $sfact=1.07; +&usage() unless GetOptions("p|p1" => \$pf, + "p0" => sub {$pf=1; $RNA::do_backtrack=0}, + "C" => \$RNA::fold_constrained, + "T=f" => \$RNA::temperature, + "4" => sub {$RNA::tetra_loop = 0}, + "d|d0" => sub {$RNA::dangles=0}, + "d2" => sub {$RNA::dangles=2},, + "noGU" => \$RNA::noGU, + "noCloseGU" => \$RNA::no_closingGU, + "noLP" => \$RNA::noLonelyPairs, + "e=i" => \$RNA::energy_set, + "P=s" => \$ParamFile, + "nsp=s" => \$ns_bases, + "S=f" => \$sfact); + + RNA::read_parameter_file($ParamFile) if ($ParamFile); + +if ($ns_bases) { + $RNA::nonstandards = ""; + foreach my $p ( split(/,/, $ns_bases) ) { + if ($p =~ s/^-//) { + $RNA::nonstandards .= reverse($p) + } + $RNA::nonstandards .= $p; + } + print "$RNA::nonstandards\n"; +} + +my $istty = (-t STDIN) && (-t STDOUT); +if (($RNA::fold_constrained)&&($istty)) { +print < : base i is paired with a base j>i +matching brackets ( ): base i pairs base j +END +} + +if ($istty) { + print "\nInput string (upper or lower case); @ to quit\n"; + for (1..8) { print "....,....$_";} + print "\n"; +} +my $fname; +while (<>) { # main loop: continue until end of file + my ($string, $cstruc, $structure, $min_en); + # skip comment lines and get filenames + if (/^>\s*(\S*)/) { + $fname = $1; + next; + } + last if (/^@/); + + if (/(\S+)/) { + $string = $1; + } else { + next; + } + + $string = uc($string); + my $length = length($string); + printf("length = %d\n", $length) if ($istty); + + if ($RNA::fold_constrained) { + $_ = <>; + $cstruc = $1 if (/(\S+)/); + die("constraint string has wrong length") + if (length($cstruc)!=$length); + ($structure, $min_en) = RNA::fold($string, $structure); + } else { + ($structure, $min_en) = RNA::fold($string); + } + print "$string\n$structure"; + if ($istty) { + printf("\n minimum free energy = %6.2f kcal/mol\n", $min_en); + } else { + printf(" (%6.2f)\n", $min_en); + } + my $ffname = ($fname) ? ($fname . '_ss.ps') : 'rna.ps'; + RNA::PS_rna_plot($string, $structure, $ffname); + + if ($pf) { + + # recompute with dangles as in pf_fold() + $RNA::dangles=2 if ($RNA::dangles); + $min_en = RNA::energy_of_struct($string, $structure); + + my $kT = ($RNA::temperature+273.15)*1.98717/1000.; # in Kcal + $RNA::pf_scale = exp(-($sfact*$min_en)/$kT/$length); + print STDERR "scaling factor $RNA::pf_scale\n" if ($length>2000); + + $structure = $cstruc if ($RNA::fold_constrained); + my $energy = RNA::pf_fold($string, $structure); + + if ($RNA::do_backtrack) { + print $structure; + printf(" [%6.2f]\n", $energy) if (!$istty); + print "\n"; + } + if (($istty)||(!$RNA::do_backtrack)) { + printf(" free energy of ensemble = %6.2f kcal/mol\n", $energy); + printf(" frequency of mfe structure in ensemble %g\n", + exp(($energy-$min_en)/$kT)); + } + + if ($RNA::do_backtrack) { + $ffname = ($fname)?($fname . "_dp.ps"):"dot.ps"; + &RNA::PS_dot_plot($string, $ffname); + } + } + undef $fname; +} + + RNA::free_pf_arrays() if ($pf); + RNA::free_arrays(); + +sub usage() +{ + die("usage: " . + "RNAfold [-p[0]] [-C] [-T temp] [-4] [-d] [-noGU] [-noCloseGU]\n" . + " [-e e_set] [-P paramfile] [-nsp pairs] [-S scale]"); +} + + +# End of file