#!/usr/local/bin/perl -w
$| =1;
use strict;
use sigtrap;
use GD;
use BPlite;
use Getopt::Std;
use vars qw($opt_x $opt_y);
getopts('x:y:');
my $usage = "
plotBlast - a simple BLAST report plotter using GD.
usage: $0 [options]
options:
-x width of image (scales)
-y height of image (does not scale)
";
die $usage unless @ARGV == 2;
my ($blast_file, $fasta_file) = @ARGV;
my ($qdef, $qseq) = getFasta($fasta_file);
# image setup
my $BEGIN = 0;
my $END = length($qseq) -1;
my $LENGTH = length($qseq);
my $WIDTH = $opt_x ? $opt_x : 1000;
my $HEIGHT = $opt_y ? $opt_y : 500;
my $STEP = 10;
my $TEXTAREA = 100;
my $MARGIN = $WIDTH - 10;
my $bpp = ($MARGIN-$TEXTAREA)/($END-$BEGIN);
my $B = $BEGIN * $bpp - $TEXTAREA;
my $img = new GD::Image($WIDTH, $HEIGHT);
# color allocation
my $black = $img->colorAllocate(0,0,0);
my $white = $img->colorAllocate(255,255,255);
my $gray = $img->colorAllocate(128,128,128);
my $red = $img->colorAllocate(255,0,0);
my $yellow = $img->colorAllocate(255,255,0);
my $green = $img->colorAllocate(0,255,0);
my $cyan = $img->colorAllocate(0,255,255);
my $blue = $img->colorAllocate(0,0,255);
# header
$img->string(gdSmallFont, 10,10, "plotBlast", $white);
# tick marks and numbers
my $MARK = $LENGTH / 10000;
if ($MARK < 1) {$MARK = 1}
else {$MARK = int($MARK)}
$MARK *= 1000;
$img->line($TEXTAREA, 10, $MARGIN, 10, $white);
for(my $i=$BEGIN;$i<=$END;$i++) {
if ($i % $MARK == 0 and $i != $BEGIN) {
$img->line($i*$bpp - $B, 10, $i*$bpp - $B, 15, $white);
$img->string(gdSmallFont, $i*$bpp-2-$B, 17, $i, $white);
}
}
# plot blast report
my $Y = 40; # spacing from top
open(BLAST, $blast_file) or die;
my $blast = new BPlite(\*BLAST);
while (my $sbjct = $blast->nextSbjct) {
my ($name) = $sbjct->name =~ /^>(\S+)/;
$img->string(gdSmallFont, 10, $Y, $name, $white);
$Y += $STEP;
while (my $hsp = $sbjct->nextHSP) {
# color ramp
my $color;
if ($hsp->bits > 200) {$color = $white}
elsif ($hsp->bits > 150) {$color = $blue}
elsif ($hsp->bits > 100) {$color = $cyan}
elsif ($hsp->bits > 80) {$color = $green}
elsif ($hsp->bits > 60) {$color = $yellow}
elsif ($hsp->bits > 40) {$color = $red}
else {$color = $gray}
my $x1 = $hsp->qb * $bpp;
my $x2 = $hsp->qe * $bpp;
my $y1 = $Y;
my $y2 = $Y + 5;
my $method = 'filledRectangle';
if ($x1 > $x2) {
($x1, $x2) = ($x2, $x1); # swqp coordiates
$method = 'rectangle'; # open rectangle on negative strand
}
$img->$method($x1 + $TEXTAREA, $y1, $x2 + $TEXTAREA, $y2, $color);
}
}
print $img->png;
sub getFasta {
my ($file) = @_;
open(FILE, $file) or die;
my ($def, @seq) = ;
chomp @seq;
my $seq = join("", @seq);
close FILE;
return ($def, $seq);
}