Skip to content

Geometry, segmentation and resampling — reference

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

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 kP0_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, decided by orientations_equivalent() / spacings_equivalent() / positions_equivalent().

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 only. Never linear-interpolate a mask.
  • B-spline overshoots at sharp boundaries → clip HU afterward (TRACE-CT: [-1000, +400]).
  • 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.
  • Linear-interpolating a mask — 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)

  • Slice normalseg_decode_resample.py compute_slice_normal().
  • Affineseg_decode_resample.py build_affine(); rider_radiomics_extraction.py affine_to_sitk_geometry().
  • Geometric orderingproject_position() (IPP projected onto the slice normal).
  • Geometry equivalenceorientations_equivalent(), spacings_equivalent(), positions_equivalent(); EXACT_GRID/RESAMPLED in classify_transfer_pair().
  • FoR — carried on every CTVolume; checked in preflight, never trusted alone.
  • SEG decode + mask transferdecode_seg_file() (per-frame functional groups, source SOP UIDs), build_source_mask_on_ct_grid().
  • Inverse NN mask resampleresample_mask_inverse_nn() (label-preserving, bitwise-verified by verify_determinism()).
  • Image resample + clipresample_image_isotropic() (B-spline to 1 mm) + clip_intensities() in rider_radiomics_extraction.py.
  • Visual QCscripts/visual_qc/* and tools/slicer/open_qc_bundle.py.

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 3).

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 resample nearest-neighbour. Same dimensions ≠ same physical space.
  6. Resampling changes features. Pin and report every preprocessing choice.