Geometry, segmentation and resampling — reference

Reference · Dense lookup material — return here later. The teaching versions live in the chapter.

Coordinate systems

System	What it is	Order / unit
Array indices	voxel position in the NumPy array	[k, j, i] = slice, row, column (integers)
Physical patient	position in the patient	(x, y, z) LPS mm (continuous)
Display	screen pixels	view-dependent

LPS: +x → Left, +y → Posterior, +z → Superior. 3D Slicer uses RAS (sign flip on x, y). NumPy has no coordinate system at all.

The index ↔ patient equation (Image Plane Module, PS3.3 C.7.6.2.1.1)

DICOM pixel indices: i = column, j = row.

• ImageOrientationPatient (0020,0037) = [rx,ry,rz, cx,cy,cz]: first triplet = row direction (increasing i); second = column direction (increasing j).
• PixelSpacing (0028,0030) = [row spacing, col spacing].

In-plane:

P(i, j) = P0 + i·(colSpacing)·[rx,ry,rz] + j·(rowSpacing)·[cx,cy,cz]

P0 = ImagePositionPatient (0020,0032) of the top-left pixel. Slice normal n = cross([rx,ry,rz],[cx,cy,cz]); slice step = slice spacing; slice index k → P0_k = ImagePositionPatient of slice k.

Index orders (the trap)

DICOM pixel       (i, j)     = (column i, row j)
SimpleITK index   (x, y, z)  = (column i, row j, slice k)   # NOT NumPy order
NumPy array       [k, j, i]  = (slice k, row j, column i)

img.GetSize() is (x,y,z); arr.shape is (z,y,x). Always print GetSpacing()/GetOrigin()/GetDirection(); convert at the boundary with GetArrayFromImage / GetImageFromArray.

Slice ordering

t_k = dot(ImagePositionPatient_k, slice_normal); sort slices by t_k ascending. Spacing truth = t_{k+1} − t_k; SpacingBetweenSlices is only a claim. InstanceNumber and filenames are not geometry.

Frame of Reference decision

same FoR      ≠  same voxel grid       (verify origin/spacing/direction)
different FoR ≠  impossible to relate  (explicit registration can map them)
different FoR ⇒  do NOT overlay without an explicit spatial relationship

Pairing classification: EXACT_GRID / RESAMPLING_REQUIRED / INCOMPATIBLE.

Segmentation representations

Representation	Stores	Carries	Loses on conversion
Binary mask	voxel labels	image geometry (if kept with image)	source references; fractional edges
DICOM SEG	per-frame labelled planes	source SOP refs; per-frame geometry; segment identity	fractional→binary; off-grid snapping
RTSTRUCT	coplanar contours	contour points + reference grid	rasterisation is lossy/tool-dependent

SEG Segmentation Type (0062,0001): BINARY / FRACTIONAL / LABELMAP.

Resampling

• Image: linear or B-spline (radiomics prefers B-spline, then clip).
• Mask: nearest-neighbour is the label-preserving/default-safe path. For an explicitly binary 0/1 ROI, IBSI also permits linear/trilinear interpolation to a common target physical geometry when declared thresholding/re-binarisation is part of the contract. Arbitrary multiclass numeric IDs are not continuous values.
• B-spline evaluates a smooth piecewise polynomial fitted through the neighbourhood → can overshoot past the input range at sharp boundaries → clip HU afterward (TRACE-CT: [-1000, +400]).
• Operation order is feature-determining: TRACE-CT freezes HU → 1 mm external resample (B-spline image / NN mask) → clip → PyRadiomics. resample → clip ≠ clip → resample (interpolation sees different inputs); the order is part of the measurement definition, not a style choice.
• No double resampling: the PyRadiomics parameter file must not set resampledPixelSpacing when resampling is done externally — two interpolation passes are a different measurement than one.
• Clipping changes intensities; resegmentation (resegmentRange) removes voxels from the mask. TRACE-CT: external clip, resegmentRange = None.
• Upsampling creates samples, not anatomy (5 mm → 1 mm grid ≠ recovered resolution).
• Report the target spacing/interpolator with every result.

resampled_img  = sitk.Resample(img,  ref, sitk.Transform(), sitk.sitkBSpline, 0.0, sitk.sitkInt16)
resampled_mask = sitk.Resample(mask, ref, sitk.Transform(), sitk.sitkNearestNeighbor, 0, sitk.sitkUInt8)

What usually goes wrong

• Index vs physical mix-up — using array indices where physical mm are needed.
• LPS vs RAS confusion (DICOM/SimpleITK LPS; Slicer RAS; NumPy neither).
• Ordering by InstanceNumber instead of slice-normal projection.
• Trusting FoR equality without verifying origin/spacing/direction.
• Leaving an explicitly binary mask fractional after linear interpolation without declared thresholding/re-binarisation — destroys label integrity.
• Comparing array shapes instead of physical geometry.
• B-spline overshoot at sharp boundaries left unclipped.
• Resampling to a target spacing but not reporting it.

TRACE-CT connection (read-only references)

The concrete functions behind this chapter — slice normal, affine, geometric ordering, geometry equivalence, FoR handling, SEG decode and mask transfer, inverse-NN mask resampling, image resampling with clip, and visual-QC tooling — are inventoried in the TRACE-CT learning map.

Go deeper

• SimpleITK fundamentals docs (origin/spacing/direction, physical vs index space, resampling) — the single best reference for the affine and resampling.
• DICOM PS3.3 Image Plane Module + Multi-frame Functional Group macros — normative IOP/IPP and SEG per-frame geometry.
• dcmqi / highdicom docs — SEG encoding/decoding beyond pydicom.
• IBSI image-processing workflow (Zwanenburg et al. 2020) — which resampling/discretisation choices a radiomics pipeline must pin (Chapter 5).

What to retain (chapter summary)

1. Three coordinate systems: array indices, physical patient (LPS), display. The affine connects the first two.
2. Geometry = origin + spacing + direction. Slice normal = cross(row, col); order slices by projection onto it — never by InstanceNumber.
3. FrameOfReferenceUID is a fast rejection / slow acceptance test; always verify geometry.
4. Decode SEGs by per-frame functional groups and source SOP UIDs; round-trips lose fractional edges and off-grid frames.
5. Images interpolate (linear/B-spline); masks use nearest-neighbour as the direct label-preserving/default-safe path. For an explicitly binary 0/1 ROI, IBSI also permits linear/trilinear interpolation to a common target physical geometry when declared thresholding/re-binarisation follows; arbitrary multiclass numeric IDs are not continuous values. Same dimensions ≠ same physical space.
6. Resampling changes features. Pin and report every preprocessing choice.