Spatial Conventions¶
ConfUSIus works with three kinds of coordinate systems:
- the voxel space, linked to the underlying array storage and indexed by integer voxel coordinates,
- the physical space embedded in every DataArray's coordinates,
- and any number of reference spaces (atlas, scanner, etc.) linked to the physical
space through affine transforms stored in
attrs["affines"].
Understanding these three spaces and the axis-ordering convention used throughout ConfUSIus makes it much easier to reason about I/O, registration, and downstream statistical analysis.
---
config:
layout: elk
---
flowchart LR
V["<b>Voxel space</b><br>(integer indices)"]
P["<b>Physical space</b><br>(probe-relative)"]
W1["<b>Scanner space</b>"]
ellipsis{{"..."}}
W2["<b>Atlas space</b>"]
V -->|".coords"| P
P -->|".attrs[affines]"| W1
P -->|".attrs[affines]"| W2
P --> ellipsis
ellipsis@{ shape: text }
Axis Ordering: (time, z, y, x)¶
Every ConfUSIus DataArray that represents a fUSI recording uses the dimension order
(time, z, y, x), where:
| Dimension | Physical axis | Typical size |
|---|---|---|
time |
Acquisition time | Thousands |
z |
Elevation (stacking direction) | 1 for 2D acquisitions |
y |
Axial / depth | Tens to hundreds |
x |
Lateral | Tens to hundreds |
Dimension ordering is mostly transparent in Xarray
Users familiar with neuroimaging may be more accustomed to spatiotemporal
conventions like (x, y, z, t). Thankfully, Xarray makes dimension ordering largely
transparent in practice: you can always refer to dimensions by name and in any
order (e.g. data.mean("time"), data.sel(x=4.54, y=-2.48, z=0.0)) rather than by
axis index, so you won't have to remember the order of the dimensions.
This ordering is motivated by several considerations.
- Equivalence with NIfTI: NIfTI stores arrays in column-major (Fortran) order as
(x, y, z, time). Transposing to the more Pythonic row-major (C) order is a zero-copy operation that yields(time, z, y, x). - Memory layout for volume-wise processing: In row-major order the last axes are
contiguous in memory, so
data[t]—a single spatial volume—is a contiguous block, which is the natural unit of work for IQ processing, motion correction, and similar operations. - Statistical analysis convention: After spatial processing, fUSI data is typically
reshaped to
(time, voxels)for GLMs and dimensionality reduction. This isdata.stack(voxels=["z", "y", "x"])in Xarray, matching the standard(samples, features)convention of scikit-learn and statsmodels. - Alignment with neuroanatomical atlases: For coronal preclinical fUSI,
(z, y, x) = (elevation, axial/depth, lateral)maps to(antero-posterior, superior-inferior, left-right), sharing the first two axes with BrainGlobe atlases (e.g. Allen CCFv3). The physical → reference affine captures any remaining orientation difference (e.g. a lateral mirror). - Visualization: Most visualization tools (e.g. napari) expect the last two axes to
be the display axes of a 2D image.
data.sel(time=t, z=z)yields a(y, x)array that plots correctly without transposing.
Coordinate Systems¶
Voxel Space¶
Voxel space has its origin at voxel (0, 0, 0) and integer indices along each
spatial axis. It is the natural indexing space of the underlying array: DataArrays can
be indexed in voxel space using the standard Xarray integer-location indexer (.isel).
Physical Space¶
The physical space is the coordinate system embedded in the DataArray's dimension
coordinates. Its axes are (z, y, x) corresponding to (elevation, axial/depth,
lateral). The unit of the coordinates is determined by the units attribute of each
coordinate array and is not fixed by ConfUSIus — millimeters are typical for fUSI, but
any consistent unit can be used.
The origin is typically the center of the probe surface, but users are free to define any physical space they find convenient. What matters is that the coordinate values are internally consistent and carry a meaningful physical scale.
Physical coordinates are set at data-loading time and depend on the source format:
- EchoFrame: Lateral and axial coordinates are read from the acquisition metadata file.
- AUTC: Lateral and axial coordinates are supplied by the user as parameters to the conversion function. If coordinates are omitted, ConfUSIus falls back to bare voxel indices and emits a warning.
- Iconeus SCAN: Coordinates are derived from the
voxelsToProbeaffine embedded in the SCAN file. The axial coordinate (y) is sign-flipped so that it is always positive and increases with depth. - NIfTI: Coordinates are derived from the translation and scale components of the "best" affine transformation found in the file header.
- Hand-constructed DataArrays: The physical space is whatever the user assigns to the dimension coordinates.
The "best" NIfTI affine
NIfTI files can store two affine transforms in their header: an sform and a
qform, each with an associated integer code indicating whether the affine is
valid (code > 0) and which space it points to. ConfUSIus follows the NIfTI
specification: if sform_code > 0 the sform is used; otherwise, if
qform_code > 0 the qform is used. If both codes are zero a warning is emitted
and coordinates fall back to a diagonal affine built from the pixdim field.
Each spatial coordinate also carries a voxdim attribute that records the native voxel
size along that dimension:
This is set at load time alongside the physical coordinates and is preserved through
any Xarray operation that propagates coordinate attributes. It is particularly useful
after downsampling: the coordinate values themselves reflect the new, coarser spacing,
but voxdim retains the original acquisition resolution. It is used wherever native
voxel dimensions are needed — for example, for aspect-ratio-correct visualization,
NIfTI export, and as a fallback when a dimension is later reduced to a single point
(e.g. after .isel(z=0)):
da_down = da.isel(y=slice(None, None, 2), x=slice(None, None, 2))
# da_down.coords["y"].attrs["voxdim"] still holds the native voxel size.
da_slice = da_down.isel(z=0)
# da_slice.coords["z"].attrs["voxdim"] is now the only record of the z voxel size,
# since the coordinate has been collapsed to a scalar.
Reference Spaces¶
ConfUSIus stores affine transformations between the DataArray's physical space and any
other reference physical space in da.attrs["affines"], a dictionary keyed by affine
name. These target spaces can be external coordinate systems defined by other tools or
standards, for example an atlas space (Allen CCFv3), a scanner space, or a user-defined
stereotactic space. Each value is a (4, 4) homogeneous matrix in (z, y, x)
convention that maps a physical-space point to the corresponding point in the reference
space:
where (pz, py, px) are the physical coordinates stored in da.coords.
Several loaders populate da.attrs["affines"] automatically:
- NIfTI: NIfTI files store
qformandsformaffines in their header that map voxel indices to reference coordinates.load_niftireads the relevant affine(s), converts them from voxel → reference to physical → reference form, and stores them under the keys"physical_to_sform"and/or"physical_to_qform"depending on which codes are valid in the header. When a file declares both forms, the physical coordinates come from the primary (usuallysform) affine, andphysical_to_{secondary}(usuallyphysical_to_qform) is anchored to that same physical space, so both keys map the one set of coordinates inda.coords, just to different reference spaces. When saving back to NIfTI,save_niftican write any named affine inattrs["affines"]to the header via itsqform=andsform=arguments, with"physical_to_qform"and"physical_to_sform"used as convenience fallbacks. - Iconeus SCAN:
load_scanstores a"physical_to_lab"affine mapping ConfUSIus physical coordinates(z, y, x)to the Iconeus lab coordinate system. For multi-pose acquisitions (3Dscan,4Dscan), one affine per pose is stored, with shape(npose, 4, 4).
Registration transforms are handled separately from attrs["affines"].
register_volume returns the estimated
transform explicitly; for linear registration it follows the SimpleITK pull convention
and maps fixed physical coordinates to moving physical coordinates. That pull transform
is useful for resampling, but it is not stored automatically in attrs["affines"].
When registration is run with resample=True, the returned DataArray already lives on
the fixed/reference grid. In that case ConfUSIus copies fixed.attrs["affines"] onto
the result, so any named physical-space relationships already attached to the
fixed/reference DataArray are preserved automatically.
Why physical → reference, not voxel → reference?
The standard NIfTI affine maps voxel indices → reference coordinates. ConfUSIus uses a physical → reference affine instead, for one practical reason: it is invariant to slicing and downsampling.
A voxel → reference affine encodes the origin as the reference position of voxel
(0,0,0) specifically. The moment you crop or subsample the DataArray—even a simple
.isel(y=slice(10, 50))—voxel (0,0,0) is no longer in the array, and the affine
silently points to the wrong location.
A physical → reference affine operates on the coordinate values already stored in
da.coords. Those values travel with the data through any Xarray operation that
preserves coordinates, so the affine remains valid without any adjustment.
save_nifti reconstructs the full voxel → reference
NIfTI affine internally by composing the physical → reference affine with the
voxel → physical map built from the dimension coordinate origin and spacing. Using
the full affine (not just its orientation) means each affine keeps its own origin, so
a file with distinct sform and qform round-trips faithfully.
Switching Physical Spaces¶
The physical space is a choice of coordinate frame, not a fixed property of the data.
When a DataArray carries an affine from its current physical space (described by its
coordinates) to a reference space in attrs["affines"], you can use
.fusi.affine.apply to re-express its
physical coordinates in that reference frame.
ConfUSIus physical coordinates are stored as three independent 1D arrays, so they can
encode scaling and translation but not rotations and shears. Consequently,
.fusi.affine.apply may perform either a
complete or partial change of physical frame:
- If the affine can be represented by the coordinate arrays (that is, it is
axis-aligned),
.fusi.affine.applyabsorbs the full transform. The requested reference frame becomes the new physical frame, and the returnedorientationis the identity. - If the affine contains axis-mixing components (rotation or shear),
.fusi.affine.applyabsorbs only the axis-aligned part into the coordinates and returns the remaining transform asorientation. In that case, the change of frame is incomplete until the residual transform is handled separately.
This operation is useful for visualization, because representable changes appear directly in coordinate ticks and world axes. It can also provide a useful initialization for registration by expressing datasets in approximately the same frame. Any returned residual orientation must still be included when comparing or registering datasets.
Take the sform/qform example from Reference Spaces. A NIfTI
file with sform_code > 0 anchors its physical space to the sform frame, absorbing
only its axis-aligned part. The remaining orientation component of sform is carried on
the physical_to_sform affine attribute, while physical_to_qform records the mapping
of the physical space to the qform frame.
The following example uses an axis-aligned sform, so physical_to_sform is the
identity. The qform is shifted by +10 along z and +5 along y relative to the
sform, so physical_to_qform contains only a translation:
da.coords["z"].values # array([0., 1., 2.])
da.coords["y"].values # array([0., 1., 2., 3.])
da.attrs["affines"]["physical_to_sform"] # identity
# array([[ 1., 0., 0., 0.],
# [ 0., 1., 0., 0.],
# [ 0., 0., 1., 0.],
# [ 0., 0., 0., 1.]])
da.attrs["affines"]["physical_to_qform"] # (+10, +5) shift in (z, y):
# array([[ 1., 0., 0., 10.],
# [ 0., 1., 0., 5.],
# [ 0., 0., 1., 0.],
# [ 0., 0., 0., 1.]])
Applying physical_to_qform absorbs the translation into the coordinate arrays. The
qform frame therefore becomes the new physical frame: physical_to_qform becomes the
identity, while every other stored affine is re-expressed relative to the new physical
coordinates.
da_q.coords["z"].values # array([10., 11., 12.]) <- shifted by +10
da_q.coords["y"].values # array([5., 6., 7., 8.]) <- shifted by +5
da_q.attrs["affines"]["physical_to_sform"] # (-10, -5) shift back to sform
# array([[ 1., 0., 0., -10.],
# [ 0., 1., 0., -5.],
# [ 0., 0., 1., 0.],
# [ 0., 0., 0., 1.]])
da_q.attrs["affines"]["physical_to_qform"] # identity, the new physical frame
# array([[ 1., 0., 0., 0.],
# [ 0., 1., 0., 0.],
# [ 0., 0., 1., 0.],
# [ 0., 0., 0., 1.]])
orientation # identity, nothing left over
Because orientation is the identity, the change of physical frame is complete.
Rotations and Shears¶
Independent one-dimensional coordinate arrays cannot represent a transform that mixes
axes. For example, consider a physical_to_qform containing a 90° rotation in the
(z, y) plane:
da.attrs["affines"]["physical_to_qform"]
# array([[ 0., -1., 0., 0.],
# [ 1., 0., 0., 0.],
# [ 0., 0., 1., 0.],
# [ 0., 0., 0., 1.]])
After calling .fusi.affine.apply, the
coordinate arrays remain unchanged because none of this rotation can be represented as
independent changes to z, y, and x:
da_q, orientation = da.fusi.affine.apply(da.attrs["affines"]["physical_to_qform"])
da_q.coords["z"].values # array([0., 1., 2.])
da_q.coords["y"].values # array([0., 1., 2., 3.])
The unabsorbed transform is returned as orientation:
orientation
# array([[ 0., -1., 0., 0.],
# [ 1., 0., 0., 0.],
# [ 0., 0., 1., 0.],
# [ 0., 0., 0., 1.]])
Because the residual is non-identity, the DataArray has not been fully re-anchored to
the qform frame. The stored affines remain valid relative to the unchanged physical
coordinates:
da_q.attrs["affines"]["physical_to_sform"] # identity
da_q.attrs["affines"]["physical_to_qform"] # unchanged 90° rotation
In general, .fusi.affine.apply absorbs
the diagonal scale and translation components into the coordinate arrays and returns any
remaining axis-mixing component, such that
Always check the returned orientation. When it is non-identity, retain or compose it
with subsequent transformations—or resample the data onto the target grid—rather than
treating the coordinate arrays alone as coordinates in the target frame.