InsightSoftwareConsortium / ITK

Insight Toolkit (ITK) -- Official Repository. ITK builds on a proven, spatially-oriented architecture for processing, segmentation, and registration of scientific images in two, three, or more dimensions.
https://itk.org
Apache License 2.0
1.38k stars 661 forks source link

ITK SpatialOrientation is opposite of DICOM and ITK annotation. #1163

Open blowekamp opened 4 years ago

blowekamp commented 4 years ago

The spatial orientation enums are defined in the following file: https://github.com/InsightSoftwareConsortium/ITK/blob/master/Modules/Core/Common/include/itkSpatialOrientation.h#L39-L40

DICOM is defined in LPS, where L indication from right-to-left. However, the enums are documented as R meaning right-to-left. And the tool converted an ITK's identity direction into RIP, where we are documented as being in LPS.

The itk::OrientImageFilter uses these enums and labeling. The filter and tools with the SpatialOrientation namespace are consistent, but they are defined in the opposite direction as DICOM, and common medial imaging standards.

😕

seanm commented 4 years ago

CC @seanm @vfonov

stale[bot] commented 4 years ago

This issue has been automatically marked as stale because it has not had recent activity. Thank you for your contributions.

dzenanz commented 2 years ago

Would moving https://github.com/InsightSoftwareConsortium/ITKSimpleITKFilters/blob/master/include/itkDICOMOrientImageFilter.h here resolve this issue?

kapsner commented 11 months ago

We face a similar issue when trying to convert MRI DWI and T1 and T2 series to nifti using SimpleITK.

We found out that especially some of the DWI-Series need to be flipped across an axis (mostly z-axis) so that the resulting niftis are displayed correctly. The axis across which needs to be flipped can be determined by some heurisitc, depending on the directions of the origins.

For example, the reference T1 is:

Spacing: (1.3671875, 1.3671875, 4.400000042385526) Origin: (-168.94673109055003, -88.77255249023402, -108.12769317626999) Direction: (1.0, 0.0, 0.0, 0.0, 1.0, 0.0, 0.0, 0.0, 1.0)

The corresponding DWI is:

Spacing: (1.3671875, 1.3671875, 4.400000042385527) Origin: (-168.94673109055, -88.772552490234, 50.272308349609) Direction: (1.0, 0.0, 0.0, 0.0, 1.0, 0.0, 0.0, 0.0, 1.0)

-> Heuristic: the z-components of the origin point into different directions thus the DWI needs to be flipped across this axis.

This is pretty robust, although we also encountered some cases where the z-component is also negative (e.g. -25) but it needs to be flipped as well.

Is there any possibility to find out more robulsty when series need to be flipped?

Using the DICOMOrientImageFilter shows "LPS" for both sequences so it doesn't help us. We have the assumption that somehow the origin needs to be taken into account in order to determine if a volume need to be flipped (according to a reference image) or not.

Would anyone have an idea how to continue on from this?

Thanks in advance, Lorenz

dzenanz commented 11 months ago

You seem to be dealing with a case of wrong metadata. There is no good solution to this. One possibility is for a human to review this, and flip axes as necessary. Another is to compute registration similarity metric between the DWI and T1/T2, both with and without axis flips and try to determine if a flip is needed.

kapsner commented 11 months ago

Hi @dzenanz , Thanks for your reply. You are right, we have now found out that the volumes corresponding to the different b-values of the DWI were read in false direction, leading to the above mentioned behaviour.

In our case, all slices of the different b-values that belong to the same DWI are stored inside one single series-folder. In order to allocate slices to a specific b-value, we implemented some logic that was kind of:


def tag_reader(slicepath):
    reader = sitk.ImageFileReader()
    reader.SetFileName(slicepath)
    reader.LoadPrivateTagsOn()
    reader.ReadImageInformation()
    return reader

dicoms = sitk.ImageSeriesReader_GetGDCMSeriesFileNames(filepath)

b_value_list = []
file_dict = {}

