diff --git a/.github/workflows/workflow.yml b/.github/workflows/workflow.yml
index b2263a460b..69156dc0e0 100644
--- a/.github/workflows/workflow.yml
+++ b/.github/workflows/workflow.yml
@@ -157,7 +157,7 @@ jobs:
       run: |
         ASPIREDIR=${{ env.WORK_DIR }} python -c \
         "import aspire; print(aspire.config['ray']['temp_dir'])"
-        ASPIREDIR=${{ env.WORK_DIR }} python -m pytest --durations=50
+        ASPIREDIR=${{ env.WORK_DIR }} PYTHONWARNINGS=error python -m pytest --durations=50
     - name: Cache Data
       run: |
         ASPIREDIR=${{ env.WORK_DIR }} python -c \
diff --git a/src/aspire/abinitio/commonline_c3_c4.py b/src/aspire/abinitio/commonline_c3_c4.py
index 68469618dd..0e0ca76565 100644
--- a/src/aspire/abinitio/commonline_c3_c4.py
+++ b/src/aspire/abinitio/commonline_c3_c4.py
@@ -850,7 +850,13 @@ def cl_angles_to_ind(cl_angles, n_theta):
         thetas = np.mod(thetas, 2 * np.pi)
 
         # linear scale from [0,2*pi) to [0,n_theta).
-        return np.mod(np.round(thetas / (2 * np.pi) * n_theta), n_theta).astype(int)
+        ind = np.mod(np.round(thetas / (2 * np.pi) * n_theta), n_theta).astype(int)
+
+        # Return scalar for single value.
+        if ind.size == 1:
+            ind = ind.flat[0]
+
+        return ind
 
     @staticmethod
     def g_sync(rots, order, rots_gt):
diff --git a/src/aspire/abinitio/commonline_sync3n.py b/src/aspire/abinitio/commonline_sync3n.py
index 3c40eb3ac5..6efb429818 100644
--- a/src/aspire/abinitio/commonline_sync3n.py
+++ b/src/aspire/abinitio/commonline_sync3n.py
@@ -1,5 +1,6 @@
 import logging
 import os.path
+import warnings
 
 import numpy as np
 from numpy.linalg import norm
@@ -727,7 +728,13 @@ def fun(x, B, P, b, x0, A=A, a=a):
 
         # Calculate probabilities
         ln_f_ind, ln_f_arb = self._pairs_probabilities(Rijs, P**2, A, a, B, b, x0)
-        Pij = 1 / (1 + (1 - P) / P * np.exp(ln_f_arb - ln_f_ind))
+
+        with warnings.catch_warnings():
+            # For large values of (ln_f_arb - ln_f_ind), numpy exponential will overflow. We still
+            # get the intended result of Pij = 0, so we capture and ignore the overflow warning.
+            warnings.filterwarnings("ignore", r".*overflow encountered in exp.*")
+
+            Pij = 1 / (1 + (1 - P) / P * np.exp(ln_f_arb - ln_f_ind))
 
         # Fix singular output
         num_nan = np.sum(np.isnan(Pij))
diff --git a/src/aspire/nufft/cufinufft.py b/src/aspire/nufft/cufinufft.py
index 218fbd5fb7..fd869aacfd 100644
--- a/src/aspire/nufft/cufinufft.py
+++ b/src/aspire/nufft/cufinufft.py
@@ -107,7 +107,8 @@ def transform(self, signal):
                 " In the future this will be an error."
             )
 
-        signal = cp.asarray(signal, dtype=self.complex_dtype)
+        # Note, if not C order, cuFINUFFT will copy-cast anyway.
+        signal = cp.asarray(signal, order="C", dtype=self.complex_dtype)
 
         sig_shape = signal.shape
         res_shape = self.num_pts
diff --git a/src/aspire/numeric/scipy.py b/src/aspire/numeric/scipy.py
index 8f2e7d86d8..c913e917a3 100644
--- a/src/aspire/numeric/scipy.py
+++ b/src/aspire/numeric/scipy.py
@@ -8,10 +8,11 @@
 
 def cg(*args, **kwargs):
     """
-    Supports scipy cg before and after 1.14.0.
+    Supports scipy cg before and after 1.12.0.
     """
 
