changed retrieval of densities

modules/Bio/EnsEMBL/DBSQL/DensityFeatureAdaptor.pm

@@ -162,6 +162,7 @@ sub fetch_all_by_Slice {
   my $density_type = undef;
   my $best_ratio_large = undef;
   my $density_type_large = undef;
+  my %dt_ratio_hash;
 
   foreach my $dt (@dtypes) {
 
@@ -175,29 +176,19 @@ sub fetch_all_by_Slice {
       my $block_size = $slice->seq_region_length() / $dt->region_features();
       $ratio = $wanted_block_size / $block_size;
     }
-
+    
     # we prefer to use a block size that's smaller than the required one
-    # (better results on interpolation). remember larger block sizes though
-    # in case there is no smaller one in the database
-    if ($ratio < 1) {
-      if(!defined($best_ratio_large) || $ratio > $best_ratio_large) {
-        $best_ratio_large = $ratio;
-        $density_type_large = $dt;
-      }
-    } else {
-      if(!defined($best_ratio) || $ratio < $best_ratio) {
-        $best_ratio = $ratio;
-        $density_type = $dt;
-      }
+    # (better results on interpolation).
+    # give larger bits a disadvantage and make them comparable
+    if( $ratio < 1 ) {
+      $ratio = 5/$ratio;
     }
-  }
-  # fall back to larger block size if there is no smaller black size 
-  # or it would require retrieving too many features
-  if( !$best_ratio || $best_ratio > 30 ) {
-    $best_ratio = $best_ratio_large;
-    $density_type = $density_type_large;
+
+    $dt_ratio_hash{ $ratio } = $dt;
   }
 
+  $best_ratio = (sort {$a<=>$b} (keys %dt_ratio_hash))[0];
+  
   #the ratio was not good enough, or this logic name was not in the
   #density type table, return an empty list
   if(!defined($best_ratio) ||
@@ -205,6 +196,8 @@ sub fetch_all_by_Slice {
     return [];
   }
 
+  $density_type = $dt_ratio_hash{$best_ratio};
+
   my $constraint = "df.density_type_id = " . $density_type->dbID();
 
   my @features =

modules/t/densityFeatureAdaptor.t

@@ -64,7 +64,6 @@ $dta->store($dt);
 ok($dt->is_stored($db));
 
 
-
 my $slice_adaptor = $db->get_SliceAdaptor();
 my $slice = $slice_adaptor->fetch_by_region('chromosome', '20', 1, ($block_size*10));
 
@@ -131,16 +130,17 @@ ok($dt->is_stored($db));
 
 @density_features = ();
 my $chr_20_slice = $slice_adaptor->fetch_by_region( "chromosome", 20 );
-my $step = $chr_20_slice->length() / 300;
+my $step = POSIX::ceil( $chr_20_slice->length() / 300);
 $start = 1;
 while( $start < $chr_20_slice->length() ) {
-  my $end = int( $start + $step );
+  my $end = $start+$step -1;
+  if( $end > $chr_20_slice->length ) { $end = $chr_20_slice->length();}
   push @density_features, Bio::EnsEMBL::DensityFeature->new(-seq_region    => $chr_20_slice,
-							     -start         => int($start),
-							     -end           => int( $start+$step),
+							     -start         => $start,
+							     -end           => $end,
 							     -density_type  => $dt,
 							     -density_value => 5 );
-  $start += $step + 1;
+  $start += $step;
 }
 
 $dfa->store( @density_features );
@@ -155,6 +155,38 @@ ok( abs( $stored_features[0]->density_value() - 15) < 0.0001 );
 debug( "Density value = ".$stored_features[0]->density_value() );
 
 
+# test the retreival of the right sized features
+# first yet another density size
+# comes to about 1000bp on chr 20
+
+$dt = Bio::EnsEMBL::DensityType->new(-analysis   => $analysis,
+				     -region_features => 30_000,
+			             -value_type => 'sum');
+$dta->store($dt);	
+
+# need some features 
+@density_features = ();
+$step = POSIX::ceil( $chr_20_slice->length() / 30_000);
+for my $arr ( [1,1000,1],[1001,2000,2],[2001,3000,3],[31_000_000,31_000_999,31.0],[31_500_000, 31_500_999,31.5] ) {
+  push @density_features, Bio::EnsEMBL::DensityFeature->new(-seq_region    => $chr_20_slice,
+							     -start         => $arr->[0],
+							     -end           => $arr->[1],
+							     -density_type  => $dt,
+							     -density_value => $arr->[2] );
+}
+$dfa->store( @density_features );
+
+# now check for retrieval
+my $sub_Slice = $chr_20_slice->sub_Slice( 1, 1000_000 );
+@stored_features = @{$dfa->fetch_all_by_Slice( $sub_Slice, 'GeneDensityTest', 2 )};
+ok( $stored_features[0]->length() > 10000 );
+@stored_features = @{$dfa->fetch_all_by_Slice( $sub_Slice, 'GeneDensityTest', 10 )};
+ok( $stored_features[0]->length() > 10000 );
+@stored_features = @{$dfa->fetch_all_by_Slice( $sub_Slice, 'GeneDensityTest', 100 )};
+ok( $stored_features[0]->length() == 1000 );
+
+
+
 
 
 $multi->restore('core', 'analysis');