Add support for converting MR images

examples/ObjectGeneration.ipynb

@@ -69,14 +69,11 @@
     "df_doses = df[df[\"modality\"] == \"RTDOSE\"]\n",
     "\n",
     "for _, dose_row in df_doses.iterrows():\n",
-    "    \n",
     "    if \"EQD2_ab\" in dose_row.hashed_uid:\n",
     "        # This is an already scaled dose\n",
     "        continue\n",
-    "    \n",
-    "    df_linked_plan = df[\n",
-    "            df[\"sop_instance_uid\"] == dose_row.referenced_sop_instance_uid\n",
-    "        ]\n",
+    "\n",
+    "    df_linked_plan = df[df[\"sop_instance_uid\"] == dose_row.referenced_sop_instance_uid]\n",
     "\n",
     "    linked_plan = df_linked_plan.iloc[0]\n",
     "    ds_plan = load_object_metadata(linked_plan)\n",
@@ -85,30 +82,26 @@
     "    fractions = int(ds_plan.FractionGroupSequence[0].NumberOfFractionsPlanned)\n",
     "\n",
     "    print(f\"{dose_row.patient_id} has {fractions} fractions\")\n",
-    "    \n",
+    "\n",
     "    # Load the dose grid\n",
     "    dose_path = Path(dose_row.path).joinpath(\"RTDOSE.nii.gz\")\n",
     "    dose = sitk.ReadImage(str(dose_path))\n",
     "    dose = sitk.Cast(dose, sitk.sitkFloat64)\n",
     "\n",
     "    dose_id = f\"{dose_row.hashed_uid}_EQD2_ab{alpha_beta}\"\n",
-    "    \n",
+    "\n",
     "    if len(df_doses[df_doses.hashed_uid == dose_id]) > 0:\n",
     "        print(f\"Already converted dose for {dose_id}\")\n",
     "        continue\n",
-    "    \n",
+    "\n",
     "    # Apply the EQD2 correction\n",
-    "    eqd2_dose = dose * (((dose/actual_frac) + alpha_beta) / (2 + alpha_beta))\n",
-    "    \n",
+    "    eqd2_dose = dose * (((dose / fractions) + alpha_beta) / (2 + alpha_beta))\n",
+    "\n",
     "    # Save off the new dose grid\n",
     "    try:\n",
     "        print(f\"Saving dose grid with ID: {dose_id}\")\n",
     "        pydicer.add_dose_object(\n",
-    "            eqd2_dose,\n",
-    "            dose_id,\n",
-    "            dose_row.patient_id,\n",
-    "            linked_plan,\n",
-    "            dose_row.for_uid\n",
+    "            eqd2_dose, dose_id, dose_row.patient_id, linked_plan, dose_row.for_uid\n",
     "        )\n",
     "    except SystemError:\n",
     "        print(f\"Dose object {dose_id} already exists!\")"

pydicer/constants.py

@@ -3,6 +3,7 @@
 RT_PLAN_STORAGE_UID = "1.2.840.10008.5.1.4.1.1.481.5"
 CT_IMAGE_STORAGE_UID = "1.2.840.10008.5.1.4.1.1.2"
 PET_IMAGE_STORAGE_UID = "1.2.840.10008.5.1.4.1.1.128"
+MR_IMAGE_STORAGE_UID = "1.2.840.10008.5.1.4.1.1.4"
 
 PYDICER_DIR_NAME = ".pydicer"
 CONVERTED_DIR_NAME = "data"

pydicer/convert/data.py

@@ -26,6 +26,7 @@
     RT_STRUCTURE_STORAGE_UID,
     CT_IMAGE_STORAGE_UID,
     PET_IMAGE_STORAGE_UID,
+    MR_IMAGE_STORAGE_UID,
 )
 from pydicer.logger import PatientLogger
 
@@ -411,6 +412,34 @@ def convert(self, patient=None, force=True):
                     entry["path"] = str(output_dir.relative_to(self.working_directory))
 
                     self.add_entry(entry)
