From 3388b5e76bcc16a96f7465a47621a1e2b62dcdc3 Mon Sep 17 00:00:00 2001 From: Florian Rupprecht Date: Tue, 19 May 2026 16:13:52 -0400 Subject: [PATCH 1/2] Resample BOLD with nitransforms instead of per-volume ANTs Closes #112. Replaces the fan-out of antsApplyTransforms per volume with a single nitransforms + scipy.ndimage.map_coordinates pass: the static chain (anat->template warp, BBR affine, optional distortion) is composed into one template-grid coord map up front, then we loop over volumes applying the per-volume motion affine and sampling. Same idea applied to the longitudinal func_transform, so split_4d goes away. Distortion warps now write the canonical ANTs 5D (X, Y, Z, 1, 3) shape so DenseFieldTransform.from_filename can load them directly. --- pyproject.toml | 1 + src/rbc/core/common.py | 28 +- src/rbc/core/functional/distortion.py | 30 +- src/rbc/core/functional/resampling.py | 366 ++++++++++++------ src/rbc/core/longitudinal/transform.py | 71 +--- src/rbc/workflows/longitudinal/functional.py | 8 +- .../integration/functional/test_resampling.py | 29 +- tests/unit/core/test_longitudinal.py | 37 +- tests/unit/test_resampling.py | 220 ++++++++++- uv.lock | 60 +++ 10 files changed, 566 insertions(+), 284 deletions(-) diff --git a/pyproject.toml b/pyproject.toml index 4ece700f..baa59881 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -9,6 +9,7 @@ requires-python = ">=3.12" dependencies = [ "bids2table>=2.1.2", "nibabel>=5.3.3", + "nitransforms>=24.0.0", "niwrap>=0.9.3", "numpy>=2.4.2", "polars>=1.38.1", diff --git a/src/rbc/core/common.py b/src/rbc/core/common.py index fd3f9d61..3c51c5dc 100644 --- a/src/rbc/core/common.py +++ b/src/rbc/core/common.py @@ -2,7 +2,7 @@ Currently provides: - Deobliquing and RPI reorientation (initial preprocessing for T1w and BOLD). -- 4D NIfTI splitting and merging utilities. +- 4D NIfTI merging utility. """ from __future__ import annotations @@ -19,9 +19,8 @@ from rbc.core.fileops import file_tmp_copy from rbc.core.nifti import strip_afni_volatile_metadata -from rbc.core.niwrap import generate_exec_folder -__all__ = ["deoblique_and_reorient", "merge_3d_to_4d", "split_4d"] +__all__ = ["deoblique_and_reorient", "merge_3d_to_4d"] def deoblique_and_reorient( @@ -55,29 +54,6 @@ def deoblique_and_reorient( ) -def split_4d(img_4d: Path) -> list[Path]: - """Split a 4D NIfTI timeseries into individual 3D volumes. - - Volumes are written as uncompressed NIfTI (.nii) to avoid gzip - overhead on float intermediates that are read back immediately. - - Args: - img_4d: Path to a 4D NIfTI image. - - Returns: - Sorted list of paths to the individual 3D volume files. - """ - img = nib.nifti1.load(img_4d) - volumes = nib.four_to_three(img) - out_dir = generate_exec_folder(suffix="split4d") - paths: list[Path] = [] - for idx, vol in enumerate(volumes): - out_path = out_dir / f"vol_{idx:04d}.nii" - nib.save(vol, out_path) - paths.append(out_path) - return paths - - def merge_3d_to_4d(volumes: Sequence[Path], output: Path) -> Path: """Merge a sequence of 3D NIfTI volumes into a single 4D timeseries. diff --git a/src/rbc/core/functional/distortion.py b/src/rbc/core/functional/distortion.py index 2e2c5a30..3c433203 100644 --- a/src/rbc/core/functional/distortion.py +++ b/src/rbc/core/functional/distortion.py @@ -23,7 +23,7 @@ from niwrap import fsl from scipy.ndimage import binary_erosion -from rbc.core.functional.resampling import merge_3d_to_4d +from rbc.core.common import merge_3d_to_4d from rbc.core.niwrap import generate_exec_folder if TYPE_CHECKING: @@ -187,7 +187,9 @@ def _shiftmap_to_itk_warp( FUGUE shift maps contain per-voxel shifts (in voxels) along the PE axis. This converts to a 3-component mm displacement field in LPS - (ANTs/ITK convention). + (ANTs/ITK convention), written in the canonical ANTs 5D shape + ``(X, Y, Z, 1, 3)`` so it loads cleanly via + :class:`nitransforms.io.itk.ITKDisplacementsField`. Args: shiftmap: FUGUE shift map NIfTI (voxel units, single component). @@ -202,17 +204,15 @@ def _shiftmap_to_itk_warp( voxel_sizes = np.array(shift_img.header.get_zooms()[:3]) axis = _pe_axis_index(pe_direction) + warp = np.zeros((*shift_data.shape[:3], 1, 3), dtype=np.float32) + warp[..., 0, axis] = shift_data * voxel_sizes[axis] - # Build 3-component displacement field (mm, RAS) - warp = np.zeros((*shift_data.shape[:3], 3), dtype=np.float32) - warp[..., axis] = shift_data * voxel_sizes[axis] - - # RAS → LPS for ANTs - warp[..., 0] *= -1 - warp[..., 1] *= -1 + # RAS -> LPS for ANTs + warp[..., 0, 0] *= -1 + warp[..., 0, 1] *= -1 header = shift_img.header.copy() - header.set_intent(1007) # NIFTI_INTENT_VECTOR + header.set_intent("vector") header.set_data_dtype(np.float32) nib.save(nib.Nifti1Image(warp, shift_img.affine, header), output) @@ -249,15 +249,15 @@ def _fieldmap_hz_to_itk_warp( axis = _pe_axis_index(pe_direction) sign = _pe_sign(pe_direction) - warp = np.zeros((*fmap_data.shape[:3], 3), dtype=np.float32) - warp[..., axis] = fmap_data * readout_time * voxel_sizes[axis] * sign + warp = np.zeros((*fmap_data.shape[:3], 1, 3), dtype=np.float32) + warp[..., 0, axis] = fmap_data * readout_time * voxel_sizes[axis] * sign # RAS -> LPS for ANTs - warp[..., 0] *= -1 - warp[..., 1] *= -1 + warp[..., 0, 0] *= -1 + warp[..., 0, 1] *= -1 header = fmap_img.header.copy() - header.set_intent(1007) + header.set_intent("vector") header.set_data_dtype(np.float32) nib.save(nib.Nifti1Image(warp, fmap_img.affine, header), output) diff --git a/src/rbc/core/functional/resampling.py b/src/rbc/core/functional/resampling.py index 7c8e39a1..652f4937 100644 --- a/src/rbc/core/functional/resampling.py +++ b/src/rbc/core/functional/resampling.py @@ -1,13 +1,22 @@ """Resampling utilities for BOLD timeseries. -Provides two functions: - -- :func:`apply_motion_transforms`: applies per-volume mcflirt affines to - STC volumes, producing desc-preproc_bold (MC + STC in native space). - This is an exported derivative only, not consumed by downstream steps. -- :func:`resample_bold_to_template`: single-step resampling of STC BOLD - to template space, applying motion + BBR + anat-to-template transforms - in one interpolation pass per volume. +Three public entry points: + +- :func:`apply_motion_transforms`: per-volume mcflirt affines applied to + an STC timeseries (desc-preproc_bold derivative; not consumed + downstream). +- :func:`resample_bold_to_template`: single-step BOLD->template + resampling that fuses motion + BBR + anat-to-template (+ optional + distortion) into one interpolation pass per volume. +- :func:`resample_image`: 3D/4D resampling through a single composite + displacement field; used by the longitudinal pipeline. + +The static (time-invariant) parts of the transform chain are composed +into a single template-grid RAS coordinate map up front; the per-volume +loop then applies only the motion affine and samples with +``scipy.ndimage.map_coordinates``. This is conceptually similar to +``nitransforms.resampling.apply`` with ``serialize_nvols=1`` but lets +us cache the static coord map once across all volumes. """ from __future__ import annotations @@ -15,31 +24,179 @@ from typing import TYPE_CHECKING import nibabel as nib +import nitransforms as nt import numpy as np +from nitransforms.base import ImageGrid +from nitransforms.nonlinear import DenseFieldTransform +from scipy import ndimage as ndi + +from rbc.core.niwrap import generate_exec_folder if TYPE_CHECKING: from pathlib import Path -from niwrap import ants -from rbc.core.common import merge_3d_to_4d, split_4d -from rbc.core.fsl2itk import mat_to_itk -from rbc.core.niwrap import generate_exec_folder +# Interpolation defaults shared by both resampling paths. Order 3 is the +# scipy.ndimage equivalent of the cubic spline used by fMRIPrep. +_INTERP_ORDER = 3 +_INTERP_MODE = "constant" +_INTERP_CVAL = 0.0 + + +def _load_ants_warp(path: Path) -> DenseFieldTransform: + """Load an ANTs/ITK NIfTI displacement field as a DenseFieldTransform. + + Delegates to nitransforms' ``ITKDisplacementsField``, which enforces + the canonical 5D ``(X, Y, Z, 1, 3)`` ANTs shape and handles the LPS + -> RAS sign flip. + """ + return DenseFieldTransform.from_filename(str(path), fmt="itk") + + +def _ras_to_voxel_grid( + coords_ras: np.ndarray, affine_inv: np.ndarray, shape: tuple[int, ...] +) -> np.ndarray: + """Convert (N, 3) RAS coordinates to a (3, *shape) voxel-coord grid.""" + n = coords_ras.shape[0] + homog = np.empty((4, n), dtype=np.float64) + homog[:3] = coords_ras.T + homog[3] = 1.0 + voxel = (affine_inv @ homog)[:3] + return voxel.reshape((3, *shape)) + + +def _resample_4d( + src_img: nib.Nifti1Image, + reference_img: nib.Nifti1Image, + static_coords_ras: np.ndarray, + out_path: Path, + motion_xfms: list[nt.linear.Affine] | None = None, + order: int = _INTERP_ORDER, +) -> Path: + """Resample a 3D or 4D image volume-by-volume into reference space. + + For each output volume, applies the (optional) per-volume motion-affine + inverse to the static RAS coordinate map and samples with + ``ndi.map_coordinates``. Output dimensionality matches the input + (3D in, 3D out; 4D in, 4D out). The source's TR is preserved for 4D + outputs; the spatial header comes from *reference_img*. + + Args: + src_img: 3D or 4D source image. + reference_img: Image whose grid defines the output space. + static_coords_ras: ``(N, 3)`` RAS coordinates obtained by pulling + the reference grid through every static transform; ``N`` must + equal ``prod(reference_img.shape[:3])``. + out_path: Destination NIfTI path. + motion_xfms: Optional per-volume forward (volume->reference) + affines. If ``None``, no per-volume affine is composed in. + order: Spline interpolation order (0-5); default cubic. + + Returns: + Path to the written NIfTI. + """ + ref_shape = reference_img.shape[:3] + is_4d = src_img.ndim == 4 + n_vols = src_img.shape[3] if is_4d else 1 + src_affine_inv = np.linalg.inv(src_img.affine).astype(np.float64) + + if motion_xfms is not None and len(motion_xfms) != n_vols: + raise ValueError( + f"Count mismatch: ({len(motion_xfms)}) motion mats, ({n_vols}) volumes" + ) + + static_voxel_grid: np.ndarray | None = None + if motion_xfms is None: + static_voxel_grid = _ras_to_voxel_grid( + static_coords_ras, src_affine_inv, ref_shape + ) + + out_data = np.zeros((*ref_shape, n_vols), dtype=np.float32) + dataobj = src_img.dataobj + + for t in range(n_vols): + if motion_xfms is not None: + vol_coords_ras = motion_xfms[t].map(static_coords_ras, inverse=True) + voxel_grid = _ras_to_voxel_grid( + vol_coords_ras, src_affine_inv, ref_shape + ) + else: + assert static_voxel_grid is not None # noqa: S101 + voxel_grid = static_voxel_grid + + vol_data = np.asanyarray( + dataobj[..., t] if is_4d else dataobj, dtype=np.float32 + ) + ndi.map_coordinates( + vol_data, + voxel_grid, + output=out_data[..., t], + order=order, + mode=_INTERP_MODE, + cval=_INTERP_CVAL, + prefilter=True, + ) + + final_data = out_data if is_4d else out_data[..., 0] + out_img = nib.Nifti1Image(final_data, reference_img.affine) + if is_4d: + zooms = ( + reference_img.header.get_zooms()[:3] + src_img.header.get_zooms()[3:4] + ) + else: + zooms = reference_img.header.get_zooms()[:3] + out_img.header.set_zooms(zooms) + nib.save(out_img, out_path) + return out_path -def _restore_tr(resampled: Path, source: Path) -> None: - """Copy pixdim[4] (TR) from *source* into *resampled* NIfTI header. +def resample_image( + src: Path, reference: Path, warp: Path, order: int = _INTERP_ORDER +) -> Path: + """Resample a 3D or 4D image into *reference* space via an ANTs warp. + + Loads *warp* as an ANTs composite displacement field (the pull mapping + from reference to source) and applies it to every volume of *src*. + + Args: + src: 3D or 4D source image to be resampled. + reference: Image whose grid defines the output space. + warp: ANTs/ITK composite displacement field NIfTI. + order: Spline interpolation order (0-5); default cubic. - antsApplyTransforms sets pixdim from the reference image, which for - template-space outputs is a 3D template with no meaningful TR. This - restores the original temporal zoom from the source BOLD. + Returns: + Path to the resampled NIfTI in *reference* space. """ - src_img = nib.nifti1.load(source) - res_img = nib.nifti1.load(resampled) - data = np.asarray(res_img.dataobj) - zooms = res_img.header.get_zooms()[:3] + src_img.header.get_zooms()[3:] - res_img.header.set_zooms(zooms) - nib.save(nib.Nifti1Image(data, res_img.affine, res_img.header), resampled) + src_img = nib.nifti1.load(src) + ref_img = nib.nifti1.load(reference) + warp_xfm = _load_ants_warp(warp) + + ref_coords_ras = ImageGrid(ref_img).ndcoords.astype(np.float64) + src_coords_ras = warp_xfm.map(ref_coords_ras) + + out_path = generate_exec_folder("resample_image") / "resampled.nii.gz" + return _resample_4d( + src_img=src_img, + reference_img=ref_img, + static_coords_ras=src_coords_ras, + out_path=out_path, + order=order, + ) + + +def _load_motion_xfms( + motion_mat_dir: Path, bold_ref_img: nib.Nifti1Image +) -> list[nt.linear.Affine]: + """Load mcflirt per-volume motion affines from a directory of .mat files.""" + motion_mats = sorted(motion_mat_dir.glob("MAT_*")) + if not motion_mats: + raise FileNotFoundError(f"No motion .mat files found in {motion_mat_dir}") + return [ + nt.linear.load( + str(m), fmt="fsl", reference=bold_ref_img, moving=bold_ref_img + ) + for m in motion_mats + ] def apply_motion_transforms( @@ -60,58 +217,30 @@ def apply_motion_transforms( Args: stc_img: Slice-timing corrected 4D BOLD timeseries. motion_mat_dir: Directory of per-volume MAT_* matrices from mcflirt. - bold_ref: BOLD reference volume (used as both source and reference - for FSL-to-ITK matrix conversion). + bold_ref: BOLD reference volume (used both as FSL reference and as + the output spatial grid). Returns: 4D BOLD (MC + STC) in native space. Raises: FileNotFoundError: No motion .mat files are found in the specified directory. - ValueError: Number of motion matrix files found does not match the number + ValueError: Number of motion matrix files does not match the number of slice-timing corrected volumes. """ - motion_mats = sorted(motion_mat_dir.glob("MAT_*")) - if not motion_mats: - raise FileNotFoundError(f"No motion .mat files found in {motion_mat_dir}") - - stc_vols = split_4d(stc_img) - - if len(motion_mats) != len(stc_vols): - raise ValueError( - f"Count mismatch: ({len(motion_mats)}) mats, ({len(stc_vols)}) volumes" - ) - - transformed_vols = [] - motion_dir = generate_exec_folder(suffix="motionTransforms") - for idx, (motion_mat, stc_vol) in enumerate( - zip(motion_mats, stc_vols, strict=True) - ): - motion_itk = mat_to_itk( - motion_mat, bold_ref, bold_ref, motion_dir / f"motion_{idx:04d}.txt" - ) - result = ants.ants_apply_transforms( - input_image=stc_vol, - reference_image=bold_ref, # bold in native space - transform=[ants.ants_apply_transforms_transform_file_name(motion_itk)], - interpolation=ants.ants_apply_transforms_lanczos_windowed_sinc(), - float_=True, - default_value=0, - dimensionality=3, - output=ants.ants_apply_transforms_warped_output( - f"vol_{idx:04d}_motion.nii" - ), - ) - transformed_vols.append(result.output.output_image_outfile) + src_img = nib.nifti1.load(stc_img) + bold_ref_img = nib.nifti1.load(bold_ref) + motion_xfms = _load_motion_xfms(motion_mat_dir, bold_ref_img) + static_coords_ras = ImageGrid(bold_ref_img).ndcoords.astype(np.float64) out_path = generate_exec_folder("preproc_bold_merge") / "preproc_bold.nii.gz" - merged = merge_3d_to_4d(transformed_vols, out_path) - - # antsApplyTransforms writes the reference image's pixdim into the output - # header, so pixdim[4] (TR) is lost. Restore it from the source BOLD. - _restore_tr(merged, stc_img) - - return merged + return _resample_4d( + src_img=src_img, + reference_img=bold_ref_img, + static_coords_ras=static_coords_ras, + motion_xfms=motion_xfms, + out_path=out_path, + ) def resample_bold_to_template( @@ -126,19 +255,20 @@ def resample_bold_to_template( ) -> Path: """Single-step resampling of STC BOLD to template space. - Applies all spatial transforms (motion + BBR + anat-to-template, and - optionally distortion correction) in a single ``antsApplyTransforms`` - call per volume. This avoids multiple interpolation passes. + Composes the static transforms (anat-to-template warp, BBR affine, + optional distortion warp) into a single template-grid RAS coordinate + map, then loops over volumes applying the per-volume motion affine + and sampling with cubic spline interpolation. Args: stc_bold: Slice-timing corrected 4D BOLD timeseries. motion_mat_dir: Directory of per-volume MAT_* matrices from mcflirt. - bold_to_anat: BOLD to T1w affine (output from BBR). - anat_to_template: T1w to template composite warp. - bold_ref: BOLD reference volume (used for ITK conversion). + bold_to_anat: BOLD to T1w FSL affine matrix (BBR output). + anat_to_template: T1w to template ANTs composite displacement field. + bold_ref: BOLD reference volume (FSL ref/moving for affines). template: Brain template in target space. - t1w_brain: Skull-stripped T1w brain. - distortion_warp: Optional distortion correction warp. + t1w_brain: Skull-stripped T1w brain (FSL reference for bold_to_anat). + distortion_warp: Optional ANTs/ITK displacement field on the BOLD grid. Returns: Resampled 4D BOLD in template space. @@ -147,66 +277,46 @@ def resample_bold_to_template( FileNotFoundError: No motion .mat files found in the directory. ValueError: Number of motion matrices does not match STC volumes. """ - motion_mats = sorted(motion_mat_dir.glob("MAT_*")) - if not motion_mats: - raise FileNotFoundError(f"No motion .mat files found in {motion_mat_dir}") - - resample_dir = generate_exec_folder(suffix="boldToTemplateResample") - bold2anat_itk = mat_to_itk( - bold_to_anat, t1w_brain, bold_ref, resample_dir / "bold2anat.txt" - ) - - stc_vols = split_4d(stc_bold) - - if len(motion_mats) != len(stc_vols): + src_img = nib.nifti1.load(stc_bold) + template_img = nib.nifti1.load(template) + bold_ref_img = nib.nifti1.load(bold_ref) + t1w_img = nib.nifti1.load(t1w_brain) + + if src_img.ndim != 4: + raise ValueError(f"Expected 4D STC BOLD, got shape {src_img.shape}") + + # Fail fast before loading the (slow) static warps. + motion_xfms = _load_motion_xfms(motion_mat_dir, bold_ref_img) + n_vols = src_img.shape[3] + if len(motion_xfms) != n_vols: raise ValueError( - f"Count mismatch: ({len(motion_mats)}) mats, ({len(stc_vols)}) volumes" + f"Count mismatch: ({len(motion_xfms)}) mats, ({n_vols}) volumes" ) - # Shared transforms (applied to every volume) - base_transforms = [anat_to_template, bold2anat_itk] - if distortion_warp: - base_transforms.append(distortion_warp) - - # Per-volume: all transforms in one antsApplyTransforms call - # Order (last applied first): motion -> distortion -> bold2anat -> anat2template - transformed_vols = [] - for idx, (motion_mat, stc_vol) in enumerate( - zip(motion_mats, stc_vols, strict=True) - ): - motion_itk = mat_to_itk( - motion_mat, bold_ref, bold_ref, resample_dir / f"motion_{idx:04d}.txt" - ) - transforms: list[ - ants.AntsApplyTransformsTransformFileNameParamsDictTagged - | ants.AntsApplyTransformsUseInverseParamsDictTagged - ] = [ - ants.ants_apply_transforms_transform_file_name(t) - for t in [*base_transforms, motion_itk] - ] - result = ants.ants_apply_transforms( - input_image=stc_vol, - reference_image=template, - transform=transforms, - interpolation=ants.ants_apply_transforms_lanczos_windowed_sinc(), - float_=True, - default_value=0, - dimensionality=3, - output=ants.ants_apply_transforms_warped_output( - f"vol_{idx:04d}_template.nii" - ), - ) - transformed_vols.append(result.output.output_image_outfile) + bold_to_anat_xfm = nt.linear.load( + str(bold_to_anat), fmt="fsl", reference=t1w_img, moving=bold_ref_img + ) + anat_to_tpl_xfm = _load_ants_warp(anat_to_template) + distortion_xfm = ( + _load_ants_warp(distortion_warp) if distortion_warp is not None else None + ) + + # Static pull chain: template_ras -> anat_ras -> bold_undistorted_ras + # (-> bold_distorted_ras if distortion). + tpl_coords_ras = ImageGrid(template_img).ndcoords.astype(np.float64) + anat_coords_ras = anat_to_tpl_xfm.map(tpl_coords_ras) + bold_coords_ras = bold_to_anat_xfm.map(anat_coords_ras, inverse=True) + if distortion_xfm is not None: + bold_coords_ras = distortion_xfm.map(bold_coords_ras) out_path = ( generate_exec_folder("bold_to_template_merge") / "bold_to_template_resampled.nii.gz" ) - merged = merge_3d_to_4d(transformed_vols, out_path) - - # antsApplyTransforms writes the reference (template) image's pixdim into - # the output header, so pixdim[4] (TR) is lost. Restore it from the - # source BOLD. - _restore_tr(merged, stc_bold) - - return merged + return _resample_4d( + src_img=src_img, + reference_img=template_img, + static_coords_ras=bold_coords_ras, + motion_xfms=motion_xfms, + out_path=out_path, + ) diff --git a/src/rbc/core/longitudinal/transform.py b/src/rbc/core/longitudinal/transform.py index 7a147271..62065f24 100644 --- a/src/rbc/core/longitudinal/transform.py +++ b/src/rbc/core/longitudinal/transform.py @@ -2,13 +2,11 @@ from __future__ import annotations -from typing import TYPE_CHECKING, Literal +from typing import TYPE_CHECKING from niwrap import ants -from rbc.core.common import merge_3d_to_4d, split_4d -from rbc.core.functional.resampling import _restore_tr -from rbc.core.niwrap import generate_exec_folder +from rbc.core.functional.resampling import resample_image if TYPE_CHECKING: from pathlib import Path @@ -73,83 +71,26 @@ def compose_transform(ref: Path, bold_to_anat_itk: Path, anat_to_tpl_xfm: Path) ).output.output_image_outfile -def func_transform( - in_file: Path, - template: Path, - xfm: Path, - strategy: Literal["single", "chunked"] = "chunked", -) -> Path: - """Apply transformation to functional data using ANTs. +def func_transform(in_file: Path, template: Path, xfm: Path) -> Path: + """Apply *xfm* to a 3D or 4D functional image with linear interpolation. Args: in_file: Input file path to apply transform to. template: Longitudinal template file path for reference. - xfm: Transformation to apply. - strategy: Transformation strategy to apply. + xfm: ANTs/ITK composite displacement field. Returns: Path to transformed file. Raises: FileNotFoundError: if input file or transformation not found. - ValueError: if strategy is unknown (not 'single' or 'chunked') """ if not in_file.exists(): raise FileNotFoundError(f"Input file not found: {in_file}") if not xfm.exists(): raise FileNotFoundError(f"Transformation not found: {xfm}") - match strategy: - case "chunked": - return _transform_4d_chunked(in_file=in_file, template=template, xfm=xfm) - case "single": - return _transform_4d(in_file=in_file, template=template, xfm=xfm) - case _: - raise ValueError( - f"Unknown strategy: {strategy!r}. Must be 'single' or 'chunked'" - ) - - -def _transform_4d(in_file: Path, template: Path, xfm: Path) -> Path: - """Apply transformation directly on 4D image.""" - return ants.ants_apply_transforms( - dimensionality=3, - reference_image=template, - input_image=in_file, - input_image_type=3, - transform=[ants.ants_apply_transforms_transform_file_name(xfm)], - output=ants.ants_apply_transforms_warped_output("subj_bold_to_template.nii.gz"), - interpolation=ants.ants_apply_transforms_linear(), - ).output.output_image_outfile - - -def _transform_4d_chunked(in_file: Path, template: Path, xfm: Path) -> Path: - """Apply transformation using a chunked strategy (split 4D image into 3D).""" - func_vols = split_4d(in_file) - - transformed_vols = [] - for idx, func_vol in enumerate(func_vols): - result = ants.ants_apply_transforms( - input_image=func_vol, - reference_image=template, - transform=[ants.ants_apply_transforms_transform_file_name(xfm)], - float_=True, - default_value=0, - dimensionality=3, - interpolation=ants.ants_apply_transforms_linear(), - output=ants.ants_apply_transforms_warped_output( - f"vol_{idx:04d}_template.nii" - ), - ) - transformed_vols.append(result.output.output_image_outfile) - - out_path = ( - generate_exec_folder("bold_to_longitudinal_merge") - / "bold_to_longitudinal.nii.gz" - ) - merged = merge_3d_to_4d(transformed_vols, out_path) - _restore_tr(merged, in_file) - return merged + return resample_image(src=in_file, reference=template, warp=xfm, order=1) def mask_transform(mask: Path, template: Path, xfm: Path) -> Path: diff --git a/src/rbc/workflows/longitudinal/functional.py b/src/rbc/workflows/longitudinal/functional.py index fdf7b669..5bae3370 100644 --- a/src/rbc/workflows/longitudinal/functional.py +++ b/src/rbc/workflows/longitudinal/functional.py @@ -88,12 +88,8 @@ def longitudinal_process( anat_to_tpl_xfm=anat_to_template_xfm, ) - long_sbref = func_transform( - in_file=sbref, template=template, xfm=bold_to_tpl_xfm, strategy="single" - ) - long_bold = func_transform( - in_file=bold, template=template, xfm=bold_to_tpl_xfm, strategy="chunked" - ) + long_sbref = func_transform(in_file=sbref, template=template, xfm=bold_to_tpl_xfm) + long_bold = func_transform(in_file=bold, template=template, xfm=bold_to_tpl_xfm) long_mask = mask_transform(mask=bold_mask, template=template, xfm=bold_to_tpl_xfm) regressed: dict[str, Path] = {} diff --git a/tests/integration/functional/test_resampling.py b/tests/integration/functional/test_resampling.py index 7b1afae4..47d0563d 100644 --- a/tests/integration/functional/test_resampling.py +++ b/tests/integration/functional/test_resampling.py @@ -54,20 +54,21 @@ def _create_synthetic_wm(t1w: Path) -> Path: return wm_file -def _create_identity_affine() -> Path: - """Create an ITK-format identity affine transform for testing.""" - out_dir = generate_exec_folder("synthetic_transform") - mat_file = out_dir / "identity_affine.txt" - - mat_file.write_text( - "#Insight Transform File V1.0\n" - "#Transform 0\n" - "Transform: MatrixOffsetTransformBase_double_3_3\n" - "Parameters: 1 0 0 0 1 0 0 0 1 0 0 0\n" - "FixedParameters: 0 0 0\n" - ) +def _create_identity_warp(reference: Path) -> Path: + """Create an identity ANTs/ITK displacement field on the reference grid. - return mat_file + The new nitransforms-based resampler expects ``anat_to_template`` to + be a composite displacement field (the format produced by + ``ants_apply_transforms`` in production), not a `.txt` affine. + """ + ref_img = nib.nifti1.load(reference) + warp = np.zeros((*ref_img.shape[:3], 1, 3), dtype=np.float32) + img = nib.Nifti1Image(warp, ref_img.affine) + img.header.set_intent("vector") + out_dir = generate_exec_folder("synthetic_transform") + out_path = out_dir / "identity_warp.nii.gz" + nib.save(img, out_path) + return out_path @pytest.mark.slow @@ -75,7 +76,7 @@ def test_resample_bold_to_template(test_subject: TestSubjectData) -> None: """Test single-step resampling of STC BOLD to template space.""" template_mni = REGISTRATION_TEMPLATES.brain_2mm synthetic_wm = _create_synthetic_wm(test_subject.t1w) - anat_to_template = _create_identity_affine() + anat_to_template = _create_identity_warp(template_mni) reoriented = deoblique_and_reorient(in_file=test_subject.bold) truncated = afni.v_3dcalc( diff --git a/tests/unit/core/test_longitudinal.py b/tests/unit/core/test_longitudinal.py index 54017e13..070537fe 100644 --- a/tests/unit/core/test_longitudinal.py +++ b/tests/unit/core/test_longitudinal.py @@ -88,43 +88,22 @@ def test_missing_xfm(self, tmp_files: tuple[Path, ...]) -> None: with pytest.raises(FileNotFoundError, match="not found"): func_transform(in_file=in_file, template=template, xfm=xfm) - @patch("rbc.core.longitudinal.transform.ants") - @patch("rbc.core.longitudinal.transform.split_4d") - @patch("rbc.core.longitudinal.transform.merge_3d_to_4d") - @patch("rbc.core.longitudinal.transform._restore_tr") - @pytest.mark.parametrize("strategy", ["chunked", "single"]) + @patch("rbc.core.longitudinal.transform.resample_image") def test_returns_output_path( self, - mock_restore_tr: MagicMock, - mock_merge_3d_to_4d: MagicMock, - mock_split_4d: MagicMock, - mock_ants: MagicMock, + mock_resample_image: MagicMock, tmp_files: tuple[Path, ...], - strategy: str, ) -> None: """Successful functional transformation to template.""" in_file, template, xfm = tmp_files - expected = Path("/out/subj_bold_to_template.nii.gz") + expected = Path("/out/resampled.nii.gz") + mock_resample_image.return_value = expected - mock_ants.ants_apply_transforms.return_value.output.output_image_outfile = ( - expected - ) - mock_ants.ants_apply_transforms_warped_output.return_value = MagicMock() - mock_ants.ants_apply_transforms_linear.return_value = MagicMock() - mock_ants.ants_apply_transforms_transform_file_name.return_value = MagicMock() - - if strategy == "chunked": - mock_split_4d.return_value = [in_file] - mock_merge_3d_to_4d.return_value = expected - mock_restore_tr.return_value = None - - result = func_transform( - in_file=in_file, - template=template, - xfm=xfm, - strategy=strategy, # type: ignore [arg-type] - ) + result = func_transform(in_file=in_file, template=template, xfm=xfm) assert result == expected + mock_resample_image.assert_called_once_with( + src=in_file, reference=template, warp=xfm, order=1 + ) class TestComposeTransform: diff --git a/tests/unit/test_resampling.py b/tests/unit/test_resampling.py index 1effad29..58b8855a 100644 --- a/tests/unit/test_resampling.py +++ b/tests/unit/test_resampling.py @@ -8,12 +8,45 @@ import numpy as np import pytest -from rbc.core.functional.resampling import merge_3d_to_4d +from rbc.core.common import merge_3d_to_4d +from rbc.core.functional.resampling import ( + apply_motion_transforms, + resample_bold_to_template, +) if TYPE_CHECKING: from pathlib import Path +def _write_identity_warp( + path: Path, shape: tuple[int, int, int], affine: np.ndarray +) -> Path: + """Write a zero-displacement ANTs/ITK warp at the given grid (5D format).""" + warp = np.zeros((*shape, 1, 3), dtype=np.float32) + img = nib.Nifti1Image(warp, affine) + img.header.set_intent(1007) # NIFTI_INTENT_VECTOR + nib.save(img, path) + return path + + +def _write_fsl_mat(path: Path, mat: np.ndarray) -> Path: + """Write a 4x4 matrix to disk in FSL .mat (whitespace-delimited) format.""" + np.savetxt(path, mat) + return path + + +def _make_bold( + path: Path, shape: tuple[int, int, int], n_vols: int, affine: np.ndarray, tr: float +) -> Path: + """Create a 4D BOLD NIfTI with reproducible random voxel values.""" + rng = np.random.default_rng(42) + data = rng.standard_normal((*shape, n_vols)).astype(np.float32) + img = nib.Nifti1Image(data, affine) + img.header.set_zooms((*affine[:3, :3].diagonal().tolist(), tr)) + nib.save(img, path) + return path + + def _make_3d_nifti(path: Path, shape: tuple = (4, 5, 6)) -> Path: """Create a minimal 3D NIfTI image with random data.""" rng = np.random.default_rng(42) @@ -67,3 +100,188 @@ def test_empty_list_raises(self, tmp_path: Path) -> None: """Passing an empty list raises ValueError.""" with pytest.raises(ValueError, match="empty"): merge_3d_to_4d([], tmp_path / "merged.nii.gz") + + +class TestApplyMotionTransforms: + """Tests for apply_motion_transforms (resamples in BOLD native space).""" + + def test_identity_motion_is_roundtrip(self, tmp_path: Path) -> None: + """Identity motion mats reproduce the input STC data within interp tolerance.""" + shape = (8, 9, 7) + n_vols = 3 + affine = np.diag([2.0, 2.0, 2.0, 1.0]) + bold_path = _make_bold(tmp_path / "stc.nii.gz", shape, n_vols, affine, tr=1.5) + + ref_path = tmp_path / "bold_ref.nii.gz" + src_data = np.asanyarray(nib.nifti1.load(bold_path).dataobj) + nib.save(nib.Nifti1Image(src_data[..., 1], affine), ref_path) + + mat_dir = tmp_path / "mats" + mat_dir.mkdir() + for i in range(n_vols): + _write_fsl_mat(mat_dir / f"MAT_{i:04d}", np.eye(4)) + + out = apply_motion_transforms( + stc_img=bold_path, motion_mat_dir=mat_dir, bold_ref=ref_path + ) + out_img = nib.nifti1.load(out) + assert out_img.shape == (*shape, n_vols) + np.testing.assert_allclose(out_img.get_fdata(), src_data, atol=1e-3) + + def test_preserves_tr(self, tmp_path: Path) -> None: + """Output keeps the source BOLD's TR in pixdim[4].""" + affine = np.diag([2.0, 2.0, 2.0, 1.0]) + bold_path = _make_bold(tmp_path / "stc.nii.gz", (4, 5, 4), 2, affine, tr=2.7) + ref_path = tmp_path / "bold_ref.nii.gz" + ref_data = np.zeros((4, 5, 4), dtype=np.float32) + nib.save(nib.Nifti1Image(ref_data, affine), ref_path) + mat_dir = tmp_path / "mats" + mat_dir.mkdir() + for i in range(2): + _write_fsl_mat(mat_dir / f"MAT_{i:04d}", np.eye(4)) + + out = apply_motion_transforms( + stc_img=bold_path, motion_mat_dir=mat_dir, bold_ref=ref_path + ) + assert nib.nifti1.load(out).header.get_zooms()[3] == pytest.approx(2.7) + + def test_missing_mats_raises(self, tmp_path: Path) -> None: + """Empty motion directory raises FileNotFoundError.""" + affine = np.diag([2.0, 2.0, 2.0, 1.0]) + bold_path = _make_bold(tmp_path / "stc.nii.gz", (4, 5, 4), 2, affine, tr=1.0) + ref_path = tmp_path / "bold_ref.nii.gz" + ref_data = np.zeros((4, 5, 4), dtype=np.float32) + nib.save(nib.Nifti1Image(ref_data, affine), ref_path) + empty = tmp_path / "empty" + empty.mkdir() + with pytest.raises(FileNotFoundError, match=r"No motion \.mat files"): + apply_motion_transforms( + stc_img=bold_path, motion_mat_dir=empty, bold_ref=ref_path + ) + + def test_mat_count_mismatch_raises(self, tmp_path: Path) -> None: + """Mismatched motion-mat / volume counts raise ValueError.""" + affine = np.diag([2.0, 2.0, 2.0, 1.0]) + bold_path = _make_bold(tmp_path / "stc.nii.gz", (4, 5, 4), 3, affine, tr=1.0) + ref_path = tmp_path / "bold_ref.nii.gz" + ref_data = np.zeros((4, 5, 4), dtype=np.float32) + nib.save(nib.Nifti1Image(ref_data, affine), ref_path) + mat_dir = tmp_path / "mats" + mat_dir.mkdir() + for i in range(2): # one too few + _write_fsl_mat(mat_dir / f"MAT_{i:04d}", np.eye(4)) + with pytest.raises(ValueError, match="Count mismatch"): + apply_motion_transforms( + stc_img=bold_path, motion_mat_dir=mat_dir, bold_ref=ref_path + ) + + +class TestResampleBoldToTemplate: + """Tests for resample_bold_to_template (single-step BOLD->template).""" + + def _set_up( + self, tmp_path: Path, n_vols: int = 3, tr: float = 2.0 + ) -> dict[str, Path]: + """Build a minimal synthetic input set.""" + bold_shape = (8, 9, 7) + bold_affine = np.diag([2.0, 2.0, 2.0, 1.0]) + bold_path = _make_bold( + tmp_path / "stc.nii.gz", bold_shape, n_vols, bold_affine, tr=tr + ) + + ref_path = tmp_path / "bold_ref.nii.gz" + src_data = np.asanyarray(nib.nifti1.load(bold_path).dataobj) + nib.save(nib.Nifti1Image(src_data[..., 0], bold_affine), ref_path) + + t1w_affine = np.diag([1.5, 1.5, 1.5, 1.0]) + t1w_path = tmp_path / "t1w.nii.gz" + nib.save( + nib.Nifti1Image(np.zeros((10, 12, 9), dtype=np.float32), t1w_affine), + t1w_path, + ) + + tpl_shape = (6, 7, 5) + tpl_affine = np.diag([3.0, 3.0, 3.0, 1.0]) + tpl_path = tmp_path / "tpl.nii.gz" + nib.save( + nib.Nifti1Image(np.zeros(tpl_shape, dtype=np.float32), tpl_affine), + tpl_path, + ) + + bold2anat = _write_fsl_mat(tmp_path / "bold2anat.mat", np.eye(4)) + anat2tpl = _write_identity_warp( + tmp_path / "anat2tpl.nii.gz", tpl_shape, tpl_affine + ) + + mat_dir = tmp_path / "mats" + mat_dir.mkdir() + for i in range(n_vols): + _write_fsl_mat(mat_dir / f"MAT_{i:04d}", np.eye(4)) + + return { + "stc_bold": bold_path, + "motion_mat_dir": mat_dir, + "bold_to_anat": bold2anat, + "anat_to_template": anat2tpl, + "bold_ref": ref_path, + "template": tpl_path, + "t1w_brain": t1w_path, + } + + def test_output_shape_and_voxel_size_match_template(self, tmp_path: Path) -> None: + """Output spatial grid matches the template, time dim matches BOLD.""" + kwargs = self._set_up(tmp_path, n_vols=4) + out = resample_bold_to_template(**kwargs) + out_img = nib.nifti1.load(out) + tpl_img = nib.nifti1.load(kwargs["template"]) + assert out_img.shape == (*tpl_img.shape, 4) + np.testing.assert_array_equal(out_img.affine, tpl_img.affine) + assert out_img.header.get_zooms()[:3] == tpl_img.header.get_zooms()[:3] + + def test_preserves_tr_from_source(self, tmp_path: Path) -> None: + """Output pixdim[4] equals the BOLD source's TR, not the template's.""" + kwargs = self._set_up(tmp_path, tr=2.7) + out = resample_bold_to_template(**kwargs) + assert nib.nifti1.load(out).header.get_zooms()[3] == pytest.approx(2.7) + + def test_output_is_finite(self, tmp_path: Path) -> None: + """Resampled output contains no NaN or Inf voxels.""" + kwargs = self._set_up(tmp_path) + out = resample_bold_to_template(**kwargs) + data = nib.nifti1.load(out).get_fdata() + assert np.isfinite(data).all() + + def test_distortion_warp_accepted(self, tmp_path: Path) -> None: + """Optional distortion warp is composed in without crashing.""" + # A zero-displacement warp does NOT exactly reproduce the no-distortion + # result: DenseFieldTransform.map cubic-spline-samples the stored RAS + # coord field, which rings near the grid boundary. So we only check + # shape + finiteness here; numerical agreement is in integration. + kwargs = self._set_up(tmp_path) + bold_affine = nib.nifti1.load(kwargs["bold_ref"]).affine + distortion = _write_identity_warp( + tmp_path / "distortion.nii.gz", (8, 9, 7), bold_affine + ) + out = resample_bold_to_template(distortion_warp=distortion, **kwargs) + out_img = nib.nifti1.load(out) + tpl_img = nib.nifti1.load(kwargs["template"]) + assert out_img.shape == (*tpl_img.shape, 3) + assert np.isfinite(out_img.get_fdata()).all() + + def test_missing_mats_raises(self, tmp_path: Path) -> None: + """Empty motion directory raises FileNotFoundError.""" + kwargs = self._set_up(tmp_path) + empty = tmp_path / "empty" + empty.mkdir() + kwargs["motion_mat_dir"] = empty + with pytest.raises(FileNotFoundError, match=r"No motion \.mat files"): + resample_bold_to_template(**kwargs) + + def test_mat_count_mismatch_raises(self, tmp_path: Path) -> None: + """Mismatched motion-mat / volume counts raise ValueError.""" + kwargs = self._set_up(tmp_path, n_vols=3) + # Drop one mat so the count no longer matches. + for m in sorted(kwargs["motion_mat_dir"].glob("MAT_*"))[-1:]: + m.unlink() + with pytest.raises(ValueError, match="Count mismatch"): + resample_bold_to_template(**kwargs) diff --git a/uv.lock b/uv.lock index 8b77ed8c..f32f999c 100644 --- a/uv.lock +++ b/uv.lock @@ -313,6 +313,49 @@ wheels = [ { url = "https://files.pythonhosted.org/packages/81/47/dd9a212ef6e343a6857485ffe25bba537304f1913bdbed446a23f7f592e1/filelock-3.29.0-py3-none-any.whl", hash = "sha256:96f5f6344709aa1572bbf631c640e4ebeeb519e08da902c39a001882f30ac258", size = 39812, upload-time = "2026-04-19T15:39:08.752Z" }, ] +[[package]] +name = "h5py" +version = "3.16.0" +source = { registry = "https://pypi.org/simple" } +dependencies = [ + { name = "numpy" }, +] +sdist = { url = "https://files.pythonhosted.org/packages/db/33/acd0ce6863b6c0d7735007df01815403f5589a21ff8c2e1ee2587a38f548/h5py-3.16.0.tar.gz", hash = "sha256:a0dbaad796840ccaa67a4c144a0d0c8080073c34c76d5a6941d6818678ef2738", size = 446526, upload-time = "2026-03-06T13:49:08.07Z" } +wheels = [ + { url = "https://files.pythonhosted.org/packages/c8/c0/5d4119dba94093bbafede500d3defd2f5eab7897732998c04b54021e530b/h5py-3.16.0-cp312-cp312-macosx_10_13_x86_64.whl", hash = "sha256:c5313566f4643121a78503a473f0fb1e6dcc541d5115c44f05e037609c565c4d", size = 3685604, upload-time = "2026-03-06T13:48:04.198Z" }, + { url = "https://files.pythonhosted.org/packages/b0/42/c84efcc1d4caebafb1ecd8be4643f39c85c47a80fe254d92b8b43b1eadaf/h5py-3.16.0-cp312-cp312-macosx_11_0_arm64.whl", hash = "sha256:42b012933a83e1a558c673176676a10ce2fd3759976a0fedee1e672d1e04fc9d", size = 3061940, upload-time = "2026-03-06T13:48:05.783Z" }, + { url = "https://files.pythonhosted.org/packages/89/84/06281c82d4d1686fde1ac6b0f307c50918f1c0151062445ab3b6fa5a921d/h5py-3.16.0-cp312-cp312-manylinux_2_28_aarch64.whl", hash = "sha256:ff24039e2573297787c3063df64b60aab0591980ac898329a08b0320e0cf2527", size = 5198852, upload-time = "2026-03-06T13:48:07.482Z" }, + { url = "https://files.pythonhosted.org/packages/9e/e9/1a19e42cd43cc1365e127db6aae85e1c671da1d9a5d746f4d34a50edb577/h5py-3.16.0-cp312-cp312-manylinux_2_28_x86_64.whl", hash = "sha256:dfc21898ff025f1e8e67e194965a95a8d4754f452f83454538f98f8a3fcb207e", size = 5405250, upload-time = "2026-03-06T13:48:09.628Z" }, + { url = "https://files.pythonhosted.org/packages/b7/8e/9790c1655eabeb85b92b1ecab7d7e62a2069e53baefd58c98f0909c7a948/h5py-3.16.0-cp312-cp312-musllinux_1_2_aarch64.whl", hash = "sha256:698dd69291272642ffda44a0ecd6cd3bda5faf9621452d255f57ce91487b9794", size = 5190108, upload-time = "2026-03-06T13:48:11.26Z" }, + { url = "https://files.pythonhosted.org/packages/51/d7/ab693274f1bd7e8c5f9fdd6c7003a88d59bedeaf8752716a55f532924fbb/h5py-3.16.0-cp312-cp312-musllinux_1_2_x86_64.whl", hash = "sha256:2b2c02b0a160faed5fb33f1ba8a264a37ee240b22e049ecc827345d0d9043074", size = 5419216, upload-time = "2026-03-06T13:48:13.322Z" }, + { url = "https://files.pythonhosted.org/packages/03/c1/0976b235cf29ead553e22f2fb6385a8252b533715e00d0ae52ed7b900582/h5py-3.16.0-cp312-cp312-win_amd64.whl", hash = "sha256:96b422019a1c8975c2d5dadcf61d4ba6f01c31f92bbde6e4649607885fe502d6", size = 3182868, upload-time = "2026-03-06T13:48:15.759Z" }, + { url = "https://files.pythonhosted.org/packages/14/d9/866b7e570b39070f92d47b0ff1800f0f8239b6f9e45f02363d7112336c1f/h5py-3.16.0-cp312-cp312-win_arm64.whl", hash = "sha256:39c2838fb1e8d97bcf1755e60ad1f3dd76a7b2a475928dc321672752678b96db", size = 2653286, upload-time = "2026-03-06T13:48:17.279Z" }, + { url = "https://files.pythonhosted.org/packages/0f/9e/6142ebfda0cb6e9349c091eae73c2e01a770b7659255248d637bec54a88b/h5py-3.16.0-cp313-cp313-macosx_10_13_x86_64.whl", hash = "sha256:370a845f432c2c9619db8eed334d1e610c6015796122b0e57aa46312c22617d9", size = 3671808, upload-time = "2026-03-06T13:48:19.737Z" }, + { url = "https://files.pythonhosted.org/packages/b0/65/5e088a45d0f43cd814bc5bec521c051d42005a472e804b1a36c48dada09b/h5py-3.16.0-cp313-cp313-macosx_11_0_arm64.whl", hash = "sha256:42108e93326c50c2810025aade9eac9d6827524cdccc7d4b75a546e5ab308edb", size = 3045837, upload-time = "2026-03-06T13:48:21.854Z" }, + { url = "https://files.pythonhosted.org/packages/da/1e/6172269e18cc5a484e2913ced33339aad588e02ba407fafd00d369e22ef3/h5py-3.16.0-cp313-cp313-manylinux_2_28_aarch64.whl", hash = "sha256:099f2525c9dcf28de366970a5fb34879aab20491589fa89ce2863a84218bb524", size = 5193860, upload-time = "2026-03-06T13:48:24.071Z" }, + { url = "https://files.pythonhosted.org/packages/bd/98/ef2b6fe2903e377cbe870c3b2800d62552f1e3dbe81ce49e1923c53d1c5c/h5py-3.16.0-cp313-cp313-manylinux_2_28_x86_64.whl", hash = "sha256:9300ad32dea9dfc5171f94d5f6948e159ed93e4701280b0f508773b3f582f402", size = 5400417, upload-time = "2026-03-06T13:48:25.728Z" }, + { url = "https://files.pythonhosted.org/packages/bc/81/5b62d760039eed64348c98129d17061fdfc7839fc9c04eaaad6dee1004e4/h5py-3.16.0-cp313-cp313-musllinux_1_2_aarch64.whl", hash = "sha256:171038f23bccddfc23f344cadabdfc9917ff554db6a0d417180d2747fe4c75a7", size = 5185214, upload-time = "2026-03-06T13:48:27.436Z" }, + { url = "https://files.pythonhosted.org/packages/28/c4/532123bcd9080e250696779c927f2cb906c8bf3447df98f5ceb8dcded539/h5py-3.16.0-cp313-cp313-musllinux_1_2_x86_64.whl", hash = "sha256:7e420b539fb6023a259a1b14d4c9f6df8cf50d7268f48e161169987a57b737ff", size = 5414598, upload-time = "2026-03-06T13:48:29.49Z" }, + { url = "https://files.pythonhosted.org/packages/c3/d9/a27997f84341fc0dfcdd1fe4179b6ba6c32a7aa880fdb8c514d4dad6fba3/h5py-3.16.0-cp313-cp313-win_amd64.whl", hash = "sha256:18f2bbcd545e6991412253b98727374c356d67caa920e68dc79eab36bf5fedad", size = 3175509, upload-time = "2026-03-06T13:48:31.131Z" }, + { url = "https://files.pythonhosted.org/packages/a5/23/bb8647521d4fd770c30a76cfc6cb6a2f5495868904054e92f2394c5a78ff/h5py-3.16.0-cp313-cp313-win_arm64.whl", hash = "sha256:656f00e4d903199a1d58df06b711cf3ca632b874b4207b7dbec86185b5c8c7d4", size = 2647362, upload-time = "2026-03-06T13:48:33.411Z" }, + { url = "https://files.pythonhosted.org/packages/48/3c/7fcd9b4c9eed82e91fb15568992561019ae7a829d1f696b2c844355d95dd/h5py-3.16.0-cp314-cp314-macosx_10_15_x86_64.whl", hash = "sha256:9c9d307c0ef862d1cd5714f72ecfafe0a5d7529c44845afa8de9f46e5ba8bd65", size = 3678608, upload-time = "2026-03-06T13:48:35.183Z" }, + { url = "https://files.pythonhosted.org/packages/6a/b7/9366ed44ced9b7ef357ab48c94205280276db9d7f064aa3012a97227e966/h5py-3.16.0-cp314-cp314-macosx_11_0_arm64.whl", hash = "sha256:8c1eff849cdd53cbc73c214c30ebdb6f1bb8b64790b4b4fc36acdb5e43570210", size = 3054773, upload-time = "2026-03-06T13:48:37.139Z" }, + { url = "https://files.pythonhosted.org/packages/58/a5/4964bc0e91e86340c2bbda83420225b2f770dcf1eb8a39464871ad769436/h5py-3.16.0-cp314-cp314-manylinux_2_28_aarch64.whl", hash = "sha256:e2c04d129f180019e216ee5f9c40b78a418634091c8782e1f723a6ca3658b965", size = 5198886, upload-time = "2026-03-06T13:48:38.879Z" }, + { url = "https://files.pythonhosted.org/packages/f1/16/d905e7f53e661ce2c24686c38048d8e2b750ffc4350009d41c4e6c6c9826/h5py-3.16.0-cp314-cp314-manylinux_2_28_x86_64.whl", hash = "sha256:e4360f15875a532bc7b98196c7592ed4fc92672a57c0a621355961cafb17a6dd", size = 5404883, upload-time = "2026-03-06T13:48:41.324Z" }, + { url = "https://files.pythonhosted.org/packages/4b/f2/58f34cb74af46d39f4cd18ea20909a8514960c5a3e5b92fd06a28161e0a8/h5py-3.16.0-cp314-cp314-musllinux_1_2_aarch64.whl", hash = "sha256:3fae9197390c325e62e0a1aa977f2f62d994aa87aab182abbea85479b791197c", size = 5192039, upload-time = "2026-03-06T13:48:43.117Z" }, + { url = "https://files.pythonhosted.org/packages/ce/ca/934a39c24ce2e2db017268c08da0537c20fa0be7e1549be3e977313fc8f5/h5py-3.16.0-cp314-cp314-musllinux_1_2_x86_64.whl", hash = "sha256:43259303989ac8adacc9986695b31e35dba6fd1e297ff9c6a04b7da5542139cc", size = 5421526, upload-time = "2026-03-06T13:48:44.838Z" }, + { url = "https://files.pythonhosted.org/packages/3e/14/615a450205e1b56d16c6783f5ccd116cde05550faad70ae077c955654a75/h5py-3.16.0-cp314-cp314-win_amd64.whl", hash = "sha256:fa48993a0b799737ba7fd21e2350fa0a60701e58180fae9f2de834bc39a147ab", size = 3183263, upload-time = "2026-03-06T13:48:47.117Z" }, + { url = "https://files.pythonhosted.org/packages/7b/48/a6faef5ed632cae0c65ac6b214a6614a0b510c3183532c521bdb0055e117/h5py-3.16.0-cp314-cp314-win_arm64.whl", hash = "sha256:1897a771a7f40d05c262fc8f37376ec37873218544b70216872876c627640f63", size = 2663450, upload-time = "2026-03-06T13:48:48.707Z" }, + { url = "https://files.pythonhosted.org/packages/5d/32/0c8bb8aedb62c772cf7c1d427c7d1951477e8c2835f872bc0a13d1f85f86/h5py-3.16.0-cp314-cp314t-macosx_10_15_x86_64.whl", hash = "sha256:15922e485844f77c0b9d275396d435db3baa58292a9c2176a386e072e0cf2491", size = 3760693, upload-time = "2026-03-06T13:48:50.453Z" }, + { url = "https://files.pythonhosted.org/packages/1d/1f/fcc5977d32d6387c5c9a694afee716a5e20658ac08b3ff24fdec79fb05f2/h5py-3.16.0-cp314-cp314t-macosx_11_0_arm64.whl", hash = "sha256:df02dd29bd247f98674634dfe41f89fd7c16ba3d7de8695ec958f58404a4e618", size = 3181305, upload-time = "2026-03-06T13:48:52.221Z" }, + { url = "https://files.pythonhosted.org/packages/f5/a1/af87f64b9f986889884243643621ebbd4ac72472ba8ec8cec891ac8e2ca1/h5py-3.16.0-cp314-cp314t-manylinux_2_28_aarch64.whl", hash = "sha256:0f456f556e4e2cebeebd9d66adf8dc321770a42593494a0b6f0af54a7567b242", size = 5074061, upload-time = "2026-03-06T13:48:54.089Z" }, + { url = "https://files.pythonhosted.org/packages/cc/d0/146f5eaff3dc246a9c7f6e5e4f42bd45cc613bce16693bcd4d1f7c958bf5/h5py-3.16.0-cp314-cp314t-manylinux_2_28_x86_64.whl", hash = "sha256:3e6cb3387c756de6a9492d601553dffea3fe11b5f22b443aac708c69f3f55e16", size = 5279216, upload-time = "2026-03-06T13:48:56.75Z" }, + { url = "https://files.pythonhosted.org/packages/a1/9d/12a13424f1e604fc7df9497b73c0356fb78c2fb206abd7465ce47226e8fd/h5py-3.16.0-cp314-cp314t-musllinux_1_2_aarch64.whl", hash = "sha256:8389e13a1fd745ad2856873e8187fd10268b2d9677877bb667b41aebd771d8b7", size = 5070068, upload-time = "2026-03-06T13:48:59.169Z" }, + { url = "https://files.pythonhosted.org/packages/41/8c/bbe98f813722b4873818a8db3e15aa3e625b59278566905ac439725e8070/h5py-3.16.0-cp314-cp314t-musllinux_1_2_x86_64.whl", hash = "sha256:346df559a0f7dcb31cf8e44805319e2ab24b8957c45e7708ce503b2ec79ba725", size = 5300253, upload-time = "2026-03-06T13:49:02.033Z" }, + { url = "https://files.pythonhosted.org/packages/32/9e/87e6705b4d6890e7cecdf876e2a7d3e40654a2ae37482d79a6f1b87f7b92/h5py-3.16.0-cp314-cp314t-win_amd64.whl", hash = "sha256:4c6ab014ab704b4feaa719ae783b86522ed0bf1f82184704ed3c9e4e3228796e", size = 3381671, upload-time = "2026-03-06T13:49:04.351Z" }, + { url = "https://files.pythonhosted.org/packages/96/91/9fad90cfc5f9b2489c7c26ad897157bce82f0e9534a986a221b99760b23b/h5py-3.16.0-cp314-cp314t-win_arm64.whl", hash = "sha256:faca8fb4e4319c09d83337adc80b2ca7d5c5a343c2d6f1b6388f32cfecca13c1", size = 2740706, upload-time = "2026-03-06T13:49:06.347Z" }, +] + [[package]] name = "identify" version = "2.6.19" @@ -542,6 +585,21 @@ wheels = [ { url = "https://files.pythonhosted.org/packages/3b/d7/601b6396b33536811668935faa790112266c70661be94555999be431f86f/nibabel-5.4.2-py3-none-any.whl", hash = "sha256:553482c5f1e1034fc312edf6fb7f32236c0056439845d1c29293b7e8c98d4854", size = 3300985, upload-time = "2026-03-11T13:31:50.028Z" }, ] +[[package]] +name = "nitransforms" +version = "25.1.0" +source = { registry = "https://pypi.org/simple" } +dependencies = [ + { name = "h5py" }, + { name = "nibabel" }, + { name = "numpy" }, + { name = "scipy" }, +] +sdist = { url = "https://files.pythonhosted.org/packages/be/dc/770284133765ea963575064ea836e6acb1e6ccf291f063daa1f4ec10a3f9/nitransforms-25.1.0.tar.gz", hash = "sha256:59cb34895fc43422e1b231fa8530e6c6802c34037dab377d9fec166e2034e1a5", size = 30093782, upload-time = "2025-09-22T15:02:05.471Z" } +wheels = [ + { url = "https://files.pythonhosted.org/packages/fa/71/e399826fda4c9f41bd3b895c7d5612d4eae76269519ea923ece9f798f79f/nitransforms-25.1.0-py3-none-any.whl", hash = "sha256:8497a12d9040efb67f3e75925ed9cc61363543b6dfb589c7916aba6e08e38771", size = 18562465, upload-time = "2025-09-22T15:02:03.178Z" }, +] + [[package]] name = "niwrap" version = "0.9.3" @@ -1071,6 +1129,7 @@ source = { editable = "." } dependencies = [ { name = "bids2table" }, { name = "nibabel" }, + { name = "nitransforms" }, { name = "niwrap" }, { name = "numpy" }, { name = "polars" }, @@ -1098,6 +1157,7 @@ docs = [ requires-dist = [ { name = "bids2table", specifier = ">=2.1.2" }, { name = "nibabel", specifier = ">=5.3.3" }, + { name = "nitransforms", specifier = ">=24.0.0" }, { name = "niwrap", specifier = ">=0.9.3" }, { name = "numpy", specifier = ">=2.4.2" }, { name = "polars", specifier = ">=1.38.1" }, From 68db8a5aa0a59c7a3aaa9ab0cd626b767b2edda3 Mon Sep 17 00:00:00 2001 From: Florian Rupprecht Date: Tue, 19 May 2026 16:20:51 -0400 Subject: [PATCH 2/2] ruff format --- src/rbc/core/functional/resampling.py | 12 +++--------- 1 file changed, 3 insertions(+), 9 deletions(-) diff --git a/src/rbc/core/functional/resampling.py b/src/rbc/core/functional/resampling.py index 652f4937..c5d4c759 100644 --- a/src/rbc/core/functional/resampling.py +++ b/src/rbc/core/functional/resampling.py @@ -117,9 +117,7 @@ def _resample_4d( for t in range(n_vols): if motion_xfms is not None: vol_coords_ras = motion_xfms[t].map(static_coords_ras, inverse=True) - voxel_grid = _ras_to_voxel_grid( - vol_coords_ras, src_affine_inv, ref_shape - ) + voxel_grid = _ras_to_voxel_grid(vol_coords_ras, src_affine_inv, ref_shape) else: assert static_voxel_grid is not None # noqa: S101 voxel_grid = static_voxel_grid @@ -140,9 +138,7 @@ def _resample_4d( final_data = out_data if is_4d else out_data[..., 0] out_img = nib.Nifti1Image(final_data, reference_img.affine) if is_4d: - zooms = ( - reference_img.header.get_zooms()[:3] + src_img.header.get_zooms()[3:4] - ) + zooms = reference_img.header.get_zooms()[:3] + src_img.header.get_zooms()[3:4] else: zooms = reference_img.header.get_zooms()[:3] out_img.header.set_zooms(zooms) @@ -192,9 +188,7 @@ def _load_motion_xfms( if not motion_mats: raise FileNotFoundError(f"No motion .mat files found in {motion_mat_dir}") return [ - nt.linear.load( - str(m), fmt="fsl", reference=bold_ref_img, moving=bold_ref_img - ) + nt.linear.load(str(m), fmt="fsl", reference=bold_ref_img, moving=bold_ref_img) for m in motion_mats ]