{ "cells": [ { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "import sys\n", "from pathlib import Path\n", "import numpy as np\n", "import matplotlib.pyplot as plt\n", "\n", "triskele_path = Path('../triskele/python/')\n", "sys.path.append(str(triskele_path.resolve()))\n", "import triskele" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "# Specific Utils\n", "\n", "def DFC_filter(raster):\n", " raster[raster > 1e4] = raster[raster < 1e4].max()\n", "\n", "def show(im, im_size=1, save=None):\n", " plt.figure(figsize=(16*im_size,3*im_size))\n", " plt.imshow(im)\n", " plt.colorbar()\n", " \n", " if save is not None:\n", " plt.savefig(save, bbox_inches='tight', pad_inches=1)\n", " \n", " plt.show()\n", "\n", "def mshow(Xs, titles=None, im_size=1, save=None):\n", " s = len(Xs)\n", "\n", " plt.figure(figsize=(16*im_size,3*im_size*s))\n", "\n", " for i in range(s):\n", " plt.subplot(s,1,i+1)\n", " plt.imshow(Xs[i])\n", " \n", " if titles is not None:\n", " plt.title(titles[i])\n", " \n", " plt.colorbar()\n", " \n", " if save is not None:\n", " plt.savefig(save, bbox_inches='tight', pad_inches=1)\n", " \n", " plt.show()" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "layers_files = [\n", " '../Data/phase1_rasters/DEM+B_C123/UH17_GEM051_TR.tif',\n", " '../Data/phase1_rasters/DEM_C123_3msr/UH17_GEG051_TR.tif',\n", " '../Data/phase1_rasters/DEM_C123_TLI/UH17_GEG05_TR.tif',\n", " '../Data/phase1_rasters/DSM_C12/UH17c_GEF051_TR.tif',\n", " '../Data/phase1_rasters/Intensity_C1/UH17_GI1F051_TR.tif',\n", " '../Data/phase1_rasters/Intensity_C2/UH17_GI2F051_TR.tif',\n", " '../Data/phase1_rasters/Intensity_C3/UH17_GI3F051_TR.tif',\n", " #'../Data/ground_truth/2018_IEEE_GRSS_DFC_GT_TR.tif'\n", "]" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Define dataset dependent raster filtering" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "def DFC_filter(raster):\n", " ## Remove extrem values\n", " #raster[raster == raster.max()] = raster[raster != raster.max()].max()\n", " raster[raster > 1e4] = raster[raster < 1e4].max()\n", " #raster[raster == np.finfo(raster.dtype).max] = raster[raster != raster.max()].max()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Load rasters data" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "layers = list()\n", "\n", "for file in layers_files:\n", " print('Loading {}'.format(file))\n", " layer = triskele.read(file)\n", " DFC_filter(layer)\n", " layers.append(layer)\n", "\n", "layers_stack = np.stack(layers, axis=2)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Display rasters" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "for i in range(layers_stack.shape[2]):\n", " plt.figure(figsize=(16*2,3*2))\n", " plt.imshow(layers_stack[:,:,i])\n", " plt.colorbar()\n", " plt.title(layers_files[i])\n", " plt.show()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Attributes filter with TRISKELE !" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "area = np.array([10, 100, 1e3, 1e4, 1e5])\n", "sd = np.array([0.5,0.9,0.99,0.999])#,1e4,1e5,5e5])\n", "moi = np.array([0.01,0.02,0.03,0.04,0.05,0.06,0.07,0.08,0.09,0.1,0.2,0.3,0.4,0.5,0.6,0.7,0.8,0.9,0.99])\n", "\n", "t = triskele.Triskele(layers_stack[:,:,:2], verbose=False)\n", "attributes_min = t.filter(tree='min-tree',\n", " area=area,\n", " standard_deviation=sd,\n", " moment_of_inertia=moi\n", " )\n", "attributes_max = t.filter(tree='max-tree',\n", " area=area,\n", " standard_deviation=sd,\n", " moment_of_inertia=moi\n", " )\n", "\n", "attributes_min_lbl = ['origin']\n", "attributes_min_lbl.extend(['Thickening area {}'.format(x) for x in area])\n", "attributes_min_lbl.extend(['Thickening σ {}'.format(x) for x in sd])\n", "attributes_min_lbl.extend(['Thickening moment of inertia {}'.format(x) for x in sd])\n", "\n", "attributes_max_lbl = ['origin']\n", "attributes_max_lbl.extend(['Thinning area {}'.format(x) for x in area])\n", "attributes_max_lbl.extend(['Thinning σ {}'.format(x) for x in sd])\n", "attributes_max_lbl.extend(['Thinning moment of inertia {}'.format(x) for x in sd])\n", "\n", "attributes_min.shape, attributes_max.shape" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "attributes = np.dstack((attributes_min[:,:,:0:-1], attributes_max))\n", "attributes_lbl = attributes_min_lbl[:0:-1]\n", "attributes_lbl.extend(attributes_max_lbl)\n", "\n", "attributes.shape, attributes_lbl" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "figs = list()\n", "\n", "attributes[0,0,:] = 255 # J'ai honte...\n", "\n", "for i in range(attributes.shape[-1]):\n", " figs.append(attributes[:,:,i])\n", "\n", "mshow(figs, attributes_lbl, 2)" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "f = plt.Figure(figsize=(16,9))\n", "f.imshow(attributes[:,:,0])\n", "f.plt.show()" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "show(attributes[:,:,17])" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "{'I_{{{}}}'.format(42)}" ] } ], "metadata": { "kernelspec": { "display_name": "Python 3", "language": "python", "name": "python3" }, "language_info": { "codemirror_mode": { "name": "ipython", "version": 3 }, "file_extension": ".py", "mimetype": "text/x-python", "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3" } }, "nbformat": 4, "nbformat_minor": 2 }