Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Dicom series loaded with incorrect origin and spacing #2117

Closed
chris-rapson-formus opened this issue May 8, 2024 · 11 comments
Closed

Dicom series loaded with incorrect origin and spacing #2117

chris-rapson-formus opened this issue May 8, 2024 · 11 comments

Comments

@chris-rapson-formus
Copy link

chris-rapson-formus commented May 8, 2024

Describe the bug
When loading a Dicom series from a folder, I am getting the origin as [0, 0, 0] and spacing as [1, 1, 1], which doesn't match the true origin and spacing.

This has been observed with multiple Dicom series, which have the same source. I haven't observed this before on any Dicom series loaded with SimpleITK. I imagine there is something unusual about the Dicoms, but at the same time, the relevant metadata seems to be ok.

There is no error or warning message.

When I inspect the metadata either with SimpleITK, pydicom or MITK Workbench, the spacing is [0.503906, 0.503906, 1]. I believe the origin should be [-208.748, -260.748, 1014.1], however MITK Workbench shows the origin as [0, 0, 0]. Here is a copy of the metadata from MITK Workbench:

GetGeometryForTimeStep(0):       SlicedGeometry3D (0x7f0d80018d90)
         IndexToWorldTransform: 
        Matrix: 
          0.503906 0 0 
          0 0.503906 0 
          0 0 1 
        Offset: [0, 0, 0]
        Center: [0, 0, 0]
        Translation: [0, 0, 0]
        Inverse: 
          1.9845 0 0 
          0 1.9845 0 
          0 0 1 
        Scale : 1 1 1 
         BoundingBox:         ( 0,512 0,512 0,527  )
         Origin: [0, 0, 0]
         ImageGeometry: 1
         Spacing: [0.503906, 0.503906, 1]
         EvenlySpaced: 1
         DirectionVector: [0, 0, 1]
         Slices: 527

To Reproduce
Steps to reproduce the behavior:

  1. Operating system, version, and architecture
  • OS: Ubuntu 20.04.6
  • Architecture: x86_64 (Intel i7-6800K)
  1. Programming language, and version Python 3.8
  2. Version of SimpleITK 2.3.1
  3. How was SimpleITK installed?
  • pip install -r requirements.txt (where requirements includes simpleitk==2.3.1)
  1. A minimal working example which causes the error.
import glob

import SimpleITK as sitk
import pydicom

dicom_folder = '/path/to/dicoms'
dicom_file_paths = glob.glob(dicom_folder + '/*')

# load as normal... origin and spacing are wrong, but metadata is correct
file_reader = sitk.ImageFileReader()
file_reader.SetFileName(dicom_file_paths[0])
file_reader.ReadImageInformation()
series_ID = file_reader.GetMetaData('0020|000e')
sorted_file_names = sitk.ImageSeriesReader.GetGDCMSeriesFileNames(dicom_folder, series_ID)
img = sitk.ReadImage(sorted_file_names)
print(f'{img.GetSpacing()=}')  # img.GetSpacing()=(1.0, 1.0, 1.0)
print(f'{img.GetOrigin()=}')  # img.GetOrigin()=(0.0, 0.0, 0.0)

ipp0 = file_reader.GetMetaData('0020|0032')  # Image Position (Patient)
iop0 = file_reader.GetMetaData('0020|0037')  # Image Orientation (Patient)
px_spacing0 = file_reader.GetMetaData('0028|0030')  # Pixel Spacing
# slice_spacing0 = file_reader.GetMetaData('0018|0088')  # Spacing Between Slices... missing for these DICOMs

file_reader_1 = sitk.ImageFileReader()
file_reader_1.SetFileName(dicom_file_paths[1])
file_reader_1.ReadImageInformation()
ipp1 = file_reader_1.GetMetaData('0020|0032')  # Image Position (Patient)
z0 = float(ipp0.split("\\")[-1])
z1 = float(ipp1.split("\\")[-1])

