BaseFeatureAdaptor.pm 30.2 KB
Newer Older
1
#
Simon Potter's avatar
Simon Potter committed
2
# EnsEMBL module for Bio::EnsEMBL::DBSQL::BaseFeatureAdaptor
3
#
Graham McVicker's avatar
Graham McVicker committed
4
# Copyright (c) 2003 Ensembl
5 6 7 8 9
#
# You may distribute this module under the same terms as perl itself

=head1 NAME

10
Bio::EnsEMBL::DBSQL::BaseFeatureAdaptor - An Abstract Base class for all
Graham McVicker's avatar
Graham McVicker committed
11
FeatureAdaptors
12 13 14

=head1 SYNOPSIS

15
Abstract class - should not be instantiated.  Implementation of
16 17 18 19 20
abstract methods must be performed by subclasses.

=head1 DESCRIPTION

This is a base adaptor for feature adaptors. This base class is simply a way
21
of eliminating code duplication through the implementation of methods
22 23
common to all feature adaptors.

24 25
=head1 CONTACT

Graham McVicker's avatar
Graham McVicker committed
26
Contact Ensembl development list for info: <ensembl-dev@ebi.ac.uk>
27 28 29

=head1 METHODS

30 31 32
=cut

package Bio::EnsEMBL::DBSQL::BaseFeatureAdaptor;
33
use vars qw(@ISA);
34 35 36
use strict;

use Bio::EnsEMBL::DBSQL::BaseAdaptor;
37
use Bio::EnsEMBL::Utils::Cache;
38 39
use Bio::EnsEMBL::Utils::Exception qw(warning throw deprecate);
use Bio::EnsEMBL::Utils::Argument qw(rearrange);
40 41 42

@ISA = qw(Bio::EnsEMBL::DBSQL::BaseAdaptor);

43 44
our $SLICE_FEATURE_CACHE_SIZE = 4;
our $MAX_SPLIT_QUERY_SEQ_REGIONS = 3;
45 46 47 48 49 50 51 52 53 54 55 56 57 58

=head2 new

  Arg [1]    : list of args @args
               Superclass constructor arguments
  Example    : none
  Description: Constructor which just initializes internal cache structures
  Returntype : Bio::EnsEMBL::BaseFeatureAdaptor
  Exceptions : none
  Caller     : implementing subclass constructors

=cut

sub new {
59
  my $caller = shift;
60

61
  my $class = ref($caller) || $caller;
62

63
  my $self = $class->SUPER::new(@_);
64

65 66 67 68
  #initialize an LRU cache
  my %cache;
  tie(%cache, 'Bio::EnsEMBL::Utils::Cache', $SLICE_FEATURE_CACHE_SIZE);
  $self->{'_slice_feature_cache'} = \%cache;
69 70 71

  return $self;
}
72

73

74

75
# _straight_join
76

77 78 79 80 81 82 83 84 85
#   Arg [1]    : (optional) boolean $new_val
#   Example    : $self->_straight_join(1);
#                $self->generic_fetch($constraint);
#                $self->_straight_join(0);
#   Description: PROTECTED Getter/Setter that turns on/off the use of 
#                a straight join in queries.
#   Returntype : boolean
#   Exceptions : none
#   Caller     : general
86 87 88 89 90 91 92 93 94 95 96

sub _straight_join {
  my $self = shift;
  if(@_) {
    $self->{'_straight_join'} = shift;
  }

  return $self->{'_straight_join'};
}


97 98 99 100
=head2 generic_fetch

  Arg [1]    : (optional) string $constraint
               An SQL query constraint (i.e. part of the WHERE clause)
101 102 103 104 105
  Arg [2]    : (optional) Bio::EnsEMBL::AssemblyMapper $mapper
               A mapper object used to remap features
               as they are retrieved from the database
  Arg [3]    : (optional) Bio::EnsEMBL::Slice $slice
               A slice that features should be remapped to
106
  Example    : $fts = $a->generic_fetch('contig_id in (1234, 1235)', 'Swall');
107 108
  Description: Performs a database fetch and returns feature objects in
               contig coordinates.
109
  Returntype : listref of Bio::EnsEMBL::SeqFeature in contig coordinates
110 111 112 113
  Exceptions : none
  Caller     : BaseFeatureAdaptor, ProxyDnaAlignFeatureAdaptor::generic_fetch

=cut
114

