← Back to branch summary

~ubuntu-branches/ubuntu/oneiric/bioperl/oneiric

« back to all changes in this revision

Viewing changes to Bio/DB/Fasta.pm

  • Committer: Bazaar Package Importer
  • Author(s): Matt Hope
  • Date: 2004-04-18 14:24:11 UTC
  • mfrom: (1.2.1 upstream) (2.1.1 warty)
  • Revision ID: james.westby@ubuntu.com-20040418142411-gr92uexquw4w8liq
Tags: 1.4-1
* New upstream release
* Examples and working code are installed by default to usr/bin,
  this has been moved to usr/share/doc/bioperl/bin

Show diffs side-by-side

added added

removed removed

Lines of Context:
Bio/DB/Fasta.pm
392
392
use File::Basename qw(basename dirname);
393
393
use Bio::DB::SeqI;
394
394
use Bio::Root::Root;
395
 
use vars qw($VERSION @ISA);
 
395
use vars qw(@ISA);
396
396
 
397
397
@ISA = qw(Bio::DB::SeqI Bio::Root::Root);
398
398
 
399
 
$VERSION = '1.03';
400
 
 
401
399
*seq = *sequence = \&subseq;
402
400
*ids = \&get_all_ids;
403
401
*get_seq_by_primary_id = *get_Seq_by_acc  = \&get_Seq_by_id;
501
499
  my %offsets;
502
500
  my $flags = $write ? O_CREAT|O_RDWR : O_RDONLY;
503
501
  my @dbmargs = $self->dbmargs;
504
 
  tie %offsets,'AnyDBM_File',$index,$flags,0644,@dbmargs or $self->throw( "Can't open cache file: $!");
 
502
  tie %offsets,'AnyDBM_File',$index,$flags,0644,@dbmargs or $self->throw( "Can't open cache file $index: $!");
505
503
  return \%offsets;
506
504
}
507
505
 
532
530
  unlink $index if $force_reindex;
533
531
 
534
532
  # get the modification time of the index
535
 
  my $indextime = (stat($index))[9] || 0;
 
533
  my $indextime   = 0;
 
534
  for my $suffix('','.pag','.dir') {
 
535
    $indextime ||= (stat("${index}${suffix}"))[9];
 
536
  }
536
537
 
537
538
  # get the most recent modification time of any of the contents
538
539
  my $modtime = 0;
573
574
sub get_Seq_by_id {
574
575
  my $self = shift;
575
576
  my $id   = shift;
576
 
  return Bio::PrimarySeqI->new($self,$id);
 
577
  return unless exists $self->{offsets}{$id};
 
578
  return Bio::PrimarySeq::Fasta->new($self,$id);
577
579
}
578
580
 
579
581
=head2 index_file
598
600
  unlink $index if $force_reindex;
599
601
 
600
602
  # get the modification time of the index
601
 
  my $indextime = (stat($index))[9];
602
 
  my $modtime   = (stat($file))[9];
 
603
  my $indextime = (stat($index))[9] || 0;
 
604
  my $modtime   = (stat($file))[9]  || 0;
603
605
 
604
606
  my $reindex = $force_reindex || $indextime < $modtime;
605
607
  my $offsets = $self->_open_index($index,$reindex) or return;
670
672
 
671
673
  my $fh = IO::File->new($file) or $self->throw( "Can't open $file: $!");
672
674
  warn "indexing $file\n" if $self->{debug};
673
 
  my ($offset,$id,$linelength,$type,$firstline,$count,%offsets);
 
675
  my ($offset,$id,$linelength,$type,$firstline,$count,$termination_length,%offsets);
