Add target resolution parameter to GenerateSamplingReference

niworkflows/interfaces/nibabel.py

@@ -385,6 +385,12 @@ class _GenerateSamplingReferenceInputSpec(BaseInterfaceInputSpec):
         'the volume extent given by fixed_image, fast forward '
         'fixed_image otherwise.',
     )
+    target_resolution = traits.Either(
+        None,
+        traits.Tuple(traits.Float, traits.Float, traits.Float),
+        desc='target resolution (mm)',
+        usedefault=True,
+    )
 
 
 class _GenerateSamplingReferenceOutputSpec(TraitedSpec):
@@ -423,6 +429,7 @@ def _run_interface(self, runtime):
             force_xform_code=self.inputs.xform_code,
             message=f'{self.__class__.__name__} (niworkflows v{__version__})',
             newpath=runtime.cwd,
+            target_resolution=self.inputs.target_resolution,
         )
         return runtime
 
@@ -569,6 +576,42 @@ def reorient_image(img: nb.spatialimages.SpatialImage, target_ornt: str):
     return r_img
 
 
+def _calculate_target_affine(base_img, target_resolution):
+    """Calculate the target affine and shape for a given base image and target resolution.
+
+    Parameters
+    ----------
+    base_img : nibabel.SpatialImage
+        The base image to calculate the target affine and shape for.
+    target_resolution : tuple of 3 floats
+        The target resolution to calculate the target affine and shape for.
+
+    Returns
+    -------
+    new_affine : 4x4 numpy.ndarray
+        The target affine.
+    new_shape : tuple of 3 ints
+        The target shape.
+    """
+    import numpy as np
+
+    if len(target_resolution) != 3:
+        raise ValueError('target_resolution must be a tuple of 3 floats')
+
+    # determine appropriate shape
+    zooms = np.array(base_img.header.get_zooms())[:3]
+    ratios = zooms / np.array(target_resolution)
+    new_shape = np.array(base_img.shape) * ratios
+    new_shape = tuple(np.round(new_shape).astype(int))
+
+    # patch in voxel sizes to affine
+    new_affine = base_img.affine.copy()
+    for i in range(3):
+        new_affine[i, i] = target_resolution[i]
+
+    return new_affine, new_shape
+
+
 def _gen_reference(
     fixed_image,
     moving_image,
@@ -577,6 +620,7 @@ def _gen_reference(
     message=None,
     force_xform_code=None,
     newpath=None,
+    target_resolution=None,
 ):
     """Generate a sampling reference, and makes sure xform matrices/codes are correct."""
     import nilearn.image as nli
@@ -586,14 +630,24 @@ def _gen_reference(
 
     # Moving images may not be RAS/LPS (more generally, transverse-longitudinal-axial)
     reoriented_moving_img = nb.as_closest_canonical(nb.load(moving_image))
-    new_zooms = reoriented_moving_img.header.get_zooms()[:3]
 
-    # Avoid small differences in reported resolution to cause changes to
-    # FOV. See https://github.com/nipreps/fmriprep/issues/512
-    # A positive diagonal affine is RAS, hence the need to reorient above.
-    new_affine = np.diag(np.round(new_zooms, 3))
-
-    resampled = nli.resample_img(fixed_image, target_affine=new_affine, interpolation='nearest')
+    if target_resolution is not None:
+        new_affine, new_shape = _calculate_target_affine(reoriented_moving_img, target_resolution)
+    else:
+        new_zooms = reoriented_moving_img.header.get_zooms()[:3]
+
+        # Avoid small differences in reported resolution to cause changes to
+        # FOV. See https://github.com/nipreps/fmriprep/issues/512
+        # A positive diagonal affine is RAS, hence the need to reorient above.
+        new_affine = np.diag(np.round(new_zooms, 3))
+        new_shape = fixed_image.shape[:3]
+
+    resampled = nli.resample_img(
+        fixed_image,
+        target_affine=new_affine,
+        target_shape=new_shape,
+        interpolation='nearest',
+    )
 
     if fov_mask is not None:
         # If we have a mask, resample again dropping (empty) samples