{"metadata":{"kernelspec":{"language":"python","display_name":"Python 3","name":"python3"},"language_info":{"name":"python","version":"3.11.13","mimetype":"text/x-python","codemirror_mode":{"name":"ipython","version":3},"pygments_lexer":"ipython3","nbconvert_exporter":"python","file_extension":".py"},"kaggle":{"accelerator":"none","dataSources":[{"sourceId":99552,"databundleVersionId":13851420,"sourceType":"competition"}],"dockerImageVersionId":31089,"isInternetEnabled":true,"language":"python","sourceType":"notebook","isGpuEnabled":false}},"nbformat_minor":4,"nbformat":4,"cells":[{"cell_type":"markdown","source":"The [Med3D](https://arxiv.org/pdf/1904.00625) study reports that medical volumes often exhibit heterogeneous voxel spacing due to differences in scanners and acquisition protocols. To mitigate this variability, the median spacing for each axis (x, y, z) is calculated within each domain, where a domain is defined by specific equipment and acquisition protocols. This approach ensures that the estimated spacing remains robust against outliers and minor data errors.\n\nIn this notebook, a domain is defined by imaging modality, resulting in four distinct domains. For each domain, the median spacing of each series is calculated. The method for determining spacing in each series depends on the type of DICOM file present:\n- For series containing a single DICOM file representing the entire volume (3D), spacing is obtained from the attribute specifying the spacing for the whole volume.\n- For series containing a collection of DICOM files from slices (2D), spacing is determined as the mode of all spacings. It is assumed that most spacings within each series are consistent. Instances with differing spacings are considered data errors and are excluded from the final median computation for each domain, by using only the mode of all spacings.\n\n- To calculate the normalized volume dimension, we compute:\n```python\nnew_size = (current_size / current_spacing) * median_spacing\n```","metadata":{}},{"cell_type":"code","source":"import os\nfrom concurrent.futures import ProcessPoolExecutor, as_completed\nimport multiprocessing\nfrom tqdm import tqdm\n\nimport pydicom\n\nimport numpy as np\nimport polars as pl\n\nimport pickle","metadata":{"_uuid":"8f2839f25d086af736a60e9eeb907d3b93b6e0e5","_cell_guid":"b1076dfc-b9ad-4769-8c92-a6c4dae69d19","trusted":true,"execution":{"iopub.status.busy":"2025-09-26T15:03:05.168472Z","iopub.execute_input":"2025-09-26T15:03:05.168742Z","iopub.status.idle":"2025-09-26T15:03:06.571963Z","shell.execute_reply.started":"2025-09-26T15:03:05.168721Z","shell.execute_reply":"2025-09-26T15:03:06.570837Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"BASE_PATH = \"/kaggle/input/rsna-intracranial-aneurysm-detection\"\nSERIES_PATH = f\"{BASE_PATH}/series\"","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-09-26T15:03:06.573739Z","iopub.execute_input":"2025-09-26T15:03:06.574113Z","iopub.status.idle":"2025-09-26T15:03:06.578881Z","shell.execute_reply.started":"2025-09-26T15:03:06.574082Z","shell.execute_reply":"2025-09-26T15:03:06.577734Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"train_df = pl.read_csv(f\"{BASE_PATH}/train.csv\")\ntrain_df.head(2)","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-09-26T15:03:06.58016Z","iopub.execute_input":"2025-09-26T15:03:06.580578Z","iopub.status.idle":"2025-09-26T15:03:06.933527Z","shell.execute_reply.started":"2025-09-26T15:03:06.58054Z","shell.execute_reply":"2025-09-26T15:03:06.932396Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"series_uid, modalities = train_df[\"SeriesInstanceUID\"].to_numpy(), train_df[\"Modality\"].to_numpy()","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-09-26T15:03:06.935979Z","iopub.execute_input":"2025-09-26T15:03:06.936311Z","iopub.status.idle":"2025-09-26T15:03:06.947569Z","shell.execute_reply.started":"2025-09-26T15:03:06.936286Z","shell.execute_reply":"2025-09-26T15:03:06.946319Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"def make_domain_spacing_dict():\n    return { modality: [] for modality in np.unique(modalities) }\n\nmake_domain_spacing_dict()","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-09-26T15:03:06.94844Z","iopub.execute_input":"2025-09-26T15:03:06.948716Z","iopub.status.idle":"2025-09-26T15:03:06.973471Z","shell.execute_reply.started":"2025-09-26T15:03:06.948693Z","shell.execute_reply":"2025-09-26T15:03:06.972618Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"def get_serie_info(serie_uid, domain, store_filenames):\n    filenames = os.listdir(f\"{SERIES_PATH}/{serie_uid}\")\n    if len(filenames) == 1:\n        return \"3d\", [serie_uid, domain, filenames if (store_filenames is not None) else None]\n    else:\n        return \"2d\", [serie_uid, domain, filenames if (store_filenames is not None) else None]\n\ndef divide_series_by_dicom_dim(series_uid, domains, store_filenames):\n\n    series_with_3d_dicoms = []\n    series_with_2d_dicoms = []\n    \n    series_uids = os.listdir(SERIES_PATH)\n    n_series = len(series_uids)\n\n    with ProcessPoolExecutor(max_workers=multiprocessing.cpu_count()) as executor:\n        futures = []\n        for i in range(n_series):\n            serie_uid = series_uid[i]\n            domain = domains[i] if (domains is not None) else None\n            futures.append(executor.submit(get_serie_info, serie_uid, domain, store_filenames))\n\n        n_jobs = len(futures)\n        for future in tqdm(as_completed(futures), total=n_jobs):\n            dim, serie_info = future.result()\n            if dim == \"3d\":\n                series_with_3d_dicoms.append(serie_info)\n            else:\n                series_with_2d_dicoms.append(serie_info)\n\n    return series_with_3d_dicoms, series_with_2d_dicoms","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-09-26T15:03:06.974365Z","iopub.execute_input":"2025-09-26T15:03:06.974601Z","iopub.status.idle":"2025-09-26T15:03:06.993005Z","shell.execute_reply.started":"2025-09-26T15:03:06.974582Z","shell.execute_reply":"2025-09-26T15:03:06.991719Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"series_with_3d_dicoms, series_with_2d_dicoms = divide_series_by_dicom_dim(series_uid, modalities, True)","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-09-26T15:03:06.994246Z","iopub.execute_input":"2025-09-26T15:03:06.994621Z","iopub.status.idle":"2025-09-26T15:03:39.726778Z","shell.execute_reply.started":"2025-09-26T15:03:06.994589Z","shell.execute_reply":"2025-09-26T15:03:39.725521Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"series_with_3d_dicoms[0], series_with_2d_dicoms[0]","metadata":{"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"print(len(series_with_3d_dicoms), len(series_with_2d_dicoms))","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-09-26T15:03:39.738882Z","iopub.execute_input":"2025-09-26T15:03:39.739208Z","iopub.status.idle":"2025-09-26T15:03:39.76248Z","shell.execute_reply.started":"2025-09-26T15:03:39.739167Z","shell.execute_reply":"2025-09-26T15:03:39.761181Z"}},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"# Get spacing from 3D DICOM","metadata":{}},{"cell_type":"code","source":"def get_dicom3d_spacing(ds):\n\n    pixel_spacing = None\n    slice_thickness = None\n    \n    shared_functional_groups_sequence = getattr(ds, \"SharedFunctionalGroupsSequence\", None)\n    if shared_functional_groups_sequence is not None:\n        pixel_measures_sequence = getattr(shared_functional_groups_sequence[0], \"PixelMeasuresSequence\", None)\n        if pixel_measures_sequence is not None:\n            pixel_spacing = getattr(pixel_measures_sequence[0], \"PixelSpacing\", None)\n            slice_thickness = getattr(pixel_measures_sequence[0], \"SliceThickness\", None)\n    \n    if pixel_spacing is None or slice_thickness is None:\n        raise Exception(\"Missing either Pixel Spacing or Slice Thickness.\")\n    \n    pixel_spacing = [float(axis_spacing) for axis_spacing in pixel_spacing]\n    slice_thickness = float(slice_thickness)\n    spacing = [*pixel_spacing, slice_thickness]\n    \n    return spacing","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-09-26T15:03:39.765727Z","iopub.execute_input":"2025-09-26T15:03:39.766005Z","iopub.status.idle":"2025-09-26T15:03:39.790475Z","shell.execute_reply.started":"2025-09-26T15:03:39.765985Z","shell.execute_reply":"2025-09-26T15:03:39.789429Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"serie_uid_temp = series_with_3d_dicoms[0][0]\ninstance_filename_temp = series_with_3d_dicoms[0][2][0]\nds_temp = pydicom.dcmread(f\"{SERIES_PATH}/{serie_uid_temp}/{instance_filename_temp}\", stop_before_pixels=True)\nprint(get_dicom3d_spacing(ds_temp))","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-09-26T15:03:39.791454Z","iopub.execute_input":"2025-09-26T15:03:39.791769Z","iopub.status.idle":"2025-09-26T15:03:39.842738Z","shell.execute_reply.started":"2025-09-26T15:03:39.791741Z","shell.execute_reply":"2025-09-26T15:03:39.841589Z"}},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"# Get spacing from 2D DICOMs","metadata":{}},{"cell_type":"code","source":"def get_dicom2d_spacing(ds):\n    \n    # Pixel Spacing (x, y)\n    pixel_spacing = getattr(ds, \"PixelSpacing\", None)\n    slice_thickness = getattr(ds, \"SliceThickness\", None)\n\n    if pixel_spacing is None or slice_thickness is None:\n        raise Exception(\"Missing either Pixel Spacing or Slice Thickness.\")\n    \n    pixel_spacing = [float(axis_spacing) for axis_spacing in pixel_spacing]\n    slice_thickness = float(slice_thickness)\n    spacing = [*pixel_spacing, slice_thickness]\n\n    return spacing","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-09-26T15:03:39.843943Z","iopub.execute_input":"2025-09-26T15:03:39.84428Z","iopub.status.idle":"2025-09-26T15:03:39.890211Z","shell.execute_reply.started":"2025-09-26T15:03:39.844243Z","shell.execute_reply":"2025-09-26T15:03:39.889Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"serie_uid_temp = series_with_2d_dicoms[0][0]\ninstance_filename_temp = series_with_2d_dicoms[0][2][0]\nds_temp = pydicom.dcmread(f\"{SERIES_PATH}/{serie_uid_temp}/{instance_filename_temp}\", stop_before_pixels=True)\nprint(get_dicom2d_spacing(ds_temp))","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-09-26T15:03:39.891284Z","iopub.execute_input":"2025-09-26T15:03:39.891667Z","iopub.status.idle":"2025-09-26T15:03:39.948516Z","shell.execute_reply.started":"2025-09-26T15:03:39.891636Z","shell.execute_reply":"2025-09-26T15:03:39.947265Z"}},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"# Get overall spacings","metadata":{}},{"cell_type":"code","source":"def get_serie_spacing(serie_uid, instances_filename_l, n_instances):\n\n    if not instances_filename_l:\n        instances_filename_l = os.listdir(f\"{SERIES_PATH}/{serie_uid}\")\n\n    if not n_instances:\n        n_instances = len(instances_filename_l)\n\n    if n_instances == 1:\n        ds = pydicom.dcmread(f\"{SERIES_PATH}/{serie_uid}/{instances_filename_l[0]}\", stop_before_pixels=True)\n        spacing = np.array(get_dicom3d_spacing(ds))\n    else:\n        spacings = np.zeros((n_instances, 3))\n        for i, instance_filename in enumerate(instances_filename_l):\n            ds = pydicom.dcmread(f\"{SERIES_PATH}/{serie_uid}/{instance_filename}\", stop_before_pixels=True)\n            spacings[i] = get_dicom2d_spacing(ds)\n        spacings, counts = np.unique(spacings, return_counts=True, axis=0)\n        spacing = spacings[np.argmax(counts)]\n\n    return spacing","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-09-26T15:03:39.949589Z","iopub.execute_input":"2025-09-26T15:03:39.949887Z","iopub.status.idle":"2025-09-26T15:03:39.956778Z","shell.execute_reply.started":"2025-09-26T15:03:39.949854Z","shell.execute_reply":"2025-09-26T15:03:39.955727Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"serie_idx = 1\n\n# With pre-fetched filenames\nprint(\"---- Pre-fetched filenames ----\\n\")\n\nserie_uid_temp = series_with_2d_dicoms[serie_idx][0]\ninstances_filename_l_temp = series_with_2d_dicoms[serie_idx][2]\nspacing = get_serie_spacing(serie_uid_temp, instances_filename_l_temp, None)\nprint(f\"2D DICOMS: {spacing}\")\n\nserie_uid_temp = series_with_3d_dicoms[serie_idx][0]\ninstances_filename_l_temp = series_with_3d_dicoms[serie_idx][2]\nspacing = get_serie_spacing(serie_uid_temp, instances_filename_l_temp, 1)\nprint(f\"3D DICOMS: {spacing}\\n\\n\")\n\n# Without pre-fetched filenames\nprint(\"---- No pre-fetched filenames ----\\n\")\n\nserie_uid_temp = series_with_2d_dicoms[serie_idx][0]\nspacing = get_serie_spacing(serie_uid_temp, None, None)\nprint(f\"2D DICOMS: {spacing}\")\n\nserie_uid_temp = series_with_3d_dicoms[serie_idx][0]\nspacing = get_serie_spacing(serie_uid_temp, None, None)\nprint(f\"3D DICOMS: {spacing}\")","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-09-26T15:03:39.957921Z","iopub.execute_input":"2025-09-26T15:03:39.958244Z","iopub.status.idle":"2025-09-26T15:03:42.057174Z","shell.execute_reply.started":"2025-09-26T15:03:39.958216Z","shell.execute_reply":"2025-09-26T15:03:42.056234Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"def get_series_spacing(series_uid, instances_filename_l, n_instances):\n\n    n_series = len(series_uid)\n    \n    if instances_filename_l is None:\n        instances_filename_l = [None] * n_series\n    if n_instances is None:\n        n_instances = [None] * n_series\n\n    assert len(instances_filename_l) == n_series and len(n_instances) == n_series, \"Length of instances_filename_l or n_instances not the same as of series_uid.\"\n\n    spacings = np.zeros((n_series, 3))\n    \n    for i in range(n_series):\n        spacings[i] = get_serie_spacing(series_uid[i], instances_filename_l[i], n_instances[i])\n\n    return spacings","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-09-26T15:03:42.058113Z","iopub.execute_input":"2025-09-26T15:03:42.05864Z","iopub.status.idle":"2025-09-26T15:03:42.065474Z","shell.execute_reply.started":"2025-09-26T15:03:42.058612Z","shell.execute_reply":"2025-09-26T15:03:42.064405Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"series_uid_temp = [serie[0] for serie in series_with_3d_dicoms[0:10]]\ninstances_filename_l_temp = [serie[2] for serie in series_with_3d_dicoms[0:10]]\nspacings = get_series_spacing(series_uid_temp, instances_filename_l_temp, None)\nprint(spacings)","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-09-26T15:03:42.066737Z","iopub.execute_input":"2025-09-26T15:03:42.06701Z","iopub.status.idle":"2025-09-26T15:03:42.31472Z","shell.execute_reply.started":"2025-09-26T15:03:42.066989Z","shell.execute_reply":"2025-09-26T15:03:42.313772Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"def get_domain_spacings(series, has_filenames, batch_size):\n    \n    domain_spacing_dict = make_domain_spacing_dict()\n\n    n_series = len(series)\n    with ProcessPoolExecutor(max_workers=multiprocessing.cpu_count()) as executor:\n        job_to_idx_dict = {}\n        for i in tqdm(range(0, n_series, batch_size)):\n            ixs = slice(i, min(n_series, i+batch_size))\n            series_i = series[ixs]\n            \n            series_uid = [serie[0] for serie in series_i]\n            if has_filenames:\n                instances_filename_l = [serie[2] for serie in series_i]\n            else:\n                instances_filename_l = None\n\n            job = executor.submit(get_series_spacing, series_uid, instances_filename_l, None)\n            job_to_idx_dict[job] = ixs\n\n        n_jobs = len(job_to_idx_dict)\n        for job in tqdm(as_completed(job_to_idx_dict), total=n_jobs):\n            try:\n                spacings = job.result()\n    \n                job_i = job_to_idx_dict[job]\n                for i, serie in enumerate(series[job_i]):\n                    domain = serie[1]\n                    domain_spacing_dict[domain].append(spacings[i])\n            except Exception as e:\n                print(e)\n\n    # Stack lists of numpy arrays\n    for domain, spacings in domain_spacing_dict.items():\n        if len(spacings) > 0:\n            domain_spacing_dict[domain] = np.vstack(spacings)\n        else:\n            domain_spacing_dict[domain] = None\n    \n    return domain_spacing_dict","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-09-26T15:03:42.315666Z","iopub.execute_input":"2025-09-26T15:03:42.316381Z","iopub.status.idle":"2025-09-26T15:03:42.325738Z","shell.execute_reply.started":"2025-09-26T15:03:42.316325Z","shell.execute_reply":"2025-09-26T15:03:42.324645Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"series = series_with_3d_dicoms + series_with_2d_dicoms\ndomain_spacing_dict_temp = get_domain_spacings(series, True, 32)","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-09-26T15:03:42.32681Z","iopub.execute_input":"2025-09-26T15:03:42.327362Z","iopub.status.idle":"2025-09-26T15:38:25.106961Z","shell.execute_reply.started":"2025-09-26T15:03:42.327311Z","shell.execute_reply":"2025-09-26T15:38:25.105446Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"pickle.dump(domain_spacing_dict_temp, open(\"/kaggle/working/domain_spacing_dict.pkl\", \"wb\"))\n#domain_spacing_dict_temp = pickle.load(open(\"/kaggle/working/domain_spacing_dict.pkl\", \"rb\"))","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-09-26T15:38:25.110241Z","iopub.execute_input":"2025-09-26T15:38:25.110593Z","iopub.status.idle":"2025-09-26T15:38:25.117212Z","shell.execute_reply.started":"2025-09-26T15:38:25.110563Z","shell.execute_reply":"2025-09-26T15:38:25.11621Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"for domain, spacings in domain_spacing_dict_temp.items():\n    print(f\"{domain}: {np.median(spacings, axis=0)}\")","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-09-26T15:38:25.118784Z","iopub.execute_input":"2025-09-26T15:38:25.119151Z","iopub.status.idle":"2025-09-26T15:38:25.145311Z","shell.execute_reply.started":"2025-09-26T15:38:25.119118Z","shell.execute_reply":"2025-09-26T15:38:25.144175Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"","metadata":{"trusted":true},"outputs":[],"execution_count":null}]}