-    # older scipy cg interface uses `tol` instead of `rtol`
-    if Version(scipy.__version__) < Version("1.14.0"):
+    # older (<1.12.0) scipy cg interface uses `tol` instead of `rtol`.
+    # `tol` will be removed in scipy 1.14.0.
+    if Version(scipy.__version__) < Version("1.12.0"):
         kwargs["tol"] = kwargs.pop("rtol", None)
     return scipy.sparse.linalg.cg(*args, **kwargs)
diff --git a/src/aspire/sinogram/sinogram.py b/src/aspire/sinogram/sinogram.py
index 34451a396d..7c7bb43662 100644
--- a/src/aspire/sinogram/sinogram.py
+++ b/src/aspire/sinogram/sinogram.py
@@ -42,7 +42,7 @@ def __init__(self, data, dtype=None):
         self.shape = self._data.shape
         self.stack_shape = self._data.shape[:-2]
         self.stack_n_dim = self._data.ndim - 2
-        self.n = np.product(self.stack_shape)
+        self.n = np.prod(self.stack_shape)
         self.n_angles = self._data.shape[-2]
         self.n_radial_points = self._data.shape[-1]
 
diff --git a/src/aspire/utils/resolution_estimation.py b/src/aspire/utils/resolution_estimation.py
index 14c142bdc5..af5a46f8a9 100644
--- a/src/aspire/utils/resolution_estimation.py
+++ b/src/aspire/utils/resolution_estimation.py
@@ -354,7 +354,8 @@ def plot(self, cutoff=None, save_to_file=False, labels=None):
         plt.ylabel("Correlation")
         plt.ylim([0, 1.1])
         for i, line in enumerate(self.correlations):
-            _label = None
+            # Set default label for single correlation (required by plt.legend() below).
+            _label = "correlation"
             if len(self.correlations) > 1:
                 _label = f"{i}"
                 if labels is not None:
diff --git a/src/aspire/utils/rotation.py b/src/aspire/utils/rotation.py
index 08bec4ca3d..07a31df9df 100644
--- a/src/aspire/utils/rotation.py
+++ b/src/aspire/utils/rotation.py
@@ -408,6 +408,10 @@ def angle_dist(r1, r2, dtype=None):
         theta = (tr_r[non_zero_dist_ind] - 1) / 2
         theta = np.maximum(np.minimum(theta, 1), -1)  # Clamp theta in [-1,1]
         dist[non_zero_dist_ind] = np.arccos(theta, dtype=dtype)
+
+        # Return scalar for single value.
+        if dist.size == 1:
+            dist = dist.flat[0]
         return dist
 
     @staticmethod
diff --git a/tests/test_FLEbasis2D.py b/tests/test_FLEbasis2D.py
index ffb1f8f7d1..0873d61bab 100644
--- a/tests/test_FLEbasis2D.py
+++ b/tests/test_FLEbasis2D.py
@@ -230,22 +230,23 @@ def testLowPass():
 def testRadialConvolution():
     # test ability to accurately convolve with a radial
     # (e.g. CTF) function via FLE coefficients
-
     L = 32
-    basis = FLEBasis2D(L, match_fb=False)
+
     # load test radial function
     x = np.load(os.path.join(DATA_DIR, "fle_radial_fn_32x32.npy")).reshape(1, 32, 32)
     x = x / np.max(np.abs(x.flatten()))
 
     # get sample images
     ims = create_images(L, 10)
+
     # convolve using coefficients
+    basis = FLEBasis2D(L, match_fb=False, dtype=ims.dtype)
     coefs = basis.evaluate_t(ims)
     coefs_convolved = basis.radial_convolve(coefs, x)
     imgs_convolved_fle = basis.evaluate(coefs_convolved).asnumpy()
 
     # convolve using FFT
