{"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":"# RSNA-2022 : [beginner] Detect breast area - OpenCV - connectedComponents","metadata":{"id":"6LF5_Ea6qRFQ"}},{"cell_type":"markdown","source":"**[Change Log]**\n- Ver.0 : 1st notebook.","metadata":{}},{"cell_type":"markdown","source":"# Objective\nThis notebook tries to describe the <font color=red>breast area detection in simple and comfortable way to read</font> for the [RSNA Screening Mammography Breast Cancer Detection](https://www.kaggle.com/competitions/rsna-breast-cancer-detection), so some codes might seem redundant but that is the concept.\n\nThis notebook is inspirated by the https://www.kaggle.com/code/vslaykovsky/rsna-cut-off-empty-space-from-images/data as an example of breast area detection. \n\nThe key contents are;\n1. Read DICOM file by **pydicom**\n2. Get breast area by **Opencv - connectedComponents**\n3. I. Show rectangle for breast area by **Opencv - rectangle & putText**\n4. Create png file\n5. II. Show rectangle for breast area by **Opencv - rectangle & putText**","metadata":{}},{"cell_type":"code","source":"# Import library\n\n## Basic ##\n%matplotlib inline\nimport pandas as pd\nimport numpy as np\nimport matplotlib.pyplot as plt\nimport os\nimport cv2\nimport glob\nfrom tqdm.notebook import tqdm\n\n## dicom file relevant\n!pip install -qU pylibjpeg\nimport pylibjpeg\nimport pydicom\n\n## Setting for matplotlib\nplt.rcParams['font.size'] = 14\nplt.rcParams['figure.figsize'] = (6,6)\nplt.rcParams['axes.grid'] = True\n","metadata":{"id":"EBjRX49eqRFd","execution":{"iopub.status.busy":"2022-12-27T03:47:32.506073Z","iopub.execute_input":"2022-12-27T03:47:32.50672Z","iopub.status.idle":"2022-12-27T03:47:45.926307Z","shell.execute_reply.started":"2022-12-27T03:47:32.506687Z","shell.execute_reply":"2022-12-27T03:47:45.925214Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Data Preparation","metadata":{"id":"obUglB5x18Lk"}},{"cell_type":"markdown","source":"## csv data","metadata":{}},{"cell_type":"code","source":"test_df = pd.read_csv('/kaggle/input/rsna-breast-cancer-detection/test.csv')\ntest_df","metadata":{"execution":{"iopub.status.busy":"2022-12-27T03:47:45.928623Z","iopub.execute_input":"2022-12-27T03:47:45.928956Z","iopub.status.idle":"2022-12-27T03:47:45.964331Z","shell.execute_reply.started":"2022-12-27T03:47:45.928924Z","shell.execute_reply":"2022-12-27T03:47:45.96321Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# train data\ntrain_df = pd.read_csv('/kaggle/input/rsna-breast-cancer-detection/train.csv')\ntrain_df.head()","metadata":{"id":"-utQupbe6dYb","execution":{"iopub.status.busy":"2022-12-27T03:47:45.96557Z","iopub.execute_input":"2022-12-27T03:47:45.965974Z","iopub.status.idle":"2022-12-27T03:47:46.084388Z","shell.execute_reply.started":"2022-12-27T03:47:45.965942Z","shell.execute_reply":"2022-12-27T03:47:46.083319Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# train data with cancer == 1\ntrain_cancer_df = train_df[train_df['cancer'] == 1]\ntrain_cancer_df.head()","metadata":{"execution":{"iopub.status.busy":"2022-12-27T03:47:46.085854Z","iopub.execute_input":"2022-12-27T03:47:46.086383Z","iopub.status.idle":"2022-12-27T03:47:46.112278Z","shell.execute_reply.started":"2022-12-27T03:47:46.086352Z","shell.execute_reply":"2022-12-27T03:47:46.111239Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# path to image data\nfname = glob.glob('/kaggle/input/rsna-breast-cancer-detection/train_images/*/*.dcm')[0]","metadata":{"execution":{"iopub.status.busy":"2022-12-27T03:47:46.115203Z","iopub.execute_input":"2022-12-27T03:47:46.115606Z","iopub.status.idle":"2022-12-27T03:48:26.291487Z","shell.execute_reply.started":"2022-12-27T03:47:46.115566Z","shell.execute_reply":"2022-12-27T03:48:26.290457Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"As described in the [EDA](https://www.kaggle.com/code/masatakaitakura/eda-for-beginner-rsna-mammography-breast-cancer), mammography image have a \"Blank area\", so we use the <font color=red>**OpenCV - connectedComponents** </font>and show the breast area by rectangles.","metadata":{}},{"cell_type":"markdown","source":"# 1. DICOM image data","metadata":{}},{"cell_type":"code","source":"# Read DICOM file\nsize=512\n\n# obtain ids\npatient_id = fname.split('/')[-2]\nimage_id = fname.split('/')[-1][:-4]\n\n# image file\ndicom = pydicom.dcmread(fname)\nimg = dicom.pixel_array\n# normalize image\nimg = (img - img.min()) / (img.max() - img.min())\nif dicom.PhotometricInterpretation == \"MONOCHROME1\":\n    img = 1 - img\nimg = cv2.resize(img, (size, size))","metadata":{"execution":{"iopub.status.busy":"2022-12-27T04:08:18.278933Z","iopub.execute_input":"2022-12-27T04:08:18.27941Z","iopub.status.idle":"2022-12-27T04:08:18.994032Z","shell.execute_reply.started":"2022-12-27T04:08:18.279373Z","shell.execute_reply":"2022-12-27T04:08:18.992992Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# show image\nfig, ax = plt.subplots(1, 1)\nax.imshow(img)\n\nimg.shape","metadata":{"execution":{"iopub.status.busy":"2022-12-27T04:08:19.748156Z","iopub.execute_input":"2022-12-27T04:08:19.748948Z","iopub.status.idle":"2022-12-27T04:08:20.0471Z","shell.execute_reply.started":"2022-12-27T04:08:19.748908Z","shell.execute_reply":"2022-12-27T04:08:20.045992Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# 2. Get breast area by Opencv - connectedComponents\nAs shown above image, the mammography image contain the \"blank\" area. To focus the breast for the training, we create the fit image by using **<font color=red>cv2.connectedcomponent</font>**.","metadata":{}},{"cell_type":"markdown","source":"The image file's contrast is described in pixel_array with 0 to 1 integer.   \n**(Dark) 0 <--> 1 (Blight)**","metadata":{}},{"cell_type":"code","source":"print(f'max: {img.max()}, min: {img.min()}')","metadata":{"execution":{"iopub.status.busy":"2022-12-27T03:48:27.584225Z","iopub.execute_input":"2022-12-27T03:48:27.584566Z","iopub.status.idle":"2022-12-27T03:48:27.591651Z","shell.execute_reply.started":"2022-12-27T03:48:27.584532Z","shell.execute_reply":"2022-12-27T03:48:27.590603Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"We can obtain the blight area by the following comparison operator and obtain the boolean equation. \n- False : darker than criteria\n- True : lighter tahn criteria","metadata":{}},{"cell_type":"code","source":"# obtain the blight area\nimg > 0.05","metadata":{"execution":{"iopub.status.busy":"2022-12-27T03:51:32.680139Z","iopub.execute_input":"2022-12-27T03:51:32.680584Z","iopub.status.idle":"2022-12-27T03:51:32.689454Z","shell.execute_reply.started":"2022-12-27T03:51:32.68055Z","shell.execute_reply":"2022-12-27T03:51:32.688202Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Next, convert the boolean to **np.array** and show the image.","metadata":{}},{"cell_type":"code","source":"# convert boolean to np\nimg_01 = (img > 0.05).astype(np.uint8)[:, :]\nprint(img_01)\n\nplt.imshow(img_01)","metadata":{"execution":{"iopub.status.busy":"2022-12-27T04:55:43.049591Z","iopub.execute_input":"2022-12-27T04:55:43.050809Z","iopub.status.idle":"2022-12-27T04:55:43.706614Z","shell.execute_reply.started":"2022-12-27T04:55:43.050745Z","shell.execute_reply":"2022-12-27T04:55:43.705343Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"The [connectedComponentsWithStats](https://docs.opencv.org/4.6.0/d3/dc0/group__imgproc__shape.html#ga5ed7784614678adccb699c70fb841075) computes \n- the connected components labeled image of boolean image\n- a statistics output for each label","metadata":{}},{"cell_type":"code","source":"# label the regions of non-empty pixels\nretval, labels, stats, centroids = cv2.connectedComponentsWithStats(image=(img > 0.05).astype(np.uint8)[:, :], connectivity=8, ltype=cv2.CV_32S)","metadata":{"execution":{"iopub.status.busy":"2022-12-27T03:56:22.381531Z","iopub.execute_input":"2022-12-27T03:56:22.382588Z","iopub.status.idle":"2022-12-27T03:56:22.39826Z","shell.execute_reply.started":"2022-12-27T03:56:22.382538Z","shell.execute_reply":"2022-12-27T03:56:22.396801Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"The followings are specific outputs.\n- retval:\nReturns the number of labels. Note that the background is also counted as one.\n- labels:\nOutput labeling image\n- stats:\nIt stores the region information for each label.\n  **[region top left x coordinate, region top left y coordinate, region width, region height, area]** Area is the number of pixels in the region.\n-centroids:\nCentroid information for each label is stored.","metadata":{}},{"cell_type":"code","source":"print(f'retval: {retval}\\n\\nlabels: {labels}\\nlabel shape={labels.shape}\\n\\nstats={stats}\\n\\ncentroids={centroids}')","metadata":{"execution":{"iopub.status.busy":"2022-12-27T03:56:23.637873Z","iopub.execute_input":"2022-12-27T03:56:23.638859Z","iopub.status.idle":"2022-12-27T03:56:23.647741Z","shell.execute_reply.started":"2022-12-27T03:56:23.638817Z","shell.execute_reply":"2022-12-27T03:56:23.646492Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"The `labels` can show the colored grouping as shown below.\n\n### Findings\n- some image includes the letter of [laterality] and [view], those area are around 10 - 80.","metadata":{}},{"cell_type":"code","source":"plt.imshow(labels)","metadata":{"execution":{"iopub.status.busy":"2022-12-27T03:59:03.642935Z","iopub.execute_input":"2022-12-27T03:59:03.643388Z","iopub.status.idle":"2022-12-27T03:59:04.074579Z","shell.execute_reply.started":"2022-12-27T03:59:03.643354Z","shell.execute_reply":"2022-12-27T03:59:04.073829Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"for i in range(0,retval): # retval[0] is background\n    x, y, width, height, area = stats[i]\n    if area > 10: # judge if the area is larger than the value\n        print(f'index-{i}: area={area}')","metadata":{"execution":{"iopub.status.busy":"2022-12-27T03:57:06.957103Z","iopub.execute_input":"2022-12-27T03:57:06.957582Z","iopub.status.idle":"2022-12-27T03:57:06.966151Z","shell.execute_reply.started":"2022-12-27T03:57:06.957546Z","shell.execute_reply":"2022-12-27T03:57:06.964732Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# 3.I. Show rectangle for breast area by Opencv - rectangle & putText\nNext, to identify the each detected areas, let's add rectangle and text to each.  \nNote: To avoid the change of original `img`, copy the image and create `img_rect`.","metadata":{}},{"cell_type":"code","source":"# set img_rect\nimport copy\nimg_rect = copy.copy(img)\n# plt.imshow(img_rect)\n\ni = 0\n\nfor i in range(1,retval): # retval[0] = background\n    x, y, width, height, area = stats[i] \n\n    if area > 10: # detect more than 10 pixcel area\n        cv2.rectangle(img_rect,\n                      pt1=(x,y),\n                      pt2=(x+width,y+height),\n                      color=(0,0,255),\n                      thickness=5)\n        cv2.putText(img_rect, f\"[{i}]:{area}\", (x, y-10), cv2.FONT_HERSHEY_PLAIN, 1, (255, 0, 0), 1, cv2.LINE_AA)\nplt.imshow(img_rect)","metadata":{"execution":{"iopub.status.busy":"2022-12-27T04:09:14.018473Z","iopub.execute_input":"2022-12-27T04:09:14.018998Z","iopub.status.idle":"2022-12-27T04:09:14.305838Z","shell.execute_reply.started":"2022-12-27T04:09:14.018954Z","shell.execute_reply":"2022-12-27T04:09:14.304907Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"I could not show the OpenCV - rectangle on the image which is directly created from DICOM file, so I create png file from DICOM file and use for OpenCV image.","metadata":{}},{"cell_type":"markdown","source":"# 4. Create png file","metadata":{}},{"cell_type":"code","source":"train_images = glob.glob(\"/kaggle/input/rsna-breast-cancer-detection/train_images/*/*.dcm\")\n\nlen(train_images)  # 54706","metadata":{"execution":{"iopub.status.busy":"2022-12-27T04:19:24.452787Z","iopub.execute_input":"2022-12-27T04:19:24.453273Z","iopub.status.idle":"2022-12-27T04:19:32.641061Z","shell.execute_reply.started":"2022-12-27T04:19:24.453236Z","shell.execute_reply":"2022-12-27T04:19:32.63973Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"SAVE_FOLDER = \"png/\"\nSIZE = 512\nEXTENSION = \"png\"\n\nos.makedirs(SAVE_FOLDER, exist_ok=True)","metadata":{"execution":{"iopub.status.busy":"2022-12-27T04:19:32.643054Z","iopub.execute_input":"2022-12-27T04:19:32.64365Z","iopub.status.idle":"2022-12-27T04:19:32.648355Z","shell.execute_reply.started":"2022-12-27T04:19:32.643616Z","shell.execute_reply":"2022-12-27T04:19:32.647439Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# create and save png file\ndef process(f, size=512, save_folder=\"png/\", extension=\"png\"):\n    patient = f.split('/')[-2]\n    image = f.split('/')[-1][:-4]\n\n    dicom = pydicom.dcmread(f)\n    img = dicom.pixel_array\n\n    img = (img - img.min()) / (img.max() - img.min())\n\n    if dicom.PhotometricInterpretation == \"MONOCHROME1\":\n        img = 1 - img\n\n    img = cv2.resize(img, (size, size))\n\n    cv2.imwrite(save_folder + f\"{patient}_{image}.{extension}\", (img * 255).astype(np.uint8))","metadata":{"execution":{"iopub.status.busy":"2022-12-27T04:19:32.649499Z","iopub.execute_input":"2022-12-27T04:19:32.649821Z","iopub.status.idle":"2022-12-27T04:19:32.660447Z","shell.execute_reply.started":"2022-12-27T04:19:32.649792Z","shell.execute_reply":"2022-12-27T04:19:32.659494Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"for f in tqdm(train_images[:6]): #if you want to process all, remove [:25]\n    process(f, save_folder=SAVE_FOLDER)","metadata":{"execution":{"iopub.status.busy":"2022-12-27T04:19:32.662478Z","iopub.execute_input":"2022-12-27T04:19:32.662821Z","iopub.status.idle":"2022-12-27T04:19:42.449345Z","shell.execute_reply.started":"2022-12-27T04:19:32.662791Z","shell.execute_reply":"2022-12-27T04:19:42.448261Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"### Confirm the created png file","metadata":{}},{"cell_type":"code","source":"png_images = glob.glob('/kaggle/working/png/*.png')\n\nlen(png_images)","metadata":{"execution":{"iopub.status.busy":"2022-12-27T04:19:44.179477Z","iopub.execute_input":"2022-12-27T04:19:44.18024Z","iopub.status.idle":"2022-12-27T04:19:44.188601Z","shell.execute_reply.started":"2022-12-27T04:19:44.180197Z","shell.execute_reply":"2022-12-27T04:19:44.187718Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Plot result\nrow = 2\ncol = len(png_images)//row\n\nfig, axes = plt.subplots(row, col, figsize=(20,20))\n    \nfor idx, image in enumerate(png_images):\n    frame = cv2.imread(image)\n    i = idx // col\n    j = idx % col\n    axes[i, j].imshow(frame)\n    axes[i, j].set_title(image.split('/')[4]) # change the number with the file path\n\nplt.subplots_adjust(wspace=0.2, hspace=-0.4)\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2022-12-27T04:20:53.224489Z","iopub.execute_input":"2022-12-27T04:20:53.225812Z","iopub.status.idle":"2022-12-27T04:20:54.392495Z","shell.execute_reply.started":"2022-12-27T04:20:53.225761Z","shell.execute_reply":"2022-12-27T04:20:54.391364Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Now we can prepare the png files by specify the train_images dcm file path.","metadata":{}},{"cell_type":"markdown","source":"# 5.II. Show rectangle for breast area by Opencv - rectangle & putText","metadata":{}},{"cell_type":"code","source":"# path to image data (already set at the begining)\nfname = glob.glob('/kaggle/input/rsna-breast-cancer-detection/train_images/*/*.dcm')[0]\nfname","metadata":{"execution":{"iopub.status.busy":"2022-12-27T04:25:03.997183Z","iopub.execute_input":"2022-12-27T04:25:03.997657Z","iopub.status.idle":"2022-12-27T04:25:10.353254Z","shell.execute_reply.started":"2022-12-27T04:25:03.997622Z","shell.execute_reply":"2022-12-27T04:25:10.352265Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# prepare png file\nprocess(fname)\n\n# read png file\npng_path = '/kaggle/working/png/' + fname.split('/')[5] + '_' + fname.split('/')[6][:-4] + '.png'\npng = cv2.imread(png_path) # png for display\npng_gray = cv2.imread(png_path, cv2.IMREAD_GRAYSCALE) # grayscale png for labelling\n\n# label the regions of non-empty pixels\nretval, labels, stats, centroids = cv2.connectedComponentsWithStats(image=(png_gray > 0.05).astype(np.uint8)[:, :], connectivity=8, ltype=cv2.CV_32S)\n\n# show rectangle\ni = 0\n\nfor i in range(1,retval): # retval[0] = background\n    x, y, width, height, area = stats[i] \n\n    if area > 10: # detect more than 10 pixcel area\n        cv2.rectangle(png,\n                      pt1=(x,y),\n                      pt2=(x+width,y+height),\n                      color=(255,0,0),\n                      thickness=5)\n        cv2.putText(png, f\"[{i}]area:{area}\", (x, y+100), cv2.FONT_HERSHEY_PLAIN, 1.5, (255, 0, 0), 2, cv2.LINE_AA)\nplt.imshow(png)","metadata":{"execution":{"iopub.status.busy":"2022-12-27T04:41:15.664187Z","iopub.execute_input":"2022-12-27T04:41:15.664718Z","iopub.status.idle":"2022-12-27T04:41:16.699175Z","shell.execute_reply.started":"2022-12-27T04:41:15.664652Z","shell.execute_reply":"2022-12-27T04:41:16.697317Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Now we can display the rectangle and text on the png image successfully.","metadata":{}},{"cell_type":"markdown","source":"Thank you very much for reading, I hope this notebook helps you!\n\n**If you enjoyed the notebook, please <font color=red>upvote!</font> 🙏 Thank you, appreciate your support!**","metadata":{}}]}