From e90108a2c83f57a6163d4a5091eada63111ec943 Mon Sep 17 00:00:00 2001 From: Karamaz0V1 Date: Fri, 31 Aug 2018 10:54:26 +0200 Subject: [PATCH 1/2] Update factory with latest tresholds --- Notebooks/Raster Factory.ipynb | 44 +++++++++++++--------------------- 1 file changed, 16 insertions(+), 28 deletions(-) diff --git a/Notebooks/Raster Factory.ipynb b/Notebooks/Raster Factory.ipynb index 25ddbc4..9587aec 100644 --- a/Notebooks/Raster Factory.ipynb +++ b/Notebooks/Raster Factory.ipynb @@ -41,16 +41,10 @@ "outputs": [], "source": [ "C1 = ra.bulk_load('../Data/lidar/C1', 'C1', filter_treshold=.01, clip_treshold=.1, dtype=np.float32)\n", - "#C2 = ra.bulk_load('../Data/lidar/C2', 'C2', filter_treshold=.5, dtype=np.float32)\n", - "#C3 = ra.bulk_load('../Data/lidar/C3', 'C3', filter_treshold=.5, dtype=np.float32)\n", - "#C123 = ra.bulk_load('../Data/lidar', 'C123', filter_treshold=.5, dtype=np.float32)" - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "# TMP" + "C2 = ra.bulk_load('../Data/lidar/C2', 'C2', filter_treshold=.01, clip_treshold=.1, dtype=np.float32)\n", + "C3 = ra.bulk_load('../Data/lidar/C3', 'C3', filter_treshold=.2, clip_treshold=.1, dtype=np.float32)\n", + "C12 = ra.bulk_load('../Data/lidar/C1', 'C1', filter_treshold=.01, clip_treshold=.1, dtype=np.float32)\n", + "C123 = ra.bulk_load('../Data/lidar', 'C123', filter_treshold=.08, clip_treshold=.1, dtype=np.float32)" ] }, { @@ -148,7 +142,7 @@ "metadata": {}, "outputs": [], "source": [ - "ra.auto_filter(dfc, treshold=0.5)" + "ra.extremum_clip(dfc, treshold=0.1)" ] }, { @@ -215,16 +209,7 @@ "cell_type": "markdown", "metadata": {}, "source": [ - "## Raster Pack" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "C123.name" + "## Raster Pool Process" ] }, { @@ -238,25 +223,28 @@ "#mp.set_start_method('spawn')\n", "\n", "def rasterize_cache_mp(data_var, field, res, method, reverse, cache):\n", + " print('FROM WORKER: Start processing {} {} with {} resolution and {} interpolation {}'.format(data_var, field, res, method, 'resversed' if reverse else ''))\n", " if data_var == 'C1':\n", - " ra.rasterize_cache(C1, field, res, method, reverse, cache)\n", + " ra.rasterize_cache(field, C1, res, method, reverse, cache) \n", " if data_var == 'C2':\n", - " ra.rasterize_cache(C2, field, res, method, reverse, cache)\n", + " ra.rasterize_cache(field, C2, res, method, reverse, cache)\n", " if data_var == 'C3':\n", - " ra.rasterize_cache(C3, field, res, method, reverse, cache)\n", + " ra.rasterize_cache(field, C3, res, method, reverse, cache)\n", " if data_var == 'C123':\n", - " ra.rasterize_cache(C123, field, res, method, reverse, cache)\n", + " ra.rasterize_cache(field, C123, res, method, reverse, cache)\n", " \n", - "pool = Pool(processes=5)\n", + "pool = Pool(processes=9)\n", "\n", "job_args = list()\n", "\n", "for res in (0.5, 1., 2., 3., 5., 10., .1):\n", " for reverse in (False, True):\n", - " for inter in ('linear', 'nearest'):\n", + " for inter in ('linear', 'nearest', 'cubic-clip'):\n", " for field in ('z', 'intensity', 'num_returns'):\n", " for data in ('C1', 'C2', 'C3', 'C123'):\n", - " job_args.append([data, field, res, inter, reverse, '../Res/HVR/'])\n", + " job_args.append([data, field, res, inter, reverse, '../Res/enrichment_rasters'])\n", + " \n", + "print(\"Job list length: {}\".format(len(job_args)))\n", " \n", "for i in pool.starmap(rasterize_cache_mp, job_args):\n", " pass" From d7f1163c06615e56a24f054bcbd2368b6bca7195 Mon Sep 17 00:00:00 2001 From: Karamaz0V1 Date: Fri, 31 Aug 2018 15:55:52 +0200 Subject: [PATCH 2/2] Try some raster combination --- Notebooks/Raster Factory.ipynb | 6 +- Notebooks/Raster combination.ipynb | 246 +++++++++++++++++++++++++++++ 2 files changed, 250 insertions(+), 2 deletions(-) create mode 100644 Notebooks/Raster combination.ipynb diff --git a/Notebooks/Raster Factory.ipynb b/Notebooks/Raster Factory.ipynb index 9587aec..46b8bca 100644 --- a/Notebooks/Raster Factory.ipynb +++ b/Notebooks/Raster Factory.ipynb @@ -43,7 +43,7 @@ "C1 = ra.bulk_load('../Data/lidar/C1', 'C1', filter_treshold=.01, clip_treshold=.1, dtype=np.float32)\n", "C2 = ra.bulk_load('../Data/lidar/C2', 'C2', filter_treshold=.01, clip_treshold=.1, dtype=np.float32)\n", "C3 = ra.bulk_load('../Data/lidar/C3', 'C3', filter_treshold=.2, clip_treshold=.1, dtype=np.float32)\n", - "C12 = ra.bulk_load('../Data/lidar/C1', 'C1', filter_treshold=.01, clip_treshold=.1, dtype=np.float32)\n", + "C12 = ra.bulk_load(['../Data/lidar/C1', '../Data/lidar/C2'], 'C12', filter_treshold=.01, clip_treshold=.1, dtype=np.float32)\n", "C123 = ra.bulk_load('../Data/lidar', 'C123', filter_treshold=.08, clip_treshold=.1, dtype=np.float32)" ] }, @@ -230,6 +230,8 @@ " ra.rasterize_cache(field, C2, res, method, reverse, cache)\n", " if data_var == 'C3':\n", " ra.rasterize_cache(field, C3, res, method, reverse, cache)\n", + " if data_var == 'C12':\n", + " ra.rasterize_cache(field, C12, res, method, reverse, cache)\n", " if data_var == 'C123':\n", " ra.rasterize_cache(field, C123, res, method, reverse, cache)\n", " \n", @@ -241,7 +243,7 @@ " for reverse in (False, True):\n", " for inter in ('linear', 'nearest', 'cubic-clip'):\n", " for field in ('z', 'intensity', 'num_returns'):\n", - " for data in ('C1', 'C2', 'C3', 'C123'):\n", + " for data in ('C1', 'C2', 'C3', 'C12', 'C123'):\n", " job_args.append([data, field, res, inter, reverse, '../Res/enrichment_rasters'])\n", " \n", "print(\"Job list length: {}\".format(len(job_args)))\n", diff --git a/Notebooks/Raster combination.ipynb b/Notebooks/Raster combination.ipynb new file mode 100644 index 0000000..fdc081b --- /dev/null +++ b/Notebooks/Raster combination.ipynb @@ -0,0 +1,246 @@ +{ + "cells": [ + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "import sys\n", + "import numpy as np\n", + "import matplotlib.pyplot as plt\n", + "\n", + "sys.path.append(\"..\")\n", + "import rasterizer\n", + "import raster_assistant as ra\n", + "\n", + "sys.path.append('../triskele/python/')\n", + "import triskele\n", + "\n", + "figsize = np.array((16, 3)) * 1.5" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "# Raster Combination\n", + "\n", + "## DSM - DTM " + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "C12 = ra.bulk_load(['../Data/lidar/C1', '../Data/lidar/C2'], 'C12', filter_treshold=.01, clip_treshold=.1, dtype=np.float32)" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "dsm = ra.rasterize_cache('z', C12, .5, 'nearest', False, '../Res/enrichment_rasters/')" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "dtm = ra.rasterize_cache('z', C12, .5, 'nearest', True, '../Res/enrichment_rasters/')" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "plt.figure(figsize=figsize)\n", + "plt.imshow(dsm)\n", + "plt.colorbar()\n", + "plt.title('DSM')\n", + "plt.show()\n", + "\n", + "plt.figure(figsize=figsize)\n", + "plt.imshow(dtm)\n", + "plt.colorbar()\n", + "plt.title('DTM')\n", + "plt.show()\n", + "\n", + "plt.figure(figsize=figsize)\n", + "plt.imshow(dsm - dtm)\n", + "plt.colorbar()\n", + "plt.title('DSM - DTM')\n", + "plt.show()\n", + "\n", + "plt.imsave('../Res/dsm-dtm.png', dsm - dtm)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## NDVI\n", + "\n", + "With \n", + "\n", + "$NDVI = \\frac{NIR - Red}{NIR + Red}$\n", + "\n", + "and the wavelenth of the Titan\n", + "\n", + "| Lazer | Wavelenght | Color |\n", + "| --- | ---: | --- |\n", + "| C1 | 1550 nm | IR? |\n", + "| C2 | 1064 nm | NIR |\n", + "| C3 | 532 nm | Green |\n", + "\n", + "we can compute a NDVI like intensity raster with\n", + "\n", + "$NDVI_{like} = \\frac{C1 - C2}{C1 + C2}$" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "C1 = ra.bulk_load('../Data/lidar/C1', 'C1', filter_treshold=.01, clip_treshold=.1, dtype=np.float32)\n", + "C2 = ra.bulk_load('../Data/lidar/C2', 'C2', filter_treshold=.01, clip_treshold=.1, dtype=np.float32)" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "C1_raster = ra.rasterize_cache('intensity', C1, .5, 'linear', False, '../Res/enrichment_rasters/')\n", + "C2_raster = ra.rasterize_cache('intensity', C2, .5, 'linear', False, '../Res/enrichment_rasters/')" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "ndvi = (C1_raster - C2_raster) / (C1_raster + C2_raster)" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "plt.figure(figsize=figsize)\n", + "plt.imshow(C1_raster)\n", + "plt.title('C1')\n", + "plt.colorbar()\n", + "plt.show()\n", + "\n", + "plt.figure(figsize=figsize)\n", + "plt.imshow(C2_raster)\n", + "plt.title('C2')\n", + "plt.colorbar()\n", + "plt.show()\n", + "\n", + "plt.figure(figsize=figsize)\n", + "plt.imshow(ndvi)\n", + "plt.title('NDVI')\n", + "plt.colorbar()\n", + "plt.show()\n", + "\n", + "plt.imsave('../Res/ndvi_linear.png', ndvi)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Normalized NDVI" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "C1_raster /= C1_raster.max()\n", + "C2_raster /= C2_raster.max()" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "ndvi = (C1_raster - C2_raster) / (C1_raster + C2_raster)" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "plt.figure(figsize=figsize)\n", + "plt.imshow(C1_raster)\n", + "plt.title('C1')\n", + "plt.colorbar()\n", + "plt.show()\n", + "\n", + "plt.figure(figsize=figsize)\n", + "plt.imshow(C2_raster)\n", + "plt.title('C2')\n", + "plt.colorbar()\n", + "plt.show()\n", + "\n", + "plt.figure(figsize=figsize)\n", + "plt.imshow(ndvi)\n", + "plt.title('NDVI')\n", + "plt.colorbar()\n", + "plt.show()\n", + "\n", + "plt.imsave('../Res/ndvi_normalized_linear.png', ndvi)" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [] + } + ], + "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 +}