fix: in reality, there's nothing we must do :D

Author: oesteban

docs/notebooks/02_afni_deoblique.ipynb

@@ -2,7 +2,7 @@
 from math import pi
 import numpy as np
 from nibabel.affines import (
-    from_matvec,
+    apply_affine,
     obliquity,
     voxel_sizes,
 )
@@ -39,25 +39,8 @@ def to_string(self, banner=True):
     @classmethod
     def from_ras(cls, ras, moving=None, reference=None):
         """Create an AFNI affine from a nitransform's RAS+ matrix."""
-        pre = LPS
-        post = LPS
-
-        if reference is not None:
-            reference = _ensure_image(reference)
-
-        if reference is not None and _is_oblique(reference.affine):
-            print("Reference affine axes are oblique.")
-            ras = ras @ _afni_warpdrive(reference.affine, ras=True, forward=False)
-
-        if moving is not None:
-            moving = _ensure_image(moving)
-
-        if moving is not None and _is_oblique(moving.affine):
-            print("Moving affine axes are oblique.")
-            ras = _afni_warpdrive(reference.affine, ras=True) @ ras
-
         # swapaxes is necessary, as axis 0 encodes series of transforms
-        parameters = np.swapaxes(post @ ras @ pre, 0, 1)
+        parameters = np.swapaxes(LPS @ ras @ LPS, 0, 1)
 
         tf = cls()
         tf.structarr["parameters"] = parameters.T
@@ -89,23 +72,8 @@ def from_string(cls, string):
 
     def to_ras(self, moving=None, reference=None):
         """Return a nitransforms internal RAS+ matrix."""
-        pre = LPS
-        post = LPS
-
-        if reference is not None:
-            reference = _ensure_image(reference)
-
-        if reference is not None and _is_oblique(reference.affine):
-            raise NotImplementedError
-
-        if moving is not None:
-            moving = _ensure_image(moving)
-
-        if moving is not None and _is_oblique(moving.affine):
-            raise NotImplementedError
-
         # swapaxes is necessary, as axis 0 encodes series of transforms
-        return post @ np.swapaxes(self.structarr["parameters"].T, 0, 1) @ pre
+        return LPS @ np.swapaxes(self.structarr["parameters"].T, 0, 1) @ LPS
 
 
 class AFNILinearTransformArray(BaseLinearTransformList):
@@ -183,7 +151,7 @@ def _is_oblique(affine, thres=OBLIQUITY_THRESHOLD_DEG):
     return (obliquity(affine).min() * 180 / pi) > thres
 
 
-def _afni_warpdrive(nii, forward=True, ras=False):
+def _afni_warpdrive(oblique, shape, forward=True, ras=False):
     """
     Calculate AFNI's ``WARPDRIVE_MATVEC_FOR_000000`` (de)obliquing affine.
 
@@ -204,15 +172,18 @@ def _afni_warpdrive(nii, forward=True, ras=False):
         to be oblique.
 
     """
-    oblique = nii.affine
-    plumb = oblique[:3, :3] / np.abs(oblique[:3, :3]).max(0)
-    plumb[np.abs(plumb) < 1.0] = 0
-    plumb *= voxel_sizes(oblique)
-
-    R = from_matvec(plumb @ np.linalg.inv(oblique[:3, :3]), (0, 0, 0))
-    plumb_orig = np.linalg.inv(R[:3, :3]) @ oblique[:3, 3]
-    print(plumb_orig)
-    R[:3, 3] = R[:3, :3] @ (plumb_orig - oblique[:3, 3])
+    shape = np.array(shape[:3])
+    plumb_r = oblique[:3, :3] / np.abs(oblique[:3, :3]).max(0)
+    plumb_r[np.abs(plumb_r) < 1.0] = 0
+    plumb_r *= voxel_sizes(oblique)
+    plumb = np.eye(4)
+    plumb[:3, :3] = plumb_r
+    obliq_o = apply_affine(oblique, 0.5 * (shape - 1))
+    plumb_c = apply_affine(plumb, 0.5 * (shape - 1))
+    plumb[:3, 3] = -plumb_c + obliq_o
+    print(obliq_o, apply_affine(plumb, 0.5 * (shape - 1)))
+
+    R = plumb @ np.linalg.inv(oblique)
     if not ras:
         # Change sign to match AFNI's warpdrive_matvec signs
         B = np.ones((2, 2))
@@ -228,5 +199,8 @@ def _afni_header(nii, field="WARPDRIVE_MATVEC_FOR_000000"):
         root.find(f".//*[@atr_name='{field}']").text,
         sep="\n",
         dtype="float32"
-    ).reshape((3, 4))
-    return np.vstack((retval, (0, 0, 0, 1)))
+    )
+    if retval.size == 12:
+        return np.vstack((retval.reshape((3, 4)), (0, 0, 0, 1)))
+
+    return retval

nitransforms/io/afni.py

@@ -0,0 +1 @@
+affine-RAS.afni

nitransforms/tests/data/affine-oblique.afni

@@ -176,9 +176,6 @@ def test_LT_conversions(data_path, fname):
 )
 @pytest.mark.parametrize("sw", ["afni", "fsl", "fs", "itk"])
 def test_Linear_common(tmpdir, data_path, sw, image_orientation, get_testdata):
-    if (image_orientation, sw) == ("oblique", "afni"):
-        pytest.skip("AFNI Deoblique unsupported.")
-
     tmpdir.chdir()
 
     moving = get_testdata[image_orientation]

nitransforms/tests/test_io.py

@@ -34,7 +34,7 @@
 }
 
 
-@pytest.mark.parametrize("matrix", [[0.0], np.ones((3, 3, 3)), np.ones((3, 4)),])
+@pytest.mark.parametrize("matrix", [[0.0], np.ones((3, 3, 3)), np.ones((3, 4)), ])
 def test_linear_typeerrors1(matrix):
     """Exercise errors in Affine creation."""
     with pytest.raises(TypeError):
@@ -141,9 +141,6 @@ def test_loadsave(tmp_path, data_path, testdata_path, fmt):
 @pytest.mark.parametrize("sw_tool", ["itk", "fsl", "afni", "fs"])
 def test_linear_save(tmpdir, data_path, get_testdata, image_orientation, sw_tool):
     """Check implementation of exporting affines to formats."""
-    if (image_orientation, sw_tool) == ("oblique", "afni"):
-        pytest.skip("AFNI Deoblique unsupported.")
-
     tmpdir.chdir()
     img = get_testdata[image_orientation]
     # Generate test transform

nitransforms/tests/test_linear.py

@@ -112,8 +112,6 @@ def test_displacements_field1(
     axis,
 ):
     """Check a translation-only field on one or more axes, different image orientations."""
-    if (image_orientation, sw_tool) == ("oblique", "afni") and axis in ((1, 2), (0, 1, 2)):
-        pytest.skip("AFNI Deoblique unsupported.")
     os.chdir(str(tmp_path))
     nii = get_testdata[image_orientation]
     msk = get_testmask[image_orientation]
@@ -176,6 +174,7 @@ def test_displacements_field1(
     sw_moved = nb.load("resampled.nii.gz")
 
     nt_moved = xfm.apply(nii, order=0)
+    nt_moved.set_data_dtype(nii.get_data_dtype())
     nt_moved.to_filename("nt_resampled.nii.gz")
     sw_moved.set_data_dtype(nt_moved.get_data_dtype())
     diff = (