sub new
{ my ($class,@args) = @_;
my $self = $class->SUPER::new(@args);
$self->program('fgenesh') if(!$self->program);
$self->matrix('hum.dat') if(!$self->matrix);
return $self;
}
} |
sub parse_results
{ my ($self, $results) = @_;
if(!$results){
$results = $self->resultsfile;
}
open(OUT, "<".$results) or throw("FAILED to open ".$results.
"Genscan:parse_results");
if (<OUT> =~ m|no reliable predictions|i ){ print STDERR "No genes predicted\n"; return;
}
while (<OUT>) {
chomp;
my $flag = 0;
if (/^\s*-[-\s]+$/) {
GENES : while (<OUT>) {
my @lines;
until (/^$/) {
if (/CDSl|CDSi|CDSf|CDSo/i) {
my @element = split;
push @lines,\@ element;
throw("Unable to parse fgenesh ouput (".@element.
") Line: $_\n") unless (scalar(@element) == 11);
}elsif (/Predicted protein/i) {
$flag = 1;
last GENES ;
}
$_ = <OUT>;
}
if ($lines[0]->[1] eq '+') {
@lines = sort {$a->[3] <=> $b->[3]} @lines;
} elsif ($lines[0]->[1] eq '-') {
@lines = reverse sort {$a->[3] <=> $b->[3]} @lines;
}
my $exon_num=1;
foreach my $line (@lines) {
my ($start, $end, $strand, $score, $pvalue, $phase, $seqname);
$seqname = $line->[0]+($exon_num/100); if ($line->[1] eq '+') {
$start = $line->[3];
$end = $line->[5];
$strand = 1;
$phase = (3-($line->[7]-$line->[3]))% 3;
} elsif ($line->[1] eq '-') {
$start = $line->[3];
$end = $line->[5];
$strand = -1;
$phase = (3-($line->[5]-$line->[9]))% 3;
}
$score = $line->[6];
my $exon = $self->feature_factory->create_prediction_exon
($start, $end, $strand, $score, '0', $phase, $seqname,
$self->query, $self->analysis);
$self->exon_groups($line->[0], $exon);
$exon_num++;
}
if ( $flag==1 ) {
last;
}
}
}
}
$self->create_transcripts();
close(OUT) or throw("FAILED to close ".$results.
"Fgenesh:parse_results");
}
1; } |