115
sub generic_fetch {
116
  my ($self, $constraint, $mapper, $slice) = @_;
117

118
  my @tabs = $self->_tables;
119
  my $columns = join(', ', $self->_columns());
120

121
  my $db = $self->db();
122 123 124 125 126

  #
  # Construct a left join statement if one was defined, and remove the
  # left-joined table from the table list
  #
127
  my @left_join_list = $self->_left_join();
128 129
  my $left_join = '';
  my @tables;
130 131
  if(@left_join_list) {
    my %left_join_hash = map { $_->[0] => $_->[1] } @left_join_list;
132
    while(my $t = shift @tabs) {
133 134 135
      if( exists $left_join_hash{ $t->[0] } ) {
        my $condition = $left_join_hash{ $t->[0] };
        my $syn = $t->[1];
Web Admin's avatar
Web Admin committed
136
        $left_join .=  "LEFT JOIN\n       ".$t->[0]." $syn ON $condition ";
137
      } else {
138
        push @tables, $t;
139 140 141 142 143
      }
    }
  } else {
    @tables = @tabs;
  }
144

145 146 147 148 149 150 151

  my $straight_join = '';

  if($self->_straight_join()) {
    $straight_join = "STRAIGHT_JOIN";
  }

152 153 154
  #construct a nice table string like 'table1 t1, table2 t2'
  my $tablenames = join(', ', map({ join(' ', @$_) } @tables));

155
  my $sql = "SELECT $straight_join $columns\n  FROM $tablenames $left_join";
156 157 158 159 160

  my $default_where = $self->_default_where_clause;
  my $final_clause = $self->_final_clause;

  #append a where clause if it was defined
161
  if($constraint) {
Web Admin's avatar
Web Admin committed
162
    $sql .= "\n WHERE $constraint ";
163
    if($default_where) {
Web Admin's avatar
Web Admin committed
164
      $sql .= " AND\n       $default_where ";
165 166
    }
  } elsif($default_where) {
Web Admin's avatar
Web Admin committed
167
    $sql .= "\n WHERE $default_where ";
168 169 170
  }

  #append additional clauses which may have been defined
Web Admin's avatar
Web Admin committed
171
  $sql .= "\n$final_clause";
172

173
  my $sth = $db->prepare($sql);
Web Admin's avatar
Web Admin committed
174
     $sth->execute;
175 176 177

  my $res = $self->_objs_from_sth($sth, $mapper, $slice);

178

179
  return $res;
180 181 182 183
}


=head2 fetch_by_dbID
184 185

  Arg [1]    : int $id
186 187 188
               The unique database identifier for the feature to be obtained
  Example    : $feat = $adaptor->fetch_by_dbID(1234));
               $feat = $feat->transform('contig');
189
  Description: Returns the feature created from the database defined by the
190 191 192 193 194 195 196 197 198
               the id $id.  The feature will be returned in its native
               coordinate system.  That is, the coordinate system in which it
               is stored in the database.  In order to convert it to a
               particular coordinate system use the transfer() or transform()
               method.  If the feature is not found in the database then
               undef is returned instead
  Returntype : Bio::EnsEMBL::Feature or undef
  Exceptions : thrown if $id arg is not provided
               does not exist
199
  Caller     : general
200 201 202 203 204 205

=cut

sub fetch_by_dbID{
  my ($self,$id) = @_;

206
  throw("id argument is required") if(!defined $id);
207

208 209
  #construct a constraint like 't1.table1_id = 123'
  my @tabs = $self->_tables;
210 211
  my ($name, $syn) = @{$tabs[0]};
  my $constraint = "${syn}.${name}_id = $id";
212

213 214 215 216 217
  #Should only be one
  my ($feat) = @{$self->generic_fetch($constraint)};

  return undef if(!$feat);

218
  return $feat;
219 220 221
}


222 223 224 225 226 227 228 229 230 231 232 233 234 235 236 237
=head2 fetch_all_by_dbID_list

  Arg [1]    : listref of ints $id_list
               The unique database identifiers for the features to be obtained
  Example    : @feats = @{$adaptor->fetch_by_dbID_list([1234, 2131, 982]))};
  Description: Returns the features created from the database defined by the
               the ids in contained in the id list $id_list.  The features 
               will be returned in their native coordinate system. That is, 
               the coordinate system in which they are stored in the database.
               In order to convert the features to a particular coordinate 
               system use the transfer() or transform() method.  If none of the
               features are found in the database a reference to an empty 
               list is returned.
  Returntype : listref of Bio::EnsEMBL::Features
  Exceptions : thrown if $id arg is not provided
               does not exist
