fix: write tests to avoid regressions

Resolves: #267.

Author: oesteban

nitransforms/tests/test_nonlinear.py

@@ -10,7 +10,6 @@
 
 import numpy as np
 import nibabel as nb
-from nitransforms.resampling import apply
 from nitransforms.base import TransformError
 from nitransforms.nonlinear import (
     BSplineFieldTransform,
@@ -21,6 +20,17 @@
 rng = np.random.default_rng()
 
 
+SOME_TEST_POINTS = np.array(
+    [
+        [0.0, 0.0, 0.0],
+        [1.0, 2.0, 3.0],
+        [10.0, -10.0, 5.0],
+        [-5.0, 7.0, -2.0],
+        [12.0, 0.0, -11.0],
+    ]
+)
+
+
 def test_displacements_init():
     identity1 = DenseFieldTransform(
         np.zeros((10, 10, 10, 3)),
@@ -65,17 +75,8 @@ def test_bsplines_references(testdata_path):
             testdata_path / "someones_bspline_coefficients.nii.gz"
         ).to_field()
 
-    with pytest.raises(TransformError):
-        apply(
-            BSplineFieldTransform(
-                testdata_path / "someones_bspline_coefficients.nii.gz"
-            ),
-            testdata_path / "someones_anatomy.nii.gz",
-        )
-
-    apply(
-        BSplineFieldTransform(testdata_path / "someones_bspline_coefficients.nii.gz"),
-        testdata_path / "someones_anatomy.nii.gz",
+    BSplineFieldTransform(
+        testdata_path / "someones_bspline_coefficients.nii.gz",
         reference=testdata_path / "someones_anatomy.nii.gz",
     )
 

nitransforms/tests/test_resampling.py

@@ -150,13 +150,13 @@ def test_apply_linear_transform(
 
 
 @pytest.mark.xfail(
-    reason="Disable while #266 is developed.",
+    reason="GH-267: disabled while debugging",
     strict=False,
 )
 @pytest.mark.parametrize("image_orientation", ["RAS", "LAS", "LPS", "oblique"])
 @pytest.mark.parametrize("sw_tool", ["itk", "afni"])
 @pytest.mark.parametrize("axis", [0, 1, 2, (0, 1), (1, 2), (0, 1, 2)])
-def test_displacements_field1(
+def test_apply_displacements_field1(
     tmp_path,
     get_testdata,
     get_testmask,
@@ -190,16 +190,17 @@ def test_displacements_field1(
     else:
         field.to_filename(xfm_fname)
 
-    xfm = nitnl.load(xfm_fname, fmt=sw_tool)
+    # xfm = nitnl.load(xfm_fname, fmt=sw_tool)
+    xfm = nitnl.DenseFieldTransform(fieldmap, reference=nii)
 
+    ants_output = tmp_path / "ants_brainmask.nii.gz"
     # Then apply the transform and cross-check with software
     cmd = APPLY_NONLINEAR_CMD[sw_tool](
         transform=os.path.abspath(xfm_fname),
         reference=tmp_path / "mask.nii.gz",
         moving=tmp_path / "mask.nii.gz",
-        output=tmp_path / "resampled_brainmask.nii.gz",
-        extra="",
-        # extra="--output-data-type uchar" if sw_tool == "itk" else "",
+        output=ants_output,
+        extra="--output-data-type uchar" if sw_tool == "itk" else "",
     )
 
     # skip test if command is not available on host
@@ -209,11 +210,13 @@ def test_displacements_field1(
 
     # resample mask
     exit_code = check_call([cmd], shell=True)
-    sw_moved_mask = nb.load("resampled_brainmask.nii.gz")
+    assert exit_code == 0
+    sw_moved_mask = nb.load(ants_output)
     nt_moved_mask = apply(xfm, msk, order=0)
 
-    # Calculate xor between both:
-    sw_mask = np.asanyarray(sw_moved_mask.dataobj, dtype=bool)
+    nt_moved_mask.to_filename(tmp_path / "nit_brainmask.nii.gz")
+
+    assert np.sqrt((diff**2).mean()) < RMSE_TOL_LINEAR
     brainmask = np.asanyarray(nt_moved_mask.dataobj, dtype=bool)
     percent_diff = (sw_mask != brainmask)[5:-5, 5:-5, 5:-5].sum() / brainmask.size
 
@@ -403,3 +406,33 @@ def test_apply_4d(serialize_4d):
     data = np.asanyarray(moved.dataobj)
     idxs = [tuple(np.argwhere(data[..., i])[0]) for i in range(nvols)]
     assert idxs == [(9 - i, 2, 2) for i in range(nvols)]
+
+
+@pytest.mark.xfail(
+    reason="GH-267: disabled while debugging",
+    strict=False,
+)
+def test_apply_bspline(tmp_path, testdata_path):
+    """Cross-check B-Splines and deformation field."""
+    os.chdir(str(tmp_path))
+
+    img_name = testdata_path / "someones_anatomy.nii.gz"
+    disp_name = testdata_path / "someones_displacement_field.nii.gz"
+    bs_name = testdata_path / "someones_bspline_coefficients.nii.gz"
+
+    bsplxfm = nitnl.BSplineFieldTransform(bs_name, reference=img_name)
+    dispxfm = nitnl.DenseFieldTransform(disp_name)
+
+    out_disp = apply(dispxfm, img_name)
+    out_bspl = apply(bsplxfm, img_name)
+
+    out_disp.to_filename("resampled_field.nii.gz")
+    out_bspl.to_filename("resampled_bsplines.nii.gz")
+
+    assert (
+        np.sqrt(
+            (out_disp.get_fdata(dtype="float32") - out_bspl.get_fdata(dtype="float32"))
+            ** 2
+        ).mean()
+        < 0.2
+    )