ENH: New image.roi module for extracting image ROI/expanding FOV

fsl/utils/image/roi.py

+#!/usr/bin/env python
+#
+# roi.py - Extract an ROI of an image.
+#
+# Author: Paul McCarthy <pauldmccarthy@gmail.com>
+#
+"""This module provides the :func:`roi` function, which can be used to extract
+a region-of-interest from, or expand the field-of-view of, an :class:`.Image`.
+"""
+
+
+import numpy as np
+
+import fsl.data.image      as fslimage
+import fsl.utils.transform as transform
+
+
+def _normaliseBounds(shape, bounds):
+    """Adjust the given ``bounds`` so that it is compatible with the given
+    ``shape``.
+
+    Bounds must be specified for at least three dimensions - if the shape
+    has more than three dimensions, additional bounds are added.
+
+    A :exc:`ValueError` is raised if the provided bounds are invalid.
+
+    :arg bounds: Sequence of ``(lo, hi)`` bounds - see :func:`roi`.
+    :returns:    An adjusted sequence of bounds.
+    """
+
+    bounds = list(bounds)
+
+    if len(bounds) < 3:
+        raise ValueError('')
+
+    if len(bounds) > len(shape):
+        raise ValueError('')
+
+    for b in bounds:
+        if len(b) != 2 or b[0] >= b[1]:
+            raise ValueError('')
+
+    if len(bounds) < len(shape):
+        for s in shape[len(bounds):]:
+            bounds.append((0, s))
+
+    return bounds
+
+
+def roi(image, bounds):
+    """Extract an ROI from the given ``image`` according to the given
+    ``bounds``.
+
+    This function can also be used to pad, or expand the field-of-view, of an
+    image, by passing in negative low bound values, or high bound values which
+    are larger than the image shape. The padded region will contain zeroes.
+
+    :arg image:  :class:`.Image` object
+
+    :arg bounds: Must be a sequence of tuples, containing the low/high bounds
+                 for each voxel dimension, where the low bound is *inclusive*,
+                 and the high bound is *exclusive*. For 4D images, the bounds
+                 for the fourth dimension are optional.
+
+    :returns:    A new :class:`.Image` object containing the region specified
+                 by the ``bounds``.
+    """
+
+    bounds = _normaliseBounds(image.shape, bounds)
+
+    newshape = [hi - lo for lo, hi in bounds]
+    oldslc   = []
+    newslc   = []
+
+    # Figure out how to slice the input image
+    # data array, and the corresponding slice
+    # in the output data array.
+    for (lo, hi), oldlen, newlen in zip(bounds, image.shape, newshape):
+
+        oldlo = max(lo, 0)
+        oldhi = min(hi, oldlen)
+        newlo = max(0, -lo)
+        newhi = newlo + (oldhi - oldlo)
+
+        oldslc.append(slice(oldlo, oldhi))
+        newslc.append(slice(newlo, newhi))
+
+    oldslc = tuple(oldslc)
+    newslc = tuple(newslc)
+
+    # Copy the ROI into the new data array
+    newdata         = np.zeros(newshape, dtype=image.dtype)
+    newdata[newslc] = image.data[oldslc]
+
+    # Create a new affine for the ROI,
+    # with an appropriate offset along
+    # each spatial dimension
+    oldaff = image.voxToWorldMat
+    offset = [lo for lo, hi in bounds[:3]]
+    offset = transform.scaleOffsetXform([1, 1, 1], offset)
+    newaff = transform.concat(oldaff, offset)
+
+    return fslimage.Image(newdata,
+                          xform=newaff,
+                          header=image.header,
+                          name=image.name + '_roi')