DnaAlignFeatureAdaptor.pm 21.2 KB
Newer Older
1
2
=head1 LICENSE

3
  Copyright (c) 1999-2013 The European Bioinformatics Institute and
4
5
6
7
8
9
10
11
12
13
  Genome Research Limited.  All rights reserved.

  This software is distributed under a modified Apache license.
  For license details, please see

    http://www.ensembl.org/info/about/code_licence.html

=head1 CONTACT

  Please email comments or questions to the public Ensembl
14
  developers list at <dev@ensembl.org>.
15
16
17
18
19

  Questions may also be sent to the Ensembl help desk at
  <helpdesk@ensembl.org>.

=cut
20
21
22
23
24
25
26

=head1 NAME

Bio::EnsEMBL::DBSQL::DnaAlignFeatureAdaptor - Adaptor for DnaAlignFeatures

=head1 SYNOPSIS

27
  $dafa = $registry->get_adaptor( 'Human', 'Core', 'DnaAlignFeature' );
28

29
  my @features = @{ $dafa->fetch_all_by_Slice($slice) };
30

31
  $dafa->store(@features);
32
33
34

=head1 DESCRIPTION

35
36
37
This is an adaptor responsible for the retrieval and storage of
DnaDnaAlignFeatures from the database. This adaptor inherits most of its
functionality from the BaseAlignFeatureAdaptor and BaseFeatureAdaptor
38
superclasses.
39

40
=head1 METHODS
41
42
43
44
45
46
47

=cut


package Bio::EnsEMBL::DBSQL::DnaAlignFeatureAdaptor;
use vars qw(@ISA);
use strict;
48
use Bio::EnsEMBL::DnaDnaAlignFeature;
49
use Bio::EnsEMBL::DBSQL::BaseAlignFeatureAdaptor;
Laura Clarke's avatar
 
Laura Clarke committed
50
use Bio::EnsEMBL::Utils::Exception qw(throw warning);
51

52
@ISA = qw(Bio::EnsEMBL::DBSQL::BaseAlignFeatureAdaptor);
53
54


55
=head2 _tables
Graham McVicker's avatar
Graham McVicker committed
56
57

  Args       : none
58
  Example    : @tabs = $self->_tables
Graham McVicker's avatar
Graham McVicker committed
59
  Description: PROTECTED implementation of the abstract method inherited from
60
61
               BaseFeatureAdaptor.  Returns list of [tablename, alias] pairs
  Returntype : list of listrefs of strings
Graham McVicker's avatar
Graham McVicker committed
62
63
  Exceptions : none
  Caller     : Bio::EnsEMBL::DBSQL::BaseFeatureAdaptor::generic_fetch
64
  Status     : Stable
Graham McVicker's avatar
Graham McVicker committed
65
66
67

=cut

68
sub _tables {
69
  my $self = shift;
70

71
  return (['dna_align_feature', 'daf'],['external_db','exdb']);
72
73
}

Graham McVicker's avatar
Graham McVicker committed
74

75
76
77
78
sub _left_join{
    return (['external_db',"exdb.external_db_id = daf.external_db_id"]);
}

Graham McVicker's avatar
Graham McVicker committed
79
80
81
82
83
84
85
86
87
=head2 _columns

  Args       : none
  Example    : @columns = $self->_columns
  Description: PROTECTED implementation of abstract superclass method.  
               Returns a list of columns that are needed for object creation.
  Returntype : list of strings
  Exceptions : none
  Caller     : Bio::EnsEMBL::DBSQL::BaseFeatureAdaptor::generic_fetch
88
  Status     : Stable
Graham McVicker's avatar
Graham McVicker committed
89
90
91

=cut

92
93
sub _columns {
  my $self = shift;
94

95
  #warning, implementation of _objs_from_sth method depends on order of list
96
97
98
99
100
101
102
103
104
105
106
107
108
  return qw(daf.dna_align_feature_id
            daf.seq_region_id
            daf.analysis_id
            daf.seq_region_start
            daf.seq_region_end
            daf.seq_region_strand
            daf.hit_start
            daf.hit_end
            daf.hit_name
            daf.hit_strand
            daf.cigar_line
            daf.evalue
            daf.perc_ident
109
110
            daf.score
            daf.external_db_id
Eugene Kulesha's avatar
typo  
Eugene Kulesha committed
111
            daf.hcoverage
112
	    daf.external_data
113
	    daf.pair_dna_align_feature_id
114
115
	    exdb.db_name
	    exdb.db_display_name);
116
}
117

