sync w/ main trunk

[bioperl-live.git] / Bio / LiveSeq / Transcript.pm

# $Id$
#
# bioperl module for Bio::LiveSeq::Transcript
#
# Please direct questions and support issues to <bioperl-l@bioperl.org> 
#
# Cared for by Joseph Insana <insana@ebi.ac.uk> <jinsana@gmx.net>
#
# Copyright Joseph Insana
#
# You may distribute this module under the same terms as perl itself
#
# POD documentation - main docs before the code

=head1 NAME

Bio::LiveSeq::Transcript - Transcript class for LiveSeq

=head1 SYNOPSIS

  # documentation needed

=head1 DESCRIPTION

This stores informations about coding sequences (CDS).
The implementation is that a Transcript object accesses a collection of
Exon objects, inferring from them the nucleotide structure and sequence.

=head1 AUTHOR - Joseph A.L. Insana

Email:  Insana@ebi.ac.uk, jinsana@gmx.net

=head1 APPENDIX

The rest of the documentation details each of the object
methods. Internal methods are usually preceded with a _

=cut

# Let the code begin...

package Bio::LiveSeq::Transcript;

use strict;
# use Carp qw(carp cluck);
use Bio::LiveSeq::Exon; # uses Exon to create new exon in case of deletion
use base qw(Bio::LiveSeq::SeqI);

=head2 new

  Title   : new
  Usage   : $transcript = Bio::LiveSeq::Transcript->new(-exons => \@obj_refs);

  Function: generates a new Bio::LiveSeq::Transcript
  Returns : reference to a new object of class Transcript
  Errorcode -1
  Args    : reference to an array of Exon object references

=cut

sub new {
  my ($thing, %args) = @_;
  my $class = ref($thing) || $thing;
  my ($obj,%transcript);

  my @exons=@{$args{-exons}};

  $obj = \%transcript;
  $obj = bless $obj, $class;

  unless (@exons) {
    $obj->warn("$class not initialised because exons array empty");
    return(-1);
  }

  # now useless, after start and end methods have been overridden here
  my $firstexon = $exons[0];
  #my $lastexon = $exons[-1];
  #my $start = $firstexon->start;
  #my $end = $lastexon->end;
  my $strand = $firstexon->strand;
  my $seq = $firstexon->{'seq'};
  $obj->alphabet('rna');

  unless (_checkexons(\@exons)) {
    $obj->warn("$class not initialised because of problems in the exon structure");
    return(-1);
  }
  $obj->{'strand'}=$strand;
  $obj->{'exons'}=\@exons;
  $obj->{'seq'}=$seq;

  # set Transcript into each Exon
  my $exon;
  foreach $exon (@exons) {
    $exon->{'transcript'}=$obj;
  }
  return $obj;
}


=head2 all_Exons

 Title   : all_Exons
 Usage   : $transcript_obj->all_Exons()
 Function: returns references to all Exon objects the Transcript is composed of
 Example : foreach $exon ($transcript->all_Exons()) { do_something }
 Returns : array of object references
 Args    : none

=cut

sub all_Exons {
  my $self=shift;
  my $exonsref=$self->{'exons'};
  my @exons=@{$exonsref};
  my @newexons;
  my $exon;
  foreach $exon (@exons) {
    unless ($exon->obj_valid) {
      $self->warn("$exon no more valid, start or end label lost, skipping....",1); # ignorable
    } else {
      push(@newexons,$exon);
    }
  }
  if ($#exons != $#newexons) {
    # update exons field
    $self->{'exons'}=\@newexons;
  }
  return (@newexons);
}

=head2 downstream_seq

 Title   : downstream_seq
 Usage   : $transcript_obj->downstream_seq()
         : $transcript_obj->downstream_seq(64)
 Function: returns a string of nucleotides downstream of the end of the
           CDS. If there is some information of the real mRNA, from features in
           an attached Gene object, it will return up to those boundaries.
           Otherwise it will return 1000 nucleotides.
           If an argument is given it will override the default 1000 number
           and return instead /that/ requested number of nucleotides.
           But if a Gene object is attached, this argument will be ignored.
 Returns : string
 Args    : an optional integer number of nucleotides to be returned instead of
           the default if no gene attached