print(f'{ipp0=}')  # ipp0='-208.748\\-260.748\\1239.7'
print(f'{iop0=}')  # iop0='1\\0\\0\\-0\\1\\-0 '
print(f'{px_spacing0=}')  # px_spacing0='0.503906\\0.503906 '
# print(f'{slice_spacing0=}')
print(f'slice spacing from IPP={z1-z0}')  # slice spacing from IPP=40.200000000000045

# try loading a single image... origin and spacing are still wrong
single_image = sitk.ReadImage(dicom_file_paths[0])
print(f'{single_image.GetOrigin()=}')  # single_image.GetOrigin()=(0.0, 0.0, 0.0)
print(f'{single_image.GetSpacing()=}')  # single_image.GetSpacing()=(1.0, 1.0, 1.0)

# try loading with pydicom. Confirm that metadata is correct.
scan = [pydicom.dcmread(f) for f in dicom_file_paths]
min_z = 9999998
second_min_z = 9999999
min_z_index = -1
second_min_z_index = -1
for ix, sc in enumerate(scan):
    sc_z = sc.ImagePositionPatient[2]
    if sc_z < min_z:
        second_min_z = min_z
        min_z = sc_z
        second_min_z_index = min_z_index
        min_z_index = ix
    elif sc_z < second_min_z:
        second_min_z = sc_z
        second_min_z_index = ix
ipp0_pydicom = scan[min_z_index].ImagePositionPatient
ipp1_pydicom = scan[second_min_z_index].ImagePositionPatient

print(f'{ipp0_pydicom=}')  # ipp0_pydicom=[-208.748, -260.748, 1014.1]
print(f'slice spacing from IPP with pydicom={ipp1_pydicom[2] - ipp0_pydicom[2]}')  # slice spacing from IPP with pydicom=0.6000000000000227
@gdevenyi
Copy link
Contributor

gdevenyi commented May 8, 2024

Have you tried generating a NIFTI with say dcm2niix, to see what it thinks the coordinate system is?

@blowekamp
Copy link
Member

Also if you read an individual file is the meta-data information correct?

Getting a list of the DICOM tags would be useful: https://simpleitk.readthedocs.io/en/master/link_DicomImagePrintTags_docs.html

Are the original images able to be shared?

@chris-rapson-formus
Copy link
Author

@gdevenyi I haven't used dcm2niix before, so hopefully this is correct. I installed it with pip, and then ran dcm2niix /path/to/dicoms. It output a JSON file and a .nii file. I can't see anything relevant in the JSON. When loading the .nii file into MITK Workbench, it now shows the correct spacing and origin:

GetGeometryForTimeStep(0):       SlicedGeometry3D (0x560bc86dc590)
         IndexToWorldTransform: 
        Matrix: 
          0.503906 -0 0 
          0 -0.503906 0 
          0 0 0.600037 
        Offset: [-208.748, -3.25201, 1014.1]
        Center: [0, 0, 0]
        Translation: [-208.748, -3.25201, 1014.1]
        Inverse: 
          1.9845 0 0 
          0 -1.9845 0 
          0 0 1.66656 
        Scale : 1 1 1 
         BoundingBox:         ( 0,512 0,512 0,527  )
         Origin: [-208.748, -3.25201, 1014.1]
         ImageGeometry: 1
         Spacing: [0.503906, 0.503906, 0.600037]

@blowekamp If I've understood how file_reader works, then I think it already reads in metadata from an individual file in my example above? The same for pydicom. I went back and checked ImagePositionPatient from all files, and they're all fine. It seems to me that either SimpleITK isn't using ImagePositionPatient, or the problem occurs when combining them into a 3d image.

I've attached the DICOM tags that were output from the script you linked. Maybe you'll be able to spot something there. At this stage I don't think the images can be shared, but I'll double check that.
dicom_tags.txt

@chris-rapson-formus
Copy link
Author

