diff --git a/gallery/tutorials/aspire_introduction.py b/gallery/tutorials/aspire_introduction.py index cffe6d544e..b0013335f5 100644 --- a/gallery/tutorials/aspire_introduction.py +++ b/gallery/tutorials/aspire_introduction.py @@ -571,7 +571,7 @@ def noise_function(x, y): # Generate several CTFs. ctf_filters = [ - RadialCTFFilter(pixel_size=vol_ds.pixel_size, defocus=d) + RadialCTFFilter(defocus=d) for d in np.linspace(defocus_min, defocus_max, defocus_ct) ] diff --git a/gallery/tutorials/pipeline_demo.py b/gallery/tutorials/pipeline_demo.py index ad6143c669..0a192ecaf1 100644 --- a/gallery/tutorials/pipeline_demo.py +++ b/gallery/tutorials/pipeline_demo.py @@ -67,7 +67,7 @@ defocus_ct = 7 ctf_filters = [ - RadialCTFFilter(pixel_size=original_vol.pixel_size, defocus=d) + RadialCTFFilter(defocus=d) for d in np.linspace(defocus_min, defocus_max, defocus_ct) ] diff --git a/gallery/tutorials/tutorials/cov2d_simulation.py b/gallery/tutorials/tutorials/cov2d_simulation.py index bf64684db4..5d9b469b04 100644 --- a/gallery/tutorials/tutorials/cov2d_simulation.py +++ b/gallery/tutorials/tutorials/cov2d_simulation.py @@ -68,7 +68,7 @@ print("Initialize simulation object and CTF filters.") # Create filters ctf_filters = [ - RadialCTFFilter(pixel_size, voltage, defocus=d, Cs=2.0, alpha=0.1) + RadialCTFFilter(voltage, defocus=d, Cs=2.0, alpha=0.1) for d in np.linspace(defocus_min, defocus_max, defocus_ct) ] @@ -94,6 +94,7 @@ amplitudes=1.0, dtype=dtype, noise_adder=noise_adder, + pixel_size=pixel_size, ) @@ -108,7 +109,9 @@ h_idx = sim.filter_indices # Evaluate CTF in the 8X8 FB basis -h_ctf_fb = [ffbbasis.filter_to_basis_mat(filt) for filt in ctf_filters] +h_ctf_fb = [ + ffbbasis.filter_to_basis_mat(filt, pixel_size=pixel_size) for filt in ctf_filters +] # Get clean images from projections of 3D map. print("Apply CTF filters to clean images.") diff --git a/gallery/tutorials/tutorials/cov3d_simulation.py b/gallery/tutorials/tutorials/cov3d_simulation.py index fdec8f601e..4e46596a36 100644 --- a/gallery/tutorials/tutorials/cov3d_simulation.py +++ b/gallery/tutorials/tutorials/cov3d_simulation.py @@ -36,6 +36,7 @@ L=img_size, C=3, dtype=dtype, + pixel_size=10, ).generate() # Create a simulation object with specified filters @@ -43,9 +44,7 @@ L=img_size, n=num_imgs, vols=vols, - unique_filters=[ - RadialCTFFilter(pixel_size=10, defocus=d) for d in np.linspace(1.5e4, 2.5e4, 7) - ], + unique_filters=[RadialCTFFilter(defocus=d) for d in np.linspace(1.5e4, 2.5e4, 7)], dtype=dtype, ) @@ -89,7 +88,7 @@ eigs_est, lambdas_est = eigs(covar_est, num_eigs) # Eigs returns column-major, so we transpose and construct a volume. -eigs_est = Volume(np.transpose(eigs_est, (3, 0, 1, 2))) +eigs_est = Volume(np.transpose(eigs_est, (3, 0, 1, 2)), pixel_size=vols.pixel_size) # Truncate the eigendecomposition. Since we know the true rank of the # covariance matrix, we enforce it here. diff --git a/gallery/tutorials/tutorials/ctf.py b/gallery/tutorials/tutorials/ctf.py index a06ce51d87..218612dc8c 100644 --- a/gallery/tutorials/tutorials/ctf.py +++ b/gallery/tutorials/tutorials/ctf.py @@ -24,6 +24,7 @@ # Image size to use throughout the demo. IMG_SIZE = 512 + # %% # Visualizing the CTF # ------------------- @@ -40,7 +41,6 @@ from aspire.operators import CTFFilter, RadialCTFFilter radial_ctf_filter = RadialCTFFilter( - pixel_size=1, # angstrom voltage=200, # kV defocus=10000, # angstrom, 10000 A = 1 um Cs=2.26, # Spherical aberration constant @@ -51,7 +51,7 @@ # The CTF filter can be visualized as an image once it is evaluated at a specific resolution. # More specifically the following code will return a transfer function as an array, # which is then plotted. -rctf_fn = radial_ctf_filter.evaluate_grid(IMG_SIZE) +rctf_fn = radial_ctf_filter.evaluate_grid(IMG_SIZE, pixel_size=1) plt.imshow(rctf_fn) plt.colorbar() plt.show() @@ -66,7 +66,6 @@ # and the values more typically differ by a few percent. ctf_filter = CTFFilter( - pixel_size=1, # angstrom voltage=200, # kV defocus_u=15000, # angstrom, 10000 A = 1 um defocus_v=10000, @@ -77,7 +76,7 @@ ) # Again, we plot it, and note the difference from the RadialCTFFilter. -plt.imshow(ctf_filter.evaluate_grid(IMG_SIZE)) +plt.imshow(ctf_filter.evaluate_grid(IMG_SIZE, pixel_size=1)) plt.colorbar() plt.show() @@ -90,7 +89,7 @@ # array returned by ASPIRE's ``CTFFilter.evaluate_grid``. -ctf_sign = np.sign(radial_ctf_filter.evaluate_grid(IMG_SIZE)) +ctf_sign = np.sign(radial_ctf_filter.evaluate_grid(IMG_SIZE, pixel_size=1)) plt.imshow(ctf_sign) plt.colorbar() plt.show() @@ -140,10 +139,9 @@ def generate_example_image(L, noise_variance=0.1): return img -img = generate_example_image(IMG_SIZE) -plt.imshow(img) -plt.colorbar() -plt.show() +img = Image(generate_example_image(IMG_SIZE), pixel_size=1) +img.show() + # %% # Apply CTF and Phase Flipping @@ -157,26 +155,19 @@ def generate_example_image(L, noise_variance=0.1): # Construct a range of CTF filters. defoci = [2500, 5000, 10000, 20000] ctf_filters = [ - RadialCTFFilter(pixel_size=1 / 2, voltage=200, defocus=d, Cs=2.26, alpha=0.07, B=0) - for d in defoci + RadialCTFFilter(voltage=200, defocus=d, Cs=2.26, alpha=0.07, B=0) for d in defoci ] -# %% -# .. note:: -# Pixel size was chosen to demonstrate effects similar to lecture notes, -# but at a higher resolution. - - # %% # Generate CTF corrupted Images # """"""""""""""""""""""""""""" # Generate images corrupted by progressively increasing defocus. # For each defocus, apply filter to the base image. -imgs = np.empty((len(defoci), IMG_SIZE, IMG_SIZE)) +imgs = Image(np.empty((len(defoci), IMG_SIZE, IMG_SIZE)), pixel_size=1) for i, ctf in enumerate(ctf_filters): - imgs[i] = Image(img).filter(ctf)[0] -Image(imgs).show() + imgs[i] = img.filter(ctf)[0] +imgs.show() # %% # Generate Phase Flipped Images @@ -191,13 +182,13 @@ def generate_example_image(L, noise_variance=0.1): # Compute the signs of this CTF # In practice, this would be an estimated CTF, # but in the demo we have the luxury of using the model CTF that was applied. - signs = np.sign(ctf.evaluate_grid(IMG_SIZE)) + signs = np.sign(ctf.evaluate_grid(IMG_SIZE, pixel_size=1)) # Apply to the image in Fourier space. phase_flipped_imgs_f[i] = signs * imgs_f[i] # Construct the centered 2D FFT of the images. phase_flipped_imgs = aspire.numeric.fft.centered_ifft2(phase_flipped_imgs_f).real -Image(phase_flipped_imgs).show() +Image(phase_flipped_imgs, pixel_size=1).show() # %% # .. warning:: @@ -218,7 +209,6 @@ def generate_example_image(L, noise_variance=0.1): # along with an erroneous CTF filter. bad_est_ctf_filter = RadialCTFFilter( - pixel_size=1, voltage=200, defocus=14000, # Modeled CTF was 10000 Cs=2.26, @@ -226,7 +216,7 @@ def generate_example_image(L, noise_variance=0.1): B=0, ) # Evaluate Filter, returning a Numpy array. -bad_ctf_fn = bad_est_ctf_filter.evaluate_grid(IMG_SIZE) +bad_ctf_fn = bad_est_ctf_filter.evaluate_grid(IMG_SIZE, pixel_size=1) c = IMG_SIZE // 2 + 1 plt.plot(rctf_fn[c, c:], label="Model CTF") # radial_ctf_filter @@ -257,10 +247,8 @@ def generate_example_image(L, noise_variance=0.1): from aspire.ctf import estimate_ctf # Using our radial_ctf_filter from earlier, corrupt an image. -test_img = Image(img).filter(radial_ctf_filter) -plt.imshow(test_img.asnumpy()[0]) -plt.colorbar() -plt.show() +test_img = img.filter(radial_ctf_filter) +test_img.show() # Create the image file in a tmp dir with TemporaryDirectory() as d: @@ -268,7 +256,7 @@ def generate_example_image(L, noise_variance=0.1): radial_ctf_est = estimate_ctf( data_folder=d, - pixel_size=radial_ctf_filter.pixel_size, + pixel_size=1, cs=radial_ctf_filter.Cs, amplitude_contrast=radial_ctf_filter.alpha, voltage=radial_ctf_filter.voltage, @@ -296,14 +284,13 @@ def generate_example_image(L, noise_variance=0.1): # Create a filter and evaluate. est_ctf = RadialCTFFilter( - pixel_size=est["pixel_size"], voltage=est["voltage"], defocus=defocus, # Modeled CTF was 10000 Cs=est["cs"], alpha=est["amplitude_contrast"], B=0, ) -est_ctf_fn = est_ctf.evaluate_grid(IMG_SIZE) +est_ctf_fn = est_ctf.evaluate_grid(IMG_SIZE, pixel_size=1) # Compare the model CTF with the estimated CTF. c = IMG_SIZE // 2 + 1 @@ -347,7 +334,7 @@ def generate_example_image(L, noise_variance=0.1): from aspire.source import Simulation # Create the Source. ``ctf_filters`` are re-used from earlier section. -src = Simulation(L=64, n=4, unique_filters=ctf_filters) +src = Simulation(L=64, n=4, unique_filters=ctf_filters, pixel_size=1) src.images[:4].show() # %% diff --git a/gallery/tutorials/tutorials/image_class.py b/gallery/tutorials/tutorials/image_class.py index ca4df61961..0a2a33eb3a 100644 --- a/gallery/tutorials/tutorials/image_class.py +++ b/gallery/tutorials/tutorials/image_class.py @@ -22,7 +22,9 @@ # Create an ASPIRE Image instance from the data # We'll tell it to convert to floating point data as well. -im = Image(img_data, dtype=np.float64) +# Adding a `pixel_size` will pass through to subsequent codes like +# filtering which may require it. +im = Image(img_data, pixel_size=1, dtype=np.float64) # %% # Plot the Image Stack @@ -57,6 +59,6 @@ # CTF Filter # ---------- -# pixel_size/defous_u/defocus_v in angstrom, voltage in kV -filter = CTFFilter(pixel_size=1, voltage=100, defocus_u=1500, defocus_v=2000) +# defous_u/defocus_v in angstrom, voltage in kV +filter = CTFFilter(voltage=100, defocus_u=1500, defocus_v=2000) im.filter(filter).show() diff --git a/gallery/tutorials/tutorials/micrograph_source.py b/gallery/tutorials/tutorials/micrograph_source.py index 65a32ecf49..62b3a78b6f 100644 --- a/gallery/tutorials/tutorials/micrograph_source.py +++ b/gallery/tutorials/tutorials/micrograph_source.py @@ -143,6 +143,7 @@ vol = AsymmetricVolume( L=100, C=1, + pixel_size=4, seed=1234, dtype=np.float32, ).generate() @@ -181,7 +182,7 @@ # Create our CTF Filter and add it to a list. # This configuration will apply the same CTF to all particles. ctfs = [ - RadialCTFFilter(pixel_size=4, voltage=200, defocus=15000, Cs=2.26, alpha=0.07, B=0), + RadialCTFFilter(voltage=200, defocus=15000, Cs=2.26, alpha=0.07, B=0), ] src = MicrographSimulation( diff --git a/gallery/tutorials/tutorials/orient3d_simulation.py b/gallery/tutorials/tutorials/orient3d_simulation.py index 085eb9ac23..142fe177f4 100644 --- a/gallery/tutorials/tutorials/orient3d_simulation.py +++ b/gallery/tutorials/tutorials/orient3d_simulation.py @@ -40,7 +40,6 @@ # Specify the CTF parameters not used for this example # but necessary for initializing the simulation object -pixel_size = 5 # Pixel size of the images (in angstroms) voltage = 200 # Voltage (in KV) defocus_min = 1.5e4 # Minimum defocus value (in angstroms) defocus_max = 2.5e4 # Maximum defocus value (in angstroms) @@ -51,7 +50,7 @@ print("Initialize simulation object and CTF filters.") # Create CTF filters filters = [ - RadialCTFFilter(pixel_size, voltage, defocus=d, Cs=2.0, alpha=0.1) + RadialCTFFilter(voltage, defocus=d, Cs=2.0, alpha=0.1) for d in np.linspace(defocus_min, defocus_max, defocus_ct) ] @@ -74,7 +73,9 @@ # Create a simulation object with specified filters and the downsampled 3D map print("Use downsampled map to creat simulation object.") -sim = Simulation(L=img_size, n=num_imgs, vols=vols, unique_filters=filters, dtype=dtype) +sim = Simulation( + L=img_size, n=num_imgs, vols=vols, unique_filters=filters, pixel_size=5, dtype=dtype +) print("Get true rotation angles generated randomly by the simulation object.") rots_true = sim.rotations diff --git a/gallery/tutorials/tutorials/preprocess_imgs_sim.py b/gallery/tutorials/tutorials/preprocess_imgs_sim.py index af6c7a4314..4e70c5f4cc 100644 --- a/gallery/tutorials/tutorials/preprocess_imgs_sim.py +++ b/gallery/tutorials/tutorials/preprocess_imgs_sim.py @@ -54,7 +54,7 @@ print("Initialize simulation object and CTF filters.") # Create CTF filters ctf_filters = [ - RadialCTFFilter(pixel_size, voltage, defocus=d, Cs=2.0, alpha=0.1) + RadialCTFFilter(voltage, defocus=d, Cs=2.0, alpha=0.1) for d in np.linspace(defocus_min, defocus_max, defocus_ct) ] @@ -75,6 +75,7 @@ vols=vols, unique_filters=ctf_filters, noise_adder=noise_adder, + pixel_size=pixel_size, ) # %% diff --git a/src/aspire/basis/ffb_2d.py b/src/aspire/basis/ffb_2d.py index f39e89f683..2f871af9f5 100644 --- a/src/aspire/basis/ffb_2d.py +++ b/src/aspire/basis/ffb_2d.py @@ -247,6 +247,8 @@ def filter_to_basis_mat(self, f, **kwargs): " Use `method=None`." ) + pixel_size = kwargs.get("pixel_size", None) + # These form a circular dependence, import locally until time to clean up. from aspire.basis.basis_utils import lgwt @@ -269,7 +271,9 @@ def filter_to_basis_mat(self, f, **kwargs): omegay = k * np.sin(theta) omega = 2 * np.pi * np.vstack((omegax.flatten("C"), omegay.flatten("C"))) - h_vals2d = h_fun(omega).reshape(n_k, n_theta).astype(self.dtype) + h_vals2d = ( + h_fun(omega, pixel_size=pixel_size).reshape(n_k, n_theta).astype(self.dtype) + ) h_vals = np.sum(h_vals2d, axis=1) / n_theta # Represent 1D function values in basis diff --git a/src/aspire/basis/fle_2d.py b/src/aspire/basis/fle_2d.py index 7e960cf471..cc6ea04672 100644 --- a/src/aspire/basis/fle_2d.py +++ b/src/aspire/basis/fle_2d.py @@ -779,6 +779,7 @@ def filter_to_basis_mat(self, f, **kwargs): "`FLEBasis2D.filter_to_basis_mat` method {method} not supported." " Use `method=None`." ) + pixel_size = kwargs.get("pixel_size", None) # Get the filter's evaluate function. h_fun = f.evaluate @@ -802,7 +803,7 @@ def filter_to_basis_mat(self, f, **kwargs): omega = 2 * xp.pi * xp.vstack((omegax.flatten("C"), omegay.flatten("C"))) h_vals2d = ( - xp.asarray(h_fun(omega)) + xp.asarray(h_fun(omega, pixel_size=pixel_size)) .reshape(n_k, n_theta) .astype(self.dtype, copy=False) ) diff --git a/src/aspire/basis/fspca.py b/src/aspire/basis/fspca.py index 07d0c65a07..0375f91833 100644 --- a/src/aspire/basis/fspca.py +++ b/src/aspire/basis/fspca.py @@ -617,7 +617,7 @@ def shift(self, coef, shifts): self.evaluate_to_image_basis(coef).shift(shifts) ) - def filter_to_basis_mat(self, f): + def filter_to_basis_mat(self, f, **kwargs): """ Convert a filter into a basis representation. diff --git a/src/aspire/basis/steerable.py b/src/aspire/basis/steerable.py index 7012c0eeba..b5051065da 100644 --- a/src/aspire/basis/steerable.py +++ b/src/aspire/basis/steerable.py @@ -481,7 +481,7 @@ def to_complex(self, coef): # implemented. This is intended to encourage future basis authors # to consider this method for their application. @abc.abstractmethod - def filter_to_basis_mat(self, f, method="evaluate_t", truncate=True): + def filter_to_basis_mat(self, f, method="evaluate_t", truncate=True, **kwargs): """ Convert a filter into a basis operator representation. @@ -504,7 +504,13 @@ def filter_to_basis_mat(self, f, method="evaluate_t", truncate=True): " Try `evaluate_t` or `expand`." ) - coef = Coef(self, np.eye(self.count, dtype=self.dtype)) + # Note this may raise at filter run time for filters requiring `pixel_size` (eg CTFFilter) + # Alternative is to make `pixel_size` required for all calls to `filter_to_basis_mat`. + coef = Coef( + self, + np.eye(self.count, dtype=self.dtype), + pixel_size=kwargs.get("pixel_size", None), + ) img = coef.evaluate() # Expansion can fail for some filters on specific basis vectors. diff --git a/src/aspire/covariance/covar2d.py b/src/aspire/covariance/covar2d.py index cb60c63fa6..f3898f3fe8 100644 --- a/src/aspire/covariance/covar2d.py +++ b/src/aspire/covariance/covar2d.py @@ -535,7 +535,10 @@ def _build(self): logger.info("Represent CTF filters in basis") unique_filters = src.unique_filters self.ctf_idx = src.filter_indices - self.ctf_basis = [self.basis.filter_to_basis_mat(f) for f in unique_filters] + self.ctf_basis = [ + self.basis.filter_to_basis_mat(f, pixel_size=self.src.pixel_size) + for f in unique_filters + ] def _calc_rhs(self): src = self.src diff --git a/src/aspire/image/image.py b/src/aspire/image/image.py index 63f6a8b45e..5280f47bad 100644 --- a/src/aspire/image/image.py +++ b/src/aspire/image/image.py @@ -594,7 +594,10 @@ def filter(self, filter): # Second note, filter and grid dtype may not match image dtype, # upcast both here for most accurate convolution. filter_values = xp.asarray( - filter.evaluate_grid(self.resolution, dtype=np.float64), dtype=np.float64 + filter.evaluate_grid( + self.resolution, dtype=np.float64, pixel_size=self.pixel_size + ), + dtype=np.float64, ) # Convolve diff --git a/src/aspire/image/xform.py b/src/aspire/image/xform.py index 316d8ca3c8..93b40228a2 100644 --- a/src/aspire/image/xform.py +++ b/src/aspire/image/xform.py @@ -423,7 +423,7 @@ def _indexed_operation(self, im, indices, which): fn_handle = getattr(xform, which) im_data[im_data_indices] = fn_handle(im[im_data_indices]).asnumpy() - return Image(im_data) + return Image(im_data, pixel_size=im.pixel_size) def _forward(self, im, indices): return self._indexed_operation(im, indices, "forward") diff --git a/src/aspire/noise/noise.py b/src/aspire/noise/noise.py index d9e136bf71..7153282ae8 100644 --- a/src/aspire/noise/noise.py +++ b/src/aspire/noise/noise.py @@ -61,7 +61,7 @@ def _forward(self, im, indices): im_s = Image(im_s).filter(self.noise_filter).asnumpy()[0] _im[i] += im_s[: im.resolution, : im.resolution] - return Image(_im) + return Image(_im, pixel_size=im.pixel_size) @abc.abstractproperty def noise_var(self): diff --git a/src/aspire/operators/filters.py b/src/aspire/operators/filters.py index d3788d0d1f..3c3f362fc3 100644 --- a/src/aspire/operators/filters.py +++ b/src/aspire/operators/filters.py @@ -35,7 +35,7 @@ def evaluate_src_filters_on_grid(src, indices=None): for i, filt in enumerate(src.unique_filters): idx_k = np.where(src.filter_indices[indices] == i)[0] if len(idx_k) > 0: - filter_values = filt.evaluate(omega) + filter_values = filt.evaluate(omega, pixel_size=src.pixel_size) h[:, idx_k] = np.column_stack((filter_values,) * len(idx_k)) h = np.reshape(h, grid2d["x"].shape + (len(indices),)) @@ -60,7 +60,7 @@ def __str__(self): """ return self.__class__.__name__ - def evaluate(self, omega): + def evaluate(self, omega, **kwargs): """ Evaluate the filter at specified frequencies. @@ -82,17 +82,17 @@ def evaluate(self, omega): omega, idx = np.unique(omega, return_inverse=True) omega = np.vstack((omega, np.zeros_like(omega))) - h = self._evaluate(omega) + h = self._evaluate(omega, **kwargs) if self.radial: h = np.take(h, idx) return h - def _evaluate(self, omega): + def _evaluate(self, omega, **kwargs): raise NotImplementedError("Subclasses should implement this method") - def basis_mat(self, basis): + def basis_mat(self, basis, **kwargs): """ Represent the filter in `basis`. @@ -100,7 +100,7 @@ def basis_mat(self, basis): :return: `basis` representation of this filter. Return type will depend on `basis`. """ - return basis.filter_to_basis_mat(self) + return basis.filter_to_basis_mat(self, **kwargs) def scale(self, c=1): """ @@ -154,8 +154,8 @@ def __init__(self, filter_in): self._filter = filter_in super().__init__() - def evaluate(self, omega): - return self._filter.evaluate(-omega) + def evaluate(self, omega, **kwargs): + return self._filter.evaluate(-omega, **kwargs) class FunctionFilter(Filter): @@ -177,7 +177,8 @@ def __init__(self, f, dim=None): # (i.e. at runtime, we will still expect the incoming omega values to have x and y components). super().__init__(dim=dim, radial=dim > n_args) - def _evaluate(self, omega): + def _evaluate(self, omega, **kwargs): + # Note kwargs are not used here, this might be trouble return self.f(*omega) @@ -201,8 +202,8 @@ def __init__(self, filter, power=1, epsilon=None): self._epsilon = epsilon super().__init__(dim=filter.dim, radial=filter.radial) - def _evaluate(self, omega): - return self._filter.evaluate(omega) ** self._power + def _evaluate(self, omega, **kwargs): + return self._filter.evaluate(omega, **kwargs) ** self._power @lru_cache(maxsize=config["cache"]["filter_cache_size"].get()) # noqa: B019 def evaluate_grid(self, L, *args, dtype=np.float32, **kwargs): @@ -245,8 +246,8 @@ def __init__(self, filter, f): self._f = f super().__init__(dim=filter.dim, radial=filter.radial) - def _evaluate(self, omega): - return self._f(self._filter.evaluate(omega)) + def _evaluate(self, omega, **kwargs): + return self._f(self._filter.evaluate(omega, **kwargs)) class MultiplicativeFilter(Filter): @@ -258,10 +259,10 @@ def __init__(self, *args): super().__init__(dim=args[0].dim, radial=all(c.radial for c in args)) self._components = args - def _evaluate(self, omega): + def _evaluate(self, omega, **kwargs): res = 1 for c in self._components: - res *= c.evaluate(omega) + res *= c.evaluate(omega, **kwargs) return res @@ -275,8 +276,8 @@ def __init__(self, filt, scale): self._scale = scale super().__init__(dim=filt.dim, radial=filt.radial) - def _evaluate(self, omega): - return self._filter.evaluate(omega / self._scale) + def _evaluate(self, omega, **kwargs): + return self._filter.evaluate(omega / self._scale, **kwargs) def __str__(self): """ @@ -323,7 +324,7 @@ def __init__(self, xfer_fn_array): self.xfer_fn_array = xfer_fn_array - def _evaluate(self, omega): + def _evaluate(self, omega, **kwargs): _input_pts = tuple(np.linspace(1, x, x) for x in self.xfer_fn_array.shape) # TODO: This part could do with some documentation - not intuitive! @@ -391,7 +392,7 @@ def __init__(self, dim=None, value=1): def __repr__(self): return f"Scalar Filter (dim={self.dim}, value={self.value})" - def _evaluate(self, omega): + def _evaluate(self, omega, **kwargs): return self.value * np.ones_like(omega) @@ -415,7 +416,6 @@ class CTFFilter(Filter): def __init__( self, - pixel_size=1, voltage=200, defocus_u=15000, defocus_v=15000, @@ -430,7 +430,6 @@ def __init__( Note if comparing to legacy MATLAB cryo_CTF_Relion, take care regarding defocus unit conversion to nm. - :param pixel_size: Pixel size in angstrom, default 1. :param voltage: Electron voltage in kV :param defocus_u: Defocus depth along the u-axis in angstrom :param defocus_v: Defocus depth along the v-axis in angstrom @@ -440,7 +439,6 @@ def __init__( :param B: Envelope decay in inverse square angstrom (default 0) """ super().__init__(dim=2, radial=defocus_u == defocus_v) - self.pixel_size = float(pixel_size) self.voltage = voltage self.wavelength = voltage_to_wavelength(self.voltage) self.defocus_u = defocus_u @@ -454,7 +452,16 @@ def __init__( self._defocus_mean_nm = 0.05 * (self.defocus_u + self.defocus_v) self._defocus_diff_nm = 0.05 * (self.defocus_u - self.defocus_v) - def _evaluate(self, omega): + def _evaluate(self, omega, **kwargs): + # Ensure we have a pixel size, + pixel_size = kwargs.get("pixel_size", None) + if pixel_size is None: + raise RuntimeError( + f"{self.__class__.__name__}.evaluate must be passed kwarg `pixel_size`." + ) + # and that it is a floating point value. + pixel_size = float(pixel_size) + # Reference MATLAB code, includes reference to paper # Mindell, J. A.; Grigorieff, N. (2003). # https://github.com/PrincetonUniversity/aspire/blob/760a43b35453e55ff2d9354339e9ffa109a25371/projections/cryo_CTF_Relion.m#L34 @@ -478,7 +485,7 @@ def _evaluate(self, omega): # Divide by 10 to make pixel size in nm. BW is the # bandwidth of the signal corresponding to the given pixel size. - BW = 1 / (self.pixel_size / 10) + BW = 1 / (pixel_size / 10) s = s * BW DFavg = self._defocus_mean_nm # (DefocusU+DefocusV)/2 @@ -498,25 +505,10 @@ def _evaluate(self, omega): return h - def scale(self, c=1): - return CTFFilter( - pixel_size=self.pixel_size * c, - voltage=self.voltage, - defocus_u=self.defocus_u, - defocus_v=self.defocus_v, - defocus_ang=self.defocus_ang, - Cs=self.Cs, - alpha=self.alpha, - B=self.B, - ) - class RadialCTFFilter(CTFFilter): - def __init__( - self, pixel_size=1, voltage=200, defocus=15000, Cs=2.26, alpha=0.07, B=0 - ): + def __init__(self, voltage=200, defocus=15000, Cs=2.26, alpha=0.07, B=0): super().__init__( - pixel_size=pixel_size, voltage=voltage, defocus_u=defocus, defocus_v=defocus, @@ -539,7 +531,7 @@ def __init__(self, dim=None, var=1): def __repr__(self): return f"BlueFilter(dim={self.dim}, var={self.var})" - def _evaluate(self, omega): + def _evaluate(self, omega, **kwargs): f = np.sqrt(omega[0]) m = np.mean(f) f = f / m @@ -559,7 +551,7 @@ def __init__(self, dim=None, var=1): def __repr__(self): return f"PinkFilter(dim={self.dim}, var={self.var})" - def _evaluate(self, omega): + def _evaluate(self, omega, **kwargs): step = np.abs(np.subtract(*omega[0][:2])) # Avoid zero division f = np.sqrt(2 * step / (omega[0] + step)) diff --git a/src/aspire/source/coordinates.py b/src/aspire/source/coordinates.py index 7c8973792c..ff8bc5bfba 100644 --- a/src/aspire/source/coordinates.py +++ b/src/aspire/source/coordinates.py @@ -431,11 +431,9 @@ def _extract_ctf(self, data_block): logger.warning( "Unable to assign source pixel_size from CTFFilters, multiple pixel_sizes found." ) - # construct filters self.unique_filters = [ CTFFilter( - pixel_size=filter_params[i, 6], voltage=filter_params[i, 0], defocus_u=filter_params[i, 1], defocus_v=filter_params[i, 2], diff --git a/src/aspire/source/image.py b/src/aspire/source/image.py index b38c667748..196f60f344 100644 --- a/src/aspire/source/image.py +++ b/src/aspire/source/image.py @@ -751,7 +751,7 @@ def cache(self, batch_size=512): for start in trange(0, len(self), batch_size): end = min(start + batch_size, len(self)) im[start:end] = self.images[start:end] - self._cached_im = Image(im) + self._cached_im = Image(im, pixel_size=self.pixel_size) self.generation_pipeline.reset() @property diff --git a/src/aspire/source/micrograph.py b/src/aspire/source/micrograph.py index 5aa35e389c..66e3057775 100644 --- a/src/aspire/source/micrograph.py +++ b/src/aspire/source/micrograph.py @@ -309,7 +309,7 @@ def __init__( """ A cryo-EM MicrographSimulation object that supplies micrographs. - `dtype` and `particle_box_size` are inferred from `volume`, where `dtype` is the data type of the micrographs and `particle_box_size` is the size of the particle images. + `pixel_size`, `dtype` and `particle_box_size` are inferred from `volume`, where `dtype` is the data type of the micrographs and `particle_box_size` is the size of the particle images. :param volume: `Volume` instance to be used in `Simulation`. An `(L,L,L)` `Volume` will generate `(L,L)` particle images. @@ -339,9 +339,11 @@ def __init__( self.seed = seed + # Note pixel_size is taken from `volume`. super().__init__( micrograph_count=micrograph_count, micrograph_size=micrograph_size, + pixel_size=self.volume.pixel_size, dtype=self.volume.dtype, ) @@ -397,6 +399,7 @@ def __init__( self.ctf_filters = ctf_filters + # Note pixel_size is taken from `volume`. self.simulation = Simulation( n=self.total_particle_count, vols=self.volume, diff --git a/src/aspire/source/relion.py b/src/aspire/source/relion.py index 09bf77b8ea..a26b8cedec 100644 --- a/src/aspire/source/relion.py +++ b/src/aspire/source/relion.py @@ -153,7 +153,6 @@ def __init__( for row in filter_params: filters.append( CTFFilter( - pixel_size=self.pixel_size, voltage=row[0], defocus_u=row[1], defocus_v=row[2], diff --git a/src/aspire/source/simulation.py b/src/aspire/source/simulation.py index 0d83024584..10add38c88 100644 --- a/src/aspire/source/simulation.py +++ b/src/aspire/source/simulation.py @@ -156,7 +156,6 @@ def __init__( if unique_filters is None: unique_filters = [] self.unique_filters = unique_filters - self._check_filter_pixel_size(unique_filters) # sim_filters must be a deep copy so that it is not changed # when unique_filters is changed self.sim_filters = copy.deepcopy(unique_filters) @@ -245,29 +244,6 @@ def _populate_ctf_metadata(self, filter_indices): filter_values, ) - def _check_filter_pixel_size(self, unique_filters): - """ - Private method to ensure user provided filters match `Simulation` pixel size. - - When `Simulation.pixel_size` is not `None`, any - `unique_filters` having a non-matching `pixel_size` attribute - will raise. - """ - - # Skip when Simulation pixel_size is not explicitly provided. - if self.pixel_size is None: - return - - for f in unique_filters: - f_pixel_size = getattr(f, "pixel_size", None) - if f_pixel_size is not None and not np.isclose( - f_pixel_size, self.pixel_size - ): - raise ValueError( - f"`Simulation.pixel_size` {self.pixel_size} does not match filter {f} pixel size {f_pixel_size}." - "Ensure provided `pixel_size` attributes match." - ) - @property def projections(self): """ @@ -338,6 +314,7 @@ def _images(self, indices, clean_images=False): if not clean_images and self.noise_adder is not None: im = self.noise_adder.forward(im, indices=indices) + # scaling pixel_size in source, scaling filter, and scaling in IMage.downsample in conflict... # Finally, apply transforms to resulting Image return self.generation_pipeline.forward(im, indices) diff --git a/tests/test_anisotropic_noise.py b/tests/test_anisotropic_noise.py index 516e1652d8..419ef82f9c 100644 --- a/tests/test_anisotropic_noise.py +++ b/tests/test_anisotropic_noise.py @@ -15,14 +15,12 @@ class SimTestCase(TestCase): def setUp(self): self.dtype = np.float32 - self.vol = LegacyVolume(L=8, dtype=self.dtype).generate() + self.vol = LegacyVolume(L=8, pixel_size=10, dtype=self.dtype).generate() self.sim = _LegacySimulation( n=1024, vols=self.vol, unique_filters=[ - # Set legacy pixel size - RadialCTFFilter(pixel_size=10, defocus=d) - for d in np.linspace(1.5e4, 2.5e4, 7) + RadialCTFFilter(defocus=d) for d in np.linspace(1.5e4, 2.5e4, 7) ], dtype=self.dtype, ) diff --git a/tests/test_batched_covar2d.py b/tests/test_batched_covar2d.py index fbc6d3d7dc..f239354700 100644 --- a/tests/test_batched_covar2d.py +++ b/tests/test_batched_covar2d.py @@ -35,6 +35,7 @@ def setUp(self): L, n, unique_filters=self.filters, + pixel_size=5, dtype=self.dtype, noise_adder=noise_adder, ) @@ -253,7 +254,7 @@ class BatchedRotCov2DTestCaseCTF(BatchedRotCov2DTestCase): @property def filters(self): return [ - RadialCTFFilter(5, 200, defocus=d, Cs=2.0, alpha=0.1) + RadialCTFFilter(200, defocus=d, Cs=2.0, alpha=0.1) for d in np.linspace(1.5e4, 2.5e4, 7) ] @@ -263,4 +264,7 @@ def ctf_idx(self): @property def ctf_basis(self): - return [self.basis.filter_to_basis_mat(f) for f in self.src.unique_filters] + return [ + self.basis.filter_to_basis_mat(f, pixel_size=self.src.pixel_size) + for f in self.src.unique_filters + ] diff --git a/tests/test_coordinate_source.py b/tests/test_coordinate_source.py index c8e8704e69..38103829cc 100644 --- a/tests/test_coordinate_source.py +++ b/tests/test_coordinate_source.py @@ -584,7 +584,6 @@ def _testCtfFilters(self, src, uniform_pixel_sizes=True): 700.0, 600.0, 500.0, - self.pixel_size, ], dtype=src.dtype, ), @@ -596,7 +595,6 @@ def _testCtfFilters(self, src, uniform_pixel_sizes=True): filter0.Cs, filter0.alpha, filter0.voltage, - filter0.pixel_size, ] ), ) @@ -615,7 +613,6 @@ def _testCtfFilters(self, src, uniform_pixel_sizes=True): 701.0, 601.0, 501.0, - pixel_size1, ], dtype=src.dtype, ), @@ -627,7 +624,6 @@ def _testCtfFilters(self, src, uniform_pixel_sizes=True): filter1.Cs, filter1.alpha, filter1.voltage, - filter1.pixel_size, ] ), ) @@ -734,28 +730,6 @@ def testCommand(self): self.assertTrue(result_star.exit_code == 0) self.assertTrue(result_preprocess.exit_code == 0) - def testPixelSizeWarning(self): - """ - Test source having a pixel size that conflicts with the CTFFilter instances. - """ - manual_pixel_size = 0.789 - src = BoxesCoordinateSource(self.files_box, pixel_size=manual_pixel_size) - # Capture and compare warning message - with pytest.warns(UserWarning, match=r".*Pixel size mismatch.*"): - src.import_relion_ctf(self.relion_ctf_file) - np.testing.assert_approx_equal(src.pixel_size, manual_pixel_size) - - def testMultiplePixelSizeWarning(self): - """ - Test source having multiple pixel sizes in CTFFilter instances. - """ - src = BoxesCoordinateSource(self.files_box) # pixel_size=None - # Capture and compare warning message - with self._caplog.at_level(logging.WARNING): - src.import_aspire_ctf(self.ctf_files) # not uniform_pixel_sizes - assert src.pixel_size is None - assert "multiple pixel_sizes found" in self._caplog.text - def testPixelSize(self): """ Test explicitly providing correct pixel_size. @@ -764,14 +738,6 @@ def testPixelSize(self): src.import_relion_ctf(self.relion_ctf_file) np.testing.assert_approx_equal(src.pixel_size, self.pixel_size) - def testPixelSizeNone(self): - """ - Test not providing pixel_size. - """ - src = BoxesCoordinateSource(self.files_box) - src.import_relion_ctf(self.relion_ctf_file) - np.testing.assert_approx_equal(src.pixel_size, self.pixel_size) - def create_test_rectangular_micrograph_and_star(tmp_path, voxel_size=(2.0, 2.0, 1.0)): # Create a rectangular micrograph (e.g., 128x256) diff --git a/tests/test_covar2d.py b/tests/test_covar2d.py index f27bc61855..b1bf41e231 100644 --- a/tests/test_covar2d.py +++ b/tests/test_covar2d.py @@ -56,7 +56,8 @@ def img_size(request): def volume(dtype, img_size): # Get a volume v = Volume( - np.load(os.path.join(DATA_DIR, "clean70SRibosome_vol_down8.npy")).astype(dtype) + np.load(os.path.join(DATA_DIR, "clean70SRibosome_vol_down8.npy")).astype(dtype), + pixel_size=5.0 * 65 / 8, ) # 1e3 is hardcoded to match legacy test files. return v * 1.0e3 @@ -84,16 +85,17 @@ def cov2d_fixture(volume, basis, ctf_enabled): # Popluate CTF if ctf_enabled: unique_filters = [ - RadialCTFFilter( - 5.0 * 65 / volume.resolution, 200, defocus=d, Cs=2.0, alpha=0.1 - ) + RadialCTFFilter(200, defocus=d, Cs=2.0, alpha=0.1) for d in np.linspace(1.5e4, 2.5e4, 7) ] # Copied from simulation defaults to match legacy test files. h_idx = randi(len(unique_filters), n, seed=0) - 1 - h_ctf_fb = [basis.filter_to_basis_mat(f) for f in unique_filters] + h_ctf_fb = [ + basis.filter_to_basis_mat(f, pixel_size=volume.pixel_size) + for f in unique_filters + ] noise_adder = WhiteNoiseAdder(var=NOISE_VAR) @@ -107,7 +109,7 @@ def cov2d_fixture(volume, basis, ctf_enabled): dtype=volume.dtype, noise_adder=noise_adder, ) - sim.cache() + sim = sim.cache() cov2d = RotCov2D(basis) coef_clean = basis.evaluate_t(sim.projections[:]) diff --git a/tests/test_covar2d_denoiser.py b/tests/test_covar2d_denoiser.py index 79828ddfe3..8a4660c5c6 100644 --- a/tests/test_covar2d_denoiser.py +++ b/tests/test_covar2d_denoiser.py @@ -16,7 +16,7 @@ noise_adder = WhiteNoiseAdder(var=noise_var) pixel_size = 5 filters = [ - RadialCTFFilter(pixel_size, 200, defocus=d, Cs=2.0, alpha=0.1) + RadialCTFFilter(200, defocus=d, Cs=2.0, alpha=0.1) for d in np.linspace(1.5e4, 2.5e4, 7) ] @@ -182,11 +182,14 @@ def test_filter_to_basis_mat_ctf(coef, basis): } # Create a RadialCTFFilter - filt = RadialCTFFilter(pixel_size=1) + filt = RadialCTFFilter() # Apply the basis filter operator. # Note transpose because `apply` expects and returns column vectors. - coef_ftbm = (basis.filter_to_basis_mat(filt, truncate=False) @ coef.asnumpy().T).T + coef_ftbm = ( + basis.filter_to_basis_mat(filt, truncate=False, pixel_size=pixel_size) + @ coef.asnumpy().T + ).T # Apply evaluate->filter->expand manually imgs = coef.evaluate() diff --git a/tests/test_covar3d.py b/tests/test_covar3d.py index 8b8df4421a..98d0fdea13 100644 --- a/tests/test_covar3d.py +++ b/tests/test_covar3d.py @@ -24,7 +24,7 @@ class Covar3DTestCase(TestCase): @classmethod def setUpClass(cls): cls.dtype = np.float32 - cls.vols = LegacyVolume(L=8, dtype=cls.dtype).generate() + cls.vols = LegacyVolume(L=8, dtype=cls.dtype, pixel_size=1).generate() cls.sim = _LegacySimulation( n=1024, vols=cls.vols, @@ -38,7 +38,8 @@ def setUpClass(cls): cls.mean_estimator = MeanEstimator(cls.sim, basis=basis) cls.mean_est = Volume( - np.load(os.path.join(DATA_DIR, "mean_8_8_8.npy")).astype(cls.dtype) + np.load(os.path.join(DATA_DIR, "mean_8_8_8.npy")).astype(cls.dtype), + pixel_size=1, ) # Passing in a mean_kernel argument to the following constructor speeds up some calculations @@ -402,7 +403,7 @@ def testClustering(self): # TODO, alter refs after RCOPT complete eigs_est_trunc = np.moveaxis(eigs_est[:, :, :, : C - 1], -1, 0) - eigs_est_trunc = Volume(eigs_est_trunc) + eigs_est_trunc = Volume(eigs_est_trunc, pixel_size=1) lambdas_est_trunc = lambdas_est[: C - 1, : C - 1] diff --git a/tests/test_downsample.py b/tests/test_downsample.py index 2879fc3dd7..d41db2e616 100644 --- a/tests/test_downsample.py +++ b/tests/test_downsample.py @@ -222,6 +222,7 @@ def test_downsample_legacy(volume, res_ds): np.testing.assert_allclose(ims_ds_legacy, ims_ds_py, atol=1e-08) +@pytest.mark.xfail(reason="Issue #1318, double application of pixel_size scaling.") def test_simulation_relion_downsample(): """ Test that Simulation.downsample corresponds to RelionSource.downsample @@ -233,7 +234,7 @@ def test_simulation_relion_downsample(): defocus_ct = 7 ctf_filters = [ - RadialCTFFilter(pixel_size=1, defocus=d) + RadialCTFFilter(defocus=d) for d in np.linspace(defocus_min, defocus_max, defocus_ct) ] @@ -244,6 +245,7 @@ def test_simulation_relion_downsample(): C=1, unique_filters=ctf_filters, noise_adder=WhiteNoiseAdder.from_snr(snr=1), + pixel_size=1, ) src_ds = src.downsample(src.L // 2) @@ -253,10 +255,11 @@ def test_simulation_relion_downsample(): src.save(starfile) # Load Simulation source as a RelionSource - rln_src = RelionSource(starfile) + rln_src = RelionSource(starfile, pixel_size=1) # Downsample and test that images and attributes correspond to src_ds rln_src_ds = rln_src.downsample(src.L // 2) + np.testing.assert_allclose( src_ds.images[:], rln_src_ds.images[:], atol=utest_tolerance(src.dtype) ) diff --git a/tests/test_filters.py b/tests/test_filters.py index 6b84e29be4..b96bc6dff3 100644 --- a/tests/test_filters.py +++ b/tests/test_filters.py @@ -23,6 +23,7 @@ class SimTestCase(TestCase): test_filter = ArrayFilter(np.random.randn(8, 8)) + filter_eval_kwargs = dict() def setUp(self): self.dtype = np.float32 @@ -34,7 +35,7 @@ def tearDown(self): def testFunctionFilter(self): filt = FunctionFilter(lambda x, y: np.exp(-(x**2 + y**2) / 2)) - result = filt.evaluate(self.omega) + result = filt.evaluate(self.omega, **self.filter_eval_kwargs) self.assertEqual(result.shape, (256,)) self.assertTrue( np.allclose( @@ -56,13 +57,13 @@ def testZeroFilter(self): self.assertTrue(np.allclose(result, np.zeros(256))) def testIdentityFilter(self): - result = IdentityFilter().evaluate(self.omega) + result = IdentityFilter().evaluate(self.omega, **self.filter_eval_kwargs) # For all filters, we should get a 1d ndarray back on evaluate self.assertEqual(result.shape, (256,)) self.assertTrue(np.allclose(result, np.ones(256))) def testScalarFilter(self): - result = ScalarFilter(value=1.5).evaluate(self.omega) + result = ScalarFilter(value=1.5).evaluate(self.omega, **self.filter_eval_kwargs) self.assertEqual(result.shape, (256,)) self.assertTrue(np.allclose(result, np.repeat(1.5, 256))) @@ -71,7 +72,7 @@ def testPowerFilter(self): filter=FunctionFilter(lambda x, y: np.exp(-(x**2 + y**2) / 2)), power=0.5, ) - result = filt.evaluate(self.omega) + result = filt.evaluate(self.omega, **self.filter_eval_kwargs) self.assertEqual(result.shape, (256,)) self.assertTrue( np.allclose( @@ -89,35 +90,30 @@ def testPowerFilter(self): ) ) - def testCTFFilter(self): - filter = CTFFilter(defocus_u=1.5e4, defocus_v=1.5e4) - result = filter.evaluate(self.omega) - self.assertEqual(result.shape, (256,)) - def testScaledFilter(self): scale_value = 2.5 - result1 = self.test_filter.evaluate(self.omega) - # ScaledFilter scales the pixel size which cancels out - # a corresponding scaling in omega + result1 = self.test_filter.evaluate(self.omega, **self.filter_eval_kwargs) + filt2 = ScaledFilter(self.test_filter, scale_value) - result2 = filt2.evaluate(self.omega * scale_value) + result2 = filt2.evaluate(self.omega * scale_value, **self.filter_eval_kwargs) self.assertTrue(np.allclose(result1, result2, atol=utest_tolerance(self.dtype))) - def testRadialCTFFilter(self): - filter = RadialCTFFilter(defocus=2.5e4) - result = filter.evaluate(self.omega) - self.assertEqual(result.shape, (256,)) - def testDualFilter(self): - result = self.test_filter.evaluate(-self.omega) + result = self.test_filter.evaluate(-self.omega, **self.filter_eval_kwargs) dual_filter = self.test_filter.dual() - dual_result = dual_filter.evaluate(self.omega) + dual_result = dual_filter.evaluate(self.omega, **self.filter_eval_kwargs) self.assertTrue(np.allclose(result, dual_result)) def testFilterSigns(self): - signs = np.sign(self.test_filter.evaluate(self.omega)) + signs = np.sign( + self.test_filter.evaluate(self.omega, **self.filter_eval_kwargs) + ) sign_filter = self.test_filter.sign - self.assertTrue(np.allclose(sign_filter.evaluate(self.omega), signs)) + self.assertTrue( + np.allclose( + sign_filter.evaluate(self.omega, **self.filter_eval_kwargs), signs + ) + ) class SimTestCaseCTFFilter(SimTestCase): @@ -126,6 +122,27 @@ class SimTestCaseCTFFilter(SimTestCase): """ test_filter = CTFFilter() + filter_eval_kwargs = dict(pixel_size=1) + + def testCTFFilter(self): + filter = CTFFilter(defocus_u=1.5e4, defocus_v=1.5e4) + result = filter.evaluate(self.omega, **self.filter_eval_kwargs) + self.assertEqual(result.shape, (256,)) + + def testRadialCTFFilter(self): + filter = RadialCTFFilter(defocus=2.5e4) + result = filter.evaluate(self.omega, **self.filter_eval_kwargs) + self.assertEqual(result.shape, (256,)) + + def testCTFScale(self): + filt = CTFFilter(defocus_u=1.5e4, defocus_v=1.5e4) + result1 = filt.evaluate(self.omega, **self.filter_eval_kwargs) + scale_value = 2.5 + filt = filt.scale(scale_value) + # scaling a CTFFilter scales the pixel size which cancels out + # a corresponding scaling in omega + result2 = filt.evaluate(self.omega * scale_value, **self.filter_eval_kwargs) + self.assertTrue(np.allclose(result1, result2, atol=utest_tolerance(self.dtype))) DTYPES = [np.float32, np.float64] @@ -198,7 +215,6 @@ def test_ctf_reference(): Test CTFFilter against a MATLAB reference. """ fltr = CTFFilter( - pixel_size=4.56, voltage=200, defocus_u=10000, defocus_v=15000, @@ -206,7 +222,7 @@ def test_ctf_reference(): Cs=2.0, alpha=0.1, ) - h = fltr.evaluate_grid(5) + h = fltr.evaluate_grid(5, pixel_size=4.56) # Compare with MATLAB. Note DF converted to nm # >> n=5; V=200; DF1=1000; DF2=1500; theta=1.23; Cs=2.0; A=0.1; pxA=4.56; diff --git a/tests/test_indexed_source.py b/tests/test_indexed_source.py index c7b7161a1b..80ada30003 100644 --- a/tests/test_indexed_source.py +++ b/tests/test_indexed_source.py @@ -84,7 +84,6 @@ def test_filter_mapping(): defoci = np.linspace(1000, 25000, N // 2) ctf_filters = [ CTFFilter( - v.pixel_size, 200, defocus_u=defoci[d], defocus_v=defoci[-d], diff --git a/tests/test_mean_estimator.py b/tests/test_mean_estimator.py index 9b5e787c9f..c7083d4ef1 100644 --- a/tests/test_mean_estimator.py +++ b/tests/test_mean_estimator.py @@ -47,18 +47,16 @@ def dtype(request): @pytest.fixture(scope="module") def sim(L, dtype): - px_sz = 1.234 sim = Simulation( L=L, n=256, C=1, # single volume unique_filters=[ - RadialCTFFilter(defocus=d, pixel_size=px_sz) - for d in np.linspace(1.5e4, 2.5e4, 7) + RadialCTFFilter(defocus=d) for d in np.linspace(1.5e4, 2.5e4, 7) ], dtype=dtype, seed=SEED, - pixel_size=px_sz, + pixel_size=1.234, ) sim = sim.cache() # precompute images diff --git a/tests/test_micrograph_simulation.py b/tests/test_micrograph_simulation.py index 14f1212e0d..c7870d96f0 100644 --- a/tests/test_micrograph_simulation.py +++ b/tests/test_micrograph_simulation.py @@ -251,12 +251,8 @@ def test_sim_save(): Specifically tests interoperability with CentersCoordinateSource """ - v = AsymmetricVolume(L=16, C=1, dtype=np.float64).generate() - ctfs = [ - RadialCTFFilter( - pixel_size=4, voltage=200, defocus=15000, Cs=2.26, alpha=0.07, B=0 - ) - ] + v = AsymmetricVolume(L=16, C=1, pixel_size=4, dtype=np.float64).generate() + ctfs = [RadialCTFFilter(voltage=200, defocus=15000, Cs=2.26, alpha=0.07, B=0)] mg_sim = MicrographSimulation( volume=v, @@ -312,12 +308,8 @@ def test_save_overwrite(caplog): Specifically tests interoperability with CentersCoordinateSource """ - v = AsymmetricVolume(L=16, C=1, dtype=np.float64).generate() - ctfs = [ - RadialCTFFilter( - pixel_size=4, voltage=200, defocus=15000, Cs=2.26, alpha=0.07, B=0 - ) - ] + v = AsymmetricVolume(L=16, C=1, pixel_size=4, dtype=np.float64).generate() + ctfs = [RadialCTFFilter(voltage=200, defocus=15000, Cs=2.26, alpha=0.07, B=0)] mg_sim = MicrographSimulation( volume=v, diff --git a/tests/test_preprocess_pipeline.py b/tests/test_preprocess_pipeline.py index ebd6e1c593..83cc070930 100644 --- a/tests/test_preprocess_pipeline.py +++ b/tests/test_preprocess_pipeline.py @@ -33,6 +33,7 @@ def get_sim_object(L, dtype): RadialCTFFilter(defocus=d) for d in np.linspace(1.5e4, 2.5e4, 7) ], noise_adder=noise_adder, + pixel_size=1, dtype=dtype, ) return sim @@ -64,6 +65,7 @@ def testEmptyPhaseFlip(caplog): sim = Simulation( L=8, n=num_images, + pixel_size=1, dtype=np.float32, ) # assert we log a warning to the user diff --git a/tests/test_simulation.py b/tests/test_simulation.py index fde8ee9f98..f0b6a1688c 100644 --- a/tests/test_simulation.py +++ b/tests/test_simulation.py @@ -124,8 +124,7 @@ def setUp(self): L=self.L, vols=self.vols, unique_filters=[ - RadialCTFFilter(pixel_size=self._pixel_size, defocus=d) - for d in np.linspace(1.5e4, 2.5e4, 7) + RadialCTFFilter(defocus=d) for d in np.linspace(1.5e4, 2.5e4, 7) ], noise_adder=WhiteNoiseAdder(var=1), dtype=self.dtype, @@ -173,9 +172,7 @@ def testSimulationCached(self): vols=self.vols, offsets=self.sim.offsets, unique_filters=[ - # Set legacy pixel size - RadialCTFFilter(pixel_size=self._pixel_size, defocus=d) - for d in np.linspace(1.5e4, 2.5e4, 7) + RadialCTFFilter(defocus=d) for d in np.linspace(1.5e4, 2.5e4, 7) ], noise_adder=WhiteNoiseAdder(var=1), dtype=self.dtype, @@ -657,10 +654,15 @@ def test_cached_image_accessors(): Test the behavior of image caching. """ # Create a CTF - ctf = [RadialCTFFilter(pixel_size=5)] + ctf = [RadialCTFFilter()] # Create a Simulation with noise and `ctf` src = Simulation( - L=32, n=3, C=1, noise_adder=WhiteNoiseAdder(var=0.123), unique_filters=ctf + L=32, + n=3, + C=1, + noise_adder=WhiteNoiseAdder(var=0.123), + unique_filters=ctf, + pixel_size=5, ) # Cache the simulation cached_src = src.cache() @@ -761,51 +763,3 @@ def check_metadata(sim_src, relion_src): np.testing.assert_allclose( v, np.array(relion_src._metadata[k]).astype(type(v[0])) ) - - -def test_pixel_size(caplog): - """ - Check pixel size is instantiated properly and warnings occur if pixel - size is overridden. - """ - vol_px_sz = 10.0 - L = 8 - data = np.ones(L**3).reshape(L, L, L) - vol = Volume(data, pixel_size=vol_px_sz) - - # Ensure vol pixel size - np.testing.assert_array_equal(vol.pixel_size, vol_px_sz) - - # Generate Simulation and check pixel_size is inhereted from vol - sim = Simulation(vols=vol) - np.testing.assert_array_equal(sim.pixel_size, vol_px_sz) - - # Generate Simulation with provided pixel_size and check - # that vol.pixel_size is overridden. - caplog.clear() - caplog.set_level(logging.WARN) - - sim_px_sz = 5.0 - msg = ( - f"Overriding volume pixel size, {vol_px_sz}, with " - f"user provided pixel size of {sim_px_sz} angstrom." - ) - - assert msg not in caplog.text - - sim = Simulation(vols=vol, pixel_size=sim_px_sz) - - assert msg in caplog.text - np.testing.assert_array_equal(sim.pixel_size, sim_px_sz) - - -def test_mismatched_pixel_size(): - """ - Confirm raises error when explicit Simulation and CTFFilter pixel sizes mismatch. - """ - # Create a CTF with a pixel_size - filts = [RadialCTFFilter(pixel_size=5)] - - # Try to create a Simulation with a different pixel_size - with raises(ValueError, match=r"pixel_size.*does not match filter.*"): - _ = Simulation(L=8, n=1, C=1, pixel_size=10, unique_filters=filts) diff --git a/tests/test_weighted_mean_estimator.py b/tests/test_weighted_mean_estimator.py index 3f58fcdd00..96f11f2cab 100644 --- a/tests/test_weighted_mean_estimator.py +++ b/tests/test_weighted_mean_estimator.py @@ -57,6 +57,7 @@ def sim(L, dtype): unique_filters=[ RadialCTFFilter(defocus=d) for d in np.linspace(1.5e4, 2.5e4, 7) ], + pixel_size=1, dtype=dtype, seed=SEED, )