=cut

sub downstream_seq {
  my ($self,$howmany)=@_;
  my $str;
  if (defined ($howmany)) {
    unless ($howmany > 0) {
      $self->throw("No sense in asking less than 1 downstream nucleotides!");
    }
  } else {
    unless ($self->{'seq'}->alphabet eq 'rna') { # if rna retrieve until the end
      #$str=$DNAobj->labelsubseq($self->end,undef,undef,"unsecuremoderequested");
      #return(substr($str,1)); # delete first nucleotide that is the last of Transcript
      if ($self->gene) { # if there is Gene object attached fetch relevant info
        $str=$self->{'seq'}->labelsubseq($self->end,undef,$self->gene->maxtranscript->end); # retrieve from end of this Transcript to end of the maxtranscript
        $str=substr($str,1); # delete first nucleotide that is the last of Transcript
        if (CORE::length($str) > 0) {
          return($str);
        } else { # if there was no downstream through the gene's maxtranscript, go the usual way
          $howmany = 1000;
        }
      } else {
        $howmany = 1000;
      }
    }
  }
  my @exons=$self->all_Exons;
  my $strand=$self->strand();
  my $lastexon=$exons[-1];
  my $lastexonlength=$lastexon->length;
  # $howmany nucs after end of last exon
  #my $downstream_seq=$lastexon->subseq($lastexonlength+1,undef,$howmany);
  my $downstream_seq;

  if ($howmany) {
      $downstream_seq=substr($lastexon->labelsubseq($self->end,$howmany,undef,"unsecuremoderequested"),1);
  } else {
    if ($strand == 1) {
      $downstream_seq=substr($lastexon->labelsubseq($self->end,undef,$self->{'seq'}->end,"unsecuremoderequested"),1);
    } else {
      $downstream_seq=substr($lastexon->labelsubseq($self->end,undef,$self->{'seq'}->start,"unsecuremoderequested"),1);
    }
  }
  return $downstream_seq;
}

=head2 upstream_seq

 Title   : upstream_seq
 Usage   : $transcript_obj->upstream_seq()
         : $transcript_obj->upstream_seq(64)
 Function: just like downstream_seq but returns nucleotides before the ATG
 Note    : the default, if no Gene information present and no nucleotides
           number given, is to return up to 400 nucleotides.

=cut

sub upstream_seq {
  my ($self,$howmany)=@_;
  if (defined ($howmany)) {
    unless ($howmany > 0) {
      $self->throw("No sense in asking less than 1 upstream nucleotides!");
    }
  } else {
    unless ($self->{'seq'}->alphabet eq 'rna') { # if rna retrieve from the start
      if ($self->gene) { # if there is Gene object attached fetch relevant info
        my $str=$self->{'seq'}->labelsubseq($self->gene->maxtranscript->start,undef,$self->start); # retrieve from start of maxtranscript to start of this Transcript
        chop $str; # delete last nucleotide that is the A of starting ATG
        if (length($str) > 0) {
          return($str);
        } else { # if there was no upstream through the gene's maxtranscript, go the usual way
          $howmany = 400;
        }
      } else {
        $howmany = 400;
      }
    }
  }
  my @exons=$self->all_Exons;
  my $firstexon=$exons[0];
  
  my $upstream_seq;
  my $strand=$self->strand();

  if ($howmany) {# $howmany nucs before begin of first exon
    my $labelbefore=$firstexon->label(-$howmany,$firstexon->start);
    if ($labelbefore < 1) {
      if ($strand == 1) {
        $labelbefore=$self->{'seq'}->start;
      } else {
        $labelbefore=$self->{'seq'}->end;
      }
    }
    $upstream_seq=$firstexon->labelsubseq($labelbefore,undef,$firstexon->start,"unsecuremoderequested");
    chop $upstream_seq;
  } else {
    if ($strand == 1) {
      $upstream_seq=$firstexon->labelsubseq($self->{'seq'}->start,undef,$self->start,"unsecuremoderequested");
      chop $upstream_seq; # delete last nucleotide that is the A of starting ATG
    } else {
      $upstream_seq=$firstexon->labelsubseq($self->{'seq'}->end,undef,$self->start,"unsecuremoderequested");
      chop $upstream_seq; # delete last nucleotide that is the A of starting ATG
    }
  }
  return $upstream_seq;
}

