use Bio::DB::Fasta;

  # create database from directory of fasta files
  my $db      = Bio::DB::Fasta->new('/path/to/fasta/files');

  # simple access (for those without Bioperl)
  my $seq     = $db->seq('CHROMOSOME_I',4_000_000 => 4_100_000);
  my $revseq  = $db->seq('CHROMOSOME_I',4_100_000 => 4_000_000);
  my @ids     = $db->ids;
  my $length  = $db->length('CHROMOSOME_I');
  my $alphabet = $db->alphabet('CHROMOSOME_I');
  my $header  = $db->header('CHROMOSOME_I');

  # Bioperl-style access
  my $db      = Bio::DB::Fasta->new('/path/to/fasta/files');

  my $obj     = $db->get_Seq_by_id('CHROMOSOME_I');
  my $seq     = $obj->seq;
  my $subseq  = $obj->subseq(4_000_000 => 4_100_000);
  my $length  = $obj->length;
  # (etc)

  # Bio::SeqIO-style access
  my $stream  = Bio::DB::Fasta->new('/path/to/fasta/files')->get_PrimarySeq_stream;
  while (my $seq = $stream->next_seq) {
    # Bio::PrimarySeqI stuff
  }

  my $fh = Bio::DB::Fasta->newFh('/path/to/fasta/files');
  while (my $seq = <$fh>) {
    # Bio::PrimarySeqI stuff
  }

  # tied hash access
  tie %sequences,'Bio::DB::Fasta','/path/to/fasta/files';
  print $sequences{'CHROMOSOME_I:1,20000'};

  @AnyDBM_File::ISA = qw(DB_File GDBM_File NDBM_File SDBM_File)

    shift->throw("Read-only database");

    shift->throw("Read-only database");

  my $self = shift;
  if ($self->{indexing}) {  # killed prematurely, so index file is no good!
    warn "indexing was interrupted, so unlinking $self->{indexing}";
    unlink $self->{indexing};
  }
}

#-------------------------------------------------------------
# Bio::PrimarySeqI compatibility
#
package Bio::PrimarySeq::Fasta;
use overload '""' => 'display_id';

use vars '@ISA';
eval {
  require Bio::PrimarySeqI;
  require Bio::Root::Root;
} && (@ISA = ('Bio::Root::Root','Bio::PrimarySeqI'));

  defined shift->offset(@_);

  shift->subseq(@_);

 tied(%{shift->{offsets}})->NEXTKEY(@_);

  my $self = shift;
  $self->next_seq;
}

1;

__END__

    shift->throw("Read-only database");

  my $class = shift;
  my $db    = shift;
  return $class->new($db);

  my $self = shift;
  return $self->new(@_);

  my $self = shift;
  my ($index,$write) = @_;
  my %offsets;
  my $flags = $write ? O_CREAT|O_RDWR : O_RDONLY;
  my @dbmargs = $self->dbmargs;
  tie %offsets,'AnyDBM_File',$index,$flags,0644,@dbmargs or $self->throw( "Can't open cache file: $!");
  return\% offsets;

  shift;
  pack STRUCT,@_;

  shift;
  local $_ = shift;
  return /^[gatcnGATCN*-]+$/   ? DNA
         : /^[gaucnGAUCN*-]+$/ ? RNA
	 : PROTEIN;

  shift;
  unpack STRUCT,shift;

  my $self = shift;
  return "unknown";

  my $self = shift;
  my $id   = shift;
  my $offset = $self->{offsets}{$id} or return;
  my $type = ($self->_unpack($offset))[4];
  return $type == DNA ? 'dna'
         : $type == RNA ? 'rna'
         : 'protein';

  my $self = shift;
  return $self->{db}->alphabet($self->{id});

  my $self = shift;
  my ($file,$offsets) = @_;
  my $base = $self->path2fileno(basename($file));

  my $fh = IO::File->new($file) or $self->throw( "Can't open $file: $!");
  warn "indexing $file\n" if $self->{debug};
  my ($offset,$id,$linelength,$type,$firstline,$count,%offsets);
  while (<$fh>) {		# don't try this at home
    if (/^>(\S+)/) {
      print STDERR "indexed $count sequences...\n" 
	if $self->{debug} && (++$count%1000) == 0;
      my $pos = tell($fh);
      if ($id) {
	my $seqlength    = $pos - $offset - length($_) - 1;
	$seqlength      -= int($seqlength/$linelength);
	$offsets->{$id}  = $self->_pack($offset,$seqlength,
					$linelength,$firstline,
					$type,$base);
      }
      $id = ref($self->{makeid}) eq 'CODE' ? $self->{makeid}->($_) : $1;
      ($offset,$firstline,$linelength) = ($pos,length($_),0);
    } else {
      $linelength ||= length($_);
      $type       ||= $self->_type($_);
    }
    }
  # deal with last entry
  if ($id) {
    my $pos = tell($fh);

#    my $seqlength   = $pos - $offset - length($_) - 1;
    # $_ is always null should not be part of this calculation
    my $seqlength   = $pos - $offset  - 1;

    if ($linelength == 0) { # yet another pesky empty chr_random.fa file
      $seqlength = 0;
    } else {
      $seqlength -= int($seqlength/$linelength);
    };
    $offsets->{$id} = $self->_pack($offset,$seqlength,
				   $linelength,$firstline,
				   $type,$base);
  }
  return\% offsets;

  my $self = shift;
  my $id   = shift;
  my $a    = shift()-1;
  my ($offset,$seqlength,$linelength,$firstline,$type,$file) = $self->_unpack($self->{offsets}{$id});
  $a = 0            if $a < 0;
  $a = $seqlength-1 if $a >= $seqlength;
  $offset + $linelength * int($a/($linelength-1)) + $a % ($linelength-1);

  my $self = shift;
  my $args = $self->{dbmargs} or return;
  return ref($args) eq 'ARRAY' ? @$args : $args;

     my $self = shift;
    return '';
}

