{"metadata":{"kernelspec":{"language":"python","display_name":"Python 3","name":"python3"},"language_info":{"name":"python","version":"3.10.12","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":6942858,"sourceType":"datasetVersion","datasetId":3987147},{"sourceId":7408599,"sourceType":"datasetVersion","datasetId":4308913},{"sourceId":7420945,"sourceType":"datasetVersion","datasetId":4317540},{"sourceId":7465479,"sourceType":"datasetVersion","datasetId":4345568},{"sourceId":7466757,"sourceType":"datasetVersion","datasetId":4346322},{"sourceId":7499683,"sourceType":"datasetVersion","datasetId":4367118},{"sourceId":7574620,"sourceType":"datasetVersion","datasetId":4410168},{"sourceId":7745413,"sourceType":"datasetVersion","datasetId":4527670},{"sourceId":7796452,"sourceType":"datasetVersion","datasetId":4564459},{"sourceId":162313529,"sourceType":"kernelVersion"}],"dockerImageVersionId":30635,"isInternetEnabled":false,"language":"python","sourceType":"notebook","isGpuEnabled":true}},"nbformat_minor":4,"nbformat":4,"cells":[{"cell_type":"code","source":"import torch as tc \nimport torch.nn as nn  \nimport numpy as np\nfrom tqdm import tqdm\nimport os,sys,cv2\nfrom torch.cuda.amp import autocast\nimport matplotlib.pyplot as plt\nimport albumentations as A\n# import segmentation_models_pytorch as smp\nfrom albumentations.pytorch import ToTensorV2\nfrom torch.utils.data import Dataset, DataLoader\nfrom torch.nn.parallel import DataParallel\nfrom glob import glob\nimport random\nimport torch\nimport pandas as pd\nimport argparse  \nimport sys\nfrom PIL import Image\nimport sys\n\nsys.path.append('/kaggle/input/sennetmetrics-no-smp-forcewithme/sennetmetrics')\nsys.path.append('/kaggle/input/sennetmetrics-no-smp-forcewithme/sennetmetrics/src')\nfrom sennet_metrices import *\nimport gc\nimport time\n\ndef rle_encode(mask):\n    pixel = mask.flatten()\n    pixel = np.concatenate([[0], pixel, [0]])\n    run = np.where(pixel[1:] != pixel[:-1])[0] + 1\n    run[1::2] -= run[::2]\n    rle = ' '.join(str(r) for r in run)\n    if rle == '':\n        rle = '1 0'\n    return rle\n\ndef seed_everything(seed):\n    random.seed(seed)\n    os.environ['PYTHONHASHSEED'] = str(seed)\n    np.random.seed(seed)\n    torch.manual_seed(seed)\n    torch.cuda.manual_seed(seed)\n    torch.backends.cudnn.deterministic = True\n    torch.backends.cudnn.benchmark = True\nseed_everything(42)\nCHOPPING_PER =1e-3","metadata":{"execution":{"iopub.status.busy":"2024-03-09T09:05:28.972069Z","iopub.execute_input":"2024-03-09T09:05:28.972371Z","iopub.status.idle":"2024-03-09T09:05:35.498513Z","shell.execute_reply.started":"2024-03-09T09:05:28.972320Z","shell.execute_reply":"2024-03-09T09:05:35.497678Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"import torch\nfrom torch.utils.data import DataLoader\nfrom transformers import ViTFeatureExtractor, ViTForImageClassification# , ViTForImageSegmentation\nimport pdb\nimport os\nimport sys\nimport torch\n","metadata":{"_uuid":"8f2839f25d086af736a60e9eeb907d3b93b6e0e5","_cell_guid":"b1076dfc-b9ad-4769-8c92-a6c4dae69d19","execution":{"iopub.status.busy":"2024-03-09T09:05:35.499970Z","iopub.execute_input":"2024-03-09T09:05:35.500413Z","iopub.status.idle":"2024-03-09T09:05:47.575260Z","shell.execute_reply.started":"2024-03-09T09:05:35.500386Z","shell.execute_reply":"2024-03-09T09:05:47.574484Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"work_dir = os.getcwd()\n# !git clone https://github.com/MIC-DKFZ/nnUNet.git\nbase_dir = '/kaggle/input/nnunet-zip-eps'\nos.chdir(base_dir)\nrespository_dir = os.path.join(base_dir,'nnUNet')\n# os.chdir(respository_dir)\nos.chdir('/kaggle/input/nnunet-zip-eps/nnUNet')\n!pip install -q -e .\n\nos.chdir('/kaggle/input/batchgenerators-helper/batchgenerators')\n!pip install -q -e .\n\nos.chdir('/kaggle/input/medpy-helper/medpy')\n!pip install -q -e .\n\n# os.chdir('/kaggle/input/segmentation-models-helper/segmentation_models.pytorch')\n# !pip install -q -e .\n\n# os.chdir('/kaggle/input/pretrained-helper/pretrained-models.pytorch')\n# !pip install -q -e .\n\nos.chdir(work_dir)\nsys.path.append(respository_dir)\nsys.path.append('/kaggle/input/batchgenerators-helper/batchgenerators')\nsys.path.append('/kaggle/input/medpy-helper/medpy')\n# sys.path.append('/kaggle/input/segmentation-models-helper/segmentation_models.pytorch')\n# sys.path.append('/kaggle/input/pretrained-helper/pretrained-models.pytorch')\n\n# !ls /kaggle/input\n# !unzip /kaggle/input/nnunet-repo-code","metadata":{"execution":{"iopub.status.busy":"2024-03-09T09:05:47.576260Z","iopub.execute_input":"2024-03-09T09:05:47.576849Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"from nnunet.training.network_training.nnUNetTrainerV2 import nnUNetTrainerV2\nimport pickle\ndevice = 'cuda' if torch.cuda.is_available() else 'cpu'\nwith open('/kaggle/input/task-008-hepaticvessel/model_final_checkpoint.model.pkl', 'rb') as f:\n    saved_checkpoint = pickle.load(f)\nsaved_model_nnUnet = torch.load('/kaggle/input/task-008-hepaticvessel/model_final_checkpoint.model', device)\n\ninit = saved_checkpoint['init']\nplans = saved_checkpoint['plans']","metadata":{"execution":{"iopub.status.busy":"2024-03-09T09:41:38.987650Z","iopub.execute_input":"2024-03-09T09:41:38.988104Z","iopub.status.idle":"2024-03-09T09:41:39.211277Z","shell.execute_reply.started":"2024-03-09T09:41:38.988064Z","shell.execute_reply":"2024-03-09T09:41:39.210468Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def build_model():\n    model = nnUNetTrainerV2(*init)\n    model.process_plans(plans)\n    model.was_initialized = True\n    model.log_file = 1\n    model.initialize_network()\n    model.initialize_optimizer_and_scheduler()\n    # best_checkpoint = torch.load('/kaggle/input/best-weights-alpha/best_epoch.bin', device)\n    best_checkpoint = torch.load('/kaggle/input/new-augmentations-weights/best_epoch_new_aug.bin', device)\n\n    model.network.load_state_dict(best_checkpoint)\n    model = model.network\n    return model","metadata":{"execution":{"iopub.status.busy":"2024-03-09T09:50:39.144811Z","iopub.execute_input":"2024-03-09T09:50:39.145877Z","iopub.status.idle":"2024-03-09T09:50:39.151546Z","shell.execute_reply.started":"2024-03-09T09:50:39.145836Z","shell.execute_reply":"2024-03-09T09:50:39.150476Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"import os\nimport random\nfrom tqdm import tqdm\nimport pandas as pd\nimport numpy as np\nfrom glob import glob\nimport gc\nimport time\nfrom collections import defaultdict\nimport  matplotlib.pyplot as plt\nfrom matplotlib.patches import Rectangle\nimport copy\nimport cv2\nimport torch.nn as nn\nfrom torch.utils.data import Dataset, DataLoader\nfrom torch.optim import lr_scheduler\nfrom torch.cuda import amp\nimport torch.optim as optim\nimport albumentations as A\nimport torch.nn.functional as F\n\n\nfrom colorama import Fore, Back, Style\nc_  = Fore.GREEN\nsr_ = Style.RESET_ALL","metadata":{"execution":{"iopub.status.busy":"2024-03-09T09:41:43.198265Z","iopub.execute_input":"2024-03-09T09:41:43.198710Z","iopub.status.idle":"2024-03-09T09:41:43.208688Z","shell.execute_reply.started":"2024-03-09T09:41:43.198671Z","shell.execute_reply":"2024-03-09T09:41:43.207734Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"from torch.nn.modules.loss import _Loss\nfrom typing import Optional, List\n\nBINARY_MODE: str = \"binary\"\n\nMULTICLASS_MODE: str = \"multiclass\"\n\nMULTILABEL_MODE: str = \"multilabel\"\ndef soft_dice_score(\n    output: torch.Tensor,\n    target: torch.Tensor,\n    smooth: float = 0.0,\n    eps: float = 1e-7,\n    dims=None,\n) -> torch.Tensor:\n    assert output.size() == target.size()\n    if dims is not None:\n        intersection = torch.sum(output * target, dim=dims)\n        cardinality = torch.sum(output + target, dim=dims)\n    else:\n        intersection = torch.sum(output * target)\n        cardinality = torch.sum(output + target)\n    dice_score = (2.0 * intersection + smooth) / (cardinality + smooth).clamp_min(eps)\n    return dice_score\n\nclass DiceLoss(_Loss):\n    def __init__(\n        self,\n        mode: str,\n        classes: Optional[List[int]] = None,\n        log_loss: bool = False,\n        from_logits: bool = True,\n        smooth: float = 0.0,\n        ignore_index: Optional[int] = None,\n        eps: float = 1e-7,\n    ):\n        \"\"\"Dice loss for image segmentation task.\n        It supports binary, multiclass and multilabel cases\n\n        Args:\n            mode: Loss mode 'binary', 'multiclass' or 'multilabel'\n            classes:  List of classes that contribute in loss computation. By default, all channels are included.\n            log_loss: If True, loss computed as `- log(dice_coeff)`, otherwise `1 - dice_coeff`\n            from_logits: If True, assumes input is raw logits\n            smooth: Smoothness constant for dice coefficient (a)\n            ignore_index: Label that indicates ignored pixels (does not contribute to loss)\n            eps: A small epsilon for numerical stability to avoid zero division error\n                (denominator will be always greater or equal to eps)\n\n        Shape\n             - **y_pred** - torch.Tensor of shape (N, C, H, W)\n             - **y_true** - torch.Tensor of shape (N, H, W) or (N, C, H, W)\n\n        Reference\n            https://github.com/BloodAxe/pytorch-toolbelt\n        \"\"\"\n        assert mode in {BINARY_MODE, MULTILABEL_MODE, MULTICLASS_MODE}\n        super(DiceLoss, self).__init__()\n        self.mode = mode\n        if classes is not None:\n            assert mode != BINARY_MODE, \"Masking classes is not supported with mode=binary\"\n            classes = to_tensor(classes, dtype=torch.long)\n\n        self.classes = classes\n        self.from_logits = from_logits\n        self.smooth = smooth\n        self.eps = eps\n        self.log_loss = log_loss\n        self.ignore_index = ignore_index\n\n    def forward(self, y_pred: torch.Tensor, y_true: torch.Tensor) -> torch.Tensor:\n\n        assert y_true.size(0) == y_pred.size(0)\n\n        if self.from_logits:\n            # Apply activations to get [0..1] class probabilities\n            # Using Log-Exp as this gives more numerically stable result and does not cause vanishing gradient on\n            # extreme values 0 and 1\n            if self.mode == MULTICLASS_MODE:\n                y_pred = y_pred.log_softmax(dim=1).exp()\n            else:\n                y_pred = F.logsigmoid(y_pred).exp()\n\n        bs = y_true.size(0)\n        num_classes = y_pred.size(1)\n        dims = (0, 2)\n\n        if self.mode == BINARY_MODE:\n            y_true = y_true.view(bs, 1, -1)\n            y_pred = y_pred.view(bs, 1, -1)\n\n            if self.ignore_index is not None:\n                mask = y_true != self.ignore_index\n                y_pred = y_pred * mask\n                y_true = y_true * mask\n\n        if self.mode == MULTICLASS_MODE:\n            y_true = y_true.view(bs, -1)\n            y_pred = y_pred.view(bs, num_classes, -1)\n\n            if self.ignore_index is not None:\n                mask = y_true != self.ignore_index\n                y_pred = y_pred * mask.unsqueeze(1)\n\n                y_true = F.one_hot((y_true * mask).to(torch.long), num_classes)  # N,H*W -> N,H*W, C\n                y_true = y_true.permute(0, 2, 1) * mask.unsqueeze(1)  # N, C, H*W\n            else:\n                y_true = F.one_hot(y_true, num_classes)  # N,H*W -> N,H*W, C\n                y_true = y_true.permute(0, 2, 1)  # N, C, H*W\n\n        if self.mode == MULTILABEL_MODE:\n            y_true = y_true.view(bs, num_classes, -1)\n            y_pred = y_pred.view(bs, num_classes, -1)\n\n            if self.ignore_index is not None:\n                mask = y_true != self.ignore_index\n                y_pred = y_pred * mask\n                y_true = y_true * mask\n\n        scores = self.compute_score(y_pred, y_true.type_as(y_pred), smooth=self.smooth, eps=self.eps, dims=dims)\n\n        if self.log_loss:\n            loss = -torch.log(scores.clamp_min(self.eps))\n        else:\n            loss = 1.0 - scores\n\n        # Dice loss is undefined for non-empty classes\n        # So we zero contribution of channel that does not have true pixels\n        # NOTE: A better workaround would be to use loss term `mean(y_pred)`\n        # for this case, however it will be a modified jaccard loss\n\n        mask = y_true.sum(dims) > 0\n        loss *= mask.to(loss.dtype)\n\n        if self.classes is not None:\n            loss = loss[self.classes]\n\n        return self.aggregate_loss(loss)\n\n    def aggregate_loss(self, loss):\n        return loss.mean()\n\n    def compute_score(self, output, target, smooth=0.0, eps=1e-7, dims=None) -> torch.Tensor:\n        return soft_dice_score(output, target, smooth, eps, dims)\n    \nclass SoftBCEWithLogitsLoss(nn.Module):\n\n    __constants__ = [\n        \"weight\",\n        \"pos_weight\",\n        \"reduction\",\n        \"ignore_index\",\n        \"smooth_factor\",\n    ]\n\n    def __init__(\n        self,\n        weight: Optional[torch.Tensor] = None,\n        ignore_index: Optional[int] = -100,\n        reduction: str = \"mean\",\n        smooth_factor: Optional[float] = None,\n        pos_weight: Optional[torch.Tensor] = None,\n    ):\n        \"\"\"Drop-in replacement for torch.nn.BCEWithLogitsLoss with few additions: ignore_index and label_smoothing\n\n        Args:\n            ignore_index: Specifies a target value that is ignored and does not contribute to the input gradient.\n            smooth_factor: Factor to smooth target (e.g. if smooth_factor=0.1 then [1, 0, 1] -> [0.9, 0.1, 0.9])\n\n        Shape\n             - **y_pred** - torch.Tensor of shape NxCxHxW\n             - **y_true** - torch.Tensor of shape NxHxW or Nx1xHxW\n\n        Reference\n            https://github.com/BloodAxe/pytorch-toolbelt\n\n        \"\"\"\n        super().__init__()\n        self.ignore_index = ignore_index\n        self.reduction = reduction\n        self.smooth_factor = smooth_factor\n        self.register_buffer(\"weight\", weight)\n        self.register_buffer(\"pos_weight\", pos_weight)\n\n    def forward(self, y_pred: torch.Tensor, y_true: torch.Tensor) -> torch.Tensor:\n        \"\"\"\n        Args:\n            y_pred: torch.Tensor of shape (N, C, H, W)\n            y_true: torch.Tensor of shape (N, H, W)  or (N, 1, H, W)\n\n        Returns:\n            loss: torch.Tensor\n        \"\"\"\n\n        if self.smooth_factor is not None:\n            soft_targets = (1 - y_true) * self.smooth_factor + y_true * (1 - self.smooth_factor)\n        else:\n            soft_targets = y_true\n\n        loss = F.binary_cross_entropy_with_logits(\n            y_pred,\n            soft_targets,\n            self.weight,\n            pos_weight=self.pos_weight,\n            reduction=\"none\",\n        )\n\n        if self.ignore_index is not None:\n            not_ignored_mask = y_true != self.ignore_index\n            loss *= not_ignored_mask.type_as(loss)\n\n        if self.reduction == \"mean\":\n            loss = loss.mean()\n\n        if self.reduction == \"sum\":\n            loss = loss.sum()\n\n        return loss\nimport gc\ndef clean_cuda():\n    torch.cuda.empty_cache()\n    gc.collect()\nclean_cuda()","metadata":{"execution":{"iopub.status.busy":"2024-03-09T09:41:43.795179Z","iopub.execute_input":"2024-03-09T09:41:43.795726Z","iopub.status.idle":"2024-03-09T09:41:44.106305Z","shell.execute_reply.started":"2024-03-09T09:41:43.795692Z","shell.execute_reply":"2024-03-09T09:41:44.105391Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"%matplotlib inline\nclass CFG:\n    seed          = 42 # *10\n    debug         = False # set debug=False for Full Training\n    exp_name      = 'nnUnet-pretrained-liver'\n    comment       = 'unet-efficientnet_b1-256x256'\n    output_dir    = './'\n    model_name    = 'Unet'\n    backbone      = 'efficientnet-b1'\n    train_bs      = 8\n    valid_bs      = 8\n    img_size      = [512, 512]\n    epochs        = 10 ## increase to 15\n    n_accumulate  = max(1, 64//train_bs)\n    lr            = 2e-3\n    scheduler     = 'CosineAnnealingLR'\n    min_lr        = 1e-6\n    T_max         = int(2279/(train_bs*n_accumulate)*epochs)+50\n    T_0           = 25\n    warmup_epochs = 0\n    wd            = 1e-6\n    n_fold        = 5\n    num_classes   = 1\n    device        = torch.device(\"cuda:0\" if torch.cuda.is_available() else \"cpu\")\n    rescale       = np.array([img_size[0]//2,img_size[0]//2,\n                              img_size[0]//2,img_size[0]//2])\n\n    gt_df = \"/kaggle/input/sennet-hoa-gt-data/gt.csv\"\n    data_root = \"/kaggle/input\"\n    train_groups = [\"kidney_1_dense\", 'kidney_1_voi', 'kidney_2', \"kidney_3_sparse\"]#, 'kidney_3_dense']\n    # train_groups = [\"kidney_1_dense\", \"kidney_3_sparse\"]#, 'kidney_3_dense']\n    valid_groups = ['kidney_3_dense'] # 501\n    # alid_groups = ['kidney_1_dense'] # 2279\n    loss_func     = \"DiceLoss\"\n\n    data_transforms = {\n        \"train\": A.Compose([\n            A.Resize(*img_size, interpolation=cv2.INTER_NEAREST),\n            A.HorizontalFlip(p=0.5),\n        ], p=1.0),\n        \n        \"valid\": A.Compose([\n            A.Resize(*img_size, interpolation=cv2.INTER_NEAREST),\n            ToTensorV2(),\n        ], p=1.0)\n    }\nvalid_aug_list = [\n    ToTensorV2(),\n]\nvalid_aug = A.Compose(valid_aug_list)","metadata":{"execution":{"iopub.status.busy":"2024-03-09T09:41:44.107942Z","iopub.execute_input":"2024-03-09T09:41:44.108291Z","iopub.status.idle":"2024-03-09T09:41:44.125392Z","shell.execute_reply.started":"2024-03-09T09:41:44.108263Z","shell.execute_reply":"2024-03-09T09:41:44.124423Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def set_seed(seed = 42):\n    '''Sets the seed of the entire notebook so results are the same every time we run.\n    This is for REPRODUCIBILITY.'''\n    np.random.seed(seed)\n    random.seed(seed)\n    torch.manual_seed(seed)\n    torch.cuda.manual_seed(seed)\n    # When running on the CuDNN backend, two further options must be set\n    torch.backends.cudnn.deterministic = True\n    torch.backends.cudnn.benchmark = False\n    # Set a fixed value for the hash seed\n    os.environ['PYTHONHASHSEED'] = str(seed)\nset_seed(CFG.seed)","metadata":{"execution":{"iopub.status.busy":"2024-03-09T09:41:44.324445Z","iopub.execute_input":"2024-03-09T09:41:44.325037Z","iopub.status.idle":"2024-03-09T09:41:44.330696Z","shell.execute_reply.started":"2024-03-09T09:41:44.325007Z","shell.execute_reply":"2024-03-09T09:41:44.329760Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def normalize_img(in_img, eps=1e-9):\n    min_ = in_img.min()\n    max_ = in_img.max()\n    return (255 * (in_img - min_) / (max_ - min_ + eps)).astype(np.uint8)\n\n\ndef add_noise(x:tc.Tensor,max_randn_rate=0.1,randn_rate=None,x_already_normed=False):\n    \"\"\"input.shape=(batch,f1,f2,...) output's var will be normalizate  \"\"\"\n    ndim=x.ndim-1\n    if x_already_normed:\n        x_std=tc.ones([x.shape[0]]+[1]*ndim,device=x.device,dtype=x.dtype)\n        x_mean=tc.zeros([x.shape[0]]+[1]*ndim,device=x.device,dtype=x.dtype)\n    else: \n        dim=list(range(1,x.ndim))\n        x_std=x.std(dim=dim,keepdim=True)\n        x_mean=x.mean(dim=dim,keepdim=True)\n    if randn_rate is None:\n        randn_rate=max_randn_rate*np.random.rand()*tc.rand(x_mean.shape,device=x.device,dtype=x.dtype)\n    cache=(x_std**2+(x_std*randn_rate)**2)**0.5\n    return (x-x_mean+tc.randn(size=x.shape,device=x.device,dtype=x.dtype)*randn_rate*x_std)/(cache+1e-7)\n\ndef filter_noise(x):\n    TH=x.reshape(-1)\n    index = -int(len(TH) * CHOPPING_PER)\n    TH:int = np.partition(TH, index)[index]\n    x[x>TH]=int(TH)\n    ########################################################################\n    TH=x.reshape(-1)\n    index = -int(len(TH) * CHOPPING_PER)\n    TH:int = np.partition(TH, -index)[-index]\n    x[x<TH]=int(TH)\n    return x","metadata":{"execution":{"iopub.status.busy":"2024-03-09T09:41:44.657242Z","iopub.execute_input":"2024-03-09T09:41:44.657568Z","iopub.status.idle":"2024-03-09T09:41:44.669392Z","shell.execute_reply.started":"2024-03-09T09:41:44.657541Z","shell.execute_reply":"2024-03-09T09:41:44.668394Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def get_indexes_along_axis0(num_slice, height, width, image_size, stride, slice_offset, nonempty_slices):\n    indexes = []\n    cur_h, cur_w = 0, 0\n    flag = 0\n    while(cur_h < height or cur_w < width):\n        if (cur_h + image_size <= height) and (cur_w + image_size <= width):\n            indexes.append((cur_h, cur_w))\n            cur_w += stride\n        elif (cur_h + image_size <= height) and (cur_w + image_size > width):\n            indexes.append((cur_h, width - image_size)) \n            cur_h += stride\n            cur_w = 0\n\n        if (cur_h + image_size > height):\n            if flag == 0:\n                cur_h = height - image_size \n                flag = 1\n            else:\n                break\n\n    indexes_3d = []\n    if nonempty_slices is None: slice_ls = range(slice_offset, num_slice-slice_offset)\n    else: slice_ls = nonempty_slices\n    for slice_idx in slice_ls:\n        for h, w in indexes:\n            indexes_3d.append((slice_idx, h, w))\n    return indexes_3d\n\ndef norm_with_clip(x:torch.Tensor,smooth=1e-5):\n    dim=list(range(1,x.ndim))\n    mean=x.mean(dim=dim,keepdim=True)\n    std=x.std(dim=dim,keepdim=True)\n    x=(x-mean)/(std+smooth)\n    x[x>5]=(x[x>5]-5)*1e-3 +5\n    x[x<-3]=(x[x<-3]+3)*1e-3-3\n    return x\n\ndef min_max_normalization(x:tc.Tensor)->tc.Tensor:\n    \"\"\"input.shape=(batch,f1,...)\"\"\"\n    shape=x.shape\n    if x.ndim>2:\n        x=x.reshape(x.shape[0],-1)\n    \n    min_=x.min(dim=-1,keepdim=True)[0]\n    max_=x.max(dim=-1,keepdim=True)[0]\n    if min_.mean()==0 and max_.mean()==1:\n        return x.reshape(shape)\n    \n    x=(x-min_)/(max_-min_+1e-9)\n    return x.reshape(shape)\n\ndef pad_hw(images, pad_h, pad_w):\n    pad_width = [(0, 0),  \n                 (pad_h, pad_h),  \n                 (pad_w, pad_w)]  \n    images = np.pad(images, pad_width=pad_width, mode='constant')\n    return images\n\nclass Dataset3D(Dataset):\n    def __init__(self, image_list, label_list=None, trans_axis=0,image_size=512, in_chans=1,stride=512,aug=False, pad_h=0, pad_w=0, nonempty_slices=None):\n        super(Dataset,self).__init__()\n\n        self.image_size=image_size\n        self.in_chans=in_chans\n        assert self.in_chans % 2 == 1\n        self.slice_offset = (self.in_chans - 1)//2\n        images = [cv2.imread(x,cv2.IMREAD_GRAYSCALE)[np.newaxis, :, :] for x in image_list]\n        images = np.concatenate(images, axis=0)## N, H, W\n\n        self.label_list = label_list\n        images = torch.tensor(filter_noise(images))\n        images = (min_max_normalization(images.to(tc.float16)[None])[0]*255).to(tc.uint8).numpy()\n        if label_list is not None:\n            labels = [cv2.imread(x,cv2.IMREAD_GRAYSCALE)[np.newaxis, :, :] for x in label_list]\n            labels = np.concatenate(labels, axis=0).astype(np.uint8)\n        else:\n            labels = None      \n\n        if trans_axis == 1:\n            images = np.transpose(images, (1,2,0))\n            if label_list is not None:\n                labels = np.transpose(labels, (1,2,0))\n        if trans_axis == 2:\n            images = np.transpose(images, (2,0,1))\n            if label_list is not None:\n                labels = np.transpose(labels, (2,0,1))\n        self.trans_axis = trans_axis       \n        images = pad_hw(images, pad_h, pad_w)\n        \n        self.images = [images]\n        self.labels = [labels]\n        num_slice, height, width = images.shape\n        slide_image_size = min(image_size, height-1, width-1)\n        slide_stride = min(stride, height-1, width-1)\n        self.coors = get_indexes_along_axis0(num_slice, height, width, slide_image_size, slide_stride, self.slice_offset, nonempty_slices)\n        self.idx_3d = [0] * len(self.coors)\n        self.transform=valid_aug\n\n    def __len__(self):\n        return len(self.coors)\n\n    def __getitem__(self,index):\n        n, h, w = self.coors[index]\n        idx_3d = self.idx_3d[index]\n        if self.slice_offset == 0:\n            image = self.images[idx_3d][n, h:h+self.image_size, w:w+self.image_size] # (1, H, W)\n        else:\n            image = self.images[idx_3d][n-self.slice_offset:n+self.slice_offset+1, h:h+self.image_size, w:w+self.image_size] # (3(5), H, W)\n            image = np.transpose(image, (1, 2, 0))            \n        data = self.transform(image=image)\n        return data['image'], idx_3d, n, h, w\n    \n    def get_labels(self):\n        return self.labels\n    \n\ndef permute_axis(pred, trans_axis=0, direction=0):\n    \"\"\"\n    direction: \n    0: from original shape to transposed shape\n    1: from transposed shape to original shape\n    \"\"\"\n    if trans_axis == 0: \n        return pred  ## N, H, W\n    elif trans_axis == 1:\n        if direction == 0:\n            pred = pred.permute(1, 2, 0) ## H, W, N\n        else:\n            pred = pred.permute(2, 0, 1)  ## N, H, W\n    elif trans_axis == 2:\n        if direction == 0:\n            pred = pred.permute(2, 0, 1) ## W, N, H\n        else:\n            pred = pred.permute(1, 2, 0) ## N, H, W\n    return pred  ## N, H, W\n\ndef get_nhw(trans_axis, NUM_SLICES, HEIGHT, WIDTH):\n    if trans_axis == 0:\n        num_slices, height, width = NUM_SLICES, HEIGHT, WIDTH\n    elif trans_axis == 1:\n        num_slices, height, width = HEIGHT, WIDTH, NUM_SLICES\n    elif trans_axis == 2:\n        num_slices, height, width = WIDTH, NUM_SLICES, HEIGHT\n    return num_slices, height, width\n\ndef get_pad_hw(height, width, image_size, edge_size=16):\n    if height <= image_size:\n        pad_h = (image_size - height) // 2 + edge_size\n    else:\n        pad_h = 0\n    if width <= image_size:\n        pad_w = (image_size - width) // 2 + edge_size\n    else:\n        pad_w = 0\n    return pad_h, pad_w\n\ndef perd_add2map(preds, preds_cnt, n, h, w, height, width, pad_h, pad_w, image_size, pred, sample_idx):\n    ### w=0, 4 是pad后在整个slice上的坐标\n    ### \n    w_start = max(w-pad_w, 0)\n    w_end = min(width, w_start+image_size)\n    h_start = max(h-pad_h, 0)\n    h_end = min(height, h_start+image_size)\n    \n    \n    if pad_h != 0:\n        h_pred_start = pad_h - h\n        h_pred_end = h_pred_start + h_end - h_start\n    else:\n        h_pred_start, h_pred_end = 0, image_size\n    if pad_w != 0:\n        w_pred_start = pad_w - w\n        w_pred_end = w_pred_start + w_end - w_start\n    else:\n        w_pred_start, w_pred_end = 0, image_size\n\n    preds[n, h_start:h_end, w_start:w_end] += pred[sample_idx][h_pred_start:h_pred_end, w_pred_start:w_pred_end]\n    preds_cnt[n, h_start:h_end, w_start:w_end] += 1\n    \n    return preds, preds_cnt\n\n\ndef search_thr(preds, labels, image_ids, width, height, min_thr=0.01, max_thr=0.5, interval=0.01):\n    thr_list = np.arange(min_thr, max_thr, interval)\n    thr_list = [round(x, 3) for x in thr_list]\n    best_dice, best_thr = 0, 0\n    for thr in tqdm(thr_list, total=len(thr_list)):\n        ## -------------- 1. Use the thre to get binary pred ---------------\n        bin_preds = preds > thr\n        ## -------------- 2. RLE Encode label and preds for each slice---------------\n        tmp_preds, tmp_labels = [], []\n        for pred, label in zip(bin_preds, labels):\n            tmp_preds.append(rle_encode(pred))\n            tmp_labels.append(rle_encode(label))\n        ## -------------- 3. Get df---------------\n        submit = pd.DataFrame({'id': image_ids, 'rle': tmp_preds, 'width':[width]*len(image_ids), 'height': [height]*len(image_ids)})  \n        label_df = pd.DataFrame({'id': image_ids, 'rle': tmp_labels, 'width':[width]*len(image_ids), 'height': [height]*len(image_ids)})\n        ## -------------- 4. Surface Dice --------------\n        surface_dice = compute_surface_dice_score(submit, label_df)\n        print(f'Surface dice at threshold {thr} is: {surface_dice}')\n        if surface_dice > best_dice:\n            best_dice, best_thr = surface_dice, thr\n\n    print(f'Best Surface dice at threshold {best_thr} is: {best_dice}')\n    \ndef filter_empty_slice(preds, threshold=50):\n    num_slice, h, w = preds.shape\n    for i in range(num_slice):\n        if np.sum(preds[i, :, :]) < threshold:\n            preds[i, :, :] = np.zeros((h, w))\n    return preds","metadata":{"execution":{"iopub.status.busy":"2024-03-09T10:25:46.348343Z","iopub.execute_input":"2024-03-09T10:25:46.348735Z","iopub.status.idle":"2024-03-09T10:25:46.392985Z","shell.execute_reply.started":"2024-03-09T10:25:46.348704Z","shell.execute_reply":"2024-03-09T10:25:46.391999Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"from glob import glob\nSPLIT = 'test'\nBASE_DIR = f'/kaggle/input/blood-vessel-segmentation/{SPLIT}'\nif SPLIT == 'test':\n    kidney_ls = os.listdir(BASE_DIR)\n    image_ls, label_ls = [], []\n    for kidney in kidney_ls:\n        tmp_img_ls = glob(os.path.join(f'{BASE_DIR}/{kidney}/images', '*.tif'))\n        tmp_img_ls.sort()\n        image_ls.append(tmp_img_ls)\n        label_ls.append(None)\nelse:\n\n    kidney = 'kidney_3'\n    if kidney == 'kidney_2':\n        tmp_img_ls = glob(os.path.join(f'{BASE_DIR}/{kidney}/images', '*.tif'))\n        tmp_img_ls.sort()\n        tmp_img_ls = tmp_img_ls[:600]\n        tmp_img_ls = glob(os.path.join(f'{BASE_DIR}/{kidney}/images', '*.tif'))\n        tmp_label_ls = glob(os.path.join(f'{BASE_DIR}/{kidney}/images', '*.tif'))\n    elif kidney == 'kidney_3':\n        path1=f\"{BASE_DIR}/kidney_3_sparse\"\n        path2=f\"{BASE_DIR}/kidney_3_dense\"\n        tmp_label_ls=glob(f\"{path2}/labels/*\")\n        tmp_img_ls=[x.replace(\"labels\",\"images\").replace(\"dense\",\"sparse\") for x in tmp_label_ls]\n        tmp_img_ls.sort()\n        tmp_label_ls.sort()        \n#         tmp_img_ls = tmp_img_ls[100:400]\n#         tmp_label_ls = tmp_label_ls[100:400]\n        \n    image_ls = [tmp_img_ls]    \n    label_ls = [tmp_label_ls]","metadata":{"execution":{"iopub.status.busy":"2024-03-09T10:26:02.610926Z","iopub.execute_input":"2024-03-09T10:26:02.611658Z","iopub.status.idle":"2024-03-09T10:26:02.631697Z","shell.execute_reply.started":"2024-03-09T10:26:02.611626Z","shell.execute_reply":"2024-03-09T10:26:02.630930Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"model_ls = {\n        \"model_name\":\"nnUnet\", \n        \"backbone\":\"nnUnet\",\n#         \"weight\":'/kaggle/input/sennet0121/UNetmaxvitlarge512_900onwards/UNetmaxvitlarge512_900onwards/valonK3_best_loss.pt', \n        \"weight\": '/kaggle/input/new-augmentations-weights/best_epoch_new_aug.bin',\n        \"image_size\":512, \n        \"stride\":512, \n        \"batch_size\":4,\n        \"in_chans\": 1,\n        \"have_large_res\":False,\n        \"tta_ls\":[[2], [2,3]],\n    },\n\nrefine_model_ls = {\n        \"model_name\":\"Unet\", \n        \"backbone\":\"tu-tf_efficientnetv2_s\",\n        \"weight\":'/kaggle/input/model4pseudolabel/NonemptyMask_UNeteffv2s_900onwards_Scale55-105_noval/NonemptyMask_UNeteffv2s_900onwards_Scale55-105_noval/valonK-1_last.pt', \n        \"image_size\":512, \n        \"stride\":512, \n        \"batch_size\":32,\n        \"in_chans\": 1\n    }","metadata":{"execution":{"iopub.status.busy":"2024-03-09T10:26:02.840841Z","iopub.execute_input":"2024-03-09T10:26:02.841202Z","iopub.status.idle":"2024-03-09T10:26:02.849265Z","shell.execute_reply.started":"2024-03-09T10:26:02.841175Z","shell.execute_reply":"2024-03-09T10:26:02.848047Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"REFINE = False\nTTA = True\n# TTA_LS = [[2], [3], [2,3]]\nTRANS_AXIS = [0,1,2]  #0,1,2\nif (len(image_ls[0]) == 3): TRANS_AXIS = [0]\nStage1_Thres = 0.2\nThreshold = 0.4# 0.15\nDiscard_thres = 30\nEMPTY_THRES = 50\n\nsample_submission = {'id':[], 'rle':[]}\n\nmodel = build_model()\n\ndef infer_kidney(NUM_SLICES, HEIGHT, WIDTH,  tmp_image_ls, tmp_label_ls, preds, preds_cnt, nonempty_dict, model_ls=model_ls, trans_axiss=TRANS_AXIS, model=model):\n    for model_idx, model_card in enumerate(model_ls):\n        TTA_LS = model_card['tta_ls']\n        for trans_axis in trans_axiss:\n                    \n            num_slices, height, width = get_nhw(trans_axis, NUM_SLICES, HEIGHT, WIDTH)   \n            if  model_card['have_large_res'] and height > model_card[\"large_size\"] and width > model_card[\"large_size\"]:\n                image_size = model_card[\"large_size\"]\n                stride = model_card[\"large_stride\"]\n                weight_path = model_card['largeRes_weight']\n            else:\n                image_size = model_card[\"image_size\"]\n                stride = model_card[\"stride\"]\n                weight_path = model_card[\"weight\"]\n            pad_h, pad_w = get_pad_hw(height, width, image_size, edge_size=2)\n            \n            val_dataset = Dataset3D(tmp_image_ls, label_list=tmp_label_ls, trans_axis=trans_axis, image_size=image_size, \n                                    stride=stride, aug=False, pad_h=pad_h, pad_w=pad_w, \n                                    in_chans=model_card[\"in_chans\"], nonempty_slices=nonempty_dict[trans_axis])\n            \n            if SPLIT == 'train':\n                if trans_axis==0 and model_idx==len(model_ls)-1:\n                    labels = val_dataset.get_labels()[0]\n            else:\n                labels = None\n                \n            val_dataset = DataLoader(val_dataset, batch_size=model_card[\"batch_size\"] ,num_workers=4, shuffle=False, drop_last=True)\n\n#             model=build_model(model_card[\"model_name\"], model_card[\"backbone\"], in_chans=model_card[\"in_chans\"])\n#             model.load_state_dict(tc.load(weight_path,\"cuda:0\"), strict=True)\n            model = build_model()\n            model.eval()\n\n            preds = permute_axis(preds, trans_axis, 0)\n            preds_cnt = permute_axis(preds_cnt, trans_axis, 0)\n            for batch_idx, (x,idx_3d, n_bs, h_bs, w_bs) in enumerate(tqdm(val_dataset, total=len(val_dataset))):\n\n                x=x.cuda().to(tc.float32)\n                x=norm_with_clip(x.reshape(-1,*x.shape[2:])).reshape(x.shape)\n\n                with autocast():\n                    with tc.no_grad():\n#                         preds = model(x)[0].mean(1)\n#                         preds = preds[0].mean(1)\n                        pred=torch.sigmoid(model(x)[0].mean(1))\n                        if TTA and (len(TTA_LS)>0):\n                            for axis in TTA_LS:  # [2],[3],\n#                                 preds = model(torch.flip(x, dims=axis))[0].mean(1)\n#                                 pdb.set_trace()\n#                                 preds = preds[0].mean(1)\n                                tmp_pred = torch.sigmoid(model(torch.flip(x, dims=axis))[0].mean(1))\n                                axis = [tmp_x-1 for tmp_x in axis]\n                                tmp_pred = torch.flip(tmp_pred, dims=axis)\n                                pred += tmp_pred\n                            pred = pred / (len(TTA_LS)+1)\n                            \n                \n                for sample_idx, (n, h, w) in enumerate(zip(n_bs, h_bs, w_bs)):\n                    if pad_h == 0 and pad_w == 0:\n#                         pdb.set_trace()\n                        preds[n, h:h+image_size, w:w+image_size] += pred[sample_idx]\n                        preds_cnt[n, h:h+image_size, w:w+image_size] += 1\n                    else:\n                        perd_add2map(preds, preds_cnt, n, h, w, height, width, pad_h, pad_w, image_size, pred, sample_idx)\n            preds = permute_axis(preds, trans_axis, 1)\n            preds_cnt = permute_axis(preds_cnt, trans_axis, 1)\n            del model, val_dataset\n            clean_cuda()\n            torch.cuda.empty_cache()\n            \n    return preds, preds_cnt, labels\n\ndef get_empty_slice_dict_with_seg(preds, preds_cnt, pos_thres=0.2, empty_thres=50):\n    tmp_preds = torch.div(preds.cpu(), preds_cnt.cpu()).numpy()\n    tmp_preds = (tmp_preds > pos_thres).astype(np.int8)\n    nonempty_slices = {0:[], 1:[], 2:[]}\n    n, h, w = tmp_preds.shape\n    for i in range(n):\n        if np.sum(tmp_preds[i, :, :]) > empty_thres:\n            nonempty_slices[0].append(i)\n\n    for i in range(h):\n        if np.sum(tmp_preds[:, i, :]) > empty_thres:\n            nonempty_slices[1].append(i)\n\n    for i in range(w):\n        if np.sum(tmp_preds[:, :, i]) > empty_thres:\n            nonempty_slices[2].append(i)\n    del tmp_preds\n    clean_cuda()\n    return nonempty_slices","metadata":{"execution":{"iopub.status.busy":"2024-03-09T10:26:03.085596Z","iopub.execute_input":"2024-03-09T10:26:03.086016Z","iopub.status.idle":"2024-03-09T10:26:03.118271Z","shell.execute_reply.started":"2024-03-09T10:26:03.085980Z","shell.execute_reply":"2024-03-09T10:26:03.117295Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"for tmp_image_ls, tmp_label_ls in zip(image_ls, label_ls):\n    image_ids = [image_path.split('/')[-3] + '_' +image_path.split('/')[-1].split('.')[0] for image_path in tmp_image_ls]\n    sample_submission['id'].extend(image_ids)\n    tmp_img = cv2.imread(tmp_image_ls[0],cv2.IMREAD_GRAYSCALE)\n    HEIGHT, WIDTH = tmp_img.shape\n    NUM_SLICES = len(image_ids)\n    preds, preds_cnt = torch.zeros((NUM_SLICES, HEIGHT, WIDTH), dtype=torch.float16).cuda(), torch.zeros((NUM_SLICES, HEIGHT, WIDTH), dtype=torch.int8).cuda()\n    del tmp_img\n    gc.collect()\n    preds, preds_cnt, labels = infer_kidney(NUM_SLICES, HEIGHT, WIDTH,  tmp_image_ls, tmp_label_ls, preds, preds_cnt, \n                                            nonempty_dict=[None, None, None], model_ls=model_ls, trans_axiss=TRANS_AXIS)\n    \n    if REFINE:\n        nonempty_dict = get_empty_slice_dict_with_seg(preds, preds_cnt, pos_thres=Stage1_Thres, empty_thres=EMPTY_THRES)\n        preds, preds_cnt, _ = infer_kidney(NUM_SLICES, HEIGHT, WIDTH,  tmp_image_ls, tmp_label_ls, preds, preds_cnt, nonempty_dict=nonempty_dict, model_ls=model_ls, trans_axiss=TRANS_AXIS)\n    \n    preds = torch.div(preds.cpu(), preds_cnt.cpu()).numpy()\n    \n    if SPLIT == 'test':\n        ## ------------------ post processing and RLE --------------------\n        preds = (preds > Threshold).astype(np.int8)\n#         preds = filter_empty_slice(preds, threshold=Discard_thres)\n        for pred in preds:\n            sample_submission['rle'].append(rle_encode(pred))\n        del preds, preds_cnt\n        gc.collect()\n    else:\n        search_thr(preds, labels, image_ids, WIDTH, HEIGHT, min_thr=0.1, max_thr=0.5, interval=0.01)\n","metadata":{"execution":{"iopub.status.busy":"2024-03-09T10:26:03.341389Z","iopub.execute_input":"2024-03-09T10:26:03.342506Z","iopub.status.idle":"2024-03-09T11:08:07.144523Z","shell.execute_reply.started":"2024-03-09T10:26:03.342469Z","shell.execute_reply":"2024-03-09T11:08:07.143481Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"if SPLIT == 'test':\n    sample_submission = pd.DataFrame(sample_submission)\n    sample_submission.to_csv('submission.csv', index=False)\n    sample_submission","metadata":{"execution":{"iopub.status.busy":"2024-03-09T11:08:07.146589Z","iopub.execute_input":"2024-03-09T11:08:07.146971Z","iopub.status.idle":"2024-03-09T11:08:07.151688Z","shell.execute_reply.started":"2024-03-09T11:08:07.146930Z","shell.execute_reply":"2024-03-09T11:08:07.150692Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# def load_img(path):\n#     img = cv2.imread(path, cv2.IMREAD_UNCHANGED)\n#     img = np.tile(img[...,None], [1, 1, 3]) # gray to rgb\n#     img = img.astype('float32') # original is uint16\n#     mx = np.max(img)\n#     if mx:\n#         img/=mx # scale image to [0, 1]\n#     return img\n\n# def load_msk(path):\n#     msk = cv2.imread(path, cv2.IMREAD_UNCHANGED)\n#     msk = msk.astype('float32')\n#     msk/=255.0\n#     return msk","metadata":{"execution":{"iopub.status.busy":"2024-03-09T11:08:07.152748Z","iopub.execute_input":"2024-03-09T11:08:07.153005Z","iopub.status.idle":"2024-03-09T11:08:07.163700Z","shell.execute_reply.started":"2024-03-09T11:08:07.152983Z","shell.execute_reply":"2024-03-09T11:08:07.162834Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# class BuildDataset(torch.utils.data.Dataset):\n#     def __init__(self, img_paths, msk_paths=[], transforms=None):\n#         self.img_paths  = img_paths\n#         self.msk_paths  = msk_paths\n#         self.transforms = transforms\n        \n#     def __len__(self):\n#         return len(self.img_paths)\n    \n#     def __getitem__(self, index):\n#         img_path  = self.img_paths[index]\n#         img = load_img(img_path)\n        \n#         if len(self.msk_paths)>0:\n#             msk_path = self.msk_paths[index]\n#             msk = load_msk(msk_path)\n#             if self.transforms:\n#                 data = self.transforms(image=img, mask=msk)\n#                 img  = data['image']\n#                 msk  = data['mask']\n#             img = np.transpose(img, (2, 0, 1))\n#             return torch.tensor(img), torch.tensor(msk)\n#         else:\n#             orig_size = img.shape\n#             if self.transforms:\n#                 data = self.transforms(image=img)\n#                 img  = data['image']\n#             img = np.transpose(img, (2, 0, 1))\n#             return torch.tensor(img), torch.tensor(np.array([orig_size[0], orig_size[1]]))\n    ","metadata":{"execution":{"iopub.status.busy":"2024-03-09T11:08:07.165952Z","iopub.execute_input":"2024-03-09T11:08:07.166242Z","iopub.status.idle":"2024-03-09T11:08:07.177730Z","shell.execute_reply.started":"2024-03-09T11:08:07.166217Z","shell.execute_reply":"2024-03-09T11:08:07.176945Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# import gc\n# def clean_cuda():\n#     torch.cuda.empty_cache()\n#     gc.collect()\n# clean_cuda()","metadata":{"execution":{"iopub.status.busy":"2024-03-09T11:08:07.178761Z","iopub.execute_input":"2024-03-09T11:08:07.179049Z","iopub.status.idle":"2024-03-09T11:08:07.188077Z","shell.execute_reply.started":"2024-03-09T11:08:07.179024Z","shell.execute_reply":"2024-03-09T11:08:07.187290Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# train_groups = CFG.train_groups\n# valid_groups = CFG.valid_groups\n# gt_df = pd.read_csv(CFG.gt_df)\n# gt_df[\"img_path\"] = gt_df[\"img_path\"].apply(lambda x: os.path.join(CFG.data_root, x))\n# gt_df[\"msk_path\"] = gt_df[\"msk_path\"].apply(lambda x: os.path.join(CFG.data_root, x))\n# train_df = gt_df.query(\"group in @train_groups\").reset_index(drop=True)\n# valid_df = gt_df.query(\"group in @valid_groups\").reset_index(drop=True)\n# train_img_paths = train_df[\"img_path\"].values.tolist()\n# train_msk_paths = train_df[\"msk_path\"].values.tolist()\n# valid_img_paths = valid_df[\"img_path\"].values.tolist()\n# valid_msk_paths = valid_df[\"msk_path\"].values.tolist()\n# if CFG.debug:\n#     train_img_paths = train_img_paths[:CFG.train_bs*5]\n#     train_msk_paths = train_msk_paths[:CFG.train_bs*5]\n#     valid_img_paths = valid_img_paths[:CFG.valid_bs*3]\n#     valid_msk_paths = valid_msk_paths[:CFG.valid_bs*3]\n# train_dataset = BuildDataset(train_img_paths, train_msk_paths, transforms=CFG.data_transforms['train'])\n# valid_dataset = BuildDataset(valid_img_paths, valid_msk_paths, transforms=CFG.data_transforms['valid'])\n# train_loader = DataLoader(train_dataset, batch_size=CFG.train_bs, num_workers=0, shuffle=True, pin_memory=True, drop_last=False)\n# valid_loader = DataLoader(valid_dataset, batch_size=CFG.valid_bs, num_workers=0, shuffle=False, pin_memory=True)","metadata":{"execution":{"iopub.status.busy":"2024-03-09T11:08:07.189149Z","iopub.execute_input":"2024-03-09T11:08:07.189445Z","iopub.status.idle":"2024-03-09T11:08:07.198092Z","shell.execute_reply.started":"2024-03-09T11:08:07.189420Z","shell.execute_reply":"2024-03-09T11:08:07.197370Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# sample_ids = [random.randint(0, len(train_img_paths)) for _ in range(5)]\n# for sample_id in sample_ids:\n#     data_name = train_df.loc[sample_id][\"id\"]\n#     img, msk = train_dataset[sample_id]\n#     # img, msk = train_dataset[sample_id]['image'], train_dataset[sample_id]['mask']\n    \n#     img = img.permute((1, 2, 0)).numpy()*255.0\n#     img = img.astype('uint8')\n#     msk = (msk*255).numpy().astype('uint8')\n#     plt.figure(figsize=(9, 4))\n#     print(data_name)\n#     plt.axis('off')\n#     plt.subplot(1,3,1)\n#     plt.imshow(img)\n#     plt.subplot(1,3,2)\n#     plt.imshow(msk)\n#     plt.subplot(1,3,3)\n#     plt.imshow(img, cmap='bone')\n#     plt.imshow(msk, alpha=0.5)\n#     plt.show()","metadata":{"execution":{"iopub.status.busy":"2024-03-09T10:25:53.353714Z","iopub.execute_input":"2024-03-09T10:25:53.353987Z","iopub.status.idle":"2024-03-09T10:25:53.365379Z","shell.execute_reply.started":"2024-03-09T10:25:53.353953Z","shell.execute_reply":"2024-03-09T10:25:53.364543Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# #clean_cuda()\n# torch.cuda.empty_cache()\n# gc.collect()\n","metadata":{"execution":{"iopub.status.busy":"2024-03-09T10:25:53.366561Z","iopub.execute_input":"2024-03-09T10:25:53.367336Z","iopub.status.idle":"2024-03-09T10:25:53.377202Z","shell.execute_reply.started":"2024-03-09T10:25:53.367300Z","shell.execute_reply":"2024-03-09T10:25:53.376266Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# clean_cuda()\n# model.train()\n# # box_np = np.array([[0 , 0, CFG.img_size[0], 224]])\n# # transfer box_np t0 1024x1024 scale\n# # box_224 = np.array([np.array([0, 0,  CFG.img_size[0],  \n# #                               CFG.img_size[0]]) for _ in range(CFG.train_bs)]) # box_np / np.array([224, 224, 224, 224]) * 224\n# # box_224 = np.array([0, 0, 224, 224])\n# # output = model(dt_img.mean(dim=1, keepdim=True))","metadata":{"execution":{"iopub.status.busy":"2024-03-09T09:49:19.809987Z","iopub.status.idle":"2024-03-09T09:49:19.810354Z","shell.execute_reply.started":"2024-03-09T09:49:19.810170Z","shell.execute_reply":"2024-03-09T09:49:19.810187Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# def show_mask(mask, ax, random_color=False):\n#     if random_color:\n#         color = np.concatenate([np.random.random(3), np.array([0.6])], axis=0)\n#         color = np.concatenate([np.array([1, 1, 0]), np.array([0.6])], axis=0)\n\n#     else:\n#         color = np.array([251/255, 252/255, 30/255, 0.6])\n#     h, w = mask.shape[-2:]\n#     mask_image = mask.reshape(h, w, 1) * color.reshape(1, 1, -1)\n#     ax.imshow(mask_image)\n","metadata":{"execution":{"iopub.status.busy":"2024-03-09T09:49:19.812155Z","iopub.status.idle":"2024-03-09T09:49:19.812623Z","shell.execute_reply.started":"2024-03-09T09:49:19.812385Z","shell.execute_reply":"2024-03-09T09:49:19.812407Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"### Only pretrained on hepatic vessel segmentation vs finetuning on the HOA dataset","metadata":{}},{"cell_type":"code","source":"# # # from segmentation_models_pytorch import Unet\n\n# # # Load the MedSAM model\n# # medsam = ViTForImageClassification.from_pretrained(\"wanglab/medsam-vit-base\")\n\n# # # Load the U-Net model\n# # unet = smp.Unet(encoder_name=medsam)# , encoder_weights=\"imagenet\")\n\n# first_iter = iter(valid_loader)\n# dt = next(first_iter)\n# dt_img = dt[0].to(device)\n# dt_mask = dt[1].to(device)\n# # dt = dt['image'].to('cuda')\n# dt_mask.shape","metadata":{"execution":{"iopub.status.busy":"2024-03-02T17:35:19.833920Z","iopub.execute_input":"2024-03-02T17:35:19.834329Z","iopub.status.idle":"2024-03-02T17:35:24.901385Z","shell.execute_reply.started":"2024-03-02T17:35:19.834299Z","shell.execute_reply":"2024-03-02T17:35:24.900276Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# i = 7\n\n# fig, ax = plt.subplots(1, 4, figsize=(10, 5))\n# # box_np = np.array([[0 , 0, 224, 224]])\n# box_np = np.array([[10,10, CFG.img_size[0]-5, CFG.img_size[0]-5]])\n# img = dt_img[i].permute((1, 2, 0)).cpu().numpy()*255.0\n# img = img.astype('uint8')\n# ax[0].imshow(img)\n# ax[0].set_title(\"Input Image\")\n\n# ### Only hepatic pretraining\n# model.load_state_dict(saved_model_nnUnet['state_dict'])\n\n# output_box_for_plot = model(dt_img[i].mean(dim=0, keepdim=True)[None,:])\n# output_box_for_plot = output_box_for_plot[0].mean(dim=1, keepdim=True)#\n# output_box_for_plot = (nn.Sigmoid()(output_box_for_plot)>0.5).double()\n# output_box_for_plot = output_box_for_plot.detach().cpu().numpy()\n\n# ax[1].imshow(img)\n# show_mask(output_box_for_plot, ax[1], True)\n# ax[1].set_title(\"Hepatic nnUnet\")\n\n# ### Fine-tuning\n# model.load_state_dict(best_checkpoint)\n\n# output_box_for_plot = model(dt_img[i].mean(dim=0, keepdim=True)[None,:])\n# output_box_for_plot = output_box_for_plot[0].mean(dim=1, keepdim=True)#\n# output_box_for_plot = (nn.Sigmoid()(output_box_for_plot)>0.5).double()\n# output_box_for_plot = output_box_for_plot.detach().cpu().numpy()\n\n# ax[2].imshow(img)\n# show_mask(output_box_for_plot, ax[2], True)\n\n# ax[2].set_title(\"Hepatic + Renal nnUnet\")\n\n# ax[3].imshow(img)\n# show_mask(dt_mask[i].cpu().numpy(), ax[3], True)\n# # show_box(box_np[0], ax[1])\n# ax[3].set_title('Ground truth mask')\n\n# plt.show()","metadata":{"execution":{"iopub.status.busy":"2024-03-02T17:35:24.902587Z","iopub.execute_input":"2024-03-02T17:35:24.902917Z","iopub.status.idle":"2024-03-02T17:35:26.654551Z","shell.execute_reply.started":"2024-03-02T17:35:24.902891Z","shell.execute_reply":"2024-03-02T17:35:26.653523Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# np.sum(output_box_for_plot)","metadata":{"execution":{"iopub.status.busy":"2024-03-02T17:35:26.655792Z","iopub.execute_input":"2024-03-02T17:35:26.656131Z","iopub.status.idle":"2024-03-02T17:35:26.663489Z","shell.execute_reply.started":"2024-03-02T17:35:26.656096Z","shell.execute_reply":"2024-03-02T17:35:26.662492Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# \"\"\"Surface Dice metric for HuBMAP 4.\"\"\"\n\n# import numpy as np\n# import pandas as pd\n# import pandas.api.types\n# from numba import jit\n# from scipy.ndimage import label, generate_binary_structure\n# from skimage.transform import resize\n# from typing import Optional, Tuple, Union\n\n\n# class ParticipantVisibleError(Exception):\n#     pass\n\n\n# def score(\n#     solution: pd.DataFrame,\n#     submission: pd.DataFrame,\n#     row_id_column_name: str,\n#     rle_column_name: str,\n#     tolerance: float = 1.0,\n#     image_id_column_name: Optional[str] = None,\n#     slice_id_column_name: Optional[str] = None,\n#     resize_fraction: float = 1.0,\n# ) -> float:\n#     \"\"\"Mean Surface Dice over collections of 2D or 3D data.\n\n#     This metric is adapted from the Google DeepMind surface dice metric as found here:\n#     https://github.com/google-deepmind/surface-distance/tree/master.\n\n#     Can be used with either 2D or 3D data. When used with 3D data, each row in the solution and\n#     submission files should contain run-length encoded masks of a width x height slice.\n\n#     Parameters\n#     ----------\n#     solution : Pandas dataframe, the ground truth values.\n\n#     submission : Pandas dataframe, the predicted values.\n\n#     row_id_column_name : str, the name of the ID column used by Kaggle's preprocessing code to\n#         align the solution and submission dataframes.\n\n#     rle_column_name : str, the name of column containing run-length encoded masks.\n\n#     tolerance : float, the distance, in millimeters, the predicted mask surfaces are allowed to\n#         vary from the ground-truth masks.\n\n#     image_id_column_name : str (optional), for 3D data, the name of the column identifying\n#         the image a slice belongs to.\n\n#     slice_id_column_name : str (optional), for 3D data, the name of the column enumerating\n#         the slices within each image.\n\n#     resize_fraction : float, the fraction by which to resize the decoded masks. Useful for memory\n#         efficiency in the case of large images.\n\n#     Returns\n#     -------\n#     mean_surface_dice : float\n\n#     Examples\n#     --------\n#     No groups (2D images).\n#     >>> solution = pd.DataFrame({\n#     ...     'id': [0, 1],\n#     ...     'rle': ['1 12 20 2', '1 6'],\n#     ...     'width': [5, 5],\n#     ...     'height': [5, 5],\n#     ... })\n\n#     Perfect submission.\n#     >>> submission = pd.DataFrame({\n#     ...     'id': [0, 1],\n#     ...     'rle': ['1 12 20 2', '1 6'],\n#     ... })\n#     >>> score(solution, submission, 'id', 'rle', 0.0)\n#     1.0\n\n#     One group with two slices.\n#     >>> solution = pd.DataFrame({\n#     ...     'id': [0, 1],\n#     ...     'rle': ['1 12 20 2', '1 6'],\n#     ...     'width': [5, 5],\n#     ...     'height': [5, 5],\n#     ...     'group': ['a', 'a'],\n#     ...     'slice': [0, 1],\n#     ... })\n\n#     Perfect submission.\n#     >>> submission = pd.DataFrame({\n#     ...     'id': [0, 1],\n#     ...     'rle': ['1 12 20 2', '1 6'],\n#     ... })\n#     >>> score(solution, submission, 'id', 'rle', 0.0, 'group', 'slice')\n#     1.0\n\n#     Null submission.\n#     >>> submission = pd.DataFrame({\n#     ...     'id': [0, 1],\n#     ...     'rle': ['', ''],\n#     ... })\n#     >>> score(solution, submission, 'id', 'rle', 0.0, 'group', 'slice')\n#     0.0\n\n#     Two groups with multiple slices.\n#     >>> solution = pd.DataFrame({\n#     ...     'id': [0, 1, 2, 3, 4],\n#     ...     'rle': ['1 12', '1 12 ', '1 12', '20 5', '22 2'],\n#     ...     'width': [5, 5, 5, 8, 8],\n#     ...     'height': [5, 5, 5, 4, 4],\n#     ...     'group': ['a', 'a', 'a', 'b', 'b'],\n#     ...     'slice': [0, 1, 2, 0, 1],\n#     ... })\n#     >>> submission = pd.DataFrame({\n#     ...     'id': [0, 1, 2, 3, 4],\n#     ...     'rle': ['1 12', '1 12 ', '1 12', '20 5', '22 2'],\n#     ... })\n#     >>> score(solution, submission, 'id', 'rle', 0.0, 'group', 'slice')\n#     1.0\n\n#     >>> submission = pd.DataFrame({\n#     ...     'id': [0, 1, 2, 3, 4],\n#     ...     'rle': ['1 11', '1 12 ', '1 13', '20 4', '22 3'],\n#     ... })\n\n#     With non-zero tolerance.\n#     >>> score(solution, submission, 'id', 'rle', 5.0, 'group', 'slice')\n#     1.0\n\n#     Cubes, adapted from https://github.com/google-deepmind/surface-distance/blob/master/surface_distance_test.py#L307\n#     >>> solution = pd.DataFrame({\n#     ...     'id': np.arange(100),\n#     ...     'rle': ['1 10000' if k <= 49 else '' for k in np.arange(100)],\n#     ...     'width': [100] * 100,\n#     ...     'height': [100] * 100,\n#     ...     'group': ['a'] * 100,\n#     ...     'slice': np.arange(100),\n#     ... })\n#     >>> submission = pd.DataFrame({\n#     ...     'id': np.arange(100),\n#     ...     'rle': ['1 10000' if k <= 50 else '' for k in np.arange(100)],\n#     ... })\n#     >>> score(solution, submission, 'id', 'rle', 0.0, 'group', 'slice')\n#     0.750877...\n\n#     \"\"\"\n#     solution = solution.set_index(row_id_column_name)\n#     submission = submission.set_index(row_id_column_name)\n\n#     # Check that both are defined or neither are\n#     if (image_id_column_name is None) != (slice_id_column_name is None):\n#         raise ValueError(\"If one of `image_id_column_name` or `slice_id_column_name` is given, the other must be given also.\")\n\n#     if tolerance < 0.0:\n#         raise ValueError(\"`tolerance` must be non-negative.\")\n\n#     # Setup spacing_mm\n#     if image_id_column_name is None:\n#         spacing_mm = (1, 1)  # (height, width)\n#     else:\n#         spacing_mm = (1, 1, 1)  # (height, width, depth)\n\n#     # Create joined dataframe to iterate over\n#     joined = solution.join(submission,\n#                            lsuffix='_sol',\n#                            rsuffix='_sub')\n\n#     # Compute surface dice for each group of slices and total\n#     total_dice = 0.0\n#     group_col = row_id_column_name if image_id_column_name is None else image_id_column_name\n#     for group, df in joined.groupby(group_col):\n#         # Make indexing easier\n#         df = df.reset_index(drop=True)\n\n#         # Check that we're stacking slices in order\n#         if image_id_column_name is not None:\n#             assert df.loc[:, slice_id_column_name].is_monotonic_increasing\n\n#         group_mask_sol, group_mask_sub = [], []\n#         height, width = df.loc[0, 'height'], df.loc[0, 'width']\n#         assert (df.loc[:, 'height'].nunique() == 1) and (df.loc[:, 'width'].nunique() == 1),\\\n#             \"Height and width must be constant within each group.\"\n\n#         # Decode slice RLEs to arrays\n#         for row in df.itertuples():\n#             # Mask transformation\n#             shape = (height, width)\n#             group_mask_sol.append(make_mask(row.rle_sol, shape, resize_fraction))\n#             group_mask_sub.append(make_mask(row.rle_sub, shape, resize_fraction))\n\n#         # Stack slices to create 3D arrays\n#         group_mask_sol = np.stack(group_mask_sol, axis=-1)\n#         assert group_mask_sol.ndim == 3\n#         group_mask_sub = np.stack(group_mask_sub, axis=-1)\n#         assert group_mask_sub.ndim == 3\n\n#         # Minimize dims: 2D or 3D\n#         group_mask_sol = group_mask_sol.squeeze()\n#         group_mask_sub = group_mask_sub.squeeze()\n\n#         # Compute surface distance\n#         surface_dist = compute_surface_distances(\n#             mask_gt=group_mask_sol,\n#             mask_pred=group_mask_sub,\n#             spacing_mm=spacing_mm,\n#         )\n#         # Compute surface dice and add to running total\n#         total_dice += compute_surface_dice_at_tolerance(\n#             surface_dist,\n#             tolerance_mm=tolerance,\n#         )\n\n#     # Compute mean from running total\n#     if image_id_column_name is None:\n#         ngroups = len(joined)\n#     else:\n#         ngroups = joined.loc[:, image_id_column_name].nunique()\n#     mean_surface_dice = total_dice / ngroups\n#     return mean_surface_dice\n\n\n# def make_mask(rle, shape, resize_fraction):\n#     if resize_fraction < 0.0 or resize_fraction > 1.0:\n#         raise ValueError(\"`resize_fraction` must be between 0.0 and 1.0, inclusive.\")\n\n#     mask = rle_decode(rle, shape)\n#     if resize_fraction < 1.0:\n#         new_shape = int(shape[0] * resize_fraction), int(shape[1] * resize_fraction)\n#         mask = voting_resize(mask, new_shape)\n#     mask = mask.astype(bool)\n#     return mask\n\n\n# def voting_resize(mask, new_shape):\n#     interpolated_mask = resize(mask, new_shape, order=1, preserve_range=True)  # Using bilinear interpolation\n#     voting_mask = (interpolated_mask > 0.5).astype(np.uint8)\n#     return voting_mask\n\n\n# # Below is code from https://www.kaggle.com/code/paulorzp/run-length-encode-and-decode\n# def rle_encode(img):\n#     '''\n#     img: numpy array, 1 - mask, 0 - background\n#     Returns run length as string formated\n#     '''\n#     pixels = img.flatten()\n#     pixels = np.concatenate([[0], pixels, [0]])\n#     runs = np.where(pixels[1:] != pixels[:-1])[0] + 1\n#     d\n#     return ' '.join(str(x) for x in runs)\n\n\n# def rle_decode(mask_rle, shape):\n#     '''\n#     mask_rle: run-length as string formated (start length)\n#     shape: (height,width) of array to return\n#     Returns numpy array, 1 - mask, 0 - background\n\n#     '''\n#     s = mask_rle.split()\n#     starts, lengths = [\n#         np.asarray(x, dtype=int) for x in (s[0:][::2], s[1:][::2])\n#     ]\n#     starts -= 1\n#     ends = starts + lengths\n#     img = np.zeros(shape[0] * shape[1], dtype=np.uint8)\n#     for lo, hi in zip(starts, ends):\n#         img[lo:hi] = 1\n#     return img.reshape(shape)\n","metadata":{"_kg_hide-input":true,"_kg_hide-output":true,"execution":{"iopub.status.busy":"2024-03-02T17:35:26.665063Z","iopub.execute_input":"2024-03-02T17:35:26.665448Z","iopub.status.idle":"2024-03-02T17:35:27.234641Z","shell.execute_reply.started":"2024-03-02T17:35:26.665412Z","shell.execute_reply":"2024-03-02T17:35:27.233835Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# def rle_encode(img):\n#     '''\n#     img: numpy array, 1 - mask, 0 - background\n#     Returns run length as string formated\n#     '''\n#     pixels = img.flatten()\n#     pixels = np.concatenate([[0], pixels, [0]])\n#     runs = np.where(pixels[1:] != pixels[:-1])[0] + 1\n#     runs[1::2] -= runs[::2]\n#     rle = ' '.join(str(x) for x in runs)\n#     if rle=='':\n#         rle = '1 0'\n#     return rle","metadata":{"execution":{"iopub.status.busy":"2024-03-02T17:35:27.235860Z","iopub.execute_input":"2024-03-02T17:35:27.236188Z","iopub.status.idle":"2024-03-02T17:35:27.243303Z","shell.execute_reply.started":"2024-03-02T17:35:27.236161Z","shell.execute_reply":"2024-03-02T17:35:27.242142Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# DATASET_FOLDER = \"/kaggle/input/blood-vessel-segmentation\"\n# ls_images = glob(os.path.join(DATASET_FOLDER, \"test\", \"*\", \"*\", \"*.tif\"))\n# test_dataset = BuildDataset(ls_images, [], transforms=CFG.data_transforms['valid'])\n# test_loader = DataLoader(test_dataset, batch_size=CFG.valid_bs, num_workers=0, shuffle=False, pin_memory=True, drop_last=False)\n# print(f\"found images: {len(ls_images)}\")","metadata":{"execution":{"iopub.status.busy":"2024-02-06T16:06:38.313363Z","iopub.execute_input":"2024-02-06T16:06:38.314353Z","iopub.status.idle":"2024-02-06T16:06:38.323386Z","shell.execute_reply.started":"2024-02-06T16:06:38.314316Z","shell.execute_reply":"2024-02-06T16:06:38.322370Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# pred[0].shape\n# output_box_for_plot[0][0].shape","metadata":{"execution":{"iopub.status.busy":"2024-02-06T14:44:40.809444Z","iopub.execute_input":"2024-02-06T14:44:40.809900Z","iopub.status.idle":"2024-02-06T14:44:40.814219Z","shell.execute_reply.started":"2024-02-06T14:44:40.809842Z","shell.execute_reply":"2024-02-06T14:44:40.813164Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# # preds.shape\n# cv2.resize(pred[0], (shape[1], shape[0]), cv2.INTER_NEAREST).shape\n# dummy = cv2.resize(output_box_for_plot[0][0], (shape[1],shape[0]), cv2.INTER_NEAREST)\n# np.sum(dummy)\n# dummy.shape\n","metadata":{"execution":{"iopub.status.busy":"2024-02-06T14:44:46.957484Z","iopub.execute_input":"2024-02-06T14:44:46.957852Z","iopub.status.idle":"2024-02-06T14:44:46.962304Z","shell.execute_reply.started":"2024-02-06T14:44:46.957823Z","shell.execute_reply":"2024-02-06T14:44:46.961385Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# ### Inference for nnUnet\n# rles = []\n# pbar = tqdm(enumerate(test_loader), total=len(test_loader), desc='Inference ')\n# # pbar = tqdm(enumerate(valid_loader), total=len(valid_loader), desc='Inference ')\n\n# for step, (images, shapes) in pbar:\n#     # pdb.set_trace()\n#     shapes = shapes.numpy()\n#     images = images.to(CFG.device, dtype=torch.float)\n#     batch_size = images.size(0)\n#     # print(images.shape)\n#     with torch.no_grad():\n#         images = images.mean(dim=1, keepdim=True)\n#         preds = model(images)\n#         preds = preds[0].mean(1)\n#         preds = (nn.Sigmoid()(preds)>0.5).double()\n#     preds = preds.cpu().numpy().astype(np.uint8)\n\n#     for pred, shape in zip(preds, shapes):\n#         # pred = cv2.resize(pred[0], (shape[1], shape[0]), cv2.INTER_NEAREST)\n#         # Target size\n# #         original_size = pred.shape\n# #         target_size = shape\n\n# #         # Calculate the scale factor while maintaining the aspect ratio\n# #         scale_factor = min(target_size[0] / original_size[0], target_size[1] / original_size[1])\n# #         new_size = (int(original_size[0] * scale_factor), int(original_size[1] * scale_factor))\n\n# #         # Resize the mask\n# #         resized_mask = cv2.resize(pred, new_size, interpolation=cv2.INTER_NEAREST)\n\n# #         # Calculate padding\n# #         pad_width = (target_size[0] - new_size[0]) // 2\n# #         pad_height = (target_size[1] - new_size[1]) // 2\n\n# #         # Pad the resized mask to match the target size\n# #         padded_mask = cv2.copyMakeBorder(resized_mask, pad_height, pad_height, pad_width, pad_width, cv2.BORDER_CONSTANT, value=0)\n\n#         #pred = cv2.resize(pred[0], (256, 256))\n# #         pred = cv2.resize(pred[0], (1303,  912), cv2.INTER_NEAREST)\n#         rle = rle_encode(pred)\n#         rles.append(rle)\n# #     for pred in preds:\n# #         # Target size\n# #         original_size = pred.shape\n# #         target_size = (1303, 912)\n\n# #         # Calculate the scale factor while maintaining the aspect ratio\n# #         scale_factor = min(target_size[0] / original_size[0], target_size[1] / original_size[1])\n# #         new_size = (int(original_size[0] * scale_factor), int(original_size[1] * scale_factor))\n\n# #         # Resize the mask\n# #         resized_mask = cv2.resize(pred, new_size, interpolation=cv2.INTER_NEAREST)\n\n# #         # Calculate padding\n# #         pad_width = (target_size[0] - new_size[0]) // 2\n# #         pad_height = (target_size[1] - new_size[1]) // 2\n\n# #         # Pad the resized mask to match the target size\n# #         padded_mask = cv2.copyMakeBorder(resized_mask, pad_height, pad_height, pad_width, pad_width, cv2.BORDER_CONSTANT, value=0)\n\n# #         #pred = cv2.resize(pred[0], (256, 256))\n# # #         pred = cv2.resize(pred[0], (1303,  912), cv2.INTER_NEAREST)\n# #         rle = rle_encode(resized_mask)\n# #         rles.append(rle)","metadata":{"execution":{"iopub.status.busy":"2024-02-06T16:05:45.663914Z","iopub.execute_input":"2024-02-06T16:05:45.664348Z","iopub.status.idle":"2024-02-06T16:05:45.943697Z","shell.execute_reply.started":"2024-02-06T16:05:45.664317Z","shell.execute_reply":"2024-02-06T16:05:45.942688Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# fig, ax = plt.subplots(1, 1, figsize=(10, 5))\n\n# val_pred = rle_decode(rles[0], (256, 256))\n# # shapes # 1303, 912\n# # shapes.shape\n# # plt.imshow(val_pred)\n# show_mask(val_pred, ax)","metadata":{"execution":{"iopub.status.busy":"2024-02-06T16:05:56.647290Z","iopub.execute_input":"2024-02-06T16:05:56.648233Z","iopub.status.idle":"2024-02-06T16:05:56.652178Z","shell.execute_reply.started":"2024-02-06T16:05:56.648197Z","shell.execute_reply":"2024-02-06T16:05:56.651251Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# fig, ax = plt.subplots(1, 1, figsize=(10, 5))\n\n# val_pred = rle_decode(rles[2], (256, 256))\n# val_pred = rle_decode(rles[0], (1303, 912))\n# # shapes # 1303, 912\n# # shapes.shape\n# # plt.imshow(val_pred)\n# show_mask(val_pred, ax)","metadata":{"execution":{"iopub.status.busy":"2024-02-06T16:05:59.906990Z","iopub.execute_input":"2024-02-06T16:05:59.907734Z","iopub.status.idle":"2024-02-06T16:05:59.911778Z","shell.execute_reply.started":"2024-02-06T16:05:59.907704Z","shell.execute_reply":"2024-02-06T16:05:59.910762Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# ids = [f'{p.split(\"/\")[-3]}_{os.path.basename(p).split(\".\")[0]}' for p in ls_images]\n# submission = pd.DataFrame.from_dict({\n#     \"id\": ids,\n#     \"rle\": rles\n# })\n# submission.head()","metadata":{"execution":{"iopub.status.busy":"2024-02-06T16:06:01.715497Z","iopub.execute_input":"2024-02-06T16:06:01.716412Z","iopub.status.idle":"2024-02-06T16:06:01.731455Z","shell.execute_reply.started":"2024-02-06T16:06:01.716374Z","shell.execute_reply":"2024-02-06T16:06:01.730489Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# test_img = next(iter(test_loader))[0]\n# fst_tst_img = test_img[0]\n# fst_tst_img.shape\n# fst_tst_img.shape","metadata":{"execution":{"iopub.status.busy":"2024-02-06T16:06:05.976262Z","iopub.execute_input":"2024-02-06T16:06:05.977391Z","iopub.status.idle":"2024-02-06T16:06:06.100139Z","shell.execute_reply.started":"2024-02-06T16:06:05.977355Z","shell.execute_reply":"2024-02-06T16:06:06.099104Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# encoded = rle_encode(output_box_for_plot)\n# decoded = rle_decode(encoded, (256, 256))\n# np.sum(decoded != output_box_for_plot)","metadata":{"execution":{"iopub.status.busy":"2024-02-06T16:06:13.835880Z","iopub.execute_input":"2024-02-06T16:06:13.836308Z","iopub.status.idle":"2024-02-06T16:06:13.840905Z","shell.execute_reply.started":"2024-02-06T16:06:13.836277Z","shell.execute_reply":"2024-02-06T16:06:13.839901Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# i = 0\n\n# fig, ax = plt.subplots(1, 3, figsize=(10, 5))\n# test_img = next(iter(valid_loader))[0]\n# test_img = test_img[i].to(device)\n# img = test_img.permute((1, 2, 0)).cpu().numpy()*255.0 # dt_img[i]\n# img = img.astype('uint8')\n# ax[0].imshow(img)\n# ax[0].set_title(\"Test Input Image\")\n\n# ### Only hepatic pretraining\n# model.load_state_dict(saved_model_nnUnet['state_dict'])\n\n# output_box_for_plot = model(test_img.mean(dim=0, keepdim=True)[None,:])\n# output_box_for_plot = output_box_for_plot[0].mean(dim=1, keepdim=True).detach().cpu().numpy()\n\n# ax[1].imshow(img)\n# show_mask(output_box_for_plot, ax[1], True)\n# ax[1].set_title(\"Hepatic nnUnet\")\n\n# ### Fine-tuning\n# model.load_state_dict(best_checkpoint)\n\n# output_box_for_plot = model(test_img.mean(dim=0, keepdim=True)[None,:])\n# output_box_for_plot = output_box_for_plot[0].mean(dim=1, keepdim=True).detach().cpu().numpy()\n\n# ax[2].imshow(img)\n# show_mask(output_box_for_plot, ax[2], True)\n\n# ax[2].set_title(\"Hepatic + Renal nnUnet\")\n\n\n# plt.show()","metadata":{"execution":{"iopub.status.busy":"2024-02-06T16:06:17.877754Z","iopub.execute_input":"2024-02-06T16:06:17.878806Z","iopub.status.idle":"2024-02-06T16:06:20.505320Z","shell.execute_reply.started":"2024-02-06T16:06:17.878767Z","shell.execute_reply":"2024-02-06T16:06:20.504249Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# val_score = score(submission, _gt_df, \"id\", \"rle\", 0.0, \"group\", \"slice\")\n# print(val_score)","metadata":{"execution":{"iopub.status.busy":"2024-01-28T17:27:18.468220Z","iopub.execute_input":"2024-01-28T17:27:18.469117Z","iopub.status.idle":"2024-01-28T17:27:18.473117Z","shell.execute_reply.started":"2024-01-28T17:27:18.469079Z","shell.execute_reply":"2024-01-28T17:27:18.472114Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# submission.to_csv(\"submission.csv\", index=False)","metadata":{"execution":{"iopub.status.busy":"2024-01-28T17:27:33.163478Z","iopub.execute_input":"2024-01-28T17:27:33.163855Z","iopub.status.idle":"2024-01-28T17:27:33.172373Z","shell.execute_reply.started":"2024-01-28T17:27:33.163827Z","shell.execute_reply":"2024-01-28T17:27:33.171417Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# IGNORE","metadata":{}},{"cell_type":"code","source":"\n\n# # Below is adapted from https://github.com/google-deepmind/surface-distance\n\n# # Copyright 2018 Google Inc. All Rights Reserved.\n# #\n# # Licensed under the Apache License, Version 2.0 (the \"License\");\n# # you may not use this file except in compliance with the License.\n# # You may obtain a copy of the License at\n# #\n# #      http://www.apache.org/licenses/LICENSE-2.0\n# #\n# # Unless required by applicable law or agreed to in writing, software\n# # distributed under the License is distributed on an \"AS-IS\" BASIS,\n# # WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.\n# # See the License for the specific language governing permissions and\n# # limitations under the License.\n\n# # surface_distance/lookup_tables.py\n# \"\"\"Lookup tables used by surface distance metrics.\"\"\"\n\n# import math\n# import numpy as np\n\n# ENCODE_NEIGHBOURHOOD_3D_KERNEL = np.array([[[128, 64], [32, 16]],\n#                                            [[8, 4], [2, 1]]])\n\n# # _NEIGHBOUR_CODE_TO_NORMALS is a lookup table.\n# # For every binary neighbour code\n# # (2x2x2 neighbourhood = 8 neighbours = 8 bits = 256 codes)\n# # it contains the surface normals of the triangles (called \"surfel\" for\n# # \"surface element\" in the following). The length of the normal\n# # vector encodes the surfel area.\n# #\n# # created using the marching_cube algorithm\n# # see e.g. https://en.wikipedia.org/wiki/Marching_cubes\n# # pylint: disable=line-too-long\n# _NEIGHBOUR_CODE_TO_NORMALS = [[[0, 0, 0]], [[0.125, 0.125, 0.125]],\n#                               [[-0.125, -0.125, 0.125]],\n#                               [[-0.25, -0.25, 0.0], [0.25, 0.25, -0.0]],\n#                               [[0.125, -0.125, 0.125]],\n#                               [[-0.25, -0.0, -0.25], [0.25, 0.0, 0.25]],\n#                               [[0.125, -0.125, 0.125], [-0.125, -0.125,\n#                                                         0.125]],\n#                               [[0.5, 0.0, -0.0], [0.25, 0.25, 0.25],\n#                                [0.125, 0.125, 0.125]], [[-0.125, 0.125,\n#                                                          0.125]],\n#                               [[0.125, 0.125, 0.125], [-0.125, 0.125, 0.125]],\n#                               [[-0.25, 0.0, 0.25], [-0.25, 0.0, 0.25]],\n#                               [[0.5, 0.0, 0.0], [-0.25, -0.25, 0.25],\n#                                [-0.125, -0.125, 0.125]],\n#                               [[0.25, -0.25, 0.0], [0.25, -0.25, 0.0]],\n#                               [[0.5, 0.0, 0.0], [0.25, -0.25, 0.25],\n#                                [-0.125, 0.125, -0.125]],\n#                               [[-0.5, 0.0, 0.0], [-0.25, 0.25, 0.25],\n#                                [-0.125, 0.125, 0.125]],\n#                               [[0.5, 0.0, 0.0], [0.5, 0.0, 0.0]],\n#                               [[0.125, -0.125, -0.125]],\n#                               [[0.0, -0.25, -0.25], [0.0, 0.25, 0.25]],\n#                               [[-0.125, -0.125, 0.125],\n#                                [0.125, -0.125, -0.125]],\n#                               [[0.0, -0.5, 0.0], [0.25, 0.25, 0.25],\n#                                [0.125, 0.125, 0.125]],\n#                               [[0.125, -0.125, 0.125], [0.125, -0.125,\n#                                                         -0.125]],\n#                               [[0.0, 0.0, -0.5], [0.25, 0.25, 0.25],\n#                                [-0.125, -0.125, -0.125]],\n#                               [[-0.125, -0.125, 0.125], [0.125, -0.125, 0.125],\n#                                [0.125, -0.125, -0.125]],\n#                               [[-0.125, -0.125, -0.125], [-0.25, -0.25, -0.25],\n#                                [0.25, 0.25, 0.25], [0.125, 0.125, 0.125]],\n#                               [[-0.125, 0.125, 0.125], [0.125, -0.125,\n#                                                         -0.125]],\n#                               [[0.0, -0.25, -0.25], [0.0, 0.25, 0.25],\n#                                [-0.125, 0.125, 0.125]],\n#                               [[-0.25, 0.0, 0.25], [-0.25, 0.0, 0.25],\n#                                [0.125, -0.125, -0.125]],\n#                               [[0.125, 0.125, 0.125], [0.375, 0.375, 0.375],\n#                                [0.0, -0.25, 0.25], [-0.25, 0.0, 0.25]],\n#                               [[0.125, -0.125, -0.125], [0.25, -0.25, 0.0],\n#                                [0.25, -0.25, 0.0]],\n#                               [[0.375, 0.375, 0.375], [0.0, 0.25, -0.25],\n#                                [-0.125, -0.125, -0.125], [-0.25, 0.25, 0.0]],\n#                               [[-0.5, 0.0, 0.0], [-0.125, -0.125, -0.125],\n#                                [-0.25, -0.25, -0.25], [0.125, 0.125, 0.125]],\n#                               [[-0.5, 0.0, 0.0], [-0.125, -0.125, -0.125],\n#                                [-0.25, -0.25, -0.25]], [[0.125, -0.125,\n#                                                          0.125]],\n#                               [[0.125, 0.125, 0.125], [0.125, -0.125, 0.125]],\n#                               [[0.0, -0.25, 0.25], [0.0, 0.25, -0.25]],\n#                               [[0.0, -0.5, 0.0], [0.125, 0.125, -0.125],\n#                                [0.25, 0.25, -0.25]],\n#                               [[0.125, -0.125, 0.125], [0.125, -0.125, 0.125]],\n#                               [[0.125, -0.125, 0.125], [-0.25, -0.0, -0.25],\n#                                [0.25, 0.0, 0.25]],\n#                               [[0.0, -0.25, 0.25], [0.0, 0.25, -0.25],\n#                                [0.125, -0.125, 0.125]],\n#                               [[-0.375, -0.375, 0.375], [-0.0, 0.25, 0.25],\n#                                [0.125, 0.125, -0.125], [-0.25, -0.0, -0.25]],\n#                               [[-0.125, 0.125, 0.125], [0.125, -0.125, 0.125]],\n#                               [[0.125, 0.125, 0.125], [0.125, -0.125, 0.125],\n#                                [-0.125, 0.125, 0.125]],\n#                               [[-0.0, 0.0, 0.5], [-0.25, -0.25, 0.25],\n#                                [-0.125, -0.125, 0.125]],\n#                               [[0.25, 0.25, -0.25], [0.25, 0.25, -0.25],\n#                                [0.125, 0.125, -0.125], [-0.125, -0.125,\n#                                                         0.125]],\n#                               [[0.125, -0.125, 0.125], [0.25, -0.25, 0.0],\n#                                [0.25, -0.25, 0.0]],\n#                               [[0.5, 0.0, 0.0], [0.25, -0.25, 0.25],\n#                                [-0.125, 0.125, -0.125], [0.125, -0.125,\n#                                                          0.125]],\n#                               [[0.0, 0.25, -0.25], [0.375, -0.375, -0.375],\n#                                [-0.125, 0.125, 0.125], [0.25, 0.25, 0.0]],\n#                               [[-0.5, 0.0, 0.0], [-0.25, -0.25, 0.25],\n#                                [-0.125, -0.125, 0.125]],\n#                               [[0.25, -0.25, 0.0], [-0.25, 0.25, 0.0]],\n#                               [[0.0, 0.5, 0.0], [-0.25, 0.25, 0.25],\n#                                [0.125, -0.125, -0.125]],\n#                               [[0.0, 0.5, 0.0], [0.125, -0.125, 0.125],\n#                                [-0.25, 0.25, -0.25]],\n#                               [[0.0, 0.5, 0.0], [0.0, -0.5, 0.0]],\n#                               [[0.25, -0.25, 0.0], [-0.25, 0.25, 0.0],\n#                                [0.125, -0.125, 0.125]],\n#                               [[-0.375, -0.375, -0.375], [-0.25, 0.0, 0.25],\n#                                [-0.125, -0.125, -0.125], [-0.25, 0.25, 0.0]],\n#                               [[0.125, 0.125, 0.125], [0.0, -0.5, 0.0],\n#                                [-0.25, -0.25, -0.25], [-0.125, -0.125,\n#                                                        -0.125]],\n#                               [[0.0, -0.5, 0.0], [-0.25, -0.25, -0.25],\n#                                [-0.125, -0.125, -0.125]],\n#                               [[-0.125, 0.125, 0.125], [0.25, -0.25, 0.0],\n#                                [-0.25, 0.25, 0.0]],\n#                               [[0.0, 0.5, 0.0], [0.25, 0.25, -0.25],\n#                                [-0.125, -0.125, 0.125],\n#                                [-0.125, -0.125, 0.125]],\n#                               [[-0.375, 0.375, -0.375], [-0.25, -0.25, 0.0],\n#                                [-0.125, 0.125, -0.125], [-0.25, 0.0, 0.25]],\n#                               [[0.0, 0.5, 0.0], [0.25, 0.25, -0.25],\n#                                [-0.125, -0.125, 0.125]],\n#                               [[0.25, -0.25, 0.0], [-0.25, 0.25, 0.0],\n#                                [0.25, -0.25, 0.0], [0.25, -0.25, 0.0]],\n#                               [[-0.25, -0.25, 0.0], [-0.25, -0.25, 0.0],\n#                                [-0.125, -0.125, 0.125]],\n#                               [[0.125, 0.125, 0.125], [-0.25, -0.25, 0.0],\n#                                [-0.25, -0.25, 0.0]],\n#                               [[-0.25, -0.25, 0.0], [-0.25, -0.25, 0.0]],\n#                               [[-0.125, -0.125, 0.125]],\n#                               [[0.125, 0.125, 0.125], [-0.125, -0.125, 0.125]],\n#                               [[-0.125, -0.125, 0.125],\n#                                [-0.125, -0.125, 0.125]],\n#                               [[-0.125, -0.125, 0.125], [-0.25, -0.25, 0.0],\n#                                [0.25, 0.25, -0.0]],\n#                               [[0.0, -0.25, 0.25], [0.0, -0.25, 0.25]],\n#                               [[0.0, 0.0, 0.5], [0.25, -0.25, 0.25],\n#                                [0.125, -0.125, 0.125]],\n#                               [[0.0, -0.25, 0.25], [0.0, -0.25, 0.25],\n#                                [-0.125, -0.125, 0.125]],\n#                               [[0.375, -0.375, 0.375], [0.0, -0.25, -0.25],\n#                                [-0.125, 0.125, -0.125], [0.25, 0.25, 0.0]],\n#                               [[-0.125, -0.125, 0.125], [-0.125, 0.125,\n#                                                          0.125]],\n#                               [[0.125, 0.125, 0.125], [-0.125, -0.125, 0.125],\n#                                [-0.125, 0.125, 0.125]],\n#                               [[-0.125, -0.125, 0.125], [-0.25, 0.0, 0.25],\n#                                [-0.25, 0.0, 0.25]],\n#                               [[0.5, 0.0, 0.0], [-0.25, -0.25, 0.25],\n#                                [-0.125, -0.125, 0.125],\n#                                [-0.125, -0.125, 0.125]],\n#                               [[-0.0, 0.5, 0.0], [-0.25, 0.25, -0.25],\n#                                [0.125, -0.125, 0.125]],\n#                               [[-0.25, 0.25, -0.25], [-0.25, 0.25, -0.25],\n#                                [-0.125, 0.125, -0.125],\n#                                [-0.125, 0.125, -0.125]],\n#                               [[-0.25, 0.0, -0.25], [0.375, -0.375, -0.375],\n#                                [0.0, 0.25, -0.25], [-0.125, 0.125, 0.125]],\n#                               [[0.5, 0.0, 0.0], [-0.25, 0.25, -0.25],\n#                                [0.125, -0.125, 0.125]],\n#                               [[-0.25, 0.0, 0.25], [0.25, 0.0, -0.25]],\n#                               [[-0.0, 0.0, 0.5], [-0.25, 0.25, 0.25],\n#                                [-0.125, 0.125, 0.125]],\n#                               [[-0.125, -0.125, 0.125], [-0.25, 0.0, 0.25],\n#                                [0.25, 0.0, -0.25]],\n#                               [[-0.25, -0.0, -0.25], [-0.375, 0.375, 0.375],\n#                                [-0.25, -0.25, 0.0], [-0.125, 0.125, 0.125]],\n#                               [[0.0, 0.0, -0.5], [0.25, 0.25, -0.25],\n#                                [-0.125, -0.125, 0.125]],\n#                               [[-0.0, 0.0, 0.5], [0.0, 0.0, 0.5]],\n#                               [[0.125, 0.125, 0.125], [0.125, 0.125, 0.125],\n#                                [0.25, 0.25, 0.25], [0.0, 0.0, 0.5]],\n#                               [[0.125, 0.125, 0.125], [0.25, 0.25, 0.25],\n#                                [0.0, 0.0, 0.5]],\n#                               [[-0.25, 0.0, 0.25], [0.25, 0.0, -0.25],\n#                                [-0.125, 0.125, 0.125]],\n#                               [[-0.0, 0.0, 0.5], [0.25, -0.25, 0.25],\n#                                [0.125, -0.125, 0.125], [0.125, -0.125, 0.125]],\n#                               [[-0.25, 0.0, 0.25], [-0.25, 0.0, 0.25],\n#                                [-0.25, 0.0, 0.25], [0.25, 0.0, -0.25]],\n#                               [[0.125, -0.125, 0.125], [0.25, 0.0, 0.25],\n#                                [0.25, 0.0, 0.25]],\n#                               [[0.25, 0.0, 0.25], [-0.375, -0.375, 0.375],\n#                                [-0.25, 0.25, 0.0], [-0.125, -0.125, 0.125]],\n#                               [[-0.0, 0.0, 0.5], [0.25, -0.25, 0.25],\n#                                [0.125, -0.125, 0.125]],\n#                               [[0.125, 0.125, 0.125], [0.25, 0.0, 0.25],\n#                                [0.25, 0.0, 0.25]],\n#                               [[0.25, 0.0, 0.25], [0.25, 0.0, 0.25]],\n#                               [[-0.125, -0.125, 0.125], [0.125, -0.125,\n#                                                          0.125]],\n#                               [[0.125, 0.125, 0.125], [-0.125, -0.125, 0.125],\n#                                [0.125, -0.125, 0.125]],\n#                               [[-0.125, -0.125, 0.125], [0.0, -0.25, 0.25],\n#                                [0.0, 0.25, -0.25]],\n#                               [[0.0, -0.5, 0.0], [0.125, 0.125, -0.125],\n#                                [0.25, 0.25, -0.25], [-0.125, -0.125, 0.125]],\n#                               [[0.0, -0.25, 0.25], [0.0, -0.25, 0.25],\n#                                [0.125, -0.125, 0.125]],\n#                               [[0.0, 0.0, 0.5], [0.25, -0.25, 0.25],\n#                                [0.125, -0.125, 0.125], [0.125, -0.125, 0.125]],\n#                               [[0.0, -0.25, 0.25], [0.0, -0.25, 0.25],\n#                                [0.0, -0.25, 0.25], [0.0, 0.25, -0.25]],\n#                               [[0.0, 0.25, 0.25], [0.0, 0.25, 0.25],\n#                                [0.125, -0.125, -0.125]],\n#                               [[-0.125, 0.125, 0.125], [0.125, -0.125, 0.125],\n#                                [-0.125, -0.125, 0.125]],\n#                               [[-0.125, 0.125, 0.125], [0.125, -0.125, 0.125],\n#                                [-0.125, -0.125, 0.125], [0.125, 0.125, 0.125]],\n#                               [[-0.0, 0.0, 0.5], [-0.25, -0.25, 0.25],\n#                                [-0.125, -0.125, 0.125],\n#                                [-0.125, -0.125, 0.125]],\n#                               [[0.125, 0.125, 0.125], [0.125, -0.125, 0.125],\n#                                [0.125, -0.125, -0.125]],\n#                               [[-0.0, 0.5, 0.0], [-0.25, 0.25, -0.25],\n#                                [0.125, -0.125, 0.125], [0.125, -0.125, 0.125]],\n#                               [[0.125, 0.125, 0.125], [-0.125, -0.125, 0.125],\n#                                [0.125, -0.125, -0.125]],\n#                               [[0.0, -0.25, -0.25], [0.0, 0.25, 0.25],\n#                                [0.125, 0.125, 0.125]],\n#                               [[0.125, 0.125, 0.125], [0.125, -0.125, -0.125]],\n#                               [[0.5, 0.0, -0.0], [0.25, -0.25, -0.25],\n#                                [0.125, -0.125, -0.125]],\n#                               [[-0.25, 0.25, 0.25], [-0.125, 0.125, 0.125],\n#                                [-0.25, 0.25, 0.25], [0.125, -0.125, -0.125]],\n#                               [[0.375, -0.375, 0.375], [0.0, 0.25, 0.25],\n#                                [-0.125, 0.125, -0.125], [-0.25, 0.0, 0.25]],\n#                               [[0.0, -0.5, 0.0], [-0.25, 0.25, 0.25],\n#                                [-0.125, 0.125, 0.125]],\n#                               [[-0.375, -0.375, 0.375], [0.25, -0.25, 0.0],\n#                                [0.0, 0.25, 0.25], [-0.125, -0.125, 0.125]],\n#                               [[-0.125, 0.125, 0.125], [-0.25, 0.25, 0.25],\n#                                [0.0, 0.0, 0.5]],\n#                               [[0.125, 0.125, 0.125], [0.0, 0.25, 0.25],\n#                                [0.0, 0.25, 0.25]],\n#                               [[0.0, 0.25, 0.25], [0.0, 0.25, 0.25]],\n#                               [[0.5, 0.0, -0.0], [0.25, 0.25, 0.25],\n#                                [0.125, 0.125, 0.125], [0.125, 0.125, 0.125]],\n#                               [[0.125, -0.125, 0.125], [-0.125, -0.125, 0.125],\n#                                [0.125, 0.125, 0.125]],\n#                               [[-0.25, -0.0, -0.25], [0.25, 0.0, 0.25],\n#                                [0.125, 0.125, 0.125]],\n#                               [[0.125, 0.125, 0.125], [0.125, -0.125, 0.125]],\n#                               [[-0.25, -0.25, 0.0], [0.25, 0.25, -0.0],\n#                                [0.125, 0.125, 0.125]],\n#                               [[0.125, 0.125, 0.125], [-0.125, -0.125, 0.125]],\n#                               [[0.125, 0.125, 0.125], [0.125, 0.125, 0.125]],\n#                               [[0.125, 0.125, 0.125]], [[0.125, 0.125, 0.125]],\n#                               [[0.125, 0.125, 0.125], [0.125, 0.125, 0.125]],\n#                               [[0.125, 0.125, 0.125], [-0.125, -0.125, 0.125]],\n#                               [[-0.25, -0.25, 0.0], [0.25, 0.25, -0.0],\n#                                [0.125, 0.125, 0.125]],\n#                               [[0.125, 0.125, 0.125], [0.125, -0.125, 0.125]],\n#                               [[-0.25, -0.0, -0.25], [0.25, 0.0, 0.25],\n#                                [0.125, 0.125, 0.125]],\n#                               [[0.125, -0.125, 0.125], [-0.125, -0.125, 0.125],\n#                                [0.125, 0.125, 0.125]],\n#                               [[0.5, 0.0, -0.0], [0.25, 0.25, 0.25],\n#                                [0.125, 0.125, 0.125], [0.125, 0.125, 0.125]],\n#                               [[0.0, 0.25, 0.25], [0.0, 0.25, 0.25]],\n#                               [[0.125, 0.125, 0.125], [0.0, 0.25, 0.25],\n#                                [0.0, 0.25, 0.25]],\n#                               [[-0.125, 0.125, 0.125], [-0.25, 0.25, 0.25],\n#                                [0.0, 0.0, 0.5]],\n#                               [[-0.375, -0.375, 0.375], [0.25, -0.25, 0.0],\n#                                [0.0, 0.25, 0.25], [-0.125, -0.125, 0.125]],\n#                               [[0.0, -0.5, 0.0], [-0.25, 0.25, 0.25],\n#                                [-0.125, 0.125, 0.125]],\n#                               [[0.375, -0.375, 0.375], [0.0, 0.25, 0.25],\n#                                [-0.125, 0.125, -0.125], [-0.25, 0.0, 0.25]],\n#                               [[-0.25, 0.25, 0.25], [-0.125, 0.125, 0.125],\n#                                [-0.25, 0.25, 0.25], [0.125, -0.125, -0.125]],\n#                               [[0.5, 0.0, -0.0], [0.25, -0.25, -0.25],\n#                                [0.125, -0.125, -0.125]],\n#                               [[0.125, 0.125, 0.125], [0.125, -0.125, -0.125]],\n#                               [[0.0, -0.25, -0.25], [0.0, 0.25, 0.25],\n#                                [0.125, 0.125, 0.125]],\n#                               [[0.125, 0.125, 0.125], [-0.125, -0.125, 0.125],\n#                                [0.125, -0.125, -0.125]],\n#                               [[-0.0, 0.5, 0.0], [-0.25, 0.25, -0.25],\n#                                [0.125, -0.125, 0.125], [0.125, -0.125, 0.125]],\n#                               [[0.125, 0.125, 0.125], [0.125, -0.125, 0.125],\n#                                [0.125, -0.125, -0.125]],\n#                               [[-0.0, 0.0, 0.5], [-0.25, -0.25, 0.25],\n#                                [-0.125, -0.125, 0.125],\n#                                [-0.125, -0.125, 0.125]],\n#                               [[-0.125, 0.125, 0.125], [0.125, -0.125, 0.125],\n#                                [-0.125, -0.125, 0.125], [0.125, 0.125, 0.125]],\n#                               [[-0.125, 0.125, 0.125], [0.125, -0.125, 0.125],\n#                                [-0.125, -0.125, 0.125]],\n#                               [[0.0, 0.25, 0.25], [0.0, 0.25, 0.25],\n#                                [0.125, -0.125, -0.125]],\n#                               [[0.0, -0.25, -0.25], [0.0, 0.25, 0.25],\n#                                [0.0, 0.25, 0.25], [0.0, 0.25, 0.25]],\n#                               [[0.0, 0.0, 0.5], [0.25, -0.25, 0.25],\n#                                [0.125, -0.125, 0.125], [0.125, -0.125, 0.125]],\n#                               [[0.0, -0.25, 0.25], [0.0, -0.25, 0.25],\n#                                [0.125, -0.125, 0.125]],\n#                               [[0.0, -0.5, 0.0], [0.125, 0.125, -0.125],\n#                                [0.25, 0.25, -0.25], [-0.125, -0.125, 0.125]],\n#                               [[-0.125, -0.125, 0.125], [0.0, -0.25, 0.25],\n#                                [0.0, 0.25, -0.25]],\n#                               [[0.125, 0.125, 0.125], [-0.125, -0.125, 0.125],\n#                                [0.125, -0.125, 0.125]],\n#                               [[-0.125, -0.125, 0.125], [0.125, -0.125,\n#                                                          0.125]],\n#                               [[0.25, 0.0, 0.25], [0.25, 0.0, 0.25]],\n#                               [[0.125, 0.125, 0.125], [0.25, 0.0, 0.25],\n#                                [0.25, 0.0, 0.25]],\n#                               [[-0.0, 0.0, 0.5], [0.25, -0.25, 0.25],\n#                                [0.125, -0.125, 0.125]],\n#                               [[0.25, 0.0, 0.25], [-0.375, -0.375, 0.375],\n#                                [-0.25, 0.25, 0.0], [-0.125, -0.125, 0.125]],\n#                               [[0.125, -0.125, 0.125], [0.25, 0.0, 0.25],\n#                                [0.25, 0.0, 0.25]],\n#                               [[-0.25, -0.0, -0.25], [0.25, 0.0, 0.25],\n#                                [0.25, 0.0, 0.25], [0.25, 0.0, 0.25]],\n#                               [[-0.0, 0.0, 0.5], [0.25, -0.25, 0.25],\n#                                [0.125, -0.125, 0.125], [0.125, -0.125, 0.125]],\n#                               [[-0.25, 0.0, 0.25], [0.25, 0.0, -0.25],\n#                                [-0.125, 0.125, 0.125]],\n#                               [[0.125, 0.125, 0.125], [0.25, 0.25, 0.25],\n#                                [0.0, 0.0, 0.5]],\n#                               [[0.125, 0.125, 0.125], [0.125, 0.125, 0.125],\n#                                [0.25, 0.25, 0.25], [0.0, 0.0, 0.5]],\n#                               [[-0.0, 0.0, 0.5], [0.0, 0.0, 0.5]],\n#                               [[0.0, 0.0, -0.5], [0.25, 0.25, -0.25],\n#                                [-0.125, -0.125, 0.125]],\n#                               [[-0.25, -0.0, -0.25], [-0.375, 0.375, 0.375],\n#                                [-0.25, -0.25, 0.0], [-0.125, 0.125, 0.125]],\n#                               [[-0.125, -0.125, 0.125], [-0.25, 0.0, 0.25],\n#                                [0.25, 0.0, -0.25]],\n#                               [[-0.0, 0.0, 0.5], [-0.25, 0.25, 0.25],\n#                                [-0.125, 0.125, 0.125]],\n#                               [[-0.25, 0.0, 0.25], [0.25, 0.0, -0.25]],\n#                               [[0.5, 0.0, 0.0], [-0.25, 0.25, -0.25],\n#                                [0.125, -0.125, 0.125]],\n#                               [[-0.25, 0.0, -0.25], [0.375, -0.375, -0.375],\n#                                [0.0, 0.25, -0.25], [-0.125, 0.125, 0.125]],\n#                               [[-0.25, 0.25, -0.25], [-0.25, 0.25, -0.25],\n#                                [-0.125, 0.125, -0.125],\n#                                [-0.125, 0.125, -0.125]],\n#                               [[-0.0, 0.5, 0.0], [-0.25, 0.25, -0.25],\n#                                [0.125, -0.125, 0.125]],\n#                               [[0.5, 0.0, 0.0], [-0.25, -0.25, 0.25],\n#                                [-0.125, -0.125, 0.125],\n#                                [-0.125, -0.125, 0.125]],\n#                               [[-0.125, -0.125, 0.125], [-0.25, 0.0, 0.25],\n#                                [-0.25, 0.0, 0.25]],\n#                               [[0.125, 0.125, 0.125], [-0.125, -0.125, 0.125],\n#                                [-0.125, 0.125, 0.125]],\n#                               [[-0.125, -0.125, 0.125], [-0.125, 0.125,\n#                                                          0.125]],\n#                               [[0.375, -0.375, 0.375], [0.0, -0.25, -0.25],\n#                                [-0.125, 0.125, -0.125], [0.25, 0.25, 0.0]],\n#                               [[0.0, -0.25, 0.25], [0.0, -0.25, 0.25],\n#                                [-0.125, -0.125, 0.125]],\n#                               [[0.0, 0.0, 0.5], [0.25, -0.25, 0.25],\n#                                [0.125, -0.125, 0.125]],\n#                               [[0.0, -0.25, 0.25], [0.0, -0.25, 0.25]],\n#                               [[-0.125, -0.125, 0.125], [-0.25, -0.25, 0.0],\n#                                [0.25, 0.25, -0.0]],\n#                               [[-0.125, -0.125, 0.125],\n#                                [-0.125, -0.125, 0.125]],\n#                               [[0.125, 0.125, 0.125], [-0.125, -0.125, 0.125]],\n#                               [[-0.125, -0.125, 0.125]],\n#                               [[-0.25, -0.25, 0.0], [-0.25, -0.25, 0.0]],\n#                               [[0.125, 0.125, 0.125], [-0.25, -0.25, 0.0],\n#                                [-0.25, -0.25, 0.0]],\n#                               [[-0.25, -0.25, 0.0], [-0.25, -0.25, 0.0],\n#                                [-0.125, -0.125, 0.125]],\n#                               [[-0.25, -0.25, 0.0], [-0.25, -0.25, 0.0],\n#                                [-0.25, -0.25, 0.0], [0.25, 0.25, -0.0]],\n#                               [[0.0, 0.5, 0.0], [0.25, 0.25, -0.25],\n#                                [-0.125, -0.125, 0.125]],\n#                               [[-0.375, 0.375, -0.375], [-0.25, -0.25, 0.0],\n#                                [-0.125, 0.125, -0.125], [-0.25, 0.0, 0.25]],\n#                               [[0.0, 0.5, 0.0], [0.25, 0.25, -0.25],\n#                                [-0.125, -0.125, 0.125],\n#                                [-0.125, -0.125, 0.125]],\n#                               [[-0.125, 0.125, 0.125], [0.25, -0.25, 0.0],\n#                                [-0.25, 0.25, 0.0]],\n#                               [[0.0, -0.5, 0.0], [-0.25, -0.25, -0.25],\n#                                [-0.125, -0.125, -0.125]],\n#                               [[0.125, 0.125, 0.125], [0.0, -0.5, 0.0],\n#                                [-0.25, -0.25, -0.25], [-0.125, -0.125,\n#                                                        -0.125]],\n#                               [[-0.375, -0.375, -0.375], [-0.25, 0.0, 0.25],\n#                                [-0.125, -0.125, -0.125], [-0.25, 0.25, 0.0]],\n#                               [[0.25, -0.25, 0.0], [-0.25, 0.25, 0.0],\n#                                [0.125, -0.125, 0.125]],\n#                               [[0.0, 0.5, 0.0], [0.0, -0.5, 0.0]],\n#                               [[0.0, 0.5, 0.0], [0.125, -0.125, 0.125],\n#                                [-0.25, 0.25, -0.25]],\n#                               [[0.0, 0.5, 0.0], [-0.25, 0.25, 0.25],\n#                                [0.125, -0.125, -0.125]],\n#                               [[0.25, -0.25, 0.0], [-0.25, 0.25, 0.0]],\n#                               [[-0.5, 0.0, 0.0], [-0.25, -0.25, 0.25],\n#                                [-0.125, -0.125, 0.125]],\n#                               [[0.0, 0.25, -0.25], [0.375, -0.375, -0.375],\n#                                [-0.125, 0.125, 0.125], [0.25, 0.25, 0.0]],\n#                               [[0.5, 0.0, 0.0], [0.25, -0.25, 0.25],\n#                                [-0.125, 0.125, -0.125], [0.125, -0.125,\n#                                                          0.125]],\n#                               [[0.125, -0.125, 0.125], [0.25, -0.25, 0.0],\n#                                [0.25, -0.25, 0.0]],\n#                               [[0.25, 0.25, -0.25], [0.25, 0.25, -0.25],\n#                                [0.125, 0.125, -0.125], [-0.125, -0.125,\n#                                                         0.125]],\n#                               [[-0.0, 0.0, 0.5], [-0.25, -0.25, 0.25],\n#                                [-0.125, -0.125, 0.125]],\n#                               [[0.125, 0.125, 0.125], [0.125, -0.125, 0.125],\n#                                [-0.125, 0.125, 0.125]],\n#                               [[-0.125, 0.125, 0.125], [0.125, -0.125, 0.125]],\n#                               [[-0.375, -0.375, 0.375], [-0.0, 0.25, 0.25],\n#                                [0.125, 0.125, -0.125], [-0.25, -0.0, -0.25]],\n#                               [[0.0, -0.25, 0.25], [0.0, 0.25, -0.25],\n#                                [0.125, -0.125, 0.125]],\n#                               [[0.125, -0.125, 0.125], [-0.25, -0.0, -0.25],\n#                                [0.25, 0.0, 0.25]],\n#                               [[0.125, -0.125, 0.125], [0.125, -0.125, 0.125]],\n#                               [[0.0, -0.5, 0.0], [0.125, 0.125, -0.125],\n#                                [0.25, 0.25, -0.25]],\n#                               [[0.0, -0.25, 0.25], [0.0, 0.25, -0.25]],\n#                               [[0.125, 0.125, 0.125], [0.125, -0.125, 0.125]],\n#                               [[0.125, -0.125, 0.125]],\n#                               [[-0.5, 0.0, 0.0], [-0.125, -0.125, -0.125],\n#                                [-0.25, -0.25, -0.25]],\n#                               [[-0.5, 0.0, 0.0], [-0.125, -0.125, -0.125],\n#                                [-0.25, -0.25, -0.25], [0.125, 0.125, 0.125]],\n#                               [[0.375, 0.375, 0.375], [0.0, 0.25, -0.25],\n#                                [-0.125, -0.125, -0.125], [-0.25, 0.25, 0.0]],\n#                               [[0.125, -0.125, -0.125], [0.25, -0.25, 0.0],\n#                                [0.25, -0.25, 0.0]],\n#                               [[0.125, 0.125, 0.125], [0.375, 0.375, 0.375],\n#                                [0.0, -0.25, 0.25], [-0.25, 0.0, 0.25]],\n#                               [[-0.25, 0.0, 0.25], [-0.25, 0.0, 0.25],\n#                                [0.125, -0.125, -0.125]],\n#                               [[0.0, -0.25, -0.25], [0.0, 0.25, 0.25],\n#                                [-0.125, 0.125, 0.125]],\n#                               [[-0.125, 0.125, 0.125], [0.125, -0.125,\n#                                                         -0.125]],\n#                               [[-0.125, -0.125, -0.125], [-0.25, -0.25, -0.25],\n#                                [0.25, 0.25, 0.25], [0.125, 0.125, 0.125]],\n#                               [[-0.125, -0.125, 0.125], [0.125, -0.125, 0.125],\n#                                [0.125, -0.125, -0.125]],\n#                               [[0.0, 0.0, -0.5], [0.25, 0.25, 0.25],\n#                                [-0.125, -0.125, -0.125]],\n#                               [[0.125, -0.125, 0.125], [0.125, -0.125,\n#                                                         -0.125]],\n#                               [[0.0, -0.5, 0.0], [0.25, 0.25, 0.25],\n#                                [0.125, 0.125, 0.125]],\n#                               [[-0.125, -0.125, 0.125],\n#                                [0.125, -0.125, -0.125]],\n#                               [[0.0, -0.25, -0.25], [0.0, 0.25, 0.25]],\n#                               [[0.125, -0.125, -0.125]],\n#                               [[0.5, 0.0, 0.0], [0.5, 0.0, 0.0]],\n#                               [[-0.5, 0.0, 0.0], [-0.25, 0.25, 0.25],\n#                                [-0.125, 0.125, 0.125]],\n#                               [[0.5, 0.0, 0.0], [0.25, -0.25, 0.25],\n#                                [-0.125, 0.125, -0.125]],\n#                               [[0.25, -0.25, 0.0], [0.25, -0.25, 0.0]],\n#                               [[0.5, 0.0, 0.0], [-0.25, -0.25, 0.25],\n#                                [-0.125, -0.125, 0.125]],\n#                               [[-0.25, 0.0, 0.25], [-0.25, 0.0, 0.25]],\n#                               [[0.125, 0.125, 0.125], [-0.125, 0.125, 0.125]],\n#                               [[-0.125, 0.125, 0.125]],\n#                               [[0.5, 0.0, -0.0], [0.25, 0.25, 0.25],\n#                                [0.125, 0.125, 0.125]],\n#                               [[0.125, -0.125, 0.125], [-0.125, -0.125,\n#                                                         0.125]],\n#                               [[-0.25, -0.0, -0.25], [0.25, 0.0, 0.25]],\n#                               [[0.125, -0.125, 0.125]],\n#                               [[-0.25, -0.25, 0.0], [0.25, 0.25, -0.0]],\n#                               [[-0.125, -0.125, 0.125]],\n#                               [[0.125, 0.125, 0.125]], [[0, 0, 0]]]\n# # pylint: enable=line-too-long\n\n\n# def create_table_neighbour_code_to_surface_area(spacing_mm):\n#     \"\"\"Returns an array mapping neighbourhood code to the surface elements area.\n\n#   Note that the normals encode the initial surface area. This function computes\n#   the area corresponding to the given `spacing_mm`.\n\n#   Args:\n#     spacing_mm: 3-element list-like structure. Voxel spacing in x0, x1 and x2\n#       direction.\n#   \"\"\"\n#     # compute the area for all 256 possible surface elements\n#     # (given a 2x2x2 neighbourhood) according to the spacing_mm\n#     neighbour_code_to_surface_area = np.zeros([256])\n#     for code in range(256):\n#         normals = np.array(_NEIGHBOUR_CODE_TO_NORMALS[code])\n#         sum_area = 0\n#         for normal_idx in range(normals.shape[0]):\n#             # normal vector\n#             n = np.zeros([3])\n#             n[0] = normals[normal_idx, 0] * spacing_mm[1] * spacing_mm[2]\n#             n[1] = normals[normal_idx, 1] * spacing_mm[0] * spacing_mm[2]\n#             n[2] = normals[normal_idx, 2] * spacing_mm[0] * spacing_mm[1]\n#             area = np.linalg.norm(n)\n#             sum_area += area\n#         neighbour_code_to_surface_area[code] = sum_area\n\n#     return neighbour_code_to_surface_area\n\n\n# # In the neighbourhood, points are ordered: top left, top right, bottom left,\n# # bottom right.\n# ENCODE_NEIGHBOURHOOD_2D_KERNEL = np.array([[8, 4], [2, 1]])\n\n\n# def create_table_neighbour_code_to_contour_length(spacing_mm):\n#     \"\"\"Returns an array mapping neighbourhood code to the contour length.\n\n#   For the list of possible cases and their figures, see page 38 from:\n#   https://nccastaff.bournemouth.ac.uk/jmacey/MastersProjects/MSc14/06/thesis.pdf\n\n#   In 2D, each point has 4 neighbors. Thus, are 16 configurations. A\n#   configuration is encoded with '1' meaning \"inside the object\" and '0' \"outside\n#   the object\". The points are ordered: top left, top right, bottom left, bottom\n#   right.\n\n#   The x0 axis is assumed vertical downward, and the x1 axis is horizontal to the\n#   right:\n#    (0, 0) --> (0, 1)\n#      |\n#    (1, 0)\n\n#   Args:\n#     spacing_mm: 2-element list-like structure. Voxel spacing in x0 and x1\n#       directions.\n#   \"\"\"\n#     neighbour_code_to_contour_length = np.zeros([16])\n\n#     vertical = spacing_mm[0]\n#     horizontal = spacing_mm[1]\n#     diag = 0.5 * math.sqrt(spacing_mm[0]**2 + spacing_mm[1]**2)\n#     # pyformat: disable\n#     neighbour_code_to_contour_length[int(\"00\" \"01\", 2)] = diag\n\n#     neighbour_code_to_contour_length[int(\"00\" \"10\", 2)] = diag\n\n#     neighbour_code_to_contour_length[int(\"00\" \"11\", 2)] = horizontal\n\n#     neighbour_code_to_contour_length[int(\"01\" \"00\", 2)] = diag\n\n#     neighbour_code_to_contour_length[int(\"01\" \"01\", 2)] = vertical\n\n#     neighbour_code_to_contour_length[int(\"01\" \"10\", 2)] = 2 * diag\n\n#     neighbour_code_to_contour_length[int(\"01\" \"11\", 2)] = diag\n\n#     neighbour_code_to_contour_length[int(\"10\" \"00\", 2)] = diag\n\n#     neighbour_code_to_contour_length[int(\"10\" \"01\", 2)] = 2 * diag\n\n#     neighbour_code_to_contour_length[int(\"10\" \"10\", 2)] = vertical\n\n#     neighbour_code_to_contour_length[int(\"10\" \"11\", 2)] = diag\n\n#     neighbour_code_to_contour_length[int(\"11\" \"00\", 2)] = horizontal\n\n#     neighbour_code_to_contour_length[int(\"11\" \"01\", 2)] = diag\n\n#     neighbour_code_to_contour_length[int(\"11\" \"10\", 2)] = diag\n#     # pyformat: enable\n\n#     return neighbour_code_to_contour_length\n\n\n# # surface_distance/metrics.py\n# \"\"\"Module exposing surface distance based measures.\"\"\"\n\n# import numpy as np\n# from scipy import ndimage\n\n\n# def _assert_is_numpy_array(name, array):\n#     \"\"\"Raises an exception if `array` is not a numpy array.\"\"\"\n#     if not isinstance(array, np.ndarray):\n#         raise ValueError(\"The argument {!r} should be a numpy array, not a \"\n#                          \"{}\".format(name, type(array)))\n\n\n# def _check_nd_numpy_array(name, array, num_dims):\n#     \"\"\"Raises an exception if `array` is not a `num_dims`-D numpy array.\"\"\"\n#     if len(array.shape) != num_dims:\n#         raise ValueError(\"The argument {!r} should be a {}D array, not of \"\n#                          \"shape {}\".format(name, num_dims, array.shape))\n\n\n# def _check_2d_numpy_array(name, array):\n#     _check_nd_numpy_array(name, array, num_dims=2)\n\n\n# def _check_3d_numpy_array(name, array):\n#     _check_nd_numpy_array(name, array, num_dims=3)\n\n\n# def _assert_is_bool_numpy_array(name, array):\n#     _assert_is_numpy_array(name, array)\n#     if array.dtype != bool:\n#         raise ValueError(\n#             \"The argument {!r} should be a numpy array of type bool, \"\n#             \"not {}\".format(name, array.dtype))\n\n\n# def _compute_bounding_box(mask):\n#     \"\"\"Computes the bounding box of the masks.\n\n#   This function generalizes to arbitrary number of dimensions great or equal\n#   to 1.\n\n#   Args:\n#     mask: The 2D or 3D numpy mask, where '0' means background and non-zero means\n#       foreground.\n\n#   Returns:\n#     A tuple:\n#      - The coordinates of the first point of the bounding box (smallest on all\n#        axes), or `None` if the mask contains only zeros.\n#      - The coordinates of the second point of the bounding box (greatest on all\n#        axes), or `None` if the mask contains only zeros.\n#   \"\"\"\n#     num_dims = len(mask.shape)\n#     bbox_min = np.zeros(num_dims, np.int64)\n#     bbox_max = np.zeros(num_dims, np.int64)\n\n#     # max projection to the x0-axis\n#     proj_0 = np.amax(mask, axis=tuple(range(num_dims))[1:])\n#     idx_nonzero_0 = np.nonzero(proj_0)[0]\n#     if len(idx_nonzero_0) == 0:  # pylint: disable=g-explicit-length-test\n#         return None, None\n\n#     bbox_min[0] = np.min(idx_nonzero_0)\n#     bbox_max[0] = np.max(idx_nonzero_0)\n\n#     # max projection to the i-th-axis for i in {1, ..., num_dims - 1}\n#     for axis in range(1, num_dims):\n#         max_over_axes = list(range(num_dims))  # Python 3 compatible\n#         max_over_axes.pop(axis)  # Remove the i-th dimension from the max\n#         max_over_axes = tuple(max_over_axes)  # numpy expects a tuple of ints\n#         proj = np.amax(mask, axis=max_over_axes)\n#         idx_nonzero = np.nonzero(proj)[0]\n#         bbox_min[axis] = np.min(idx_nonzero)\n#         bbox_max[axis] = np.max(idx_nonzero)\n\n#     return bbox_min, bbox_max\n\n\n# def _crop_to_bounding_box(mask, bbox_min, bbox_max):\n#     \"\"\"Crops a 2D or 3D mask to the bounding box specified by `bbox_{min,max}`.\"\"\"\n#     # we need to zeropad the cropped region with 1 voxel at the lower,\n#     # the right (and the back on 3D) sides. This is required to obtain the\n#     # \"full\" convolution result with the 2x2 (or 2x2x2 in 3D) kernel.\n#     # TODO:  This is correct only if the object is interior to the\n#     # bounding box.\n#     cropmask = np.zeros((bbox_max - bbox_min) + 2, np.uint8)\n\n#     num_dims = len(mask.shape)\n#     # pyformat: disable\n#     if num_dims == 2:\n#         cropmask[0:-1, 0:-1] = mask[bbox_min[0]:bbox_max[0] + 1,\n#                                     bbox_min[1]:bbox_max[1] + 1]\n#     elif num_dims == 3:\n#         cropmask[0:-1, 0:-1, 0:-1] = mask[bbox_min[0]:bbox_max[0] + 1,\n#                                           bbox_min[1]:bbox_max[1] + 1,\n#                                           bbox_min[2]:bbox_max[2] + 1]\n#     # pyformat: enable\n#     else:\n#         assert False\n\n#     return cropmask\n\n\n# def _sort_distances_surfels(distances, surfel_areas):\n#     \"\"\"Sorts the two list with respect to the tuple of (distance, surfel_area).\n\n#   Args:\n#     distances: The distances from A to B (e.g. `distances_gt_to_pred`).\n#     surfel_areas: The surfel areas for A (e.g. `surfel_areas_gt`).\n\n#   Returns:\n#     A tuple of the sorted (distances, surfel_areas).\n#   \"\"\"\n#     sorted_surfels = np.array(sorted(zip(distances, surfel_areas)))\n#     return sorted_surfels[:, 0], sorted_surfels[:, 1]\n\n\n# def compute_surface_distances(mask_gt, mask_pred, spacing_mm):\n#     \"\"\"Computes closest distances from all surface points to the other surface.\n\n#   This function can be applied to 2D or 3D tensors. For 2D, both masks must be\n#   2D and `spacing_mm` must be a 2-element list. For 3D, both masks must be 3D\n#   and `spacing_mm` must be a 3-element list. The description is done for the 2D\n#   case, and the formulation for the 3D case is present is parenthesis,\n#   introduced by \"resp.\".\n\n#   Finds all contour elements (resp surface elements \"surfels\" in 3D) in the\n#   ground truth mask `mask_gt` and the predicted mask `mask_pred`, computes their\n#   length in mm (resp. area in mm^2) and the distance to the closest point on the\n#   other contour (resp. surface). It returns two sorted lists of distances\n#   together with the corresponding contour lengths (resp. surfel areas). If one\n#   of the masks is empty, the corresponding lists are empty and all distances in\n#   the other list are `inf`.\n\n#   Args:\n#     mask_gt: 2-dim (resp. 3-dim) bool Numpy array. The ground truth mask.\n#     mask_pred: 2-dim (resp. 3-dim) bool Numpy array. The predicted mask.\n#     spacing_mm: 2-element (resp. 3-element) list-like structure. Voxel spacing\n#       in x0 anx x1 (resp. x0, x1 and x2) directions.\n\n#   Returns:\n#     A dict with:\n#     \"distances_gt_to_pred\": 1-dim numpy array of type float. The distances in mm\n#         from all ground truth surface elements to the predicted surface,\n#         sorted from smallest to largest.\n#     \"distances_pred_to_gt\": 1-dim numpy array of type float. The distances in mm\n#         from all predicted surface elements to the ground truth surface,\n#         sorted from smallest to largest.\n#     \"surfel_areas_gt\": 1-dim numpy array of type float. The length of the\n#       of the ground truth contours in mm (resp. the surface elements area in\n#       mm^2) in the same order as distances_gt_to_pred.\n#     \"surfel_areas_pred\": 1-dim numpy array of type float. The length of the\n#       of the predicted contours in mm (resp. the surface elements area in\n#       mm^2) in the same order as distances_gt_to_pred.\n\n#   Raises:\n#     ValueError: If the masks and the `spacing_mm` arguments are of incompatible\n#       shape or type. Or if the masks are not 2D or 3D.\n#   \"\"\"\n#     # The terms used in this function are for the 3D case. In particular, surface\n#     # in 2D stands for contours in 3D. The surface elements in 3D correspond to\n#     # the line elements in 2D.\n\n#     _assert_is_bool_numpy_array(\"mask_gt\", mask_gt)\n#     _assert_is_bool_numpy_array(\"mask_pred\", mask_pred)\n\n#     if not len(mask_gt.shape) == len(mask_pred.shape) == len(spacing_mm):\n#         raise ValueError(\n#             \"The arguments must be of compatible shape. Got mask_gt \"\n#             \"with {} dimensions ({}) and mask_pred with {} dimensions \"\n#             \"({}), while the spacing_mm was {} elements.\".format(\n#                 len(mask_gt.shape), mask_gt.shape, len(mask_pred.shape),\n#                 mask_pred.shape, len(spacing_mm)))\n\n#     num_dims = len(spacing_mm)\n#     if num_dims == 2:\n#         _check_2d_numpy_array(\"mask_gt\", mask_gt)\n#         _check_2d_numpy_array(\"mask_pred\", mask_pred)\n\n#         # compute the area for all 16 possible surface elements\n#         # (given a 2x2 neighbourhood) according to the spacing_mm\n#         neighbour_code_to_surface_area = (\n#             create_table_neighbour_code_to_contour_length(spacing_mm))\n#         kernel = ENCODE_NEIGHBOURHOOD_2D_KERNEL\n#         full_true_neighbours = 0b1111\n#     elif num_dims == 3:\n#         _check_3d_numpy_array(\"mask_gt\", mask_gt)\n#         _check_3d_numpy_array(\"mask_pred\", mask_pred)\n\n#         # compute the area for all 256 possible surface elements\n#         # (given a 2x2x2 neighbourhood) according to the spacing_mm\n#         neighbour_code_to_surface_area = (\n#             create_table_neighbour_code_to_surface_area(spacing_mm))\n#         kernel = ENCODE_NEIGHBOURHOOD_3D_KERNEL\n#         full_true_neighbours = 0b11111111\n#     else:\n#         raise ValueError(\"Only 2D and 3D masks are supported, not \"\n#                          \"{}D.\".format(num_dims))\n\n#     # compute the bounding box of the masks to trim the volume to the smallest\n#     # possible processing subvolume\n#     bbox_min, bbox_max = _compute_bounding_box(mask_gt | mask_pred)\n#     # Both the min/max bbox are None at the same time, so we only check one.\n#     if bbox_min is None:\n#         return {\n#             \"distances_gt_to_pred\": np.array([]),\n#             \"distances_pred_to_gt\": np.array([]),\n#             \"surfel_areas_gt\": np.array([]),\n#             \"surfel_areas_pred\": np.array([]),\n#         }\n\n#     # crop the processing subvolume.\n#     cropmask_gt = _crop_to_bounding_box(mask_gt, bbox_min, bbox_max)\n#     cropmask_pred = _crop_to_bounding_box(mask_pred, bbox_min, bbox_max)\n\n#     # compute the neighbour code (local binary pattern) for each voxel\n#     # the resulting arrays are spacially shifted by minus half a voxel in each\n#     # axis.\n#     # i.e. the points are located at the corners of the original voxels\n#     neighbour_code_map_gt = ndimage.correlate(cropmask_gt.astype(np.uint8),\n#                                               kernel,\n#                                               mode=\"constant\",\n#                                               cval=0)\n#     neighbour_code_map_pred = ndimage.correlate(cropmask_pred.astype(np.uint8),\n#                                                 kernel,\n#                                                 mode=\"constant\",\n#                                                 cval=0)\n\n#     # create masks with the surface voxels\n#     borders_gt = ((neighbour_code_map_gt != 0) &\n#                   (neighbour_code_map_gt != full_true_neighbours))\n#     borders_pred = ((neighbour_code_map_pred != 0) &\n#                     (neighbour_code_map_pred != full_true_neighbours))\n\n#     # compute the distance transform (closest distance of each voxel to the\n#     # surface voxels)\n#     if borders_gt.any():\n#         distmap_gt = distance_transform_edt(~borders_gt)\n\n#     else:\n#         distmap_gt = np.Inf * np.ones(borders_gt.shape)\n\n#     distances_pred_to_gt = distmap_gt[borders_pred]\n#     del distmap_gt\n\n#     if borders_pred.any():\n#         distmap_pred = distance_transform_edt(~borders_pred)\n\n#     else:\n#         distmap_pred = np.Inf * np.ones(borders_pred.shape)\n\n#     distances_gt_to_pred = distmap_pred[borders_gt]\n#     del distmap_pred\n\n#     # compute the area of each surface element\n#     surface_area_map_gt = neighbour_code_to_surface_area[neighbour_code_map_gt]\n#     surfel_areas_gt = surface_area_map_gt[borders_gt]\n#     del neighbour_code_map_gt, surface_area_map_gt, borders_gt\n\n#     surface_area_map_pred = neighbour_code_to_surface_area[neighbour_code_map_pred]\n#     surfel_areas_pred = surface_area_map_pred[borders_pred]\n#     del neighbour_code_map_pred, surface_area_map_pred, borders_pred\n\n#     # sort them by distance\n#     if distances_gt_to_pred.shape != (0, ):\n#         distances_gt_to_pred, surfel_areas_gt = _sort_distances_surfels(\n#             distances_gt_to_pred, surfel_areas_gt)\n\n#     if distances_pred_to_gt.shape != (0, ):\n#         distances_pred_to_gt, surfel_areas_pred = _sort_distances_surfels(\n#             distances_pred_to_gt, surfel_areas_pred)\n\n#     return {\n#         \"distances_gt_to_pred\": distances_gt_to_pred,\n#         \"distances_pred_to_gt\": distances_pred_to_gt,\n#         \"surfel_areas_gt\": surfel_areas_gt,\n#         \"surfel_areas_pred\": surfel_areas_pred,\n#     }\n\n\n# def compute_average_surface_distance(surface_distances):\n#     \"\"\"Returns the average surface distance.\n\n#   Computes the average surface distances by correctly taking the area of each\n#   surface element into account. Call compute_surface_distances(...) before, to\n#   obtain the `surface_distances` dict.\n\n#   Args:\n#     surface_distances: dict with \"distances_gt_to_pred\", \"distances_pred_to_gt\"\n#     \"surfel_areas_gt\", \"surfel_areas_pred\" created by\n#     compute_surface_distances()\n\n#   Returns:\n#     A tuple with two float values:\n#       - the average distance (in mm) from the ground truth surface to the\n#         predicted surface\n#       - the average distance from the predicted surface to the ground truth\n#         surface.\n#   \"\"\"\n#     distances_gt_to_pred = surface_distances[\"distances_gt_to_pred\"]\n#     distances_pred_to_gt = surface_distances[\"distances_pred_to_gt\"]\n#     surfel_areas_gt = surface_distances[\"surfel_areas_gt\"]\n#     surfel_areas_pred = surface_distances[\"surfel_areas_pred\"]\n#     average_distance_gt_to_pred = (\n#         np.sum(distances_gt_to_pred * surfel_areas_gt) /\n#         np.sum(surfel_areas_gt))\n#     average_distance_pred_to_gt = (\n#         np.sum(distances_pred_to_gt * surfel_areas_pred) /\n#         np.sum(surfel_areas_pred))\n#     return (average_distance_gt_to_pred, average_distance_pred_to_gt)\n\n\n# def compute_robust_hausdorff(surface_distances, percent):\n#     \"\"\"Computes the robust Hausdorff distance.\n\n#   Computes the robust Hausdorff distance. \"Robust\", because it uses the\n#   `percent` percentile of the distances instead of the maximum distance. The\n#   percentage is computed by correctly taking the area of each surface element\n#   into account.\n\n#   Args:\n#     surface_distances: dict with \"distances_gt_to_pred\", \"distances_pred_to_gt\"\n#       \"surfel_areas_gt\", \"surfel_areas_pred\" created by\n#       compute_surface_distances()\n#     percent: a float value between 0 and 100.\n\n#   Returns:\n#     a float value. The robust Hausdorff distance in mm.\n#   \"\"\"\n#     distances_gt_to_pred = surface_distances[\"distances_gt_to_pred\"]\n#     distances_pred_to_gt = surface_distances[\"distances_pred_to_gt\"]\n#     surfel_areas_gt = surface_distances[\"surfel_areas_gt\"]\n#     surfel_areas_pred = surface_distances[\"surfel_areas_pred\"]\n#     if len(distances_gt_to_pred) > 0:  # pylint: disable=g-explicit-length-test\n#         surfel_areas_cum_gt = np.cumsum(surfel_areas_gt) / np.sum(\n#             surfel_areas_gt)\n#         idx = np.searchsorted(surfel_areas_cum_gt, percent / 100.0)\n#         perc_distance_gt_to_pred = distances_gt_to_pred[min(\n#             idx,\n#             len(distances_gt_to_pred) - 1)]\n#     else:\n#         perc_distance_gt_to_pred = np.Inf\n\n#     if len(distances_pred_to_gt) > 0:  # pylint: disable=g-explicit-length-test\n#         surfel_areas_cum_pred = (np.cumsum(surfel_areas_pred) /\n#                                  np.sum(surfel_areas_pred))\n#         idx = np.searchsorted(surfel_areas_cum_pred, percent / 100.0)\n#         perc_distance_pred_to_gt = distances_pred_to_gt[min(\n#             idx,\n#             len(distances_pred_to_gt) - 1)]\n#     else:\n#         perc_distance_pred_to_gt = np.Inf\n\n#     return max(perc_distance_gt_to_pred, perc_distance_pred_to_gt)\n\n\n# def compute_surface_overlap_at_tolerance(surface_distances, tolerance_mm):\n#     \"\"\"Computes the overlap of the surfaces at a specified tolerance.\n\n#   Computes the overlap of the ground truth surface with the predicted surface\n#   and vice versa allowing a specified tolerance (maximum surface-to-surface\n#   distance that is regarded as overlapping). The overlapping fraction is\n#   computed by correctly taking the area of each surface element into account.\n\n#   Args:\n#     surface_distances: dict with \"distances_gt_to_pred\", \"distances_pred_to_gt\"\n#       \"surfel_areas_gt\", \"surfel_areas_pred\" created by\n#       compute_surface_distances()\n#     tolerance_mm: a float value. The tolerance in mm\n\n#   Returns:\n#     A tuple of two float values. The overlap fraction in [0.0, 1.0] of the\n#     ground truth surface with the predicted surface and vice versa.\n#   \"\"\"\n#     distances_gt_to_pred = surface_distances[\"distances_gt_to_pred\"]\n#     distances_pred_to_gt = surface_distances[\"distances_pred_to_gt\"]\n#     surfel_areas_gt = surface_distances[\"surfel_areas_gt\"]\n#     surfel_areas_pred = surface_distances[\"surfel_areas_pred\"]\n#     rel_overlap_gt = (\n#         np.sum(surfel_areas_gt[distances_gt_to_pred <= tolerance_mm]) /\n#         np.sum(surfel_areas_gt))\n#     rel_overlap_pred = (\n#         np.sum(surfel_areas_pred[distances_pred_to_gt <= tolerance_mm]) /\n#         np.sum(surfel_areas_pred))\n#     return (rel_overlap_gt, rel_overlap_pred)\n\n\n# def compute_surface_dice_at_tolerance(surface_distances, tolerance_mm):\n#     \"\"\"Computes the _surface_ DICE coefficient at a specified tolerance.\n\n#   Computes the _surface_ DICE coefficient at a specified tolerance. Not to be\n#   confused with the standard _volumetric_ DICE coefficient. The surface DICE\n#   measures the overlap of two surfaces instead of two volumes. A surface\n#   element is counted as overlapping (or touching), when the closest distance to\n#   the other surface is less or equal to the specified tolerance. The DICE\n#   coefficient is in the range between 0.0 (no overlap) to 1.0 (perfect overlap).\n\n#   Args:\n#     surface_distances: dict with \"distances_gt_to_pred\", \"distances_pred_to_gt\"\n#       \"surfel_areas_gt\", \"surfel_areas_pred\" created by\n#       compute_surface_distances()\n#     tolerance_mm: a float value. The tolerance in mm\n\n#   Returns:\n#     A float value. The surface DICE coefficient in [0.0, 1.0].\n#   \"\"\"\n#     distances_gt_to_pred = surface_distances[\"distances_gt_to_pred\"]\n#     distances_pred_to_gt = surface_distances[\"distances_pred_to_gt\"]\n#     surfel_areas_gt = surface_distances[\"surfel_areas_gt\"]\n#     surfel_areas_pred = surface_distances[\"surfel_areas_pred\"]\n#     overlap_gt = np.sum(surfel_areas_gt[distances_gt_to_pred <= tolerance_mm])\n#     overlap_pred = np.sum(\n#         surfel_areas_pred[distances_pred_to_gt <= tolerance_mm])\n#     surface_dice = (overlap_gt + overlap_pred) / (np.sum(surfel_areas_gt) +\n#                                                   np.sum(surfel_areas_pred))\n#     return surface_dice\n\n\n# def compute_dice_coefficient(mask_gt, mask_pred):\n#     \"\"\"Computes soerensen-dice coefficient.\n\n#   compute the soerensen-dice coefficient between the ground truth mask `mask_gt`\n#   and the predicted mask `mask_pred`.\n\n#   Args:\n#     mask_gt: 3-dim Numpy array of type bool. The ground truth mask.\n#     mask_pred: 3-dim Numpy array of type bool. The predicted mask.\n\n#   Returns:\n#     the dice coeffcient as float. If both masks are empty, the result is NaN.\n#   \"\"\"\n#     volume_sum = mask_gt.sum() + mask_pred.sum()\n#     if volume_sum == 0:\n#         return np.NaN\n#     volume_intersect = (mask_gt & mask_pred).sum()\n#     return 2 * volume_intersect / volume_sum\n\n\n# # Below is from https://github.com/Image-Py/imagepy/blob/master/imagepy/ipyalg/hydrology/edt.py\n# # Copyright (c) <2016>, <ImagePy>\n# # All rights reserved.\n\n# # Redistribution and use in source and binary forms, with or without\n# # modification, are permitted provided that the following conditions are met:\n# # 1. Redistributions of source code must retain the above copyright\n# #    notice, this list of conditions and the following disclaimer.\n# # 2. Redistributions in binary form must reproduce the above copyright\n# #    notice, this list of conditions and the following disclaimer in the\n# #    documentation and/or other materials provided with the distribution.\n# # 3. All advertising materials mentioning features or use of this software\n# #    must display the following acknowledgement:\n# #    This product includes software developed by the <organization>.\n# # 4. Neither the name of the <organization> nor the\n# #    names of its contributors may be used to endorse or promote products\n# #    derived from this software without specific prior written permission.\n\n# # THIS SOFTWARE IS PROVIDED BY <ImagePy> ''AS IS'' AND ANY\n# # EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED\n# # WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE\n# # DISCLAIMED. IN NO EVENT SHALL <ImagePy> BE LIABLE FOR ANY\n# # DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES\n# # (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES;\n# # LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND\n# # ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT\n# # (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS\n# # SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.\n\n# def neighbors(shape):\n#     dim = len(shape)\n#     block = generate_binary_structure(dim, 1)\n#     block[tuple([1]*dim)] = 0\n#     idx = np.where(block>0)\n#     idx = np.array(idx, dtype=np.uint8).T\n#     idx = np.array(idx-[1]*dim)\n#     acc = np.cumprod((1,)+shape[::-1][:-1])\n#     return np.dot(idx, acc[::-1])\n\n\n# @jit(nopython=True)\n# def dist(idx1, idx2, acc):\n#     dis = 0\n#     for i in range(len(acc)):\n#         c1 = idx1//acc[i]\n#         c2 = idx2//acc[i]\n#         dis += (c1-c2)**2\n#         idx1 -= c1*acc[i]\n#         idx2 -= c2*acc[i]\n#     return dis\n\n\n# @jit(nopython=True)\n# def _step(dis, pts, roots, s, level, nbs, acc, scale):\n#     cur = 0\n#     while cur<s:\n#         p = pts[cur]\n#         rp = roots[cur]\n#         if dis[p] == 0xffff or dis[p]>level:\n#             cur += 1\n#             continue\n#         for dp in nbs:\n#             cp = p+dp\n#             if dis[cp]<level*scale+2e-10:continue\n#             if dis[cp]==0xffff:continue\n#             tdist = dist(cp, rp, acc)\n#             if tdist<dis[cp]**2-1e-10:\n#                 dis[cp] = (tdist**0.5)*scale\n#                 pts[s] = cp\n#                 roots[s] = rp\n#                 if s == len(pts):\n#                     s, cur = clear(pts, roots, s, cur)\n#                 s+=1\n#         pts[cur] = -1\n#         cur+=1\n#     return cur\n\n\n# @jit(nopython=True)\n# def clear(pts, roots, s, cur):\n#     ns = 0; nc=0;\n#     for c in range(s):\n#         if pts[c]!=-1:\n#             pts[ns] = pts[c]\n#             roots[ns] = roots[c]\n#             ns += 1\n#             if c<cur:nc += 1\n#     return ns, nc\n\n        \n# @jit(nopython=True)\n# def collect(dis, nbs, pts, root):\n#     cur = 0\n#     for p in range(len(dis)):\n#         if dis[p]>=0xffff-1: continue # edge or back\n#         for dp in nbs:\n#             if dis[p+dp]==0xffff-1:\n#                 pts[cur] = p\n#                 root[cur] = p\n#                 cur += 1\n#                 break\n#     return cur\n\n\n# @jit(nopython=True)\n# def bufjit(line):\n#     for i in range(len(line)):\n#         line[i] = 1 if line[i]==0 else 0xffff\n\n\n# def buffer(img, dtype):\n#     buf = np.ones(tuple(np.array(img.shape)+2), dtype=dtype)\n#     buf[tuple([slice(1,-1)]*buf.ndim)] = img\n#     bufjit(buf.ravel())\n#     buf[tuple([slice(1,-1)]*buf.ndim)] -= 1\n#     return buf\n\n\n# def distance_transform_edt(img, output=np.float32, scale=1):\n#     dis = buffer(img, output)\n#     nbs = neighbors(dis.shape)\n#     acc = np.cumprod((1,)+dis.shape[::-1][:-1])[::-1]\n#     line = dis.ravel()\n#     pts = np.zeros(max(line.size//4, 1024**2), dtype=np.int64)\n#     roots = np.zeros(max(line.size//4, 1024**2), dtype=np.int64)\n#     s = collect(line, nbs, pts, roots)\n#     for level in range(10000):\n#         s, c = clear(pts, roots, s, 0)\n#         s = _step(line, pts, roots, s, level, nbs, acc, scale)\n#         if s==0:break\n#     return dis[(slice(1,-1),)*img.ndim]","metadata":{},"execution_count":null,"outputs":[]}]}