# These get redefined here, overriding the SeqI one because they draw their
# information from the Exons a Transcript is built of
# optional argument: firstlabel. If not given, it checks coordinate_start
#                                This is useful when called by Translation
#                                also used by _delete
sub label {
  my ($self,$position,$firstlabel)=@_;
  unless ($position) {  # if position = 0 complain ?
    $self->warn("Position not given or position 0");
    return (-1);
  }
  my ($start,$end,$strand)=($self->start(),$self->end(),$self->strand());
  my ($label,@labels,$length,$arraypos);
  unless (defined ($firstlabel)) {
    $firstlabel=$self->coordinate_start; # this is inside Transcript obj
  }
  my $coord_pos=$self->_inside_position($firstlabel);
  $length=$self->length;
  #if ($strand == 1) {
    if ($position < 1) {
      $position++; # to account for missing of 0 position
    }
    $arraypos=$position+$coord_pos-2;
    #print "\n=-=-=-=-DEBUG: arraypos $arraypos, pos $position, coordpos: $coord_pos";
    if ($arraypos < 0) {
      $label=$self->{'seq'}->label($arraypos,$start,$strand); #?
    } elsif ($arraypos >= $length) {
      $label=$self->{'seq'}->label($arraypos-$length+2,$end,$strand); #?
    } else { # inside the Transcript
      @labels=$self->all_labels;
      $label=$labels[$arraypos];
    }
  #}
}

# argument: label
# returns: position of label according to coord_start
# errorcode: 0 label not found
# optional argument: firstlabel. If not given, it checks coordinate_start
#                                This is useful when called by Translation
sub position {
  my ($self,$label,$firstlabel)=@_;
  unless ($self->{'seq'}->valid($label)) {
    $self->warn("label is not valid");
    return (0);
  }
  unless (defined ($firstlabel)) {
    $firstlabel=$self->coordinate_start; # this is inside Transcript obj
  }
  if ($label == $firstlabel) {
    return (1);
  }
  my ($start,$end,$strand)=($self->start(),$self->end(),$self->strand());
  my ($position,$in_pos,$out_pos,$coord_pos);
  my $length=$self->length;
  $coord_pos=$self->_inside_position($firstlabel);
  if ($self->valid($label)) { # if label is inside the Transcript
    $in_pos=$self->_inside_position($label);
    $position=$in_pos-$coord_pos+1;
    if ($position <= 0) {
      return ($position-1); # accounts for the missing of the 0 position
    }
  } else {
    if ($self->follows($end,$label)) { # label after end of transcript
      $out_pos=$self->{'seq'}->position($label,$end,$strand);
      #print "\n+++++++++DEBUG label $label FOLLOWS end $end outpos $out_pos coordpos $coord_pos";
      $position=$out_pos+$length-$coord_pos;
    } elsif ($self->follows($label,$start)) { # label before begin of transcript
      #print "\n+++++++++DEBUG label $label BEFORE start $start outpos $out_pos coordpos $coord_pos";
      $out_pos=$self->{'seq'}->position($label,$start,$strand);
      $position=$out_pos-$coord_pos+1;
    } else { # label is in intron (not valid, not after, not before)!
      $self->warn("Cannot give position of label pointing to intron according to CDS numbering!",1);
      return (0);
    }
  }
  return ($position);
}

