scripts/intervention_planning.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Intervention planning\n",
    "=====================\n",
    "\n",
    "Overview\n",
    "--------\n",
    "\n",
    "This notebook facilitates using the RiverScape Python package (Straatsma\n",
    "and Kleinhans, 2018) to parameterize and position landscaping measures\n",
    "and update the input data for the two-dimensional (2D) flow model\n",
    "Delft3D Flexible Mesh (DFM).\n",
    "\n",
    "For the current notebook version, we use the River Waal, which is the\n",
    "main distributary of the River Rhine in the Netherlands. For general\n",
    "concepts and detailed description of the approach used here in these\n",
    "notebooks we refer to the publications [Straatsma et\n",
    "al. (2017)](https://advances.sciencemag.org/content/3/11/e1602762) and\n",
    "[Straatsma et al. (2019)](https://doi.org/10.5194/nhess-19-1167-2019).\n",
    "\n",
    "How to start\n",
    "------------\n",
    "\n",
    "### Setting up the environment\n",
    "\n",
    "To run this notebook, please import the following Python modules."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Import standard modules\n",
    "import os\n",
    "import sys\n",
    "import string\n",
    "import subprocess\n",
    "import time\n",
    "import math\n",
    "import pprint as pp\n",
    "\n",
    "# Import required modules/packages\n",
    "import numpy as np\n",
    "import pandas as pd\n",
    "import geopandas as gpd\n",
    "from scipy.spatial import Delaunay\n",
    "from shapely.geometry import MultiLineString\n",
    "from shapely.ops import cascaded_union, polygonize\n",
    "\n",
    "import geoviews\n",
    "geoviews.extension('bokeh')\n",
    "\n",
    "from collections import OrderedDict\n",
    "import pcraster as pcr"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Please also make sure that this notebook file is in the same folder as\n",
    "the RiverScape Python module files that must be loaded. You can then\n",
    "import the RiverScape Python module:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import riverscape\n",
    "from riverscape import pcrRecipes\n",
    "from riverscape import msr\n",
    "\n",
    "%reload_ext autoreload\n",
    "%autoreload 2"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "%autosave 0"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Input and output folders\n",
    "\n",
    "In the following please define the input and output folders."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Default locations\n",
    "input_dir = riverscape.input_data_path()\n",
    "output_folder_name = 'output_intervent'\n",
    "output_dir = os.path.join(os.getcwd(), output_folder_name)\n",
    "if not os.path.exists(output_dir):\n",
    "    os.mkdir(output_dir)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "You may also want to set the folder interactively. You need to uncomment\n",
    "the following lines:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# input_dir  = riverscape.select_directory()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# output_dir = riverscape.select_directory()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Finally, some temporary folder for calculations will be created:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Create scratch directory and go to this folder\n",
    "scratch_dir  = os.path.join(output_dir, \"tmp\")\n",
    "pcrRecipes.make_dir(scratch_dir)\n",
    "os.chdir(scratch_dir)\n",
    "\n",
    "# print some statements about the folder locations:\n",
    "msg = \"The input folder is set on:   {}\".format(input_dir)\n",
    "print(msg)\n",
    "msg = \"The output folder is set on:  {}\".format(output_dir)\n",
    "print(msg)\n",
    "msg = \"The scratch folder is set on: {}\".format(scratch_dir)\n",
    "print(msg)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Start processing/calculation\n",
    "----------------------------\n",
    "\n",
    "To start processing, please load the following cells in order to set and\n",
    "test some basic configuration."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# go to the scratch folder\n",
    "pcrRecipes.make_dir(scratch_dir)\n",
    "os.chdir(scratch_dir)\n",
    "\n",
    "# set global option for PCRaster such that length of cells is computed in true length of cells\n",
    "pcr.setglobaloption('unittrue')\n",
    "\n",
    "# set the pcraster clone map\n",
    "current_dir = os.path.join(input_dir, 'reference_maps')\n",
    "pcr.setclone(os.path.join(current_dir, 'clone.map'))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Loading input files\n",
    "\n",
    "By running the following cells, the input files would be read. These\n",
    "input files consist of the following attributes of current conditions:\n",
    "\n",
    "-   main\\_dike: current/existing river embankment properties,\n",
    "    e.g. length, volume and height\n",
    "-   minemb: minor embankment properties, e.g. length, volume and height\n",
    "-   groynse: groyne properties, e.g. length, volume and height\n",
    "-   hydro: hydrodynamics (delft3d-fm) attributes, e.g. chezy, nikuradse,\n",
    "    specific discharge, velocity, water depth, water level, etc.\n",
    "-   mesh: delft3d mesh\n",
    "-   axis: river attributes, e.g. location, radius, turning direction,\n",
    "    velocity, water depth, water level, etc.\n",
    "-   geom: river geometry attributes, e.g. clone, dem,\n",
    "    dist\\_to\\_main\\_dike, dist\\_to\\_groyne\\_field,\n",
    "    dist\\_to\\_main\\_channel, flpl\\_width, flpl\\_narrow, flpl\\_wide,\n",
    "    main\\_channel\\_width, river\\_side, shore\\_line\n",
    "-   lulc: land use and land cover attributes, e.g. backwaters, ecotopos,\n",
    "    floodplain, groyne\\_field, main\\_channel, trachytopes, sections,\n",
    "    winter\\_bed, real\\_estate\\_value\n",
    "\n",
    "For further information about them, please check [Straatsma and\n",
    "Kleinhans (2018)](https://doi.org/10.1016/j.envsoft.2017.12.010)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# change to the 'current_dir' (input data) for reading/importing input data\n",
    "os.chdir(current_dir)\n",
    "\n",
    "\n",
    "# reading current/existing river embankment properties for main dikes, minor embankments and groynes\n",
    "# - for each, this will return location, length, volume, and height\n",
    "main_dike = msr.read_dike_maps(current_dir, 'main_dike')\n",
    "minemb    = msr.read_dike_maps(current_dir, 'minemb')\n",
    "groynes   = msr.read_dike_maps(current_dir, 'groyne')\n",
    "\n",
    "# reading Hydrodynamics (delft3d-fm) attributes\n",
    "# - chezy, nikuradse, specific discharge, velocity, water depth, water level, etc\n",
    "hydro = msr.read_hydro_maps(current_dir)\n",
    "# - hydrodynamic mesh\n",
    "mesh = msr.read_mesh_maps(current_dir)\n",
    "\n",
    "# reading RiverScape attributes\n",
    "# - axis: location, radius, turning_direction, rkm, rkm_point, rkm_line, rkm_full\n",
    "axis = msr.read_axis_maps(current_dir)\n",
    "# - geometry: clone, dem, distance to main dike, distance to groyne, distance to main channel, floodplain widths, contiguous narrow floodplain, contiguous wide floodplain\n",
    "geom = msr.read_geom_maps(current_dir)\n",
    "# - land use land cover attributes\n",
    "lulc = msr.read_lulc_maps(current_dir)\n",
    "\n",
    "msg = \"Files should be already read from: \" + current_dir\n",
    "print(\"\\n\")\n",
    "print(msg)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Listing attributes and variables\n",
    "\n",
    "After succesfully reading the files/maps, you should be able to list the\n",
    "attributes/variables of main\\_dike, minemb, groynes, hydro, mesh, axis,\n",
    "geom, and lulc in the following."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "for obj in [\"main_dike\", \"minemb\", \"groynes\", \"hydro\", \"mesh\", \"axis\", \"geom\", \"lulc\"]:\n",
    "    obj_vars = vars()[obj].__dict__.keys()\n",
    "    print('{}'.format(obj))\n",
    "\n",
    "    for var in obj_vars:\n",
    "        print('\\t{}'.format(var))\n",
    "\n",
    "    print()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Visualizing existing attributes\n",
    "\n",
    "You may want to inspect existing attributes by plotting the\n",
    "corresponding raster maps. You can use the plot function for that, for\n",
    "example to plot the digital elevation model type:\n",
    "\n",
    "``` python\n",
    "riverscape.plot(geom.dem)\n",
    "```"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# plot river main channel width\n",
    "riverscape.plot(geom.main_channel_width)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# plot height of main_dike\n",
    "riverscape.plot(main_dike.height)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# plot flpl_wide IDs\n",
    "riverscape.plot(geom.flpl_wide)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Initiating the River and its Measures\n",
    "-------------------------------------\n",
    "\n",
    "Given the aforementioned attributes, the River Waal and its current\n",
    "measured would be initiated by executing the following cells."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# The River Wall is initiated based on the aforementioned given attributes.\n",
    "waal = msr.River('Waal', axis, main_dike, minemb, groynes, hydro, mesh, lulc, geom)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Initiate the Measures for the River.\n",
    "waal_msr = msr.RiverMeasures(waal)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Specifing the different measures\n",
    "================================\n",
    "\n",
    "First, you can have a look at the current specification of the ecotope\n",
    "and trachytope classes present in the area. The floodplain sections are\n",
    "depicted as well. Ecotopes are homogeneous ecological landscape units\n",
    "w.r.t. vegetation structure or succession stage. Trachytopes are\n",
    "spatially-distributed roughness values for the channel.\n",
    "\n",
    "Open the maps with the following command. Note that generating the plots\n",
    "may take a moment."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "riverscape.plot_eco_trachy_sec()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Side channel measure:\n",
    "---------------------\n",
    "\n",
    "### Specify your own side channel properties:\n",
    "\n",
    "The measures are configured with certain properties. You can inspect and\n",
    "change a few of them:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "settings = riverscape.measures_settings()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "You may want to modify the side channel properties using the following\n",
    "interactive cell. Note that height and width are not true to scale."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "channel_values = riverscape.channel_properties()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Load the following cell to use your configuration."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "settings['channel_width'] = channel_values.kwargs['width']\n",
    "settings['channel_depth'] = channel_values.kwargs['depth']\n",
    "settings['channel_slope'] = channel_values.kwargs['slope']\n",
    "waal_msr.settings = settings"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "You can also print the current settings to check whether they are\n",
    "suitable:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "for cur in settings.items():\n",
    "    print('{:21s}: {}'.format(cur[0], cur[1]))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Selecting the region mask for measure area:\n",
    "\n",
    "Please select the areas where you want to perform this measure.\n",
    "\n",
    "Please also give a label for this measure."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "label = 'custom_label'"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "selection = riverscape.select_area(waal.geom.flpl_wide)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "mask = riverscape.generate_mask(waal.geom.flpl_wide, selection)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "riverscape.plot(mask)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "chosen_flpl_wide = pcr.ifthen(mask, waal.geom.flpl_wide)\n",
    "riverscape.plot(chosen_flpl_wide)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Implementing the measure\n",
    "\n",
    "By running the following cell, the measure will be implemented. Note\n",
    "that this step can consume noticeable computing time, depending on the\n",
    "number of areas you selected."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# measure by side channel constrcution\n",
    "# - this includes looping over floodplain IDs\n",
    "chan_msr = waal_msr.side_channel_measure(settings, mask = mask, ID = label)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Exploring the measure\n",
    "\n",
    "You can explore the measure implemented by running the following\n",
    "plotting cell. Plotting may take a while before completed."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# plot/explore side channel measure\n",
    "chan_msr.plot()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Floodplain lowering measure\n",
    "---------------------------\n",
    "\n",
    "For the floodplain lowering measure you need to select new areas where\n",
    "you want to introduce this measure. First, give a new identifier:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "ID = 'everywhere'"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Then specify the sections that will form the new mask:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "sections = pcr.readmap('flpl_sections.map')\n",
    "selection = riverscape.select_area(sections)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "and generate the new mask"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "mask = riverscape.generate_mask(sections, selection)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "For checking purposes you could plot the mask map."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "As before, you can now perform the measure and visualise the results."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# floodplain lowering measure\n",
    "lowering_msr = waal_msr.lowering_measure(settings, mask=mask, ID=ID)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# plot/explore floodplain lowering measure\n",
    "lowering_msr.plot()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Groyne lowering measure:\n",
    "------------------------\n",
    "\n",
    "Please set the ID/label for this measure and set the mask where you want\n",
    "to introduce this measure."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "ID = 'everywhere'"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "For the groyne lowering measure you can specify a specific area, defined\n",
    "by the distances from . The minimum and maximum values should be between\n",
    "867 and 960."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "min_value = 890\n",
    "max_value = 910\n",
    "rkm = pcr.readmap('rkm_full.map')\n",
    "mask = pcr.ifthen((rkm >= min_value) & (rkm <= max_value), pcr.boolean(1))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Plot the new mask:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# groyne lowering measure\n",
    "groyne_low_msr = waal_msr.groyne_lowering_msr(settings, mask=mask, ID=ID)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# plot/explore groyne lowering measure\n",
    "groyne_low_msr.plot()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Minor embankment lowering measure:\n",
    "----------------------------------\n",
    "\n",
    "Please set the ID/label for this measure and set the mask where you want\n",
    "to introduce this measure."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "mask = pcr.boolean(1)\n",
    "ID = 'everywhere'"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# measure by minor embankment lowering\n",
    "minemb_low_msr = waal_msr.minemb_lowering_msr(settings, mask=mask, ID=ID)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# plot/explore minor embankment lowering\n",
    "minemb_low_msr.plot()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Main dike raising measure:\n",
    "--------------------------\n",
    "\n",
    "Please set the ID/label for this measure and set the mask where you want\n",
    "to introduce this measure."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "mask = pcr.boolean(1)\n",
    "ID = 'everywhere'"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# measure by main dike raising\n",
    "main_dike_raise_msr = waal_msr.main_dike_raising_msr(settings, mask=mask, ID=ID)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# plot/explore main dike raising measure\n",
    "main_dike_raise_msr.plot()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Roughness smoothing measure:\n",
    "----------------------------\n",
    "\n",
    "Please set the ID/label for this measure and set the mask where you want\n",
    "to introduce this measure."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "mask = pcr.boolean(1)\n",
    "ID = 'everywhere'"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# measure by roughness smoothing\n",
    "smooth_msr = waal_msr.smoothing_measure(settings, mask=mask, ID=ID)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# plot/explore measure by roughness smoothing\n",
    "smooth_msr.plot()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Saving measures to the disk\n",
    "===========================\n",
    "\n",
    "You can store the PCRaster output maps and setting files of your\n",
    "measures to disk:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# list of measures\n",
    "msr_list = [groyne_low_msr, minemb_low_msr,\n",
    "            main_dike_raise_msr, lowering_msr, chan_msr, smooth_msr]\n",
    "# - preparing the directory\n",
    "msr_root_dir = os.path.join(output_dir, 'measures_ensemble', 'maps')\n",
    "pcrRecipes.make_dir(msr_root_dir)\n",
    "for measure in msr_list:\n",
    "    msr.write_measure(measure, msr_root_dir)"
   ]
  }
 ],
 "nbformat": 4,
 "nbformat_minor": 5,
 "metadata": {}
}