238
  Caller     : general
239 240 241

=cut

242 243
sub fetch_all_by_dbID_list {
  my ($self,$id_list) = @_;
244

245 246
  if(!defined($id_list) || ref($id_list) ne 'ARRAY') {
    throw("id_list list reference argument is required") if(!defined $id_list);
247 248
  }

249
  return [] if(!@$id_list);
250

251 252
  my @out;
  #construct a constraint like 't1.table1_id = 123'
253
  my @tabs = $self->_tables;
254
  my ($name, $syn) = @{$tabs[0]};
255

256 257
  # mysql is faster and we ensure that we do not exceed the max query size by
  # splitting large queries into smaller queries of 200 ids
258
  my $max_size = 200;
259

260 261 262 263 264 265 266
  while(@$id_list) {
    my @ids;
    if(@$id_list > $max_size) {
      @ids = splice(@$id_list, 0, $max_size);
    } else {
      @ids = splice(@$id_list, 0);
    }
267

268 269 270 271 272 273
    my $id_str;
    if(@ids > 1)  {
      $id_str = " IN (" . join(',', @ids). ")";
    } else {
      $id_str = " = " . $ids[0];
    }
274

275
    my $constraint = "${syn}.${name}_id $id_str";
276

277 278
    push @out, @{$self->generic_fetch($constraint)};
  }
Arne Stabenau's avatar
Arne Stabenau committed
279

280
  return \@out;
281 282 283 284
}



285
=head2 fetch_all_by_Slice
286 287 288 289 290

  Arg [1]    : Bio::EnsEMBL::Slice $slice
               the slice from which to obtain features
  Arg [2]    : (optional) string $logic_name
               the logic name of the type of features to obtain
291 292 293 294 295 296
  Example    : $fts = $a->fetch_all_by_Slice($slice, 'Swall');
  Description: Returns a listref of features created from the database 
               which are on the Slice defined by $slice. If $logic_name is 
               defined only features with an analysis of type $logic_name 
               will be returned. 
  Returntype : listref of Bio::EnsEMBL::SeqFeatures in Slice coordinates
297 298 299 300 301
  Exceptions : none
  Caller     : Bio::EnsEMBL::Slice

=cut

302
sub fetch_all_by_Slice {
303
  my ($self, $slice, $logic_name) = @_;
304

305
  #fetch by constraint with empty constraint
306
  return $self->fetch_all_by_Slice_constraint($slice, '', $logic_name);
307 308 309
}


310
=head2 fetch_all_by_Slice_and_score
311 312 313

  Arg [1]    : Bio::EnsEMBL::Slice $slice
               the slice from which to obtain features
314
  Arg [2]    : (optional) float $score
315 316 317
               lower bound of the the score of the features retrieved
  Arg [3]    : (optional) string $logic_name
               the logic name of the type of features to obtain
318
  Example    : $fts = $a->fetch_all_by_Slice($slice, 'Swall');
319 320 321 322 323
  Description: Returns a list of features created from the database which are 
               are on the Slice defined by $slice and which have a score 
               greated than $score. If $logic_name is defined, 
               only features with an analysis of type $logic_name will be 
               returned. 
324
  Returntype : listref of Bio::EnsEMBL::SeqFeatures in Slice coordinates
325 326 327 328 329
  Exceptions : none
  Caller     : Bio::EnsEMBL::Slice

=cut

330
sub fetch_all_by_Slice_and_score {
331 332 333
  my ($self, $slice, $score, $logic_name) = @_;
  my $constraint;

334
  if(defined $score) {
335 336 337 338
    #get the synonym of the primary_table
    my @tabs = $self->_tables;
    my $syn = $tabs[0]->[1];
    $constraint = "${syn}.score > $score";
339
  }
340 341
  return $self->fetch_all_by_Slice_constraint($slice, $constraint, 
					      $logic_name);
342
}
343 344


345
=head2 fetch_all_by_Slice_constraint
346

347 348 349 350 351
  Arg [1]    : Bio::EnsEMBL::Slice $slice
               the slice from which to obtain features
  Arg [2]    : (optional) string $constraint
               An SQL query constraint (i.e. part of the WHERE clause)
  Arg [3]    : (optional) string $logic_name
352
               the logic name of the type of features to obtain
353 354
  Example    : $fs = $a->fetch_all_by_Slice_constraint($slc, 'perc_ident > 5');
  Description: Returns a listref of features created from the database which 
