{"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_minor":4,"nbformat":4,"cells":[{"cell_type":"markdown","source":"<div class='alert alert-info' style='text-align: center'><h1>Exploring rib suppression on CXR</h1>\n    - yet another chest x-ray processing notebook -\n</div>\n\n#### It occured to me that I should try to create a rib suppression tool. I've not seen any open source software for it, although there are proprietarty systems with filters such as this.\n\n- The concept is to create masks that can be used to filter out ribs to enable better 'visualization' of underlying lung tissue. \n- Tissue density under the ribs varies greatly depending on underlying anatomy, so we'll use multiple masks to deal with this.\n- In this case, I used two masks to demonstrate the thought process.\n- This notebook relies on a dataset of 2 binary masks I made by hand specifically for the DICOM image I'm using.","metadata":{}},{"cell_type":"code","source":"import numpy as np\nimport pydicom\nimport cv2\nimport sys\nimport matplotlib.pyplot as plt","metadata":{"execution":{"iopub.status.busy":"2021-07-06T13:46:07.902878Z","iopub.execute_input":"2021-07-06T13:46:07.903302Z","iopub.status.idle":"2021-07-06T13:46:08.320780Z","shell.execute_reply.started":"2021-07-06T13:46:07.903213Z","shell.execute_reply":"2021-07-06T13:46:08.319788Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"- First we'll grab a random DICOM image from the SIIM Covid-19 Detection set, and squeeze the pixels to a range of 0-255.","metadata":{}},{"cell_type":"code","source":"# Load the image\nimage = pydicom.dcmread('../input/siim-covid19-detection/train/00a76543ed93/4a223cccbe04/ad8d4a5ba8f0.dcm')\n\n# Normalize it to 8 bit\nimage = cv2.normalize(image.pixel_array, None, 0,255, cv2.NORM_MINMAX, cv2.CV_8U)\nplt.figure(figsize= (10,5))\nplt.imshow(image, cmap='gray');","metadata":{"execution":{"iopub.status.busy":"2021-07-06T13:46:08.322519Z","iopub.execute_input":"2021-07-06T13:46:08.322942Z","iopub.status.idle":"2021-07-06T13:46:10.957992Z","shell.execute_reply.started":"2021-07-06T13:46:08.322876Z","shell.execute_reply":"2021-07-06T13:46:10.954093Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"- I made two binary masks using GIMP. They could be created via segmentation or other means (edge detection, or bit-plane slicing maybe?)\n- One of the masks is for the ribs, the other represents the lung fields without the mediastinum.\n- From these two masks, we'll create a third mask .. which is the 'difference' of the other two.","metadata":{}},{"cell_type":"code","source":"# Load the masks\nmask_rib = cv2.imread('../input/cxr-mask/ad8d4a5ba8f0_mask.jpg', 0)\nret,mask_rib = cv2.threshold(mask_rib,127,255,cv2.THRESH_BINARY)\n\nmask_lung = cv2.imread('../input/cxr-mask/ad8d4a5ba8f0_mask_2.jpg', 0)\nret2,mask_lung = cv2.threshold(mask_lung,127,255,cv2.THRESH_BINARY)\n\n# Create the third mask from the first two using CV2 BitwiseOR\nmask_diff = cv2.bitwise_or(mask_rib, mask_rib, mask = mask_lung)\n\n# Make inverted versions of the masks\nmask_rib_inv = cv2.bitwise_not(mask_rib)\nmask_lung_inv = cv2.bitwise_not(mask_lung)\nmask_diff_inv = cv2.bitwise_not(mask_diff)","metadata":{"execution":{"iopub.status.busy":"2021-07-06T13:46:10.959917Z","iopub.execute_input":"2021-07-06T13:46:10.960247Z","iopub.status.idle":"2021-07-06T13:46:11.404009Z","shell.execute_reply.started":"2021-07-06T13:46:10.960217Z","shell.execute_reply":"2021-07-06T13:46:11.403058Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"plt.figure(figsize= (20,10))\nplt.subplot(1, 3, 1)\nplt.imshow(mask_rib, cmap='gray');\nplt.subplot(1, 3, 2)\nplt.imshow(mask_lung, cmap='gray');\nplt.subplot(1, 3, 3)\nplt.imshow(mask_diff, cmap='gray');","metadata":{"execution":{"iopub.status.busy":"2021-07-06T13:46:11.405207Z","iopub.execute_input":"2021-07-06T13:46:11.405452Z","iopub.status.idle":"2021-07-06T13:46:15.855840Z","shell.execute_reply.started":"2021-07-06T13:46:11.405428Z","shell.execute_reply":"2021-07-06T13:46:15.854962Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"- Now we'll get the pixels from the original image using the rib mask and darken them, then put them back onto the image using the inverse mask.","metadata":{}},{"cell_type":"code","source":"# Get pixels from the original image using the first mask\npix = cv2.bitwise_or(image, image, mask = mask_rib)\n\n# Darken them slightly\npix = pix * .95\npix = np.where(pix < 1, 0, pix)\n \n# Make an inverse mask of the original image using the inverted rib mask\ninv = cv2.bitwise_or(image, image, mask = mask_rib_inv)\n\n# Add the two images together\noutput = pix + inv","metadata":{"execution":{"iopub.status.busy":"2021-07-06T13:46:15.856855Z","iopub.execute_input":"2021-07-06T13:46:15.857197Z","iopub.status.idle":"2021-07-06T13:46:16.020466Z","shell.execute_reply.started":"2021-07-06T13:46:15.857167Z","shell.execute_reply":"2021-07-06T13:46:16.019468Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Display the image steps\nplt.figure(figsize= (20,5))\nplt.subplot(1, 3, 1)\nplt.imshow(pix, cmap='gray');\nplt.subplot(1, 3, 2)\nplt.imshow(inv, cmap='gray');\nplt.subplot(1, 3, 3)\nplt.imshow(output, cmap='gray');","metadata":{"execution":{"iopub.status.busy":"2021-07-06T13:46:16.021538Z","iopub.execute_input":"2021-07-06T13:46:16.021801Z","iopub.status.idle":"2021-07-06T13:46:20.301799Z","shell.execute_reply.started":"2021-07-06T13:46:16.021774Z","shell.execute_reply":"2021-07-06T13:46:20.300882Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"- I only darkened the masked pixels slightly, but you can see that the ribs over the heart are slightly less visible. (hard to see in this small image).\n- Next, we'll make a mask of the masked pixels from the last step, and darken them a little bit.\n- This will darken only the pixels in the lung area and exclude the heart and mediastinum.","metadata":{}},{"cell_type":"code","source":"# Get the pixels from the first mapped set using the second mask\npix2 = cv2.bitwise_or(pix, pix, mask = mask_lung)\n\n# Darken them a bit, I randomly picked this value for demonstration purposes\npix2 = pix2 * .82\npix2 = np.where(pix2 < 1, 0, pix2)\n \n# Make an inverse mask of the output of the first masking operation\ninv2 = cv2.bitwise_or(output, output, mask = mask_diff_inv)\n\n# Add the third mask mapped pixels and the inverse mask back together\noutput2 = pix2 + inv2","metadata":{"execution":{"iopub.status.busy":"2021-07-06T13:46:20.302946Z","iopub.execute_input":"2021-07-06T13:46:20.303333Z","iopub.status.idle":"2021-07-06T13:46:20.577649Z","shell.execute_reply.started":"2021-07-06T13:46:20.303302Z","shell.execute_reply":"2021-07-06T13:46:20.576693Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Display the steps\nplt.figure(figsize= (20,5))\nplt.subplot(1, 3, 1)\nplt.imshow(pix2, cmap='gray');\nplt.subplot(1, 3, 2)\nplt.imshow(inv2, cmap='gray');\nplt.subplot(1, 3, 3)\nplt.imshow(output2, cmap='gray');","metadata":{"execution":{"iopub.status.busy":"2021-07-06T13:46:20.579634Z","iopub.execute_input":"2021-07-06T13:46:20.580190Z","iopub.status.idle":"2021-07-06T13:46:24.759134Z","shell.execute_reply.started":"2021-07-06T13:46:20.580145Z","shell.execute_reply":"2021-07-06T13:46:24.758244Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Display the original image and the final image\nplt.figure(figsize= (20,10))\nplt.subplot(1, 2, 1)\nplt.imshow(image, cmap='gray');\nplt.subplot(1, 2, 2)\nplt.imshow(output2, cmap='gray');","metadata":{"execution":{"iopub.status.busy":"2021-07-06T13:46:24.760330Z","iopub.execute_input":"2021-07-06T13:46:24.760579Z","iopub.status.idle":"2021-07-06T13:46:27.713537Z","shell.execute_reply.started":"2021-07-06T13:46:24.760553Z","shell.execute_reply":"2021-07-06T13:46:27.712908Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"#### Conclusion\n- Of course this is a quick idea and my results aren't very good but, you can see in the final image on the right that the ribs are much less prominent than they are in the original image.\n- The soft tissue of the lung has a little less noise in it and the underlying pixel values haven't been disturbed.\n- With a few more masks and some filtering (maybe edge detection for the denser periosteum on the rib edges?), this could be useful.","metadata":{}},{"cell_type":"markdown","source":"**Here are some other processing notebooks I made:**\n- Lung segmentation without CNN -> https://www.kaggle.com/davidbroberts/lung-segmentation-without-cnn\n- Applying filters to x-rays -> https://www.kaggle.com/davidbroberts/applying-filters-to-chest-x-rays\n- Manual DICOM VOI LUT -> https://www.kaggle.com/davidbroberts/manual-dicom-voi-lut\n- Apply Unsharp Mask to Chest X-Rays -> https://www.kaggle.com/davidbroberts/unsharp-masking-chest-x-rays\n- Cropping Chest X-Rays -> https://www.kaggle.com/davidbroberts/cropping-chest-x-rays\n- Bounding Boxes on Cropped Images -> https://www.kaggle.com/davidbroberts/bounding-boxes-on-cropped-images\n- Visualizing Chest X-Ray bit planes -> https://www.kaggle.com/davidbroberts/visualizing-chest-x-ray-bitplanes\n- DICOM full range pixels as CNN input -> https://www.kaggle.com/davidbroberts/dicom-full-range-pixels-as-cnn-input\n- Standardizing Chest X-Ray Dataset Exports -> https://www.kaggle.com/davidbroberts/standardizing-cxr-datasets","metadata":{}}]}