Why is the first axis swapped in ITK versus SimpleITK?

Community
, ,
1

I want to read a DICOM series with ITK (5.1rc) and compare it with the SimpleITK version:

def read_dcm_series_sitk(files):
    file_reader = sitk.ImageFileReader()
    file_reader.SetFileName(str(files[0]))
    file_reader.ReadImageInformation()
    series_ids = [file_reader.GetMetaData('0020|000e')]
    with tempfile.TemporaryDirectory() as tmpdir_name:
        for f in files:
            os.symlink(os.path.abspath(f), os.path.join(tmpdir_name, os.path.basename(f)))
        sorted_filenames = sitk.ImageSeriesReader.GetGDCMSeriesFileNames(tmpdir_name, series_ids[0])
        sitk_image = sitk.ReadImage(sorted_filenames)

    
    return sitk_image

def read_dcm_series_itk(files):
    PixelType = itk.ctype('double') # is int32 available?
    ImageType = itk.Image[PixelType, 3]

    reader = itk.ImageSeriesReader[ImageType].New()
    dicomIO = itk.GDCMImageIO.New()
    reader.SetImageIO(dicomIO)
    reader.SetFileNames(list(files))

    reader.Update()
    itk_image = reader.GetOutput()

    return itk_image

I now have a list of filenames which constitute a complete CT volume (after applying slope and intercept which ITK internally does becomes int32)

So I read as follows:

itk_image = read_dcm_series_itk(fns)
sitk_image = read_dcm_series_sitk(fns)

Nevertheless, the first axis is flipped:

np.all(itk.array_from_image(itk_image)[::-1] == sitk.GetArrayFromImage(sitk_image))

True

Is this an oversight of mine, and expected behaviour?

2

Hello,

It looks like the (arbitrarily ordered) file list is directly being used in native ITK python series reader. While in SimpleITK the series ID and the GDCM series sorting functionality are being used.

You can see SimpleITK’s implementation of GetGDCMSeriesFileNames for details of the code that needs to be added to the read_dcm_series_itk.

You can inspect the different between sorted_filenames and files to see if that accounts for the difference.

3 Likes
3

If the files are sorted correctly, a shorter call is:

itk_image = itk.imread(list(files))
1 Like
4

Thanks for the replies! I presort the images based on the dicom headers (InstanceNumber) as I understood how that works, but you are right, @blowekamp that is an easy check. Will write down my solution here when I figured it out.

By the way, the reason that I am looking into ITK versus SimpleITK is the new xarray ability in ITK and the back-and-forth passing to ITK images and the ability to store metadata in the xarray (for instance DCE-MRI you can store a lot of the acquisition parameters in the xarray and set chunking when saving). For most of my applications SimpleITK is by far enough but this xarray is quite useful!

1 Like