Reading Legacy Converted Enhanced Images

[hackermd]

I've converted a DICOM Series of Single-Frame MR Images to a Legacy Converted Enhanced MR Image using MultiFrameImageFactory. When I try to load the generated multi-frame image using ITK's imread function, I get the following error:

----> 1 image = itk.imread(input_filename)

~/.pyenvs/imglib/lib/python3.6/site-packages/itkExtras.py in imread(fileName, pixelType)
    443         reader = itk.ImageFileReader[ImageType].New(FileName=fileName)
    444     else:
--> 445         reader = itk.ImageFileReader.New(FileName=fileName)
    446     reader.Update()
    447     return reader.GetOutput()

~/.pyenvs/imglib/lib/python3.6/site-packages/itkTemplate.py in New(self, *args, **kwargs)
    381         cur = itk.auto_pipeline.current
    382         if self.__name__ == "itk::ImageFileReader":
--> 383             return self._NewImageFileReader(*args, **kwargs)
    384         primary_input_methods = ('Input', 'InputImage', 'Input1')
    385         if len(args) != 0:

~/.pyenvs/imglib/lib/python3.6/site-packages/itkTemplate.py in _NewImageFileReader(self, *args, **kwargs)
    418         imageIO = itk.ImageIOFactory.CreateImageIO( inputFileName, itk.ImageIOFactory.ReadMode )
    419         if not imageIO:
--> 420             raise RuntimeError("No ImageIO is registered to handle the given file.")
    421         componentTypeDic= {"float": itk.F, "double": itk.D,
    422         "unsigned_char": itk.UC, "unsigned_short": itk.US, "unsigned_int": itk.UI,

RuntimeError: No ImageIO is registered to handle the given file.

Support for Legacy Converted Enhanced Image Storage classes has only recently been added to GDCM: https://github.com/malaterre/GDCM/blob/513a75c2f4dc704549223b831879c4a5fab00076/Source/DataStructureAndEncodingDefinition/gdcmMediaStorage.cxx#L144

This functionality is not yet in the GDCM code that comes with ITK: https://github.com/InsightSoftwareConsortium/ITK/blob/6dffe29eb245f77ac27a17d0398fa3073cb6bf4c/Modules/ThirdParty/GDCM/src/gdcm/Source/DataStructureAndEncodingDefinition/gdcmMediaStorage.cxx#L143

Is there a change that the GDCM source code will be updated for the 5.0 release?

See also related SimpleITK issue

[blowekamp]

I have created a patch which updates ITK's GDCM to the latest: http://review.source.kitware.com/#/c/23486/

Would you be able to build ITK with this patch and check to see if it adds support for your particular file?

xref: SimpleITK/SimpleITK#440

[hackermd]

Thanks @blowekamp! It worked!

These are the steps to reproduce:

1) Install ITK from source using patch 23485:

git clone git://itk.org/ITK.git /tmp/itk
cd /tmp/itk
git fetch http://review.source.kitware.com/ITK refs/changes/85/23485/1 && git checkout FETCH_HEAD
mkdir -p /tmp/local/share
cat << EOF > /tmp/local/share/config.site                                                                                                                                       
CPPFLAGS=-I/tmp/local/include                                                                                                                                                
LDFLAGS=-L/tmp/local/lib                                                                                                                                                     
EOF
mkdir /tmp/itk/build
cd /tmp/itk/build
cmake ../ -DCMAKE_INSTALL_PREFIX=/tmp/local -DITK_WRAP_PYTHON:BOOL=ON
make
make install

2) Test reading Legacy Converted Enhanced MR Image:

cd /tmp/local/lib/python3.6/site-packages
python test_legacy_converted.py /tmp/testfile.dcm

Where the content of test_legacy_converted.py looks as follows

import itk
import sys

if __name__ == '__main__':

    input_filename = sys.args[1]
    image = itk.imread(input_filename)
    metadata = image.GetMetaDataDictionary()
    print(metadata['0008|0016'])  # expected: '1.2.840.10008.5.1.4.1.1.4.4'

