diff --git a/nipype/interfaces/image.py b/nipype/interfaces/image.py new file mode 100644 index 0000000000..f5279c9bc6 --- /dev/null +++ b/nipype/interfaces/image.py @@ -0,0 +1,145 @@ +# -*- coding: utf-8 -*- +# emacs: -*- mode: python; py-indent-offset: 4; indent-tabs-mode: nil -*- +# vi: set ft=python sts=4 ts=4 sw=4 et: + +import numpy as np +import nibabel as nb + +from ..utils.filemanip import fname_presuffix +from .base import (SimpleInterface, TraitedSpec, BaseInterfaceInputSpec, + traits, File) + +_axes = ('RL', 'AP', 'SI') +_orientations = tuple( + ''.join((x[i], y[j], z[k])) + for x in _axes for y in _axes for z in _axes + if x != y != z != x + for i in (0, 1) for j in (0, 1) for k in (0, 1)) + + +class ReorientInputSpec(BaseInterfaceInputSpec): + in_file = File(exists=True, mandatory=True, desc='Input image') + orientation = traits.Enum(_orientations, usedefault=True, + desc='Target axis orientation') + + +class ReorientOutputSpec(TraitedSpec): + out_file = File(exists=True, desc='Reoriented image') + transform = File(exists=True, + desc='Affine transform from input orientation to output') + + +class Reorient(SimpleInterface): + """Conform an image to a given orientation + + Flips and reorder the image data array so that the axes match the + directions indicated in ``orientation``. + The default ``RAS`` orientation corresponds to the first axis being ordered + from left to right, the second axis from posterior to anterior, and the + third axis from inferior to superior. + + For oblique images, the original orientation is considered to be the + closest plumb orientation. + + No resampling is performed, and thus the output image is not de-obliqued + or registered to any other image or template. + + The effective transform is calculated from the original affine matrix to + the reoriented affine matrix. + + Examples + -------- + + If an image is not reoriented, the original file is not modified + + >>> import numpy as np + >>> from nipype.interfaces.image import Reorient + >>> reorient = Reorient(orientation='LPS') + >>> reorient.inputs.in_file = 'segmentation0.nii.gz' + >>> res = reorient.run() + >>> res.outputs.out_file + 'segmentation0.nii.gz' + + >>> print(np.loadtxt(res.outputs.transform)) + [[1. 0. 0. 0.] + [0. 1. 0. 0.] + [0. 0. 1. 0.] + [0. 0. 0. 1.]] + + >>> reorient.inputs.orientation = 'RAS' + >>> res = reorient.run() + >>> res.outputs.out_file # doctest: +ELLIPSIS + '.../segmentation0_ras.nii.gz' + + >>> print(np.loadtxt(res.outputs.transform)) + [[-1. 0. 0. 60.] + [ 0. -1. 0. 72.] + [ 0. 0. 1. 0.] + [ 0. 0. 0. 1.]] + + """ + input_spec = ReorientInputSpec + output_spec = ReorientOutputSpec + + def _run_interface(self, runtime): + from nibabel.orientations import ( + axcodes2ornt, ornt_transform, inv_ornt_aff) + + fname = self.inputs.in_file + orig_img = nb.load(fname) + + # Find transform from current (approximate) orientation to + # target, in nibabel orientation matrix and affine forms + orig_ornt = nb.io_orientation(orig_img.affine) + targ_ornt = axcodes2ornt(self.inputs.orientation) + transform = ornt_transform(orig_ornt, targ_ornt) + affine_xfm = inv_ornt_aff(transform, orig_img.shape) + + # Check can be eliminated when minimum nibabel version >= 2.2 + if hasattr(orig_img, 'as_reoriented'): + reoriented = orig_img.as_reoriented(transform) + else: + reoriented = _as_reoriented_backport(orig_img, transform) + + # Image may be reoriented + if reoriented is not orig_img: + suffix = '_' + self.inputs.orientation.lower() + out_name = fname_presuffix(fname, suffix=suffix, + newpath=runtime.cwd) + reoriented.to_filename(out_name) + else: + out_name = fname + + mat_name = fname_presuffix(fname, suffix='.mat', + newpath=runtime.cwd, use_ext=False) + np.savetxt(mat_name, affine_xfm, fmt='%.08f') + + self._results['out_file'] = out_name + self._results['transform'] = mat_name + + return runtime + + +def _as_reoriented_backport(img, ornt): + """Backport of img.as_reoriented as of nibabel 2.2.0""" + from nibabel.orientations import inv_ornt_aff + if np.array_equal(ornt, [[0, 1], [1, 1], [2, 1]]): + return img + + t_arr = nb.apply_orientation(img.get_data(), ornt) + new_aff = img.affine.dot(inv_ornt_aff(ornt, img.shape)) + reoriented = img.__class__(t_arr, new_aff, img.header) + + if isinstance(reoriented, nb.Nifti1Pair): + # Also apply the transform to the dim_info fields + new_dim = list(reoriented.header.get_dim_info()) + for idx, value in enumerate(new_dim): + # For each value, leave as None if it was that way, + # otherwise check where we have mapped it to + if value is None: + continue + new_dim[idx] = np.where(ornt[:, 0] == idx)[0] + + reoriented.header.set_dim_info(*new_dim) + + return reoriented diff --git a/nipype/interfaces/tests/test_auto_Reorient.py b/nipype/interfaces/tests/test_auto_Reorient.py new file mode 100644 index 0000000000..19d322df4c --- /dev/null +++ b/nipype/interfaces/tests/test_auto_Reorient.py @@ -0,0 +1,30 @@ +# AUTO-GENERATED by tools/checkspecs.py - DO NOT EDIT +from __future__ import unicode_literals +from ..image import Reorient + + +def test_Reorient_inputs(): + input_map = dict( + ignore_exception=dict( + deprecated='1.0.0', + nohash=True, + usedefault=True, + ), + in_file=dict(mandatory=True, ), + orientation=dict(usedefault=True, ), + ) + inputs = Reorient.input_spec() + + for key, metadata in list(input_map.items()): + for metakey, value in list(metadata.items()): + assert getattr(inputs.traits()[key], metakey) == value +def test_Reorient_outputs(): + output_map = dict( + out_file=dict(), + transform=dict(), + ) + outputs = Reorient.output_spec() + + for key, metadata in list(output_map.items()): + for metakey, value in list(metadata.items()): + assert getattr(outputs.traits()[key], metakey) == value diff --git a/nipype/interfaces/tests/test_image.py b/nipype/interfaces/tests/test_image.py new file mode 100644 index 0000000000..bb4adf1d01 --- /dev/null +++ b/nipype/interfaces/tests/test_image.py @@ -0,0 +1,64 @@ +# emacs: -*- mode: python; py-indent-offset: 4; indent-tabs-mode: nil -*- +# vi: set ft=python sts=4 ts=4 sw=4 et: +import numpy as np +import nibabel as nb +import pytest + +from nibabel.orientations import axcodes2ornt, ornt_transform + +from ..image import _as_reoriented_backport, _orientations +from ... import LooseVersion + +nibabel22 = LooseVersion(nb.__version__) >= LooseVersion('2.2.0') + + +@pytest.mark.skipif(not nibabel22, + reason="Old nibabel - can't directly compare") +def test_reorientation_backport(): + pixdims = ((1, 1, 1), (2, 2, 3)) + data = np.random.normal(size=(17, 18, 19, 2)) + + for pixdim in pixdims: + # Generate a randomly rotated affine + angles = np.random.uniform(-np.pi, np.pi, 3) * [1, 0.5, 1] + rot = nb.eulerangles.euler2mat(*angles) + scale = np.diag(pixdim) + translation = np.array((17, 18, 19)) / 2 + affine = nb.affines.from_matvec(rot.dot(scale), translation) + + # Create image + img = nb.Nifti1Image(data, affine) + dim_info = {'freq': 0, 'phase': 1, 'slice': 2} + img.header.set_dim_info(**dim_info) + + # Find a random, non-identity transform + targ_ornt = orig_ornt = nb.io_orientation(affine) + while np.array_equal(targ_ornt, orig_ornt): + new_code = np.random.choice(_orientations) + targ_ornt = axcodes2ornt(new_code) + + identity = ornt_transform(orig_ornt, orig_ornt) + transform = ornt_transform(orig_ornt, targ_ornt) + + # Identity transform returns exact image + assert img.as_reoriented(identity) is img + assert _as_reoriented_backport(img, identity) is img + + reoriented_a = img.as_reoriented(transform) + reoriented_b = _as_reoriented_backport(img, transform) + + flips_only = img.shape == reoriented_a.shape + + # Reorientation changes affine and data array + assert not np.allclose(img.affine, reoriented_a.affine) + assert not (flips_only and + np.allclose(img.get_data(), reoriented_a.get_data())) + # Dimension info changes iff axes are reordered + assert flips_only == np.array_equal(img.header.get_dim_info(), + reoriented_a.header.get_dim_info()) + + # Both approaches produce equivalent images + assert np.allclose(reoriented_a.affine, reoriented_b.affine) + assert np.array_equal(reoriented_a.get_data(), reoriented_b.get_data()) + assert np.array_equal(reoriented_a.header.get_dim_info(), + reoriented_b.header.get_dim_info())