355 356 357 358
               are on the Slice defined by $slice and fulfill the SQL 
               constraint defined by $constraint. If logic name is defined, 
               only features with an analysis of type $logic_name will be 
               returned. 
359
  Returntype : listref of Bio::EnsEMBL::SeqFeatures in Slice coordinates
360 361
  Exceptions : thrown if $slice is not defined
  Caller     : Bio::EnsEMBL::Slice
362 363 364

=cut

365
sub fetch_all_by_Slice_constraint {
366 367 368
  my($self, $slice, $constraint, $logic_name) = @_;

  my @result;
369

370
  if(!ref($slice) || !$slice->isa("Bio::EnsEMBL::Slice")) {
Graham McVicker's avatar
Graham McVicker committed
371
    throw("Bio::EnsEMBL::Slice argument expected.");
372 373
  }

374 375
  $constraint ||= '';
  $constraint = $self->_logic_name_to_constraint($constraint, $logic_name);
376 377

  #if the logic name was invalid, undef was returned
378
  return [] if(!defined($constraint));
379

380
  #check the cache and return if we have already done this query
381
  my $key = uc(join(':', $slice->name, $constraint));
382 383 384 385 386

  if(exists($self->{'_slice_feature_cache'}->{$key})) {
    return $self->{'_slice_feature_cache'}->{$key};
  }

387
  my $sa = $slice->adaptor();
388

389 390
  # Hap/PAR support: retrieve normalized 'non-symlinked' slices
  my @proj = @{$sa->fetch_normalized_slice_projection($slice)};
391

392
  if(@proj == 0) {
393 394 395 396
    throw('Could not retrieve normalized Slices. Database contains ' .
          'incorrect assembly_exception information.');
  }

397 398
  # Want to get features on the FULL original slice
  # as well as any symlinked slices
399

400 401
  # Filter out partial slices from projection that are on
  # same seq_region as original slice
402

403
  my $sr_id = $slice->get_seq_region_id();
404

405
  @proj = grep { $_->to_Slice->get_seq_region_id() != $sr_id } @proj;
406

407 408 409
  my $segment = bless([1,$slice->length(),$slice ],
                      'Bio::EnsEMBL::ProjectionSegment');
  push( @proj, $segment );
410

411 412 413 414
  # fetch features for the primary slice AND all symlinked slices
  foreach my $segment (@proj) {
    my $offset = $segment->from_start();
    my $seg_slice  = $segment->to_Slice();
415

416
    my $features = $self->_slice_fetch($seg_slice, $constraint);
417

418 419 420 421 422 423 424 425 426 427 428 429 430 431 432 433
    # if this was a symlinked slice offset the feature coordinates as needed
    if($seg_slice->name() ne $slice->name()) {
      foreach my $f (@$features) {
        if($offset != 1) {
          $f->{'start'} += $offset-1;
          $f->{'end'}   += $offset-1;
        }
        $f->{'slice'} = $slice;
        push @result, $f;
      }
    } else {
      push @result, @$features;
    }
  }

  $self->{'_slice_feature_cache'}->{$key} = \@result;
434

435 436
  return \@result;
}
437 438