[blowekamp]

You can build SimpleITK against this ITK to get this support too, if you prefer that interface.

[hackermd]

I've run a few tests to compare ITK Image objects resulting from reading the original series of legacy MR Image and Legacy Converted Enhanced MR Image instances. The actual pixel data are identical:

import itk
import numpy as np

ImageType = itk.Image[itk.ctype('signed short'), 3]

series_reader = itk.ImageSeriesReader[ImageType].New(FileNames=legacy_filenames)
legacy_image = series_reader.GetOutput()
legacy_array = itk.GetArrayFromImage(legacy_image)

file_reader = itk.ImageFileReader[ImageType].New(FileName=enhanced_filename)
enhanced_image = file_reader.GetOutput()
enhanced_array = itk.GetArrayFromImage(enhanced_image)

np.testing.assert_array_equal(legacy_array, enhanced_array)

However, the spacing is different:

assert legacy_image.GetImageDimension() == enhanced_image.GetImageDimension()
for i in range(legacy_image.GetImageDimensions()):
    assert legacy_image.GetSpacing().GetElement(i) == enhanced_image.GetSpacing().GetElement(i), \
        "spacing along image dimension {} is not identical".format(i)

This is because the ImageFileReader does not interpret the DICOM metadata correctly:

legacy_metadata = series_reader.GetMetaDataDictionaryArray()
legacy_metadata[0]['0028|0030']  # Pixel Spacing
legacy_metadata[0]['0018|0088']  # Spacing Between Slices

enhanced_metadata = enhanced_image.GetMetaDataDictionary()
enhanced_metadata['5200|9229']  # Shared Functional Group Sequence

The problem is that in case of the Legacy Converted Enhanced MR Image information object definition class the Pixel Spacing and Spacing Between Slices attributes are contained in the Pixel Measures Sequence attribute within the Shared Functional Group Sequence attribute (see PS3.3 A.71.4 Legacy Converted Enhanced MR Image Functional Group Macros and PS3.3 C.7.6.16.2.1 Pixel Measures Macro).

