Bio/SearchIO/hmmer2.pm

   1 #
   2 # BioPerl module for Bio::SearchIO::hmmer2
   3 #
   4 # Please direct questions and support issues to <bioperl-l@bioperl.org>
   5 #
   6 # Cared for by Jason Stajich <jason@bioperl.org>
   7 #
   8 # Copyright Jason Stajich
   9 #
  10 # You may distribute this module under the same terms as perl itself
  11
  12 # POD documentation - main docs before the code
  13
  14 =head1 NAME
  15
  16 Bio::SearchIO::hmmer2 - A parser for HMMER output (hmmpfam, hmmsearch)
  17
  18 =head1 SYNOPSIS
  19
  20     # do not use this class directly it is available through Bio::SearchIO
  21     use Bio::SearchIO;
  22     my $in = Bio::SearchIO->new(-format => 'hmmer2',
  23                                -file   => 't/data/L77119.hmmer');
  24     while( my $result = $in->next_result ) {
  25         # this is a Bio::Search::Result::HMMERResult object
  26         print $result->query_name(), " for HMM ", $result->hmm_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";
  31             }
  32         }
  33     }
  34
  35 =head1 DESCRIPTION
  36
  37 This object implements a parser for HMMER output.
  38
  39 =head1 FEEDBACK
  40
  41 =head2 Mailing Lists
  42
  43 User feedback is an integral part of the evolution of this and other
  44 Bioperl modules. Send your comments and suggestions preferably to
  45 the Bioperl mailing list.  Your participation is much appreciated.
  46
  47   bioperl-l@bioperl.org                  - General discussion
  48   http://bioperl.org/wiki/Mailing_lists  - About the mailing lists
  49
  50 =head2 Support
  51
  52 Please direct usage questions or support issues to the mailing list:
  53
  54 I<bioperl-l@bioperl.org>
  55
  56 rather than to the module maintainer directly. Many experienced and
  57 reponsive experts will be able look at the problem and quickly
  58 address it. Please include a thorough description of the problem
  59 with code and data examples if at all possible.
  60
  61 =head2 Reporting Bugs
  62
  63 Report bugs to the Bioperl bug tracking system to help us keep track
  64 of the bugs and their resolution. Bug reports can be submitted via the
  65 web:
  66
  67   https://redmine.open-bio.org/projects/bioperl/
  68
  69 =head1 AUTHOR - Jason Stajich
  70
  71 Email jason-at-bioperl.org
  72
  73 =head1 APPENDIX
  74
  75 The rest of the documentation details each of the object methods.
  76 Internal methods are usually preceded with a _
  77
  78 =cut
  79
  80 # Let the code begin...
  81
  82 package Bio::SearchIO::hmmer2;
  83
  84 use strict;
  85
  86 use Bio::Factory::ObjectFactory;
  87
  88 use vars qw(%MAPPING %MODEMAP
  89 );
  90
  91 use base qw(Bio::SearchIO::hmmer);
  92
  93 BEGIN {
  94
  95     # mapping of HMMER items to Bioperl hash keys
  96     %MODEMAP = (
  97         'HMMER_Output' => 'result',
  98         'Hit'          => 'hit',
  99         'Hsp'          => 'hsp'
 100     );
 101
 102     %MAPPING = (
 103         'Hsp_bit-score'   => 'HSP-bits',
 104         'Hsp_score'       => 'HSP-score',
 105         'Hsp_evalue'      => 'HSP-evalue',
 106         'Hsp_query-from'  => 'HSP-query_start',
 107         'Hsp_query-to'    => 'HSP-query_end',
 108         'Hsp_hit-from'    => 'HSP-hit_start',
 109         'Hsp_hit-to'      => 'HSP-hit_end',
 110         'Hsp_positive'    => 'HSP-conserved',
 111         'Hsp_identity'    => 'HSP-identical',
 112         'Hsp_gaps'        => 'HSP-hsp_gaps',
 113         'Hsp_hitgaps'     => 'HSP-hit_gaps',
 114         'Hsp_querygaps'   => 'HSP-query_gaps',
 115         'Hsp_qseq'        => 'HSP-query_seq',
 116         'Hsp_hseq'        => 'HSP-hit_seq',
 117         'Hsp_midline'     => 'HSP-homology_seq',
 118         'Hsp_align-len'   => 'HSP-hsp_length',
 119         'Hsp_query-frame' => 'HSP-query_frame',
 120         'Hsp_hit-frame'   => 'HSP-hit_frame',
 121
 122         'Hit_id'        => 'HIT-name',
 123         'Hit_len'       => 'HIT-length',
 124         'Hit_accession' => 'HIT-accession',
 125         'Hit_desc'      => 'HIT-description',
 126         'Hit_signif'    => 'HIT-significance',
 127         'Hit_score'     => 'HIT-score',
 128
 129         'HMMER_program'   => 'RESULT-algorithm_name',
 130         'HMMER_version'   => 'RESULT-algorithm_version',
 131         'HMMER_query-def' => 'RESULT-query_name',
 132         'HMMER_query-len' => 'RESULT-query_length',
 133         'HMMER_query-acc' => 'RESULT-query_accession',
 134         'HMMER_querydesc' => 'RESULT-query_description',
 135         'HMMER_hmm'       => 'RESULT-hmm_name',
 136         'HMMER_seqfile'   => 'RESULT-sequence_file',
 137         'HMMER_db'        => 'RESULT-database_name',
 138     );
 139 }
 140
 141 =head2 next_result
 142
 143  Title   : next_result
 144  Usage   : my $hit = $searchio->next_result;
 145  Function: Returns the next Result from a search
 146  Returns : Bio::Search::Result::ResultI object
 147  Args    : none
 148
 149 =cut
 150
 151 sub next_result {
 152     my ($self) = @_;
 153     my $seentop = 0;
 154     my $reporttype;
 155     my ( $last, @hitinfo, @hspinfo, %hspinfo, %hitinfo );
 156     local $/ = "\n";
 157     local $_;
 158
 159     my $verbose = $self->verbose;    # cache for speed?
 160     $self->start_document();
 161     local ($_);
 162     while ( defined( $_ = $self->_readline ) ) {
 163         my $lineorig = $_;
 164         chomp;
 165         if (/^HMMER\s+(\S+)\s+\((.+)\)/o) {
 166             my ( $prog, $version ) = split;
 167             if ($seentop) {
 168                 $self->_pushback($_);
 169                 $self->end_element( { 'Name' => 'HMMER_Output' } );
 170                 return $self->end_document();
 171             }
 172             $self->{'_hmmidline'} = $_;
 173             $self->start_element( { 'Name' => 'HMMER_Output' } );
 174             $self->{'_result_count'}++;
 175             $seentop = 1;
 176             if ( defined $last ) {
 177                 ($reporttype) = split( /\s+/, $last );
 178                 $reporttype = uc($reporttype) if defined $reporttype;
 179                 $self->element(
 180                     {
 181                         'Name' => 'HMMER_program',
 182                         'Data' => $reporttype
 183                     }
 184                 );
 185             }
 186             $self->element(
 187                 {
 188                     'Name' => 'HMMER_version',
 189                     'Data' => $version
 190                 }
 191             );
 192         }
 193         elsif (s/^HMM file:\s+//o) {
 194             $self->{'_hmmfileline'} = $lineorig;
 195             $self->element(
 196                 {
 197                     'Name' => 'HMMER_hmm',
 198                     'Data' => $_
 199                 }
 200             );
 201         }
 202         elsif (s/^Sequence\s+(file|database):\s+//o) {
 203             $self->{'_hmmseqline'} = $lineorig;
 204             if ( $1 eq 'database' ) {
 205                 $self->element(
 206                     {
 207                         'Name' => 'HMMER_db',
 208                         'Data' => $_
 209                     }
 210                 );
 211             }
 212             $self->element(
 213                 {
 214                     'Name' => 'HMMER_seqfile',
 215                     'Data' => $_
 216                 }
 217             );
 218         }
 219         elsif (s/^Query(\s+(sequence|HMM))?(?:\s+\d+)?:\s+//o) {
 220             if ( !$seentop ) {
 221
 222                 # we're in a multi-query report
 223                 $self->_pushback($lineorig);
 224                 $self->_pushback( $self->{'_hmmseqline'} );
 225                 $self->_pushback( $self->{'_hmmfileline'} );
 226                 $self->_pushback( $self->{'_hmmidline'} );
 227                 next;
 228             }
 229             s/\s+$//;
 230             $self->element(
 231                 {
 232                     'Name' => 'HMMER_query-def',
 233                     'Data' => $_
 234                 }
 235             );
 236         }
 237         elsif (s/^Accession:\s+//o) {
 238             s/\s+$//;
 239             $self->element(
 240                 {
 241                     'Name' => 'HMMER_query-acc',
 242                     'Data' => $_
 243                 }
 244             );
 245         }
 246         elsif (s/^Description:\s+//o) {
 247             s/\s+$//;
 248             $self->element(
 249                 {
 250                     'Name' => 'HMMER_querydesc',
 251                     'Data' => $_
 252                 }
 253             );
 254         }
 255         elsif ( defined $self->{'_reporttype'}
 256             && $self->{'_reporttype'} eq 'HMMSEARCH' )
 257         {
 258
 259             # PROCESS HMMSEARCH RESULTS HERE
 260             if (/^Scores for complete sequences/o) {
 261                 while ( defined( $_ = $self->_readline ) ) {
 262                     last if (/^\s+$/);
 263                     next if ( /^Sequence\s+Description/o || /^\-\-\-/o );
 264                     my @line = split;
 265                     my ( $name, $n, $evalue, $score ) =
 266                       ( shift @line, pop @line, pop @line, pop @line );
 267                     my $desc = join( ' ', @line );
 268                     push @hitinfo, [ $name, $desc, $evalue, $score ];
 269                     $hitinfo{$name} = $#hitinfo;
 270                 }
 271             }
 272             elsif (/^Parsed for domains:/o) {
 273                 @hspinfo = ();
 274
 275                 while ( defined( $_ = $self->_readline ) ) {
 276                     last if (/^\s+$/);
 277                     if (m!^//!) {
 278                         $self->_pushback($_);
 279                         last;
 280                     }
 281                     next if ( /^(Model|Sequence)\s+Domain/ || /^\-\-\-/ );
 282
 283                     chomp;
 284                     if (
 285                         my ( $n, $domainnum, $domainct, @vals ) = (
 286                             m!^(\S+)\s+      # host name
 287                         (\d+)/(\d+)\s+   # num/num (ie 1 of 2)
 288                         (\d+)\s+(\d+).+? # sequence start and end
 289                         (\d+)\s+(\d+)\s+ # hmm start and end
 290                         \S+\s+           # []
 291                         (\S+)\s+         # score
 292                         (\S+)            # evalue
 293                         \s*$!ox
 294                         )
 295                       )
 296                     {
 297
 298                         # array lookup so that we can get rid of things
 299                         # when they've been processed
 300                         my $info = $hitinfo[ $hitinfo{$n} ];
 301                         if ( !defined $info ) {
 302                             $self->warn(
 303 "Incomplete Sequence information, can't find $n hitinfo says $hitinfo{$n}"
 304                             );
 305                             next;
 306                         }
 307                         push @hspinfo, [ $n, @vals ];
 308                     }
 309                 }
 310             }
 311             elsif (/^Alignments of top/o) {
 312                 my ( $prelength, $lastdomain, $count, $width );
 313                 $count = 0;
 314                 my %domaincounter;
 315                 my $second_tier = 0;
 316                 while ( defined( $_ = $self->_readline ) ) {
 317                     next if ( /^Align/o
 318                         || /^\s+RF\s+[x\s]+$/o );
 319                     if ( /^Histogram/o || m!^//!o ) {
 320                         if ( $self->in_element('hsp') ) {
 321                             $self->end_element( { 'Name' => 'Hsp' } );
 322                         }
 323                         if ( $self->within_element('hit') ) {
 324                             $self->end_element( { 'Name' => 'Hit' } );
 325                         }
 326                         last;
 327                     }
 328                     chomp;
 329
 330                     if (
 331                         m/^\s*(.+):\s+domain\s+(\d+)\s+of\s+(\d+)\,\s+
 332                         from\s+(\d+)\s+to\s+(\d+)/x
 333                       )
 334                     {
 335                         my ( $name, $domainct, $domaintotal, $from, $to ) =
 336                           ( $1, $2, $3, $4, $5 );
 337                         $domaincounter{$name}++;
 338                         if ( $self->within_element('hit') ) {
 339                             if ( $self->within_element('hsp') ) {
 340                                 $self->end_element( { 'Name' => 'Hsp' } );
 341                             }
 342                             $self->end_element( { 'Name' => 'Hit' } );
 343                         }
 344
 345                         $self->start_element( { 'Name' => 'Hit' } );
 346                         my $info = [
 347                             @{
 348                                 $hitinfo[ $hitinfo{$name} ] || $self->throw(
 349 "Could not find hit info for $name: Insure that your database contains only unique sequence names"
 350                                 )
 351                               }
 352                         ];
 353                         if ( $info->[0] ne $name ) {
 354                             $self->throw(
 355 "Somehow the Model table order does not match the order in the domains (got "
 356                                   . $info->[0]
 357                                   . ", expected $name)" );
 358                         }
 359                         $self->element(
 360                             {
 361                                 'Name' => 'Hit_id',
 362                                 'Data' => shift @{$info}
 363                             }
 364                         );
 365                         $self->element(
 366                             {
 367                                 'Name' => 'Hit_desc',
 368                                 'Data' => shift @{$info}
 369                             }
 370                         );
 371                         $self->element(
 372                             {
 373                                 'Name' => 'Hit_signif',
 374                                 'Data' => shift @{$info}
 375                             }
 376                         );
 377                         $self->element(
 378                             {
 379                                 'Name' => 'Hit_score',
 380                                 'Data' => shift @{$info}
 381                             }
 382                         );
 383
 384                         $self->start_element( { 'Name' => 'Hsp' } );
 385                         $self->element(
 386                             {
 387                                 'Name' => 'Hsp_identity',
 388                                 'Data' => 0
 389                             }
 390                         );
 391                         $self->element(
 392                             {
 393                                 'Name' => 'Hsp_positive',
 394                                 'Data' => 0
 395                             }
 396                         );
 397                         my $HSPinfo = shift @hspinfo;
 398                         my $id      = shift @$HSPinfo;
 399
 400                         if ( $id ne $name ) {
 401                             $self->throw(
 402 "Somehow the domain list details do not match the table (got $id, expected $name)"
 403                             );
 404                         }
 405                         if ( $domaincounter{$name} == $domaintotal ) {
 406                             $hitinfo[ $hitinfo{$name} ] = undef;
 407                         }
 408                         $self->element(
 409                             {
 410                                 'Name' => 'Hsp_hit-from',
 411                                 'Data' => shift @$HSPinfo
 412                             }
 413                         );
 414                         $self->element(
 415                             {
 416                                 'Name' => 'Hsp_hit-to',
 417                                 'Data' => shift @$HSPinfo
 418                             }
 419                         );
 420                         $self->element(
 421                             {
 422                                 'Name' => 'Hsp_query-from',
 423                                 'Data' => shift @$HSPinfo
 424                             }
 425                         );
 426                         $self->element(
 427                             {
 428                                 'Name' => 'Hsp_query-to',
 429                                 'Data' => shift @$HSPinfo
 430                             }
 431                         );
 432                         $self->element(
 433                             {
 434                                 'Name' => 'Hsp_score',
 435                                 'Data' => shift @$HSPinfo
 436                             }
 437                         );
 438                         $self->element(
 439                             {
 440                                 'Name' => 'Hsp_evalue',
 441                                 'Data' => shift @$HSPinfo
 442                             }
 443                         );
 444                         $lastdomain = $name;
 445                     }
 446                     else {
 447
 448                         # Might want to change this so that it
 449                         # accumulates all the of the alignment lines into
 450                         # three array slots and then tests for the
 451                         # end of the line
 452                         if (/^(\s+\*\-\>)(\S+)/o) {    # start of domain
 453                             $prelength = CORE::length($1);
 454                             $width     = 0;
 455
 456                             # deal with fact that start en stop is on same line
 457                             my $data = $2;
 458                             if ($data =~ s/\<\-?\*?\s*$//)
 459                             {
 460                                 $width = CORE::length($data);
 461                             }
 462
 463                             $self->element(
 464                                 {
 465                                     'Name' => 'Hsp_qseq',
 466                                     'Data' => $data
 467                                 }
 468                             );
 469                             $count       = 0;
 470                             $second_tier = 0;
 471                         }
 472                         elsif (/^(\s+)(\S+)\<\-\*\s*$/o) {    #end of domain
 473                             $self->element(
 474                                 {
 475                                     'Name' => 'Hsp_qseq',
 476                                     'Data' => $2
 477                                 }
 478                             );
 479                             $width = CORE::length($2);
 480                             $count = 0;
 481                         }
 482                         elsif (( $count != 1 && /^\s+$/o )
 483                             || CORE::length($_) == 0
 484                             || /^\s+\-?\*\s*$/ )
 485                         {
 486                             next;
 487                         }
 488                         elsif ( $count == 0 ) {
 489                             $prelength -= 3 unless ( $second_tier++ );
 490                             unless ( defined $prelength ) {
 491
 492                                 # $self->warn("prelength not set");
 493                                 next;
 494                             }
 495                             $self->element(
 496                                 {
 497                                     'Name' => 'Hsp_qseq',
 498                                     'Data' => substr( $_, $prelength )
 499                                 }
 500                             );
 501                         }
 502                         elsif ( $count == 1 ) {
 503                             if ( !defined $prelength ) {
 504                                 $self->warn("prelength not set");
 505                             }
 506                             if ($width) {
 507                                 $self->element(
 508                                     {
 509                                         'Name' => 'Hsp_midline',
 510                                         'Data' =>
 511                                           substr( $_, $prelength, $width )
 512                                     }
 513                                 );
 514                             }
 515                             else {
 516                                 $self->element(
 517                                     {
 518                                         'Name' => 'Hsp_midline',
 519                                         'Data' => substr( $_, $prelength )
 520                                     }
 521                                 );
 522                             }
 523                         }
 524                         elsif ( $count == 2 ) {
 525                             if (/^\s+(\S+)\s+(\d+|\-)\s+(\S*)\s+(\d+|\-)/o) {
 526                                 $self->element(
 527                                     {
 528                                         'Name' => 'Hsp_hseq',
 529                                         'Data' => $3
 530                                     }
 531                                 );
 532                             }
 533                             else {
 534                                 $self->warn("unrecognized line: $_\n");
 535                             }
 536                         }
 537                         $count = 0 if $count++ >= 2;
 538                     }
 539                 }
 540             }
 541             elsif ( /^Histogram/o || m!^//!o ) {
 542                 while ( my $HSPinfo = shift @hspinfo ) {
 543                     my $id   = shift @$HSPinfo;
 544                     my $info = [ @{ $hitinfo[ $hitinfo{$id} ] } ];
 545                     next unless defined $info;
 546                     $self->start_element( { 'Name' => 'Hit' } );
 547                     $self->element(
 548                         {
 549                             'Name' => 'Hit_id',
 550                             'Data' => shift @{$info}
 551                         }
 552                     );
 553                     $self->element(
 554                         {
 555                             'Name' => 'Hit_desc',
 556                             'Data' => shift @{$info}
 557                         }
 558                     );
 559                     $self->element(
 560                         {
 561                             'Name' => 'Hit_signif',
 562                             'Data' => shift @{$info}
 563                         }
 564                     );
 565                     $self->element(
 566                         {
 567                             'Name' => 'Hit_score',
 568                             'Data' => shift @{$info}
 569                         }
 570                     );
 571                     $self->start_element( { 'Name' => 'Hsp' } );
 572                     $self->element(
 573                         {
 574                             'Name' => 'Hsp_query-from',
 575                             'Data' => shift @$HSPinfo
 576                         }
 577                     );
 578                     $self->element(
 579                         {
 580                             'Name' => 'Hsp_query-to',
 581                             'Data' => shift @$HSPinfo
 582                         }
 583                     );
 584                     $self->element(
 585                         {
 586                             'Name' => 'Hsp_hit-from',
 587                             'Data' => shift @$HSPinfo
 588                         }
 589                     );
 590                     $self->element(
 591                         {
 592                             'Name' => 'Hsp_hit-to',
 593                             'Data' => shift @$HSPinfo
 594                         }
 595                     );
 596                     $self->element(
 597                         {
 598                             'Name' => 'Hsp_score',
 599                             'Data' => shift @$HSPinfo
 600                         }
 601                     );
 602                     $self->element(
 603                         {
 604                             'Name' => 'Hsp_evalue',
 605                             'Data' => shift @$HSPinfo
 606                         }
 607                     );
 608                     $self->element(
 609                         {
 610                             'Name' => 'Hsp_identity',
 611                             'Data' => 0
 612                         }
 613                     );
 614                     $self->element(
 615                         {
 616                             'Name' => 'Hsp_positive',
 617                             'Data' => 0
 618                         }
 619                     );
 620                     $self->element(
 621                         {
 622                             'Name' => 'Hsp_positive',
 623                             'Data' => 0
 624                         }
 625                     );
 626                     $self->end_element( { 'Name' => 'Hsp' } );
 627                     $self->end_element( { 'Name' => 'Hit' } );
 628                 }
 629                 @hitinfo = ();
 630                 %hitinfo = ();
 631                 last;
 632             }
 633         }
 634         elsif ( defined $self->{'_reporttype'}
 635             && $self->{'_reporttype'} eq 'HMMPFAM' )
 636         {
 637             # process HMMPFAM results here
 638             if (/^Scores for sequence family/o) {
 639                 while ( defined( $_ = $self->_readline ) ) {
 640                     last if (/^\s+$/);
 641                     next if ( /^Model\s+Description/o || /^\-\-\-/o );
 642                     chomp;
 643                     my @line = split;
 644                     my ( $model, $n, $evalue, $score ) =
 645                       ( shift @line, pop @line, pop @line, pop @line );
 646                     my $desc = join( ' ', @line );
 647                     push @hitinfo, [ $model, $desc, $score, $evalue, $n ];
 648                     $hitinfo{$model} = $#hitinfo;
 649                 }
 650             }
 651             elsif (/^Parsed for domains:/o) {
 652                 @hspinfo = ();
 653                 while ( defined( $_ = $self->_readline ) ) {
 654                     last if (/^\s+$/);
 655                     if (m!^//!) {
 656                         $self->_pushback($_);
 657                         last;
 658                     }
 659                     next if ( /^Model\s+Domain/o || /^\-\-\-/o );
 660                     chomp;
 661                     if (
 662                         my ( $n, $domainnum, $domainct, @vals ) = (
 663                             m!^(\S+)\s+         # domain name
 664                             (\d+)/(\d+)\s+      # domain num out of num
 665                             (\d+)\s+(\d+).+?    # seq start, end
 666                             (\d+)\s+(\d+)\s+    # hmm start, end
 667                             \S+\s+              # []
 668                             (\S+)\s+            # score
 669                             (\S+)               # evalue
 670                             \s*$!ox
 671                         )
 672                       )
 673                     {
 674                         my $hindex = $hitinfo{$n};
 675                         if ( !defined $hindex ) {
 676                             push @hitinfo,
 677                               [ $n, '', $vals[5], $vals[6], $domainct ];
 678                             $hitinfo{$n} = $#hitinfo;
 679                             $hindex = $#hitinfo;
 680                         }
 681                         my $info = $hitinfo[$hindex];
 682                         if ( !defined $info ) {
 683                             $self->warn(
 684 "incomplete Domain information, can't find $n hitinfo says $hitinfo{$n}"
 685                             );
 686                             next;
 687                         }
 688                         push @hspinfo, [ $n, @vals ];
 689                     }
 690                 }
 691             }
 692             elsif (/^Alignments of top/o) {
 693                 my ( $prelength, $lastdomain, $count, $width );
 694                 $count = 0;
 695                 my $second_tier = 0;
 696                 while ( defined( $_ = $self->_readline ) ) {
 697                     next
 698                       if (
 699                         /^Align/o
 700                         || ( $count != 1
 701                             && /^\s+RF\s+[x\s]+$/o )
 702                       );
 703                     # fix for bug 2632
 704                     next if ($_ =~ m/^\s+CS\s+/o && $count == 0);
 705                     if ( /^Histogram/o || m!^//!o || /^Query sequence/o ) {
 706                         if ( $self->in_element('hsp') ) {
 707                             $self->end_element( { 'Name' => 'Hsp' } );
 708                         }
 709                         if ( $self->in_element('hit') ) {
 710                             $self->end_element( { 'Name' => 'Hit' } );
 711                         }
 712                         $self->_pushback($_);
 713                         last;
 714                     }
 715                     chomp;
 716                     if (m/(\S+):.*from\s+(\d+)\s+to\s+(\d+)/o) {
 717                         my ( $name, $from, $to ) = ( $1, $2, $3 );
 718
 719                         if ( $self->within_element('hit') ) {
 720                             if ( $self->in_element('hsp') ) {
 721                                 $self->end_element( { 'Name' => 'Hsp' } );
 722                             }
 723                             $self->end_element( { 'Name' => 'Hit' } );
 724                         }
 725                         my $info = [ @{ $hitinfo[ $hitinfo{$name} ] } ];
 726                         if ( !defined $info
 727                             || $info->[0] ne $name )
 728                         {
 729                             $self->warn(
 730 "Somehow the Model table order does not match the order in the domains (got "
 731                                   . $info->[0]
 732                                   . ", expected $name). We're back loading this from the alignment information instead"
 733                             );
 734                             $info = [
 735                                 $name, '',
 736                                 /score\s+([^,\s]+),\s+E\s+=\s+(\S+)/ox
 737                             ];
 738                             push @hitinfo, $info;
 739                             $hitinfo{$name} = $#hitinfo;
 740                         }
 741                         $self->start_element( { 'Name' => 'Hit' } );
 742
 743                         $self->element(
 744                             {
 745                                 'Name' => 'Hit_id',
 746                                 'Data' => shift @{$info}
 747                             }
 748                         );
 749                         $self->element(
 750                             {
 751                                 'Name' => 'Hit_desc',
 752                                 'Data' => shift @{$info}
 753                             }
 754                         );
 755                         $self->element(
 756                             {
 757                                 'Name' => 'Hit_score',
 758                                 'Data' => shift @{$info}
 759                             }
 760                         );
 761                         $self->element(
 762                             {
 763                                 'Name' => 'Hit_signif',
 764                                 'Data' => shift @{$info}
 765                             }
 766                         );
 767
 768                         $self->start_element( { 'Name' => 'Hsp' } );
 769                         $self->element(
 770                             {
 771                                 'Name' => 'Hsp_identity',
 772                                 'Data' => 0
 773                             }
 774                         );
 775                         $self->element(
 776                             {
 777                                 'Name' => 'Hsp_positive',
 778                                 'Data' => 0
 779                             }
 780                         );
 781                         my $HSPinfo = shift @hspinfo;
 782                         my $id      = shift @$HSPinfo;
 783
 784                         if ( $id ne $name ) {
 785                             $self->throw(
 786 "Somehow the domain list details do not match the table (got $id, expected $name)"
 787                             );
 788                         }
 789                         $self->element(
 790                             {
 791                                 'Name' => 'Hsp_query-from',
 792                                 'Data' => shift @$HSPinfo
 793                             }
 794                         );
 795                         $self->element(
 796                             {
 797                                 'Name' => 'Hsp_query-to',
 798                                 'Data' => shift @$HSPinfo
 799                             }
 800                         );
 801                         $self->element(
 802                             {
 803                                 'Name' => 'Hsp_hit-from',
 804                                 'Data' => shift @$HSPinfo
 805                             }
 806                         );
 807                         $self->element(
 808                             {
 809                                 'Name' => 'Hsp_hit-to',
 810                                 'Data' => shift @$HSPinfo
 811                             }
 812                         );
 813                         $self->element(
 814                             {
 815                                 'Name' => 'Hsp_score',
 816                                 'Data' => shift @$HSPinfo
 817                             }
 818                         );
 819                         $self->element(
 820                             {
 821                                 'Name' => 'Hsp_evalue',
 822                                 'Data' => shift @$HSPinfo
 823                             }
 824                         );
 825                         $lastdomain = $name;
 826                     }
 827                     else {
 828                         if (/^(\s+\*\-\>)(\S+)/o) {
 829
 830                             # start of domain
 831                             $prelength = CORE::length($1);
 832                             $width     = 0;
 833
 834                             # deal with fact that start en stop is on same line
 835                             my $data = $2;
 836                             if ($data =~ s/\<\-?\*?\s*$//)
 837                             {
 838                               $width = CORE::length($data);
 839                             }
 840
 841                             $self->element(
 842                                 {
 843                                     'Name' => 'Hsp_hseq',
 844                                     'Data' => $data
 845                                 }
 846                             );
 847                             $count       = 0;
 848                             $second_tier = 0;
 849
 850                         }
 851                         elsif (/^(\s+)(\S+)\<\-?\*?\s*$/o) {
 852
 853                             #end of domain
 854                             $prelength -= 3 unless ( $second_tier++ );
 855                             $self->element(
 856                                 {
 857                                     'Name' => 'Hsp_hseq',
 858                                     'Data' => $2
 859                                 }
 860                             );
 861                             $width = CORE::length($2);
 862                             $count = 0;
 863                         }
 864                         elsif (CORE::length($_) == 0
 865                             || ( $count != 1 && /^\s+$/o )
 866                             || /^\s+\-?\*\s*$/
 867                             || /^\s+\S+\s+\-\s+\-\s*$/ )
 868                         {
 869                             next;
 870                         }
 871                         elsif ( $count == 0 ) {
 872                             $prelength -= 3 unless ( $second_tier++ );
 873                             unless ( defined $prelength ) {
 874
 875                                 # $self->warn("prelength not set");
 876                                 next;
 877                             }
 878                             $self->element(
 879                                 {
 880                                     'Name' => 'Hsp_hseq',
 881                                     'Data' => substr( $_, $prelength )
 882                                 }
 883                             );
 884                         }
 885                         elsif ( $count == 1 ) {
 886                             if ( !defined $prelength ) {
 887                                 $self->warn("prelength not set");
 888                             }
 889                             if ($width) {
 890                                 $self->element(
 891                                     {
 892                                         'Name' => 'Hsp_midline',
 893                                         'Data' =>
 894                                           substr( $_, $prelength, $width )
 895                                     }
 896                                 );
 897                             }
 898                             else {
 899                                 $self->element(
 900                                     {
 901                                         'Name' => 'Hsp_midline',
 902                                         'Data' => substr( $_, $prelength )
 903                                     }
 904                                 );
 905                             }
 906                         }
 907                         elsif ( $count == 2 ) {
 908                             if (   /^\s+(\S+)\s+(\d+)\s+(\S+)\s+(\d+)/o
 909                                 || /^\s+(\S+)\s+(\-)\s+(\S*)\s+(\-)/o )
 910                             {
 911                                 $self->element(
 912                                     {
 913                                         'Name' => 'Hsp_qseq',
 914                                         'Data' => $3
 915                                     }
 916                                 );
 917                             }
 918                             else {
 919                                 $self->throw(
 920                                     "unrecognized line ($count): $_\n");
 921                             }
 922                         }
 923                         $count = 0 if $count++ >= 2;
 924                     }
 925                 }
 926             }
 927             elsif ( /^Histogram/o || m!^//!o ) {
 928
 929                 while ( my $HSPinfo = shift @hspinfo ) {
 930                     my $id   = shift @$HSPinfo;
 931                     my $info = [ @{ $hitinfo[ $hitinfo{$id} ] } ];
 932                     next unless defined $info;
 933                     $self->start_element( { 'Name' => 'Hit' } );
 934                     $self->element(
 935                         {
 936                             'Name' => 'Hit_id',
 937                             'Data' => shift @{$info}
 938                         }
 939                     );
 940                     $self->element(
 941                         {
 942                             'Name' => 'Hit_desc',
 943                             'Data' => shift @{$info}
 944                         }
 945                     );
 946                     $self->element(
 947                         {
 948                             'Name' => 'Hit_signif',
 949                             'Data' => shift @{$info}
 950                         }
 951                     );
 952                     $self->element(
 953                         {
 954                             'Name' => 'Hit_score',
 955                             'Data' => shift @{$info}
 956                         }
 957                     );
 958                     $self->start_element( { 'Name' => 'Hsp' } );
 959                     $self->element(
 960                         {
 961                             'Name' => 'Hsp_query-from',
 962                             'Data' => shift @$HSPinfo
 963                         }
 964                     );
 965                     $self->element(
 966                         {
 967                             'Name' => 'Hsp_query-to',
 968                             'Data' => shift @$HSPinfo
 969                         }
 970                     );
 971                     $self->element(
 972                         {
 973                             'Name' => 'Hsp_hit-from',
 974                             'Data' => shift @$HSPinfo
 975                         }
 976                     );
 977                     $self->element(
 978                         {
 979                             'Name' => 'Hsp_hit-to',
 980                             'Data' => shift @$HSPinfo
 981                         }
 982                     );
 983                     $self->element(
 984                         {
 985                             'Name' => 'Hsp_score',
 986                             'Data' => shift @$HSPinfo
 987                         }
 988                     );
 989                     $self->element(
 990                         {
 991                             'Name' => 'Hsp_evalue',
 992                             'Data' => shift @$HSPinfo
 993                         }
 994                     );
 995                     $self->element(
 996                         {
 997                             'Name' => 'Hsp_identity',
 998                             'Data' => 0
 999                         }
1000                     );
1001                     $self->element(
1002                         {
1003                             'Name' => 'Hsp_positive',
1004                             'Data' => 0
1005                         }
1006                     );
1007                     $self->element(
1008                         {
1009                             'Name' => 'Hsp_positive',
1010                             'Data' => 0
1011                         }
1012                     );
1013                     $self->end_element( { 'Name' => 'Hsp' } );
1014                     $self->end_element( { 'Name' => 'Hit' } );
1015                 }
1016                 @hitinfo = ();
1017                 %hitinfo = ();
1018                 last;
1019             }
1020             # uncomment to see missed lines with verbose on
1021             #else {
1022             #    $self->debug($_);
1023             #}
1024         }
1025         $last = $_;
1026     }
1027     $self->end_element( { 'Name' => 'HMMER_Output' } ) unless !$seentop;
1028     return $self->end_document();
1029 }
1030
1031 =head2 start_element
1032
1033  Title   : start_element
1034  Usage   : $eventgenerator->start_element
1035  Function: Handles a start element event
1036  Returns : none
1037  Args    : hashref with at least 2 keys 'Data' and 'Name'
1038
1039
1040 =cut
1041
1042 sub start_element {
1043     my ( $self, $data ) = @_;
1044
1045     # we currently don't care about attributes
1046     my $nm   = $data->{'Name'};
1047     my $type = $MODEMAP{$nm};
1048     if ($type) {
1049         if ( $self->_eventHandler->will_handle($type) ) {
1050             my $func = sprintf( "start_%s", lc $type );
1051             $self->_eventHandler->$func( $data->{'Attributes'} );
1052         }
1053         unshift @{ $self->{'_elements'} }, $type;
1054     }
1055     if ( defined $type
1056         && $type eq 'result' )
1057     {
1058         $self->{'_values'} = {};
1059         $self->{'_result'} = undef;
1060     }
1061 }
1062
1063 =head2 end_element
1064
1065  Title   : start_element
1066  Usage   : $eventgenerator->end_element
1067  Function: Handles an end element event
1068  Returns : none
1069  Args    : hashref with at least 2 keys 'Data' and 'Name'
1070
1071
1072 =cut
1073
1074 sub end_element {
1075     my ( $self, $data ) = @_;
1076     my $nm   = $data->{'Name'};
1077     my $type = $MODEMAP{$nm};
1078     my $rc;
1079
1080     if ( $nm eq 'HMMER_program' ) {
1081         if ( $self->{'_last_data'} =~ /(HMM\S+)/i ) {
1082             $self->{'_reporttype'} = uc $1;
1083         }
1084     }
1085
1086     # Hsp are sort of weird, in that they end when another
1087     # object begins so have to detect this in end_element for now
1088     if ( $nm eq 'Hsp' ) {
1089         foreach (qw(Hsp_qseq Hsp_midline Hsp_hseq)) {
1090             my $data = $self->{'_last_hspdata'}->{$_};
1091             if ($data && $_ eq 'Hsp_hseq') {
1092                 # replace hmm '.' gap symbol by '-'
1093                 $data =~ s/\./-/g;
1094             }
1095             $self->element(
1096                 {
1097                     'Name' => $_,
1098                     'Data' => $data
1099                 }
1100             );
1101         }
1102         $self->{'_last_hspdata'} = {};
1103     }
1104     if ($type) {
1105         if ( $self->_eventHandler->will_handle($type) ) {
1106             my $func = sprintf( "end_%s", lc $type );
1107             $rc = $self->_eventHandler->$func( $self->{'_reporttype'},
1108                 $self->{'_values'} );
1109         }
1110         my $lastelem = shift @{ $self->{'_elements'} };
1111     }
1112     elsif ( $MAPPING{$nm} ) {
1113         if ( ref( $MAPPING{$nm} ) =~ /hash/i ) {
1114             my $key = ( keys %{ $MAPPING{$nm} } )[0];
1115             $self->{'_values'}->{$key}->{ $MAPPING{$nm}->{$key} } =
1116               $self->{'_last_data'};
1117         }
1118         else {
1119             $self->{'_values'}->{ $MAPPING{$nm} } = $self->{'_last_data'};
1120         }
1121     }
1122     else {
1123         $self->debug("unknown nm $nm, ignoring\n");
1124     }
1125     $self->{'_last_data'} = '';    # remove read data if we are at
1126                                    # end of an element
1127     $self->{'_result'} = $rc if ( defined $type && $type eq 'result' );
1128     return $rc;
1129 }
1130
1131 =head2 element
1132
1133  Title   : element
1134  Usage   : $eventhandler->element({'Name' => $name, 'Data' => $str});
1135  Function: Convience method that calls start_element, characters, end_element
1136  Returns : none
1137  Args    : Hash ref with the keys 'Name' and 'Data'
1138
1139
1140 =cut
1141
1142 sub element {
1143     my ( $self, $data ) = @_;
1144     $self->start_element($data);
1145     $self->characters($data);
1146     $self->end_element($data);
1147 }
1148
1149 =head2 characters
1150
1151  Title   : characters
1152  Usage   : $eventgenerator->characters($str)
1153  Function: Send a character events
1154  Returns : none
1155  Args    : string
1156
1157
1158 =cut
1159
1160 sub characters {
1161     my ( $self, $data ) = @_;
1162
1163     if (   $self->in_element('hsp')
1164         && $data->{'Name'} =~ /Hsp\_(qseq|hseq|midline)/o
1165         && defined $data->{'Data'} )
1166     {
1167         $self->{'_last_hspdata'}->{ $data->{'Name'} } .= $data->{'Data'};
1168     }
1169     return unless ( defined $data->{'Data'} && $data->{'Data'} !~ /^\s+$/o );
1170
1171     $self->{'_last_data'} = $data->{'Data'};
1172 }
1173
1174 =head2 within_element
1175
1176  Title   : within_element
1177  Usage   : if( $eventgenerator->within_element($element) ) {}
1178  Function: Test if we are within a particular element
1179            This is different than 'in' because within can be tested
1180            for a whole block.
1181  Returns : boolean
1182  Args    : string element name
1183
1184
1185 =cut
1186
1187 sub within_element {
1188     my ( $self, $name ) = @_;
1189     return 0
1190       if ( !defined $name
1191         || !defined $self->{'_elements'}
1192         || scalar @{ $self->{'_elements'} } == 0 );
1193     foreach ( @{ $self->{'_elements'} } ) {
1194         return 1 if ( $_ eq $name );
1195     }
1196     return 0;
1197 }
1198
1199 =head2 in_element
1200
1201  Title   : in_element
1202  Usage   : if( $eventgenerator->in_element($element) ) {}
1203  Function: Test if we are in a particular element
1204            This is different than 'within' because 'in' only
1205            tests its immediete parent.
1206  Returns : boolean
1207  Args    : string element name
1208
1209
1210 =cut
1211
1212 sub in_element {
1213     my ( $self, $name ) = @_;
1214     return 0 if !defined $self->{'_elements'}->[0];
1215     return ( $self->{'_elements'}->[0] eq $name );
1216 }
1217
1218 =head2 start_document
1219
1220  Title   : start_document
1221  Usage   : $eventgenerator->start_document
1222  Function: Handle a start document event
1223  Returns : none
1224  Args    : none
1225
1226
1227 =cut
1228
1229 sub start_document {
1230     my ($self) = @_;
1231     $self->{'_lasttype'} = '';
1232     $self->{'_values'}   = {};
1233     $self->{'_result'}   = undef;
1234     $self->{'_elements'} = [];
1235 }
1236
1237 =head2 end_document
1238
1239  Title   : end_document
1240  Usage   : $eventgenerator->end_document
1241  Function: Handles an end document event
1242  Returns : Bio::Search::Result::ResultI object
1243  Args    : none
1244
1245
1246 =cut
1247
1248 sub end_document {
1249     my ($self) = @_;
1250     return $self->{'_result'};
1251 }
1252
1253 =head2 result_count
1254
1255  Title   : result_count
1256  Usage   : my $count = $searchio->result_count
1257  Function: Returns the number of results we have processed
1258  Returns : integer
1259  Args    : none
1260
1261
1262 =cut
1263
1264 sub result_count {
1265     my $self = shift;
1266     return $self->{'_result_count'};
1267 }
1268
1269 1;