Getting slice index from background layer

Support
1

Hi,
We want to use this code to get the slice number/index (found it on script repository but now the code is gone and we are unable to refer back to it as script repository documentation is updated) and wanted to confirm few things:

if (d->sliceLogic != nullptr) {
		vtkMRMLSliceLayerLogic* BL = d->sliceLogic->GetBackgroundLayer();
		vtkGeneralTransform* xyToIJK = BL->GetXYToIJKTransform();
		const float inCoord[3] = { 0,0,0 };
		float outCoord[3] = { 0,0,0 };
		xyToIJK->InternalTransformPoint(inCoord, outCoord);
		k = outCoord[c];
	}
	return k;
  1. What is the inCoord for, as per my understanding, is it defining the origin of 3D volume?
  2. Does the k, which is the z-axis always corresponds to Axial view?
  3. We are storing ‘k’ as int type but the coordinate is float type. Should we round the outCoord to nearest integer before storing into ‘k’? because the way it is right now, it will always floor.

Regards
SI

1 Like
2

DICOM does not allow any kind of indexing based on slice position, so no geometric conversion exists between physical coordinates to “slice index”.

The number that you see in 2D DICOM viewers when you browse the images is the Instance number, which is usually assigned to image slices based on where the image slice is located and/or when it was acquired. It is usually offset by one compared to the third voxel coordinate value (as instance numbers are typically 1-based, while voxel coordinates are 0-based) and often also inverted (as voxel coordinate increases, instance number can increase or decrease), and sometimes (e.g., when not all originally frames are available) then the offset may change throughout the series. Since a slice view is not always aligned with the original image slices (reformatted in orthogonal or oblique orientations), the instance number also often changes within a slice (different regions of the displayed image is reconstructed from different original slices).

Therefore, the only reliable way to get the instance number for a physical point position is the following:

Implementation in Python (after the IJK coordinates are computed):

ijk = [12, 34, 56]
volumeNode = slicer.mrmlScene.GetFirstNodeByClass('vtkMRMLVolumeNode')
instanceUID = volumeNode.GetAttribute('DICOM.instanceUIDs').split(' ')[ijk[2]]
instanceNumber = slicer.dicomDatabase.instanceValue(instanceUID, '0020,0013')
print(f"Instance Number: {instanceNumber}")

You can add this code snippet to DataProbe.py to make the instance number show up in the Data Probe in the bottom-left of the screen. If you implement it and works well for you then send a pull request, we may consider adding it to the Slicer core, as it can be useful for cross-referencing positions with slices shown in other DICOM viewers.

2 Likes
3

Hi Andras,

Thanks for the detailed information. I have a follow up question:

For scans that we know for sure were acquired as a series of axially scanned slices each stored as DICOM file, is the code I posted valid, for determining the slice index? Basically, the math corresponding to the code I posted seems to be as follows:

20210816_191803
20210816_1918031920×727 293 KB

so it seems like the a12 element always encodes the “slice index” is that correct?

1 Like
4

As I wrote above, you cannot determine the DICOM instance number from the slice number (a12) because there are several unknown factors, out of your control that influences this mapping. Even if you find that for some images, slice index + 1 = DICOM instance number, it does not mean that will be true for other images as well. It is enough that the imaging technologist slightly adjusts an acquisition parameter and you get a slightly or completely different mapping.

Fortunately, there is a reliable solution that I described above.

1 Like
5

Hi Andras, thanks got it about the slice number/index aka K NOT having a simple relation with the InstanceNumber, and hence we have to use K to query the DICOM database tor retrieve the InstanceNumber. But, our main question is, can we ALWAYS trust a12 to be the same as the slice number/index aka K? , where a12 is the (3, 4)th element of the transformation matrix obtained via:
xyToIJK = BL->GetXYToIJKTransform()

Would greatly appreciate your confirmation of this. Thanks.

1 Like
6

Yes, when you load an image using DICOM module using DICOMScalarVolumePlugin then third voxel coordinate value (K) is always the same as the index in the instanceUIDs list (because the same list is used to populate the instanceUIDs list that is passed to the DICOM reader).

2 Likes
7

Thank you very much @lassoan for your support, closing this thread.

1 Like