This is the output of the GDCM program dcmdump:

    (0028,9110) SQ (Sequence with undefined length #=1)     # u/l, 1 PixelMeasuresSequence
      (fffe,e000) na (Item with undefined length #=2)         # u/l, 1 Item
        (0018,0050) DS [5]                                      #   2, 1 SliceThickness
        (0028,0030) DS [0.9375\0.9375]                          #  14, 2 PixelSpacing
      (fffe,e00d) na (ItemDelimitationItem)                   #   0, 0 ItemDelimitationItem
    (fffe,e0dd) na (SequenceDelimitationItem)               #   0, 0 SequenceDelimitationItem

Key "5200|9229" of the Shared Functional Group Sequence data element is not even contained in the MetaDataDictionary.

[hackermd]

The relevant attributes may also end up in the Unassigned Per-Frame Converted Attributes Sequence attributes (see PS3.3 C.7.6.16.2.25.2 Unassigned Per-Frame Converted Attributes Macro).

[hackermd]

The return value of legacy_image.GetSpacing() is itkVectorD3 ([0.9375, 0.9375, 5.99997]), while enhanced_image.GetSpacing() returns itkVectorD3 ([1, 1, 1]).

To me, it seems this problem originates in GDCM. It appears GDCM looks for the Pixel Spacing and Spacing Between Slices attributes only at the root level of the dataset, but not in the Pixel Measures Sequence of the Shared Functional Groups Sequence (see gdcmSpacing.h). If GDCM cannot determine the spacing, it returns ones as default values (see gdcmImage.h).

@blowekamp Do you think this evaluation is accurate?

[hackermd]

This is worrying, because it probably also affects multi-frame images of other SOP classes, such as Enhanced MR Image (see Enhanced MR Image Functional Group Macros) and Enhanced CT Image (see Enhanced CT Image Functional Group Macros).

[malaterre]

@hackermd I am one of the main author of GDCM and I can guarantee that GDCM is reading the Pixel Spacing from the correct location (esp. Per-Frame or Shared Function Groups Sequence). Please add a link to your generated DICOM DataSet?

[malaterre]

@hackermd If you cannot redistribute your DICOM DataSet, please include the output from gdcminfo my_enh.dcm also specify the version of GDCM you are using.

[hackermd]

@malaterre thanks for your feedback. Please find below the output of gdcminfo (GDCM version 2.8.6) for the DICOM file.

MediaStorage is 1.2.840.10008.5.1.4.1.1.4.4 [Legacy Converted Enhanced MR Image Storage]
TransferSyntax is 1.2.840.10008.1.2.1 [Explicit VR Little Endian]
NumberOfDimensions: 3
Dimensions: (256,256,26)
SamplesPerPixel    :1
BitsAllocated      :16
BitsStored         :16
HighBit            :15
PixelRepresentation:1
ScalarType found   :INT16
PhotometricInterpretation: MONOCHROME2
PlanarConfiguration: 0
TransferSyntax: 1.2.840.10008.1.2.1
Origin: (0,0,0)
Spacing: (1,1,1)
DirectionCosines: (1,0,0,0,1,0)
Rescale Intercept/Slope: (0,1)
Orientation Label: AXIAL

I will determine whether I can share the DICOM object.

[hackermd]

@malaterre See this comment on the corresponding SimpleITK github issue on how to generate the enhanced image object.

[hackermd]

Reading multi-frame images now works with the 5.0 Beta 1 release. However, the values of the spacing, origin and direction attributes are not set correctly but represent the defaults.

[thewtex]

@hackermd please try the 5.0 Beta 3 release where this should be addressed per:

https://github.com/InsightSoftwareConsortium/ITK/pull/315 https://github.com/InsightSoftwareConsortium/ITK/pull/208

python -m pip install --upgrade --pre itk

[hackermd]

@thewtex thanks for your follow up.

I just tested with version 5.0b3. The enhanced image is now read correctly via the imread() function:

import itk

enhanced_image = itk.imread(enhanced_filename)
enhanced_image.GetSpacing()

This returns itkVectorD3 ([0.9375, 0.9375, 6]).

However, it does not work with the ImageFileReader class:

import itk

file_reader = itk.ImageFileReader.New(FileName=enhanced_filename)
enhanced_image = file_reader.GetOutput()
enhanced_image.GetSpacing()

This returns itkVectorD3 ([1, 1, 1]).

Notably, the output of ImageSeriesReader is also not correct for the series of legacy image instances:

import itk

series_reader = itk.ImageSeriesReader.New(FileNames=legacy_filenames)
legacy_image = series_reader.GetOutput()
legacy_image.GetSpacing()

This returns itkVectorD3 ([1, 1, 1]) as well.

The series is read correctly with SimpleITK version 1.1.0:

import SimpleITK as sitk

series_reader = sitk.ImageSeriesReader()
series_reader.SetFileNames(legacy_filenames)
legacy_image = series_reader.Execute()
legacy_image.GetSpacing()

This returns (0.9375, 0.9375, 5.999974727630615).

[thewtex]

@hackermd thanks for the update. I am glad to hear it is working with itk.imread.

Is there a a series that we can use to reproduce the legacy filenames issue and create a regression test?

[hackermd]

@thewtex please see this comment on issue #440 of SimpleITK repository on how to reproduce

[hackermd]

I tested with series 1.3.6.1.4.1.9328.50.50.104552199813796461463367647403696421653 of study 1.3.6.1.4.1.9328.50.50.78410421624738603252174560245217471974.

Expected spacing: itkVectorD3 ([0.859375, 0.859375, 7.5])

[thewtex]

::jumps from Tardis ☎️ ::

@hackermd Can you please try the recent itk 5.4rc3 python packages:

pip install --upgrade -pre itk

and let us know how it goes? These include: https://github.com/InsightSoftwareConsortium/ITK/pull/4521

[thewtex]

More improvements made in ITK 5.4 RC 4. Available in the itk 5.4.0 wheels.