Skip to content

GitLab

O
OTB Guided Tour
Commit a425da5a authored by David Youssefi's avatar David Youssefi
Browse files
Options

ENH: add step 3

parent 0edd69ff
Hide whitespace changes
Inline Side-by-side
Showing with 202 additions and 47 deletions
otb_workshop.ipynb
View file @ a425da5a
......@@ -137,7 +137,7 @@
"cell_type": "markdown",
"metadata": {},
"source": [
"### Go to the [Step 1](http://localhost:8888/notebooks/otb-guided-tour/step1_compute_radiometric_indices.ipynb)"
"### Go to the [Step 1](./step1_compute_radiometric_indices.ipynb)"
]
}
],
......
step3_segment_watermask_threshold_NDWI.ipynb 0 → 100644
View file @ a425da5a
{
"cells": [
{
"cell_type": "markdown",
"metadata": {},
"source": [
"<img style=\"float: left; margin:0px 15px 15px 0px; width:120px\" src=\"https://www.orfeo-toolbox.org/wp-content/uploads/2016/03/logo-orfeo-toolbox.png\">\n",
"\n",
"# OTB Guided Tour - JURSE 2019 Vannes France - May 21th\n",
"## Yannick TANGUY and David YOUSSEFI (CNES, French Space Agency)\n",
"\n",
"<br>\n",
"\n",
"<b> Press <span style=\"color:black;background:yellow\">SHIFT+ENTER</span> to execute the notebook interactively cell by cell </b></div>\n",
"\n",
"## Step 3 : Segment watermask thresholding NDWI"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### Functions of the full OTB Pipeline"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"import otbApplication\n",
"\n",
"def app_radiometricindices(images):\n",
" index = 0\n",
" for im in images:\n",
" # Here we create an application\n",
" app = otbApplication.Registry.CreateApplication(\"BandMath\")\n",
" # BandMath takes a list of images as input,\n",
" # so we have to give a Python list [ image ], or [ image1, image2, .., imageN]\n",
" app.SetParameterStringList(\"il\",[im])\n",
" app.SetParameterString(\"out\",im.replace(\"RVBPIR.tif\",\"_indices.tif\"))\n",
" # the main parameter is the mathematical expression\n",
" # here, we compute NDVI : (nir - red) / (nir + red)\n",
" # The corresponding bands for NIR and Red are respectively the 4th and the 3rd bands (b4, b3)\n",
" # of the first image (im1)\n",
" app.SetParameterString(\"exp\",\"(im1b2+im1b4) > 0 ? (im1b2-im1b4)/(im1b2+im1b4) : -1 \")\n",
" app.Execute()\n",
" if (index == 0):\n",
" tab_App = [app]\n",
" else:\n",
" tab_App.append(app)\n",
" index += 1\n",
" return tab_App\n",
"\n",
"def compute_watermask(list_app,outputfilename):\n",
" app = otbApplication.Registry.CreateApplication(\"BandMath\")\n",
"\n",
" index = 1\n",
" input_list = []\n",
" exp = \"max(\"\n",
" for appli in list_app:\n",
" app.AddImageToParameterInputImageList(\"il\",appli.GetParameterOutputImage(\"out\"))\n",
" if index==1:\n",
" exp = exp+\"im\"+str(index)+\"b1\"\n",
" else:\n",
" exp = exp+\", \"+\"im\"+str(index)+\"b1\"\n",
" index += 1\n",
"\n",
" exp = exp +\")\"\n",
" print(exp)\n",
" print(input_list)\n",
"\n",
" app.SetParameterString(\"out\",outputfilename)\n",
" app.SetParameterString(\"exp\",exp)\n",
" app.Execute()\n",
" return app\n",
"\n",
"def thresholdNDWI(appWatermask,outputfilename,threshold):\n",
" app = otbApplication.Registry.CreateApplication(\"BandMath\")\n",
" app.AddImageToParameterInputImageList(\"il\", appWatermask.GetParameterOutputImage(\"out\"))\n",
" app.SetParameterString(\"exp\", \"im1b1 > \"+str(threshold)+\" ? 1 : 0\")\n",
" app.SetParameterString(\"out\", outputfilename)\n",
" app.ExecuteAndWriteOutput()\n",
" return app\n",
"\n",
"def smoothResult(thresholdNDWI, radius, outputfilename):\n",
" app = otbApplication.Registry.CreateApplication(\"BinaryMorphologicalOperation\")\n",
" app.SetParameterInputImage(\"in\",thresholdNDWI.GetParameterOutputImage(\"out\"))\n",
" app.SetParameterString(\"out\",outputfilename)\n",
" app.SetParameterString(\"structype\",\"ball\")\n",
" app.SetParameterFloat(\"structype.ball.xradius\", radius)\n",
" app.SetParameterFloat(\"structype.ball.yradius\", radius)\n",
" app.SetParameterString(\"filter\",\"closing\")\n",
" app.ExecuteAndWriteOutput()\n"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### Full OTB Pipeline"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"import os\n",
"\n",
"import utils\n",
"\n",
"# Threshold\n",
"threshold = 0.3\n",
"\n",
"# Data directory\n",
"DATA_DIR = \"data\"\n",
"\n",
"# Output directory\n",
"OUTPUT_DIR = \"output\"\n",
"\n",
"# Input / Output filenames\n",
"images_list = utils.list_images(DATA_DIR,\"_RVBPIR.tif\")\n",
"ndwi_max_fname = os.path.join(OUTPUT_DIR, \"max_NDWI2.tif\")\n",
"ndwi_thres_fname = os.path.join(OUTPUT_DIR, \"ndwi_threshold%d.tif\" % (threshold*100) )\n",
"\n",
"# Compute radiometric indices on a stack of images\n",
"list_app = app_radiometricindices(images_list)\n",
"\n",
"# Compute a water mask from NDWI2\n",
"wmask = compute_watermask(list_app, ndwi_max_fname)\n",
"\n",
"# Use a threshold to extract the mask\n",
"ndwi = thresholdNDWI(wmask, ndwi_thres_fname,threshold)\n",
"\n",
"# Smooth the result by applying a morphological operation\n",
"# smoothResult(ndwi,2.0,wdir+\"/result_smoothed.tif\")\n",
"\n",
"print (\"Full OTB Pipeline completed\")"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### Display OTB result"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"import rasterio\n",
"import display_api\n",
"\n",
"raster = rasterio.open(ndwi_thres_fname)\n",
"m, dc = display_api.rasters_on_map([raster], OUTPUT_DIR, [\"OTB result\"])\n",
"m"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### Go to the [Step 4](step4_filter_shapefile.ipynb)"
]
}
],
"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",
"version": "3.5.6"
}
},
"nbformat": 4,
"nbformat_minor": 2
}
step4_filter_shapefile.ipynb
View file @ a425da5a
......@@ -25,17 +25,9 @@
},
{
"cell_type": "code",
"execution_count": 1,
"execution_count": null,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"Converting output/ndwi_threshold30.tif into a shapefile output/watermask.sqlite\n"
]
}
],
"outputs": [],
"source": [
"import os\n",
"import subprocess\n",
......@@ -66,7 +58,7 @@
},
{
"cell_type": "code",
"execution_count": 2,
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
......@@ -86,18 +78,9 @@
},
{
"cell_type": "code",
"execution_count": 3,
"execution_count": null,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"Nb features (before filtering): 344\n",
"Nb islands (after filtering): 45\n"
]
}
],
"outputs": [],
"source": [
"from osgeo import ogr, osr\n",
"\n",
......@@ -155,24 +138,9 @@
},
{
"cell_type": "code",
"execution_count": 7,
"execution_count": null,
"metadata": {},
"outputs": [
{
"data": {
"application/vnd.jupyter.widget-view+json": {
"model_id": "17fd0284fe9f4b5c8e7a7861779afc87",
"version_major": 2,
"version_minor": 0
},
"text/plain": [
"Map(basemap={'url': 'https://{s}.tile.openstreetmap.org/{z}/{x}/{y}.png', 'attribution': 'Map data (c) <a href…"
]
},
"metadata": {},
"output_type": "display_data"
}
],
"outputs": [],
"source": [
"import rasterio\n",
"from glob import glob\n",
......@@ -195,13 +163,6 @@
"m, dc = display_api.rasters_on_map([raster], OUTPUT_DIR, [DATE], geojson_data=fc)\n",
"m"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": []
}
],
"metadata": {
......
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment