Hi all,
What is the most accurate method to take common DICOM tags and calculate voxel resolution? For instance, how do FOV and percent sampling impact the calculation of pixel spacing?
Many thanks!
Maggie
You can determine in plane resolution from the DICOM tag 0028,0030 Pixel Spacing Attribute, which is specified by a numeric pair - adjacent row spacing (delimiter) adjacent column spacing in mm (note the order, in NIfTI the order of the three spatial dimensions are column [i], row [j] and slice [k]).
In theory, one can determine distance between voxel centers in the slice direction using DICOM tag 0018,0088 Spacing Between Slices - which should be the sum of the SliceThickness (0018,0050) and the slice gap (slice gaps are common in MR to avoid interference between slices). Unfortunately, not all vendors used 0018,0088 consistently and therefore the robust method to determine distance between voxel centers in the slice direction is to calculate the distance between the 0020,0032 ImagePositionPatient of two neighboring slices. While this value is virtually always the same for all slices in a MRI sequence, in CT it is common to have variable slice thickness (e.g. for the head, more slices in the critical-for-life brain stem than the cortex).
Here is an example for the critical tags from a typical MRI scan:
(0028,0030) DS [1.75\1.75] # 10, 2 PixelSpacing
(0018,0050) DS [2.5] # 4, 1 SliceThickness
(0018,0088) DS [2.5] # 4, 1 SpacingBetweenSlices
(0020,0032) DS [-114.58726620674-104.09609556198-74.025070190429] # 50, 3 ImagePositionPatient
1 Like
Thank you for the information. Very helpful!
Follow up, how does percent sampling impact voxel size? Or are the values from the pixel spacing dicom tag equivalent to voxel size (non-slice dimensions)?
I have a bit of code based on an equation for calculating voxel resolution that I can no longer find and I’m looking to determine the accuracy of the approach.
fov_frequency = round(dicom_header.PixelSpacing[0] * dicom_header.Columns)
# FOV_phase = FOV * percent phase FOV
# Percent phase FOV (also known as RFOV) is the ratio of FOV dimension in phase
# direction to the FOV dimension in the frequency direction.
rfov = dicom_header.PercentPhaseFieldOfView / 100
fov_phase = round(dicom_header.PixelSpacing[1] * dicom_header.Rows * rfov)
# Percent sampling is the fraction of acquisition matrix lines acquired.
percent_sampling = dicom_header.PercentSampling / 100
if (
dicom_header.InPlanePhaseEncodingDirection == "ROW"
or dicom_header.InPlanePhaseEncodingDirection == "OTHER"
):
acquisition_matrix_frequency = dicom_header.AcquisitionMatrix[1]
acquisition_matrix_phase = round(
dicom_header.AcquisitionMatrix[2] * rfov * percent_sampling
)
elif dicom_header.InPlanePhaseEncodingDirection == "COL":
acquisition_matrix_frequency = dicom_header.AcquisitionMatrix[0]
acquisition_matrix_phase = round(
dicom_header.AcquisitionMatrix[3] * rfov * percent_sampling
)
pixel_size_frequency = fov_frequency / acquisition_matrix_frequency
pixel_size_phase = fov_phase / acquisition_matrix_phase
voxel_x = float("{:03.3f}".format(pixel_size_frequency))
voxel_y = float("{:03.3f}".format(pixel_size_phase))
voxel_z = float("{:03.3f}".format(dicom_header.SliceThickness))
Do you have any suggestions? Or stick to the pixel spacing for reporting voxel resolution?
Many thanks for the advice!
Maggie
@Maggie
I would suggest you test your algorithm on some DICOM images. Here are the results for DICOMs where interpolation, partial FOV and mosaics are acquired.
I would suggest the method to use would be:
pixdim_1 = dicom_header.PixelSpacing[1]
pixdim_2 = dicom_header.PixelSpacing[0]
pixdim_3 = dicom_header.SliceThickness
as long as you can assume that the vendor is accurately reporting the distance between voxel centers in SliceThickness (compare distance between the ImagePositionPatient of two adjacent slices for a more robust solution).
Your example does not appear robust in my sample data (and a sample dataset with rectangular voxels will reveal that you need to swap the order of pixel spacing):
>./pixdim.py PF_PAT_POS_BW_INTERP_TEST.MR.HEAD_HARMS.0029.0001.2017.09.20.16.16.44.39766.120172126.IMA @margaretmahan pixdim 19.2 19.2 2.4 @Chris_Rorden pixdim 1.2000000476837 1.2000000476837 2.4000000953674 >./pixdim.py PF_PAT_POS_BW_INTERP_TEST.MR.HEAD_HARMS.0033.0001.2017.09.20.16.16.44.39766.120190570.IMA @margaretmahan pixdim 19.2 19.2 2.4 @Chris_Rorden pixdim 2.4000000953674 2.4000000953674 2.4000000953674 >./pixdim.py PF_PAT_POS_BW_INTERP_TEST.MR.HEAD_HARMS.0035.0001.2017.09.20.16.16.44.39766.120198766.IMA @margaretmahan pixdim 19.2 19.2 2.4 @Chris_Rorden pixdim 2.4000000953674 2.4000000953674 2.4000000953674
Here is the code to replicate this:
#!/usr/bin/env python3
# -*- coding: utf-8 -*-
# python3 d_speed_size.py : test speed/compression for folder 'corpus'
# python3 d_speed_size.py indir : test speed/compression for folder 'indir'
import os
import sys
import pydicom as dcm
if __name__ == '__main__':
"""Compare speed and size for different compression tools
Parameters
----------
fnm : str
name of DICOM file to inspect
"""
if len(sys.argv) > 1:
fnm = sys.argv[1]
else:
print('No filename provided')
sys.exit()
if not os.path.isfile(fnm):
print('Unable to find "' + fnm +'"')
sys.exit()
dicom_header = dcm.read_file(fnm)
# FOV_phase = FOV * percent phase FOV
# Percent phase FOV (also known as RFOV) is the ratio of FOV dimension in phase
# direction to the FOV dimension in the frequency direction.
rfov = dicom_header.PercentPhaseFieldOfView / 100
fov_phase = round(dicom_header.PixelSpacing[1] * dicom_header.Rows * rfov)
fov_frequency = round(dicom_header.PixelSpacing[0] * dicom_header.Columns)
# Percent sampling is the fraction of acquisition matrix lines acquired.
percent_sampling = dicom_header.PercentSampling / 100
if (
dicom_header.InPlanePhaseEncodingDirection == "ROW"
or dicom_header.InPlanePhaseEncodingDirection == "OTHER"
):
acquisition_matrix_frequency = dicom_header.AcquisitionMatrix[1]
acquisition_matrix_phase = round(dicom_header.AcquisitionMatrix[2] * rfov * percent_sampling)
elif dicom_header.InPlanePhaseEncodingDirection == "COL":
acquisition_matrix_frequency = dicom_header.AcquisitionMatrix[0]
acquisition_matrix_phase = round(
dicom_header.AcquisitionMatrix[3] * rfov * percent_sampling
)
pixel_size_frequency = fov_frequency / acquisition_matrix_frequency
pixel_size_phase = fov_phase / acquisition_matrix_phase
voxel_x = float("{:03.3f}".format(pixel_size_frequency))
voxel_y = float("{:03.3f}".format(pixel_size_phase))
voxel_z = float("{:03.3f}".format(dicom_header.SliceThickness))
print('@margaretmahan pixdim {0} {1} {2}'.format(voxel_x, voxel_y, voxel_z))
pixdim_1 = dicom_header.PixelSpacing[1]
pixdim_2 = dicom_header.PixelSpacing[0]
pixdim_3 = dicom_header.SliceThickness
print('@Chris_Rorden pixdim {0} {1} {2}'.format(pixdim_1, pixdim_2, pixdim_3))@margaretmahan
Hello,
i just see your answers related to dicom resolution.
I have a similar question.
actually i have a CT dicom slices images, i want to calibrate the 3d image as per the x-ray image dpi.
For that i need to find the resolution of each slice or whole 3d image. I have no idea.
can ypu suggest what to do and how i can find the resolution using head information of a dicom slices.
@Payal_Maken the solution is the same for CT and MRI scans. You will want to be aware of two issues.
- CT scans can have gantry tilt. This shear does not influent the volume, but you need to be aware that the voxels might be rhomboidal.
- MRI scans virtually always have even spacing between slices. In contrast, some CT sequences choose to have thin slices near the life-critical brain stem and thicker slices for the cortex. In this case the volume of an individual voxel varies across the image. Formats like NIfTI require that all slices are equidistant, so in these cases tools like dcm2niix will create an image with interpolated equidistant slices (these volumes have the
_Eqpost fix in their file name).
Thank you @neurolabusc
for the issues i have to do the image calibration.
For that i will be making the dpi or resolution of the 3d image equal to the x-ray image dpi which i am usng for 3d reconstruction.
i am stuck here in which i have to extract informaation from the header files of the dicom slices about the resolution of the 3d image.
can you tell me… which tag data i can use.
i mean how pixel spacing, thickness and spcing between slices contribute for finding resolution of the 2d dicom slices or whole 3d image?
also if you can please suggest how to work with shear factor?
You can reliably determine the distance along the normal to the common orientation (parallel slices) using the Image Position (Patient) 0020,0032 origin offset and the Image Orientation (Patient) 0020,0037 row and column vectors. There are several formulas to compute shear. Given the image position of two neighboring slices and the slice normal you can use sohcahtoa. Alternatively, if you create the 4x4 spatial transformation use AFNI method.
I’ve only looked at gantry tilt in CT data once many years ago, but there is a record here of how I handled that particular case.
@dglen you can get validation DICOM CT scans with gantry tilt here. The In folder includes the raw DICOMs, while the Ref folder is dcm2niix attempt to save as NIfTI. For images with gantry tilt, dcm2niix will convert the raw data with the shear stored in the sform (e.g. 1F_STEREOTAXIS_Head_201.nii as well as an interpolated volume which corrects for a tilt 1F_STEREOTAXIS_Head_201_Tilt_1.nii. Not all tools support shear in datasets (e.g. any VTK based tool) and only the NIfTI sform can represent shears, the quaternion can not store this transformation.
1 Like