{"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":"nvidiaTeslaT4","dataSources":[{"sourceId":61446,"databundleVersionId":6962461,"sourceType":"competition"},{"sourceId":7183197,"sourceType":"datasetVersion","datasetId":4152203}],"dockerImageVersionId":30683,"isInternetEnabled":true,"language":"python","sourceType":"notebook","isGpuEnabled":true}},"nbformat_minor":4,"nbformat":4,"cells":[{"cell_type":"code","source":"import torch\nimport torch.nn as nn\nimport torch.nn.functional as F\nimport numpy as np\nimport cv2\nimport math\n\nimport timm\nprint('timm', timm.__version__)\nfrom timm.models.resnet import *\n\nimport matplotlib.pyplot as plt\nfrom IPython.display import Image, display\n\n###########################################33\n#helper\ndef image_show_norm(name, image, min=None, max=None, type='bgr', resize=1):\n\tif max is None: max = image.max()\n\tif min is None: min = image.min()\n\tif type == 'rgb': image = np.ascontiguousarray(image[:, :, ::-1])\n\n\tH, W = image.shape[0:2]\n\tcv2.namedWindow(name, cv2.WINDOW_GUI_NORMAL)  # WINDOW_NORMAL\n\tcv2.imshow(name, (np.clip((image - min) / (max - min), 0, 1) * 255).astype(np.uint8))\n\tcv2.resizeWindow(name, round(resize * W), round(resize * H))\n\n\n\ndef np_sigmoid(x):\n\treturn 1 / (1 + np.exp(-x))\n\nprint('IMPORT OK')","metadata":{"_uuid":"8f2839f25d086af736a60e9eeb907d3b93b6e0e5","_cell_guid":"b1076dfc-b9ad-4769-8c92-a6c4dae69d19","execution":{"iopub.status.busy":"2024-04-18T16:48:03.212973Z","iopub.execute_input":"2024-04-18T16:48:03.213321Z","iopub.status.idle":"2024-04-18T16:48:11.061024Z","shell.execute_reply.started":"2024-04-18T16:48:03.213294Z","shell.execute_reply":"2024-04-18T16:48:11.060105Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"!pip install torch torchvision torchaudio -f https://download.pytorch.org/whl/cu11.2/torch_stable.html\n","metadata":{"execution":{"iopub.status.busy":"2024-04-18T16:48:11.062797Z","iopub.execute_input":"2024-04-18T16:48:11.063081Z","iopub.status.idle":"2024-04-18T16:48:24.378234Z","shell.execute_reply.started":"2024-04-18T16:48:11.063057Z","shell.execute_reply":"2024-04-18T16:48:24.377231Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"root_dir = \\\n\t'/kaggle/input/blood-vessel-segmentation'\n\t#'/home/user/share1/kaggle/2023/blood-vessel-segmentation'\n\n\nIMAGE_SIZE  = 640\nIMAGE_DEPTH = 32\nclass MyDecoderBlock3d(nn.Module):\n\tdef __init__(\n\t\tself,\n\t\tin_channel,\n\t\tskip_channel,\n\t\tout_channel,\n\t):\n\t\tsuper().__init__()\n\t\tself.conv1 = nn.Sequential(\n\t\t\tnn.Conv3d(in_channel + skip_channel, out_channel, kernel_size=3, padding=1, bias=False),\n\t\t\tnn.BatchNorm3d(out_channel, eps=1e-4),\n\t\t\tnn.ReLU(inplace=True),\n\t\t)\n\t\tself.attention1 = nn.Identity()\n\t\tself.conv2 = nn.Sequential(\n\t\t\tnn.Conv3d(out_channel,out_channel,kernel_size=3, padding=1, bias=False),\n\t\t\tnn.BatchNorm3d(out_channel, eps=1e-4),\n\t\t\tnn.ReLU(inplace=True),\n\t\t)\n\t\tself.attention2 = nn.Identity()\n\n\tdef forward(self, x, skip=None):\n\t\tx = F.interpolate(x, scale_factor=2, mode='nearest')\n\t\tif skip is not None:\n\t\t\tx = torch.cat([x, skip], dim=1)\n\t\t\tx = self.attention1(x)\n\t\tx = self.conv1(x)\n\t\tx = self.conv2(x)\n\t\tx = self.attention2(x)\n\t\treturn x\n    \n    \nroot_dir = \\\n\t'/kaggle/input/blood-vessel-segmentation'\n\t#'/home/user/share1/kaggle/2023/blood-vessel-segmentation'\n\n\nIMAGE_SIZE  = 640\nIMAGE_DEPTH = 32\n\n# Net\nclass MyDecoderBlock3d(nn.Module):\n\tdef __init__(\n\t\tself,\n\t\tin_channel,\n\t\tskip_channel,\n\t\tout_channel,\n\t):\n\t\tsuper().__init__()\n\t\tself.conv1 = nn.Sequential(\n\t\t\tnn.Conv3d(in_channel + skip_channel, out_channel, kernel_size=3, padding=1, bias=False),\n\t\t\tnn.BatchNorm3d(out_channel, eps=1e-4),\n\t\t\tnn.ReLU(inplace=True),\n\t\t)\n\t\tself.attention1 = nn.Identity()\n\t\tself.conv2 = nn.Sequential(\n\t\t\tnn.Conv3d(out_channel,out_channel,kernel_size=3, padding=1, bias=False),\n\t\t\tnn.BatchNorm3d(out_channel, eps=1e-4),\n\t\t\tnn.ReLU(inplace=True),\n\t\t)\n\t\tself.attention2 = nn.Identity()\n\n\tdef forward(self, x, skip=None):\n\t\tx = F.interpolate(x, scale_factor=2, mode='nearest')\n\t\tif skip is not None:\n\t\t\tx = torch.cat([x, skip], dim=1)\n\t\t\tx = self.attention1(x)\n\t\tx = self.conv1(x)\n\t\tx = self.conv2(x)\n\t\tx = self.attention2(x)\n\t\treturn x\n\nclass MyUnetDecoder3d(nn.Module):\n\tdef __init__(self,\n\t\t\t\t in_channel,\n\t\t\t\t skip_channel,\n\t\t\t\t out_channel,\n\t\t\t\t ):\n\t\tsuper().__init__()\n\t\tself.center = nn.Identity()\n\n\t\ti_channel = [in_channel, ] + out_channel[:-1]\n\t\ts_channel = skip_channel\n\t\to_channel = out_channel\n\t\tblock = [\n\t\t\tMyDecoderBlock3d(i, s, o)\n\t\t\tfor i, s, o in zip(i_channel, s_channel, o_channel)\n\t\t]\n\t\tself.block = nn.ModuleList(block)\n\n\tdef forward(self, feature, skip):\n\t\td = self.center(feature)\n\t\tdecode = []\n\t\tfor i, block in enumerate(self.block):\n\t\t\t#print(i, d.shape, skip[i].shape if skip[i] is not None else 'none')\n\t\t\t#print(block.conv1[0])\n\t\t\t#print('')\n\n\t\t\ts = skip[i]\n\t\t\td = block(d, s)\n\t\t\tdecode.append(d)\n\t\tlast = d\n\t\treturn last, decode\n    \nclass Net(nn.Module):\n\tdef __init__(self, ):\n\t\tsuper().__init__()\n\t\tencoder_dim = [64, 256, 512, 1024, 2048]\n\t\tdecoder_dim = [256, 128, 128, 64, 32 ]\n\n\t\t#self.encoder = seresnext26d_32x4d(pretrained=True, in_chans=3)\n\t\tself.encoder = resnet50d(pretrained=False, in_chans=3)\n\t\tself.add_conv1 = nn.Sequential(\n\t\t\tnn.Conv3d(1, 32, kernel_size=3, padding=1, bias=False),\n\t\t\tnn.BatchNorm3d(32, eps=1e-4),\n\t\t\tnn.ReLU(inplace=True),\n\t\t\tnn.Conv3d(32, 32, kernel_size=3, padding=1, bias=False),\n\t\t\tnn.BatchNorm3d(32, eps=1e-4),\n\t\t\tnn.ReLU(inplace=True),\n\t\t\tnn.Conv3d(32, 32, kernel_size=3, padding=1, bias=False),\n\t\t\tnn.BatchNorm3d(32, eps=1e-4),\n\t\t\tnn.ReLU(inplace=True),\n\t\t)\n\n\n\t\tself.decoder = MyUnetDecoder3d(\n\t\t\tin_channel  = encoder_dim[-1],\n\t\t\tskip_channel= encoder_dim[:-1][::-1]+[32],\n\t\t\tout_channel = decoder_dim,\n\t\t)\n\t\tself.vessel = nn.Conv3d(decoder_dim[-1], 1, kernel_size=1)\n\n\t\t#just a simple demo. please improve this\n\t\tself.convert_3d = torch.nn.ModuleList([\n\t\t\t\tnn.Identity(),\n\t\t\t\tnn.Sequential(nn.Conv3d(  64,  64,kernel_size=( 2,1,1),stride=( 2,1,1),padding=(0,0,0)),  nn.BatchNorm3d(  64, eps=1e-4), nn.ReLU(inplace=True), ),\n\t\t\t\tnn.Sequential(nn.Conv3d( 256, 256,kernel_size=( 4,1,1),stride=( 4,1,1),padding=(0,0,0)),  nn.BatchNorm3d( 256, eps=1e-4), nn.ReLU(inplace=True), ),\n\t\t\t\tnn.Sequential(nn.Conv3d( 512, 512,kernel_size=( 8,1,1),stride=( 8,1,1),padding=(0,0,0)),  nn.BatchNorm3d( 512, eps=1e-4), nn.ReLU(inplace=True), ),\n\t\t\t\tnn.Sequential(nn.Conv3d(1024,1024,kernel_size=(16,1,1),stride=(16,1,1),padding=(0,0,0)),  nn.BatchNorm3d(1024, eps=1e-4), nn.ReLU(inplace=True), ),\n\t\t\t\tnn.Sequential(nn.Conv3d(2048,2048,kernel_size=(32,1,1),stride=(32,1,1),padding=(0,0,0)),  nn.BatchNorm3d(2048, eps=1e-4), nn.ReLU(inplace=True), ),\n\t\t])\n\n\tdef forward(self, subvolume):\n\t\txx = subvolume\n\t\tB, D, H, W = xx.shape\n\n\t\tx = xx.reshape(B * D, 1, H, W)\n\t\tx = x.expand(-1, 3, -1, -1)\n\n\t\tencode = []\n\t\txx = self.add_conv1(xx.unsqueeze(1)); encode.append(xx)\n\n\t\te = self.encoder\n\t\tx = e.conv1(x)\n\t\tx = e.bn1(x)\n\t\tx = e.act1(x);   encode.append(x)\n\t\tx = F.avg_pool2d(x, kernel_size=2, stride=2)\n\t\tx = e.layer1(x); encode.append(x)\n\t\tx = e.layer2(x); encode.append(x)\n\t\tx = e.layer3(x); encode.append(x)\n\t\tx = e.layer4(x); encode.append(x)\n\n\t\tfor i, x in enumerate(encode):\n\t\t\tif i == 0: continue\n\t\t\tx = encode[i]\n\t\t\tBD, c, h, w = x.shape\n\t\t\tx = x.reshape(B, D, c, h, w)\n\t\t\tx = x.transpose(1, 2)\n\t\t\tx = self.convert_3d[i](x)\n\t\t\tencode[i] = x\n\n\t\tlast, decode = self.decoder(\n\t\t\tfeature=encode[-1], skip=encode[:-1][::-1]\n\t\t)\n\n\t\tvessel = self.vessel(last).squeeze(1)\n\t\tvessel = torch.sigmoid(vessel.float())\n\t\treturn vessel\n\n#############################################################################33\n\ndef run_check_net():\n\theight, width = 480, 480\n\tdepth = 32\n\tbatch_size = 2\n\n\tsubvolume = torch.from_numpy(np.random.uniform(0, 1, (batch_size, depth, height, width))).float().cuda()\n\tnet = Net().cuda()\n\n\n\twith torch.no_grad():\n\t\twith torch.cuda.amp.autocast(enabled=True):\n\t\t\tvessel = net(subvolume)\n\tprint(subvolume.shape)\n\tprint(vessel.shape)\n\nrun_check_net()\nprint('NET OK!!!')","metadata":{"execution":{"iopub.status.busy":"2024-04-18T16:48:24.380233Z","iopub.execute_input":"2024-04-18T16:48:24.380620Z","iopub.status.idle":"2024-04-18T16:48:29.507060Z","shell.execute_reply.started":"2024-04-18T16:48:24.380585Z","shell.execute_reply":"2024-04-18T16:48:29.506139Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def norm_by_min_max(x, xmin, xmax, alpha=0.01):\n\txmin = float(xmin)\n\txmax = float(xmax)\n\tx = (x-xmin)/(xmax-xmin)\n\tif 1:\n\t\tx[x>1]=(x[x>1]-1)*alpha +1\n\t\tx[x<0]=(x[x<0])*alpha\n\t#x = np.clip(x,0,1)\n\treturn x\n\ndef load_dummy_data():\n\t#x,y,z = 740,470, 389\n\tname = 'kidney_3_sparse'\n\tx,y,z = 420, 150, 356\n\n\tsubvolume = []\n\ttruth = []\n\tfor i in range(z,z+IMAGE_DEPTH):\n\t\timage_file = f'{root_dir}/train/{name}/images/{i:04d}.tif'\n\t\tm = cv2.imread(image_file, cv2.IMREAD_UNCHANGED)\n\t\tsubvolume.append(m)\n\n\t\timage_file = f'{root_dir}/train/{name}/labels/{i:04d}.tif'\n\t\tv = cv2.imread(image_file, cv2.IMREAD_GRAYSCALE)\n\t\ttruth.append(v)\n\n\tsubvolume = np.stack(subvolume)\n\tsubvolume = subvolume[:,y:y+IMAGE_SIZE,x:x+IMAGE_SIZE]\n\txmin, xmax = (18806, 21903) #precompute low and high percentile\n\tsubvolume = norm_by_min_max(subvolume, xmin, xmax)\n\n\ttruth = np.stack(truth)\n\ttruth = truth[:,y:y+IMAGE_SIZE,x:x+IMAGE_SIZE]\n\ttruth = truth//255\n\n\treturn subvolume, truth","metadata":{"execution":{"iopub.status.busy":"2024-04-18T16:48:29.509361Z","iopub.execute_input":"2024-04-18T16:48:29.509660Z","iopub.status.idle":"2024-04-18T16:48:29.518841Z","shell.execute_reply.started":"2024-04-18T16:48:29.509635Z","shell.execute_reply":"2024-04-18T16:48:29.517999Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def run_demo():\n\n\tsubvolume, truth = load_dummy_data()\n\n\tnet = Net()\n\tcheckpoint_file ='/kaggle/input/2d-to-3d-demo-data/00000828.pth'\n\tstate_dict = torch.load(checkpoint_file, map_location=lambda storage, loc: storage)['state_dict']\n\tprint(net.load_state_dict(state_dict, strict=True))  # True\n\tnet = net.eval()\n\tnet = net.cuda()\n\n\ttensor = torch.from_numpy(subvolume).float().cuda()\n\twith torch.no_grad():\n\t\twith torch.cuda.amp.autocast(enabled=True):\n\t\t\tvessel = net(tensor.unsqueeze(0)).squeeze(0)\n\tvessel = vessel.float().data.cpu().numpy()\n\n\n\t#--------------------------------------------------------------\n\tif 0: #visualisation (run offline. does not work in kaggle\n\n\t\tD,H,W = subvolume.shape\n\t\tsubvolume = np.clip(subvolume,0,1)\n\t\tx_mean = subvolume.mean(0)\n\t\tp_mean = vessel.mean(0)\n\t\ty_mean = truth.mean(0)\n\n\t\t#---\n\t\tdef add_contrast_for_x(x_mean):\n\t\t\tx_mean = np_sigmoid((1-x_mean)*5)\n\t\t\tx_mean = 1- (x_mean-x_mean.min())/(x_mean.max()-x_mean.min())\n\t\t\treturn x_mean\n\n\t\tdef add_contrast_for_py(p_mean, y_mean):\n\t\t\tzero = np.zeros((H,W))\n\t\t\tpy_mean = np.dstack([zero,p_mean,y_mean])\n\t\t\tpy_mean = ((py_mean-py_mean.min())/(py_mean.max()-py_mean.min()))**0.25\n\t\t\treturn py_mean\n\t\t#---\n\n\t\tx_mean  = add_contrast_for_x(x_mean)\n\t\tpy_mean = add_contrast_for_py(p_mean, y_mean)\n\t\timage_show_norm('x_mean',x_mean, min=0,max=1,resize=1)\n\t\timage_show_norm('py_mean',py_mean, min=0,max=1,resize=1)\n\t\tcv2.waitKey(0)\n\telse:\n\t\tdisplay(Image(filename='../input/2d-to-3d-demo-data/Selection_999(4367).png'))\n\t\tdisplay(Image(filename='../input/2d-to-3d-demo-data/Selection_999(4368).png'))\n \n\n\t#--------------------------------------------------------------\n\tif 0: #3d visualisation (run offline. does not work in kaggle)\n\t\timport pyvista as pv\n\t\ty = truth>0\n\t\tp = vessel>0.5\n\t\thit = y*p\n\t\tfp = (1-y)*p\n\t\tmiss = y*(1-p)\n\n\t\tpl = pv.Plotter()\n\n\t\tmhit  = pv.PolyData(np.stack(np.where(hit > 0.1)).T).glyph(geom=pv.Cube())\n\t\tmfp   = pv.PolyData(np.stack(np.where(fp  > 0.1)).T).glyph(geom=pv.Cube())\n\t\tmmiss = pv.PolyData(np.stack(np.where(miss > 0.1)).T).glyph(geom=pv.Cube())\n\t\tpl.add_mesh(mhit,  color='yellow' )\n\t\tpl.add_mesh(mfp,   color='green' )\n\t\tpl.add_mesh(mmiss, color='red' )\n\t\tpl.show()\n\telse:\n\t\tprint('why cannot play GIF ???????')        \n\t\tdisplay(Image(filename='../input/2d-to-3d-demo-data/Peek 2023-12-12 02-38.gif', format='png'))\n\t\tdisplay(Image(filename='../input/2d-to-3d-demo-data/Peek 2023-12-12 02-36.gif', format='png'))\n    \n\n#run_demo()\nprint('DEMO OK!!!')","metadata":{"execution":{"iopub.status.busy":"2024-04-18T16:48:29.519961Z","iopub.execute_input":"2024-04-18T16:48:29.520234Z","iopub.status.idle":"2024-04-18T16:48:29.538945Z","shell.execute_reply.started":"2024-04-18T16:48:29.520211Z","shell.execute_reply":"2024-04-18T16:48:29.538115Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"import torch\nprint(torch.__version__)\n","metadata":{"execution":{"iopub.status.busy":"2024-04-18T16:48:29.539904Z","iopub.execute_input":"2024-04-18T16:48:29.540182Z","iopub.status.idle":"2024-04-18T16:48:29.554288Z","shell.execute_reply.started":"2024-04-18T16:48:29.540139Z","shell.execute_reply":"2024-04-18T16:48:29.553477Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"pip install --upgrade torch torchvision\n","metadata":{"execution":{"iopub.status.busy":"2024-04-18T16:48:29.555320Z","iopub.execute_input":"2024-04-18T16:48:29.555637Z","iopub.status.idle":"2024-04-18T16:50:45.799412Z","shell.execute_reply.started":"2024-04-18T16:48:29.555606Z","shell.execute_reply":"2024-04-18T16:50:45.798238Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"pip install connected-components-3d","metadata":{"execution":{"iopub.status.busy":"2024-04-18T16:50:45.801030Z","iopub.execute_input":"2024-04-18T16:50:45.801757Z","iopub.status.idle":"2024-04-18T16:50:58.723085Z","shell.execute_reply.started":"2024-04-18T16:50:45.801719Z","shell.execute_reply":"2024-04-18T16:50:58.722011Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"try:\n    import cc3d\nexcept:\n    #https://pypi.org/project/connected-components-3d/\n    #!pip install connected-components-3d\n\n    !ls /kaggle/input/installation-connected-components-3d\n    !pip install connected-components-3d --no-index --find-links=file:///kaggle/input/installation-connected-components-3d/\n\nimport cc3d\n\n###### \nimport sys, os\nsys.path.append('/kaggle/input/blood-vessel-segmentation-third-party')\nsys.path.append('/kaggle/input/blood-vessel-segmentation-01')\n\nfrom helper import *\n\nimport cv2\nimport pandas as pd\nfrom glob import glob\nimport numpy as np\nfrom skimage.filters import apply_hysteresis_threshold\n\nfrom timeit import default_timer as timer\nimport gc\n\nimport torch\nimport torch.nn as nn\nimport torch.nn.functional as F\n\nimport matplotlib\nimport matplotlib.pyplot as plt\n\nprint('IMPORT OK  !!!!')","metadata":{"execution":{"iopub.status.busy":"2024-04-18T16:50:58.724565Z","iopub.execute_input":"2024-04-18T16:50:58.724890Z","iopub.status.idle":"2024-04-18T16:51:00.032977Z","shell.execute_reply.started":"2024-04-18T16:50:58.724862Z","shell.execute_reply":"2024-04-18T16:51:00.031986Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"class DotDict(dict):\n    \"\"\"Dot notation access to dictionary attributes.\"\"\"\n    __getattr__ = dict.get\n    __setattr__ = dict.__setitem__\n    __delattr__ = dict.__delitem__\n\n# Now you can use DotDict to create your dictionary with dot notation access\ndata_dir = '/kaggle/input/blood-vessel-segmentation'\nDATA_META = DotDict({\n    'kidney_1_dense': DotDict(\n        name='kidney_1_dense',\n        image_no=(1000, 1000 + 1000),\n        image_dir=f'{data_dir}/train'\n    ),\n    \n    'kidney_3_dense': DotDict(\n        name='kidney_3_dense',\n        image_no=(496, 996 + 1),\n        image_dir=f'{data_dir}/train'\n    ),\n    'kidney_2': DotDict(\n        name='kidney_2',\n        image_no=None,\n        image_dir=f'{data_dir}/train'\n    ),\n    \n})\n\ncfg = DotDict({\n    'p_threshold': 0.5,\n    'cc_threshold': -1,\n    'use_tta': False,\n})\nmode = 'local' #'local'  # 'submit' #\n\nif 'local' in mode:\n\tvalid_meta = [ DATA_META['kidney_1_dense'],DATA_META['kidney_3_dense'], ]\n\t#valid_meta = [DATA_META['kidney_3_dense'], ]\nif 'submit' in mode:\n\t#valid_meta = [DATA_META['kidney_2'], DATA_META['kidney_6']]\n\tvalid_meta = [DATA_META['kidney_2']]\n\n\n## io input, etc function\ndef build_file_list(d):\n\tif d.image_no is not None:\n\t\td.file = [f'{d.image_dir}/{d.name.replace(\"kidney_3_dense\",\"kidney_3_sparse\")}/images/{i:04d}.tif' for i in range(*d.image_no)]\n\telse:\n\t\td.file = sorted(glob(f'{d.image_dir}/{d.name}/images/*.tif'))\n\n\nfor d in valid_meta:\n\tbuild_file_list(d)\nprint(valid_meta[0].name, '\\n', valid_meta[0].file[:5]) \nprint('MODE OK  !!!!')","metadata":{"execution":{"iopub.status.busy":"2024-04-18T16:51:00.035630Z","iopub.execute_input":"2024-04-18T16:51:00.036385Z","iopub.status.idle":"2024-04-18T16:51:00.050256Z","shell.execute_reply.started":"2024-04-18T16:51:00.036356Z","shell.execute_reply":"2024-04-18T16:51:00.049389Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def file_to_id(f):\n\ts = f.split('/')\n\treturn s[-3] + '_' + s[-1][:-4]\n\n\ndef load_volume(d):\n\tvolume = [\n\t\tcv2.imread(f, cv2.IMREAD_UNCHANGED) for f in d.file\n\t]\n\tvolume = np.stack(volume)\n\treturn volume\n\n\ndef load_truth(d):\n\ttruth = [\n\t\tcv2.imread(f.replace('/images/', '/labels/'), cv2.IMREAD_GRAYSCALE) for f in d.file\n\t]\n\ttruth = np.stack(truth)\n\ttruth = truth // 255\n\treturn truth\n\n\ndef norm_by_percentile(x, low=10, high=99.8, alpha=0.01):\n\txmin = np.percentile(x, low)\n\txmax = np.percentile(x, high)\n\tx = (x - xmin) / (xmax - xmin)\n\tif 1:\n\t\tx[x > 1] = (x[x > 1] - 1) * alpha + 1\n\t\tx[x < 0] = (x[x < 0]) * alpha\n\t# x = np.clip(x,0,1)\n\treturn x\ndef rle_encode(mask):\n\tpixel = mask.flatten()\n\tpixel = np.concatenate([[0], pixel, [0]])\n\trun = np.where(pixel[1:] != pixel[:-1])[0] + 1\n\trun[1::2] -= run[::2]\n\trle = ' '.join(str(r) for r in run)\n\tif rle == '':\n\t\trle = '1 0'\n\treturn rle\n\n\n# memory efficient\ndef apply_hysteresis_threshold_in_chuck(\n\t\tx, low, high, chunk_size=32\n):\n\tD, H, W = x.shape\n\tpredict = np.zeros((D, H, W), np.uint8)\n\tfor i in range(0, D, chunk_size // 2):\n\t\tprint(i)\n\t\tpredict[i:i + chunk_size] = np.maximum(\n\t\t\tapply_hysteresis_threshold(x[i:i + chunk_size], low, high),\n\t\t\tpredict[i:i + chunk_size]\n\t\t)\n\t\tpass\n\treturn predict\n\ndef make_dummy_submission():\n    submission_df = []\n    for d in valid_meta:\n        submission_df.append(\n            pd.DataFrame(data={\n                'id': [file_to_id(f) for f in d.file],\n                'rle': ['1 0'] * len(d.file),\n            })\n        )\n    submission_df = pd.concat(submission_df).reset_index(drop=True)\n    submission_df.to_csv('submission.csv', index=False)\n    return submission_df\n\n\n\nprint('DATASET OK  !!!!')","metadata":{"execution":{"iopub.status.busy":"2024-04-18T16:51:00.051648Z","iopub.execute_input":"2024-04-18T16:51:00.052009Z","iopub.status.idle":"2024-04-18T16:51:00.068596Z","shell.execute_reply.started":"2024-04-18T16:51:00.051977Z","shell.execute_reply":"2024-04-18T16:51:00.067727Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#from model_resnet50d_2d_3d import Net  # as Net0\nimport time\n\ndef time_to_str(t, units='sec'):\n    if units == 'sec':\n        if t < 1e-3:\n            return f\"{t * 1e6:.0f} µs\"\n        elif t < 1:\n            return f\"{t * 1e3:.0f} ms\"\n        else:\n            return f\"{t:.2f} s\"\n    elif units == 'min':\n        return f\"{t / 60:.2f} min\"\n    else:\n        raise ValueError(f\"Unsupported time unit: {units}\")\n\ncheckpoint_file = \\\n    '/kaggle/input/blood-vessel-segmentation-01/100-resnet50-2d-to-3d-03-00000828.pth'\n \nnet = Net()\nstate_dict = torch.load(checkpoint_file, map_location=lambda storage, loc: storage)['state_dict']\nprint(net.load_state_dict(state_dict, strict=False))  # True\nnet = net.eval()\nnet = net.cuda()\n\n\ndef do_submit():\n    submission_df =[]\n    #submission_df = make_dummy_submission()\n    #submission_df = submission_df.set_index('id')\n    #print(submission_df)\n    \n    for d in valid_meta:\n        #if d.name in ['kidney_1_dense','kidney_6']: continue\n\n        print('load_volume() ...')\n        volume = load_volume(d)\n        volume = norm_by_percentile(volume)\n        volume = volume.astype(np.float16)\n        D, H, W = volume.shape\n        print(volume.shape)\n\n        # ----\n        depth = 32\n        size  = 640 #640\n        #calculating no. of subvolumes\n        num_D = int(np.ceil((D-8)/(depth-8)))\n        num_H = int(np.ceil((H-32)/(size-32)))\n        num_W = int(np.ceil((W-32)/(size-32)))\n        num_subvolume = num_D*num_H*num_W\n        print('num_D,num_H,num_W', (num_D,num_H,num_W))\n        #generating subvolume coordinates\n        zz = np.linspace(0,D-depth,num_D).astype(int).tolist()\n        yy = np.linspace(0,H-size,num_H).astype(int).tolist()\n        xx = np.linspace(0,W-size,num_W).astype(int).tolist()\n\n\n        prob = np.zeros((D, H, W), dtype=np.float16)\n        prob_count = np.full((D, H, W), fill_value=0.01, dtype=np.float16) \n        start_timer = timer()\n        t=0\n        for z in zz:\n            for y in yy:\n                for x in xx:\n                    print('\\r', f'{t}/{num_subvolume} : {z, y, x} ', time_to_str(timer() - start_timer, 'min'), end='')\n                    t=t+1\n                    image = torch.from_numpy(\n                        volume[z:z + depth, y:y + size, x:x + size]).cuda().unsqueeze(0)\n\n                    vessel = 0\n                    counter = 0\n                    def tta1(x, k):\n                            if k== 0: return x\n                            if k== 1: return torch.flip(x, dims=[2,])  #tta\n                            if k==-1: return torch.flip(x, dims=[2,])  #undo tta\n                    def tta2(x, k):\n                            if k== 0: return x\n                            if k== 1: return torch.flip(x, dims=[3,]) \n                            if k==-1: return torch.flip(x, dims=[3,]) \n                    def tta3(x, k):\n                            if k== 0: return x\n                            if k== 1: return torch.rot90(x, k= 1, dims=[2,3])\n                            if k==-1: return torch.rot90(x, k=-1, dims=[2,3])\n\n                    with torch.cuda.amp.autocast(enabled=True):\n                        with torch.no_grad():\n                                for t1,t2,t3 in [\n                                    [0,0,0],\n                                    [1,0,0],\n                                    [0,1,0],\n                                    [1,1,0],\n                                    [0,0,1],\n                                    [1,0,1],\n                                    [0,1,1],\n                                    [1,1,1], ]:\n                                #    v, k = net(tta1(tta2(tta3(image,t3),t2),t1))\n                                    v = net(tta1(tta2(tta3(image,t3),t2),t1))\n                                    vessel += tta3(tta2(tta1(v,-t1),-t2),-t3)\n                                #    kidney += tta3(tta2(tta1(k,-t1),-t2),-t3)\n                                    counter += 1\n                    vessel = vessel/counter\n                    vessel = vessel.half().data.cpu().numpy() #probably memory leak?\n                    vessel = vessel.squeeze(0)\n\n                    prob [z:z + depth, y:y + size, x:x + size] += vessel\n                    prob_count[z:z + depth, y:y + size, x:x + size] += 1\n\n                    #--for debug\n                    if (t<=2) and (mode=='local'):\n                        image = image.float().data.cpu().numpy()\n                        image = image.squeeze(0)\n\n                        m = image.mean(0)\n                        v = vessel.mean(0)\n                        v = np.clip(v*3,0,1)\n                        \n                        \n                        #image_show_norm('m,p', np.hstack([m, v]), min=0, max=1, resize=1)\n                        #cv2.waitKey(1) \n                        plt.figure(figsize=(12,12))\n                        plt.imshow(np.hstack([m, v]),cmap='gray')\n                        plt.show()\n                    #--for debug\n                    \n                    #del image \n                    del vessel \n                    gc.collect()\n        print('')\n        #end of all subvolume \n        prob /= prob_count\n        if (mode == 'local'):\n            np.savez_compressed(f'prob.xyz{d.name}.npz', prob=(prob*255).astype(np.uint8))\n            v = prob.mean(0)\n            v = np.clip(v*10,0,1)**0.5\n            #image_show_norm('v', v, min=0, max=1, resize=1)\n            #cv2.waitKey(0)\n            plt.figure(figsize=(12,12))\n            plt.imshow(v,cmap='gray')\n            plt.show()\n          \n        del volume \n        gc.collect()\n        \n        predict = (prob > cfg.p_threshold).astype(np.uint8) \n        # post processing ---\n        if cfg.cc_threshold > 0:\n            predict = cc3d.dust(\n                predict,\n                connectivity=26,\n                threshold=cfg.cc_threshold,\n                in_place=False\n            )\n        # ---\n        #submission_df.loc[[file_to_id(f) for f in d.file],'rle']=  [rle_encode(p) for p in predict]\n        submission_df.append(\n            pd.DataFrame(data={\n                'id' : [file_to_id(f) for f in d.file],\n                'rle': [rle_encode(p) for p in predict],\n            })\n        )\n\n        del predict \n        del prob\n        del prob_count\n        gc.collect()\n\n    submission_df = pd.concat(submission_df)\n    submission_df = submission_df.reset_index(drop=True)\n    submission_df.to_csv('submission.csv', index=False)\n    print(submission_df)\n\n\n\n#glob_file = glob(f'{data_dir}/test/kidney_5/images/*.tif')\nglob_file = glob(f'{data_dir}/train/kidney_2/images/*.tif')\nif (mode == 'submit') and (len(glob_file) == 3):  # cannot do 3d cnn because too few test files\n    submission_df = make_dummy_submission()\n    print(submission_df)\nelse:\n    do_submit()\n\n\nprint('SUBMIT OK!!!')\n\ndef np_hit_fp_metric(predict, truth):\n    \"\"\"\n    Compute hit, false positive, true sum, and predicted sum for binary segmentation.\n\n    Parameters:\n        predict (np.ndarray): Binary prediction mask.\n        truth (np.ndarray): Ground truth mask.\n\n    Returns:\n        tuple: (hit, fp, t_sum, p_sum)\n    \"\"\"\n    hit = np.sum(np.logical_and(predict, truth))\n    fp = np.sum(np.logical_and(predict, np.logical_not(truth)))\n    t_sum = np.sum(truth)\n    p_sum = np.sum(predict)\n\n    return hit, fp, t_sum, p_sum\n\n\n\n\n\nif (mode == 'local'):\n    for d in valid_meta: \n        print(d.name)\n\n        print('load_truth truth ...')\n        truth = load_truth(d)\n        prob  = np.load(f'prob.xyz{d.name}.npz')['prob']/255\n\n        for th in [0.5, 0.4, ]:\n            predict = (prob > th)\n\n            hit, fp, t_sum, p_sum = np_hit_fp_metric(predict, truth)\n            lb_score = fast_compute_surface_dice_score_from_tensor(predict, truth)\n\n            print(checkpoint_file)\n            print('th=',th)\n            print('hit      :', hit / t_sum)\n            print('fp       :', fp / p_sum)\n            print('lb_score :', lb_score)\n            print('')","metadata":{"execution":{"iopub.status.busy":"2024-04-18T16:51:00.069982Z","iopub.execute_input":"2024-04-18T16:51:00.070381Z"},"trusted":true},"execution_count":null,"outputs":[]}]}