+
+                elif sop_class_uid == MR_IMAGE_STORAGE_UID:
+                    # TODO Handle inconsistent slice spacing
+                    if config.get_config("interp_missing_slices"):
+                        series_files = handle_missing_slice(df_files)
+                    else:
+                        raise ValueError("Slice Locations are not evenly spaced")
+
+                    series_files = [str(f) for f in series_files]
+                    series = sitk.ReadImage(series_files)
+
+                    output_dir = patient_directory.joinpath("images", sop_instance_hash)
+                    output_dir.mkdir(exist_ok=True, parents=True)
+
+                    nifti_file = output_dir.joinpath("MR.nii.gz")
+                    sitk.WriteImage(series, str(nifti_file))
+                    logger.debug("Writing MR Image Series to: %s", nifti_file)
+
+                    json_file = output_dir.joinpath("metadata.json")
+                    convert_dicom_headers(
+                        series_files[0],
+                        str(nifti_file.relative_to(self.output_directory)),
+                        json_file,
+                    )
+
+                    entry["path"] = str(output_dir.relative_to(self.working_directory))
+
+                    self.add_entry(entry)
                     patient_logger.eval_module_process("convert", sop_instance_hash)
 
                 elif sop_class_uid == RT_STRUCTURE_STORAGE_UID:

pydicer/preprocess/data.py

@@ -14,6 +14,7 @@
     RT_PLAN_STORAGE_UID,
     RT_STRUCTURE_STORAGE_UID,
     CT_IMAGE_STORAGE_UID,
+    MR_IMAGE_STORAGE_UID,
 )
 from pydicer.quarantine import copy_file_to_quarantine
 from pydicer.utils import read_preprocessed_data, get_iterator
@@ -109,7 +110,11 @@ def scan_file(self, file):
                 except AttributeError:
                     logger.warning("Unable to determine Reference Series UID")
 
-            elif dicom_type_uid in (CT_IMAGE_STORAGE_UID, PET_IMAGE_STORAGE_UID):
+            elif dicom_type_uid in (
+                CT_IMAGE_STORAGE_UID,
+                PET_IMAGE_STORAGE_UID,
+                MR_IMAGE_STORAGE_UID,
+            ):
                 image_position = np.array(ds.ImagePositionPatient, dtype=float)
                 image_orientation = np.array(ds.ImageOrientationPatient, dtype=float)
 
@@ -213,8 +218,6 @@ def preprocess(self, input_directory, force=True):
 
             files = [f for f in files if str(f) not in files_already_scanned]
 
-        logger.info("Found %d files to scan", len(files))
-
         result_list = []
 
         for f in get_iterator(files, unit="files", name="preprocess"):

pydicer/visualise/data.py

@@ -52,7 +52,7 @@ def visualise(self, dataset_name=CONVERTED_DIR_NAME, patient=None, force=True):
             join_working_directory=True,
         )
 
-        visualise_modalities = ["CT", "RTSTRUCT", "RTDOSE", "PT"]
+        visualise_modalities = ["CT", "MR", "RTSTRUCT", "RTDOSE", "PT"]
         df_process = df_process[df_process.modality.isin(visualise_modalities)]
 
         for _, row in get_iterator(
@@ -139,17 +139,16 @@ def visualise(self, dataset_name=CONVERTED_DIR_NAME, patient=None, force=True):
                 patient_logger.eval_module_process("visualise", row.hashed_uid)
                 logger.debug("Created CT visualisation: %s", vis_filename)
 
-            if row.modality == "PT":
+            if row.modality in ("MR", "PT"):
                 img_path = Path(row.path)
-                vis_filename = img_path.joinpath("PT.png")
+                vis_filename = img_path.joinpath(f"{row.modality}.png")
 
                 if vis_filename.exists() and not force:
                     logger.info("Visualisation already exists at %s", vis_filename)
                     continue
 
                 img = sitk.ReadImage(str(img_path.joinpath(f"{row.modality}.nii.gz")))
 
-                # TODO find linked CT and render PT on top of that
                 vis = ImageVisualiser(img)
                 fig = vis.show()
 
@@ -160,7 +159,7 @@ def visualise(self, dataset_name=CONVERTED_DIR_NAME, patient=None, force=True):
                 plt.close(fig)
 
                 patient_logger.eval_module_process("visualise", row.hashed_uid)
-                logger.debug("Created PT visualisation: %s", vis_filename)
+                logger.debug("Created %s visualisation: %s", row.modality, vis_filename)
 
             # Visualise the structures on top of their linked image
             if row.modality == "RTSTRUCT":