#-------------------------------------------------------------
# stream-based access to the database
#
package Bio::DB::Fasta::Stream;
use Tie::Handle;
use vars qw(@ISA);
@ISA = qw(Tie::Handle);
eval {
  require Bio::DB::SeqI;
} && (push @ISA,'Bio::DB::SeqI');

  my $self = shift;
  return $self->{id};

  my $self = shift;
  my $id   = shift;
  my $file = $self->file($id) or return;
  $self->fhcache("$self->{dirname}/$file") or $self->throw( "Can't open file $file");

  my $self = shift;
  my $path = shift;
  if (!$self->{fhcache}{$path}) {
    if ($self->{curopen} >= $self->{maxopen}) {
      my @lru = sort {$self->{cacheseq}{$a} <=> $self->{cacheseq}{$b};} keys %{$self->{fhcache}};
      splice(@lru, $self->{maxopen} / 3);
      $self->{curopen} -= @lru;
      for (@lru) { delete $self->{fhcache}{$_} }
    }
    $self->{fhcache}{$path} = IO::File->new($path) or return;
    $self->{curopen}++;
  }
  $self->{cacheseq}{$path}++;
  $self->{fhcache}{$path}

  my $self = shift;
  my $id   = shift;
  my $offset = $self->{offsets}{$id} or return;
  $self->fileno2path(($self->_unpack($offset))[5]);

  my $self = shift;
  my $no   = shift;
  return $self->{offsets}{"__file_$no"};

  my $self = shift;
  return Bio::DB::Fasta::Stream->new($self);

  my $self = shift;
  my $id   = shift;
  return Bio::PrimarySeq::Fasta->new($self,$id);

 grep {!/^__/} keys %{shift->{offsets}}

  my $self = shift;
  my $id   = shift;
  my ($offset,$seqlength,$linelength,$firstline,$type,$file) 
    = $self->_unpack($self->{offsets}{$id}) or return;
  $offset -= $firstline;
  my $data;
  my $fh = $self->fh($id) or return;
  seek($fh,$offset,0);
  read($fh,$data,$firstline);
  chomp $data;
  substr($data,0,1) = '';
  $data;

    my $self = shift;
    my $id   = shift;
    return unless $self->{offsets}{$id};
    return $self->offset($id) - $self->headerlen($id);

  my $self = shift;
  my $id   = shift;
  my $offset = $self->{offsets}{$id} or return;
  ($self->_unpack($offset))[3];

  my $self = shift;
  my $dir  = shift;
  my $force_reindex = shift;

  # find all fasta files
  my @files = glob("$dir/$self->{glob}");
  $self->throw( "no fasta files in $dir") unless @files;

  # get name of index
  my $index = $self->index_name($dir,1);

  # if caller has requested reindexing, then unlink
  # the index file.
  unlink $index if $force_reindex;

  # get the modification time of the index
  my $indextime = (stat($index))[9] || 0;

  # get the most recent modification time of any of the contents
  my $modtime = 0;
  my %modtime;
  foreach (@files) {
    my $m = (stat($_))[9];
    $modtime{$_} = $m;
    $modtime = $m if $modtime < $m;
  }

  my $reindex = $force_reindex || $indextime < $modtime;
  my $offsets = $self->_open_index($index,$reindex) or return;
  $self->{offsets} = $offsets;

  # no indexing needed
  return $offsets unless $reindex;

  # otherwise reindex contents of changed files
  $self->{indexing} = $index;
  foreach (@files) {
    next if( defined $indextime && $modtime{$_} <= $indextime);
    $self->calculate_offsets($_,$offsets);
  }
  delete $self->{indexing};
  return $self->{offsets};

  my $self = shift;
  my $file = shift;
  my $force_reindex = shift;

  my $index = $self->index_name($file);
  # if caller has requested reindexing, then unlink the index
  unlink $index if $force_reindex;

  # get the modification time of the index
  my $indextime = (stat($index))[9];
  my $modtime   = (stat($file))[9];

  my $reindex = $force_reindex || $indextime < $modtime;
  my $offsets = $self->_open_index($index,$reindex) or return;
  $self->{offsets} = $offsets;

  return $self->{offsets} unless $reindex;

  $self->{indexing} = $index;
  $self->calculate_offsets($file,$offsets);
  delete $self->{indexing};
  return $self->{offsets};

  my $self  = shift;
  my ($path,$isdir) = @_;
  unless ($path) {
    my $dir = $self->{dirname} or return;
    return $self->index_name($dir,-d $dir);
  } 
  return "$path/directory.index" if $isdir;
  return "$path.index";

  my $self = shift;
  my $id   = shift;
  my $offset = $self->{offsets}{$id} or return;
  ($self->_unpack($offset))[1];

  my $self = shift;
  return $self->{db}->length($self->{id});

  my $self = shift;
  my $id   = shift;
  my $offset = $self->{offsets}{$id} or return;
  ($self->_unpack($offset))[2];

  my $class = shift;
  my $path  = shift;
  my %opts  = @_;

  my $self = bless { debug      => $opts{-debug},
		     makeid     => $opts{-makeid},
		     glob       => $opts{-glob}    || '*.{fa,fasta,FA,FASTA,fast,FAST,dna,fsa}',
		     maxopen    => $opts{-maxfh}   || 32,
		     dbmargs    => $opts{-dbmargs} || undef,
		     fhcache    => {},
		     cacheseq   => {},
		     curopen    => 0,
		     openseq    => 1,
		     dirname    => undef,
		     offsets    => undef,
		   }, $class;
  my ($offsets,$dirname);

  if (-d $path) {
    $offsets = $self->index_dir($path,$opts{-reindex});
    $dirname = $path;
  } elsif (-f _) {
    $offsets = $self->index_file($path,$opts{-reindex});
    $dirname = dirname($path);
  } else {
    $self->throw( "$path: Invalid file or dirname");
  }
  @{$self}{qw(dirname offsets)} = ($dirname,$offsets);

  $self;

  my $class = shift;
  $class = ref($class) if ref $class;
  my ($db,$id,$start,$stop) = @_;
  return bless { db    => $db,
		 id    => $id,
		 start => $start || 1,
		 stop  => $stop  || $db->length($id)
	       },$class;

  my $class = shift;
  my $db    = shift;
  my $key = $db->FIRSTKEY;
  return bless { db=>$db,key=>$key },$class;

  my $class = shift;
  my $self  = $class->new(@_);
  require Symbol;
  my $fh = Symbol::gensym or return;
  tie $$fh,'Bio::DB::Fasta::Stream',$self or return;
  $fh;

  my $self = shift;
  my ($key,$db) = @{$self}{'key','db'};
  my $value = $db->get_Seq_by_id($key);
  $self->{key} = $db->NEXTKEY($key);
  $value;

  my $self = shift;
  my $id   = shift;
  my $offset = $self->{offsets}{$id} or return;
  ($self->_unpack($offset))[0];

  my $self = shift;
  my $path = shift;
  if ( !defined $self->{offsets}{"__path_$path"} ) {
    my $fileno  = ($self->{offsets}{"__path_$path"} = 0+ $self->{fileno}++);
    $self->{offsets}{"__file_$fileno"} = $path;
  }
  return $self->{offsets}{"__path_$path"}

  my $self = shift;
  return overload::StrVal($self);

  my $self = shift;
  return $self->new(@{$self}{'db','id','stop','start'});

  my $self = shift;
  return $self->{db}->seq($self->{id},$self->{start},$self->{stop});

  my ($self,$id,$start,$stop) = @_;
  if ($id =~ /^(.+):([\d_]+)[,-]([\d_]+)$/) {
    ($id,$start,$stop) = ($1,$2,$3);
    $start =~ s/_//g;
    $stop =~ s/_//g;
  }
  $start ||= 1;
  $stop  ||= $self->length($id);

  my $reversed;
  if ($start > $stop) {
    ($start,$stop) = ($stop,$start);
    $reversed++;
  }

  my $data;

  my $fh = $self->fh($id) or return;
  my $filestart = $self->caloffset($id,$start);
  my $filestop  = $self->caloffset($id,$stop);

  seek($fh,$filestart,0);
  read($fh,$data,$filestop-$filestart+1);
  $data =~ s/\n//g;
  if ($reversed) {
    $data = reverse $data;
    $data =~ tr/gatcGATC/ctagCTAG/;
  }
  $data;

  my $self = shift;
  my ($start,$stop) = @_;
  $self->throw("Stop cannot be smaller than start")  unless $start <= $stop;
  return $self->{start} <= $self->{stop} ?  $self->new($self->{db},
						       $self->{id},
						       $self->{start}+$start-1,
						       $self->{start}+$stop-1)
                                         :  $self->new($self->{db},
						       $self->{id},
						       $self->{start}-($start-1),
						       $self->{start}-($stop-1)
						      );

 Option Name   Description               Default
 -----------   -----------               -------

 -glob         Glob expression to use    *.{fa,fasta,fast,FA,FASTA,FAST,dna}
               for searching for Fasta
	       files in directories. 

 -makeid       A code subroutine for     None
	       transforming Fasta IDs.

 -maxopen      Maximum size of		 32
	       filehandle cache.

 -debug        Turn on status		 0
	       messages.

 -reindex      Force the index to be     0
               rebuilt.

 -dbmargs      Additional arguments      none
               to pass to the DBM
               routines when tied
               (scalar or array ref).

  $db = Bio::DB::Fasta->new("file.fa",-makeid=>\&make_my_id);

  sub make_my_id {
    my $description_line = shift;
    # get a new id from the fasta header
    return $new_id;
  }

 >A12345.3 Predicted C. elegans protein egl-2

   $db->seq("$id:$start,$stop")

  my $seq = $db{CHROMOSOME_I};

  my $seq = $db{'CHROMOSOME_I:1,1000'};

 while (my ($id,$sequence) = each %db) {
    print "$id => $sequence\n";
 }

 while (my $id = each %db) {
    print "$id => $db{$sequence:1,100}\n";
    print "$id => ",tied(%db)->length($id),"\n";
 }

 my $db = Bio::DB::Fasta->new('/path/to/fasta/files');
 for (my $id=$db->FIRSTKEY; $id; $id=$db->NEXTKEY($id)) {
    # do something with sequence
 }

 my $fh = Bio::DB::Fasta->newFh('/path/to/fasta/files');
 while (my $seq = <$fh>) {
   print $seq->id,' => ',$seq->length,"\n";
 }

 Title   : new
 Usage   : my $db = new Bio::DB::Fasta( $path, @options);
 Function: initialize a new Bio::DB::Fasta object
 Returns : new Bio::DB::Fasta object
 Args    : path to dir of fasta files or a single filename

 -glob         Glob expression to use    *.{fa,fasta,fast,FA,FASTA,FAST}
               for searching for Fasta
	       files in directories. 

 -makeid       A code subroutine for     None
	       transforming Fasta IDs.

 -maxopen      Maximum size of		 32
	       filehandle cache.

 -debug        Turn on status		 0
	       messages.

 -reindex      Force the index to be     0
               rebuilt.

 -dbmargs      Additional arguments      none
               to pass to the DBM
               routines when tied
               (scalar or array ref).

 Title   : subseq
 Usage   : $seqdb->subseq($id,$start,$stop);
 Function: returns a subseq of a sequence in the db
 Returns : subsequence data
 Args    : id of sequence, starting point, ending point