Hmmer3: Final fix for Bug #3369. Eliminated the removal of
[bioperl-live.git] / Bio / SearchIO / blast_pull.pm
blobc754745a73bfb0cc1a3a9da8cca0cad2cc9a938a
2 # BioPerl module for Bio::SearchIO::blast_pull
4 # Please direct questions and support issues to <bioperl-l@bioperl.org>
6 # Cared for by Sendu Bala <bix@sendu.me.uk>
8 # Copyright Sendu Bala
10 # You may distribute this module under the same terms as perl itself
12 # POD documentation - main docs before the code
14 =head1 NAME
16 Bio::SearchIO::blast_pull - A parser for BLAST output
18 =head1 SYNOPSIS
20 # do not use this class directly it is available through Bio::SearchIO
21 use Bio::SearchIO;
22 my $in = Bio::SearchIO->new(-format => 'blast_pull',
23 -file => 't/data/new_blastn.txt');
24 while (my $result = $in->next_result) {
25 # this is a Bio::Search::Result::BlastPullResult object
26 print "Results for ", $result->query_name(), "\n";
27 while (my $hit = $result->next_hit) {
28 print $hit->name(), "\n";
29 while (my $hsp = $hit->next_hsp) {
30 print "length is ", $hsp->length(), "\n";
35 =head1 DESCRIPTION
37 This object implements a pull-parser for BLAST output. It is fast since it
38 only does work on request (hence 'pull').
40 Currently only NCBI BLASTN and BLASTP are supported.
42 =head1 FEEDBACK
44 =head2 Mailing Lists
46 User feedback is an integral part of the evolution of this and other
47 Bioperl modules. Send your comments and suggestions preferably to
48 the Bioperl mailing list. Your participation is much appreciated.
50 bioperl-l@bioperl.org - General discussion
51 http://bioperl.org/wiki/Mailing_lists - About the mailing lists
53 =head2 Support
55 Please direct usage questions or support issues to the mailing list:
57 I<bioperl-l@bioperl.org>
59 rather than to the module maintainer directly. Many experienced and
60 reponsive experts will be able look at the problem and quickly
61 address it. Please include a thorough description of the problem
62 with code and data examples if at all possible.
64 =head2 Reporting Bugs
66 Report bugs to the Bioperl bug tracking system to help us keep track
67 of the bugs and their resolution. Bug reports can be submitted via the
68 web:
70 https://redmine.open-bio.org/projects/bioperl/
72 =head1 AUTHOR - Sendu Bala
74 Email bix@sendu.me.uk
76 =head1 APPENDIX
78 The rest of the documentation details each of the object methods.
79 Internal methods are usually preceded with a _
81 =cut
83 # Let the code begin...
85 package Bio::SearchIO::blast_pull;
87 use strict;
88 use Bio::Search::Result::BlastPullResult;
90 use base qw(Bio::SearchIO Bio::PullParserI);
92 =head2 new
94 Title : new
95 Usage : my $obj = Bio::SearchIO::blast_pull->new();
96 Function: Builds a new Bio::SearchIO::blast_pull object
97 Returns : Bio::SearchIO::blast_pull
98 Args : -fh/-file => BLAST output filename
99 -format => 'blast_pull'
100 -evalue => float or scientific notation number to be used
101 as an evalue cutoff for hits
102 -score => integer or scientific notation number to be used
103 as a score value cutoff for hits
104 -piped_behaviour => 'temp_file'|'memory'|'sequential_read'
106 -piped_behaviour defines what the parser should do if the input is
107 an unseekable filehandle (eg. piped input), see
108 Bio::PullParserI::chunk for details. Default is 'memory'.
110 =cut
112 sub _initialize {
113 my ($self, @args) = @_;
115 # don't do normal SearchIO initialization
117 my ($writer, $file, $fh, $piped_behaviour, $evalue, $score) =
118 $self->_rearrange([qw(WRITER
119 FILE FH
120 PIPED_BEHAVIOUR
121 EVALUE
122 SCORE)], @args);
123 $self->writer($writer) if $writer;
125 $self->_fields( { ( header => undef,
126 algorithm => undef,
127 algorithm_version => undef,
128 algorithm_reference => '',
129 database_name => undef,
130 database_letters => undef,
131 database_entries => undef,
132 next_result => undef,
133 evalue_cutoff => '[unset]',
134 score_cutoff => '[unset]' ) } );
136 $self->_fields->{evalue_cutoff} = $evalue if $evalue;
137 $self->_fields->{score_cutoff} = $score if $score;
139 $self->_dependencies( { ( algorithm => 'header',
140 algorithm_version => 'header',
141 database_name => 'header',
142 database_letters => 'header',
143 database_entries => 'header' ) } );
145 $self->chunk($file || $fh || $self->throw("-file or -fh must be supplied"),
146 -piped_behaviour => $piped_behaviour || 'memory');
149 sub _discover_header {
150 my $self = shift;
151 $self->_chunk_seek(0);
152 my $header = $self->_get_chunk_by_end("\nQuery=");
153 $self->{_after_header} = $self->_chunk_tell;
155 #*** won't catch all types? only support blastn/p now anyway
156 $header =~ /^(\S+) (\S+\s+\S+)/;
157 $self->_fields->{algorithm} = $1;
158 $self->_fields->{algorithm_version} = $2;
160 my ($database) = $header =~ /^Database: (.+)/sm;
162 unless ($database) {
163 # earlier versions put query before database?
164 my $header2 = $self->_get_chunk_by_end(".done\n");
165 ($database) = $header2 =~ /^Database: (.+)/sm;
168 $database =~ s/\s+(\d\S+) sequences; (\d\S+) total letters.*//s;
169 my $entries = $1;
170 my $letters = $2;
171 $database =~ s/\n//g;
172 $entries =~ s/,//g;
173 $letters =~ s/,//g;
174 $self->_fields->{database_name} = $database;
175 $self->_fields->{database_entries} = $entries;
176 $self->_fields->{database_letters} = $letters;
178 $self->_fields->{header} = 1;
181 sub _discover_next_result {
182 my $self = shift;
183 return if $self->{_after_results};
184 my $type = $self->get_field('algorithm'); # also sets _after_header if not set
186 if ($type eq 'BLASTN' || $type eq 'BLASTP') {
187 unless ($self->_sequential) {
188 $self->_chunk_seek($self->{_end_of_previous_result} || $self->{_after_header});
190 my ($start, $end) = $self->_find_chunk_by_end("\nQuery=");
191 return if ($start == $end);
193 unless ($end) {
194 $start = $self->{_end_of_previous_result} || $self->{_after_header};
195 $end = undef;
198 $self->_fields->{next_result} = Bio::Search::Result::BlastPullResult->new(-chunk => [($self->chunk, $start, $end)],
199 -parent => $self);
201 $self->{_end_of_previous_result} = $end;
203 else {
204 #*** doesn't work for the last result, needs fixing - try getting the database end chunk on failure?...
205 $self->throw("sequential mode not yet implemented");
206 my $chunk = $self->_get_chunk_by_end("\nQuery=");
207 $chunk || return;
208 $self->_fields->{next_result} = Bio::Search::Result::BlastPullResult->new(-chunk => [$chunk],
209 -parent => $self);
212 else {
213 $self->throw("Can only handle NCBI BLASTN and BLASTP right now");
217 =head2 next_result
219 Title : next_result
220 Usage : my $hit = $searchio->next_result;
221 Function: Returns the next Result from a search
222 Returns : Bio::Search::Result::ResultI object
223 Args : none
225 =cut
227 sub next_result {
228 my $self = shift;
229 my $result = $self->get_field('next_result') || return;
231 undef $self->_fields->{next_result};
233 $self->{'_result_count'}++;
234 return $result;
237 =head2 result_count
239 Title : result_count
240 Usage : my $count = $searchio->result_count
241 Function: Returns the number of results we have processed.
242 Returns : integer
243 Args : none
245 =cut
247 sub result_count {
248 my $self = shift;
249 return $self->{'_result_count'};
252 =head2 rewind
254 Title : rewind
255 Usage : $searchio->rewind;
256 Function: Allow one to reset the Result iterator to the beginning, so that
257 next_result() will subsequently return the first result and so on.
259 NB: result objects are not cached, so you will get new result objects
260 each time you rewind. Also, note that result_count() counts the
261 number of times you have called next_result(), so will not be able
262 tell you how many results there were in the file if you use rewind().
264 Returns : n/a
265 Args : none
267 =cut
269 sub rewind {
270 my $self = shift;
271 delete $self->{_end_of_previous_result};