{"cells":[{"metadata":{},"cell_type":"markdown","source":"# Let's See How to Select Sequences and Sort DICOM data \n\n\n### After inspecting the data in the RSNA STR Pulmonary Embolism Detection challenge I occasionally found studies with two sequences with the same series number. This simple notebook shows an example of how to handle this edge case and uses StudyInstanceUID 759a5963508b as an example input.\n\n### For 3D approachs to CTA pulmonary embolus prediction an accurate and consistent sorting mechanism and proper sequence selection is essential."},{"metadata":{},"cell_type":"markdown","source":"# Table of Contents\n\n### 1. [Load a study](#load)\n### 2. [Utility functions for dicom display](#utility)\n### 3. [Isolate dicom series](#series)\n### 4. [Display all series unsorted](#unsorted)\n### 5. [Display all series sorted](#sorted)"},{"metadata":{"_uuid":"8f2839f25d086af736a60e9eeb907d3b93b6e0e5","_cell_guid":"b1076dfc-b9ad-4769-8c92-a6c4dae69d19","trusted":true},"cell_type":"code","source":"##### PACKAGES\nimport os\nfrom pathlib import Path\nimport pydicom as dcm  # great library for working with dicom\nfrom collections import defaultdict\nimport matplotlib.pyplot as plt","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"# Load a dicom study <a class=\"anchor\" id=\"load\"></a>"},{"metadata":{"_uuid":"d629ff2d2480ee46fbb7e2d37f6b5fab8052498a","_cell_guid":"79c7e3d0-c299-4dcb-8224-4455121ee9b0","trusted":true},"cell_type":"code","source":"data_path = Path('../input/rsna-str-pulmonary-embolism-detection/train/759a5963508b')\ndcm_paths = list(data_path.glob('**/*.dcm'))\n\nprint(f\"First 10 dicom paths: \\n\\n{dcm_paths[:10]}\\n\\n... out of {len(dcm_paths)} total dcms\")","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"# Utility functions for display <a class=\"anchor\" id=\"utility\"></a>"},{"metadata":{"trusted":true},"cell_type":"code","source":"\"\"\"Utiliy functions for displaying dcms clearly\"\"\"\n\n\ndef get_first_of_dicom_field_as_int(x):\n    #get x[0] as in int is x is a 'pydicom.multival.MultiValue', otherwise get int(x)\n    if type(x) == dcm.multival.MultiValue: return int(x[0])\n    else: return int(x)\n    \ndef get_windowing(data):\n    dicom_fields = [data[('0028','1050')].value, #window center\n                    data[('0028','1051')].value, #window width\n                    data[('0028','1052')].value, #intercept\n                    data[('0028','1053')].value] #slope\n    return [get_first_of_dicom_field_as_int(x) for x in dicom_fields]\n\ndef metadata_window(img, print_ranges=False):\n    # Get data from dcm\n    window_center, window_width, intercept, slope = get_windowing(img)\n    img = img.pixel_array\n    \n    # Window based on dcm metadata\n    img = img * slope + intercept\n    img_min = window_center - window_width // 2\n    img_max = window_center + window_width // 2\n    if print_ranges:\n        print(f'Window Level: {window_center}, Window Width:{window_width}, Range:[{img_min} {img_max}]')\n    img[img < img_min] = img_min\n    img[img > img_max] = img_max\n    \n    # Normalize\n    img = (img - img_min) / (img_max - img_min)\n    return img\n\ndef show_image_slices(dcms, title=\"dcms\"):\n    fig, axes = plt.subplots(1, len(dcms), figsize=(25,5))\n    fig.suptitle(title, fontsize=16)\n    for i, dcm in enumerate(dcms):\n        axes[i].imshow(metadata_window(dcm), cmap=\"gray\")\n        ","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"# Isolate a dicom series <a class=\"anchor\" id=\"series\"></a>"},{"metadata":{"trusted":true},"cell_type":"code","source":"\"\"\" \nNote: slices are most reliably sorted by patient position. For axial slices this will be z index\nas there is no gaurantee that the sequence numbers are different in a study so we can generate \na unique key for each sequence (series number + position-axis[0] + position-axis[1]).\nThe first two positions (x, y) in ImagePositionPatient usually vary between Axial sequences \n\"\"\"\n\n# dictionary of keys to series\ndcm_series = defaultdict(list)\nfor dcm_path in dcm_paths:\n    d = dcm.read_file(dcm_path) \n    # key = series number + first axis, most likely unique across axial series.\n    dcm_series[f\"Series #:{d.SeriesNumber}, x:{d.ImagePositionPatient[0]}, y:{d.ImagePositionPatient[1]}\"].append(d)\n    \n    \nprint(\"printing series by series-keys with lengths...\")\nfor key, series in dcm_series.items():\n    print(f\"series key [{key}] has:  {len(series)} dcms\")\n","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"# Display series data unsorted <a class=\"anchor\" id=\"unsorted\"></a>"},{"metadata":{"trusted":true},"cell_type":"code","source":"for key, series in dcm_series.items():\n    \n    show_image_slices(series[:30:5], title=key+\" UNSORTED\")\n    ","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"# Display series data sorted <a class=\"anchor\" id=\"sorted\"></a>"},{"metadata":{"trusted":true},"cell_type":"code","source":"for key, series in dcm_series.items():\n    \n    show_image_slices(sorted(series[:30:5], key=lambda dcm:dcm.ImagePositionPatient[2]), title=key+\" SORTED\")   # ImagePositionPatient[2] is the Z axis\n    ","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"## Thanks for checking out this notebook. I hope some of you find it helpful.\n\n## Have fun working with dicoms and medical imaging!\n"}],"metadata":{"kernelspec":{"language":"python","display_name":"Python 3","name":"python3"},"language_info":{"pygments_lexer":"ipython3","nbconvert_exporter":"python","version":"3.6.4","file_extension":".py","codemirror_mode":{"name":"ipython","version":3},"name":"python","mimetype":"text/x-python"}},"nbformat":4,"nbformat_minor":4}