118

119
=head2 store
120

121
  Arg [1]    : list of Bio::EnsEMBL::DnaAlignFeatures @feats
Graham McVicker's avatar
Graham McVicker committed
122
               the features to store in the database
123
  Example    : $dna_align_feature_adaptor->store(@features);
Graham McVicker's avatar
Graham McVicker committed
124
125
  Description: Stores a list of DnaAlignFeatures in the database
  Returntype : none
126
127
128
129
130
131
132
133
134
  Exceptions : throw if any of the provided features cannot be stored
               which may occur if:
                 * The feature does not have an associate Slice
                 * The feature does not have an associated analysis
                 * The Slice the feature is associated with is on a seq_region
                   unknown to this database
               A warning is given if:
                 * The feature has already been stored in this db
  Caller     : Pipeline
135
  Status     : Stable
136
137
138

=cut

139
sub store {
140
  my ( $self, @feats ) = @_;
141

142
  throw("Must call store with features") if ( scalar(@feats) == 0 );
143
144

  my @tabs = $self->_tables;
145
  my ($tablename) = @{ $tabs[0] };
146

147
  my $db               = $self->db();
148
149
150
  my $analysis_adaptor = $db->get_AnalysisAdaptor();

  my $sth = $self->prepare(
151
152
153
154
155
156
157
158
159
160
161
162
163
    "INSERT INTO $tablename (seq_region_id, seq_region_start,
                             seq_region_end, seq_region_strand,
                             hit_start, hit_end, hit_strand, hit_name,
                             cigar_line, analysis_id, score, evalue,
                             perc_ident, external_db_id, hcoverage,
                             pair_dna_align_feature_id)
     VALUES (?,?,?,?,?,?,?,?,?,?,?,?,?,?,?,?)"    # 16 arguments
  );

FEATURE:
  foreach my $feat (@feats) {
    if ( !ref $feat || !$feat->isa("Bio::EnsEMBL::DnaDnaAlignFeature") )
    {
164
      throw("feature must be a Bio::EnsEMBL::DnaDnaAlignFeature,"
165
166
167
          . " not a ["
          . ref($feat)
          . "]." );
168
169
    }

170
171
172
173
    if ( $feat->is_stored($db) ) {
      warning( "DnaDnaAlignFeature ["
          . $feat->dbID()
          . "] is already stored in this database." );
174
175
176
      next FEATURE;
    }

177
178
    my $hstart  = $feat->hstart();
    my $hend    = $feat->hend();
179
    my $hstrand = $feat->hstrand();
180
    $self->_check_start_end_strand( $hstart, $hend, $hstrand );
181
182

    my $cigar_string = $feat->cigar_string();
183
    if ( !$cigar_string ) {
184
      $cigar_string = $feat->length() . 'M';
185
186
      warning( "DnaDnaAlignFeature does not define a cigar_string.\n"
          . "Assuming ungapped block with cigar_line=$cigar_string ." );
187
    }
188
189

    my $hseqname = $feat->hseqname();
190
    if ( !$hseqname ) {
191
192
193
      throw("DnaDnaAlignFeature must define an hseqname.");
    }

194
195
196
    if ( !defined( $feat->analysis ) ) {
      throw(
        "An analysis must be attached to the features to be stored.");
197
198
    }

199
    #store the analysis if it has not been stored yet
200
201
    if ( !$feat->analysis->is_stored($db) ) {
      $analysis_adaptor->store( $feat->analysis() );
202
    }
203

204
205
    my $original = $feat;
    my $seq_region_id;
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
    ( $feat, $seq_region_id ) = $self->_pre_store($feat);

    $sth->bind_param( 1,  $seq_region_id,        SQL_INTEGER );
    $sth->bind_param( 2,  $feat->start,          SQL_INTEGER );
    $sth->bind_param( 3,  $feat->end,            SQL_INTEGER );
    $sth->bind_param( 4,  $feat->strand,         SQL_TINYINT );
    $sth->bind_param( 5,  $hstart,               SQL_INTEGER );
    $sth->bind_param( 6,  $hend,                 SQL_INTEGER );
    $sth->bind_param( 7,  $hstrand,              SQL_TINYINT );
    $sth->bind_param( 8,  $hseqname,             SQL_VARCHAR );
    $sth->bind_param( 9,  $cigar_string,         SQL_LONGVARCHAR );
    $sth->bind_param( 10, $feat->analysis->dbID, SQL_INTEGER );
    $sth->bind_param( 11, $feat->score,          SQL_DOUBLE );
    $sth->bind_param( 12, $feat->p_value,        SQL_DOUBLE );
    $sth->bind_param( 13, $feat->percent_id,     SQL_FLOAT );
    $sth->bind_param( 14, $feat->external_db_id, SQL_INTEGER );
    $sth->bind_param( 15, $feat->hcoverage,      SQL_DOUBLE );
    $sth->bind_param( 16, $feat->pair_dna_align_feature_id,
      SQL_INTEGER );
225
226

    $sth->execute();
227
228

    $original->dbID( $sth->{'mysql_insertid'} );
229
    $original->adaptor($self);
230
  } ## end foreach my $feat (@feats)
