[Bioperl-l] extracting data from SCF
Linus Taejoon Kwon
linusben@bawi.org
Sat, 4 Nov 2000 13:40:32 +0900
--xHFwDpU9dbj6ez1V
Content-Type: text/plain; charset=us-ascii
Content-Disposition: inline
> Hi,
> I'm trying to extract signal intensity from
> Standard Chrommatogram File . Is there a module
> wich is able to do that ?
> Thanks.
>
> Olivier.
I announced my SCF module in this mailing list..
It is not finished yet, however, you can extract
some part of my module to deal with SCF data.
Linus Taejoon Kwon
linusben@bawi.org
PS. I cannot understand how to deposit my modules to
bioperl project yet.. -.-a
--xHFwDpU9dbj6ez1V
Content-Type: text/plain; charset=us-ascii
Content-Disposition: attachment; filename="SCF.pm"
=head1 NAME
Bioinfo::SCF - Modules for SCF(Standard Chromatogram File)
=head1 SYNOPSIS
This provides the method to extract the real data from SCF binary format,
and to draw chromatogram graph via GD.
(example 1)
#!/usr/bin/perl -w
use Bioinfo::SCF;
my $scf = new Bioinfo::SCF;
my $scf_file = "./testfile.scf";
$scf->init_file($scf_file);
my $pos_list = "3 5 10 14 23 34";
my $begin = 2;
my $end = 5;
my $position = 3;
print "Content-type: image/png\n\n";
print $scf->get_chromat_raw(
-pos_list=>$pos_list,
-pos=>$position,
-begin=>$begin,
-end=>$end,
);
(example 2)
#!/usr/bin/perl -w
use Bioinfo::SCF;
my $scf = new Bioinfo::SCF;
my $scf_file = "./testfile.scf";
$scf->init_file($scf_file);
print "TOTAL SAMPLE NUMBER : ", $scf->get_sample, "\n";
print "CALLED BASE : ", $scf->get_called_base, "\n";
print "Accuracy of A : ", $scf->get_base_accuracy->{A}, "\n";
(example 3)
#!/usr/bin/perl -w
use Bioinfo::SCF;
use Bioinfo::Phred::DB;
use DBI;
use vars qw($DBH_PHRED);
my $scf = new Bioinfo::SCF;
my $phred_db = new Bioinfo::Phred::DB;
$DBH_PHRED = DBI->connect("DBI:mysql:trace_db", "userid", "password");
unless $DBH_PHRED;
my $rv = $phred_db->get_record(-no=>23);
my $pos_list = $rv->{POS}
$scf->init_db(-no=>23);
my $begin = 2;
my $end = 35;
my $position = 16;
print "Content-type: image/png\n\n";
print $scf->get_chromat(
-pos_list=>$pos_list,
-position=>$position,
-begin=>$begin,
-end=>$end,
);
$DBH_PHRED->disconnect;
=head1 DESCRIPTION
This module is based on scf.pm included on bioperl-0.6.1, made by
Aaron Mackey and Lincorn Stein. Previous scf.pm is focused on only
for sequence information, so I modified(almost new one) it to use
the entire informations on SCF file.
In addition, this module provides the chromatogram image using
gd library. If you think the java applet viewer is so heavy to
see some part of SCF file, it will be the better solution.
(example 1) Basic Usage to draw chromatogram
Like any other perl script, you should define the route for perl
interpriter. -w option is recommended to debug easily.
#!/usr/bin/perl -w
Define the SCF modules.
use Bioinfo::SCF;
Construct the instance for this modules, and define the SCF file name.
my $scf = new Bioinfo::SCF;
my $scf_file = "./testfile.scf";
If you are willing to use the SCF file directly(not from DBMS), you can
use init_file(FILENAME) method for initiation.
$scf->init_file($scf_file);
-pos_list means the base-calling position. It can be obtained by "phred" or
any other base-calling programs. If you give this value to get_chromat method,
the chromatogram contains the vertical grids for base-calling position and
the called-base labeling. If you are working with "phred" data, Bioinfo::Phred
module will help you to get this information via get_info_from_phd1 method.
-begin is the start position of chromatogram, and -end is the end position of that
(negative -end parameter means 'till the end of file').
If you'd like to remark a specific position on chromatogram, you can use -position
parameters.
my $pos_list = "3 5 10 14 23 34";
my $begin = 2;
my $end = 5;
my $position = 3;
print "Content-type: image/png\n\n";
print $scf->get_chromat_raw(
-pos_list=>$pos_list,
-pos=>$position,
-begin=>$begin,
-end=>$end,
);
NOTES: If you want to show this chromatogram to web browser, there are some
restrictions.
(1) IE(test on ver 5.0 and 5.5) can display "Content-type: image/png" directly.
However, NS(test on ver 4.72) cannot. So it will be better to use with HTML
wrapper like this:
print "<HTML><BODY>";
print "<IMG SRC=\"chromat.cgi\">\n
print "</BODY></HTML>";
(2) In general, the whole SCF file generates very large image. Sometimes it
causes an error, which does not display any image even though the file is
recieved (guess it from web browser status). So using -begin and -end parameter,
it will be better to seperate a large image to several small ones.
(example 2) Basic Usage to get sequence information
These lines are the same as example 1.
#!/usr/bin/perl -w
use Bioinfo::SCF;
my $scf = new Bioinfo::SCF;
my $scf_file = "./testfile.scf";
$scf->init_file($scf_file);
you can get the total sample number using get_sample method.
print "TOTAL SAMPLE NUMBER : ", $scf->get_sample, "\n";
get_called_base method provides the called base sequence.
print "CALLED BASE : ", $scf->get_called_base, "\n";
get_base_accuracy method provides the base quality
(for SCF ver 3.0 or higher).
print "Accuracy of A : ", $scf->get_base_accuracy->{A}, "\n";
(example 3) Integrated with database
In this module, the data parsing subroutines for each SCF data field are
seperated. So It can be used with DBMS system to store the SCF trace file
very easily. In this example, I present this with Bioinfo::Phred::DB module,
which can support the "phred" data management with MySQL.
If you like to use Bioinfo::Phred::DB, you should define $DBH_PHRED to
global variable, and use DBI module.
#!/usr/bin/perl -w
use Bioinfo::SCF;
use Bioinfo::Phred::DB;
use DBI;
use vars qw($DBH_PHRED);
Make the instance for SCF module and Phred::DB module.
my $scf = new Bioinfo::SCF;
my $phred_db = new Bioinfo::Phred::DB;
Connect to the database which has the SCF data. The table specification is
like this:
TABLE TBL_TRACE
NO INT NOT NULL AUTO_INCREMENT,
MARKER VARCHAR(64) NOT NULL,
TEMPLATE VARCHAR(64) NOT NULL,
PRIMER CHAR(1) NOT NULL,
SEQ TEXT NOT NULL,
SCORE TEXT NOT NULL,
POS TEXT NOT NULL,
POLY MEDIUMTEXT NOT NULL,
SCF_SAMPLE_NO INT NOT NULL,
SCF_BASE_NO INT NOT NULL,
SCF_VERSION VARCHAR(8) NOT NULL,
SCF_SAMPLE_SIZE INT NOT NULL,
SCF_COMMENTS_SIZE INT NOT NULL,
SCF_PRIVATES_SIZE INT NOT NULL,
SCF_HEADER BLOB NOT NULL,
SCF_SAMPLE BLOB NOT NULL,
SCF_BASE BLOB NOT NULL,
SCF_COMMENTS BLOB NOT NULL,
SCF_PRIVATES BLOB NOT NULL,
MODIFIED INT NOT NULL,
PRIMARY KEY (NO),
INDEX idx_snp (MARKER,TEMPLATE,PRIMER)
Values for SCF_* columns are added by SCF module, and the others are
added by other module (ex. Bioinfo::Phred module with .phd.1 or .poly
file after "phred" operation).
$DBH_PHRED = DBI->connect("DBI:mysql:trace_db", "userid", "password");
unless $DBH_PHRED;
Based on NO (serial number of record), you can get pos_list data.
my $rv = $phred_db->get_record(-no=>23);
my $pos_list = $rv->{POS}
init_db method with NO value, you can initiate the $scf instance like
init_file method with SCF file name.
$scf->init_db(-no=>23);
Other things are the same as example 1.
my $begin = 2;
my $end = 35;
my $position = 16;
print "Content-type: image/png\n\n";
print $scf->get_chromat(
-pos_list=>$pos_list,
-position=>$position,
-begin=>$begin,
-end=>$end,
);
Finally, the connection with DBMS is discarded.
$DBH_PHRED->disconnect;
=head1 Requirements
The GD.pm, version 1.30 or higher
The gd graphics library, version 1.8.3 or higher
The PNG graphics library
The zlib compression library
=head1 FEEDBACK
If you have some opinions or find some bugs, please let me know that
via email to Linus Taejoon Kwon(See AUTHORS).
=head1 AUTHORS
=head2 Linus Taejoon Kwon
Email linusben@bawi.org
=head2 Aaron Mackey
=head2 Lincorn Stein
=head1 Change Logs
=head2 ver 0.1 (29/JUL/2000)
Initiation (Linus)
=head2 ver 0.2 (06/OCT/2000)
Include init_db subroutine (Linus)
Include get_chromat_text for debugging (Linus)
Change almost all subroutine name (Linus)
=head1 Copyright
This module is open source software ; you can redistribute it and/or
modify it under the same terms as Perl itself.
=head1 See Also
=head2 web-phred project
http://white.bawi.net/~bioinfo/snp/
=head1 Miscellaneous
=head2 SCF-RFC
If you want to get more information about SCF file format, refer
http://www.mrc-lmb.cam.ac.uk/pubseq/scf-rfc.html.
=head2 Java Applet
There are some jave applet viewer for SCF file. One of them is
found at http://www.embl-heidelberg.de/~toldo/scf/SCF_viewer_applet.html.
=cut
package Bioinfo::SCF;
use Bioinfo::Phred::DB;
use GD;
use strict;
use vars qw($DBH_PHRED);
BEGIN {
use Exporter ();
use vars qw($VERSION @ISA @EXPORT @EXPORT_OK %EXPORT_TAGS);
$VERSION = do { my @r=("0.2"=~/\d+/g); sprintf "%d."."%02d"x$#r,@r };
@ISA = qw(Exporter);
@EXPORT = qw();
@EXPORT_OK = qw($DBH_PHRED);
%EXPORT_TAGS = ();
}
=head2 new
Title : new
Usage : $scf = new Bioinfo::SCF;
Function: Construct the instance
Returns : None.
Args : None.
=cut
sub new
{
my ($class) = @_;
bless {
_file => $_[1],
_header => $_[2],
_sample => $_[3], # number of sample
_base => $_[4], # number of bases
_version => $_[5], # version of SCF format
_sample_size => $_[6], # sample size (precision)
_comments_size => $_[7], # comments size
_privates_size => $_[8], # privates size
_struct_sample => $_[9],
_struct_base => $_[10],
_comments => $_[11],
_privates => $_[12],
}, $class;
}
sub get_file { return $_[0]->{_file}}
sub get_header { return $_[0]->{_header}}
sub get_sample { return $_[0]->{_sample}}
sub get_base { return $_[0]->{_base}}
sub get_version { return $_[0]->{_version}}
sub get_sample_size { return $_[0]->{_sample_size}}
sub get_comments_size { return $_[0]->{_comments_size}}
sub get_privates_size { return $_[0]->{_privates_size}}
sub get_struct_sample { return $_[0]->{_struct_sample}}
sub get_struct_base { return $_[0]->{_struct_base}}
sub get_comments { return $_[0]->{_comments}}
sub get_privates { return $_[0]->{_privates}}
=head2 init_file
Title : init_file
Usage : $scf->init_file("test.scf");
Function: Initialization for SCF data via file. All the SCF data is
stored as internal parameter of module.
Returns : None.
Args : filename SCF filename
=cut
sub init_file
{
my ($self, $file) = @_;
my $header;
my @header_spare;
my $tmp;
if($file) {
$self->{_file} = $file;
open(FILE, $file);
read FILE, $self->{_header}, 128;
( $header->{scf},
$header->{sample},
$header->{sample_offset},
$header->{base},
$header->{base_left_clip},
$header->{base_right_clip},
$header->{base_offset},
$header->{comments_size},
$header->{comments_offset},
$header->{version},
$header->{sample_size},
$header->{code_set},
$header->{privates_size},
$header->{privates_offset},
@header_spare ) = unpack "a4 NNNNNNNN a4 NNNN N18", $self->get_header;
$self->{_sample} = $header->{sample};
$self->{_base} = $header->{base};
$self->{_version} = $header->{version};
$self->{_sample_size} = $header->{sample_size};
$self->{_comments_size} = $header->{comments_size};
$self->{_privates_size} = $header->{privates_size};
read FILE, $self->{_struct_sample}, ($header->{sample} * $header->{sample_size} * 4);
read FILE, $self->{_struct_base}, ($header->{base} * 12);
read FILE, $self->{_comments}, $header->{comments_size};
read FILE, $self->{_privates}, $header->{privates_size};
close(FILE);
}
}
=head2 init_db
Title : init_db
Usage : $scf->init_db(-no=>23);
Function: Initialization for SCF data via DB. All the SCF data is
stored as internal parameter of module. For more information
about database integration, see example 3 above and "See Also".
Returns : None.
Args : -no Serial Number for trace data
=cut
sub init_db
{
my ($self, %param) = @_;
my $header;
my @header_spare;
if($param{"-no"}) {
my $phred_db = new Bioinfo::Phred::DB;
my $rv = $phred_db->get_record(-no=>$param{"-no"});
$self->{_header} = $rv->{'SCF_HEADER'};
( $header->{scf},
$header->{sample},
$header->{sample_offset},
$header->{base},
$header->{base_left_clip},
$header->{base_right_clip},
$header->{base_offset},
$header->{comments_size},
$header->{comments_offset},
$header->{version},
$header->{sample_size},
$header->{code_set},
$header->{privates_size},
$header->{privates_offset},
@header_spare ) = unpack "a4 NNNNNNNN a4 NNNN N18", $self->get_header;
$self->{_sample} = $header->{sample};
$self->{_base} = $header->{base};
$self->{_version} = $header->{version};
$self->{_sample_size} = $header->{sample_size};
$self->{_comments_size} = $header->{comments_size};
$self->{_privates_size} = $header->{privates_size};
$self->{_struct_sample} = $rv->{'SCF_SAMPLE'};
$self->{_struct_base} = $rv->{'SCF_BASE'};
$self->{_comments} = $rv->{'SCF_COMMENTS'};
$self->{_privates} = $rv->{'SCF_PRIVATES'};
}
}
=head2 get_sample_info
Title : get_sample_info
Usage : print $scf->get_sample_info->{T};
Function: Get the trace information about each nucleic acids with
white space as seperator.
Returns : hashref for each nucleic acid's trace information with
white space as seperator. Based on these values, the
chromatogram image is made.
Args : None.
=cut
sub get_sample_info
{
my $self = shift;
my $rv;
my $sample = $self->get_sample;
my $sample_size = $self->get_sample_size;
my $struct_sample = $self->get_struct_sample;
my @tmp_base = ("A","C","G","T");
if( $self->get_version < 3.0 ) {
# For version 1.0 and 2.0
for(my $i=0;$i<4;$i++)
{
$rv->{$tmp_base[$i]} = "";
my @tmp_sample;
for(my $j=0;$j<$sample;$j++)
{
my $tmp = substr($struct_sample, ($i*$sample + $j), $sample_size);
$tmp = unpack "H2", $tmp;
$tmp = hex $tmp;
$rv->{$tmp_base[$i]} .= " $tmp";
}
}
} else {
# For version 3.0 and over
for(my $i=0;$i<4;$i++)
{
my @tmp_sample;
for(my $j=0;$j<$sample;$j++)
{
my $tmp = substr($struct_sample, ($i*$sample + $j), $sample_size);
$tmp = unpack "H2", $tmp;
$tmp = hex $tmp;
push(@tmp_sample, $tmp);
}
my $p_sample = 0;
for(my $j = 0; $j < $sample ; $j++)
{
$tmp_sample[$j] += $p_sample;
while($tmp_sample[$j] >= 256) {
$tmp_sample[$j] -= 256;
}
$p_sample = $tmp_sample[$j];
}
$p_sample = 0;
for(my $j = 0; $j < $sample ; $j++)
{
$tmp_sample[$j] += $p_sample;
while($tmp_sample[$j] >= 256) {
$tmp_sample[$j] -= 256;
}
$p_sample = $tmp_sample[$j];
}
$rv->{$tmp_base[$i]} = join(" ",@tmp_sample);
}
}
return $rv;
}
=head2 get_peak_offset
Title : get_peak_offset
Usage : print $scf->get_peak_offset
Function: Get peak offset information on SCF file
Returns : raw peak_offset information on SCF file
Args : None.
=cut
sub get_peak_offset {
my $self = shift;
my $struct_base = $self->get_struct_base;
my $base = $self->get_base;
my $tmp = substr($struct_base, 0, $base*4);
my $rv = unpack "I*", $tmp;
return $rv;
}
=head2 get_base_accuracy
Title : get_base_accuracy
Usage : print $scf->get_base_accuracy->{A};
Function: Get the base-calling accuracy information on SCF file.
It is guaranteed for SCF ver 3.0 or higher.
Returns : hashref for each nucleic acid's accuracy sequence with
white space as seperator.
Args : None.
=cut
sub get_base_accuracy {
my $self = shift;
my $rv;
my @tmp_base = ('A','C','G','T');
my $struct_base = $self->get_struct_base;
my $base = $self->get_base;
for(my $i=0;$i<4;$i++)
{
my $tmp = substr($struct_base, $base*($i+4), $base);
my @tmp_accuracy = unpack "c*", $tmp;
$rv->{$tmp_base[$i]} = join(" ",@tmp_accuracy);
}
return $rv;
}
=head2 get_called_base
Title : get_called_base
Usage : print $scf->get_called_base
Function: Get the sequence information on SCF file.
Returns : plain text for called base sequence(without white space).
Args : None.
=cut
sub get_called_base {
my $self = shift;
my $rv;
my $struct_base = $self->get_struct_base;
my $base = $self->get_base;
for(my $i=0;$i<$base;$i++)
{
my $tmp = substr($struct_base, $base*8+$i, 1);
$tmp = unpack "a", $tmp;
$rv .= $tmp;
}
return $rv;
}
=head2 get_reserved
Title : get_reserved
Usage : Not yet.
Function: Get the reserved information on SCF Header
Returns : Not yet.
Args : Not yet.
=cut
sub get_reserved {
# Under Constrcution...
# No necessary... yet...
}
=head2 get_comments_info
Title : get_comments_info
Usage : print $scf->get_comments_info
Function: Get comments files on SCF file
Returns : Raw "comments" field information
Args : None.
=cut
sub get_comments_info {
my $self = shift;
my $comments_size = $self->get_comments_size;
my $comments = $self->get_comments;
my $tmp = substr($comments, 0, $comments_size);
my $rv = unpack "a*", $tmp;
return $rv;
}
=head2 get_privates_info
Title : get_privates_info
Usage : print $scf->get_privates_info
Function: Get privates files on SCF file
Returns : Raw "privates" field information
Args : None.
=cut
sub get_privates_info {
my $self = shift;
my $privates_size = $self->get_privates_size;
my $privates = $self->get_privates;
my $tmp = substr($privates, 0, $privates_size);
my $rv = unpack "a*", $tmp;
return $rv;
}
=head2 get_chromat
Title : get_chromat
Usage : print $scf->get_chormat(-begin=> 5,
-end=> 300,
-pos_list=> $pos_list,
-position=> 125);
Function: It provides the png image which contains the SCF chromatogram
information.
Returns : PNG image file.
Args : -begin start position of chromatogram to be shown
-end end position of chromatogram to be shown
-pos_list base_calling position, seperated by white space
-position specific position willing to be marked on
chromatogram
=cut
sub get_chromat {
my ($self, %param) = @_;
# Get Sample Information
my $sample_info = $self->get_sample_info;
my @sample_A = split(/\s+/,$sample_info->{A});
my @sample_C = split(/\s+/,$sample_info->{C});
my @sample_G = split(/\s+/,$sample_info->{G});
my @sample_T = split(/\s+/,$sample_info->{T});
my $sample = $self->get_sample;
# Define Image Constants
my $set;
if($param{'-begin'}) {
$set->{begin} = $param{'-begin'};
} else {
$set->{begin} = 0;
}
if($set->{begin} < 0) {
$set->{begin} = 0;
}
if($param{'-end'} ) {
$set->{end} = $param{'-end'};
} else {
$set->{end} = $sample;
}
if(($set->{end} > $sample) || ($set->{end} <= 0) ) {
$set->{end} = $sample;
}
$set->{height} = 120;
$set->{space} = 5;
$set->{left_margin} = 5;
$set->{right_margin} = 5;
$set->{top_margin} = 5;
$set->{bottom_margin} = 37;
$set->{label} = 10;
$set->{width} = ($set->{end} - $set->{begin}) * $set->{space};
my $pic_height = $set->{height}+$set->{top_margin}+$set->{bottom_margin};
my $pic_width = $set->{width}+$set->{left_margin}+$set->{right_margin};
my $bottom = $pic_height - $set->{bottom_margin};
# If there is calling position information for each base, it will be drawn.
my $is_pos_list = 0;
my $idx = 1;
my @called_base;
my @pos;
if($param{'-pos_list'}) {
$is_pos_list = 1;
@called_base = split("",$self->get_called_base);
@pos = split(/\s+/,$param{'-pos_list'});
while($pos[$idx] < $set->{begin}) {
$idx++;
}
}
# Create Image Instance
my $im = new GD::Image($pic_width, $pic_height);
# Allocate some colors
my $white = $im->colorAllocate(255,255,255);
my $black = $im->colorAllocate(0,0,0);
my $grey = $im->colorAllocate(204,204,204);
my $red = $im->colorAllocate(255,0,0);
my $green = $im->colorAllocate(0,204,0);
my $yellow = $im->colorAllocate(255,204,0);
my $blue = $im->colorAllocate(0,0,255);
my %color;
$color{A} = $green;
$color{C} = $blue;
$color{G} = $grey;
$color{T} = $red;
$color{N} = $yellow;
# make the background transparent and interlaced
$im->transparent($white);
$im->interlaced('true');
# Put a black frame around the picture
$im->filledRectangle( 0, 0, $pic_width, $pic_height, $black);
for(my $i=$set->{begin};$i<$set->{end};$i++)
{
my $new_i = $i - $set->{begin};
# Set New Scale
$sample_A[$i] = $sample_A[$i] * $set->{height} / 256;
$sample_C[$i] = $sample_C[$i] * $set->{height} / 256;
$sample_G[$i] = $sample_G[$i] * $set->{height} / 256;
$sample_T[$i] = $sample_T[$i] * $set->{height} / 256;
# Position Labeling
if($i % 20 == 0) {
$im->string( gdSmallFont,
$new_i*$set->{space}, $bottom,
$i, $grey);
}
if($is_pos_list == 1 && $pos[$idx] == $i) {
# If two or more bases are called at the same position.
if( ($idx < $#pos - 1) && ($pos[$idx+1] eq $pos[$idx]) ) {
$idx++;
}
# Individual position check ; vertical line and labeling
if((defined $param{'-pos'}) && ($i eq $param{'-pos'}))
{
$im->line( $new_i*$set->{space}+$set->{left_margin}-1,
$set->{top_margin},
$new_i*$set->{space}+$set->{left_margin}-1,
$bottom,
$white);
$im->line( $new_i*$set->{space}+$set->{left_margin},
$set->{top_margin},
$new_i*$set->{space}+$set->{left_margin},
$bottom,
$white);
$im->line( $new_i*$set->{space}+$set->{left_margin}+1,
$set->{top_margin},
$new_i*$set->{space}+$set->{left_margin}+1,
$bottom,
$white);
$im->string( gdGiantFont,
$new_i*$set->{space}, $bottom+15,
$called_base[$idx], $color{$called_base[$idx]});
} else {
$im->line( $new_i*$set->{space}+$set->{left_margin},
$set->{top_margin},
$new_i*$set->{space}+$set->{left_margin},
$bottom,
$grey);
$im->string( gdMediumBoldFont,
$new_i*$set->{space}, $bottom+15,
$called_base[$idx-1], $color{$called_base[$idx-1]});
}
# increase idx
if($idx < $#pos) {
$idx++;
}
}
if($new_i == 0) {
$im->setPixel($set->{left_margin},$bottom,$color{A});
$im->setPixel($set->{left_margin},$bottom,$color{C});
$im->setPixel($set->{left_margin},$bottom,$color{G});
$im->setPixel($set->{left_margin},$bottom,$color{T});
} else {
$im->line( ($new_i-1)*$set->{space}+$set->{left_margin},
$bottom-$sample_A[$i-1],
$new_i*$set->{space}+$set->{left_margin},
$bottom - $sample_A[$i],
$color{A});
$im->line( ($new_i-1)*$set->{space}+$set->{left_margin},
$bottom-$sample_C[$i-1],
$new_i*$set->{space}+$set->{left_margin},
$bottom-$sample_C[$i],
$color{C});
$im->line( ($new_i-1)*$set->{space}+$set->{left_margin},
$bottom-$sample_G[$i-1],
$new_i*$set->{space}+$set->{left_margin},
$bottom-$sample_G[$i],
$color{G});
$im->line( ($new_i-1)*$set->{space}+$set->{left_margin},
$bottom-$sample_T[$i-1],
$new_i*$set->{space}+$set->{left_margin},
$bottom-$sample_T[$i],
$color{T});
}
}
return $im->png;
}
=head2 get_chromat_text
Title : get_chromat_text
Usage : $scf->get_chormat_text(-begin=> 5,
-end=> 300,
-pos_list=> $pos_list,
-position=> 125);
Function: For debugging, it provides the informations which will
be shown on chromatogram image.
Returns : NO. Just print out the values.
Args : -begin start position of chromatogram to be shown
-end end position of chromatogram to be shown
-pos_list base_calling position, seperated by white space
-position specific position willing to be marked on
chromatogram
=cut
sub get_chromat_text {
my ($self, %param) = @_;
# Get Sample Information
my $sample_info = $self->get_sample_info;
my @sample_A = split(/\s+/,$sample_info->{A});
my @sample_C = split(/\s+/,$sample_info->{C});
my @sample_G = split(/\s+/,$sample_info->{G});
my @sample_T = split(/\s+/,$sample_info->{T});
my $sample = $self->get_sample;
# Define Image Constants
my $set;
if($param{'-begin'}) {
$set->{begin} = $param{'-begin'};
} else {
$set->{begin} = 0;
}
if($set->{begin} < 0) {
$set->{begin} = 0;
}
if($param{'-end'} ) {
$set->{end} = $param{'-end'};
} else {
$set->{end} = $sample;
}
if(($set->{end} > $sample) || ($set->{end} <= 0) ) {
$set->{end} = $sample;
}
$set->{height} = 120;
$set->{space} = 5;
$set->{left_margin} = 5;
$set->{right_margin} = 5;
$set->{top_margin} = 5;
$set->{bottom_margin} = 37;
$set->{label} = 10;
$set->{width} = ($set->{end} - $set->{begin}) * $set->{space};
my $pic_height = $set->{height}+$set->{top_margin}+$set->{bottom_margin};
my $pic_width = $set->{width}+$set->{left_margin}+$set->{right_margin};
my $bottom = $pic_height - $set->{bottom_margin};
print "WIDTH : $pic_width\n";
print "HEIGHT : $pic_height\n";
# If there is calling position information for each base, it will be drawn.
my $is_pos_list = 0;
my $idx = 1;
my @called_base;
my @pos;
if($param{'-pos_list'}) {
$is_pos_list = 1;
@called_base = split("",$self->get_called_base);
@pos = split(/\s+/,$param{'-pos_list'});
while($pos[$idx] < $set->{begin}) {
$idx++;
}
}
for(my $i=$set->{begin};$i<$set->{end};$i++)
{
my $new_i = $i - $set->{begin};
# Set New Scale
$sample_A[$i] = $sample_A[$i] * $set->{height} / 256;
$sample_C[$i] = $sample_C[$i] * $set->{height} / 256;
$sample_G[$i] = $sample_G[$i] * $set->{height} / 256;
$sample_T[$i] = $sample_T[$i] * $set->{height} / 256;
# Position Labeling
if($i % 20 == 0) {
print "LABEL : $i <BR>\n";
}
if($is_pos_list == 1 && $pos[$idx] == $i) {
# If two or more bases are called at the same position.
if( ($idx < $#pos - 1) && ($pos[$idx+1] eq $pos[$idx]) ) {
$idx++;
}
# Individual position check ; vertical line and labeling
if((defined $param{'-pos'}) && ($i eq $param{'-pos'}))
{
print "POSITION SITE : ",$new_i*$set->{space} + $set->{left_margin}, "[ $called_base[$idx-1] ]\n";
} else {
print "SITE : ", $new_i*$set->{space} + $set->{left_margin}
, "[ $called_base[$idx-1] ]\n";
}
# increase idx
if($idx < $#pos) {
$idx++;
}
}
if($new_i == 0) {
print "A : Start \t";
print "C : Start \t";
print "G : Start \t";
print "T : Start \t";
} else {
print "FROM : ",($new_i-1)*$set->{space}+$set->{left_margin},"\n";
print "TO : ",$new_i*$set->{space}+$set->{left_margin},"\n";
print "HEIGHT [ $i ] \t";
print "A : ", $sample_A[$i-1],":",$sample_A[$i],"\t";
print "C : ", $sample_C[$i-1],":",$sample_C[$i],"\t";
print "G : ", $sample_G[$i-1],":",$sample_G[$i],"\t";
print "T : ", $sample_T[$i-1],":",$sample_T[$i],"\n";
}
}
}
1;
--xHFwDpU9dbj6ez1V--