I've now been told that the images were anonymised using Materialise Mimics' standard anonymisation function. The tags that were modified were:

  • 0008,1110 Referenced Study Sequence
  • 0008,1140 Referenced Image Sequence
  • 0008,2112 Source Image Sequence
  • 0010,0010 Patient's Name
  • 0010,0020 Patient ID
  • 0020,000D Study Instance UID
  • 0020,000E Series Instance UID
  • 0020,0010 Study ID
  • 0020,0052 Frame of Reference UID

But again, I don't think the problem is in the DICOMs, because other programs can process the metadata correctly. My understanding was that origin and spacing depended entirely on Image Position Patient and Image Orientation Patient, and that images were sorted by the z coordinate, i.e. the third entry of Image Position Patient.

If I was going to look into the source code myself to see if I can spot anything, where should I look?

@gdevenyi
Copy link
Contributor

@chris-rapson-formus this code that does this may be the upstream https://github.com/InsightSoftwareConsortium/ITK code, perhaps a parallel issue should be opened there

@chris-rapson-formus
Copy link
Author

chris-rapson-formus commented May 16, 2024

By comparing the metadata with another file that loads successfully, I realised that the problematic file is missing three relevant tags:
-0018,0050 Slice Thickness (Required, Empty if Unknown (2))
-0018,0088 Spacing Between Slices (optional)
-0020,1041 Slice Location (optional)

I assume ITK works fine when optional tags are missing. (But that assumption could be wrong?)

I also assumed that the Slice Thickness shouldn't be relevant for calculating the spacing or the origin. Maybe this is a result of missing a required tag? Although if that's the case, I would have hoped that

  • ITK would treat a missing tag the same as an empty one
  • there would be a warning that the DICOMs are missing required information, rather than failing silently.

I'll try manually adding these tags and see what happens...

  • no difference. GetSpacing() still returns (1.0, 1.0, 1.0) and GetOrigin() still returns (0.0, 0.0, 0.0)

@blowekamp
Copy link
Member

By comparing the metadata with another file that loads successfully, I realised that the problematic file is missing three relevant tags:
-0018,0050 Slice Thickness (Required, Empty if Unknown (2))
-0018,0088 Spacing Between Slices (optional)
-0020,1041 Slice Location (optional)

How did you determine this was the problematic file? If this file is omitted from the series it load? If the series is truncated before this file it loads correctly?

Is it a viable options to load the series as individual slices, then join the series ( maybe the joint series or tile filter ) into the 3d volume with a computed spacing ( and other metadata ) manually set?

@chris-rapson-formus
Copy link
Author

Sorry, I should have said "problematic files". All of the files in this CT scan are missing those tags.

Although, it seems like those missing tags weren't actually the root cause, since I still get the same spacing and origin after adding them in.

Yes, manually calculating the spacing is a viable option, but that doesn't seem like a good way to use SimpleITK. Since there are no error/warning messages, I would have to manually calculate the spacing for every case. I'd essentially be rolling my own functionality that might have other bugs or become obsolete.

@chris-rapson-formus
Copy link
Author

The ITK maintainers were able to solve this for me: the issue is
0008,0016 SOP Class UID: 1.2.840.10008.5.1.4.1.1.7
which indicates "Secondary Capture Image Storage". This format is not supported by a released version of ITK yet. It is supported in the most up-to-date release candidate. So the bug has been fixed, I just have to wait for the release. Closing now.

@blowekamp
Copy link
Member

@chris-rapson-formus The development version of SimpleITK is built against the ITK v5.4rc04. You can download it for the GitHub release pages:
https://github.com/SimpleITK/SimpleITK/releases/tag/latest

Please let us know if that addresses the issue.

@chris-rapson-formus
Copy link
Author

chris-rapson-formus commented May 21, 2024

Thanks @blowekamp yes that dev version of SimpleITK gives the correct spacing and origin using the script from my first comment in this thread :-)

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Projects
None yet
Development

No branches or pull requests

3 participants