sub seq {
  my $self=shift;
  my ($exon,$str);
  my @exons=$self->all_Exons();
  foreach $exon (@exons) {
    $str .= $exon->seq();
  }
  return $str;
}

sub length {
  my $self=shift;
  my ($exon,$length);
  my @exons=$self->all_Exons();
  foreach $exon (@exons) {
    $length += $exon->length();
  }
  return $length;
}

sub all_labels {
  my $self=shift;
  my ($exon,@labels);
  my @exons=$self->all_Exons();
  foreach $exon (@exons) {
    push (@labels,$exon->all_labels());
  }
  return @labels;
}

# redefined here so that it will retrieve effective subseq without introns
# otherwise it would have retrieved an underlying DNA (possibly with introns)
# subsequence
# Drawback: this is really bulky, label->position and then a call to
# subseq that will do the opposite position-> label
#
# one day this can be rewritten as the main one so that the normal subseq
# will rely on this one and hence avoid this double (useless and lengthy)
# conversion between labels and positions
sub old_labelsubseq {
  my ($self,$start,$length,$end)=@_;
  my ($pos1,$pos2);
  if ($start) {
    unless ($self->valid($start)) {
      $self->warn("Start label not valid"); return (-1);
    }
    $pos1=$self->position($start);
  }
  if ($end) {
    if ($end == $start) {
      $length=1;
    } else {
      unless ($self->valid($end)) {
        $self->warn("End label not valid"); return (-1);
      }
      unless ($self->follows($start,$end) == 1) {
        $self->warn("End label does not follow Start label!"); return (-1);
      }
      $pos2=$self->position($end);
      undef $length;
    }
  }
  return ($self->subseq($pos1,$pos2,$length));
}

# rewritten, eventually

sub labelsubseq {
  my ($self,$start,$length,$end,$unsecuremode)=@_;
  unless (defined $unsecuremode && 
          $unsecuremode eq "unsecuremoderequested") 
  { # to skip security checks (faster)
    if ($start) {
      unless ($self->valid($start)) {
        $self->warn("Start label not valid"); return (-1);
      }
    } else {
      $start=$self->start;
    }
    if ($end) {
      if ($end == $start) {
        $length=1;
        undef $end;
      } else {
        undef $length; # end argument overrides length argument
        unless ($self->valid($end)) {
          $self->warn("End label not valid"); return (-1);
        }
        unless ($self->follows($start,$end) == 1) {
          $self->warn("End label does not follow Start label!"); return (-1);
        }
      }
    } else {
      $end=$self->end;
    }
  }
  my ($seq,$exon,$startexon,$endexon); my @exonlabels;
  my @exons=$self->all_Exons;
  EXONCHECK:
  foreach $exon (@exons) {
    if ((!(defined($startexon)))&&($exon->valid($start))) { # checks only if not yet found
      $startexon=$exon;
    }
    if ($exon->valid($end)) {
      $endexon=$exon;
    }
    if ((!(defined($seq)) && (defined($startexon)))) { # initializes only once
      if ((defined($endexon)) && ($endexon eq $startexon)) { # then perfect, we are finished
        if ($length) {
          $seq = $startexon->labelsubseq($start,$length,undef,"unsecuremoderequested");


          last EXONCHECK;
        } else {
          $seq = $startexon->labelsubseq($start,undef,$end,"unsecuremoderequested");
        }
        last EXONCHECK;
      } else { # get up to the end of the exon
        $seq = $startexon->labelsubseq($start,undef,undef,"unsecuremoderequested");
      }
    }
    if (($startexon)&&($exon ne $startexon)) {
      if (defined($endexon)) { # we arrived to the last exon
        $seq .= $endexon->labelsubseq(undef,undef,$end,"unsecuremoderequested"); # get from the start of the exon
        last EXONCHECK;

      } elsif (defined($startexon)) { # we are in a whole-exon-in-the-middle case
          $seq .= $exon->seq; # we add it completely to the seq
      } # else, we still have to reach the start point, exon useless, we move on
      if ($length) { # if length argument specified
        if (($seq && (CORE::length($seq) >= $length))) {
          last EXONCHECK;
        }
      }
    }
  }
  if ($length) {
    return (substr($seq,0,$length));
  } else {
    return ($seq);
  }
}