674
676
  while (<$fh>) {               # don't try this at home
 
677
    $termination_length ||= /\r\n$/ ? 2 : 1;  # account for crlf-terminated Windows files
675
678
    if (/^>(\S+)/) {
676
679
      print STDERR "indexed $count sequences...\n" 
677
680
        if $self->{debug} && (++$count%1000) == 0;
678
681
      my $pos = tell($fh);
679
682
      if ($id) {
680
683
        my $seqlength    = $pos - $offset - length($_) - 1;
681
 
        $seqlength      -= int($seqlength/$linelength);
 
684
        $seqlength      -= $termination_length * int($seqlength/$linelength);
682
685
        $offsets->{$id}  = $self->_pack($offset,$seqlength,
683
686
                                        $linelength,$firstline,
684
687
                                        $type,$base);
693
696
  # deal with last entry
694
697
  if ($id) {
695
698
    my $pos = tell($fh);
696
 
 
697
 
#    my $seqlength   = $pos - $offset - length($_) - 1;
698
 
    # $_ is always null should not be part of this calculation
699
699
    my $seqlength   = $pos - $offset  - 1;
700
700
 
701
701
    if ($linelength == 0) { # yet another pesky empty chr_random.fa file
702
702
      $seqlength = 0;
703
703
    } else {
704
 
      $seqlength -= int($seqlength/$linelength);
 
704
      $seqlength -= $termination_length * int($seqlength/$linelength);
705
705
    };
706
706
    $offsets->{$id} = $self->_pack($offset,$seqlength,
707
707
                                   $linelength,$firstline,
708
708
                                   $type,$base);
709
709
  }
 
710
  $offsets->{__termination_length} = $termination_length;
710
711
  return \%offsets;
711
712
}
712
713
 
786
787
sub path2fileno {
787
788
  my $self = shift;
788
789
  my $path = shift;
789
 
  unless (exists $self->{offsets}{"__path_$path"}) {
 
790
  if ( !defined $self->{offsets}{"__path_$path"} ) {
790
791
    my $fileno  = ($self->{offsets}{"__path_$path"} = 0+ $self->{fileno}++);
791
792
    $self->{offsets}{"__file_$fileno"} = $path;
792
793
  }
828
829
  seek($fh,$filestart,0);
829
830
  read($fh,$data,$filestop-$filestart+1);
830
831
  $data =~ s/\n//g;
 
832
  $data =~ s/\r//g;
831
833
  if ($reversed) {
832
834
    $data = reverse $data;
833
835
    $data =~ tr/gatcGATC/ctagCTAG/;
864
866
  my ($offset,$seqlength,$linelength,$firstline,$type,$file) = $self->_unpack($self->{offsets}{$id});
865
867
  $a = 0            if $a < 0;
866
868
  $a = $seqlength-1 if $a >= $seqlength;
867
 
  $offset + $linelength * int($a/($linelength-1)) + $a % ($linelength-1);
 
869
  my $tl = $self->{offsets}{__termination_length};
 
870
  $offset + $linelength * int($a/($linelength-$tl)) + $a % ($linelength-$tl);
868
871
}
869
872
 
870
873
sub fhcache {
960
963
eval {
961
964
  require Bio::PrimarySeqI;
962
965
  require Bio::Root::Root;
963
 
} && (@ISA = ('Bio::PrimarySeqI','Bio::Root::Root'));
 
966
} && (@ISA = ('Bio::Root::Root','Bio::PrimarySeqI'));
964
967
 
965
968
sub new {
966
969
  my $class = shift;
980
983
 
981
984
sub subseq {
982
985
  my $self = shift;
 
986
  $self->trunc(@_)->seq();      
 
987
}
 
988
 
 
989
sub trunc {
 
990
  my $self = shift;
983
991
  my ($start,$stop) = @_;
984
992
  $self->throw("Stop cannot be smaller than start")  unless $start <= $stop;
985
993
  return $self->{start} <= $self->{stop} ?  $self->new($self->{db},
990
998
                                                       $self->{id},
991
999
                                                       $self->{start}-($start-1),
992
1000
                                                       $self->{start}-($stop-1)
993
 
                                                      );
 
1001
                                                      );  
994
1002
        
995
1003
}
996
1004
 
1026
1034
  return $self->{db}->length($self->{id});
1027
1035
}
1028
1036
 
 
1037
sub description  { 
 
1038
    my $self = shift;
 
1039
    return '';
 
1040
}
 
1041
 
 
1042
*desc = \&description;
 
1043
 
1029
1044
#-------------------------------------------------------------
1030
1045
# stream-based access to the database
1031
1046
#
1048
1063
sub next_seq {
1049
1064
  my $self = shift;
1050
1065
  my ($key,$db) = @{$self}{'key','db'};
 
1066
  while ($key =~ /^__/) {
 
1067
    $key = $db->NEXTKEY($key);
 
1068
    return unless defined $key;
 
1069
  }
1051
1070
  my $value = $db->get_Seq_by_id($key);
1052
1071
  $self->{key} = $db->NEXTKEY($key);
1053
1072
  $value;