439 440 441 442 443 444 445 446 447 448 449 450 451 452 453 454 455 456 457 458 459 460 461 462 463 464 465 466 467 468 469 470 471 472 473 474 475 476 477
#
# helper function used by fetch_all_by_Slice_constraint method
#
sub _slice_fetch {
  my $self = shift;
  my $slice = shift;
  my $orig_constraint = shift;

  my $slice_start  = $slice->start();
  my $slice_end    = $slice->end();
  my $slice_strand = $slice->strand();
  my $slice_cs     = $slice->coord_system();
  my $slice_seq_region = $slice->seq_region_name();

  #get the synonym and name of the primary_table
  my @tabs = $self->_tables;
  my ($tab_name, $tab_syn) = @{$tabs[0]};

  #find out what coordinate systems the features are in
  my $mcc = $self->db->get_MetaCoordContainer();
  my @feat_css = @{$mcc->fetch_all_CoordSystems_by_feature_type($tab_name)};

  my $asma = $self->db->get_AssemblyMapperAdaptor();
  my @features;

  # fetch the features from each coordinate system they are stored in
 COORD_SYSTEM: foreach my $feat_cs (@feat_css) {
    my $mapper;
    my @coords;
    my @ids;

    if($feat_cs->equals($slice_cs)) {
      # no mapping is required if this is the same coord system
      my $constraint = $orig_constraint;

      my $sr_id = $self->db->get_SliceAdaptor->get_seq_region_id($slice);
      $constraint .= " AND " if($constraint);
      $constraint .=
          "${tab_syn}.seq_region_id = $sr_id AND " .
478 479
          "${tab_syn}.seq_region_start <= $slice_end AND " .
          "${tab_syn}.seq_region_end >= $slice_start";
480
      my $fs = $self->generic_fetch($constraint,undef,$slice);
481

482 483 484
      # features may still have to have coordinates made relative to slice
      # start
      $fs = $self->_remap($fs, $mapper, $slice);
485

486 487 488
      push @features, @$fs;
    } else {
      $mapper = $asma->fetch_by_CoordSystems($slice_cs, $feat_cs);
489

490 491 492 493
      # Get list of coordinates and corresponding internal ids for
      # regions the slice spans
      @coords = $mapper->map($slice_seq_region, $slice_start, $slice_end,
                             $slice_strand, $slice_cs);
494

495
      @coords = grep {!$_->isa('Bio::EnsEMBL::Mapper::Gap')} @coords;
496

497 498 499 500 501 502 503 504 505 506 507 508 509 510 511 512 513
      next COORD_SYSTEM if(!@coords);

      @ids = map {$_->id()} @coords;
      @ids = @{$asma->seq_regions_to_ids($feat_cs, \@ids)};

      # When regions are large and only partially spanned by slice
      # it is faster to to limit the query with start and end constraints.
      # Take simple approach: use regional constraints if there are less
      # than a specific number of regions covered.

      if(@coords > $MAX_SPLIT_QUERY_SEQ_REGIONS) {
        my $constraint = $orig_constraint;
        my $id_str = join(',', @ids);
        $constraint .= " AND " if($constraint);
        $constraint .= "${tab_syn}.seq_region_id IN ($id_str)";

        my $fs = $self->generic_fetch($constraint, $mapper, $slice);
514

515 516 517 518 519 520 521 522 523 524 525
        $fs = $self->_remap($fs, $mapper, $slice);

        push @features, @$fs;

      } else {
        # do multiple split queries using start / end constraints
        my $len = @coords;
        for(my $i = 0; $i < $len; $i++) {
          my $constraint = $orig_constraint;
          $constraint .= " AND " if($constraint);
          $constraint .=
526 527 528
              "${tab_syn}.seq_region_id = "     . $ids[$i] . " AND " .
              "${tab_syn}.seq_region_start <= " . $coords[$i]->end() . " AND ".
              "${tab_syn}.seq_region_end >= "   . $coords[$i]->start();
529
          my $fs = $self->generic_fetch($constraint,$mapper,$slice);
530

531
          $fs = $self->_remap($fs, $mapper, $slice);
532

533
          push @features, @$fs;
534 535 536
        }
      }
    }
537
  } #COORD system loop
538

539
  return \@features;
540
}
541

542

543

544 545 546 547 548 549 550 551 552 553 554 555 556 557 558 559 560 561
#
# Helper function containing some common feature storing functionality
#
# Given a Feature this will return a copy (or the same feature if no changes 
# to the feature are needed) of the feature which is relative to the start
# of the seq_region it is on. The seq_region_id of the seq_region it is on
# is also returned.
#
# This method will also ensure that the database knows which coordinate
# systems that this feature is stored in.
#

sub _pre_store {
  my $self    = shift;
  my $feature = shift;

  if(!ref($feature) || !$feature->isa('Bio::EnsEMBL::Feature')) {
    throw('Expected Feature argument.');
562
  }
563

564 565 566
  $self->_check_start_end_strand($feature->start(),$feature->end(),
                                 $feature->strand());

567 568 569 570

  my $db = $self->db();

  my $slice_adaptor = $db->get_SliceAdaptor();
571 572 573 574 575 576
  my $slice = $feature->slice();

  if(!ref($slice) || !$slice->isa('Bio::EnsEMBL::Slice')) {
    throw('Feature must be attached to Slice to be stored.');
  }

577
  # make sure feature coords are relative to start of entire seq_region
578
  if($slice->start != 1 || $slice->strand != 1) {
579
    #move feature onto a slice of the entire seq_region
580 581 582 583 584 585 586 587 588 589 590 591 592 593 594
    $slice = $slice_adaptor->fetch_by_region($slice->coord_system->name(),
                                             $slice->seq_region_name(),
                                             undef, #start
                                             undef, #end
                                             undef, #strand
                                             $slice->coord_system->version());

    $feature = $feature->transfer($slice);

    if(!$feature) {
      throw('Could not transfer Feature to slice of ' .
            'entire seq_region prior to storing');
    }
  }

595
  # Ensure this type of feature is known to be stored in this coord system.
596 597 598 599 600
  my $cs = $slice->coord_system;

  my ($tab) = $self->_tables();
  my $tabname = $tab->[0];

601 602 603
  my $mcc = $db->get_MetaCoordContainer();

  $mcc->add_feature_type($cs, $tabname);
604 605 606 607 608 609 610 611 612 613 614 615 616 617 618 619 620 621 622 623 624 625 626 627 628 629 630 631 632 633 634 635 636 637 638 639 640 641 642 643 644 645 646 647 648 649 650 651 652 653 654 655 656 657

  my $seq_region_id = $slice_adaptor->get_seq_region_id($slice);

  if(!$seq_region_id) {
    throw('Feature is associated with seq_region which is not in this DB.');
  }

  return ($feature, $seq_region_id);
}


