--- a/bin/barrnap
+++ b/bin/barrnap
@@ -4,8 +4,9 @@
use Time::Piece;
use List::Util qw(max);
use FindBin;
+use File::Temp qw/ tempfile tempdir /;
-# . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
+# . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
# global variables
my $EXE = $FindBin::RealScript;
@@ -29,13 +30,13 @@
);
my $MAXLEN = int( 1.2 * max(values %LENG) );
-# . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
+# . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
# command line options
my(@Options, $quiet, $kingdom, $threads, $evalue, $lencutoff, $reject, $incseq);
setOptions();
-# . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
+# . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
# check all is well
msg("This is $EXE $VERSION");
@@ -64,23 +65,38 @@
my $hmmdb = "$DBDIR/$kdom.hmm";
err("Can't find database: $hmmdb") unless -r $hmmdb;
msg("Using database: $hmmdb");
+my $nonfreehmmdb = "$DBDIR/nonfree/$kdom.hmm";
-# . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
+# . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
# run the external command
my $fasta = shift @ARGV;
$fasta && -r $fasta or err("Usage: $EXE <file.fasta>");
+if (-e $nonfreehmmdb) {
+ # prepare full pHMM file (with nonfree if it exists)
+ open ( HMMDB, "<", $hmmdb ) or die "Could not open file $hmmdb: $!";
+ open ( NONFREE, "<", $nonfreehmmdb );
+ my ($tmpfile, $tmpfilename) = tempfile();
+ while ( my $line = <HMMDB> ) {
+ print $tmpfile $line;
+ }
+ while ( my $line = <NONFREE> ) {
+ print $tmpfile $line;
+ }
+ $hmmdb = $tmpfilename
+}
+
msg("Scanning $fasta for $kdom rRNA genes... please wait");
my $cmd = "$NHMMER --cpu $threads -E $evalue --w_length $MAXLEN -o /dev/null --tblout /dev/stdout \Q$hmmdb\E \Q$fasta\E";
msg("Command: $cmd");
my @hits = qx($cmd 2>&1);
-# . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
+# . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
# process the output
my @feat;
-HIT:
+HIT:
foreach (@hits) {
chomp;
err("nhmmer failed to run - $_") if m/fail|error|core dump|bus error/i;
@@ -102,7 +118,7 @@
# check for incomplete alignments
my $note = '';
my $len = $end-$begin+1;
-
+
if ( $len < int($reject * $LENG{$gene}) ) {
msg("Rejecting short $len nt predicted $gene. Adjust via --reject option.");
next HIT;
@@ -120,7 +136,7 @@
];
}
-# . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
+# . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
# output a sorted GFF3
sub gff_sort {
@@ -168,13 +184,13 @@
sub show_citation {
print STDERR << "EOCITE";
-
+
If you use Barrnap in your work, please cite:
Seemann T (2013)
$EXE $VERSION : $DESC
$URL
-
+
Thank you.
EOCITE
@@ -193,7 +209,7 @@
{OPT=>"help", VAR=>\&usage, DESC=>"This help"},
{OPT=>"version", VAR=>\&version, DESC=>"Print version and exit"},
{OPT=>"citation",VAR=>\&show_citation, DESC=>"Print citation for referencing $EXE"},
- {OPT=>"kingdom=s", VAR=>\$kingdom, DEFAULT=>'bac',
+ {OPT=>"kingdom=s", VAR=>\$kingdom, DEFAULT=>'bac',
DESC=>"Kingdom: ".join(' ', values %KINGDOM) },
{OPT=>"quiet!", VAR=>\$quiet, DEFAULT=>0, DESC=>"No screen output"},
{OPT=>"threads=i", VAR=>\$threads, DEFAULT=>8, DESC=>"Number of threads/cores/CPUs to use"},
@@ -226,15 +242,15 @@
my $def = defined($_->{DEFAULT}) ? " (default '$_->{DEFAULT}')" : "";
$def = ($def ? ' (default OFF)' : '(default ON)') if $_->{OPT} =~ m/!$/;
my $opt = $_->{OPT};
- $opt =~ s/!$//;
- $opt =~ s/=s$/ [X]/;
+ $opt =~ s/!$//;
+ $opt =~ s/=s$/ [X]/;
$opt =~ s/=i$/ [N]/;
$opt =~ s/=f$/ [n.n]/;
printf STDERR " --%-15s %s%s\n", $opt, $_->{DESC}, $def;
}
else {
print STDERR "$_\n";
- }
+ }
}
exit(1);
}