{"metadata":{"kernelspec":{"language":"python","display_name":"Python 3","name":"python3"},"language_info":{"name":"python","version":"3.10.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":39272,"databundleVersionId":4629629,"sourceType":"competition"}],"dockerImageVersionId":30746,"isInternetEnabled":true,"language":"python","sourceType":"notebook","isGpuEnabled":false}},"nbformat_minor":4,"nbformat":4,"cells":[{"cell_type":"code","source":"# This Python 3 environment comes with many helpful analytics libraries installed\n# It is defined by the kaggle/python Docker image: https://github.com/kaggle/docker-python\n# For example, here's several helpful packages to load\n\n# import numpy as np # linear algebra\n# import pandas as pd # data processing, CSV file I/O (e.g. pd.read_csv)\n\n# # Input data files are available in the read-only \"../input/\" directory\n# # For example, running this (by clicking run or pressing Shift+Enter) will list all files under the input directory\n\n# import os\n# for dirname, _, filenames in os.walk('/kaggle/input'):\n#     for filename in filenames:\n#         print(os.path.join(dirname, filename))\n\n# You can write up to 20GB to the current directory (/kaggle/working/) that gets preserved as output when you create a version using \"Save & Run All\" \n# You can also write temporary files to /kaggle/temp/, but they won't be saved outside of the current session","metadata":{"_uuid":"8f2839f25d086af736a60e9eeb907d3b93b6e0e5","_cell_guid":"b1076dfc-b9ad-4769-8c92-a6c4dae69d19","execution":{"iopub.status.busy":"2024-08-02T18:41:16.633032Z","iopub.execute_input":"2024-08-02T18:41:16.633494Z","iopub.status.idle":"2024-08-02T18:41:16.639425Z","shell.execute_reply.started":"2024-08-02T18:41:16.633459Z","shell.execute_reply":"2024-08-02T18:41:16.638056Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# !pip install pydicom","metadata":{"execution":{"iopub.status.busy":"2024-08-02T18:41:16.641873Z","iopub.execute_input":"2024-08-02T18:41:16.64262Z","iopub.status.idle":"2024-08-02T18:41:16.656703Z","shell.execute_reply.started":"2024-08-02T18:41:16.642559Z","shell.execute_reply":"2024-08-02T18:41:16.655529Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"import cv2\nimport pydicom\nfrom pydicom.data import get_testdata_file\nfrom gzip import GzipFile\nimport matplotlib.pyplot as plt\nimport pandas as pd","metadata":{"execution":{"iopub.status.busy":"2024-08-03T16:10:09.837797Z","iopub.execute_input":"2024-08-03T16:10:09.838228Z","iopub.status.idle":"2024-08-03T16:10:10.890839Z","shell.execute_reply.started":"2024-08-03T16:10:09.838196Z","shell.execute_reply":"2024-08-03T16:10:10.889337Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Đọc file CSV\nfile_path = '/kaggle/input/rsna-breast-cancer-detection/train.csv'\ndf = pd.read_csv(file_path)\n\n# Lấy ra 3 dòng có dữ liệu là 1 trong cột 'cancer'\ncancer_1 = df[df['cancer'] == 1].head(10)\n\n# Lấy ra 3 dòng có dữ liệu là 0 trong cột 'cancer'\ncancer_0 = df[df['cancer'] == 0].head(10)\n# print(cancer_1)\n# print(cancer_0)\n# Kết hợp lại hai DataFrame\nresult = pd.concat([cancer_1, cancer_0])\n\n# In kết quả\nprint(result)\n\n# Nếu muốn lưu kết quả vào file CSV mới\n# result.to_csv('output.csv', index=False)","metadata":{"execution":{"iopub.status.busy":"2024-08-03T16:10:14.717835Z","iopub.execute_input":"2024-08-03T16:10:14.719483Z","iopub.status.idle":"2024-08-03T16:10:14.913974Z","shell.execute_reply.started":"2024-08-03T16:10:14.719433Z","shell.execute_reply":"2024-08-03T16:10:14.912451Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# !pip install pylibjpeg pylibjpeg-libjpeg pylibjpeg-openjpeg","metadata":{"execution":{"iopub.status.busy":"2024-08-02T18:41:16.76378Z","iopub.execute_input":"2024-08-02T18:41:16.7642Z","iopub.status.idle":"2024-08-02T18:41:16.768921Z","shell.execute_reply.started":"2024-08-02T18:41:16.764168Z","shell.execute_reply":"2024-08-02T18:41:16.767669Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"cancerImg = f\"/kaggle/input/rsna-breast-cancer-detection/train_images/{result.iloc[2]['patient_id']}/{result.iloc[2]['image_id']}.dcm\"\nnoncancerImg = f\"/kaggle/input/rsna-breast-cancer-detection/train_images/{result.iloc[15]['patient_id']}/{result.iloc[15]['image_id']}.dcm\"\n\n# dicomCancer = pydicom.dcmread(cancerImg)\n# dicomNoneCancer = pydicom.dcmread(noncancerImg)\n# In đường dẫn\n# print(noncancerImg)\n# print(cancerImg)","metadata":{"execution":{"iopub.status.busy":"2024-08-03T16:12:03.912358Z","iopub.execute_input":"2024-08-03T16:12:03.912842Z","iopub.status.idle":"2024-08-03T16:12:03.920963Z","shell.execute_reply.started":"2024-08-03T16:12:03.912805Z","shell.execute_reply":"2024-08-03T16:12:03.919354Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# !pip install SimpleITK","metadata":{"execution":{"iopub.status.busy":"2024-08-02T18:41:16.784806Z","iopub.execute_input":"2024-08-02T18:41:16.785378Z","iopub.status.idle":"2024-08-02T18:41:16.798874Z","shell.execute_reply.started":"2024-08-02T18:41:16.785346Z","shell.execute_reply":"2024-08-02T18:41:16.797657Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"import SimpleITK as sitk\nimport matplotlib.pyplot as plt\n\n\n# Đọc file DICOM\nimage1 = sitk.ReadImage(cancerImg)\nimage_arrayCancer = sitk.GetArrayFromImage(image1)[0]\n\nimage2 = sitk.ReadImage(noncancerImg)\nimage_arrayNonCancer = sitk.GetArrayFromImage(image2)[0]\n\n# Tạo hình ảnh và hiển thị chúng trên cùng một hàng\nfig, axes = plt.subplots(1, 2, figsize=(10, 5))\n\n# Hiển thị hình ảnh thứ nhất\naxes[0].imshow(image_arrayCancer, cmap='gray')\naxes[0].set_title('Cancer')\naxes[0].axis('off')\n\n# Hiển thị hình ảnh thứ hai\naxes[1].imshow(image_arrayNonCancer, cmap='gray')\naxes[1].set_title('No Cancer')\naxes[1].axis('off')\n\n# Hiển thị hình ảnh\nplt.tight_layout()\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2024-08-03T16:12:05.746653Z","iopub.execute_input":"2024-08-03T16:12:05.747157Z","iopub.status.idle":"2024-08-03T16:12:08.302972Z","shell.execute_reply.started":"2024-08-03T16:12:05.747115Z","shell.execute_reply":"2024-08-03T16:12:08.301485Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# 1. Meidan Filtering","metadata":{}},{"cell_type":"code","source":"import cv2\nimport matplotlib.pyplot as plt\n\n# Apply median filtering\nkernel_size = 3  # Adjust the kernel size as needed\nfiltered_image = cv2.medianBlur(image_arrayCancer, kernel_size)\nfiltered_image2 = cv2.medianBlur(image_arrayNonCancer, kernel_size)\n\n# Display the original and filtered images\nfig, ax = plt.subplots(1, 4, figsize=(12, 9))\n\n# Display original image YES\nax[0].imshow(image_arrayCancer, cmap='gray')\nax[0].set_title('Original Image YES')\nax[0].axis('off')\n\n# Display filtered image YES\nax[1].imshow(filtered_image, cmap='gray')\nax[1].set_title('Median Filtered YES')\nax[1].axis('off')\n\n# Display original image NO\nax[2].imshow(image_arrayNonCancer, cmap='gray')\nax[2].set_title('Original Image NO')\nax[2].axis('off')\n\n# Display filtered image NO\nax[3].imshow(filtered_image2, cmap='gray')\nax[3].set_title('Median Filtered NO')\nax[3].axis('off')\n\n# Show images\nplt.tight_layout()\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2024-08-03T16:12:49.906913Z","iopub.execute_input":"2024-08-03T16:12:49.907463Z","iopub.status.idle":"2024-08-03T16:12:52.73558Z","shell.execute_reply.started":"2024-08-03T16:12:49.907426Z","shell.execute_reply":"2024-08-03T16:12:52.734201Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# 2. Radiopaque artifacts suppression","metadata":{}},{"cell_type":"markdown","source":"## Binarization","metadata":{}},{"cell_type":"code","source":"import numpy as np\n\nif image_arrayCancer.dtype != np.uint8:\n    image_arrayCancer = cv2.normalize(image_arrayCancer, None, 0, 255, cv2.NORM_MINMAX).astype(np.uint8)\n\n# Apply binarization (Global Thresholding)\n_, binary_imageCancer = cv2.threshold(image_arrayCancer, 140, 255, cv2.THRESH_BINARY)\n\nif image_arrayNonCancer.dtype != np.uint8:\n    image_arrayNonCancer = cv2.normalize(image_arrayNonCancer, None, 0, 255, cv2.NORM_MINMAX).astype(np.uint8)\n\n# Apply binarization (Global Thresholding)\n_, binary_imageNonCancer = cv2.threshold(image_arrayNonCancer, 140, 255, cv2.THRESH_BINARY)\n\n# Alternatively, use Otsu's Thresholding for automatic threshold determination\n# _, binary_image = cv2.threshold(pixel_array, 0, 255, cv2.THRESH_BINARY + cv2.THRESH_OTSU)\n\n# Display the original and binarized images\nfig, ax = plt.subplots(1, 4, figsize=(20, 6))\n\n# Display original image\nax[0].imshow(image_arrayCancer, cmap='gray')\nax[0].set_title('Original Image YES')\nax[0].axis('off')\n\n# Display binarized image\nax[1].imshow(binary_imageCancer, cmap='gray')\nax[1].set_title('Binarized Image YES')\nax[1].axis('off')\n\n# Display original image\nax[2].imshow(image_arrayNonCancer, cmap='gray')\nax[2].set_title('Original Image NO')\nax[2].axis('off')\n\n# Display binarized image\nax[3].imshow(binary_imageNonCancer, cmap='gray')\nax[3].set_title('Binarized Image NO')\nax[3].axis('off')\n\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2024-08-03T16:13:38.787036Z","iopub.execute_input":"2024-08-03T16:13:38.787553Z","iopub.status.idle":"2024-08-03T16:13:41.750737Z","shell.execute_reply.started":"2024-08-03T16:13:38.787516Z","shell.execute_reply":"2024-08-03T16:13:41.749175Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## Artifacts suppression","metadata":{}},{"cell_type":"code","source":"import cv2\nimport numpy as np\nimport matplotlib.pyplot as plt\n\ndef select_largest_obj(img_bin, lab_val=255, fill_holes=False, \n                       smooth_boundary=False, kernel_size=15):\n    '''Select the largest object from a binary image and optionally\n    fill holes inside it and smooth its boundary.\n    Args:\n        img_bin(2D array): 2D numpy array of binary image.\n        lab_val([int]): integer value used for the label of the largest \n                        object. Default is 255.\n        fill_holes([boolean]): whether fill the holes inside the largest \n                               object or not. Default is false.\n        smooth_boundary([boolean]): whether smooth the boundary of the \n                                    largest object using morphological \n                                    opening or not. Default is false.\n        kernel_size([int]): the size of the kernel used for morphological \n                            operation.\n    '''\n    n_labels, img_labeled, lab_stats, _ = cv2.connectedComponentsWithStats(\n        img_bin, connectivity=8, ltype=cv2.CV_32S)\n    largest_obj_lab = np.argmax(lab_stats[1:, 4]) + 1\n    largest_mask = np.zeros(img_bin.shape, dtype=np.uint8)\n    largest_mask[img_labeled == largest_obj_lab] = lab_val\n    if fill_holes:\n        bkg_locs = np.where(img_labeled == 0)\n        bkg_seed = (bkg_locs[0][0], bkg_locs[1][0])\n        img_floodfill = largest_mask.copy()\n        h_, w_ = largest_mask.shape\n        mask_ = np.zeros((h_ + 2, w_ + 2), dtype=np.uint8)\n        cv2.floodFill(img_floodfill, mask_, seedPoint=bkg_seed, newVal=lab_val)\n        holes_mask = cv2.bitwise_not(img_floodfill)  # mask of the holes.\n        largest_mask = largest_mask + holes_mask\n    if smooth_boundary:\n        kernel_ = np.ones((kernel_size, kernel_size), dtype=np.uint8)\n        largest_mask = cv2.morphologyEx(largest_mask, cv2.MORPH_OPEN, kernel_)\n        \n    return largest_mask\n\ndef process_image(image_array, threshold=128):\n    '''Process an image array to apply thresholding, select the largest object, \n       and return the results.\n    Args:\n        image_array(2D array): The image array to process.\n        threshold(int): The threshold value for binarization.\n    '''\n    if image_array.dtype != np.uint8:\n        image_array = cv2.normalize(image_array, None, 0, 255, cv2.NORM_MINMAX).astype(np.uint8)\n    _, binary_image = cv2.threshold(image_array, threshold, 255, cv2.THRESH_BINARY)\n    largest_object_mask = select_largest_obj(binary_image, lab_val=255, fill_holes=True, smooth_boundary=True, kernel_size=15)\n    result_image = cv2.bitwise_and(image_array, largest_object_mask)\n    return image_array, binary_image, largest_object_mask, result_image\n\n# Assume image_arrayCancer and image_arrayNonCancer are loaded as numpy arrays\n\n# Process both images\nimage_arrayCancer, binary_imageCancer, largest_object_maskCancer, result_imageCancer = process_image(image_arrayCancer)\nimage_arrayNonCancer, binary_imageNonCancer, largest_object_maskNonCancer, result_imageNonCancer = process_image(image_arrayNonCancer)\n\n# Display the images in a single row\nfig, axes = plt.subplots(2, 4, figsize=(20, 12))\n\n# Display Original Image Cancer\naxes[0, 0].imshow(image_arrayCancer, cmap='gray')\naxes[0, 0].set_title('Original Image Cancer')\naxes[0, 0].axis('off')\n\n# Display Binarized Image Cancer\naxes[0, 1].imshow(binary_imageCancer, cmap='gray')\naxes[0, 1].set_title('Binarized Image Cancer')\naxes[0, 1].axis('off')\n\n# Display Largest Object Mask Cancer\naxes[0, 2].imshow(largest_object_maskCancer, cmap='gray')\naxes[0, 2].set_title('Largest Object Mask Cancer')\naxes[0, 2].axis('off')\n\n# Display Artifact Suppressed Image Cancer\naxes[0, 3].imshow(result_imageCancer, cmap='gray')\naxes[0, 3].set_title('Artifact Suppressed Image Cancer')\naxes[0, 3].axis('off')\n\n# Display Original Image NonCancer\naxes[1, 0].imshow(image_arrayNonCancer, cmap='gray')\naxes[1, 0].set_title('Original Image NonCancer')\naxes[1, 0].axis('off')\n\n# Display Binarized Image NonCancer\naxes[1, 1].imshow(binary_imageNonCancer, cmap='gray')\naxes[1, 1].set_title('Binarized Image NonCancer')\naxes[1, 1].axis('off')\n\n# Display Largest Object Mask NonCancer\naxes[1, 2].imshow(largest_object_maskNonCancer, cmap='gray')\naxes[1, 2].set_title('Largest Object Mask NonCancer')\naxes[1, 2].axis('off')\n\n# Display Artifact Suppressed Image NonCancer\naxes[1, 3].imshow(result_imageNonCancer, cmap='gray')\naxes[1, 3].set_title('Artifact Suppressed Image NonCancer')\naxes[1, 3].axis('off')\n\nplt.tight_layout()\nplt.show()\n","metadata":{"execution":{"iopub.status.busy":"2024-08-02T18:41:32.871349Z","iopub.execute_input":"2024-08-02T18:41:32.871703Z","iopub.status.idle":"2024-08-02T18:41:44.286384Z","shell.execute_reply.started":"2024-08-02T18:41:32.871672Z","shell.execute_reply":"2024-08-02T18:41:44.285249Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# 3. Pectoral muscle suppression","metadata":{}},{"cell_type":"markdown","source":"## Contrast enhancement","metadata":{}},{"cell_type":"code","source":"import cv2\nimport numpy as np\nimport matplotlib.pyplot as plt\n\ndef contrast_stretching(img):\n    min_intensity = np.min(img)\n    max_intensity = np.max(img)\n    stretched = cv2.normalize(img, None, 0, 255, cv2.NORM_MINMAX)\n    return stretched\n\ndef enhance_image(pixel_array):\n    if pixel_array.dtype != np.uint8:\n        pixel_array = cv2.normalize(pixel_array, None, 0, 255, cv2.NORM_MINMAX).astype(np.uint8)\n\n    # Apply CLAHE\n    clahe = cv2.createCLAHE(clipLimit=2.0, tileGridSize=(8, 8))\n    clahe_enhanced = clahe.apply(pixel_array)\n    \n    # Apply linear contrast stretching\n    contrast_stretched = contrast_stretching(pixel_array)\n    \n    return clahe_enhanced, contrast_stretched\n\n# Assume image_arrayCancer and image_arrayNonCancer are loaded as numpy arrays\n\n# Enhance both images\nclahe_enhanced_Cancer, contrast_stretched_Cancer = enhance_image(image_arrayCancer)\nclahe_enhanced_NonCancer, contrast_stretched_NonCancer = enhance_image(image_arrayNonCancer)\n\n# Display the images\nfig, axes = plt.subplots(2, 3, figsize=(20, 12))\n\n# Display Original Image Cancer\naxes[0, 0].imshow(image_arrayCancer, cmap='gray')\naxes[0, 0].set_title('Original Image Cancer')\naxes[0, 0].axis('off')\n\n# Display CLAHE Enhanced Image Cancer\naxes[0, 1].imshow(clahe_enhanced_Cancer, cmap='gray')\naxes[0, 1].set_title('CLAHE Enhanced Cancer')\naxes[0, 1].axis('off')\n\n# Display Contrast Stretched Image Cancer\naxes[0, 2].imshow(contrast_stretched_Cancer, cmap='gray')\naxes[0, 2].set_title('Contrast Stretched Cancer')\naxes[0, 2].axis('off')\n\n# Display Original Image NonCancer\naxes[1, 0].imshow(image_arrayNonCancer, cmap='gray')\naxes[1, 0].set_title('Original Image NonCancer')\naxes[1, 0].axis('off')\n\n# Display CLAHE Enhanced Image NonCancer\naxes[1, 1].imshow(clahe_enhanced_NonCancer, cmap='gray')\naxes[1, 1].set_title('CLAHE Enhanced NonCancer')\naxes[1, 1].axis('off')\n\n# Display Contrast Stretched Image NonCancer\naxes[1, 2].imshow(contrast_stretched_NonCancer, cmap='gray')\naxes[1, 2].set_title('Contrast Stretched NonCancer')\naxes[1, 2].axis('off')\n\nplt.tight_layout()\nplt.show()\n","metadata":{"execution":{"iopub.status.busy":"2024-08-02T18:41:44.287782Z","iopub.execute_input":"2024-08-02T18:41:44.288091Z","iopub.status.idle":"2024-08-02T18:41:52.412835Z","shell.execute_reply.started":"2024-08-02T18:41:44.288063Z","shell.execute_reply":"2024-08-02T18:41:52.41177Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## Seeded flood fill (not work)","metadata":{}},{"cell_type":"code","source":"import cv2\nimport numpy as np\nimport matplotlib.pyplot as plt\n\ndef convert_to_8bit(pixel_array):\n    if pixel_array.dtype != np.uint8:\n        return cv2.normalize(pixel_array, None, 0, 255, cv2.NORM_MINMAX).astype(np.uint8)\n    return pixel_array\n\ndef apply_seeded_flood_fill(binary_image, seed_point, new_value, lo_diff=10, up_diff=10):\n    filled_image = binary_image.copy()\n    h, w = filled_image.shape\n    mask = np.zeros((h + 2, w + 2), dtype=np.uint8)\n    cv2.floodFill(filled_image, mask, seedPoint=seed_point, newVal=new_value,\n                  loDiff=(lo_diff, lo_diff, lo_diff), upDiff=(up_diff, up_diff, up_diff))\n    return filled_image\n\n# Convert both images to 8-bit if necessary\nimage_arrayCancer_8bit = convert_to_8bit(image_arrayCancer)\nimage_arrayNonCancer_8bit = convert_to_8bit(image_arrayNonCancer)\n\n# Convert images to binary\n_, binary_imageCancer = cv2.threshold(image_arrayCancer_8bit, 128, 255, cv2.THRESH_BINARY)\n_, binary_imageNonCancer = cv2.threshold(image_arrayNonCancer_8bit, 128, 255, cv2.THRESH_BINARY)\n\n# Define seed points (adjust as needed)\nseed_point_Cancer = (50, 50)  # Example seed point for Cancer image\nseed_point_NonCancer = (50, 50)  # Example seed point for NonCancer image\n\n# Apply seeded flood fill\nfilled_imageCancer = apply_seeded_flood_fill(binary_imageCancer, seed_point_Cancer, new_value=255)\nfilled_imageNonCancer = apply_seeded_flood_fill(binary_imageNonCancer, seed_point_NonCancer, new_value=255)\n\n# Visualize the results\nfig, axes = plt.subplots(2, 3, figsize=(18, 12))\n\n# Display Binary Image Cancer\naxes[0, 0].imshow(binary_imageCancer, cmap='gray')\naxes[0, 0].set_title('Binary Image Cancer')\naxes[0, 0].axis('off')\n\n# Display Seeded Flood Fill Cancer\naxes[0, 1].imshow(filled_imageCancer, cmap='gray')\naxes[0, 1].set_title('Seeded Flood Fill Cancer')\naxes[0, 1].axis('off')\n\n# Display Original Image Cancer\naxes[0, 2].imshow(image_arrayCancer_8bit, cmap='gray')\naxes[0, 2].set_title('Original Image Cancer')\naxes[0, 2].axis('off')\n\n# Display Binary Image NonCancer\naxes[1, 0].imshow(binary_imageNonCancer, cmap='gray')\naxes[1, 0].set_title('Binary Image NonCancer')\naxes[1, 0].axis('off')\n\n# Display Seeded Flood Fill NonCancer\naxes[1, 1].imshow(filled_imageNonCancer, cmap='gray')\naxes[1, 1].set_title('Seeded Flood Fill NonCancer')\naxes[1, 1].axis('off')\n\n# Display Original Image NonCancer\naxes[1, 2].imshow(image_arrayNonCancer_8bit, cmap='gray')\naxes[1, 2].set_title('Original Image NonCancer')\naxes[1, 2].axis('off')\n\nplt.tight_layout()\nplt.show()\n","metadata":{"execution":{"iopub.status.busy":"2024-08-02T18:41:52.414638Z","iopub.execute_input":"2024-08-02T18:41:52.415046Z","iopub.status.idle":"2024-08-02T18:42:00.404216Z","shell.execute_reply.started":"2024-08-02T18:41:52.415008Z","shell.execute_reply":"2024-08-02T18:42:00.403057Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## Watershed segmentation","metadata":{}},{"cell_type":"code","source":"import cv2\nimport numpy as np\nimport matplotlib.pyplot as plt\n\ndef select_largest_obj(img_bin, lab_val=255, fill_holes=False, smooth_boundary=False, kernel_size=15):\n    '''Select the largest object from a binary image and optionally\n    fill holes inside it and smooth its boundary.'''\n    n_labels, img_labeled, lab_stats, _ = cv2.connectedComponentsWithStats(img_bin, connectivity=8, ltype=cv2.CV_32S)\n    largest_obj_lab = np.argmax(lab_stats[1:, 4]) + 1\n    largest_mask = np.zeros(img_bin.shape, dtype=np.uint8)\n    largest_mask[img_labeled == largest_obj_lab] = lab_val\n    if fill_holes:\n        bkg_locs = np.where(img_labeled == 0)\n        bkg_seed = (bkg_locs[0][0], bkg_locs[1][0])\n        img_floodfill = largest_mask.copy()\n        h_, w_ = largest_mask.shape\n        mask_ = np.zeros((h_ + 2, w_ + 2), dtype=np.uint8)\n        cv2.floodFill(img_floodfill, mask_, seedPoint=bkg_seed, newVal=lab_val)\n        holes_mask = cv2.bitwise_not(img_floodfill)  # Mask of the holes.\n        largest_mask = largest_mask + holes_mask\n    if smooth_boundary:\n        kernel_ = np.ones((kernel_size, kernel_size), dtype=np.uint8)\n        largest_mask = cv2.morphologyEx(largest_mask, cv2.MORPH_OPEN, kernel_)\n    return largest_mask\n\ndef process_image(pixel_array):\n    '''Process a single image array'''\n    # Convert to 8-bit image if necessary\n    if pixel_array.dtype != np.uint8:\n        pixel_array = cv2.normalize(pixel_array, None, 0, 255, cv2.NORM_MINMAX).astype(np.uint8)\n\n    # Apply histogram equalization for better contrast\n    if len(pixel_array.shape) == 2:  # For grayscale images\n        mammo_breast_equ = cv2.equalizeHist(pixel_array)\n    else:\n        mammo_breast_equ = cv2.equalizeHist(cv2.cvtColor(pixel_array, cv2.COLOR_RGB2GRAY))\n\n    # Apply binary thresholding\n    pect_high_inten_thres = 200\n    _, pect_binary_thres = cv2.threshold(mammo_breast_equ, pect_high_inten_thres, 255, cv2.THRESH_BINARY)\n\n    # Markers image for watershed algorithm\n    pect_marker_img = np.zeros(pect_binary_thres.shape, dtype=np.int32)\n\n    # Sure foreground\n    pect_mask_init = select_largest_obj(pect_binary_thres, lab_val=255, fill_holes=True, smooth_boundary=False)\n    kernel_ = np.ones((3, 3), dtype=np.uint8)  # Parameter to tune\n    n_erosions = 7  # Number of erosions to tune\n    pect_mask_eroded = cv2.erode(pect_mask_init, kernel_, iterations=n_erosions)\n    pect_marker_img[pect_mask_eroded > 0] = 255\n\n    # Sure background - breast\n    n_dilations = 7  # Number of dilations to tune\n    pect_mask_dilated = cv2.dilate(pect_mask_init, kernel_, iterations=n_dilations)\n    pect_marker_img[pect_mask_dilated == 0] = 128\n\n    # Sure background - general background\n    mammo_breast_mask = select_largest_obj(pect_binary_thres, lab_val=64, fill_holes=False, smooth_boundary=False)  # Adjust if needed\n    pect_marker_img[mammo_breast_mask == 0] = 64\n\n    return mammo_breast_equ, pect_binary_thres, pect_marker_img\n\n# Example pixel arrays for two images\npixel_array1 = ...  # Your first image array\npixel_array2 = ...  # Your second image array\n\n# Process both images\nmammo_breast_equ1, pect_binary_thres1, pect_marker_img1 = process_image(image_arrayCancer)\nmammo_breast_equ2, pect_binary_thres2, pect_marker_img2 = process_image(image_arrayNonCancer)\n\n# Plot results for both images\nfig, axes = plt.subplots(2, 3, figsize=(18, 12))\n\n# Image 1\naxes[0, 0].imshow(mammo_breast_equ1, cmap='gray')\naxes[0, 0].set_title('Histogram Equalized Image 1')\naxes[0, 0].axis('off')\n\naxes[0, 1].imshow(pect_binary_thres1, cmap='gray')\naxes[0, 1].set_title('Thresholded Image 1')\naxes[0, 1].axis('off')\n\naxes[0, 2].imshow(pect_marker_img1, cmap='gray')\naxes[0, 2].set_title('Watershed Markers Image 1')\naxes[0, 2].axis('off')\n\n# Image 2\naxes[1, 0].imshow(mammo_breast_equ2, cmap='gray')\naxes[1, 0].set_title('Histogram Equalized Image 2')\naxes[1, 0].axis('off')\n\naxes[1, 1].imshow(pect_binary_thres2, cmap='gray')\naxes[1, 1].set_title('Thresholded Image 2')\naxes[1, 1].axis('off')\n\naxes[1, 2].imshow(pect_marker_img2, cmap='gray')\naxes[1, 2].set_title('Watershed Markers Image 2')\naxes[1, 2].axis('off')\n\nplt.tight_layout()\nplt.show()\n","metadata":{"execution":{"iopub.status.busy":"2024-08-02T18:42:00.40901Z","iopub.execute_input":"2024-08-02T18:42:00.409492Z","iopub.status.idle":"2024-08-02T18:42:09.673569Z","shell.execute_reply.started":"2024-08-02T18:42:00.409451Z","shell.execute_reply":"2024-08-02T18:42:09.67246Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"import numpy as np\nimport cv2\nimport matplotlib.pyplot as plt\nimport pydicom\n\n# Define the helper function for selecting the largest object\ndef select_largest_obj(img_bin, lab_val=255, fill_holes=False, smooth_boundary=False, kernel_size=15):\n    '''Select the largest object from a binary image and optionally fill holes and smooth its boundary.'''\n    n_labels, img_labeled, lab_stats, _ = cv2.connectedComponentsWithStats(img_bin, connectivity=8, ltype=cv2.CV_32S)\n    largest_obj_lab = np.argmax(lab_stats[1:, 4]) + 1\n    largest_mask = np.zeros(img_bin.shape, dtype=np.uint8)\n    largest_mask[img_labeled == largest_obj_lab] = lab_val\n    \n    if fill_holes:\n        bkg_locs = np.where(img_labeled == 0)\n        bkg_seed = (bkg_locs[0][0], bkg_locs[1][0])\n        img_floodfill = largest_mask.copy()\n        h_, w_ = largest_mask.shape\n        mask_ = np.zeros((h_ + 2, w_ + 2), dtype=np.uint8)\n        cv2.floodFill(img_floodfill, mask_, seedPoint=bkg_seed, newVal=lab_val)\n        holes_mask = cv2.bitwise_not(img_floodfill)\n        largest_mask = largest_mask + holes_mask\n    \n    if smooth_boundary:\n        kernel_ = np.ones((kernel_size, kernel_size), dtype=np.uint8)\n        largest_mask = cv2.morphologyEx(largest_mask, cv2.MORPH_OPEN, kernel_)\n    \n    return largest_mask\n\ndef watershed_segmentation(img, high_inten_thres=200):\n    '''Apply watershed segmentation to the image.'''\n    # Convert to 8-bit if necessary\n    if img.dtype != np.uint8:\n        img = cv2.normalize(img, None, 0, 255, cv2.NORM_MINMAX).astype(np.uint8)\n    \n    # Histogram equalization\n    if len(img.shape) == 2:  # For grayscale images\n        img_eq = cv2.equalizeHist(img)\n    else:\n        img_eq = cv2.equalizeHist(cv2.cvtColor(img, cv2.COLOR_RGB2GRAY))\n    \n    # Apply binary thresholding\n    _, binary_thres = cv2.threshold(img_eq, high_inten_thres, 255, cv2.THRESH_BINARY)\n    \n    # Marker image for watershed algorithm\n    marker_img = np.zeros(binary_thres.shape, dtype=np.int32)\n    \n    # Sure foreground\n    mask_init = select_largest_obj(binary_thres, lab_val=255, fill_holes=True, smooth_boundary=False)\n    kernel_ = np.ones((3, 3), dtype=np.uint8)\n    mask_eroded = cv2.erode(mask_init, kernel_, iterations=7)\n    marker_img[mask_eroded > 0] = 255\n    \n    # Sure background\n    mask_dilated = cv2.dilate(mask_init, kernel_, iterations=7)\n    marker_img[mask_dilated == 0] = 128\n    \n    # General background\n    breast_mask = select_largest_obj(binary_thres, lab_val=64, fill_holes=False, smooth_boundary=False)\n    marker_img[breast_mask == 0] = 64\n    \n    # Convert the image to BGR for watershed algorithm\n    img_eq_3c = cv2.cvtColor(img_eq, cv2.COLOR_GRAY2BGR)\n    \n    # Apply watershed segmentation\n    cv2.watershed(img_eq_3c, marker_img)\n    mask_watershed = marker_img.copy()\n    img_eq_3c[mask_watershed == -1] = (0, 0, 255)  # Mark boundaries\n    \n    return mask_watershed, img_eq_3c\n\ndef process_and_plot_images(image_arrays, titles, high_inten_thres=200):\n    '''Process a list of images and plot the results.'''\n    fig, axes = plt.subplots(len(image_arrays), 2, figsize=(12, 9 * len(image_arrays)))\n    \n    for i, img_array in enumerate(image_arrays):\n        mask_watershed, segmented_img = watershed_segmentation(img_array, high_inten_thres)\n        \n        axes[i, 0].imshow(mask_watershed, cmap='gray')\n        axes[i, 0].set_title(f'{titles[i]} - Watershed Mask')\n        axes[i, 0].axis('off')\n        \n        axes[i, 1].imshow(segmented_img)\n        axes[i, 1].set_title(f'{titles[i]} - Watershed Segmentation')\n        axes[i, 1].axis('off')\n    \n    plt.show()\n\n# Example usage\n# Load DICOM files for cancer and non-cancer (replace with actual file paths)\n# cancer_dicom = pydicom.dcmread('path/to/cancer.dcm')\n# non_cancer_dicom = pydicom.dcmread('path/to/non_cancer.dcm')\n\n# Extract pixel arrays\n# cancer_pixel_array = cancer_dicom.pixel_array\n# non_cancer_pixel_array = non_cancer_dicom.pixel_array\n\n# List of images and their titles\nimage_arrays = [image_arrayCancer, image_arrayNonCancer]\ntitles = ['Cancer Image', 'Non-Cancer Image']\n\n# Process and plot images\nprocess_and_plot_images(image_arrays, titles, high_inten_thres=200)\n","metadata":{"execution":{"iopub.status.busy":"2024-08-02T18:42:09.675059Z","iopub.execute_input":"2024-08-02T18:42:09.675807Z","iopub.status.idle":"2024-08-02T18:42:19.486392Z","shell.execute_reply.started":"2024-08-02T18:42:09.675768Z","shell.execute_reply":"2024-08-02T18:42:19.485251Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# 4. A preprocessing class","metadata":{}},{"cell_type":"code","source":"import numpy as np\nimport cv2\nimport matplotlib.pyplot as plt\n\nclass DMImagePreprocessor:\n    def __init__(self):\n        pass\n\n    def process(self, img):\n        '''Process the image: applies normalization, histogram equalization, and other preprocessing steps.'''\n        # Convert to 8-bit image if necessary\n        if img.dtype != np.uint8:\n            img = cv2.normalize(img, None, 0, 255, cv2.NORM_MINMAX).astype(np.uint8)\n        \n        # Apply histogram equalization\n        if len(img.shape) == 2:  # For grayscale images\n            processed_img = cv2.equalizeHist(img)\n        else:\n            processed_img = cv2.equalizeHist(cv2.cvtColor(img, cv2.COLOR_RGB2GRAY))\n        \n        # Apply color map for visualization\n        color_img = cv2.applyColorMap(processed_img, cv2.COLORMAP_JET)\n        \n        return processed_img, color_img\n\n# Initialize the preprocessor\ndm_img_pproc = DMImagePreprocessor()\n\n# Replace `cancer_image_array` and `non_cancer_image_array` with your actual image arrays\n# cancer_image_array = np.array(cancer_dicom_file.pixel_array)\n# non_cancer_image_array = np.array(non_cancer_dicom_file.pixel_array)\n\n# Process the DICOM images\ncancer_pproc, cancer_col = dm_img_pproc.process(image_arrayCancer)\nnon_cancer_pproc, non_cancer_col = dm_img_pproc.process(image_arrayNonCancer)\n\n# Plot the results\nfig, axes = plt.subplots(2, 3, figsize=(18, 18))\n\n# Cancer image\naxes[0, 0].imshow(image_arrayCancer, cmap='gray')\naxes[0, 0].set_title('Cancer Original Image')\naxes[0, 0].axis('off')\n\naxes[0, 1].imshow(cancer_pproc, cmap='gray')\naxes[0, 1].set_title('Cancer Processed Image')\naxes[0, 1].axis('off')\n\naxes[0, 2].imshow(cancer_col)\naxes[0, 2].set_title('Cancer Colorized Image')\naxes[0, 2].axis('off')\n\n# Non-cancer image\naxes[1, 0].imshow(image_arrayNonCancer, cmap='gray')\naxes[1, 0].set_title('Non-Cancer Original Image')\naxes[1, 0].axis('off')\n\naxes[1, 1].imshow(non_cancer_pproc, cmap='gray')\naxes[1, 1].set_title('Non-Cancer Processed Image')\naxes[1, 1].axis('off')\n\naxes[1, 2].imshow(non_cancer_col)\naxes[1, 2].set_title('Non-Cancer Colorized Image')\naxes[1, 2].axis('off')\n\nplt.show()\n","metadata":{"execution":{"iopub.status.busy":"2024-08-02T18:42:19.48798Z","iopub.execute_input":"2024-08-02T18:42:19.488381Z","iopub.status.idle":"2024-08-02T18:42:30.600485Z","shell.execute_reply.started":"2024-08-02T18:42:19.488347Z","shell.execute_reply":"2024-08-02T18:42:30.599278Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"import numpy as np\nimport cv2\nimport matplotlib.pyplot as plt\nimport pydicom\n\nclass DMImagePreprocessor:\n    def __init__(self):\n        pass\n\n    def process(self, img):\n        '''Process the image: applies normalization, histogram equalization, and other preprocessing steps.'''\n        if img.dtype != np.uint8:\n            img = cv2.normalize(img, None, 0, 255, cv2.NORM_MINMAX).astype(np.uint8)\n        \n        if len(img.shape) == 2:  # For grayscale images\n            processed_img = cv2.equalizeHist(img)\n        else:\n            processed_img = cv2.equalizeHist(cv2.cvtColor(img, cv2.COLOR_RGB2GRAY))\n        \n        color_img = cv2.applyColorMap(processed_img, cv2.COLORMAP_JET)\n        \n        return processed_img, color_img\n\n    def suppress_artifacts(self, img):\n        '''Suppress artifacts in the image. This is a placeholder function.'''\n        # Example: Applying a Gaussian blur to remove noise\n        artifact_removed = cv2.GaussianBlur(img, (5, 5), 0)\n        # Create a binary mask for the breast region (example method)\n        _, binary_mask = cv2.threshold(artifact_removed, 128, 255, cv2.THRESH_BINARY)\n        return artifact_removed, binary_mask\n\n    def remove_pectoral(self, img, mask):\n        '''Remove pectoral muscle from the image using a binary mask.'''\n        # Ensure mask is uint8\n        if mask.dtype != np.uint8:\n            mask = mask.astype(np.uint8)\n        \n        # Check dimensions\n        if img.shape != mask.shape:\n            raise ValueError(\"Image and mask dimensions do not match.\")\n        \n        # Invert the mask to remove pectoral\n        mask_inv = cv2.bitwise_not(mask)\n        pectoral_removed = cv2.bitwise_and(img, img, mask=mask_inv)\n        return pectoral_removed, mask_inv\n\ndef load_dicom_image(file_path):\n    '''Load a DICOM image from a file path.'''\n    dicom_file = pydicom.dcmread(file_path)\n    return dicom_file.pixel_array\n\n# Load the DICOM files\n# cancer_image_array = load_dicom_image('path/to/cancer_image.dcm')\n# non_cancer_image_array = load_dicom_image('path/to/non_cancer_image.dcm')\n\n# Initialize the preprocessor\ndm_img_pproc = DMImagePreprocessor()\n\n# Process the DICOM images for cancer\ncancer_pproc, cancer_col = dm_img_pproc.process(image_arrayCancer)\ncancer_artif_removed, cancer_breast_mask = dm_img_pproc.suppress_artifacts(image_arrayCancer)\ncancer_pectoral_removed, _ = dm_img_pproc.remove_pectoral(cancer_artif_removed, cancer_breast_mask)\n\n# Process the DICOM images for non-cancer\nnon_cancer_pproc, non_cancer_col = dm_img_pproc.process(image_arrayNonCancer)\nnon_cancer_artif_removed, non_cancer_breast_mask = dm_img_pproc.suppress_artifacts(image_arrayNonCancer)\nnon_cancer_pectoral_removed, _ = dm_img_pproc.remove_pectoral(non_cancer_artif_removed, non_cancer_breast_mask)\n\n# Plot the results for cancer\nfig, axes = plt.subplots(2, 3, figsize=(18, 18))\n\naxes[0, 0].imshow(image_arrayCancer, cmap='gray')\naxes[0, 0].set_title('Cancer Original Image')\naxes[0, 0].axis('off')\n\naxes[0, 1].imshow(cancer_artif_removed, cmap='gray')\naxes[0, 1].set_title('Cancer Artifact Removed')\naxes[0, 1].axis('off')\n\naxes[0, 2].imshow(cancer_pectoral_removed, cmap='gray')\naxes[0, 2].set_title('Cancer Pectoral Muscle Removed')\naxes[0, 2].axis('off')\n\n# Plot the results for non-cancer\naxes[1, 0].imshow(image_arrayNonCancer, cmap='gray')\naxes[1, 0].set_title('Non-Cancer Original Image')\naxes[1, 0].axis('off')\n\naxes[1, 1].imshow(non_cancer_artif_removed, cmap='gray')\naxes[1, 1].set_title('Non-Cancer Artifact Removed')\naxes[1, 1].axis('off')\n\naxes[1, 2].imshow(non_cancer_pectoral_removed, cmap='gray')\naxes[1, 2].set_title('Non-Cancer Pectoral Muscle Removed')\naxes[1, 2].axis('off')\n\nplt.show()\n","metadata":{"execution":{"iopub.status.busy":"2024-08-02T18:42:30.602558Z","iopub.execute_input":"2024-08-02T18:42:30.603082Z","iopub.status.idle":"2024-08-02T18:42:38.619216Z","shell.execute_reply.started":"2024-08-02T18:42:30.603037Z","shell.execute_reply":"2024-08-02T18:42:38.617932Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"import numpy as np\nimport cv2\nimport matplotlib.pyplot as plt\nimport pydicom\n\nclass DMImagePreprocessor:\n    def __init__(self):\n        pass\n\n    def process(self, img, pect_removal=False, high_int_threshold=128):\n        '''Process the image: applies normalization, histogram equalization, and optionally remove pectoral muscle.'''\n        # Convert to 8-bit if necessary\n        if img.dtype != np.uint8:\n            img = cv2.normalize(img, None, 0, 255, cv2.NORM_MINMAX).astype(np.uint8)\n        \n        # Histogram equalization\n        if len(img.shape) == 2:  # For grayscale images\n            processed_img = cv2.equalizeHist(img)\n        else:\n            processed_img = cv2.equalizeHist(cv2.cvtColor(img, cv2.COLOR_RGB2GRAY))\n        \n        # Optionally remove pectoral muscle\n        if pect_removal:\n            artifact_removed, breast_mask = self.suppress_artifacts(processed_img)\n            processed_img, _ = self.remove_pectoral(artifact_removed, breast_mask)\n        \n        # Apply binary thresholding based on the high intensity threshold parameter\n        _, binary_img = cv2.threshold(processed_img, high_int_threshold, 255, cv2.THRESH_BINARY)\n        \n        return processed_img, binary_img\n\n    def suppress_artifacts(self, img):\n        '''Suppress artifacts in the image.'''\n        artifact_removed = cv2.GaussianBlur(img, (5, 5), 0)\n        _, binary_mask = cv2.threshold(artifact_removed, 128, 255, cv2.THRESH_BINARY)\n        return artifact_removed, binary_mask\n\n    def remove_pectoral(self, img, mask):\n        '''Remove pectoral muscle from the image using a binary mask.'''\n        if mask.dtype != np.uint8:\n            mask = mask.astype(np.uint8)\n        \n        if img.shape != mask.shape:\n            raise ValueError(\"Image and mask dimensions do not match.\")\n        \n        mask_inv = cv2.bitwise_not(mask)\n        pectoral_removed = cv2.bitwise_and(img, img, mask=mask_inv)\n        return pectoral_removed, mask_inv\n\n# Function to process and plot images\ndef process_and_plot_images(image_arrays, titles, pect_removal=False, high_int_threshold=128):\n    '''Process a list of images and plot the results.'''\n    dm_img_pproc = DMImagePreprocessor()\n    \n    fig, axes = plt.subplots(len(image_arrays), 2, figsize=(12, 9 * len(image_arrays)))\n    \n    for i, img_array in enumerate(image_arrays):\n        processed_img, _ = dm_img_pproc.process(img_array, pect_removal=pect_removal, high_int_threshold=high_int_threshold)\n        \n        axes[i, 0].imshow(img_array, cmap='gray')\n        axes[i, 0].set_title(f'{titles[i]} - Original')\n        axes[i, 0].axis('off')\n        \n        axes[i, 1].imshow(processed_img, cmap='gray')\n        axes[i, 1].set_title(f'{titles[i]} - Processed')\n        axes[i, 1].axis('off')\n    \n    plt.show()\n\n# Example usage\n# Load DICOM files for cancer and non-cancer (replace with actual file paths)\n# cancer_dicom = pydicom.dcmread('path/to/cancer.dcm')\n# non_cancer_dicom = pydicom.dcmread('path/to/non_cancer.dcm')\n\n# Extract pixel arrays\n# cancer_pixel_array = cancer_dicom.pixel_array\n# non_cancer_pixel_array = non_cancer_dicom.pixel_array\n\n# List of images and their titles\nimage_arrays = [image_arrayCancer, image_arrayNonCancer]\ntitles = ['Cancer Image', 'Non-Cancer Image']\n\n# Process and plot images\nprocess_and_plot_images(image_arrays, titles, pect_removal=True, high_int_threshold=100)\n","metadata":{"execution":{"iopub.status.busy":"2024-08-02T18:42:38.620833Z","iopub.execute_input":"2024-08-02T18:42:38.621199Z","iopub.status.idle":"2024-08-02T18:42:43.954031Z","shell.execute_reply.started":"2024-08-02T18:42:38.621167Z","shell.execute_reply":"2024-08-02T18:42:43.952913Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"","metadata":{},"execution_count":null,"outputs":[]}]}