view release on metacpan or  search on metacpan

Ace/Sequence/Feature.pm  view on Meta::CPAN

use Ace qw(:DEFAULT rearrange);
use Ace::Object;
use Ace::Sequence::Homol;
use Carp;
use AutoLoader 'AUTOLOAD';
use vars '@ISA','%REV';
@ISA = 'Ace::Sequence';  # for convenience sake only
%REV = ('+1' => '-1',
	'-1' => '+1');  # war is peace, &c.

use overload 
  '""' => 'asString',
  ;

# parse a line from a sequence list
sub new {
  my $pack = shift;
  my ($parent,$ref,$r_offset,$r_strand,$abs,$gff_line,$db) = @_;
  my ($sourceseq,$method,$type,$start,$end,$score,$strand,$frame,$group) = split "\t",$gff_line;
  if (defined($strand)) {
    $strand = $strand eq '-' ? '-1' : '+1';
  } else {
    $strand = 0;
  }

  # for efficiency/performance, we don't use superclass new() method, but modify directly
  # handling coordinates.  See SCRAPS below for what should be in here
  $strand = '+1' if $strand < 0 && $r_strand < 0;  # two wrongs do make a right
  ($start,$end) = ($end,$start) if $strand < 0;
  my $offset = $start - 1;
  my $length = ($end > $start) ? $end - $offset : $end - $offset - 2;

  # handle negative strands
  $offset ||= 0;
  $offset *= -1 if $r_strand < 0 && $strand != $r_strand;

  my $self= bless {
		   obj      => $ref,
		   offset   => $offset,
		   length   => $length,
		   parent   => $parent,
		   p_offset => $r_offset,
		   refseq   => [$ref,$r_offset,$r_strand],
		   strand   => $r_strand,
		   fstrand  => $strand,
		   absolute => $abs,
		   info     => {
				seqname=> $sourceseq,
				method => $method,
				type   => $type,
				score  => $score,
				frame  => $frame,
				group  => $group,
				db     => $db,
			       }
		  },$pack;
  return $self;
}

sub smapped { 1; }

# $_[0] is field name, $_[1] is self, $_[2] is optional replacement value
sub _field {
  my $self = shift;
  my $field = shift;
  my $v = $self->{info}{$field};
  $self->{info}{$field} = shift if @_;
  return if defined $v && $v eq '.';
  return $v;
}

sub strand { return $_[0]->{fstrand} }

sub seqname   { 
  my $self = shift;
  my $seq = $self->_field('seqname');
  $self->db->fetch(Sequence=>$seq); 
}

sub method    { shift->_field('method',@_) }  # ... I prefer "method"
sub subtype   { shift->_field('method',@_) }  # ... or even "subtype"
sub type      { shift->_field('type',@_)   }  # ... I prefer "type"
sub score     { shift->_field('score',@_)  }  # float indicating some sort of score
sub frame     { shift->_field('frame',@_)  }  # one of 1, 2, 3 or undef
sub info      {                  # returns Ace::Object(s) with info about the feature
  my $self = shift;
  unless ($self->{group}) {
    my $info = $self->{info}{group} || 'Method "'.$self->method.'"';
    $info =~ s/(\"[^\"]*);([^\"]*\")/$1$;$2/g;
    my @data = split(/\s*;\s*/,$info);
    foreach (@data) { s/$;/;/g }
    $self->{group} = [map {$self->toAce($_)} @data];
  }
  return wantarray ? @{$self->{group}} : $self->{group}->[0];
}

# bioperl compatibility
sub primary_tag { shift->type(@_)    }
sub source_tag  { shift->subtype(@_) }

sub db { # database identifier (from Ace::Sequence::Multi)
  my $self = shift;
  my $db = $self->_field('db',@_);
  return $db || $self->SUPER::db;
}

sub group  { $_[0]->info; }
sub target { $_[0]->info; }

sub asString {
  my $self = shift;
  my $name = $self->SUPER::asString;
  my $type = $self->type;
  return "$type:$name";
}

# unique ID
sub id {
  my $self = shift;
  my $source = $self->source->name;
  my $start = $self->start;
  my $end = $self->end;
  return "$source/$start,$end";
}

# map info into a reasonable set of ace objects
sub toAce {
    my $self = shift;
    my $thing = shift;
    my ($tag,@values) = $thing=~/(\"[^\"]+?\"|\S+)/g;
    foreach (@values) { # strip the damn quotes
      s/^\"(.*)\"$/$1/;  # get rid of leading and trailing quotes
    }
    return $self->tag2ace($tag,@values);
}

# synthesize an artificial Ace object based on the tag
sub tag2ace {
    my $self = shift;
    my ($tag,@data) = @_;

    # Special cases, hardcoded in Ace GFF code...
    my $db = $self->db;;
    my $class = $db->class;

    # for Notes we just return a text, no database associated
    return $class->new(Text=>$data[0]) if $tag eq 'Note';

    # for homols, we create the indicated Protein or Sequence object
    # then generate a bogus Homology object (for future compatability??)
    if ($tag eq 'Target') {
	my ($objname,$start,$end) = @data;
	my ($classe,$name) = $objname =~ /^(\w+):(.+)/;
	return Ace::Sequence::Homol->new_homol($classe,$name,$db,$start,$end);
    }

    # General case:
    my $obj = $class->new($tag=>$data[0],$self->db);

    return $obj if defined $obj;

    # Last resort, return a Text
    return $class->new(Text=>$data[0]);
}

sub sub_SeqFeature {
  return wantarray ? () : 0;
}

1;

=head1 NAME

Ace::Sequence::Feature - Examine Sequence Feature Tables

=head1 SYNOPSIS

    # open database connection and get an Ace::Object sequence
    use Ace::Sequence;

    # get a megabase from the middle of chromosome I
    $seq = Ace::Sequence->new(-name   => 'CHROMOSOME_I,
                              -db     => $db,
			      -offset => 3_000_000,
			      -length => 1_000_000);