test: added dtype checks for radon

Author: mrava87

pytests/test_radon.py

@@ -132,15 +132,16 @@ def test_unknown_engine():
 @pytest.mark.parametrize(
     "par", [(par1), (par2), (par3), (par4), (par5), (par6), (par7), (par8)]
 )
-def test_Radon2D(par):
+@pytest.mark.parametrize("dtype", [np.float32, np.float64])
+def test_Radon2D(par, dtype):
     """Dot-test, forward and adjoint consistency check
     (for onthefly parameter), and sparse inverse for Radon2D operator
     """
     dt, dh = 0.005, 1
-    t = np.arange(par["nt"]) * dt
-    h = np.arange(par["nhx"]) * dh
-    px = np.linspace(0, par["pxmax"], par["npx"])
-    x = np.zeros((par["npx"], par["nt"]))
+    t = np.arange(par["nt"], dtype=dtype) * dt
+    h = np.arange(par["nhx"], dtype=dtype) * dh
+    px = np.linspace(0, par["pxmax"], par["npx"], dtype=dtype)
+    x = np.zeros((par["npx"], par["nt"]), dtype=dtype)
     x[2, par["nt"] // 2] = 1
 
     Rop = Radon2D(
@@ -152,7 +153,7 @@ def test_Radon2D(par):
         kind=par["kind"],
         onthefly=False,
         engine=par["engine"],
-        dtype="float64",
+        dtype=dtype,
     )
     R1op = Radon2D(
         t,
@@ -163,82 +164,96 @@ def test_Radon2D(par):
         kind=par["kind"],
         onthefly=True,
         engine=par["engine"],
-        dtype="float64",
+        dtype=dtype,
+    )
+    assert dottest(
+        Rop,
+        par["nhx"] * par["nt"],
+        par["npx"] * par["nt"],
+        rtol=1e-3 if dtype == np.float32 else 1e-6,
     )
-    assert dottest(Rop, par["nhx"] * par["nt"], par["npx"] * par["nt"], rtol=1e-3)
 
     y = Rop * x.ravel()
     y1 = R1op * x.ravel()
+    assert y.dtype == dtype
+    assert y1.dtype == dtype
     assert_array_almost_equal(y, y1, decimal=4)
 
     xadj = Rop.H * y
     xadj1 = R1op.H * y
+    assert xadj.dtype == dtype
+    assert xadj1.dtype == dtype
     assert_array_almost_equal(xadj, xadj1, decimal=4)
 
     xinv, _, _ = fista(Rop, y, niter=30, eps=1e0)
     assert_array_almost_equal(x.ravel(), xinv, decimal=1)
 
 
-# @pytest.mark.skipif(
-#     int(os.environ.get("TEST_CUPY_PYLOPS", 0)) == 1, reason="Not CuPy enabled"
-# )
-# @pytest.mark.parametrize(
-#     "par", [(par1), (par2), (par3), (par4), (par5), (par6), (par7), (par8)]
-# )
-# def test_Radon3D(par):
-#     """Dot-test, forward and adjoint consistency check
-#     (for onthefly parameter), and sparse inverse for Radon3D operator
-#     """
-#     dt, dhy, dhx = 0.005, 1, 1
-#     t = np.arange(par["nt"]) * dt
-#     hy = np.arange(par["nhy"]) * dhy
-#     hx = np.arange(par["nhx"]) * dhx
-#     py = np.linspace(0, par["pymax"], par["npy"])
-#     px = np.linspace(0, par["pxmax"], par["npx"])
-#     x = np.zeros((par["npy"], par["npx"], par["nt"]))
-#     x[3, 2, par["nt"] // 2] = 1
-
-#     Rop = Radon3D(
-#         t,
-#         hy,
-#         hx,
-#         py,
-#         px,
-#         centeredh=par["centeredh"],
-#         interp=par["interp"],
-#         kind=par["kind"],
-#         onthefly=False,
-#         engine=par["engine"],
-#         dtype="float64",
-#     )
-#     R1op = Radon3D(
-#         t,
-#         hy,
-#         hx,
-#         py,
-#         px,
-#         centeredh=par["centeredh"],
-#         interp=par["interp"],
-#         kind=par["kind"],
-#         onthefly=True,
-#         engine=par["engine"],
-#         dtype="float64",
-#     )
-
-#     assert dottest(
-#         Rop,
-#         par["nhy"] * par["nhx"] * par["nt"],
-#         par["npy"] * par["npx"] * par["nt"],
-#         rtol=1e-3,
-#     )
-#     y = Rop * x.ravel()
-#     y1 = R1op * x.ravel()
-#     assert_array_almost_equal(y, y1, decimal=4)
-
-#     xadj = Rop.H * y
-#     xadj1 = R1op.H * y
-#     assert_array_almost_equal(xadj, xadj1, decimal=4)
-
-#     if Rop.engine == "numba":  # as numpy is too slow here...
-#         xinv, _, _ = fista(Rop, y, niter=200, eps=3e0)
-#         assert_array_almost_equal(x.ravel(), xinv, decimal=1)
+@pytest.mark.skipif(
+    int(os.environ.get("TEST_CUPY_PYLOPS", 0)) == 1, reason="Not CuPy enabled"
+)
+@pytest.mark.parametrize(
+    "par", [(par1), (par2), (par3), (par4), (par5), (par6), (par7), (par8)]
+)
+@pytest.mark.parametrize("dtype", [np.float32, np.float64])
+def test_Radon3D(par, dtype):
+    """Dot-test, forward and adjoint consistency check
+    (for onthefly parameter), and sparse inverse for Radon3D operator
+    """
+    dt, dhy, dhx = 0.005, 1, 1
+    t = np.arange(par["nt"], dtype=dtype) * dt
+    hy = np.arange(par["nhy"], dtype=dtype) * dhy
+    hx = np.arange(par["nhx"], dtype=dtype) * dhx
+    py = np.linspace(0, par["pymax"], par["npy"], dtype=dtype)
+    px = np.linspace(0, par["pxmax"], par["npx"], dtype=dtype)
+    x = np.zeros((par["npy"], par["npx"], par["nt"]), dtype=dtype)
+    x[3, 2, par["nt"] // 2] = 1
+
+    Rop = Radon3D(
+        t,
+        hy,
+        hx,
+        py,
+        px,
+        centeredh=par["centeredh"],
+        interp=par["interp"],
+        kind=par["kind"],
+        onthefly=False,
+        engine=par["engine"],
+        dtype=dtype,
+    )
+    R1op = Radon3D(
+        t,
+        hy,
+        hx,
+        py,
+        px,
+        centeredh=par["centeredh"],
+        interp=par["interp"],
+        kind=par["kind"],
+        onthefly=True,
+        engine=par["engine"],
+        dtype=dtype,
+    )
+
+    assert dottest(
+        Rop,
+        par["nhy"] * par["nhx"] * par["nt"],
+        par["npy"] * par["npx"] * par["nt"],
+        rtol=1e-3 if dtype == np.float32 else 1e-6,
+    )
+    y = Rop * x.ravel()
+    y1 = R1op * x.ravel()
+    assert y.dtype == dtype
+    assert y1.dtype == dtype
+    assert_array_almost_equal(y, y1, decimal=4)
+
+    xadj = Rop.H * y
+    xadj1 = R1op.H * y
+    assert xadj.dtype == dtype
+    assert xadj1.dtype == dtype
+    assert_array_almost_equal(xadj, xadj1, decimal=4)
+
+    if Rop.engine == "numba":  # as numpy is too slow here...
+        xinv, _, _ = fista(Rop, y, niter=200, eps=3e0)
+        assert_array_almost_equal(x.ravel(), xinv, decimal=1)