#
# helper function used to validate start/end/strand and 
# hstart/hend/hstrand etc.
#
sub _check_start_end_strand {
  my $self = shift;
  my $start = shift;
  my $end   = shift;
  my $strand = shift;

  #
  # Make sure that the start, end, strand are valid
  #
  if(int($start) != $start) {
    throw("Invalid Feature start [$start].  Must be integer.");
  }
  if(int($end) != $end) {
    throw("Invalid Feature end [$end]. Must be integer.");
  }
  if(int($strand) != $strand || $strand < -1 || $strand > 1) {
    throw("Invalid Feature strand [$strand]. Must be -1, 0 or 1.");
  }
  if($end < $start) {
    throw("Invalid Feature start/end [$start/$end]. Start must be less " .
          "than or equal to end.");
  }

  return 1;
}


#
# Given a list of features checks if they are in the correct coord system
# by looking at the first features slice.  If they are not then they are
# converted and placed on the slice.
#
sub _remap {
  my ($self, $features, $mapper, $slice) = @_;

  #check if any remapping is actually needed
  if(@$features && (!$features->[0]->isa('Bio::EnsEMBL::Feature') ||
                    $features->[0]->slice == $slice)) {
    return $features;
658
  }
659

660
  #remapping has not been done, we have to do our own conversion from
661 662 663 664 665 666 667 668 669 670
  #to slice coords

  my @out;

  my $slice_start = $slice->start();
  my $slice_end   = $slice->end();
  my $slice_strand = $slice->strand();
  my $slice_cs    = $slice->coord_system();

  my ($seq_region, $start, $end, $strand);
671

672
  my $slice_seq_region = $slice->seq_region_name();
673

674
  foreach my $f (@$features) {
675
    #since feats were obtained in contig coords, attached seq is a contig
676 677 678 679 680 681
    my $fslice = $f->slice();
    if(!$fslice) {
      throw("Feature does not have attached slice.\n");
    }
    my $fseq_region = $fslice->seq_region_name();
    my $fcs = $fslice->coord_system();
682

683 684
    if(!$slice_cs->equals($fcs)) {
      #slice of feature in different coord system, mapping required
685

686 687 688 689 690 691 692 693 694 695 696 697 698
      ($seq_region, $start, $end, $strand) =
        $mapper->fastmap($fseq_region,$f->start(),$f->end(),$f->strand(),$fcs);

      # undefined start means gap
      next if(!defined $start);
    } else {
      $start      = $f->start();
      $end        = $f->end();
      $strand     = $f->strand();
      $seq_region = $f->slice->seq_region_name();
    }

    # maps to region outside desired area
699 700
    next if ($start > $slice_end) || ($end < $slice_start) || 
      ($slice_seq_region ne $seq_region);
701

702
    #shift the feature start, end and strand in one call
703
    if($slice_strand == -1) {
704
      $f->move( $slice_end - $end + 1, $slice_end - $start + 1, $strand * -1 );
705 706 707
    } else {
      $f->move( $start - $slice_start + 1, $end - $slice_start + 1, $strand );
    }
708 709 710

    $f->slice($slice);

Web Admin's avatar
Web Admin committed
711
    push @out,$f;
712
  }
713 714 715 716 717 718 719 720 721 722 723 724 725 726 727 728 729 730 731 732 733 734 735 736 737 738 739 740 741 742 743 744 745 746 747 748 749 750 751 752 753 754 755 756 757 758 759 760 761 762

  return \@out;
}


