Skip to content
GitLab
Commit 2f963d24 authored by Paul McCarthy's avatar Paul McCarthy 🚵
Browse files

Merge branch 'enh/resample' into 'master' 

Enh/resample

See merge request fsl/fslpy!122
parents 05b35efa 73388cb5
Pipeline #3691 canceled with stages
in 9 minutes and 26 seconds
Hide whitespace changes
Inline Side-by-side
.ci/test_template.sh
View file @ 2f963d24
......@@ -48,11 +48,11 @@ xvfb-run python setup.py test --addopts="$TEST_OPTS tests/test_platform.py"
# unintuitively, includes nobody)
chmod -R a+w `pwd`
cmd="source /test.venv/bin/activate && python setup.py test"
cmd="$cmd --addopts='$TEST_OPTS tests/test_immv_imcp.py'"
cmd="$cmd --addopts='$TEST_OPTS tests/test_scripts/test_immv_imcp.py tests/test_immv_imcp.py'"
su -s /bin/bash -c "$cmd" nobody
# All other tests can be run as normal.
python setup.py test --addopts="$TEST_OPTS -m 'not longtest' --ignore=tests/test_idle.py --ignore=tests/test_platform.py --ignore=tests/test_immv_imcp.py"
python setup.py test --addopts="$TEST_OPTS -m 'not longtest' --ignore=tests/test_idle.py --ignore=tests/test_platform.py --ignore=tests/test_immv_imcp.py --ignore=tests/test_scripts/test_immv_imcp.py"
# Long tests are only run on release branches
if [[ $CI_COMMIT_REF_NAME == v* ]]; then
......
CHANGELOG.rst
View file @ 2f963d24
......@@ -2,6 +2,49 @@ This document contains the ``fslpy`` release history in reverse chronological
order.
2.2.0 (Wednesday May 8th 2019)
------------------------------
Added
^^^^^
* New :mod:`.resample_image` script.
* New :mod:`.resample` module (replacing the :func:`.Image.resample` method),
containing functions to resample an :class:`.Image`.
* New :func:`.resample.resampleToPixdim` and
:func:`.resample.resampleToReference` functions, convenience wrappers around
:func:`.resample.resample`.
* New :func:`.idle.block` function.
Changed
^^^^^^^
* The :func:`.resample` function (formerly :meth:`.Image.resample`) now
accepts ``origin`` and ``matrix`` parameters, which can be used to adjust
the alignment of the voxel grids of the input and output images.
* The :func:`.transform.decompose` function now accepts both ``(3, 3)``
and ``(4, 4)`` matrices.
Fixed
^^^^^
* Minor fixes to some :mod:`.filetree` tree definitions.
Deprecated
^^^^^^^^^^
* The :meth:`.Image.resample` method has been deprecated in favour of the
:func:`.resample.resample` function.
2.1.0 (Saturday April 13th 2019)
--------------------------------
......
fsl/data/atlases.py
View file @ 2f963d24
......@@ -53,6 +53,7 @@ import numpy as np
import fsl.data.image as fslimage
import fsl.data.constants as constants
from fsl.utils.platform import platform as platform
import fsl.utils.image.resample as resample
import fsl.utils.transform as transform
import fsl.utils.notifier as notifier
import fsl.utils.settings as fslsettings
......@@ -695,9 +696,8 @@ class Atlas(fslimage.Image):
# for resampling, as it is most likely
# that the mask is binary.
try:
mask, xform = mask.resample(self.shape[:3],
dtype=np.float32,
order=0)
mask, xform = resample.resample(
mask, self.shape[:3], dtype=np.float32, order=0)
except ValueError:
raise MaskError('Mask has wrong number of dimensions')
......
fsl/data/image.py
View file @ 2f963d24
......@@ -41,7 +41,6 @@ import warnings
import six
import numpy as np
import scipy.ndimage as ndimage
import nibabel as nib
import nibabel.fileslice as fileslice
......@@ -1166,96 +1165,12 @@ class Image(Nifti):
self.notify(topic='saveState')
def resample(self,
newShape,
sliceobj=None,
dtype=None,
order=1,
smooth=True):
"""Returns a copy of the data in this ``Image``, resampled to the
specified ``newShape``.
:arg newShape: Desired shape. May containg floating point values,
in which case the resampled image will have shape
``round(newShape)``, but the voxel sizes will
have scales ``self.shape / newShape``.
:arg sliceobj: Slice into this ``Image``. If ``None``, the whole
image is resampled, and it is assumed that it has the
same number of dimensions as ``newShape``. A
:exc:`ValueError` is raised if this is not the case.
:arg dtype: ``numpy`` data type of the resampled data. If ``None``,
the :meth:`dtype` of this ``Image`` is used.
:arg order: Spline interpolation order, passed through to the
``scipy.ndimage.affine_transform`` function - ``0``
corresponds to nearest neighbour interpolation, ``1``
(the default) to linear interpolation, and ``3`` to
cubic interpolation.
:arg smooth: If ``True`` (the default), the data is smoothed before
being resampled, but only along axes which are being
down-sampled (i.e. where
``newShape[i] < self.shape[i]``).
:returns: A tuple containing:
- A ``numpy`` array of shape ``newShape``, containing
an interpolated copy of the data in this ``Image``.
- A ``numpy`` array of shape ``(4, 4)``, containing the
adjusted voxel-to-world transformation for the spatial
dimensions of the resampled data.
"""
if sliceobj is None: sliceobj = slice(None)
if dtype is None: dtype = self.dtype
data = self[sliceobj]
data = np.array(data, dtype=dtype, copy=False)
oldShape = np.array(data.shape, dtype=np.float)
newShape = np.array(newShape, dtype=np.float)
if len(oldShape) != len(newShape):
raise ValueError('Shapes don\'t match')
if not np.all(np.isclose(oldShape, newShape)):
ratio = oldShape / newShape
newShape = np.array(np.round(newShape), dtype=np.int)
scale = np.diag(ratio)
# If interpolating and smoothing, we apply a
# gaussian filter along axes with a resampling
# ratio greater than 1.1. We do this so that
# interpolation has an effect when down-sampling
# to a resolution where the voxel centres are
# aligned (as otherwise any interpolation regime
# will be equivalent to nearest neighbour). This
# more-or-less mimics the behaviour of FLIRT.
if order > 0 and smooth:
sigma = np.array(ratio)
sigma[ratio < 1.1] = 0
sigma[ratio >= 1.1] *= 0.425
data = ndimage.gaussian_filter(data, sigma)
data = ndimage.affine_transform(data,
scale,
output_shape=newShape,
order=order)
# Construct an affine transform which
# puts the resampled image into the
# same world coordinate system as this
# image.
scale = transform.scaleOffsetXform(ratio[:3], 0)
xform = transform.concat(self.voxToWorldMat, scale)
else:
xform = self.voxToWorldMat
return data, xform
@deprecated.deprecated('2.2.0', '3.0.0',
'Use fsl.utils.image.resample instead.')
def resample(self, *args, **kwargs):
"""Deprecated - use :func:`.image.resample` instead. """
from fsl.utils.image.resample import resample
return resample(self, *args, **kwargs)
def __getitem__(self, sliceobj):
......
fsl/scripts/resample_image.py 0 → 100644
View file @ 2f963d24
#!/usr/bin/env python
#
# resample_image.py - Script to resample an image
#
# Author: Paul McCarthy <pauldmccarthy@gmail.com>
#
"""This module defines the ``resample_image`` script, for resampling
a NIfTI image.
"""
import textwrap as tw
import sys
import argparse
import numpy as np
import fsl.utils.parse_data as parse_data
import fsl.utils.image.resample as resample
import fsl.data.image as fslimage
ARGS = {
'input' : ('input',),
'output' : ('output',),
'shape' : ('-s', '--shape'),
'dim' : ('-d', '--dim'),
'reference' : ('-r', '--reference'),
'interp' : ('-i', '--interp'),
'origin' : ('-o', '--origin'),
'dtype' : ('-dt', '--dtype'),
'smooth' : ('-n', '--nosmooth')}
OPTS = {
'input' : dict(type=parse_data.Image),
'output' : dict(type=parse_data.ImageOut),
'reference' : dict(type=parse_data.Image, metavar='IMAGE'),
'shape' : dict(type=int, nargs=3, metavar=('X', 'Y', 'Z')),
'dim' : dict(type=float, nargs=3, metavar=('X', 'Y', 'Z')),
'interp' : dict(choices=('nearest', 'linear', 'cubic'),
default='linear'),
'origin' : dict(choices=('centre', 'corner'), default='centre'),
'dtype' : dict(choices=('char', 'short', 'int', 'float', 'double')),
'smooth' : dict(dest='smooth', action='store_false')}
HELPS = {
'input' : 'Input image',
'output' : 'Output image',
'shape' : 'Output shape',
'dim' : 'Output voxel dimensions',
'reference' : 'Resample input to the space of this reference image'
'(overrides --origin)',
'interp' : 'Interpolation (default: linear)',
'origin' : 'Resampling origin (default: centre)',
'dtype' : 'Data type (default: data type of input image)',
'smooth' : 'Do not smooth image when downsampling'}
DESC = tw.dedent("""
Resample an image to different dimensions.
""").strip()
DEST_DESC = tw.dedent("""
Specify the resampling destination space using one of the following
options. Note that the --reference option will cause the field-of-view
of the input image to be changed to that of the reference image.
""").strip()
USAGE = 'resample_image (--shape|--dim|--reference) [options] input output'
INTERPS = {'nearest' : 0,
'linear' : 1,
'cubic' : 3}
DTYPES = {'char' : np.uint8,
'short' : np.int16,
'int' : np.int32,
'float' : np.float32,
'double' : np.float64}
def parseArgs(argv):
"""Parses command-line arguments.
:arg argv: Sequence of command-line arguments
:returns: An ``argparse.Namespace`` object containing parsed arguments.
"""
parser = argparse.ArgumentParser(prog='resample_image',
usage=USAGE,
description=DESC)
dest = parser.add_argument_group('Resampling destination', DEST_DESC)
dest = dest.add_mutually_exclusive_group(required=True)
for a in ('input', 'output', 'interp', 'origin',
'dtype', 'smooth'):
parser.add_argument(*ARGS[a], help=HELPS[a], **OPTS[a])
for a in ('shape', 'dim', 'reference'):
dest.add_argument(*ARGS[a], help=HELPS[a], **OPTS[a])
if len(argv) == 0:
parser.print_help()
sys.exit(0)
args = parser.parse_args(argv)
args.interp = INTERPS[ args.interp]
args.dtype = DTYPES.get(args.dtype, args.input.dtype)
return args
def main(argv=None):
"""Entry point for ``resample_image``. Parses arguments, resamples the
input image, and saves it to the specified output file.
:arg argv: Sequence of command-line arguments. If not provided, taken
from ``sys.argv``.
"""
if argv is None:
argv = sys.argv[1:]
args = parseArgs(argv)
reskwargs = {
'dtype' : args.dtype,
'order' : args.interp,
'smooth' : args.smooth,
'origin' : args.origin}
# One of these is guaranteed to be set
if args.shape is not None:
func = resample.resample
resargs = (args.input, args.shape)
elif args.dim is not None:
func = resample.resampleToPixdims
resargs = (args.input, args.dim)
elif args.reference is not None:
func = resample.resampleToReference
resargs = (args.input, args.reference)
resampled, xform = func(*resargs, **reskwargs)
if args.reference is None:
hdr = args.input.header
else:
hdr = args.reference.header
xform = None
resampled = fslimage.Image(resampled, xform=xform, header=hdr)
resampled.save(args.output)
return 0
if __name__ == '__main__':
sys.exit(main())
fsl/utils/idle.py
View file @ 2f963d24
......@@ -19,6 +19,7 @@ Idle tasks
.. autosummary::
:nosignatures:
block
idle
idleWhen
inIdle
......@@ -377,6 +378,29 @@ def cancelIdle(taskName):
_idleQueueDict[taskName].timeout = -1
def block(secs, delta=0.01):
"""Blocks for the specified number of seconds, yielding to the main ``wx``
loop.
If ``wx`` is not available, or a ``wx`` application is not running, this
function is equivalent to ``time.sleep(secs)``.
:arg secs: Time in seconds to block
:arg delta: Time in seconds to sleep between successive yields to ``wx``.
"""
from fsl.utils.platform import platform as fslplatform
if not fslplatform.haveGui:
time.sleep(secs)
else:
import wx
start = time.time()
while (time.time() - start) < secs:
wx.YieldIfNeeded()
time.sleep(delta)
def idle(task, *args, **kwargs):
"""Run the given task on a ``wx.EVT_IDLE`` event.
......@@ -556,7 +580,7 @@ def idleWhen(func, condition, *args, **kwargs):
def wait(threads, task, *args, **kwargs):
"""Creates and starts a new ``Thread`` which waits for all of the ``Thread``
instances to finsih (by ``join``ing them), and then runs the given
instances to finish (by ``join``ing them), and then runs the given
``task`` via :func:`idle`.
If the ``direct`` parameter is ``True``, or a ``wx.App`` is not running,
......
fsl/utils/image/__init__.py 0 → 100644
View file @ 2f963d24
#!/usr/bin/env python
#
# __init__.py - The fsl.utils.image package
#
# Author: Paul McCarthy <pauldmccarthy@gmail.com>
#
"""The :mod:`fsl.utils.image` oackage contains algorithms and utilities for
manipulating and working with :class:`.Image` objects.
The following modules are available:
.. autosumary::
:nosignature
.image.resample
"""
fsl/utils/image/resample.py 0 → 100644
View file @ 2f963d24
#!/usr/bin/env python
#
# resample.py - The resample functino
#
# Author: Paul McCarthy <pauldmccarthy@gmail.com>
#
"""This module defines the :func:`resample` function, which can be used
to resample an :class:`.Image` object to a different resolution.
The :func:`resampleToPixdims` and :func:`resampleToReference` functions
are convenience wrappers around :func:`resample`.
The :func:`applySmoothing` and :func:`calculateMatrix` functions are
sub-functions of :func:`resample`.
"""
import collections.abc as abc
import numpy as np
import scipy.ndimage as ndimage
import fsl.utils.transform as transform
def resampleToPixdims(image, newPixdims, **kwargs):
"""Resample ``image`` so that it has the specified voxel dimensions.
This is a wrapper around :func:`resample` - refer to its documenttion
for details on the other arguments and the return values.
:arg image: :class:`.Image` to resample
:arg pixdims: New voxel dimensions to resample ``image`` to.
"""
newPixdims = np.array(newPixdims)
oldShape = np.array(image.shape)
oldPixdims = np.array(image.pixdim)
newShape = oldShape * (oldPixdims / newPixdims)
return resample(image, newShape, **kwargs)
def resampleToReference(image, reference, **kwargs):
"""Resample ``image`` into the space of the ``reference``.
This is a wrapper around :func:`resample` - refer to its documenttion
for details on the other arguments and the return values.
:arg image: :class:`.Image` to resample
:arg reference: :class:`.Nifti` defining the space to resample ``image``
into
"""
kwargs['mode'] = kwargs.get('mode', 'constant')
kwargs['newShape'] = reference.shape
kwargs['matrix'] = transform.concat(image.worldToVoxMat,
reference.voxToWorldMat)
return resample(image, **kwargs)
def resample(image,
newShape,
sliceobj=None,
dtype=None,
order=1,
smooth=True,
origin='centre',
matrix=None,
mode='nearest',
cval=0):
"""Returns a copy of the data in the ``image``, resampled to the specified
``newShape``.
The space that the image is resampled into can be defined in one of the
following ways, in decreasing order of precedence:
1. If a ``matrix`` is provided, it is applied to the voxel coordinates
when retrieving values from the ``image``
2. Otherwise the image is simply scaled according to the ratio calculated
by ``image.shape / newShape``. In this case the ``origin`` argument
may be used to adjust the alignemnt of the original and resampled
voxel grids.
See the ``scipy.ndimage.affine_transform`` function for more details,
particularly on the ``order``, ``matrix``, ``mode`` and
``cval`` arguments.
:arg newShape: Desired shape. May containg floating point values, in which
case the resampled image will have shape
``round(newShape)``, but the voxel sizes will have scales
``self.shape / newShape`` (unless ``matrix`` is specified).
:arg sliceobj: Slice into this ``Image``. If ``None``, the whole
image is resampled, and it is assumed that it has the
same number of dimensions as ``newShape``. A
:exc:`ValueError` is raised if this is not the case.
:arg dtype: ``numpy`` data type of the resampled data. If ``None``,
the :meth:`dtype` of this ``Image`` is used.
:arg order: Spline interpolation order, passed through to the
``scipy.ndimage.affine_transform`` function - ``0``
corresponds to nearest neighbour interpolation, ``1``
(the default) to linear interpolation, and ``3`` to
cubic interpolation.
:arg smooth: If ``True`` (the default), the data is smoothed before
being resampled, but only along axes which are being
down-sampled (i.e. where ``newShape[i] < self.shape[i]``).
:arg origin: ``'centre'`` (the default) or ``'corner'``. ``'centre'``
resamples the image such that the centre of the corner
voxels of this image and the resampled data are
aligned. ``'corner'`` resamples the image such that
the corner of the corner voxels are aligned (and
therefore the voxel grids are aligned).
Ignored if ``offset`` or ``matrix`` is specified.
:arg matrix: Arbitrary affine transformation matrix to apply to the
voxel coordinates of ``image`` when resampling.
:arg mode: How to handle regions which are outside of the image FOV.
Defaults to `''nearest'``.
:arg cval: Constant value to use when ``mode='constant'``.
:returns: A tuple containing:
- A ``numpy`` array of shape ``newShape``, containing
an interpolated copy of the data in this ``Image``.
- A ``numpy`` array of shape ``(4, 4)``, containing the
adjusted voxel-to-world transformation for the spatial
dimensions of the resampled data.
"""
if sliceobj is None: sliceobj = slice(None)
if dtype is None: dtype = image.dtype
if origin == 'center': origin = 'centre'
if origin not in ('centre', 'corner'):
raise ValueError('Invalid value for origin: {}'.format(origin))
data = np.array(image[sliceobj], dtype=dtype, copy=False)
if len(data.shape) != len(newShape):
raise ValueError('Data dimensions do not match new shape: '
'len({}) != len({})'.format(data.shape, newShape))
# If matrix not provided, calculate
# a scaling/offset matrix from the
# old/new shape ratio and the origin
# setting.
if matrix is None:
matrix = calculateMatrix(data.shape, newShape, origin)
# calculateMatrix will return None
# if it decides that the image
# doesn't need to be resampled
if matrix is None:
return data, image.voxToWorldMat
newShape = np.array(np.round(newShape), dtype=np.int)
# Apply smoothing if requested,
# and if not using nn interp
if order > 0 and smooth:
data = applySmoothing(data, matrix, newShape)
# Do the resample thing
data = ndimage.affine_transform(data,
matrix,
output_shape=newShape,
order=order,
mode=mode,
cval=cval)
# Construct an affine transform which
# puts the resampled image into the
# same world coordinate system as this
# image. The calculateMatrix function
# might not return a 4x4 matrix, so we
# make sure it is valid.
if matrix.shape != (4, 4):
matrix = np.vstack((matrix[:3, :4], [0, 0, 0, 1]))
matrix = transform.concat(image.voxToWorldMat, matrix)
return data, matrix
def applySmoothing(data, matrix, newShape):
"""Called by the :func:`resample` function.
If interpolating and smoothing, we apply a gaussian filter along axes with
a resampling ratio greater than 1.1. We do this so that interpolation has
an effect when down-sampling to a resolution where the voxel centres are
aligned (as otherwise any interpolation regime will be equivalent to