How to locate DICOM in space

I’ve worked in medical domain from time to time since 2010, especially to display DICOM images, and if I had a dollar every time I heard “This picture is upside down”, I think I could have bought a sandwich. Yes, that’s not so frequent but when it happens, it quickly becomes a nightmare. So let’s try to demystify DICOM image position.

Some Convention in 3D

I noticed that some convention that was obvious to me was not obvious to everyone, so let’s remind them:

The vast majority of 3D software always use the same colors to draw axes:

When an axis is normal to the image, a circle with a dot indicates this axis goes toward you:

And a circle with a cross indicate this axis is moving away from you:

Coordinate system

A very important thing to understand is that coordinates of a point always refer to a reference system. In the example below, the coordinates of point P are in both case (2, 1), but depending on which coordinate system we use, its final position won’t be the same:

When talking about Image Position (Patient) Attribute – (0020,0032) for example, one should understand in which space are expressed these coordinates. Same for Image Orientation (Patient) Attribute – (0020,0037). So, how to know in which space these coordinates are expressed ?

RCS – Reference Coordinates System

DICOM define two reference coordinates system depending on if the patient is a biped or a quadruped. It corresponds to the natural way of putting specimen in device (scanner, MRI, …) so that patient space correspond to device space. This information can be retrieved with the Anatomical Orientation Type Attribute – (0010,2210). As states in the doc: “If this Attribute is not present, the default human standard anatomical position is used”, that is to say Biped.

BIPED

X: from right to left
Y: from front to back
Z: from foot to head

QUADRUPED

X: from left to right
Y: from dorsal to ventral
Z: from tail to head

Image space

DICOM also define an image space with origin in center of first pixel. Note that most render engine consider origin in corner of texture, which cause a lot of “half pixel shift” issues in medical applications.

X: along rows, unit of 1mm
Y: along column, unit of 1mm
Z: cross product of X and Y to form a unit direct coordinate system.
Origin: Center of first pixel of image

When talking about 1st pixel of image, I’m not referring to top left one or bottom or whatever. I’m referring about the first value in pixel buffer. For example, with a 256*256 image:

If you do a simple DICOM image viewer, you can just locate 1st pixel in top left and that’s it. But if you do a 3D viewer and have to deal with patient orientation, oblique slices and so on, thing will becomes much more complicated, and you should better not focussing on displaying 1st pixel in top left.

Locating image in patient space

Now serious things begin. Let’s remind the 3 tags involved here:

Image Position (Patient) Attribute – (0020,0032) gives Image space origin in patient space.
Image Orientation (Patient) Attribute – (0020, 0037) gives Image space x and y axis direction in patient space. z can be retrieved with cross product between x and y.
Pixel Spacing Attribute – (0028, 0030) gives length of Image space x and y axes

We can then build a transformation matrix from image space to patient space:

Applying this transformation matrix to image properly locate it relative to patient:

Before taking care of image orientation on screen, you first have to ensure the image is properly located into patient space. Else, you will just do shotgun debugging by blindly flipping and rotating the image until that particular one is well oriented. And you will have a bug for the next one.

A few words about Anatomical planes

Before talking about image orientation on screen, we should introduce Anatomical Planes. Medical domains defines 3 cutting planes, relative to each of 3 patient main axes:

Sagittal: Normal to X axis
Coronal: Normal to Y axis
Axial (Transverse): Normal to Z axis

Sagittal Coronal Axial anatomical planes

Defining if an image is along Sagittal, Coronal or Axial plane is easy in case of axis aligned images, but in case of oblique slicing, this may be more difficult. To do so, we have to find which of the Sagittal, Coronal or Axial plane is the closer to slice plane.

Let’s first compute image normal, expressed in patient space. It is the cross product of image’s X and Y axes, returned by Image Orientation (Patient) Attribute – (0020, 0037).

Now, let’s find which of the Sagittal, Coronal or Axial axis is closest to slice’s axis. This can be done by computing dot product between slice normal and X, Y, Z patient space axes. But as we have slice normal directly expressed in patient space, dot product between X, Y, Z patient axes are directly the coordinates of slice’s normal:

Max Zimg coordinates gives anatomical plane

If the highest value is Zimg.x, it means the slice is closest to Sagittal plane, and thus can be considered as a sagittal slice.

Displaying a slice on screen

Until now, I didn’t talk about slice orientation on screen. I never talked about Up or Down. This is because these concept is a matter of Camera, not position in space.

Let’s assume you have a Coronal slice of the patient’s head. How would you want to see it on screen ? Like on the left or like on the right ?

Note that both orientations are valid. In both cases, Z axis goes from foot to head and X axis goes from patient’s left to right.

To decide which is the better, a convention defines how Sagittal, Coronal and Axial images should be displayed on screen. I didn’t find any literature affirming this convention, but it seems the most widely used among medical software.

Note that the convention is not the same for BIPED and for QUADRUPED.

Sagittal Biped

Patient Z axis along UP screen vector
Patient Y axis along RIGHT screen vector