-    x = basis.evaluate(basis.evaluate_t(x)).asnumpy()
+    x = basis.evaluate(basis.evaluate_t(Image(x))).asnumpy()
     ims = basis.evaluate(coefs).asnumpy()
 
     imgs_convolved_slow = np.zeros((10, L, L))
diff --git a/tests/test_covar2d_denoiser.py b/tests/test_covar2d_denoiser.py
index 7b4da5511e..ea5410fb34 100644
--- a/tests/test_covar2d_denoiser.py
+++ b/tests/test_covar2d_denoiser.py
@@ -18,12 +18,28 @@
     RadialCTFFilter(5, 200, defocus=d, Cs=2.0, alpha=0.1)
     for d in np.linspace(1.5e4, 2.5e4, 7)
 ]
+
+# For (F)PSWFBasis2D we get off-block entries which are truncated
+# when converting to block-diagonal. We filter these warnings.
 BASIS = [
     pytest.param(FBBasis2D, marks=pytest.mark.expensive),
     FFBBasis2D,
     FLEBasis2D,
-    pytest.param(PSWFBasis2D, marks=pytest.mark.expensive),
-    FPSWFBasis2D,
+    pytest.param(
+        PSWFBasis2D,
+        marks=[
+            pytest.mark.expensive,
+            pytest.mark.filterwarnings(
+                "ignore:BlkDiagMatrix.from_dense truncating values*"
+            ),
+        ],
+    ),
+    pytest.param(
+        FPSWFBasis2D,
+        marks=pytest.mark.filterwarnings(
+            "ignore:BlkDiagMatrix.from_dense truncating values*"
+        ),
+    ),
 ]
 
 
diff --git a/tests/test_diag_matrix.py b/tests/test_diag_matrix.py
index ecce899105..05805912c9 100644
--- a/tests/test_diag_matrix.py
+++ b/tests/test_diag_matrix.py
@@ -77,7 +77,7 @@ def test_repr():
     Test accessing the `repr` does not crash.
     """
 
-    d = DiagMatrix(np.empty((10, 8)))
+    d = DiagMatrix(np.ones((10, 8)))
     assert repr(d).startswith("DiagMatrix(")
 
 
@@ -86,7 +86,7 @@ def test_str():
     Test accessing the `str` does not crash.
     """
 
-    d = DiagMatrix(np.empty((10, 8)))
+    d = DiagMatrix(np.ones((10, 8)))
     assert str(d).startswith("DiagMatrix(")
 
 
@@ -104,13 +104,13 @@ def test_len():
     """
     Test the `len`.
     """
-    d = DiagMatrix(np.empty((10, 8)))
+    d = DiagMatrix(np.ones((10, 8)))
 
     assert d.size == 10
     assert d.count == 8
     assert len(d) == 10
 
-    d = DiagMatrix(np.empty((2, 5, 8)))
+    d = DiagMatrix(np.ones((2, 5, 8)))
 
     assert d.size == 10
     assert d.count == 8
@@ -121,8 +121,8 @@ def test_size_mismatch():
     """
     Test we raise operating on `DiagMatrix` having different counts.
     """
-    d1 = DiagMatrix(np.empty((10, 8)))
-    d2 = DiagMatrix(np.empty((10, 7)))
+    d1 = DiagMatrix(np.ones((10, 8)))
+    d2 = DiagMatrix(np.ones((10, 7)))
 
     with pytest.raises(RuntimeError, match=r".*not same dimension.*"):
         _ = d1 + d2
@@ -132,8 +132,8 @@ def test_dtype_mismatch():
     """
     Test we raise operating on `DiagMatrix` having different dtypes.
     """
-    d1 = DiagMatrix(np.empty((10, 8)), dtype=np.float32)
-    d2 = DiagMatrix(np.empty((10, 8)), dtype=np.float64)
+    d1 = DiagMatrix(np.ones((10, 8)), dtype=np.float32)
+    d2 = DiagMatrix(np.ones((10, 8)), dtype=np.float64)
 
     with pytest.raises(RuntimeError, match=r".*received different types.*"):
         _ = d1 + d2
