{"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":"none","dataSources":[{"sourceId":61446,"databundleVersionId":6962461,"sourceType":"competition"}],"dockerImageVersionId":30626,"isInternetEnabled":true,"language":"python","sourceType":"notebook","isGpuEnabled":false}},"nbformat_minor":4,"nbformat":4,"cells":[{"cell_type":"markdown","source":"You can make any volume slice with any rotation and translation.\nYou can also make a 3D slice this way.","metadata":{}},{"cell_type":"code","source":"import typing as tp\nfrom pathlib import Path\n\nimport numpy as np\nimport cv2\n\nfrom tqdm import tqdm\n\n\ndef load_volume(folder_path: tp.Union[Path, str], is_mask: bool) -> np.ndarray:\n    folder_path = Path(folder_path)\n    \n    assert folder_path.is_dir()\n    \n    n_slices = 0\n    \n    slice_paths: tp.Dict[int, Path] = {}\n    for slice_path in folder_path.iterdir():\n        if slice_path.is_dir():\n            continue\n        if slice_path.suffix != '.tif':\n            continue\n        slice_idx = int(slice_path.stem)\n        if slice_idx >= n_slices:\n            n_slices = slice_idx + 1\n        slice_paths[slice_idx] = slice_path\n\n    slices: tp.List[np.ndarray] = []\n    for z in tqdm(range(n_slices)):\n        if z not in slice_paths:\n            raise RuntimeError(f'Missing file: {z}.tif')\n        image = cv2.imread(str(slice_paths[z]), cv2.IMREAD_ANYDEPTH)\n        image = image.astype(np.uint8 if is_mask else np.uint16)\n        slices.append(image)\n    volume = np.stack(slices)\n    \n    if not is_mask:\n        volume //= 255\n        volume[volume > 255] = 255\n        volume = volume.astype(np.uint8)\n    \n    return volume\n\ndef downsample_volume(volume: np.ndarray, maxval: float = 255) -> np.ndarray:\n    dtype = volume.dtype\n    if volume.shape[2] % 2 > 0:\n        volume = np.concatenate([volume, volume[:, :, -1:]], axis=2)\n    if volume.shape[1] % 2 > 0:\n        volume = np.concatenate([volume, volume[:, -1:, :]], axis=1)\n    if volume.shape[0] % 2 > 0:\n        volume = np.concatenate([volume, volume[-1:, :, :]], axis=0)\n    parts = [volume[0::2, 0::2, 0::2].astype(np.float32),\n             volume[1::2, 0::2, 0::2].astype(np.float32),\n             volume[0::2, 1::2, 0::2].astype(np.float32),\n             volume[1::2, 1::2, 0::2].astype(np.float32),\n             volume[0::2, 0::2, 1::2].astype(np.float32),\n             volume[1::2, 0::2, 1::2].astype(np.float32),\n             volume[0::2, 1::2, 1::2].astype(np.float32),\n             volume[1::2, 1::2, 1::2].astype(np.float32)]\n    new_volume = sum(parts[1:], parts[0]) * (1 / 8)\n    new_volume[new_volume > maxval] = maxval\n    new_volume = new_volume.astype(dtype)\n    return new_volume","metadata":{"_uuid":"8f2839f25d086af736a60e9eeb907d3b93b6e0e5","_cell_guid":"b1076dfc-b9ad-4769-8c92-a6c4dae69d19","execution":{"iopub.status.busy":"2023-12-17T12:58:39.353410Z","iopub.execute_input":"2023-12-17T12:58:39.353774Z","iopub.status.idle":"2023-12-17T12:58:39.372188Z","shell.execute_reply.started":"2023-12-17T12:58:39.353746Z","shell.execute_reply":"2023-12-17T12:58:39.371347Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"volume = load_volume('/kaggle/input/blood-vessel-segmentation/train/kidney_1_dense/images', False)\nvolume = downsample_volume(volume)\nvolume_mask = load_volume('/kaggle/input/blood-vessel-segmentation/train/kidney_1_dense/labels', True)\nvolume_mask = downsample_volume(volume_mask)","metadata":{"execution":{"iopub.status.busy":"2023-12-17T12:58:39.374360Z","iopub.execute_input":"2023-12-17T12:58:39.374806Z","iopub.status.idle":"2023-12-17T13:01:45.279762Z","shell.execute_reply.started":"2023-12-17T12:58:39.374754Z","shell.execute_reply":"2023-12-17T13:01:45.278315Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"from scipy.spatial.transform import Rotation as R\n\n\ndef make_3d_slice(volume: np.ndarray,\n                  volume_mask: np.ndarray,\n                  out_shape: tp.Tuple[int, int], \n                  euler_angles: tp.Tuple[float, float, float],\n                  scale: float,\n                  maxval: float = 255) -> tp.Tuple[np.ndarray, np.ndarray]:\n    assert volume.shape == volume_mask.shape\n    \n    offset = [0, 0, 0]\n    \n    volume_center = np.array([*volume.shape], dtype=np.float32) * 0.5\n    grid_center = np.array([out_shape[1], out_shape[0], 0], dtype=np.float32) * 0.5\n\n    x, y = np.meshgrid(np.linspace(0, out_shape[1], out_shape[1]),\n                       np.linspace(0, out_shape[0], out_shape[0]),\n                       indexing='ij', sparse=False)\n    z = np.zeros((out_shape[1], out_shape[0]), dtype=np.float32)\n    z[:, :] = grid_center[2]\n    grid = np.stack([x, y, z], axis=-1).reshape(-1, 3)\n\n    r = R.from_euler('xyz', euler_angles, degrees=True)\n    transform_matrix = r.as_matrix() * scale\n\n    proj_grid = (grid - grid_center) @ transform_matrix.T\n    proj_grid[:, 0] += offset[0] + volume_center[0]\n    proj_grid[:, 1] += offset[1] + volume_center[1]\n    proj_grid[:, 2] += offset[2] + volume_center[2]\n    proj_grid[proj_grid[:, 0] < 0, 0] = 0\n    proj_grid[proj_grid[:, 0] > (volume.shape[0] - 1), 0] = volume.shape[0] - 1\n    proj_grid[proj_grid[:, 1] < 0, 1] = 0\n    proj_grid[proj_grid[:, 1] > (volume.shape[1] - 1), 1] = volume.shape[1] - 1\n    proj_grid[proj_grid[:, 2] < 0, 2] = 0\n    proj_grid[proj_grid[:, 2] > (volume.shape[2] - 1), 2] = volume.shape[2] - 1\n    \n    proj_grid_i: np.ndarray = np.floor(proj_grid)\n    w_proj_grid = proj_grid - proj_grid_i\n    proj_grid_i = proj_grid_i.astype(np.int32)\n    iw_proj_grid = 1 - w_proj_grid\n    \n    vertex_coords: tp.List[np.ndarray] = [proj_grid_i]\n    vertex_weights: tp.List[np.ndarray] = [iw_proj_grid[:, 0] * iw_proj_grid[:, 1] * iw_proj_grid[:, 2]]\n    \n    proj_grid_i_ = proj_grid_i.copy()\n    proj_grid_i_[:, 0] = np.minimum(proj_grid_i_[:, 0] + 1, volume.shape[0] - 1)\n    vertex_coords.append(proj_grid_i_)\n    vertex_weights.append(w_proj_grid[:, 0] * iw_proj_grid[:, 1] * iw_proj_grid[:, 2])\n    \n    proj_grid_i_ = proj_grid_i.copy()\n    proj_grid_i_[:, 0] = np.minimum(proj_grid_i_[:, 0] + 1, volume.shape[0] - 1)\n    proj_grid_i_[:, 1] = np.minimum(proj_grid_i_[:, 1] + 1, volume.shape[1] - 1)\n    vertex_coords.append(proj_grid_i_)\n    vertex_weights.append(w_proj_grid[:, 0] * w_proj_grid[:, 1] * iw_proj_grid[:, 2])\n    \n    proj_grid_i_ = proj_grid_i.copy()\n    proj_grid_i_[:, 0] = np.minimum(proj_grid_i_[:, 0] + 1, volume.shape[0] - 1)\n    proj_grid_i_[:, 2] = np.minimum(proj_grid_i_[:, 2] + 1, volume.shape[2] - 1)\n    vertex_coords.append(proj_grid_i_)\n    vertex_weights.append(w_proj_grid[:, 0] * iw_proj_grid[:, 1] * w_proj_grid[:, 2])\n    \n    proj_grid_i_ = proj_grid_i.copy()\n    proj_grid_i_[:, 0] = np.minimum(proj_grid_i_[:, 0] + 1, volume.shape[0] - 1)\n    proj_grid_i_[:, 1] = np.minimum(proj_grid_i_[:, 1] + 1, volume.shape[1] - 1)\n    proj_grid_i_[:, 2] = np.minimum(proj_grid_i_[:, 2] + 1, volume.shape[2] - 1)\n    vertex_coords.append(proj_grid_i_)\n    vertex_weights.append(w_proj_grid[:, 0] * w_proj_grid[:, 1] * w_proj_grid[:, 2])\n    \n    proj_grid_i_ = proj_grid_i.copy()\n    proj_grid_i_[:, 1] = np.minimum(proj_grid_i_[:, 1] + 1, volume.shape[1] - 1)\n    proj_grid_i_[:, 2] = np.minimum(proj_grid_i_[:, 2] + 1, volume.shape[2] - 1)\n    vertex_coords.append(proj_grid_i_)\n    vertex_weights.append(iw_proj_grid[:, 0] * w_proj_grid[:, 1] * w_proj_grid[:, 2])\n    \n    proj_grid_i_ = proj_grid_i.copy()\n    proj_grid_i_[:, 1] = np.minimum(proj_grid_i_[:, 1] + 1, volume.shape[1] - 1)\n    vertex_coords.append(proj_grid_i_)\n    vertex_weights.append(iw_proj_grid[:, 0] * w_proj_grid[:, 1] * iw_proj_grid[:, 2])\n    \n    proj_grid_i_ = proj_grid_i.copy()\n    proj_grid_i_[:, 2] = np.minimum(proj_grid_i_[:, 2] + 1, volume.shape[2] - 1)\n    vertex_coords.append(proj_grid_i_)\n    vertex_weights.append(iw_proj_grid[:, 0] * iw_proj_grid[:, 1] * w_proj_grid[:, 2])\n    \n    slice_3d = np.zeros((out_shape[1], out_shape[0]), np.float32)\n    mask_slice_3d = np.zeros((out_shape[1], out_shape[0]), np.float32)\n    flat_slice_3d = slice_3d.reshape(-1)\n    mask_flat_slice_3d = mask_slice_3d.reshape(-1)\n    for v_coords, v_weights in zip(vertex_coords, vertex_weights):\n        flat_slice_3d += volume[v_coords[:, 0], v_coords[:, 1], v_coords[:, 2]] * v_weights\n        mask_flat_slice_3d += volume_mask[v_coords[:, 0], v_coords[:, 1], v_coords[:, 2]] * v_weights\n    slice_3d[slice_3d > maxval] = maxval\n    slice_3d = slice_3d.astype(volume.dtype)\n    mask_slice_3d[mask_slice_3d > 255] = 255\n    mask_slice_3d = mask_slice_3d.astype(np.uint8)\n    \n    return slice_3d, mask_slice_3d","metadata":{"execution":{"iopub.status.busy":"2023-12-17T13:01:45.282235Z","iopub.execute_input":"2023-12-17T13:01:45.282634Z","iopub.status.idle":"2023-12-17T13:01:45.479177Z","shell.execute_reply.started":"2023-12-17T13:01:45.282598Z","shell.execute_reply":"2023-12-17T13:01:45.478133Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"%matplotlib inline\nimport matplotlib.pyplot as plt\nimport matplotlib.animation as animation\nfrom IPython.display import HTML\n\n\ndef create_animation(images, masks, figsize = (8, 4), cmap_mask='cividis', alpha_mask=0.5):\n    fig, ax = plt.subplots(figsize=figsize)\n    artists = []\n    for img, mask in zip(images, masks):\n        im1 = ax.imshow(img, cmap='gray', animated=True)\n        im2 = ax.imshow(mask, cmap=cmap_mask, alpha=alpha_mask, animated=True)\n        plt.axis(\"off\")\n        plt.close()\n        artists.append([im1, im2])\n\n    ani = animation.ArtistAnimation(fig=fig,\n                                    artists=artists,\n                                    interval=100,\n                                    blit=True,\n                                    repeat=True)\n    return ani\n    \nslices_3d = []\nmasks = []\n\nfor angle in np.linspace(0, 360, 36):\n    slice_3d, mask = make_3d_slice(volume, volume_mask, (256, 256),\n                                   euler_angles=(0, angle, 0), scale=2.0)\n    slices_3d.append(slice_3d)\n    masks.append(mask)\n\nani = create_animation(slices_3d, masks)\nHTML(ani.to_jshtml())","metadata":{"execution":{"iopub.status.busy":"2023-12-17T13:01:45.480882Z","iopub.execute_input":"2023-12-17T13:01:45.481510Z","iopub.status.idle":"2023-12-17T13:01:54.711510Z","shell.execute_reply.started":"2023-12-17T13:01:45.481477Z","shell.execute_reply":"2023-12-17T13:01:54.710050Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"slices_3d = []\nmasks = []\n\nfor angle in np.linspace(0, 360, 36):\n    slice_3d, mask = make_3d_slice(volume, volume_mask, (256, 256),\n                                   euler_angles=(angle, 0, 0), scale=2.0)\n    slices_3d.append(slice_3d)\n    masks.append(mask)\n\nani = create_animation(slices_3d, masks)\nHTML(ani.to_jshtml())","metadata":{"execution":{"iopub.status.busy":"2023-12-17T13:01:54.715006Z","iopub.execute_input":"2023-12-17T13:01:54.715583Z","iopub.status.idle":"2023-12-17T13:02:04.099181Z","shell.execute_reply.started":"2023-12-17T13:01:54.715550Z","shell.execute_reply":"2023-12-17T13:02:04.097246Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"slices_3d = []\nmasks = []\n\nfor angle in np.linspace(0, 360, 36):\n    slice_3d, mask = make_3d_slice(volume, volume_mask, (256, 256),\n                                   euler_angles=(0, 0, angle), scale=2.0)\n    slices_3d.append(slice_3d)\n    masks.append(mask)\n\nani = create_animation(slices_3d, masks)\nHTML(ani.to_jshtml())","metadata":{"execution":{"iopub.status.busy":"2023-12-17T13:02:04.101235Z","iopub.execute_input":"2023-12-17T13:02:04.101794Z","iopub.status.idle":"2023-12-17T13:02:13.576416Z","shell.execute_reply.started":"2023-12-17T13:02:04.101741Z","shell.execute_reply":"2023-12-17T13:02:13.574763Z"},"trusted":true},"execution_count":null,"outputs":[]}]}