289 lines
6.8 KiB
Plaintext
289 lines
6.8 KiB
Plaintext
{
|
|
"cells": [
|
|
{
|
|
"cell_type": "markdown",
|
|
"metadata": {},
|
|
"source": [
|
|
"# Custom Raster Gen for LD2DAPs"
|
|
]
|
|
},
|
|
{
|
|
"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": [
|
|
"## Load LAS data"
|
|
]
|
|
},
|
|
{
|
|
"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)\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', '../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)"
|
|
]
|
|
},
|
|
{
|
|
"cell_type": "code",
|
|
"execution_count": null,
|
|
"metadata": {},
|
|
"outputs": [],
|
|
"source": [
|
|
"raster = ra.rasterize_cache('num_returns', C1, .5, 'nearest', False, cache_dir='../Res/enrichment_rasters')\n",
|
|
"\n",
|
|
"plt.figure(figsize=figsize)\n",
|
|
"plt.imshow(raster, origin='upper')\n",
|
|
"plt.show()\n",
|
|
"plt.imsave('../Res/raster_validation.png', raster)"
|
|
]
|
|
},
|
|
{
|
|
"cell_type": "markdown",
|
|
"metadata": {},
|
|
"source": [
|
|
"## Raster Validation\n",
|
|
"\n",
|
|
"### Rasterize some data"
|
|
]
|
|
},
|
|
{
|
|
"cell_type": "code",
|
|
"execution_count": null,
|
|
"metadata": {},
|
|
"outputs": [],
|
|
"source": [
|
|
"raster = rasterizer.rasterize(C1.spatial, C1.intensity, 0.5, dtype=np.float32)"
|
|
]
|
|
},
|
|
{
|
|
"cell_type": "code",
|
|
"execution_count": null,
|
|
"metadata": {},
|
|
"outputs": [],
|
|
"source": [
|
|
"plt.figure(figsize=figsize)\n",
|
|
"plt.imshow(raster, origin='upper')\n",
|
|
"plt.show()\n",
|
|
"plt.imsave('../Res/raster_validation.png', raster)"
|
|
]
|
|
},
|
|
{
|
|
"cell_type": "code",
|
|
"execution_count": null,
|
|
"metadata": {},
|
|
"outputs": [],
|
|
"source": [
|
|
"plt.hist(raster.flatten(), bins=1000)\n",
|
|
"plt.show()"
|
|
]
|
|
},
|
|
{
|
|
"cell_type": "markdown",
|
|
"metadata": {},
|
|
"source": [
|
|
"### Write TIFF file"
|
|
]
|
|
},
|
|
{
|
|
"cell_type": "code",
|
|
"execution_count": null,
|
|
"metadata": {},
|
|
"outputs": [],
|
|
"source": [
|
|
"triskele.write('../Res/validation.tiff', raster)"
|
|
]
|
|
},
|
|
{
|
|
"cell_type": "markdown",
|
|
"metadata": {},
|
|
"source": [
|
|
"### Compare with DFC dataset"
|
|
]
|
|
},
|
|
{
|
|
"cell_type": "code",
|
|
"execution_count": null,
|
|
"metadata": {},
|
|
"outputs": [],
|
|
"source": [
|
|
"dfc = triskele.read('../Data/phase1_rasters/Intensity_C1/UH17_GI1F051_TR.tif')\n",
|
|
"our = triskele.read('../Res/validation.tiff')"
|
|
]
|
|
},
|
|
{
|
|
"cell_type": "markdown",
|
|
"metadata": {},
|
|
"source": [
|
|
"#### Filter DFC with same parameters"
|
|
]
|
|
},
|
|
{
|
|
"cell_type": "code",
|
|
"execution_count": null,
|
|
"metadata": {},
|
|
"outputs": [],
|
|
"source": [
|
|
"ra.extremum_clip(dfc, treshold=0.1)"
|
|
]
|
|
},
|
|
{
|
|
"cell_type": "markdown",
|
|
"metadata": {},
|
|
"source": [
|
|
"#### Display Stats"
|
|
]
|
|
},
|
|
{
|
|
"cell_type": "code",
|
|
"execution_count": null,
|
|
"metadata": {},
|
|
"outputs": [],
|
|
"source": [
|
|
"dfc.shape, our.shape"
|
|
]
|
|
},
|
|
{
|
|
"cell_type": "code",
|
|
"execution_count": null,
|
|
"metadata": {},
|
|
"outputs": [],
|
|
"source": [
|
|
"dfc.dtype, our.dtype"
|
|
]
|
|
},
|
|
{
|
|
"cell_type": "code",
|
|
"execution_count": null,
|
|
"metadata": {},
|
|
"outputs": [],
|
|
"source": [
|
|
"plt.hist(dfc.flatten(), bins=1000, label='DFC')\n",
|
|
"plt.hist(our.flatten(), bins=1000, label='Our', alpha=.8)\n",
|
|
"plt.legend()\n",
|
|
"plt.show()"
|
|
]
|
|
},
|
|
{
|
|
"cell_type": "markdown",
|
|
"metadata": {},
|
|
"source": [
|
|
"#### Display Rasters"
|
|
]
|
|
},
|
|
{
|
|
"cell_type": "code",
|
|
"execution_count": null,
|
|
"metadata": {},
|
|
"outputs": [],
|
|
"source": [
|
|
"f, axs = plt.subplots(2, figsize=figsize * 2)\n",
|
|
"\n",
|
|
"axs[0].imshow(dfc)\n",
|
|
"axs[0].set_title('DFC')\n",
|
|
"axs[1].imshow(our)\n",
|
|
"axs[1].set_title('Our')\n",
|
|
"\n",
|
|
"plt.show()"
|
|
]
|
|
},
|
|
{
|
|
"cell_type": "markdown",
|
|
"metadata": {},
|
|
"source": [
|
|
"## Raster Pool Process"
|
|
]
|
|
},
|
|
{
|
|
"cell_type": "code",
|
|
"execution_count": null,
|
|
"metadata": {},
|
|
"outputs": [],
|
|
"source": [
|
|
"from multiprocessing import Pool, Process, Queue\n",
|
|
"import multiprocessing as mp\n",
|
|
"#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(field, C1, res, method, reverse, cache) \n",
|
|
" if data_var == 'C2':\n",
|
|
" 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",
|
|
"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', 'cubic-clip'):\n",
|
|
" for field in ('z', 'intensity', 'num_returns'):\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",
|
|
" \n",
|
|
"for i in pool.starmap(rasterize_cache_mp, job_args):\n",
|
|
" pass"
|
|
]
|
|
},
|
|
{
|
|
"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
|
|
}
|