231
232

  $sth->finish();
233
} ## end sub store
234
235
236
237
238
239
240
241
242
243
244
245
246
247


sub save {
  my ($self, $features) = @_;

  my @feats = @$features;
  throw("Must call store with features") if( scalar(@feats) == 0 );

  my @tabs = $self->_tables;
  my ($tablename) = @{$tabs[0]};

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

248
  my $sql = qq{INSERT INTO $tablename (seq_region_id, seq_region_start, seq_region_end, seq_region_strand, hit_start, hit_end, hit_strand, hit_name, cigar_line, analysis_id, score, evalue, perc_ident, external_db_id, hcoverage, pair_dna_align_feature_id, external_data) VALUES (?, ?, ?, ?, ?, ?, ?, ?, ?, ?, ?, ?, ?, ?, ?, ?, ?)};
249

250
251
  my %analyses = ();

252
253
254
255
256
257
258
259
260
261
262
263
264
265
  my $sth = $self->prepare($sql);
     
 FEATURE: foreach my $feat ( @feats ) {
    if( !ref $feat || !$feat->isa("Bio::EnsEMBL::DnaDnaAlignFeature") ) {
      throw("feature must be a Bio::EnsEMBL::DnaDnaAlignFeature,"
            . " not a [".ref($feat)."].");
    }

    if($feat->is_stored($db)) {
      warning("DnaDnaAlignFeature [".$feat->dbID."] is already stored" .
              " in this database.");
      next FEATURE;
    }

266
267
268
269
270
271
272
273
    my $hstart  = $feat->hstart || 0; # defined $feat->hstart  ? $feat->hstart : $feat->start ;
    my $hend    = $feat->hend   || 0; # defined $feat->hend    ? $feat->hend : $feat->end;
    my $hstrand = $feat->hstrand|| 0; # defined $feat->hstrand ? $feat->hstrand : $feat->strand;
    if( $hstart && $hend ) {
      if($hend < $hstart) {
        throw("Invalid Feature start/end [$hstart/$hend]. Start must be less than or equal to end.");
      }
    }
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
    my $cigar_string = $feat->cigar_string();
    if(!$cigar_string) {
      $cigar_string = $feat->length() . 'M';
      warning("DnaDnaAlignFeature does not define a cigar_string.\n" .
              "Assuming ungapped block with cigar_line=$cigar_string .");
    }

    my $hseqname = $feat->hseqname();
    if(!$hseqname) {
      throw("DnaDnaAlignFeature must define an hseqname.");
    }

    if(!defined($feat->analysis)) {
      throw("An analysis must be attached to the features to be stored.");
    }

    #store the analysis if it has not been stored yet
    if(!$feat->analysis->is_stored($db)) {
      $analysis_adaptor->store($feat->analysis());
    }

295
296
    $analyses{ $feat->analysis->dbID }++;

297
298
299
300
    my $original = $feat;
    my $seq_region_id;
    ($feat, $seq_region_id) = $self->_pre_store_userdata($feat);

301
302
    my $extra_data = $feat->extra_data ? $self->dump_data($feat->extra_data) : '';

303
304
305
306
307
308
309
310
311
312
    $sth->bind_param(1,$seq_region_id,SQL_INTEGER);
    $sth->bind_param(2,$feat->start,SQL_INTEGER);
    $sth->bind_param(3,$feat->end,SQL_INTEGER);
    $sth->bind_param(4,$feat->strand,SQL_TINYINT);
    $sth->bind_param(5,$hstart,SQL_INTEGER);
    $sth->bind_param(6,$hend,SQL_INTEGER);
    $sth->bind_param(7,$hstrand,SQL_TINYINT);
    $sth->bind_param(8,$hseqname,SQL_VARCHAR);
    $sth->bind_param(9,$cigar_string,SQL_LONGVARCHAR);
    $sth->bind_param(10,$feat->analysis->dbID,SQL_INTEGER);
313
    $sth->bind_param(11,$feat->score,SQL_DOUBLE);
314
#    $sth->bind_param(11,$feat->score); # if the above statement does not work it means you need to upgrade DBD::mysql, meantime you can replace it with this line
315
316
    $sth->bind_param(12,$feat->p_value,SQL_DOUBLE);
    $sth->bind_param(13,$feat->percent_id,SQL_FLOAT);
317
318
    $sth->bind_param(14,$feat->external_db_id,SQL_INTEGER);
    $sth->bind_param(15,$feat->hcoverage,SQL_DOUBLE);
319
320
321
    $sth->bind_param(16,$feat->pair_dna_align_feature_id,SQL_INTEGER);
    $sth->bind_param(17,$extra_data,SQL_LONGVARCHAR);

322
323

    $sth->execute();
324
325
    $original->dbID($sth->{'mysql_insertid'});
    $original->adaptor($self);
326
  }
327
328

  $sth->finish();
329
330
331

## js5 hack to update meta_coord table... 
  if( keys %analyses ) {
332

333
334
    my $sth = $self->prepare( 'select sr.coord_system_id, max(daf.seq_region_end-daf.seq_region_start) from seq_region as sr, dna_align_feature as daf where daf.seq_region_id=sr.seq_region_id and analysis_id in ('.join(',',keys %analyses).') group by coord_system_id' );
    $sth->execute;
335

336
337
338
339
    foreach( @{ $sth->fetchall_arrayref } ) {
      my $sth2 = $self->prepare( qq(insert ignore into meta_coord values("dna_align_feature",$_->[0],$_->[1])) );
      $sth2->execute;
      $sth2->finish;
340
341

      $sth2 = $self->prepare( qq(update meta_coord set max_length = $_->[1] where coord_system_id = $_->[0] and table_name="dna_align_feature" and max_length < $_->[1]) );
342
343
344
      $sth2->execute;
      $sth2->finish;
    }
345

346
347
348
    $sth->finish;
  }

349
350
351
}


