Skip to content

The DICOM row/column trap

Tonight · ~25 min · read · energy: medium · setup: none

You wrote arr[x, y] to index column x, row y of a CT — and the anatomy came out transposed. Or you read PixelSpacing, assumed the first number was the x spacing, and your resampling stretched the image the wrong way. Both are the same bug: DICOM’s row/column naming crosses itself, and it crosses again against NumPy. This lesson is the one place to nail down the conventions so you never transpose by accident again.

The anchor: the same voxel, three naming systems

One physical voxel is described by three different tools using three different index orders. They all mean the same voxel, but they index it differently:

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)

NumPy and SimpleITK index the same voxel in opposite orders. That is the bite.

The cross-coupling, spelled out

The confusion has two layers. Layer 1 — the names cross:

  • ImageOrientationPatient first triplet = row direction = direction of increasing column index i.
  • PixelSpacing first value = row spacing = step when row index j increases by 1.

So the row direction cosine pairs with the column spacing (it is the step along i), and the column direction cosine pairs with the row spacing (the step along j). The word “row” appears in both but refers to different axes. This is not pedantry — getting it backwards resamples or registers along the wrong axis.

Layer 2 — the index orders cross between libraries:

A DICOM pixel at (i = column, j = row) lives in a NumPy array (slices-first, image order) as numpy[k, j, i], i.e. numpy[slice, row, column]. So arr[0] is slice 0, and within a slice arr[k, j, i] indexes row j, column i — the transpose of the naive arr[x, y].

flowchart TD
    D["DICOM: pixel (i=col, j=row)"]
    N["NumPy: arr[k, j, i]<br/>= slice, row, column"]
    S["SimpleITK: GetPixel(x=i, y=j, z=k)<br/>= column, row, slice"]
    D -- "same voxel" --> N
    D -- "same voxel" --> S
    N -. "opposite index order" .-> S

The one-line rule

When in doubt, hold this in your head:

DICOM i is column, j is row; PixelSpacing is [row, col]; NumPy is [k, j, i]; SimpleITK is (x=i, y=j, z=k).

And the operational version: always print GetSpacing()/GetOrigin()/GetDirection() rather than assuming, and convert explicitly with sitk.GetArrayFromImage / GetImageFromArray at the NumPy↔SimpleITK boundary.

A concrete check you can run

import SimpleITK as sitk
img = sitk.ReadImage("series_dir") # or build from DICOM
arr = sitk.GetArrayFromImage(img) # -> NumPy [k, j, i]
print(img.GetSize()) # (x=cols, y=rows, z=slices) -- SimpleITK order
print(arr.shape) # (slices, rows, cols) -- NumPy order (reversed)
print(img.GetSpacing()) # (sx, sy, sz) in mm, LPS

Note the reversal: img.GetSize() is (x, y, z) while arr.shape is (z, y, x). The voxel at SimpleITK (x=10, y=20, z=5) is arr[5, 20, 10]. If you forget the reversal you have silently transposed the volume.

Recognising orientations from the cosines

The two IOP triplets tell you which way the row/column axes point in the patient. For a standard axial series the row points left→right and the column points anterior→posterior, giving the axial IOP we used in lesson 2:

ImageOrientationPatient = [1, 0, 0, 0, 1, 0] # row=+x(L), col=+y(P), axial plane

A coronal or sagittal series has different cosines; an oblique series (acquired at a tilt) has non-axis cosines. Never assume “axial” — read the cosines, because the direction matrix is part of the geometry that decides whether two objects share space (lesson 5).

Stop and think — then reveal

You read PixelSpacing = [0.8, 0.7] and need the step in millimetres when the column index i increases by 1. Which value do you use — 0.8 or 0.7?

0.7 — the second value (column spacing). The column index i increases along a row, so the step is the column spacing, which is the second element of PixelSpacing. The first value (0.8, row spacing) is the step when the row index j increases. Mixing these up stretches or compresses your resample along the wrong axis — a silent error that looks fine until you overlay a mask.

What to retain

  1. The same voxel has three index orders: DICOM (i=col, j=row), SimpleITK (x=i, y=j, z=k), NumPy [k, j, i]. NumPy and SimpleITK are reversed.
  2. PixelSpacing is [row spacing, col spacing]; the row direction cosine pairs with the column spacing. The naming crosses itself.
  3. img.GetSize() is (x,y,z); arr.shape is (z,y,x). Convert at the boundary and print the geometry rather than assuming.
  4. Never assume “axial” — read ImageOrientationPatient; oblique volumes have non-identity direction matrices.

Next: with index/physical space sorted, the question that makes or breaks a 3D volume — how to order slices without trusting filenames.