@@ -144,7 +144,7 @@ def test_dtype_passthrough():
     Test that the datatype is inferred correctly.
     """
     for dtype in (int, np.float32, np.float64, np.complex64, np.complex128):
-        d_np = np.empty(42, dtype=dtype)
+        d_np = np.ones(42, dtype=dtype)
         d = DiagMatrix(d_np)
         assert d.dtype == dtype
 
@@ -154,7 +154,7 @@ def test_dtype_cast():
     Test that a datatype is cast when overridden.
     """
     for dtype in (int, np.float32, np.float64, np.complex64, np.complex128):
-        d_np = np.empty(42, dtype=np.float16)
+        d_np = np.ones(42, dtype=np.float16)
         d = DiagMatrix(d_np, dtype)
         assert d.dtype == dtype
 
@@ -444,7 +444,7 @@ def test_diag_badtype_matmul():
     """
     Test matrix multiply of `DiagMatrix` with incompatible type raises.
     """
-    d1 = DiagMatrix(np.empty(8))
+    d1 = DiagMatrix(np.ones(8))
 
     # matmul
     with pytest.raises(RuntimeError, match=r".*not implemented for.*"):
@@ -576,7 +576,7 @@ def test_bad_as_blk_diag(matrix_size, blk_diag):
     """
     with pytest.raises(RuntimeError, match=r".*only implemented for singletons.*"):
         # Construct via Numpy.
-        d_np = np.empty((2, matrix_size), dtype=blk_diag.dtype)
+        d_np = np.ones((2, matrix_size), dtype=blk_diag.dtype)
 
         # Create DiagMatrix then convert to BlkDiagMatrix
         d = DiagMatrix(d_np)
@@ -654,7 +654,7 @@ def test_diag_blk_mul():
     """
     Test mixing `BlkDiagMatrix` with `DiagMatrix` element-wise multiplication raises.
     """
-    d = DiagMatrix(np.empty(8))
+    d = DiagMatrix(np.ones(8))
 
     partition = [(4, 4), (4, 4)]
     b = BlkDiagMatrix.ones(partition, dtype=d.dtype)
@@ -672,7 +672,7 @@ def test_non_square_as_blk_diag():
     """
     Test non square partition blocks raise an error in as_blk_diag.
     """
-    d = DiagMatrix(np.empty(8))
+    d = DiagMatrix(np.ones(8))
 
     partition = [(4, 5), (4, 3)]
     with pytest.raises(RuntimeError, match=r".*not square.*"):
@@ -683,8 +683,8 @@ def test_bad_broadcast():
     """
     Test incompatible stack shapes raise appropriate error.
     """
-    d1 = DiagMatrix(np.empty((2, 3, 8)))
-    d2 = DiagMatrix(np.empty((2, 2, 8)))
+    d1 = DiagMatrix(np.ones((2, 3, 8)))
+    d2 = DiagMatrix(np.ones((2, 2, 8)))
 
     with pytest.raises(ValueError, match=r".*incompatible shapes.*"):
         _ = d1 + d2
diff --git a/tests/test_image.py b/tests/test_image.py
index 887e726c0d..89fbde4a84 100644
--- a/tests/test_image.py
+++ b/tests/test_image.py
@@ -353,34 +353,34 @@ def test_asnumpy_readonly():
         vw[0, 0, 0] = 123
 
 
-@pytest.mark.xfail(reason="Ray logging issue ray#37711", strict=False)
 def test_corrupt_mrc_load(caplog):
     """
     Test that corrupt mrc files are logged as expected.
     """
 
-    caplog.set_level(logging.WARNING)
-
     # Create a tmp dir for this test output
     with tempfile.TemporaryDirectory() as tmpdir_name:
         # tmp filename
         mrc_path = os.path.join(tmpdir_name, "bad.mrc")
 
         # Create and save image
-        Image(np.empty((1, 8, 8), dtype=np.float32)).save(mrc_path)
+        Image(np.ones((1, 8, 8), dtype=np.float32)).save(mrc_path)
 
         # Open mrc file and soft corrupt it
         with mrcfile.open(mrc_path, "r+") as fh:
             fh.header.map = -1
 
         # Check that we get a WARNING
-        _ = Image.load(mrc_path)
+        with caplog.at_level(logging.WARNING):
+            _ = Image.load(mrc_path)
+
+            # Check the message prefix
+            assert f"Image.load of {mrc_path} reporting 1 corruptions" in caplog.text
 
-    # Check the message prefix
-    assert f"Image.load of {mrc_path} reporting 1 corruptions" in caplog.text
+            # Check the message contains the file path
+            assert mrc_path in caplog.text
 
-    # Check the message contains the file path
-    assert mrc_path in caplog.text
+            caplog.clear()
 
 
 def test_load_bad_ext():
diff --git a/tests/test_micrograph_source.py b/tests/test_micrograph_source.py
index 06c5874a07..d4793cf61f 100644
--- a/tests/test_micrograph_source.py
+++ b/tests/test_micrograph_source.py
@@ -285,7 +285,7 @@ def test_rectangular_micrograph_source_files():
     """
 
     # Test inconsistent mrc files