#
# Given a logic name and an existing constraint this will
# add an analysis table constraint to the feature.  Note that if no
# analysis_id exists in the columns of the primary table then no
# constraint is added at all
#
sub _logic_name_to_constraint {
  my $self = shift;
  my $constraint = shift;
  my $logic_name = shift;

  return $constraint if(!$logic_name);

  #make sure that an analysis_id exists in the primary table
  my ($prim_tab) = $self->_tables();
  my $prim_synonym = $prim_tab->[1];

  my $found_analysis=0;
  foreach my $col ($self->_columns) {
    my ($syn,$col_name) = split(/\./,$col);
    next if($syn ne $prim_synonym);
    if($col_name eq 'analysis_id') {
      $found_analysis = 1;
      last;
    }
  }

  if(!$found_analysis) {
    warning("This feature is not associated with an analysis.\n" .
            "Ignoring logic_name argument = [$logic_name].\n");
    return $constraint;
  }

  my $aa = $self->db->get_AnalysisAdaptor();
  my $an = $aa->fetch_by_logic_name($logic_name);

  if(!$an) {
    return undef;
  }

  my $an_id = $an->dbID();

  $constraint .= ' AND' if($constraint);
  $constraint .= " ${prim_synonym}.analysis_id = $an_id";
  return $constraint;
763 764 765
}


766
=head2 store
767

768
  Arg [1]    : list of Bio::EnsEMBL::SeqFeature
769 770 771 772 773 774 775
  Example    : $adaptor->store(@feats);
  Description: ABSTRACT  Subclasses are responsible for implementing this 
               method.  It should take a list of features and store them in 
               the database.
  Returntype : none
  Exceptions : thrown method is not implemented by subclass
  Caller     : general
776 777 778 779 780 781

=cut

sub store{
  my $self = @_;

782
  throw("Abstract method store not defined by implementing subclass\n");
783 784 785
}


786 787 788 789 790 791 792 793 794
=head2 remove

  Arg [1]    : A feature $feature 
  Example    : $feature_adaptor->remove($feature);
  Description: This removes a feature from the database.  The table the
               feature is removed from is defined by the abstract method
               _tablename, and the primary key of the table is assumed
               to be _tablename() . '_id'.  The feature argument must 
               be an object implementing the dbID method, and for the
795
               feature to be removed from the database a dbID value must
796 797
               be returned.
  Returntype : none
798 799
  Exceptions : thrown if $feature arg does not implement dbID(), or if
               $feature->dbID is not a true value
800 801 802 803
  Caller     : general

=cut

804

805 806 807
sub remove {
  my ($self, $feature) = @_;

808 809
  if(!$feature || !ref($feature) || !$feature->isa('Bio::EnsEMBL::Feature')) {
    throw('Feature argument is required');
810 811
  }

812 813
  if(!$feature->is_stored($self->db)) {
    throw("This feature is not stored in this database");
814 815
  }

816 817
  my @tabs = $self->_tables;
  my ($table) = @{$tabs[0]};
818 819 820 821

  my $sth = $self->prepare("DELETE FROM $table WHERE ${table}_id = ?");
  $sth->execute($feature->dbID());

822 823 824 825
  #unset the feature dbID ad adaptor
  $feature->dbID(undef);
  $feature->adaptor(undef);

826 827 828 829
  return;
}


830
=head2 remove_by_Slice
831

832 833 834 835 836 837 838 839
  Arg [1]    : Bio::Ensembl::Slice $slice
  Example    : $feature_adaptor->remove_by_RawContig($slice);
  Description: This removes features from the database which lie on a region
               represented by the passed in slice.  Only features which are
               fully contained by the slice are deleted; features which overlap
               the edge of the slice are not removed.
               The table the features are removed from is defined by
               the abstract method_tablename.
840
  Returntype : none
841
  Exceptions : thrown if no slice is supplied
842 843 844 845
  Caller     : general

=cut

846 847
sub remove_by_Slice {
  my ($self, $slice) = @_;
848

849 850
  if(!$slice || !ref($slice) || !$slice->isa('Bio::EnsEMBL::Slice')) {
    throw("Slice argument is required");
851 852
  }

853 854
  my @tabs = $self->_tables;
  my ($table_name) = @{$tabs[0]};
855

856 857 858
  my $seq_region_id = $self->db->get_SliceAdaptor->get_seq_region_id($slice);
  my $start = $slice->start();
  my $end   = $slice->end();
859

860 861 862 863 864 865 866
  #
  # Delete only features fully on the slice, not overlapping ones
  #
  my $sth = $self->prepare("DELETE FROM $table_name " .
                           "WHERE seq_region_id = ? " .
                           "AND   seq_region_start >= ? " .
                           "AND   seq_region_end <= ?");
867

868 869
  $sth->execute($seq_region_id, $start, $end);
  $sth->finish();
870 871 872 873
}



