{"metadata":{"kernelspec":{"language":"python","display_name":"Python 3","name":"python3"},"language_info":{"pygments_lexer":"ipython3","nbconvert_exporter":"python","version":"3.6.4","file_extension":".py","codemirror_mode":{"name":"ipython","version":3},"name":"python","mimetype":"text/x-python"}},"nbformat_minor":4,"nbformat":4,"cells":[{"cell_type":"code","source":"import numpy as np\nimport pandas as pd\nimport matplotlib.pyplot as plt\nimport tifffile as tiff\nimport cv2\nimport os\nfrom tqdm.notebook import tqdm\nimport zipfile\nimport rasterio\nfrom rasterio.windows import Window","metadata":{"_uuid":"8f2839f25d086af736a60e9eeb907d3b93b6e0e5","_cell_guid":"b1076dfc-b9ad-4769-8c92-a6c4dae69d19","trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"sz = 480   #the size of tiles\nWINDOW = 1024\nOVERLAP = 204\nMASKS = '../input/hubmap-kidney-segmentation/train.csv'\nDATA = '../input/hubmap-kidney-segmentation/train/'\nOUT_TRAIN = 'train.zip'\nOUT_MASKS = 'masks.zip'","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#functions to convert encoding to mask and mask to encoding\ndef enc2mask(encs, shape):\n    img = np.zeros(shape[0]*shape[1], dtype=np.uint8)\n    for m,enc in enumerate(encs):\n        if isinstance(enc,np.float) and np.isnan(enc): continue\n        s = enc.split()\n        for i in range(len(s)//2):\n            start = int(s[2*i]) - 1\n            length = int(s[2*i+1])\n            img[start:start+length] = 1 + m\n    return img.reshape(shape).T\n\ndef mask2enc(mask, n=1):\n    pixels = mask.T.flatten()\n    encs = []\n    for i in range(1,n+1):\n        p = (pixels == i).astype(np.int8)\n        if p.sum() == 0: encs.append(np.nan)\n        else:\n            p = np.concatenate([[0], p, [0]])\n            runs = np.where(p[1:] != p[:-1])[0] + 1\n            runs[1::2] -= runs[::2]\n            encs.append(' '.join(str(x) for x in runs))\n    return encs\n\ndf_masks = pd.read_csv(MASKS).set_index('id')\ndf_masks.head()","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def make_grid(shape, window=256, min_overlap=32):\n    \"\"\"\n        Return Array of size (N,4), where N - number of tiles,\n        2nd axis represente slices: x1,x2,y1,y2 \n    \"\"\"\n    x, y = shape\n    nx = x // (window - min_overlap) + 1\n    x1 = np.linspace(0, x, num=nx, endpoint=False, dtype=np.int64)\n    x1[-1] = x - window\n    x2 = (x1 + window).clip(0, x)\n    assert np.all(x2-x1 == window), \"Row or height not equal to window. All tiles must be window_x_window\"\n    ny = y // (window - min_overlap) + 1\n    y1 = np.linspace(0, y, num=ny, endpoint=False, dtype=np.int64)\n    y1[-1] = y - window\n    y2 = (y1 + window).clip(0, y)\n    assert np.all(y2-y1 == window), \"Column or width not equal to window. All tiles must be window_x_window\"\n    slices = np.zeros((nx,ny, 4), dtype=np.int64)\n    \n    for i in range(nx):\n        for j in range(ny):\n            slices[i,j] = x1[i], x2[i], y1[j], y2[j]    \n    return slices.reshape(nx*ny,4)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"std_th = 7\nx_means,x2_means = [],[]\n\nwith zipfile.ZipFile(OUT_TRAIN, 'w') as img_out,\\\n zipfile.ZipFile(OUT_MASKS, 'w') as mask_out:\n    for index, encs in tqdm(df_masks.iterrows(),total=len(df_masks)):\n        img_dataset = rasterio.open(os.path.join(DATA,index+'.tiff'))\n        if img_dataset.count != 3:\n            img_layers = []\n            if len(img_dataset.subdatasets) > 0:\n                for i, subdataset in enumerate(img_dataset.subdatasets, 0):\n                    img_layers.append(rasterio.open(subdataset))\n        slice_coords = make_grid(img_dataset.shape, window=WINDOW, min_overlap=OVERLAP) # np.array; num_slices,4; np.int64 \n        mask = enc2mask(encs,(img_dataset.shape[1],img_dataset.shape[0]))  # np.array; H,W; np.uint8\n        \n        for i, coords in enumerate(slice_coords):\n            x1,x2,y1,y2 = coords\n            if img_dataset.count == 3:\n                img_tile = img_dataset.read(window=Window.from_slices((x1,x2),(y1,y2))) # shape: [C,H,W]\n            else:\n                img_tile = np.zeros((3,WINDOW,WINDOW),np.uint8)\n                for j, layer in enumerate(img_layers):\n                    img_tile[j,:,:] = layer.read(window=Window.from_slices((x1,x2),(y1,y2)))[0]\n            mask_tile = mask[x1:x2,y1:y2]\n            img_tile = img_tile.transpose((1,2,0)) # shape: [H,W,C], np.uint8\n            \n            #remove image slices with no real content\n            if np.all(np.array([np.std(img_tile[:,:,i]) for i in range(3)]) <= std_th):\n                continue\n            \n            img_tile = cv2.resize(img_tile, (sz,sz), interpolation=cv2.INTER_AREA)\n            mask_tile = cv2.resize(mask_tile, (sz,sz), interpolation=cv2.INTER_NEAREST)\n            \n            x_means.append((img_tile/255.0).reshape(-1,3).mean(0))\n            x2_means.append(((img_tile/255.0)**2).reshape(-1,3).mean(0))\n            \n            img_tile = cv2.imencode('.png',cv2.cvtColor(img_tile, cv2.COLOR_RGB2BGR))[1]\n            img_out.writestr(f'{index}_{i}.png', img_tile)\n            mask_tile = cv2.imencode('.png',mask_tile)[1]\n            mask_out.writestr(f'{index}_{i}.png', mask_tile)\n                \n#image stats\nimg_mean =  np.stack(x_means).mean(0)\nimg_std =  np.sqrt(np.stack(x2_means).mean(0) - img_mean**2)\nprint(f'mean:{img_mean}, std:{img_std}')","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"columns, rows = 4,4\nidx0 = 763\nfig=plt.figure(figsize=(columns*4, rows*4))\nwith zipfile.ZipFile(OUT_TRAIN, 'r') as img_arch, \\\n     zipfile.ZipFile(OUT_MASKS, 'r') as msk_arch:\n    fnames = img_arch.namelist()\n    for i in range(rows):\n        for j in range(columns):\n            idx = i*columns+j\n            img = cv2.imdecode(np.frombuffer(img_arch.read(fnames[idx0+idx]), \n                                             np.uint8), cv2.IMREAD_COLOR)\n            img = cv2.cvtColor(img, cv2.COLOR_BGR2RGB)\n            mask = cv2.imdecode(np.frombuffer(msk_arch.read(fnames[idx0+idx]), \n                                              np.uint8), cv2.IMREAD_GRAYSCALE)\n    \n            fig.add_subplot(rows, columns, idx+1)\n            plt.axis('off')\n            plt.imshow(img)\n            #plt.imshow(Image.fromarray(mask), alpha=0.2)\nplt.show()","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"columns, rows = 4,4\nidx0 = 763\nfig=plt.figure(figsize=(columns*4, rows*4))\nwith zipfile.ZipFile(OUT_TRAIN, 'r') as img_arch, \\\n     zipfile.ZipFile(OUT_MASKS, 'r') as msk_arch:\n    fnames = img_arch.namelist()\n    for i in range(rows):\n        for j in range(columns):\n            idx = i*columns+j\n            img = cv2.imdecode(np.frombuffer(img_arch.read(fnames[idx0+idx]), \n                                             np.uint8), cv2.IMREAD_COLOR)\n            img = cv2.cvtColor(img, cv2.COLOR_BGR2RGB)\n            mask = cv2.imdecode(np.frombuffer(msk_arch.read(fnames[idx0+idx]), \n                                              np.uint8), cv2.IMREAD_GRAYSCALE)\n    \n            fig.add_subplot(rows, columns, idx+1)\n            plt.axis('off')\n            plt.imshow(img)\n            plt.imshow(mask, alpha=0.2)\nplt.show()","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"! ls -alrh ","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"","metadata":{"trusted":true},"execution_count":null,"outputs":[]}]}