confusius.multipose¶
multipose ¶
Multi-pose data processing utilities.
This module provides functions for processing multi-pose fUSI data, including assembling independently loaded poses into one array, consolidating multiple poses into a single volume, slice timing correction, and other multi-pose specific operations.
Modules:
-
consolidate–Multi-pose volume consolidation.
-
slice_timing–Slice timing correction for multi-pose fUSI data.
-
stack–Assembling single-pose VoxelData arrays into one pose-dependent DataArray.
-
timing–Timing helpers shared across the multipose module and other pose-aware consumers.
Functions:
-
consolidate_poses–Merge the
posedimension into the swept voxel dimension, ordered by position. -
correct_slice_timings–Resample each sweep position to the volume's reference time.
-
stack_poses–Stack independently loaded single-pose VoxelData arrays into one pose-dependent array.
consolidate_poses ¶
Merge the pose dimension into the swept voxel dimension, ordered by position.
Per-(pose, sweep_dim) world positions are read directly from da's own world
coordinates, which requires da's primary voxel-to-world geometry to itself be
pose-dependent (a (npose, 4, 4) affine stack — see
[VoxelToWorldIndex.is_pose_dependent][confusius._utils.geometry.VoxelToWorldIndex.is_pose_dependent]).
This is the case for SCAN data, where lab space (a fixed scanner frame shared by
every pose) is the canonical world frame — see
load_scan — and for
stack_poses output. If da's primary
geometry is instead driven by a secondary, named affine in
da.attrs["affines"] (e.g. a stack of NIfTI DataArrays with their
world_to_qform affines stacked, no pose-dependent primary geometry of their
own), rebase onto it first with
.fusi.affine.apply, e.g.
da.fusi.affine.apply("world_to_qform"), before calling this function.
The swept voxel dimension is detected from da's own geometry: the per-pose
translation's dominant direction (via SVD) is matched against each voxel
dimension's world-space direction, and the best-aligned dimension is used. The
primary sweep direction is then found via singular value decomposition of all
positions along that dimension. Each voxel is projected onto that axis, the
positions are checked for regularity, then the data is reindexed in ascending
order along the consolidated sweep axis.
This function is primarily intended for consolidating multi-pose fUSI volumes
acquired with an Iconeus system using a purely translational probe sweep. In that
workflow, each pose corresponds to one probe position along the elevation axis
(k/world z), and the VoxelData array is produced by
load_scan:
scan_3d = load_scan("recording.scan") # dims: (pose, k, j, i)
volume = consolidate_poses(scan_3d) # dims: (k, j, i)
scan_4d = load_scan("recording_4d.scan") # dims: (time, pose, k, j, i)
volume = consolidate_poses(scan_4d) # dims: (time, k, j, i)
Parameters:
-
(da¶DataArray) –VoxelData array with a
posedimension and pose-dependent primary voxel-to-world geometry. Typically produced byload_scanfor3Dscanor4Dscanfiles, or bystack_poses. -
(rtol¶float, default:0.01) –Relative tolerance for the regularity check (fraction of mean spacing).
Returns:
-
DataArray–VoxelData array with
posemerged into the swept voxel dimension, sorted by world position. Every voxel keeps its world position: the output voxel-to-world affine carries over the input's non-swept columns (including any rotation), uses one regular step along the detected sweep axis for the swept column, and is anchored at the first sorted voxel. For inputs whosetimecoordinate is itself pose-dependent ((time, pose)-shaped -- see stack_poses), a consolidatedslice_timewith dims("time", <sweep_dim>)is included: each slice inherits the timestamp of the pose it came from.
Raises:
-
ValueError–If
dahas noposedimension, ifda's primary geometry is not pose-dependent, if the rotation block of the affine is not constant across poses (non-translation sweep), if the swept voxel dimension cannot be detected (fewer than two poses, identical pose positions, or a degenerate voxel axis), or if the consolidated positions are not regularly spaced withinrtol-- which also rejects a sweep that steps along more than one voxel axis at once, since no single voxel dimension can span it.
Warns:
-
UserWarning–If the sweep is not purely 1D (secondary/primary singular value ratio > 0.01).
correct_slice_timings ¶
correct_slice_timings(
da: DataArray,
method: Literal[
"linear",
"nearest",
"nearest-up",
"zero",
"slinear",
"quadratic",
"cubic",
"previous",
"next",
] = "linear",
fill_value: float
| tuple[float, float]
| Literal["extrapolate", "nan"] = "extrapolate",
) -> DataArray
Resample each sweep position to the volume's reference time.
In multi-pose fUSI acquisitions, each sweep position is acquired at a different time
within the volume period. This function resamples each position's time series so
that all positions appear to have been acquired at the time stored in the time
coordinate.
This function works on both:
- Consolidated data: dims
(time, <sweep_dim>, ...)with aslice_timecoordinate with dims(time, <sweep_dim>), typically produced byconsolidate_poses. - Unconsolidated data: dims
(time, pose, ...)with a pose-dependent(time, pose)-shapedtimecoordinate (see stack_poses), holding each pose's own real timestamp directly. The result'stimebecomes a genuine 1D coordinate, computed the same wayconsolidate_posesderives its whole-arraytimefrom per-pose timestamps -- after correction, every pose really is simultaneous, so there is no more reason fortimeto stay pose-dependent.
The sweep dimension is inferred from the second dim of whichever timing coordinate is present.
If the input is Dask-backed, the function stays lazy: computation is deferred until
.compute() is called. The time dimension must not be chunked; spatial dimensions
may be freely chunked.
Parameters:
-
(da¶DataArray) –VoxelData array with a
slice_timecoordinate, or a pose-dependent(time, pose)-shapedtimecoordinate, with dims(time, <sweep_dim>). -
(method¶(linear, nearest, nearest - up, zero, slinear, quadratic, cubic, previous, next), default:"linear") –Interpolation method passed to
scipy.interpolate.interp1d:"linear": linear interpolation."nearest": nearest-neighbour interpolation; rounds down at half-integers."nearest-up": nearest-neighbour interpolation; rounds up at half-integers."zero": zeroth-order spline (step function)."slinear": first-order spline."quadratic": second-order spline."cubic": third-order spline."previous": use previous point's value."next": use next point's value.
-
(fill_value¶float or tuple[float, float] or {extrapolate, nan}, default:"extrapolate") –How to handle target times that fall outside the range of a position's acquisition times.
"extrapolate"allows linear extrapolation."nan"inserts NaNs out of bounds. Use a float for a constant fill value, or a tuple(left, right)for different values on each side.
Returns:
-
DataArray–New VoxelData array with the same dims as the input, resampled so every sweep position appears simultaneous. For already-consolidated input,
timeis unchanged andslice_timeis dropped (avoiding accidental double-correction). For pose-dependent input,timebecomes a genuine 1D coordinate (see stack_poses for what it was before correction).
Raises:
-
ValueError–If
dahas notimedimension or only one time point, ifdahas neither aslice_timecoordinate nor a pose-dependenttimecoordinate, if the timing coordinate does not have dims(time, <sweep_dim>), or if thetimedimension is chunked.
Warns:
-
UserWarning–If a spline method fails due to too few points and falls back to
"linear".
stack_poses ¶
Stack independently loaded single-pose VoxelData arrays into one pose-dependent array.
xr.concat cannot combine N single-grid VoxelToWorldIndex objects into one
joint pose-dependent index by itself: xarray's own pre-concat alignment step only
excludes a coordinate from its equality check when that coordinate's existing
index already spans the concat dimension, which a single-grid array's z/y/x
index does not (it only spans k/j/i). This function closes that gap by first
promoting each input to a genuinely pose-dependent index of length 1 (see
[VoxelToWorldIndex][confusius._utils.geometry.VoxelToWorldIndex]) via
[attach_voxel_to_world_index][confusius._utils.geometry.attach_voxel_to_world_index]
-- once every input's z/y/x index already spans pose, alignment correctly
excludes it and dispatches to VoxelToWorldIndex.concat, which merges the pose
labels and affine stacks in order exactly as
xr.concat(..., dim="pose") expects.
Parameters:
-
(poses¶sequence[DataArray]) –VoxelData arrays to stack, one per pose, in pose order. Each must have no existing
posedimension (see ensure_voxeldata'sallow_poseparameter). Voxel dimensions, shape, voxel-space (k/j/i) coordinate values, and any non-core dimensions must otherwise agree across poses, exactly as required to merge non-concatenated variables in anyxr.concatcall. -
(pose¶sequence[hashable], default:None) –Pose coordinate labels, one per entry of
poses. If not provided, defaults to0, 1, ..., len(poses) - 1.
Returns:
-
DataArray–Stacked DataArray with a new
posedimension and pose-dependent voxel-to-world geometry (an(npose, 4, 4)affine stack, one affine per input). If atimedimension is present and every pose shares identicaltimevalues,timestays an ordinary 1D dimension coordinate. If per-posetimevalues differ (poses acquired sequentially rather than simultaneously),timeinstead becomes a genuine(time, pose)-shaped coordinate holding each pose's own real timestamp directly -- there is no single answer for "the" time of a(pose, k, j, i)voxel any more than there is a single answer for itsz/y/xposition, sotimerequires a scalarposeselection first, exactly like world coordinates already do. A 2Dtimeis not itself an index (xarray dimension coordinates must be 1D), so.sel(time=...)is unavailable until a pose is selected; after that,.set_xindex("time")promotes the resulting 1Dtimeback into a real, selectable index.
Raises:
-
ValueError–If
posesis empty, ifposedoes not have one label per entry ofposes, if any pose already has aposedimension, or if poses have mismatchedtimelengths.