352
=head2 _objs_from_sth
Graham McVicker's avatar
Graham McVicker committed
353

354
355
356
357
  Arg [1]    : DBI statement handle $sth
               an exectuted DBI statement handle generated by selecting 
               the columns specified by _columns() from the table specified 
               by _table()
Graham McVicker's avatar
Graham McVicker committed
358
359
360
  Example    : @dna_dna_align_feats = $self->_obj_from_hashref
  Description: PROTECTED implementation of superclass abstract method. 
               Creates DnaDnaAlignFeature objects from a DBI hashref
361
  Returntype : listref of Bio::EnsEMBL::DnaDnaAlignFeatures
Graham McVicker's avatar
Graham McVicker committed
362
363
  Exceptions : none
  Caller     : Bio::EnsEMBL::BaseFeatureAdaptor::generic_fetch
364
  Status     : Stable
Graham McVicker's avatar
Graham McVicker committed
365
366

=cut
367

368
sub _objs_from_sth {
369
  my ( $self, $sth, $mapper, $dest_slice ) = @_;
370

371
372
373
374
375
  #
  # This code is ugly because an attempt has been made to remove as many
  # function calls as possible for speed purposes.  Thus many caches and
  # a fair bit of gymnastics is used.
  #
376

377
378
379
380
381
  # In case of userdata we need the features on the dest_slice.  In case
  # of get_all_supporting_features dest_slice is not provided.
  my $sa = (   $dest_slice
             ? $dest_slice->adaptor()
             : $self->db()->get_SliceAdaptor() );
382
  my $aa = $self->db->get_AnalysisAdaptor();
383

384
385
386
387
388
389
  my @features;
  my %analysis_hash;
  my %slice_hash;
  my %sr_name_hash;
  my %sr_cs_hash;

390
391
392
393
394
395
396
397
398
399
400
401
402
403
404
405
406
407
408
409
410
411
  my ( $dna_align_feature_id, $seq_region_id,
       $analysis_id,          $seq_region_start,
       $seq_region_end,       $seq_region_strand,
       $hit_start,            $hit_end,
       $hit_name,             $hit_strand,
       $cigar_line,           $evalue,
       $perc_ident,           $score,
       $external_db_id,       $hcoverage,
       $extra_data,           $pair_dna_align_feature_id,
       $external_db_name,     $external_display_db_name );

  $sth->bind_columns( \( $dna_align_feature_id, $seq_region_id,
                         $analysis_id,          $seq_region_start,
                         $seq_region_end,       $seq_region_strand,
                         $hit_start,            $hit_end,
                         $hit_name,             $hit_strand,
                         $cigar_line,           $evalue,
                         $perc_ident,           $score,
                         $external_db_id,       $hcoverage,
                         $extra_data,       $pair_dna_align_feature_id,
                         $external_db_name, $external_display_db_name )
  );
412
413
414
415
416
417
418

  my $asm_cs;
  my $cmp_cs;
  my $asm_cs_vers;
  my $asm_cs_name;
  my $cmp_cs_vers;
  my $cmp_cs_name;
419
420
421
422

  if ( defined($mapper) ) {
    $asm_cs      = $mapper->assembled_CoordSystem();
    $cmp_cs      = $mapper->component_CoordSystem();
423
424
425
426
427
    $asm_cs_name = $asm_cs->name();
    $asm_cs_vers = $asm_cs->version();
    $cmp_cs_name = $cmp_cs->name();
    $cmp_cs_vers = $cmp_cs->version();
  }
428

429
430
431
432
  my $dest_slice_start;
  my $dest_slice_end;
  my $dest_slice_strand;
  my $dest_slice_length;
433
  my $dest_slice_sr_name;
434
  my $dest_slice_seq_region_id;
435

436
437
438
439
440
441
  if ( defined($dest_slice) ) {
    $dest_slice_start         = $dest_slice->start();
    $dest_slice_end           = $dest_slice->end();
    $dest_slice_strand        = $dest_slice->strand();
    $dest_slice_length        = $dest_slice->length();
    $dest_slice_sr_name       = $dest_slice->seq_region_name();
442
    $dest_slice_seq_region_id = $dest_slice->get_seq_region_id();
443
  }
444

445
446
447
FEATURE:
  while ( $sth->fetch() ) {
    # Get the analysis object.
448
449
    my $analysis = $analysis_hash{$analysis_id} ||=
      $aa->fetch_by_dbID($analysis_id);
450

451
452
    # Get the slice object.
    my $slice = $slice_hash{ "ID:" . $seq_region_id };
453

454
    if ( !defined($slice) ) {
455
      $slice = $sa->fetch_by_seq_region_id($seq_region_id);
456
457
      if ( defined($slice) ) {
        $slice_hash{ "ID:" . $seq_region_id } = $slice;
Anne Parker's avatar
 
Anne Parker committed
458
        $sr_name_hash{$seq_region_id} = $slice->seq_region_name();
459
        $sr_cs_hash{$seq_region_id}   = $slice->coord_system();
Anne Parker's avatar
 
Anne Parker committed
460
      }
461
    }
462

463
464
    my $sr_name = $sr_name_hash{$seq_region_id};
    my $sr_cs   = $sr_cs_hash{$seq_region_id};
465
466
467
468

    # Remap the feature coordinates to another coord system
    # if a mapper was provided.
    if ( defined($mapper) ) {
469
470
471
472
473
474
475
476
477
478
479
480
481
482

	if (defined $dest_slice && $mapper->isa('Bio::EnsEMBL::ChainedAssemblyMapper')  ) {
	    ( $seq_region_id,  $seq_region_start,
	      $seq_region_end, $seq_region_strand )
		=
		$mapper->map( $sr_name, $seq_region_start, $seq_region_end,
                          $seq_region_strand, $sr_cs, 1, $dest_slice);

	} else {

	    ( $seq_region_id,  $seq_region_start,
	      $seq_region_end, $seq_region_strand )
		=
		$mapper->fastmap( $sr_name, $seq_region_start, $seq_region_end,
483
                          $seq_region_strand, $sr_cs );
484
	}
485
486
487
488
489
490
491
492
493

      # Skip features that map to gaps or coord system boundaries.
      if ( !defined($seq_region_id) ) { next FEATURE }

      # Get a slice in the coord system we just mapped to.
      if ( $asm_cs == $sr_cs
           || ( $cmp_cs != $sr_cs && $asm_cs->equals($sr_cs) ) )
      {
        $slice = $slice_hash{ "ID:" . $seq_region_id } ||=
494
          $sa->fetch_by_seq_region_id($seq_region_id);
495
      } else {
496
        $slice = $slice_hash{ "ID:" . $seq_region_id } ||=
497
          $sa->fetch_by_seq_region_id($seq_region_id);
498
499
      }
    }
500

501
    # If a destination slice was provided, convert the coords.  If the
502
    # dest_slice starts at 1 and is forward strand, nothing needs doing.
503
    if ( defined($dest_slice) ) {
504
505
506
507
508
509
510
511
512
513
514
515
516
517
518
519
520
521
522
523
524
525
526
527
528
529
530
531
532
533
534
535
536
537
538
539
540
541
542
543
544
545
546
547
548
549
550
551
552
553
554
555
      my $seq_region_len = $dest_slice->seq_region_length();

      if ($dest_slice_strand == 1) { # Positive strand
		
	$seq_region_start = $seq_region_start - $dest_slice_start + 1;
	$seq_region_end   = $seq_region_end - $dest_slice_start + 1;

	if ($dest_slice->is_circular()) {
	  # Handle cicular chromosomes.

	  if ($seq_region_start > $seq_region_end) {
	    # Looking at a feature overlapping the chromsome origin.

	    if ($seq_region_end > $dest_slice_start) {

	      # Looking at the region in the beginning of the
	      # chromosome.
	      $seq_region_start -= $seq_region_len;
	    }

	    if ($seq_region_end < 0) {
	      $seq_region_end += $seq_region_len;
	    }

	  } else {

	    if (   $dest_slice_start > $dest_slice_end
		   && $seq_region_end < 0) {
	      # Looking at the region overlapping the chromosome
	      # origin and a feature which is at the beginning of the
	      # chromosome.
	      $seq_region_start += $seq_region_len;
	      $seq_region_end   += $seq_region_len;
	    }
	  }

	}		       ## end if ($dest_slice->is_circular...)

      } else {			# Negative strand

	my $start = $dest_slice_end - $seq_region_end + 1;
	my $end = $dest_slice_end - $seq_region_start + 1;

	if ($dest_slice->is_circular()) {

	  if ($dest_slice_start > $dest_slice_end) { 
	    # slice spans origin or replication

	    if ($seq_region_start >= $dest_slice_start) {
	      $end += $seq_region_len;
	      $start += $seq_region_len 
		if $seq_region_end > $dest_slice_start;
556

557
558
559
560
561
562
563
564
565
566
567
568
569
570
571
572
573
574
575
576
577
578
579
580
581
582
583
584
585
586
587
588
589
590
591
592
593
594
595
596
597
598
599
600
601
602
603
604
605
	    } elsif ($seq_region_start <= $dest_slice_end) {
	      # do nothing
	    } elsif ($seq_region_end >= $dest_slice_start) {
	      $start += $seq_region_len;
	      $end += $seq_region_len;

	    } elsif ($seq_region_end <= $dest_slice_end) {

	      $end += $seq_region_len
		if $end < 0;

	    } elsif ($seq_region_start > $seq_region_end) {
		  
	      $end += $seq_region_len;

	    } else {
		  
	    }
      
	  } else {

	    if ($seq_region_start <= $dest_slice_end and $seq_region_end >= $dest_slice_start) {
	      # do nothing
	    } elsif ($seq_region_start > $seq_region_end) {
	      if ($seq_region_start <= $dest_slice_end) {
	  
		$start -= $seq_region_len;

	      } elsif ($seq_region_end >= $dest_slice_start) {
		$end += $seq_region_len;

	      } else {
		    
	      }
	    }
	  }

	}

	$seq_region_start = $start;
	$seq_region_end = $end;
	$seq_region_strand *= -1;

      }	## end else [ if ($dest_slice_strand...)]

      # Throw away features off the end of the requested slice.
      if (    $seq_region_end < 1
	      || $seq_region_start > $dest_slice_length
	      || ( $dest_slice_seq_region_id ne $seq_region_id ) )
606
        {
607
608
          next FEATURE;
        }
609

610
      $slice = $dest_slice;
611
    }
612

613
614
615
    # Inlining the following in the hash causes major issues with 5.16 and messes up the hash 
    my $evalled_extra_data = $extra_data ? $self->get_dumped_data($extra_data) : '';

616
617
618
    # Finally, create the new DnaAlignFeature.
    push( @features,
          $self->_create_feature_fast(
619
620
621
622
623
624
625
626
627
628
629
630
631
632
633
634
635
636
             'Bio::EnsEMBL::DnaDnaAlignFeature', {
               'slice'          => $slice,
               'start'          => $seq_region_start,
               'end'            => $seq_region_end,
               'strand'         => $seq_region_strand,
               'hseqname'       => $hit_name,
               'hstart'         => $hit_start,
               'hend'           => $hit_end,
               'hstrand'        => $hit_strand,
               'score'          => $score,
               'p_value'        => $evalue,
               'percent_id'     => $perc_ident,
               'cigar_string'   => $cigar_line,
               'analysis'       => $analysis,
               'adaptor'        => $self,
               'dbID'           => $dna_align_feature_id,
               'external_db_id' => $external_db_id,
               'hcoverage'      => $hcoverage,
637
               'extra_data'     => $evalled_extra_data,
638
639
640
641
642
643
               'dbname'                    => $external_db_name,
               'db_display_name'           => $external_display_db_name,
               'pair_dna_align_feature_id' => $pair_dna_align_feature_id
             } ) );

  } ## end while ( $sth->fetch() )
644

645
  return \@features;
646
} ## end sub _objs_from_sth
647

Glenn Proctor's avatar
Glenn Proctor committed
648
649
650
651
=head2 list_dbIDs

  Arg [1]    : none
  Example    : @feature_ids = @{$dna_align_feature_adaptor->list_dbIDs()};
652
653
  Description: Gets an array of internal ids for all dna align features in 
               the current db
654
  Arg[1]     : <optional> int. not 0 for the ids to be sorted by the seq_region.
Glenn Proctor's avatar
Glenn Proctor committed
655
656
657
  Returntype : list of ints
  Exceptions : none
  Caller     : ?
658
  Status     : Stable
Glenn Proctor's avatar
Glenn Proctor committed
659
660
661
662

=cut

sub list_dbIDs {
663
   my ($self, $ordered) = @_;
Glenn Proctor's avatar
Glenn Proctor committed
664

665
   return $self->_list_dbIDs("dna_align_feature",undef, $ordered);
Glenn Proctor's avatar
Glenn Proctor committed
666
667
}

668
669
670
1;