{"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":"## About Notebook\n\nThe notebook encapsulates all DALI and ROI related code in the separate notebook. Run this code directly from your inference notebook to convert all dicom test images to png.\n\n**This notebook is a fork of Outwrest's code [Yolov5 ROI Batch DALI Preprocessing Pipeline](https://www.kaggle.com/code/outwrest/yolov5-roi-batch-dali-preprocessing-pipeline) with minor changes. Which in turn is based on David Austin's [Decode jpeg2000 dicom with DALI](https://www.kaggle.com/code/tivfrvqhs5/decode-jpeg2000-dicom-with-dali) and Theo Viel's [RSNA Breast Baseline - Faster Inference with Dali](https://www.kaggle.com/code/theoviel/rsna-breast-baseline-faster-inference-with-dali) notebooks. All credits go to all of them.**\n\n**For ROI extraction task was used [yolov5 model](https://www.kaggle.com/datasets/remekkinas/rsna-breast-cancer-detection-roi-model) proposed in Remek Kinas's [⭐️⭐️ Breast Cancer - ROI (brest) extractor ⭐️⭐️](https://www.kaggle.com/code/remekkinas/breast-cancer-roi-brest-extractor) notebook.**\n\n\n### How to use it\n\n1. Download .ipynb file and add it to dataset. Or you can use my dataset [Easy to Use DALI + ROI ipynb](https://www.kaggle.com/olegbaryshnikov/easy-to-use-dali-roi-ipynb).\n2. Add all of dataset dependencies to your inference notebook.\n3. Add `%run ../input/notebook-dataset-name/easy-to-use-dali-roi.ipynb` statement to execute the notebook.\n4. Use `%reset -f` to release memory.\n5. Use png images saved in **oimg_dir** path (default '/tmp/imgs') for inference.\n\n```\n%run ../input/easy-to-use-dali-roi-ipynb/easy-to-use-dali-roi.ipynb\n%reset -f\n\nimgs_dir = '/tmp/imgs'\n#inference code\n```\n\nUsage example: [[RSNA] CoaT [TF][Inference]](https://www.kaggle.com/code/olegbaryshnikov/rsna-coat-tf-inference)","metadata":{}},{"cell_type":"markdown","source":"## Imports","metadata":{}},{"cell_type":"code","source":"!pip install -I --no-deps --no-index --find-links /kaggle/input/offline-rsna-dali-pipeline-pips python-gdcm pydicom pylibjpeg dicomsdl nvidia-dali-nightly-cuda110","metadata":{"_kg_hide-input":true,"_kg_hide-output":true,"execution":{"iopub.status.busy":"2023-01-13T12:39:01.864971Z","iopub.execute_input":"2023-01-13T12:39:01.865433Z","iopub.status.idle":"2023-01-13T12:39:20.176139Z","shell.execute_reply.started":"2023-01-13T12:39:01.865341Z","shell.execute_reply":"2023-01-13T12:39:20.174961Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"import os\nimport gc\nimport sys\nimport glob\n\nimport cupy as cp\nimport numpy as np\nimport pandas as pd\n\nimport cv2\nimport torch\nimport pydicom\nimport dicomsdl\nimport matplotlib\n\nfrom tqdm.notebook import tqdm\nimport matplotlib.pyplot as plt\nfrom builtins import StopIteration\nfrom joblib import Parallel, delayed\nfrom pydicom.filebase import DicomBytesIO\n\nimport nvidia.dali.fn as fn\nimport nvidia.dali.ops as ops\nimport nvidia.dali.types as types\nfrom nvidia.dali import pipeline_def\nfrom nvidia.dali.pipeline import Pipeline\nfrom nvidia.dali.types import DALIDataType\nfrom nvidia.dali.plugin.pytorch import DALIGenericIterator\nfrom nvidia.dali.plugin.numba.fn.experimental import numba_function\n\nsys.path.append('/kaggle/input/yolov5-github-repo-tracker/yolov5')\n\nfrom utils.general import non_max_suppression","metadata":{"_kg_hide-input":true,"execution":{"iopub.status.busy":"2023-01-13T12:39:20.178549Z","iopub.execute_input":"2023-01-13T12:39:20.179278Z","iopub.status.idle":"2023-01-13T12:39:24.642553Z","shell.execute_reply.started":"2023-01-13T12:39:20.179228Z","shell.execute_reply":"2023-01-13T12:39:24.641636Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## Config","metadata":{}},{"cell_type":"code","source":"class Config:\n    # dirs\n    img_dir = '../input/rsna-breast-cancer-detection/test_images'\n    oimg_dir = '/tmp/imgs' # output image dir\n    \n    preserve_folder_structure = False\n    save_notebook = True\n    \n    # df\n    df_path = '/kaggle/input/rsna-breast-cancer-detection/test.csv'\n    \n    # roi pipeline \n    batch_size = 32\n    resize_height_to = 640 # yolov5 image size\n    factor = 2 # image downscale factor - saves memory\n    out_img_size_x = 512 # output image size\n    out_img_size_y = 1024 # output image size\n    \n    # dali pipeline args\n    pipeline_args = {\n        'batch_size': 32,\n        'num_threads': 2,\n        'device_id': 0,\n        'py_num_workers': 4,\n        'exec_async': False,\n        'exec_pipelined': False\n    }\n    \n    # yolov5\n    yolov5_model = '/kaggle/input/rsna-breast-cancer-detection-roi-model/rsna-roi-003.pt'\n    yolov5_path = '/kaggle/input/yolov5-github-repo-tracker/yolov5'\n    \n    debug = False","metadata":{"execution":{"iopub.status.busy":"2023-01-13T12:39:24.644156Z","iopub.execute_input":"2023-01-13T12:39:24.644741Z","iopub.status.idle":"2023-01-13T12:39:24.652385Z","shell.execute_reply.started":"2023-01-13T12:39:24.644704Z","shell.execute_reply":"2023-01-13T12:39:24.651487Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## Image processing","metadata":{}},{"cell_type":"code","source":"def __read_metadata(fn):\n    f = open(fn, 'rb')\n    md = f.read(400)\n    f.close()\n    offset = md.index(b'1.2.840.10008.1.2.4.') + 20\n    return md[offset] == 57","metadata":{"_kg_hide-input":true,"execution":{"iopub.status.busy":"2023-01-13T12:39:24.655064Z","iopub.execute_input":"2023-01-13T12:39:24.655689Z","iopub.status.idle":"2023-01-13T12:39:24.672025Z","shell.execute_reply.started":"2023-01-13T12:39:24.655652Z","shell.execute_reply":"2023-01-13T12:39:24.671172Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#get yolov5 and preserve plt backend\ndef get_yolo():\n    b = plt.get_backend()\n    model = torch.hub.load(Config.yolov5_path, 'custom', path = Config.yolov5_model, source = 'local', force_reload = True)\n    matplotlib.use(b)\n    return model","metadata":{"execution":{"iopub.status.busy":"2023-01-13T12:39:24.673547Z","iopub.execute_input":"2023-01-13T12:39:24.673883Z","iopub.status.idle":"2023-01-13T12:39:24.684146Z","shell.execute_reply.started":"2023-01-13T12:39:24.673846Z","shell.execute_reply":"2023-01-13T12:39:24.683216Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df = pd.read_csv(Config.df_path)\n\nif Config.debug:\n    df = df.sample(1000).reset_index(drop = True)\n\ndf['pis'] = Config.img_dir + '/' + df.patient_id.astype(str) + '/' + df.image_id.astype(str) + '.dcm'\ndf['path']= Config.oimg_dir + '/' + df.patient_id.astype(str) + '_' + df.image_id.astype(str) + '.png'\n# low overhead so we can increase n_jobs to max cpu usage\nis_j2k = np.array(Parallel(n_jobs = 24)(delayed(__read_metadata)(ipath) for ipath in tqdm(df.pis.values.tolist())))\ndf_j2k = df.loc[is_j2k]\ndf_jll = df.loc[~is_j2k]\n\nos.makedirs(Config.oimg_dir)\n\nprint('j2k:', len(df_j2k), 'jll:', len(df_jll))\n\nmodel = get_yolo().to('cuda')","metadata":{"execution":{"iopub.status.busy":"2023-01-13T12:39:24.685346Z","iopub.execute_input":"2023-01-13T12:39:24.685798Z","iopub.status.idle":"2023-01-13T12:39:35.506663Z","shell.execute_reply.started":"2023-01-13T12:39:24.685762Z","shell.execute_reply":"2023-01-13T12:39:35.505718Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"class J2KIterator(object):\n    def __init__(self, df, batch_size):\n        self.df = df\n        self.batch_size = batch_size\n       \n    @staticmethod\n    def convert_dicom_to_j2k(patient_id, image_id):\n        dcmfile = pydicom.dcmread(f'{Config.img_dir}/{patient_id}/{image_id}.dcm')\n        with open(f'{Config.img_dir}/{patient_id}/{image_id}.dcm', 'rb') as fp:\n            raw = DicomBytesIO(fp.read())\n            ds = pydicom.dcmread(raw)\n        offset = ds.PixelData.find(b\"\\x00\\x00\\x00\\x0C\")\n\n        return np.frombuffer(ds.PixelData[offset:], dtype=np.uint8), dcmfile.PhotometricInterpretation == 'MONOCHROME1'\n        \n    def __len__(self):\n        return int(np.ceil(len(self.df) / self.batch_size))\n\n    def __iter__(self):\n        self.i = 0\n        return self\n\n    def __next__(self):\n        if self.i > len(self.df):\n            raise StopIteration\n            \n        compressed_imgs = []\n        is_monochrome_imgs = []\n        p_ids = []\n        i_ids = []\n        \n        df = self.df.iloc[range(self.i, min(self.i + self.batch_size, len(self.df)))]\n        \n        for p_id, i_id in df[['patient_id', 'image_id']].values:\n            compressed_img, is_monochrome_img = self.convert_dicom_to_j2k(p_id, i_id)\n            \n            compressed_imgs.append(compressed_img)\n            is_monochrome_imgs.append(np.array([is_monochrome_img], dtype = np.bool_))\n            p_ids.append(np.array([p_id], dtype = np.int64)) \n            i_ids.append(np.array([i_id], dtype = np.int64))\n            \n        self.i += self.batch_size\n            \n        return compressed_imgs, is_monochrome_imgs, p_ids, i_ids\n    \nclass JLLIterator(object):\n    def __init__(self, df, batch_size):\n        self.df = df\n        self.batch_size = batch_size\n        \n    @staticmethod\n    def norm_img(img):\n        i_max = img.max()\n        i_min = img.min()\n        return (255 * ((m - i_min) / (i_max - i_min))).astype(np.uint8)\n       \n    @staticmethod\n    def decompress_jll(dm):\n        info = dm.getPixelDataInfo()\n        img = np.empty((info['Rows'], info['Cols']), dtype = info['dtype'])\n        dm.copyFrameData(0, img)\n        return img\n    \n    def process_img(self, patient_id, image_id):\n        dm = dicomsdl.open(f'{Config.img_dir}/{patient_id}/{image_id}.dcm')\n        \n#         img = norm_img(decompress_jll(dm))\n        \n#         if dm.PhotometricInterpretation== 'MONOCHROME1':\n#             img = 255 - img\n        \n#         return img\n\n        return self.decompress_jll(dm).astype(np.uint16), dm.PhotometricInterpretation == 'MONOCHROME1'\n        \n    def __len__(self):\n        return int(np.ceil(len(self.df) / self.batch_size))\n\n    def __iter__(self):\n        self.i = 0\n        return self\n\n    def __next__(self):\n        if self.i > len(self.df):\n            raise StopIteration\n        \n        decompressed_imgs = []\n        is_monochrome_imgs = []\n        p_ids = []\n        i_ids = []\n        \n        df = self.df.iloc[range(self.i, min(self.i + self.batch_size, len(self.df)))]\n        \n        for p_id, i_id in df[['patient_id', 'image_id']].values:\n            decompressed_img, is_monochrome_img = self.process_img(p_id, i_id)\n            \n            decompressed_imgs.append(np.expand_dims(decompressed_img, axis = 2))\n            is_monochrome_imgs.append(np.array([is_monochrome_img], dtype = np.bool_))\n            p_ids.append(np.array([p_id], dtype = np.int64))\n            i_ids.append(np.array([i_id], dtype = np.int64))\n            \n        self.i += self.batch_size\n                    \n        return decompressed_imgs, is_monochrome_imgs, p_ids, i_ids","metadata":{"_kg_hide-input":true,"execution":{"iopub.status.busy":"2023-01-13T12:39:35.508909Z","iopub.execute_input":"2023-01-13T12:39:35.50967Z","iopub.status.idle":"2023-01-13T12:39:35.529931Z","shell.execute_reply.started":"2023-01-13T12:39:35.509618Z","shell.execute_reply":"2023-01-13T12:39:35.528861Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"J2Ki = J2KIterator(df_j2k, Config.batch_size)\nJLLi = JLLIterator(df_jll, Config.batch_size)","metadata":{"_kg_hide-input":true,"execution":{"iopub.status.busy":"2023-01-13T12:39:35.531141Z","iopub.execute_input":"2023-01-13T12:39:35.531522Z","iopub.status.idle":"2023-01-13T12:39:35.548405Z","shell.execute_reply.started":"2023-01-13T12:39:35.531467Z","shell.execute_reply":"2023-01-13T12:39:35.547482Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def normalize_img(in_img):\n    return (255 * (in_img - in_img.min()) / (in_img.max() - in_img.min())).astype(cp.uint8)\n\ndef mux(condition, true_case, false_case):\n    neg_condition = condition ^ True\n    return condition * true_case + neg_condition * false_case\n\ndef common_pipeline(imgs, is_monochromes, p_ids, i_ids):\n    # to lower information loss\n    imgs = fn.resize(imgs, resize_y = int(Config.resize_height_to * Config.factor), interp_type = types.INTERP_NN)\n    \n    # normalize\n    imgs = fn.cast(imgs, dtype = types.INT32) # weird bug with python_function and unsigned ints\n    imgs = fn.python_function(imgs, device = 'gpu', function = normalize_img, batch_processing = False)\n\n    # monochrome invert using mux\n    imgs_inverted = 255 - imgs\n    imgs = mux(is_monochromes, imgs_inverted, imgs)\n    imgs = fn.cast(imgs, dtype = types.UINT8)\n    \n    # pad to the right to be square\n    imgs = fn.pad(imgs, axes = (1,), shape = (int(Config.resize_height_to * Config.factor)), )\n    \n    # resize for model\n    imgs_2roi = fn.resize(imgs, resize_y = Config.resize_height_to) / 255\n    \n    return imgs_2roi, imgs, p_ids, i_ids","metadata":{"_kg_hide-input":true,"execution":{"iopub.status.busy":"2023-01-13T12:39:35.549781Z","iopub.execute_input":"2023-01-13T12:39:35.55019Z","iopub.status.idle":"2023-01-13T12:39:35.5607Z","shell.execute_reply.started":"2023-01-13T12:39:35.550156Z","shell.execute_reply":"2023-01-13T12:39:35.559801Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"@pipeline_def\ndef j2k_decode_pipeline():\n    imgs, is_monochromes, p_ids, i_ids = fn.external_source(\n        source = J2Ki, num_outputs = 4, device = \"cpu\", dtype = [types.UINT8, types.BOOL, types.INT64, types.INT64],\n    )\n    \n    imgs = fn.experimental.decoders.image(\n        imgs, device = 'mixed', output_type = types.ANY_DATA, dtype = DALIDataType.UINT16\n    )\n    \n    return common_pipeline(imgs, is_monochromes, p_ids, i_ids)\n\n@pipeline_def\ndef jll_process_pipeline():\n    imgs, is_monochromes, p_ids, i_ids = fn.external_source(\n        source = JLLi, num_outputs = 4, device = \"gpu\", dtype = [types.UINT16, types.BOOL, types.INT64, types.INT64],\n    )\n    \n    imgs = fn.reinterpret(imgs, layout = \"HWC\")\n    \n    return common_pipeline(imgs, is_monochromes, p_ids, i_ids)","metadata":{"_kg_hide-input":true,"execution":{"iopub.status.busy":"2023-01-13T12:39:35.564417Z","iopub.execute_input":"2023-01-13T12:39:35.564721Z","iopub.status.idle":"2023-01-13T12:39:35.574747Z","shell.execute_reply.started":"2023-01-13T12:39:35.564697Z","shell.execute_reply":"2023-01-13T12:39:35.573872Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"pipes = [\n    (j2k_decode_pipeline, len(J2Ki)),\n    (jll_process_pipeline, len(JLLi))\n]","metadata":{"execution":{"iopub.status.busy":"2023-01-13T12:39:35.576068Z","iopub.execute_input":"2023-01-13T12:39:35.5765Z","iopub.status.idle":"2023-01-13T12:39:35.584871Z","shell.execute_reply.started":"2023-01-13T12:39:35.576441Z","shell.execute_reply":"2023-01-13T12:39:35.583819Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def _p(d):\n    if len(d) > 0 and len(d[0]) >= 1:\n        return d[0].numpy()\n    return np.array([0, 0, Config.resize_height_to, Config.resize_height_to, 1, 0])\n\nfor pipe_l, p_len in pipes:\n    print(\"Starting pipeline:\", pipe_l.__name__)\n    \n    if p_len == 0:\n        continue\n    \n    pipe = pipe_l(**Config.pipeline_args)\n    pipe.start_py_workers()\n    dali_iter = DALIGenericIterator(pipe, ['imgs_2roi', 'imgs', 'p_ids', 'i_ids'])\n    \n    for i, data in tqdm(enumerate(dali_iter), total = p_len):\n        imgs_2roi = data[0]['imgs_2roi']\n        imgs = data[0]['imgs'].cpu().numpy()\n        p_ids = data[0]['p_ids'].cpu().numpy()\n        i_ids = data[0]['i_ids'].cpu().numpy()\n\n        with torch.no_grad():\n            y = model(imgs_2roi.expand(-1, -1, -1, 3).permute(0, 3, 1, 2)).cpu()\n            dets = non_max_suppression(y, conf_thres = 0.1, classes = [0, 1])\n            # dets -> 32 x (x y x y conf cls)\n\n            # keep best \n            dets = np.array([\n                _p(d) for d in dets\n            ])\n\n        cls_ = dets[:, 4].astype(np.bool_)\n        xyxy = dets[:, [0, 1, 2, 3]] / 640\n        xyxy = np.clip(xyxy, 0, 1)\n        \n        for img, (xmin, ymin, xmax, ymax), lr, p_id, i_id in zip(imgs, xyxy, cls_, p_ids, i_ids):\n            h, w, _ = img.shape\n            img = np.squeeze(img)\n\n            xmin = (xmin * w).astype(np.uint16)\n            xmax = (xmax * w).astype(np.uint16)\n            ymin = (ymin * h).astype(np.uint16)\n            ymax = (ymax * h).astype(np.uint16)\n\n            if lr:\n                img = np.fliplr(img)\n                xmin = w - xmin\n                xmax = w - xmax\n                xmin, xmax = xmax, xmin\n\n            resulting_crop_h = ymax - ymin\n            resulting_crop_w = xmax - xmin\n\n            if (2 * resulting_crop_w) < resulting_crop_h: # increase width\n                # how much more width do I need to add?\n                needed_width = (resulting_crop_h // 2) - resulting_crop_w\n\n                # can I meet width without going past image?\n                if needed_width + xmax > w:\n                    # padding/resizing is required (I use padding, in this case)\n                    difference = xmax + needed_width - w\n\n                    img = cv2.copyMakeBorder(img, 0, 0, difference, 0, cv2.BORDER_CONSTANT, value = 0)\n\n                    # no need to offset bbox since origin is top left\n\n                # expand bbox by needed width\n                xmax += needed_width\n\n            elif (2 * resulting_crop_w) > resulting_crop_h: # increase height\n                # how much more height do I need to add?\n                needed_height = (resulting_crop_w * 2) - resulting_crop_h\n\n                # can I meet height without going past image?\n                if ymin - needed_height < 0:\n                    # padding/resizing is required (I use padding in this case)\n                    difference = needed_height - ymin\n\n                    img = cv2.copyMakeBorder(img, difference, 0, 0, 0, cv2.BORDER_CONSTANT, value = 0)\n\n                    # offset bbox values to new origin\n                    ymin += difference\n                    ymax += difference\n\n                ymin -= needed_height\n\n            if lr:\n                img = np.fliplr(img)\n                h, w = img.shape\n\n                xmin = w - xmin\n                xmax = w - xmax\n                xmin, xmax = xmax, xmin\n\n            img = img[ymin:ymax, xmin:xmax]\n\n            h, w = img.shape\n\n            img = cv2.resize(\n                img, dsize = (Config.out_img_size_x, Config.out_img_size_y), interpolation=cv2.INTER_NEAREST\n            )\n            \n            if(Config.preserve_folder_structure):\n                os.makedirs(f'{Config.oimg_dir}/{p_id[0]}/', exist_ok = True)\n                cv2.imwrite(f'{Config.oimg_dir}/{p_id[0]}/{i_id[0]}.png', img) # {img_dir}\n            else:\n                cv2.imwrite(f'{Config.oimg_dir}/{p_id[0]}_{i_id[0]}.png', img) # {img_dir}\n            \n    del img, imgs, xyxy, cls_, p_ids, i_ids, y, dets, data, dali_iter, pipe\n    gc.collect()\n    \ndel pipes, model\n_ = gc.collect()","metadata":{"execution":{"iopub.status.busy":"2023-01-13T12:39:35.586428Z","iopub.execute_input":"2023-01-13T12:39:35.586752Z","iopub.status.idle":"2023-01-13T12:39:44.352038Z","shell.execute_reply.started":"2023-01-13T12:39:35.586726Z","shell.execute_reply":"2023-01-13T12:39:44.350996Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"!ls /tmp/imgs","metadata":{"execution":{"iopub.status.busy":"2023-01-13T12:39:44.353877Z","iopub.execute_input":"2023-01-13T12:39:44.354602Z","iopub.status.idle":"2023-01-13T12:39:45.387572Z","shell.execute_reply.started":"2023-01-13T12:39:44.354564Z","shell.execute_reply":"2023-01-13T12:39:45.386348Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#image = cv2.imread('/tmp/imgs/10008_1591370361.png')\n\n#plt.imshow(image)\n#plt.show()","metadata":{"execution":{"iopub.status.busy":"2023-01-13T12:39:45.389434Z","iopub.execute_input":"2023-01-13T12:39:45.389789Z","iopub.status.idle":"2023-01-13T12:39:45.595604Z","shell.execute_reply.started":"2023-01-13T12:39:45.38975Z","shell.execute_reply":"2023-01-13T12:39:45.594608Z"},"trusted":true},"execution_count":null,"outputs":[]}]}