874 875


876 877 878 879 880 881 882 883 884 885 886 887 888 889 890 891
#_tables
#
#  Args       : none
#  Example    : $tablename = $self->_table_name()
#  Description: ABSTRACT PROTECTED Subclasses are responsible for implementing
#               this method.  It should list of [tablename, alias] pairs.  
#               Additionally the primary table (with the dbID, analysis_id, and
#               score) should be the first table in the list.
#               e.g:
#               ( ['repeat_feature',   'rf'],
#                 ['repeat_consensus', 'rc']);
#               used to obtain features.  
#  Returntype : list of [tablename, alias] pairs
#  Exceptions : thrown if not implemented by subclass
#  Caller     : BaseFeatureAdaptor::generic_fetch
#
892

893
sub _tables {
894 895
  my $self = shift;

896
  throw("abstract method _tables not defined by implementing" .
Simon Potter's avatar
Simon Potter committed
897
               " subclass of BaseFeatureAdaptor");
898 899 900
  return undef;
}

901

902 903 904 905 906 907 908 909 910 911 912
#_columns
#
#  Args       : none
#  Example    : $tablename = $self->_columns()
#  Description: ABSTRACT PROTECTED Subclasses are responsible for implementing
#               this method.  It should return a list of columns to be used
#               for feature creation
#  Returntype : list of strings
#  Exceptions : thrown if not implemented by subclass
#  Caller     : BaseFeatureAdaptor::generic_fetch
#
913 914 915 916

sub _columns {
  my $self = shift;

917
  throw("abstract method _columns not defined by implementing" .
Simon Potter's avatar
Simon Potter committed
918
               " subclass of BaseFeatureAdaptor");
919 920 921
}


922

923 924 925 926 927 928 929 930 931 932 933 934
# _default_where_clause
#
#  Arg [1]    : none
#  Example    : none
#  Description: May be overridden to provide an additional where constraint to 
#               the SQL query which is generated to fetch feature records.
#               This constraint is always appended to the end of the generated
#               where clause
#  Returntype : string
#  Exceptions : none
#  Caller     : generic_fetch
#
935

936 937
sub _default_where_clause {
  my $self = shift;
938

939 940
  return '';
}
941 942 943



944
# _left_join
945

946 947 948 949 950 951 952 953 954
#  Arg [1]    : none
#  Example    : none
#  Description: Can be overridden by a subclass to specify any left joins
#               which should occur. The table name specigfied in the join
#               must still be present in the return values of
#  Returntype : a {'tablename' => 'join condition'} pair
#  Exceptions : none
#  Caller     : general
#
955

956 957
sub _left_join {
  my $self = shift;
958

959 960
  return ();
}
961 962


963 964 965 966 967 968 969 970 971 972 973 974 975 976

#_final_clause

#  Arg [1]    : none
#  Example    : none
#  Description: May be overriden to provide an additional clause to the end
#               of the SQL query used to fetch feature records.  
#               This is useful to add a required ORDER BY clause to the 
#               query for example.
#  Returntype : string
#  Exceptions : none
#  Caller     : generic_fetch

sub _final_clause {
977 978 979 980 981 982
  my $self = shift;

  return '';
}


983
#_objs_from_sth
984

985 986 987 988 989 990 991 992 993
#  Arg [1]    : DBI::row_hashref $hashref containing key-value pairs 
#               for each of the columns specified by the _columns method
#  Example    : my @feats = $self->_obj_from_hashref
#  Description: ABSTRACT PROTECTED The subclass is responsible for implementing
#               this method.  It should take in a DBI row hash reference and
#               return a list of created features in contig coordinates.
#  Returntype : list of Bio::EnsEMBL::*Features in contig coordinates
#  Exceptions : thrown if not implemented by subclass
#  Caller     : BaseFeatureAdaptor::generic_fetch
994

995 996
sub _objs_from_sth {
  my $self = shift;
997

998 999 1000
  throw("abstract method _obj_from_sth not defined by implementing"
             . " subclass of BaseFeatureAdaptor");
}
1001

1002 1003 1004 1005 1006 1007 1008 1009 1010