-    imgs = [np.empty((7, 7)), np.empty((8, 8))]
+    imgs = [np.zeros((7, 7)), np.zeros((8, 8))]
     with tempfile.TemporaryDirectory() as tmp_output_dir:
         # Save the files
         for i, img in enumerate(imgs):
diff --git a/tests/test_rotation.py b/tests/test_rotation.py
index 9e0dba4ec6..e02e650bd5 100644
--- a/tests/test_rotation.py
+++ b/tests/test_rotation.py
@@ -173,6 +173,9 @@ def test_angle_dist(dtype):
     with pytest.raises(ValueError, match=r"r1 and r2 are not broadcastable*"):
         _ = Rotation.angle_dist(rots[:3], rots[:5])
 
+    # Test that single value returns as 0-dim.
+    assert Rotation.angle_dist(rots[0], rots[1], dtype).ndim == 0
+
 
 def test_mean_angular_distance(dtype):
     rots_z = Rotation.about_axis("z", [0, np.pi / 4, np.pi / 2], dtype=dtype).matrices
diff --git a/tests/test_sinogram.py b/tests/test_sinogram.py
index bbf448db97..dc173449b9 100644
--- a/tests/test_sinogram.py
+++ b/tests/test_sinogram.py
@@ -113,7 +113,7 @@ def test_project_multidim(num_ang):
 
     # Generate a mask
     g = grid_2d(L, normalized=True, shifted=True)
-    mask = g["r"] < 1
+    mask = g["r"] < 0.99
 
     # Generate images
     imgs = Image(np.random.random((m, n, L, L))) * mask
diff --git a/tests/test_utils.py b/tests/test_utils.py
index ffad5bc9f6..c7617d28a8 100644
--- a/tests/test_utils.py
+++ b/tests/test_utils.py
@@ -399,6 +399,9 @@ def matplotlib_no_gui():
 
         yield
 
+    # Explicitly close all figures before making backend changes.
+    matplotlib.pyplot.close("all")
+
     # Restore backend
     matplotlib.use(backend)
 
diff --git a/tests/test_volume.py b/tests/test_volume.py
index ac86c4096b..1f55645e10 100644
--- a/tests/test_volume.py
+++ b/tests/test_volume.py
@@ -846,12 +846,6 @@ def test_aglebraic_ops_symmetry_warnings(symmetric_vols):
     # Should have 4 warnings on record.
     assert len(record) == 4
 
-    # Check that warning occurs only once per line.
-    with warnings.catch_warnings(record=True) as record:
-        for _ in range(5):
-            vol_c3 + vol_c4
-    assert len(record) == 1
-
 
 def test_volume_load_with_symmetry():
     # Check we can load a Volume with symmetry_group.