# argument: label
# returns: the objref and progressive number of the Exon containing that label
# errorcode: -1
sub in_which_Exon {
  my ($self,$label)=@_;
  my ($count,$exon);
  my @exons=$self->all_Exons;
  foreach $exon (@exons) {
    $count++; # 1st exon is numbered "1"
    if ($exon->valid($label)) {
      return ($exon,$count)
    }
  }
  return (-1); # if nothing found
}

# recoded to exploit the new fast labelsubseq()
# valid only inside Transcript
sub subseq {
  my ($self,$pos1,$pos2,$length) = @_;
  my ($str,$startlabel,$endlabel);
  if (defined ($pos1)) {
    if ($pos1 == 0) {  # if position = 0 complain
      $self->warn("Position cannot be 0!"); return (-1);
    }
    if ((defined ($pos2))&&($pos1>$pos2)) {
      $self->warn("1st position($pos1) cannot be > 2nd position($pos2)!"); return (-1);
    }
    $startlabel=$self->label($pos1);
    unless ($self->valid($startlabel)) {
      $self->warn("Start label not valid"); return (-1);
    }
    if ($startlabel < 1) {
      $self->warn("position $pos1 not valid as start of subseq!"); return (-1);
    }
  } else {
    $startlabel=$self->start;
  }
  if (defined ($pos2)) {
    if ($pos2 == 0) {  # if position = 0 complain
      $self->warn("Position cannot be 0!"); return (-1);
    }
    undef $length;
    if ((defined ($pos1))&&($pos1>$pos2)) {
      $self->warn("1st position($pos1) cannot be > 2nd position($pos2)!"); return (-1);
    }
    $endlabel=$self->label($pos2);
    unless ($self->valid($endlabel)) {
      $self->warn("End label not valid"); return (-1);
    }
    if ($endlabel < 1) {
      $self->warn("position $pos2 not valid as end of subseq!"); return (-1);
    }
  } else {
    unless (defined ($length)) {
      $endlabel=$self->end;
    }
  }
  return ($self->labelsubseq($startlabel,$length,$endlabel,"unsecuremoderequested"));
}

# works only inside the transcript, complains if asked outside
sub old_subseq {
  my ($self,$pos1,$pos2,$length) = @_;
  my ($str,$startcount,$endcount,$seq,$seqlength);
  if (defined ($length)) {
    if ($length < 1) {
      $self->warn("No sense asking for a subseq of length < 1");
      return (-1);
    }
  }
  my $firstlabel=$self->coordinate_start; # this is inside Transcript obj
  my $coord_pos=$self->_inside_position($firstlabel); # TESTME old
  $seq=$self->seq;
  $seqlength=CORE::length($seq);
  unless (defined ($pos1)) {
    $startcount=1+$coord_pos-1; # i.e. coord_pos
  } else {
    if ($pos1 == 0) {  # if position = 0 complain
      $self->warn("Position cannot be 0!"); return (-1);
    } elsif ($pos1 < 0) {
      $pos1++;
    }
    if ((defined ($pos2))&&($pos1>$pos2)) {
      $self->warn("1st position ($pos1) cannot be > 2nd position ($pos2)!");
      return (-1);
    }
    $startcount=$pos1+$coord_pos-1;
  }
  unless (defined ($pos2)) {
  ;
  } else {
    if ($pos2 == 0) {  # if position = 0 complain
      $self->warn("Position cannot be 0!"); return (-1);
    } elsif ($pos2 < 0) {
      $pos2++;
    }
    if ((defined ($pos1))&&($pos1>$pos2)) {
      $self->warn("1st position ($pos1) cannot be > 2nd position ($pos2)!");
      return (-1);
    }
    $endcount=$pos2+$coord_pos-1;
    if ($endcount > $seqlength) {
      #print "\n###DEBUG###: pos1 $pos1 pos2 $pos2 coordpos $coord_pos endcount $endcount seqln $seqlength\n";
      $self->warn("Cannot access end position after the end of Transcript");
      return (-1);
    }
    $length=$endcount-$startcount+1;
  }
  #print "\n###DEBUG pos1 $pos1 pos2 $pos2 start $startcount end $endcount length $length coordpos $coord_pos\n";
  my $offset=$startcount-1;
  if ($offset < 0) {
    $self->warn("Cannot access startposition before the beginning of Transcript, returning from start",1); # ignorable
    return (substr($seq,0,$length));
  } elsif ($offset >= $seqlength) {
    $self->warn("Cannot access startposition after the end of Transcript");
    return (-1);
  } else {
    $str=substr($seq,$offset,$length);
    if (CORE::length($str) < $length) {
      $self->warn("Attention, cannot return the length requested ".
                  "for subseq",1) if $self->verbose > 0; # ignorable
    }
    return $str;
  }
}

