PowerFilter/whiten safeguard default value

gallery/experiments/simulated_abinitio_pipeline.py

@@ -54,7 +54,7 @@
 # ---------------
 # Start with the hi-res volume map EMDB-2660 sourced from EMDB,
 # https://www.ebi.ac.uk/emdb/EMD-2660, and dowloaded via ASPIRE's downloader utility.
-og_v = emdb_2660().astype(np.float64)
+og_v = emdb_2660()
 logger.info("Original volume map data" f" shape: {og_v.shape} dtype:{og_v.dtype}")
 
 
@@ -94,7 +94,6 @@ def noise_function(x, y):
     vols=og_v,
     noise_adder=custom_noise,
     unique_filters=ctf_filters,
-    dtype=np.float64,
 )
 
 # Downsample

src/aspire/operators/filters.py

@@ -217,6 +217,7 @@ def evaluate_grid(self, L, *args, dtype=np.float32, **kwargs):
             eps = self._epsilon
             if eps is None:
                 eps = np.finfo(filter_vals.dtype).eps
+                eps = (100 * eps) ** (-1 / self._power)
             condition = abs(filter_vals) < eps
             num_less_eps = np.count_nonzero(condition)
             if num_less_eps > 0:

src/aspire/source/image.py

@@ -808,6 +808,7 @@ def whiten(self, noise_estimate=None, epsilon=None):
 
         if epsilon is None:
             epsilon = np.finfo(self.dtype).eps
+            epsilon = (100 * epsilon) ** 2
 
         logger.info("Whitening source object")
         whiten_filter = PowerFilter(noise_filter, power=-0.5, epsilon=epsilon)

tests/test_filters.py

@@ -353,18 +353,19 @@ def epsilon(request):
 def test_power_filter_safeguard(dtype, epsilon, caplog):
     L = 25
     arr = np.ones((L, L), dtype=dtype)
+    power = -0.5
 
-    # Set a few values below machine epsilon.
+    # Set a few values below default safeguard.
     num_eps = 3
     eps = epsilon
     if eps is None:
-        eps = np.finfo(dtype).eps
+        eps = (100 * np.finfo(dtype).eps) ** (-1 / power)
     arr[L // 2, L // 2 : L // 2 + num_eps] = eps / 2
 
     # For negative powers, values below machine eps will be set to zero.
     filt = PowerFilter(
         filter=ArrayFilter(arr),
-        power=-0.5,
+        power=power,
         epsilon=epsilon,
     )
 

tests/test_preprocess_pipeline.py

@@ -155,6 +155,33 @@ def test_whiten_safeguard(dtype):
     np.testing.assert_allclose(whiten_filt[ind], 0.0)
 
 
+@pytest.mark.parametrize("dtype", [np.float32, np.float64])
+def test_whiten_safeguard_default(dtype):
+    """
+    This test catches the case found in the simulated_abinitio_pipeline.py
+    of images being zeroed out by to strict of a safeguard on the whiten filter.
+    """
+    L = 25
+    sim = get_sim_object(L, dtype)
+    noise_estimator = AnisotropicNoiseEstimator(sim)
+
+    # Alter noise_estimator filter values to be below machine eps.
+    # These are values comparable to those in
+    # gallery/experiments/simulated_abinitio_pipeline.py.
+    mach_eps = np.finfo(dtype).eps
+    noise_estimator.filter.xfer_fn_array *= mach_eps
+    assert noise_estimator.filter.xfer_fn_array.min() < mach_eps
+
+    # Whiten images
+    sim = sim.whiten(noise_estimator.filter)
+
+    # Get whitening_filter from generation pipeline.
+    whiten_filt = sim.generation_pipeline.xforms[0].filter.evaluate_grid(sim.L)
+
+    # Check that no values in the whiten filter have been zeroed out by safeguard.
+    assert np.count_nonzero(whiten_filt == 0) == 0
+
+
 @pytest.mark.parametrize("L, dtype", params)
 def testInvertContrast(L, dtype):
     sim1 = get_sim_object(L, dtype)