# This is an interface and cannot be instantiated

    # generally we use Bio::SearchIO to build these objects
    use Bio::SearchIO;
    my $in = new Bio::SearchIO(-format => 'hmmer_pull',
							   -file   => 'result.hmmer');

    while (my $result = $in->next_result) {
		while (my $hit = $result->next_hit) {
			while (my $hsp = $hit->next_hsp) {
                $r_type = $hsp->algorithm;
                $pvalue = $hsp->p();
                $evalue = $hsp->evalue();
                $frac_id = $hsp->frac_identical( ['query'|'hit'|'total'] );
                $frac_cons = $hsp->frac_conserved( ['query'|'hit'|'total'] );
                $gaps = $hsp->gaps( ['query'|'hit'|'total'] );
                $qseq = $hsp->query_string;
                $hseq = $hsp->hit_string;
                $homo_string = $hsp->homology_string;
                $len = $hsp->length( ['query'|'hit'|'total'] );
                $len = $hsp->length( ['query'|'hit'|'total'] );
                $rank = $hsp->rank;
            }

    my ($self, @args) = @_;
	
	# fields most subclasses probably will want
	$self->_fields( { ( hsp_length => undef,
                        identical => undef,
                        percent_identity => undef,
                        conserved => undef,
                        hsp_gaps => undef,
                        query_gaps => undef,
                        hit_gaps => undef,
						evalue => undef,
						pvalue => undef,
						score => undef,
						query_start => undef,
						query_end => undef,
						query_string => undef,
						hit_start => undef,
						hit_end => undef,
						hit_string => undef,
						homology_string => undef,
						rank => undef,
                        seq_inds => undef,
                        hit_identical_inds => undef,
                        hit_conserved_inds => undef,
                        hit_nomatch_inds => undef,
                        hit_gap_inds => undef,
                        query_identical_inds => undef,
                        query_conserved_inds => undef,
                        query_nomatch_inds => undef,
                        query_gap_inds => undef ) } );
	
	my ($parent, $chunk, $hsp_data) = $self->_rearrange([qw(PARENT
														    CHUNK
															HSP_DATA)], @args);
	
    $self->throw("Need -parent or -chunk to be defined") unless defined $parent || $chunk;
    
	$self->parent($parent) if $parent;
    
    if ($chunk) {
        my ($io, $start, $end) = (undef, 0, undef);
        if (ref($chunk) eq 'ARRAY') {
            ($io, $start, $end) = @{$chunk};
        }
        else {
            $io = $chunk;
        }
        $self->chunk($io, -start => $start, -end => $end);
    }
    
	$self->_raw_hsp_data($hsp_data) if $hsp_data;
	
    return $self;

	my $self = shift;
	if (@_) {
		$self->{_raw_hsp_data} = shift;
	}
	return $self->{_raw_hsp_data};

	return shift->get_field('algorithm');

	return shift->get_field('pvalue');

	return shift->get_field('evalue');

	my ($self, $type) = @_;
	
	$type = lc $type if defined $type;
	$type = 'hit' if (defined $type && $type =~ /subject|sbjct/);
	$type = 'total' if (! defined $type || $type eq 'hsp' || $type !~ /query|hit|subject|sbjct|total/);
	
	my $ratio = $self->num_identical($type) / $self->length($type);
    return sprintf( "%.3f", $ratio);

    my ($self, $type) = @_;
	
	$type = lc $type if defined $type;
	$type = 'hit' if (defined $type && $type =~ /subject|sbjct/);
	$type = 'total' if (! defined $type || $type eq 'hsp' || $type !~ /query|hit|subject|sbjct|total/);
	
	my $ratio = $self->num_conserved($type) / $self->length($type);
    return sprintf( "%.3f", $ratio);

    my $self = shift;
	return scalar($self->seq_inds('hit', 'identical'));

    my $self = shift;
	return scalar($self->seq_inds('hit', 'conserved-not-identical'));

    my ($self, $type) = @_;
    $type = lc $type if defined $type;
    $type = 'total' if (! defined $type || $type eq 'hsp' || $type !~ /query|hit|subject|sbjct|total/); 
    $type = 'hit' if $type =~ /sbjct|subject/;
	
	if ($type eq 'total') {
		return scalar($self->seq_inds('hit', 'gap')) + scalar($self->seq_inds('query', 'gap'));
	}
	return scalar($self->seq_inds($type, 'gap'));

	return shift->get_field('query_string');

	return shift->get_field('hit_string');

	return shift->get_field('homology_string');

    my ($self, $type) = @_;
    $type = 'total' unless defined $type;
    $type = lc $type;

    if ($type =~ /^q/i) {
        return $self->query->length;
    }
	elsif ($type =~ /^(hit|subject|sbjct)/) {
        return $self->hit->length;
    }
	else { 
        return $self->hit->length + $self->gaps('hit');
	}

	return shift->length('total');

	my ($self) = @_;
	return $self->frac_identical('hsp') * 100;

	my $self = shift;
	
    require Bio::LocatableSeq;
    require Bio::SimpleAlign;
    my $aln = new Bio::SimpleAlign;
    my $hs = $self->seq('hit');
    my $qs = $self->seq('query');
	if ($hs && $qs) {
		$aln->add_seq($hs);
		$aln->add_seq($qs);
		return $aln;
	}
	return;

    my ($self, $seqType, $class, $collapse) = @_;
    
    $seqType ||= 'query';
    $class ||= 'identical';
    $collapse ||= 0;
    $seqType = lc($seqType);
    $class = lc($class);
    $seqType = 'hit' if $seqType eq 'sbjct';
    
    my $t = substr($seqType,0,1);
    if ($t eq 'q') {
        $seqType = 'query';
    }
    elsif ($t eq 's' || $t eq 'h') {
        $seqType = 'hit';
    }
    else { 
        $self->warn("unknown seqtype $seqType using 'query'");
        $seqType = 'query';
    }
    
    $t = substr($class,0,1);
    if ($t eq 'c') {
        if ($class eq 'conserved-not-identical') {
            $class = 'conserved';
        }
        else { 
            $class = 'conservedall';
        }
    }
    elsif ($t eq 'i') {
        $class = 'identical';
    }
    elsif ($t eq 'n') {
        $class = 'nomatch';
    }
    elsif ($t eq 'g') {
        $class = 'gap';
    }
    else { 
        $self->warn("unknown sequence class $class using 'identical'");
        $class = 'identical';
    }
    
    $seqType .= '_';
    $class .= '_inds';
    
    my @ary;
    if ($class eq 'conservedall_inds') {
		my %tmp = map { $_, 1 } @{$self->get_field($seqType.'conserved_inds')},
								@{$self->get_field($seqType.'identical_inds')};
		@ary = sort {$a <=> $b} keys %tmp;
    }
    else { 
        @ary = @{$self->get_field($seqType.$class)};
    }
    
    return $collapse ? &Bio::Search::SearchUtils::collapse_nums(@ary) : @ary;

	return shift->get_field('evalue');

	return shift->get_field('score');

	return shift->get_field('bits');

    my $self = shift;
    my $val = shift;
    $val = 'query' unless defined $val;
    $val =~ s/^\s+//;

    if ($val =~ /^q/i) {
        return $self->query->strand(@_);
    }
    elsif ($val =~ /^hi|^s/i) {
        return $self->hit->strand(@_);
    }
    elsif ($val =~ /^list|array/i) {
        return ($self->query->strand(@_), $self->hit->strand(@_) );
    }
    else { 
        $self->warn("unrecognized component '$val' requested\n");
    }
    return 0;

    my $self = shift;
    my $val = shift;
    $val = 'query' unless defined $val;
    $val =~ s/^\s+//;

    if ($val =~ /^q/i) { 
        return $self->query->start(@_);
    }
    elsif ($val =~ /^(hi|s)/i) {
        return $self->hit->start(@_);
    }
    elsif ($val =~ /^list|array/i) {
        return ($self->query->start(@_), $self->hit->start(@_) );
    }
    else { 
        $self->warn("unrecognized component '$val' requested\n");
    }
    return 0;

    my $self = shift;
    my $val = shift;
    $val = 'query' unless defined $val;
    $val =~ s/^\s+//;

    if ($val =~ /^q/i) { 
        return $self->query->end(@_);
    }
    elsif ($val =~ /^(hi|s)/i) {
        return $self->hit->end(@_);
    }
    elsif ($val =~ /^list|array/i) {
        return ($self->query->end(@_), $self->hit->end(@_) );
    }
    else {
        $self->warn("unrecognized end component '$val' requested\n");
    }
    return 0;

    my ($self, $seqType) = @_; 
    $seqType ||= 'query';
    $seqType = 'hit' if $seqType eq 'sbjct';
    if ($seqType =~ /^(m|ho)/i ) {
        $self->throw("cannot call seq on the homology match string, it isn't really a sequence, use get_aln to convert the HSP to a Bio::AlignIO and generate a consensus from that.");
    }
	
    my $str = $self->seq_str($seqType) || return;
    require Bio::LocatableSeq;
    my $id = $seqType =~ /^q/i ? $self->query->seq_id : $self->hit->seq_id;
    return new Bio::LocatableSeq (  -ID    => $id,
                                    -SEQ   => $str,
                                    -START => $self->start($seqType),
                                    -END   => $self->end($seqType),
                                    -STRAND=> $self->strand($seqType),
                                    -DESC  => "$seqType sequence " );

      my $self = shift;
    my $type = shift || 'query';

    if ($type =~ /^q/i) {
        return $self->query_string(@_);
    }
    elsif ($type =~ /^(s)|(hi)/i) {
        return $self->hit_string(@_);
    }
    elsif ($type =~ /^(ho)|(ma)/i) {
        return $self->homology_string(@_);
    }
    else { 
        $self->warn("unknown sequence type $type");
    }
    return '';

    return shift->get_field('rank');

    my ($self, @args) = @_;
    my($seqType, $beg, $end) = $self->_rearrange([qw(SEQ START STOP)], @args);
    $seqType ||= 'query';
    $seqType = 'hit' if $seqType eq 'sbjct';
	
    my @data;
    if ((!defined $beg && !defined $end) || ! $self->seq_str('match')) {
        push @data, ($self->num_identical, $self->num_conserved);
    }
	else {
        $beg ||= 0;
        $end ||= 0;
        my ($start, $stop) = $self->range($seqType);
		
        if ($beg == 0) {
			$beg = $start;
			$end = $beg+$end;
		}
        elsif ($end == 0) {
			$end = $stop;
			$beg = $end-$beg;
		}
		
        if ($end >= $stop) {
			$end = $stop;
		}
        else {
			$end += 1;
		}
        if ($beg < $start) {
			$beg = $start;
		}
        
        my $seq = substr($self->seq_str('homology'), $beg-$start, ($end-$beg));
        
        if (!CORE::length $seq) {
            $self->throw("Undefined sub-sequence ($beg,$end). Valid range = $start - $stop");
        }
        ## Get data for a substring.
        $seq =~ s/ //g;  # remove space (no info).
        my $len_cons = CORE::length $seq;
        $seq =~ s/\+//g;  # remove '+' characters (conservative substitutions)
        my $len_id = CORE::length $seq;
        push @data, ($len_id, $len_cons);
    }
	
    return @data;

    return shift->get_field('num_hsps');

    my ($self, $seqType) = @_;
    $seqType ||= 'query';
	
    my ($start, $end);
    if ($seqType eq 'query') {
        $start = $self->query->start;
        $end = $self->query->end;
    }
    else {
        $start = $self->hit->start;
        $end = $self->hit->end;
    }
    return ($start, $end);

  bioperl-l@bioperl.org                  - General discussion
  http://bioperl.org/wiki/Mailing_lists  - About the mailing lists

  http://bugzilla.open-bio.org/

 Title   : query
 Usage   : my $query = $hsp->query
 Function: Returns a SeqFeature representing the query in the HSP
 Returns : Bio::SeqFeature::Similarity
 Args    : [optional] new value to set

 Title   : hit
 Usage   : my $hit = $hsp->hit
 Function: Returns a SeqFeature representing the hit in the HSP
 Returns : Bio::SeqFeature::Similarity
 Args    : [optional] new value to set