# redefined so that it doesn't require other methods (after deletions) to
# reset it.
sub start {
  my $self = shift;
  my $exonsref=$self->{'exons'};
  my @exons=@{$exonsref};
  return ($exons[0]->start);
}

sub end {
  my $self = shift;
  my $exonsref=$self->{'exons'};
  my @exons=@{$exonsref};
  return ($exons[-1]->end);
}


# internal methods begin here

# returns: position of label in transcript's all_labels
#          with STARTlabel == 1
# errorcode 0 -> label not found
# argument: label
sub _inside_position {
  my ($self,$label)=@_;
  my ($start,$end,$strand)=($self->start(),$self->end(),$self->strand());
  my ($position,$checkme);
  my @labels=$self->all_labels;
  foreach $checkme (@labels) {
    $position++;
    if ($label == $checkme) {
      return ($position);
    }
  }
  return (0);
}

# returns 1 OK or 0 ERROR
# arguments: reference to array of Exon object references
sub _checkexons {
  my ($exon,$thisstart);
  my $self=$exon;
  my $exonsref=$_[0];
  my @exons=@{$exonsref};

  my $firstexon = $exons[0];

  unless (ref($firstexon) eq "Bio::LiveSeq::Exon") {
    $self->warn("Object not of class Exon");
    return (0);
  }
  my $strand = $firstexon->strand;

  my $prevend = $firstexon->end;
  shift @exons; # skip first one
  foreach $exon (@exons) {
    unless (ref($exon) eq "Bio::LiveSeq::Exon") { # object class check
      $self->warn("Object not of class Exon");
      return (0);
    }
    if ($exon->strand != $strand) { # strand consistency check
      $self->warn("Exons' strands not consistent when trying to create Transcript");
      return (0);
    }
    $thisstart = $exon->start;
    unless ($exon->{'seq'}->follows($prevend,$thisstart,$strand)) {
      $self->warn("Exons not in correct order when trying to create Transcript");
      return (0);
    }
    $prevend = $exon->end;
  }
  return (1);
}

=head2 get_Translation

  Title   : valid
  Usage   : $translation = $obj->get_Translation()
  Function: retrieves the reference to the object of class Translation (if any)
            attached to a LiveSeq object
  Returns : object reference
  Args    : none

=cut

sub get_Translation {
  my $self=shift;
  return ($self->{'translation'}); # this is set when Translation->new is called
}

