From 528ca5a6cbe40a289b14efa2f4d6c61ebf6c412c Mon Sep 17 00:00:00 2001 From: AKA Math Date: Tue, 12 Dec 2023 21:48:11 +0100 Subject: [PATCH 1/2] Dirty version of reading in RTDose files - need to clean up and add tests next. --- pyradise/data/__init__.py | 2 +- pyradise/data/image.py | 26 ++++++++++++- pyradise/data/subject.py | 4 +- pyradise/fileio/crawling.py | 57 ++++++++++++++++++++++++++--- pyradise/fileio/dicom_conversion.py | 14 +++++-- pyradise/fileio/series_info.py | 14 +++++++ 6 files changed, 104 insertions(+), 13 deletions(-) diff --git a/pyradise/data/__init__.py b/pyradise/data/__init__.py index 3d60ed9..410bb0a 100644 --- a/pyradise/data/__init__.py +++ b/pyradise/data/__init__.py @@ -1,5 +1,5 @@ from .annotator import Annotator -from .image import Image, ImageProperties, IntensityImage, SegmentationImage +from .image import Image, ImageProperties, IntensityImage, DoseImage, SegmentationImage from .modality import Modality from .organ import Organ, OrganAnnotatorCombination from .subject import Subject diff --git a/pyradise/data/image.py b/pyradise/data/image.py index a8d449e..40a0451 100644 --- a/pyradise/data/image.py +++ b/pyradise/data/image.py @@ -16,7 +16,7 @@ TransformInfo = TypeVar("TransformInfo") -__all__ = ["Image", "IntensityImage", "SegmentationImage", "ImageProperties"] +__all__ = ["Image", "IntensityImage", "SegmentationImage", "DoseImage", "ImageProperties"] class ImageProperties: @@ -744,6 +744,30 @@ def __str__(self) -> str: return f"SegmentationImage: {self.organ.get_name()} / {self.annotator.get_name()}" +class DoseImage(IntensityImage): + """A dose image class to specialize for properties of RTDose volumes. + + Args: + image (Union[sitk.Image, itk.Image]): The image data as :class:`itk.Image` or :class:`SimpleITK.Image`. + modality (Union[Modality, str]): The image :class:`~pyradise.data.modality.Modality` or the modality's name. + """ + + def __init__(self, image: Union[sitk.Image, itk.Image], modality: Union[Modality, str]) -> None: + + # Handle situation where RTDose intensity images are 4D - with a singleton dimensions. + if image.GetDimension() == 4: + if image.GetSize()[0] == 1: + image = image[0, :, :, :] + elif image.GetSize()[1] == 1: + image = image[:, 0, :, :] + elif image.GetSize()[2] == 1: + image = image[:, :, 0, :] + elif image.GetSize()[3] == 1: + image = image[:, :, :, 0] + + super().__init__(image, modality) + + # Preparation for next release # class DoseImage(Image): # """A dose image class including a :class:`~pyradise.data.taping.TransformTape`. diff --git a/pyradise/data/subject.py b/pyradise/data/subject.py index 2f14965..acb95c8 100644 --- a/pyradise/data/subject.py +++ b/pyradise/data/subject.py @@ -3,7 +3,7 @@ from warnings import warn from .annotator import Annotator -from .image import Image, IntensityImage, SegmentationImage +from .image import Image, IntensityImage, DoseImage, SegmentationImage from .modality import Modality from .organ import Organ @@ -385,7 +385,7 @@ def get_images_by_type(self, image_type: type) -> List[Image]: Returns: List[Image]: A list of all images of the specified type. """ - if image_type == IntensityImage: + if image_type == IntensityImage or image_type == DoseImage: return self.intensity_images elif image_type == SegmentationImage: return self.segmentation_images diff --git a/pyradise/fileio/crawling.py b/pyradise/fileio/crawling.py index 5b74f14..af612b4 100644 --- a/pyradise/fileio/crawling.py +++ b/pyradise/fileio/crawling.py @@ -14,8 +14,8 @@ from .modality_config import ModalityConfiguration from .series_info import (DicomSeriesImageInfo, DicomSeriesInfo, DicomSeriesRegistrationInfo, DicomSeriesRTSSInfo, - FileSeriesInfo, IntensityFileSeriesInfo, - SegmentationFileSeriesInfo) + DicomSeriesDoseInfo, FileSeriesInfo, + IntensityFileSeriesInfo, SegmentationFileSeriesInfo) __all__ = ["Crawler", "SubjectFileCrawler", "DatasetFileCrawler", "SubjectDicomCrawler", "DatasetDicomCrawler"] @@ -451,6 +451,28 @@ def _get_rtss_files(paths: Tuple[str, ...]) -> Tuple[str, ...]: return tuple(rtss_files) + + @staticmethod + def _get_rtdose_files(paths: Tuple[str, ...]) -> Tuple[str, ...]: + """Get all DICOM RTDOSE files in the subject directory. + + Args: + paths (Tuple[str, ...]): The DICOM file paths to check if they specify a DICOM RTDOSE file. + + Returns: + Tuple[str, ...]: The DICOM RTDOSE file paths. + """ + valid_sop_class_uid = "1.2.840.10008.5.1.4.1.1.481.2" # RT Structure Set Storage + + rtdose_files = [] + for path in paths: + dataset = load_dataset_tag(path, (Tag(0x0008, 0x0016),)) + + if dataset.get("SOPClassUID", None) == valid_sop_class_uid: + rtdose_files.append(path) + + return tuple(rtdose_files) + @staticmethod def _generate_image_infos(image_paths: Tuple[Tuple[str, ...], ...]) -> Tuple[DicomSeriesImageInfo]: """Generate the :class:`~pyradise.fileio.series_info.DicomSeriesImageInfo` entries for the DICOM file paths @@ -504,7 +526,7 @@ def _generate_rtss_info(rtss_paths: Tuple[str, ...]) -> Tuple[DicomSeriesRTSSInf rtss_paths (Tuple[str, ...]): The DICOM RTSS file paths. Returns: - Tuple[DicomSeriesRTStructureSetInfo, ...]: AThe retrieved + Tuple[DicomSeriesRTStructureSetInfo, ...]: The retrieved :class:`~pyradise.fileio.series_info.DicomSeriesRTStructureSetInfo` entries. """ infos = [] @@ -515,6 +537,27 @@ def _generate_rtss_info(rtss_paths: Tuple[str, ...]) -> Tuple[DicomSeriesRTSSInf return tuple(infos) + + @staticmethod + def _generate_rtdose_info(rtdose_paths: Tuple[str, ...]) -> Tuple[DicomSeriesImageInfo]: + """Generate the :class:`~pyradise.fileio.series_info.DicomSeriesImageInfo` entries for the DICOM file + paths specified. + + Args: + rtdose_paths (Tuple[str, ...]): The DICOM RTDOSE file paths. + + Returns: + Tuple[DicomSeriesImageInfo, ...]: The retrieved + :class:`~pyradise.fileio.series_info.DicomSeriesImageInfo` entries. + """ + infos = [] + + for path in rtdose_paths: + rtdose_info = DicomSeriesDoseInfo(path) + infos.append(rtdose_info) + + return tuple(infos) + def _export_modality_config(self, infos: Tuple[DicomSeriesInfo, ...]) -> None: """Export the retrieved :class:`~pyradise.fileio.modality_config.ModalityConfiguration` to a file. @@ -647,16 +690,20 @@ def execute(self) -> Tuple[DicomSeriesInfo, ...]: remaining_paths = tuple(set(remaining_paths) - set(registration_paths)) rtss_paths = self._get_rtss_files(remaining_paths) + remaining_paths = tuple(set(remaining_paths) - set(rtss_paths)) + + rtdose_paths = self._get_rtdose_files(remaining_paths) # generate the series infos image_infos = self._generate_image_infos(image_paths) registration_infos = self._generate_registration_infos(registration_paths, image_infos) rtss_infos = self._generate_rtss_info(rtss_paths) + rtdose_infos = self._generate_rtdose_info(rtdose_paths) # apply the modality config and write it to disk if requested - self._apply_modality_config(image_infos) + self._apply_modality_config(image_infos + rtdose_infos) - return image_infos + registration_infos + rtss_infos + return image_infos + registration_infos + rtss_infos + rtdose_infos class DatasetDicomCrawler(Crawler): diff --git a/pyradise/fileio/dicom_conversion.py b/pyradise/fileio/dicom_conversion.py index ddc60fa..7cfd30e 100644 --- a/pyradise/fileio/dicom_conversion.py +++ b/pyradise/fileio/dicom_conversion.py @@ -24,14 +24,14 @@ from pydicom.uid import (PYDICOM_IMPLEMENTATION_UID, ImplicitVRLittleEndian, generate_uid) -from pyradise.data import (IntensityImage, Modality, Organ, SegmentationImage, +from pyradise.data import (IntensityImage, Modality, Organ, SegmentationImage, DoseImage, Subject, str_to_modality) from pyradise.utils import (chunkify, convert_to_itk_image, get_slice_direction, get_slice_position, get_spacing_between_slices, load_dataset, load_dataset_tag, load_datasets) -from .series_info import (DicomSeriesImageInfo, DicomSeriesRegistrationInfo, +from .series_info import (DicomSeriesImageInfo, DicomSeriesDoseInfo, DicomSeriesRegistrationInfo, DicomSeriesRTSSInfo, RegistrationInfo, SeriesInfo) __all__ = [ @@ -2594,7 +2594,10 @@ def convert(self) -> Tuple[IntensityImage, ...]: # if no registration info is available, the image is added as is if reg_info is None: - image_ = IntensityImage(image, info.modality) + if isinstance(info, DicomSeriesDoseInfo): + image_ = DoseImage(image, info.modality) + else: + image_ = IntensityImage(image, info.modality) image_.add_data({"SeriesInstanceUID": info.series_instance_uid}) images.append(image_) @@ -2610,7 +2613,10 @@ def convert(self) -> Tuple[IntensityImage, ...]: ) image = self._transform_image(image, reg_info.transform, is_intensity=True) - image_ = IntensityImage(image, info.modality) + if isinstance(info, DicomSeriesDoseInfo): + image_ = DoseImage(image, info.modality) + else: + image_ = IntensityImage(image, info.modality) image_.add_data({"SeriesInstanceUID": info.series_instance_uid}) images.append(image_) diff --git a/pyradise/fileio/series_info.py b/pyradise/fileio/series_info.py index abc6285..ff72676 100644 --- a/pyradise/fileio/series_info.py +++ b/pyradise/fileio/series_info.py @@ -21,6 +21,7 @@ "SegmentationFileSeriesInfo", "DicomSeriesInfo", "DicomSeriesImageInfo", + "DicomSeriesDoseInfo", "DicomSeriesRegistrationInfo", "DicomSeriesRTSSInfo", "ReferenceInfo", @@ -417,6 +418,19 @@ def update(self) -> None: self._is_updated = True +class DicomSeriesDoseInfo(DicomSeriesImageInfo): + """A :class:`DicomSeriesDoseInfo` class for DICOM Dose images. In addition to the information provided by the + :class:`DicomSeriesImageInfo` class, this class contains a flag to indicate the image is a Dose volume. + + Args: + paths (Tuple[str, ...]): The paths to the DICOM image files to load. + """ + + def __init__(self, paths: Tuple[str, ...]) -> None: + super().__init__(paths) + self.is_dose_image = True + + # noinspection PyUnresolvedReferences @dataclass class ReferenceInfo: From 99c976d28abc833004a80e6852a696f406b274ba Mon Sep 17 00:00:00 2001 From: AKA Math Date: Wed, 13 Dec 2023 16:13:57 +0100 Subject: [PATCH 2/2] Fixing bug where dose volumes were int32, and unscaled. Now they have meaningful numbers. --- pyradise/data/image.py | 6 +++++- pyradise/fileio/dicom_conversion.py | 4 ++-- pyradise/fileio/series_info.py | 4 ++++ 3 files changed, 11 insertions(+), 3 deletions(-) diff --git a/pyradise/data/image.py b/pyradise/data/image.py index 40a0451..c4432c5 100644 --- a/pyradise/data/image.py +++ b/pyradise/data/image.py @@ -752,7 +752,7 @@ class DoseImage(IntensityImage): modality (Union[Modality, str]): The image :class:`~pyradise.data.modality.Modality` or the modality's name. """ - def __init__(self, image: Union[sitk.Image, itk.Image], modality: Union[Modality, str]) -> None: + def __init__(self, image: Union[sitk.Image, itk.Image], modality: Union[Modality, str], scaling_value: float) -> None: # Handle situation where RTDose intensity images are 4D - with a singleton dimensions. if image.GetDimension() == 4: @@ -765,6 +765,10 @@ def __init__(self, image: Union[sitk.Image, itk.Image], modality: Union[Modality elif image.GetSize()[3] == 1: image = image[:, :, :, 0] + # See here for why scaling is needed: https://dicom.innolitics.com/ciods/rt-dose/rt-dose/3004000e + image = sitk.Cast(image, sitk.sitkFloat32) + image = image * scaling_value + super().__init__(image, modality) diff --git a/pyradise/fileio/dicom_conversion.py b/pyradise/fileio/dicom_conversion.py index 7cfd30e..1510d84 100644 --- a/pyradise/fileio/dicom_conversion.py +++ b/pyradise/fileio/dicom_conversion.py @@ -2595,7 +2595,7 @@ def convert(self) -> Tuple[IntensityImage, ...]: # if no registration info is available, the image is added as is if reg_info is None: if isinstance(info, DicomSeriesDoseInfo): - image_ = DoseImage(image, info.modality) + image_ = DoseImage(image, info.modality, info.scaling_value) else: image_ = IntensityImage(image, info.modality) image_.add_data({"SeriesInstanceUID": info.series_instance_uid}) @@ -2614,7 +2614,7 @@ def convert(self) -> Tuple[IntensityImage, ...]: image = self._transform_image(image, reg_info.transform, is_intensity=True) if isinstance(info, DicomSeriesDoseInfo): - image_ = DoseImage(image, info.modality) + image_ = DoseImage(image, info.modality, info.scaling_value) else: image_ = IntensityImage(image, info.modality) image_.add_data({"SeriesInstanceUID": info.series_instance_uid}) diff --git a/pyradise/fileio/series_info.py b/pyradise/fileio/series_info.py index ff72676..ae946f3 100644 --- a/pyradise/fileio/series_info.py +++ b/pyradise/fileio/series_info.py @@ -428,6 +428,10 @@ class DicomSeriesDoseInfo(DicomSeriesImageInfo): def __init__(self, paths: Tuple[str, ...]) -> None: super().__init__(paths) + scaling_tag = [Tag(0x3004,0x000E)] + dataset = load_dataset_tag(self.path[0], scaling_tag) + + self.scaling_value = str(dataset.get("DoseGridScaling", 1.0)) self.is_dose_image = True