for dcm in dicoms:
    metadata_reader = tag_reader(slicepath=dcm)
    # set dicom tag for b-values
    bval_dicom_tag = "0019|100c"  # Diffusion b-value, at least for siemens scanners

    try:
        bval = int(metadata_reader.GetMetaData(bval_dicom_tag).strip())
        b_value_list.append(bval)
    except Exception:
        raise RuntimeError("no b-values found for this series")

    if bval in file_dict.keys():
        file_dict[bval].append(dcm)
    else:
        file_dict[bval] = [dcm]

This file_dict now contains the paths that are associated with each b-value.

When providing the values of this dict-items directly to the imageseriesreader, the resulting sitk-image (and nifti) are flipped across the z-axis:

for _bval, _dcmlist in file_dict.items():
    reader = sitk.ImageSeriesReader()
    reader.SetFileNames(_dcmlist)
    img = reader.Execute()

However, what we have then identified is that the resulting dcm-path-lists in the file_dict are somehow ordered exactly in opposite direction as compared to taking the slice-location (tag "0020|1041") and instance-number (tag "0020|0013") into account for sorting the dcm-file paths.

The latter leads also to the correct / expected results, thus the sorting of the list of filenames provided to sitk.ImageSeriesReader().SetFileNames() plays an important role and in the case of a 4D DWI-series, this might somewhat not be detected correctly with sitk.ImageSeriesReader_GetGDCMSeriesFileNames (or it is messed up when iterating over the resulting list and slices corresponding to specific b-values are allocated manually).

For us, this is solved now; however, the reason why I described our approach in this detail is that we were wondering, if this is something that can be solved more easily (maybe this "correct" sorting is already implemented somewhere) or if this finding can somehow help to improve existing code.

dzenanz commented 11 months ago

@blowekamp or @issakomi might have some opinions.

zivy commented 10 months ago

Hello @kapsner,

This is a feature not a bug 😉 .

The ImageSeriesReader is a generic class which has no idea about DICOM. It assumes/expects the images it is given to be ordered correctly. The ImageSeriesReader_GetGDCMSeriesFileNames function identifies all the files that belong to the first series in a given directory based on the study+series UIDs (if the seriesID parameter value is provided it will list the files in that series). Sometimes this isn't enough to identify a series and you need to turn the useSeriesDetails on (sometimes this isn't enough too).

Possibly take a look at the SimpleITK code in the characterize_data.py script, specifically the inspect_series function. It traverses a directory structure and aggregates all the files that belong to each series into a list. You could do something similar with the key being {study}:{sid}:{bval}. In any case, nicer to avoid repeatedly creating an ImageFileReader as done in the tag_reader function.

Also look at how the inspect_single_series reads each series (creating a temp dir and softlinking or copying the original files) - this may be overkill for your case. It does address the situation where slices from the same series are in different directories and have the same file name (e.g. a/1.dcm b/1.dcm).

Welcome to the wonderful world of DICOM.

issakomi commented 10 months ago

And the tool converted an ITK's identity direction into RIP, where we are documented as being in LPS.

@blowekamp Following only applies to the first post in the thread:

ITK identity is RAI, not RIP. I am sure that you know the rest of this post, just for the record.

The "from" and "to" (Right-to-left vs. right-to-Left) terminology diverged long ago. The ITK space is the same as the DICOM space. But the ITK's identity is ITK_COORDINATE_ORIENTATION_RAI and the DICOM space is often called "DICOM LPS". I don't think that the DICOM standard defines any three-letter codes anywhere, it's kind of slang that somehow got established, probably later than the ITK orientation code was written. I personally always mention RAI-Code for ITK notations for clarity (like ITK-Snap does).


RAI-Code RAI is

1  0  0
0  1  0
0  0  1

the same as DICOM image with direction cosines 1,0,0,0,1,0 (0x0020,0x0037 "Image Orientation (Patient)")


RAI-Code LPS is

-1   0   0
 0  -1   0
 0   0  -1

is not defined in DICOM space, because DICOM does not provide 3rd vector, 3rd vector is defined as right-hand cross product of the 1st and 2nd