# this checks so that deletion spanning multiple exons is
# handled accordingly and correctly
# arguments: begin and end label of a deletion
# this is called BEFORE any deletion in the chain
sub _deletecheck {
  my ($self,$startlabel,$endlabel)=@_;
  my $exonsref=$self->{'exons'};
  my @exons=@{$exonsref};
  my ($startexon,$endexon,$exon);
  $startexon=$endexon=0;
  foreach $exon (@exons) {
    if (($startexon == 0)&&($exon->valid($startlabel))) {
      $startexon=$exon; # exon containing start of deletion
    }
    if (($endexon == 0)&&($exon->valid($endlabel))) {
      $endexon=$exon; # exon containing end of deletion
    }
    if (($startexon)&&($endexon)) {
      last; # don't check further
    }
  }
  my $nextend=$self->label(2,$endlabel); # retrieve the next label
  my $prevstart=$self->label(-1,$startlabel); # retrieve the prev label

  if ($startexon eq $endexon) { # intra-exon deletion
    if (($startexon->start eq $startlabel) && ($startexon->end eq $endlabel)) {
      # let's delete the entire exon
      my @newexons;
      foreach $exon (@exons) {
        unless ($exon eq $startexon) {
          push(@newexons,$exon);
        }
      }
      $self->{'exons'}=\@newexons;
    } elsif ($startexon->start eq $startlabel) { # special cases
      $startexon->{'start'}=$nextend; # set a new start of exon
    } elsif ($startexon->end eq $endlabel) {
      $startexon->{'end'}=$prevstart; # set a new end of exon
    } else {
      return; # no problem
    }
  } else { # two new exons to be created, inter-exons deletion
    my @newexons;
    my $exonobj;
    my $dna=$self->{'seq'};
    my $strand=$self->strand;
    my $notmiddle=1; # flag for skipping exons in the middle of deletion
    foreach $exon (@exons) {
      if ($exon eq $startexon) {
        $exonobj=Bio::LiveSeq::Exon->new('-seq'=>$dna,'-start'=>$exon->start,'-end'=>$prevstart,'-strand'=>$strand); # new partial exon
        push(@newexons,$exonobj);
        $notmiddle=0; # now we enter totally deleted exons
      } elsif ($exon eq $endexon) {
        $exonobj=Bio::LiveSeq::Exon->new('-seq'=>$dna,'-start'=>$nextend,'-end'=>$exon->end,'-strand'=>$strand); # new partial exon
        push(@newexons,$exonobj);
        $notmiddle=1; # exiting totally deleted exons
      } else {
        if ($notmiddle) { # if before or after exons with deletion
          push(@newexons,$exon); 
        }# else skip them
      }
    }
    $self->{'exons'}=\@newexons;
  }
}

=head2 translation_table

 Title   : translation_table
 Usage   : $name = $obj->translation_table;
         : $name = $obj->translation_table(11);
 Function: Returns or sets the translation_table used for translating the
           transcript.
           If it has never been set, it will return undef.
 Returns : an integer

=cut

sub translation_table {
  my ($self,$value) = @_;
  if (defined $value) {
    $self->{'translation_table'} = $value;
  }
  unless (exists $self->{'translation_table'}) {
    return;
  } else {
    return $self->{'translation_table'};
  }
}

=head2 frame

 Title   : frame
 Usage   : $frame = $transcript->frame($label);
 Function: Returns the frame of a particular nucleotide.
           Frame can be 0 1 or 2 and means the position in the codon triplet
           of the particulat nucleotide. 0 is the first codon_position.
           Codon_position (1 2 3) is simply frame+1.
           If the label asked for is not inside the Transcript, -1 will be
           returned.
 Args    : a label
 Returns : 0 1 or 2
 Errorcode -1

=cut

# args: label
# returns: frame of nucleotide (0 1 2)
# errorcode: -1
sub frame {
  my ($self,$inputlabel)=@_;
  my @labels=$self->all_labels;
  my ($label,$frame,$count);
  foreach $label (@labels) {
    if ($inputlabel == $label) {
      return ($count % 3);
    }
    $count++; # 0 1 2 3 4....
  }
  return (-1); # label not found amid Transcript labels
}

1;