NF: Conformation function and CLI tool to apply shape, orientation and zooms

nibabel/processing.py

@@ -22,8 +22,9 @@
 spnd, _, _ = optional_package('scipy.ndimage')
 
 from .affines import AffineError, to_matvec, from_matvec, append_diag
+from .funcs import as_closest_canonical
 from .spaces import vox2out_vox
-from .nifti1 import Nifti1Image
+from .nifti1 import Nifti1Header, Nifti1Image
 from .imageclasses import spatial_axes_first
 
 SIGMA2FWHM = np.sqrt(8 * np.log(2))
@@ -310,3 +311,96 @@ def smooth_image(img,
                                    mode=mode,
                                    cval=cval)
     return out_class(sm_data, img.affine, img.header)
+
+
+def _transform_range(x, new_min, new_max):
+    """ Transform data to a new range, while maintaining ratios.
+
+    Parameters
+    ----------
+    x : array-like
+        The data to transform.
+    new_min, new_max : scalar
+        The minimum and maximum of the output array.
+
+    Returns
+    -------
+    transformed : array-like
+        A copy of the transformed data.
+
+    Examples
+    --------
+    >>> _transform_range([2, 4, 6], -1, 1)
+    array([-1.,  0.,  1.])
+    """
+    x = np.asarray(x)
+    x_min, x_max = x.min(), x.max()
+    return (x - x_min) * (new_max - new_min) / (x_max - x_min) + new_min
+
+
+def conform(from_img,
+            out_shape=(256, 256, 256),
+            voxel_size=(1.0, 1.0, 1.0),
+            order=3,
+            cval=0.0,
+            out_class=Nifti1Image):
+    """ Resample image to ``out_shape`` with voxels of size ``voxel_size``.
+
+    Using the default arguments, this function is meant to replicate most parts
+    of FreeSurfer's ``mri_convert --conform`` command. Specifically, this
+    function:
+        - Resamples data to ``output_shape``
+        - Resamples voxel sizes to ``voxel_size``
+        - Transforms data to range [0, 255] (while maintaining ratios)
+        - Casts to unsigned eight-bit integer
+        - Reorients to RAS (``mri_convert --conform`` reorients to LIA)
+
+    Parameters
+    ----------
+    from_img : object
+        Object having attributes ``dataobj``, ``affine``, ``header`` and
+        ``shape``. If `out_class` is not None, ``img.__class__`` should be able
+        to construct an image from data, affine and header.
+    out_shape : sequence, optional
+        The shape of the output volume. Default is (256, 256, 256).
+    voxel_size : sequence, optional
+        The size in millimeters of the voxels in the resampled output. Default
+        is 1mm isotropic.
+    order : int, optional
+        The order of the spline interpolation, default is 3.  The order has to
+        be in the range 0-5 (see ``scipy.ndimage.affine_transform``)
+    cval : scalar, optional
+        Value used for points outside the boundaries of the input if
+        ``mode='constant'``. Default is 0.0 (see
+        ``scipy.ndimage.affine_transform``)
+    out_class : None or SpatialImage class, optional
+        Class of output image.  If None, use ``from_img.__class__``.
+
+    Returns
+    -------
+    out_img : object
+        Image of instance specified by `out_class`, containing data output from
+        resampling `from_img` into axes aligned to the output space of
+        ``from_img.affine``
+    """
+    if from_img.ndim != 3:
+        raise ValueError("Only 3D images are supported.")
+    # Create fake image of the image we want to resample to.
+    hdr = Nifti1Header()
+    hdr.set_data_shape(out_shape)
+    hdr.set_zooms(voxel_size)
+    dst_aff = hdr.get_best_affine()
+    to_img = Nifti1Image(np.empty(out_shape), affine=dst_aff, header=hdr)
+    # Resample input image.
+    out_img = resample_from_to(
+        from_img=from_img, to_vox_map=to_img, order=order, mode="constant",
+        cval=cval, out_class=out_class)
+    # Cast to uint8.
+    data = out_img.get_fdata()
+    data = _transform_range(data, new_min=0.0, new_max=255.0)
+    data = data.round(out=data).astype(np.uint8)
+    out_img._dataobj = data
+    out_img.set_data_dtype(data.dtype)
+    # Reorient to RAS.
+    out_img = as_closest_canonical(out_img)
+    return out_img

nibabel/tests/test_processing.py

@@ -19,10 +19,11 @@
 
 import nibabel as nib
 from nibabel.processing import (sigma2fwhm, fwhm2sigma, adapt_affine,
-                                resample_from_to, resample_to_output, smooth_image)
+                                resample_from_to, resample_to_output, smooth_image,
+                                _transform_range, conform)
 from nibabel.nifti1 import Nifti1Image
 from nibabel.nifti2 import Nifti2Image
-from nibabel.orientations import flip_axis, inv_ornt_aff
+from nibabel.orientations import aff2axcodes, flip_axis, inv_ornt_aff
 from nibabel.affines import (AffineError, from_matvec, to_matvec, apply_affine,
                              voxel_sizes)
 from nibabel.eulerangles import euler2mat
@@ -426,3 +427,33 @@ def test_against_spm_resample():
     moved2output = resample_to_output(moved_anat, 4, order=1, cval=np.nan)
     spm2output = nib.load(pjoin(DATA_DIR, 'reoriented_anat_moved.nii'))
     assert_spm_resampling_close(moved_anat, moved2output, spm2output);
+
+
+def test__transform_range():
+    assert_array_equal(_transform_range([2, 4, 6], -1, 1), [-1, 0, 1])
+    assert_array_equal(_transform_range([-1, 0, 1], 2, 6), [2, 4, 6])
+    assert_array_equal(_transform_range(np.arange(11), 0, 5),
+                       np.arange(0, 5.5, 0.5))
+    assert_array_equal(_transform_range(np.arange(-100, 101), 0, 200),
+                       np.arange(201))
+
+
+@needs_scipy
+def test_conform():
+    anat = nib.load(pjoin(DATA_DIR, 'anatomical.nii'))
+
+    c = conform(anat)
+    assert c.shape == (256, 256, 256)
+    assert c.header.get_zooms() == (1, 1, 1)
+    assert c.dataobj.dtype == np.dtype(np.uint8)
+    assert aff2axcodes(c.affine) == ('R', 'A', 'S')
+
+    c = conform(anat, out_shape=(100, 100, 200), voxel_size=(2, 2, 1.5))
+    assert c.shape == (100, 100, 200)
+    assert c.header.get_zooms() == (2, 2, 1.5)
+    assert c.dataobj.dtype == np.dtype(np.uint8)
+    assert aff2axcodes(c.affine) == ('R', 'A', 'S')
+
+    # Error on non-3D images.
+    func = nib.load(pjoin(DATA_DIR, 'functional.nii'))
+    assert_raises(ValueError, conform, func)