The patient is standing up and look to the left of the screen.

Sagittal Quadruped

Patient -Y axis along UP screen vector
Patient -Z axis along RIGHT screen vector

The patient is lying down is standing up and look to the left of the screen.

Coronal (Biped and Quadruped)

Patient Z axis along UP screen vector
Patient X axis along RIGHT screen vector

The patient is standing up and face you. Its left side is on the right of the screen.

Axial (Biped and Quadruped)

Patient -Y axis along UP screen vector
Patient X axis along RIGHT screen vector

Patient is lying down and look to the top of the screen.

Let’s write some code

Let’s wrap this into some pseudo-code:

vec3 xImageInPatientSpace, yImageInPatientSpace = readtag("(0020,0037)");

vec3 sliceNormalInPatientSpace = xImageInPatientSpace.cross(yImageInPatientSpace).normalize();

Orientation orientation;
float maxCoord = max(abs(sliceNormalInPatientSpace.x),
                     abs(sliceNormalInPatientSpace.y),
                     abs(sliceNormalInPatientSpace.z));

if (maxCoord == abs(SliceNormalInPatientSpace.x))
    orientation = Sagittal;
if (maxCoord == abs(SliceNormalInPatientSpace.y))
    orientation = Coronal;
if (maxCoord == abs(SliceNormalInPatientSpace.z))
    orientation = Transverse;


if (orientation == Coronal && biped) {
  right = yPatientInWorldSpace;
  up = zPatientInWorldSpace;
} else if (orientation == Coronal && quadruped) {
  right = -zPatientInWorldSpace;
  up = -yPatientInWorldSpace;
} else if (orientation == Sagittal) {
  right = xPatientInWorldSpace;
  up = zPatientInWorldSpace;
} else if (orientation == Transverse) {
  right = xPatientInWorldSpace;
  up = -yPatientInWorldSpace;
}
camera->lookAt(target, right, up)

This will really depend on your application, which lib you use to read DICOM, which render engine you use for your display, but the main idea is here. From Image Orientation, you can compute an UP and a RIGHT vector, and then ask your camera to look at some point using the given UP and RIGHT orientation. Most render engines have this kind of feature.

Not axis-aligned image

The code above ensure having patient space always aligned with screen axes. But what if you have an image which is not aligned with patient axes, like the case below ?

Patient axes aligned on screen or image axes aligned on screen ?

We have a perfectly Sagittal acquisition. The slice normal is along X axis, but image rows are not aligned with patient Y axis. In this case, what do you want ?

Left: Patient aligned with screen. Image not aligned with screen.
Right: Image aligned with screen. Patient not aligned with screen.

Well… it depends on what you want, but you have to choose. You cannot have both image and patient aligned on screen.

The code presented previously will give you the Left result. Patient axes aligned on screen. If you want the Right one, you should add an extra step to ensure left and up vector are aligned with image X and Y vectors:

vec3 xImageInPatientSpace, yImageInPatientSpace = readtag("(0020,0037)");

vec3 sliceNormalInPatientSpace = xImageInPatientSpace.cross(yImageInPatientSpace).normalize();

Orientation orientation;
float maxCoord = max(abs(sliceNormalInPatientSpace.x),
                     abs(sliceNormalInPatientSpace.y),
                     abs(sliceNormalInPatientSpace.z));

if (maxCoord == abs(SliceNormalInPatientSpace.x))
    orientation = Sagittal;
if (maxCoord == abs(SliceNormalInPatientSpace.y))
    orientation = Coronal;
if (maxCoord == abs(SliceNormalInPatientSpace.z))
    orientation = Transverse;


if (orientation == Coronal && biped) {
  right = yPatientInWorldSpace;
  up = zPatientInWorldSpace;
} else if (orientation == Coronal && quadruped) {
  right = -zPatientInWorldSpace;
  up = -yPatientInWorldSpace;
} else if (orientation == Sagittal) {
  right = xPatientInWorldSpace;
  up = zPatientInWorldSpace;
} else if (orientation == Transverse) {
  right = xPatientInWorldSpace;
  up = -yPatientInWorldSpace;
}

// Extra step to handle Image screen aligned {
right = findClosestVector(right, [xImageInPatientSpace, yImageInPatientSpace])
up = findClosestVector(up, [xImageInPatientSpace, yImageInPatientSpace])
// }

camera->lookAt(target, right, up)

With findClosestVector:

// Find among vectors, the vector which is the closest to refVector.
// All vectors must be normalized.
vec3 findClosest(vec3 refVector, vec3[] vectors) {
  float closestDot = 0;
  vec3 closestVector = vec3(0);
  for(v in vectors)
  {
    float dot = abs(refVector.dot(v));
    if(dot > closestDot)
      closestVector = v;
  }

  if(closestVector.dot(refVector) > 0)
    return closestVector;
  else
    return -closestVector;
}

Conclusion

I think That video perfectly summarize the kind of argument you may have when your user require “having images always aligned on screen and patient looking to the left”.