Merge pull request #186 from nipy/fix/optimize-afni-toras

FIX: Inefficient iterative reloading of reference and moving images

Author: oesteban

nitransforms/io/afni.py

@@ -130,9 +130,17 @@ class AFNILinearTransformArray(BaseLinearTransformList):
 
     def to_ras(self, moving=None, reference=None):
         """Return a nitransforms' internal RAS matrix."""
-        return np.stack(
-            [xfm.to_ras(moving=moving, reference=reference) for xfm in self.xforms]
-        )
+
+        pre_rotation = post_rotation = np.eye(4)
+        if reference is not None and _is_oblique(ref_aff := _ensure_image(reference).affine):
+            pre_rotation = _cardinal_rotation(ref_aff, True)
+        if moving is not None and _is_oblique(mov_aff := _ensure_image(moving).affine):
+            post_rotation = _cardinal_rotation(mov_aff, False)
+
+        return np.stack([
+            post_rotation @ (xfm.to_ras() @ pre_rotation)
+            for xfm in self.xforms
+        ])
 
     def to_string(self):
         """Convert to a string directly writeable to file."""
@@ -144,14 +152,22 @@ def to_string(self):
                 if line.strip()
             ]
             strings += lines
-        return "\n".join(strings)
+        return "\n".join(strings + [""])
 
     @classmethod
     def from_ras(cls, ras, moving=None, reference=None):
         """Create an ITK affine from a nitransform's RAS+ matrix."""
         _self = cls()
+
+        pre_rotation = post_rotation = np.eye(4)
+
+        if reference is not None and _is_oblique(ref_aff := _ensure_image(reference).affine):
+            pre_rotation = _cardinal_rotation(ref_aff, False)
+        if moving is not None and _is_oblique(mov_aff := _ensure_image(moving).affine):
+            post_rotation = _cardinal_rotation(mov_aff, True)
+
         _self.xforms = [
-            cls._inner_type.from_ras(ras[i, ...], moving=moving, reference=reference)
+            cls._inner_type.from_ras(post_rotation @ ras[i, ...] @ pre_rotation)
             for i in range(ras.shape[0])
         ]
         return _self

nitransforms/tests/data/affine-LAS.afni-array

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

nitransforms/tests/data/affine-LPS.afni-array

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

nitransforms/tests/data/affine-RAS.afni-array

@@ -0,0 +1,3 @@
+# 3dvolreg matrices (DICOM-to-DICOM, row-by-row):
+0.999999	-0.000999999	-0.001	-4	0.00140494	0.621609	0.783327	-2	-0.000161717	-0.783327	0.62161	-1
+0.999999	-0.000999999	-0.001	-4	0.00140494	0.621609	0.783327	-2	-0.000161717	-0.783327	0.62161	-1

nitransforms/tests/data/affine-oblique.afni-array

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

nitransforms/tests/test_io.py

@@ -180,7 +180,7 @@ def test_LT_conversions(data_path, fname):
         "oblique",
     ],
 )
-@pytest.mark.parametrize("sw", ["afni", "fsl", "fs", "itk"])
+@pytest.mark.parametrize("sw", ["afni", "fsl", "fs", "itk", "afni-array"])
 def test_Linear_common(tmpdir, data_path, sw, image_orientation, get_testdata):
     tmpdir.chdir()
 
@@ -190,6 +190,8 @@ def test_Linear_common(tmpdir, data_path, sw, image_orientation, get_testdata):
     ext = ""
     if sw == "afni":
         factory = afni.AFNILinearTransform
+    elif sw == "afni-array":
+        factory = afni.AFNILinearTransformArray
     elif sw == "fsl":
         factory = fsl.FSLLinearTransform
     elif sw == "itk":
@@ -222,6 +224,9 @@ def test_Linear_common(tmpdir, data_path, sw, image_orientation, get_testdata):
 
     # Test from_ras
     RAS = from_matvec(euler2mat(x=0.9, y=0.001, z=0.001), [4.0, 2.0, -1.0])
+    if sw == "afni-array":
+        RAS = np.array([RAS, RAS])
+
     xfm = factory.from_ras(RAS, reference=reference, moving=moving)
     assert np.allclose(xfm.to_ras(reference=reference, moving=moving), RAS)