GeoFlood_local_demo.ipynb

{
  "cells": [
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "1LXi4paalFZa"
      },
      "source": [
        "# **GeoFlood demo Jupyter notebook (local)**\n",
        "\n",
        "**Purpose of this notebook**\n",
        " - Perform all operations in the GeoFlood workflow to obtain flood inundation extent and depth based on a hydrological forecast of discharge."
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "IS1ubCZ_fuXt"
      },
      "source": [
        "---\n",
        "# **(1) Prepare the GeoFlood software environment**\n",
        " - Install GRASS GIS\n",
        "   - Refer to included GRASS_ubuntu_installation.txt if running Ubuntu 20.04\n",
        " - Create GeoFlood conda environment and use as kernel for this notebook\n",
        "   - See included geoflood-env.yml\n",
        " - Install TauDEM from GitHub\n",
        "   - https://github.com/amoodie/TauDEM (build with `cmake` and `amoodie-inunmap` branch)\n",
        "   - https://github.com/dtarb/TauDEM (original)"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "E74Dff5h6001"
      },
      "source": [
        "---\n",
        "# **(2) Configuration of folders for storing and accessing inputs and outputs**"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "pEYMHNYga4Xl"
      },
      "source": [
        "1. Make a folder called `geoflood_demo` to store input and output data for this demo notebook\n",
        "2. Extract example input data to `geoflood_demo`"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "LC5El6lRHWZJ"
      },
      "outputs": [],
      "source": [
        "%%bash\n",
        "\n",
        "tar xzf INPUT.tar.gz INPUT 2> /dev/null\n",
        "find INPUT -name '.*' -type f -delete\n",
        "mkdir geoflood_demo\n",
        "mv INPUT geoflood_demo"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "Clone GeoFlood repo in `geoflood_demo` directory"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "%%bash\n",
        "\n",
        "cd geoflood_demo\n",
        "git clone https://github.com/passaH2O/GeoFlood.git"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "uosJ0fLFHY9a"
      },
      "source": [
        "This configuration script sets up the GeoFlood file structure in our `geoflood_demo` folder"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "dqxHjnrp7t9V"
      },
      "outputs": [],
      "source": [
        "!python geoflood_demo/GeoFlood/GeoNet/pygeonet_configure.py \\\n",
        "  -dir geoflood_demo \\\n",
        "  -p OnionCreek -n OC1mTest --no_chunk \\\n",
        "  --input_dir INPUT --output_dir OUTPUT\n"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "DVTEgn8Q_Sgg"
      },
      "source": [
        "---\n",
        "# **(3) Run nonlinear filtering on a portion of Onion Creek**"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "IkN9HqME_Y56"
      },
      "outputs": [],
      "source": [
        "!python geoflood_demo/GeoFlood/GeoNet/pygeonet_nonlinear_filter.py"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "7XK3A_UYzCh7"
      },
      "source": [
        "---\n",
        "# **(4) Compute flow directions and flow accumulation**\n",
        "We are now going to use GRASS GIS functions to compute flow directions, flow accumulation, and identify subbasins within our domain based on this information. All these operations are contained within the script `pygeonet_grass_py3.py`. As usual, you can open it from your cloned repo on the left to have a look. "
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "Sesc_YK3zKIv"
      },
      "outputs": [],
      "source": [
        "!python geoflood_demo/GeoFlood/GeoNet/pygeonet_grass_py3.py"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "IgMWcYh7TmsN"
      },
      "source": [
        "---\n",
        "# **(5) Compute slope and geometric curvature on the filtered terrain**\n",
        "\n",
        "We are now going to call as input our filtered landscape and we will compute the slope and geometric curvature, which is the default definition of curvature used in GeoNet and GeoFlood. All the operations are contained within `pygeonet_slope_curvature.py`. You can open the actual script by double clicking on it from the repository on the left."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "eJ14Z9mqUH4m"
      },
      "outputs": [],
      "source": [
        "!python geoflood_demo/GeoFlood/GeoNet/pygeonet_slope_curvature.py"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "KA3pxUj6UM-b"
      },
      "source": [
        "---\n",
        "# **(6) Compute skeletons based on curvature, flow accumulation, and curvature and flow accumulation combined**\n",
        "\n",
        "We will now call a new operation in this notebook, which allows us to use curvature only, flow accumulation only, and curvature and flow accumulation combined to identify which parts of the landscape are likely channelized. All these operations are contained within `pygeonet_skeleton_definition.py`."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "rp1uML3NU6BO"
      },
      "outputs": [],
      "source": [
        "!python geoflood_demo/GeoFlood/GeoNet/pygeonet_skeleton_definition.py"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "tZpmK1TITBHX"
      },
      "source": [
        "---\n",
        "# **(7) Identify end points of NHD MR channels**\n",
        "\n",
        "We are going to now use `Network_Node_reading.py` to automatically identify the end points of the input flowline - the NHD MR information for this watershed. This operation is contained within the GeoFlood portion of the code and will generate a new output in the Output folder called `OC1mTest_endPoints`. This file contains the coordinates (start and end point) of the end points of each flowline. "
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "-OEA0gD5Tb5B"
      },
      "outputs": [],
      "source": [
        "!python geoflood_demo/GeoFlood/GeoFlood/Network_Node_Reading.py"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "CkptfiyJlBbW"
      },
      "source": [
        "---\n",
        "# **(8) Network extraction with cost function including negative height**\n",
        "\n",
        "We are now going to extract the geodesic curve between end points using the default cost function that includes negative height.  \n",
        "\n",
        "**Compute negative height information**\n",
        "\n",
        "We are going to use `Relative_Height_Estimation.py`, which will return a binary raster/array with values of 1 given to pixels at a lower elevation than the elevation associated with NHD MR Flowline pixels and a value of zero to all other pixels, i.e. pixels at a higher elevation than the NHD MR Flowlines. \n",
        "\n",
        "This operation will produce three more outputs; the `_NegaHand` is the relative height we are looking for to be added to the cost function, while the other two (`_nhdflowline` and `_Allocation`) are the rasterized NHD flowline used for calculation and another term used for identifying the corridor of values below the elevation the NHD.\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "96m3Nu1_m_NT"
      },
      "outputs": [],
      "source": [
        "!python geoflood_demo/GeoFlood/GeoFlood/Relative_Height_Estimation.py"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "r8yVVs0XnliP"
      },
      "source": [
        "---\n",
        "# **Network extraction with negative height**\n",
        "\n",
        "We can now run our `Network_Extraction.py` to obtain the geodesic curve between the end points. \n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "84mvlMtFn1yn"
      },
      "outputs": [],
      "source": [
        "!python geoflood_demo/GeoFlood/GeoFlood/Network_Extraction.py"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "fSG90g96wgS_"
      },
      "source": [
        "---\n",
        "# **(9) TauDEM preparation**"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "ttGOROEDuV5x"
      },
      "source": [
        "TauDEM can leverage a computer's resources by running parallel processes.\n",
        "We can see from the command below the number of cores available. "
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "!lscpu | grep 'per socket'"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "Let's save the number of cores in order to include it in our TauDEM operations"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "OzTAHcEuw1Nu"
      },
      "outputs": [],
      "source": [
        "numcores = !lscpu | grep 'per socket'\n",
        "numcores = int(numcores[0][-1])\n",
        "numcores"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "M_Dxsnv7oJ64"
      },
      "source": [
        "---\n",
        "# **(10) TauDEM operations to prepare for the HAND raster calculation**\n",
        "\n",
        "We will now compute the HAND raster. As it is based on Dinf flow directions, we will need to perform pit filling on the elevation data and compute flow directions. Note we are running with 8 processes."
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "X3o5hkvo6vFA"
      },
      "source": [
        "---\n",
        " - **Pit filling**"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "JdY4OlNe6Ryv"
      },
      "outputs": [],
      "source": [
        "!mpiexec -n $numcores pitremove \\\n",
        "    -z geoflood_demo/INPUT/GIS/OnionCreek/OC1mTest.tif \\\n",
        "    -fel geoflood_demo/OUTPUT/GIS/OnionCreek/OC1mTest_fel.tif"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "ZiVKbCId61Wj"
      },
      "source": [
        "---\n",
        " - **Dinf flow directions**"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "zijBtlKr6WK4"
      },
      "outputs": [],
      "source": [
        "!mpiexec -n $numcores dinfflowdir \\\n",
        "    -fel geoflood_demo/OUTPUT/GIS/OnionCreek/OC1mTest_fel.tif \\\n",
        "    -ang geoflood_demo/OUTPUT/GIS/OnionCreek/OC1mTest_ang.tif \\\n",
        "    -slp geoflood_demo/OUTPUT/GIS/OnionCreek/OC1mTest_slp.tif"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "dgMiRiz_67LH"
      },
      "source": [
        "---\n",
        "# **(11) HAND calculation with GeoFlood flowline**\n",
        "\n",
        "We will now use the previous outputs from TauDEM to compute the HAND raster based on the GeoFlood flowline. "
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "5uIJe3n26-EF"
      },
      "outputs": [],
      "source": [
        "!mpiexec -n $numcores dinfdistdown \\\n",
        "    -ang geoflood_demo/OUTPUT/GIS/OnionCreek/OC1mTest_ang.tif \\\n",
        "    -fel geoflood_demo/OUTPUT/GIS/OnionCreek/OC1mTest_fel.tif \\\n",
        "    -slp geoflood_demo/OUTPUT/GIS/OnionCreek/OC1mTest_slp.tif \\\n",
        "    -src geoflood_demo/OUTPUT/GIS/OnionCreek/OC1mTest_path.tif \\\n",
        "    -dd geoflood_demo/OUTPUT/GIS/OnionCreek/OC1mTest_hand_GeoFlood.tif \\\n",
        "    -m ave v"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "_xsU-Y0H4AV7"
      },
      "source": [
        "---\n",
        "# **(12) Calculation of hydraulic geometry for the stream reach**\n",
        "\n",
        "Now that we have the HAND raster, we will go back to GeoFlood functions to compute the hydraulic geometry properties for this river reach, which will then be used to translate the discharge forecast into a depth prediction for inundation mapping (in the next lab).\n",
        "\n",
        "The first operation we are going to perform is the segmentation of the flowline which is done within `Streamline_Segmentation.py`. This operation checks the length of the flowline; if the length is longer than a predefined maximum segment length, the script breaks the original flowline into equal-length segments with a length shorter than the predefined threshold. The input of this script is the extracted river network shapefile (`OC1mTest_channelNetwork.shp`), and the output is a river segment shapefile (`OC1mTest_channelSegment.shp`)."
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "eTj6OX6j4Vde"
      },
      "source": [
        "\n",
        "\n",
        "---\n",
        "\n",
        "# **Channel segmentation**"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "oj1rKSF44a2R"
      },
      "outputs": [],
      "source": [
        "!python geoflood_demo/GeoFlood/GeoFlood/Streamline_Segmentation.py"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "lPwYMy1g4gOE"
      },
      "source": [
        "\n",
        "\n",
        "---\n",
        "# **Catchment delineation**\n",
        "\n",
        "Now we will delineate the catchment for the segmented flowline. Similarly to what we have done in the past, the script uses a GRASS function called from `Segment_Catchment_Delineation.py`.\n",
        "\n",
        "The river segment shapefile (`OC1mTest_channelSegment.shp`) is one of the script inputs. Since this script relies on the GRASS GIS catchment delineation routine, the flow direction raster (`OC1mTest_fdr.tif`) we derived with GRASS is another input. The output of this script is a segment catchment geotiff file (`OC1mTest_segmentCatchment.tif`)."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "tcB80-F14lUM"
      },
      "outputs": [],
      "source": [
        "!python geoflood_demo/GeoFlood/GeoFlood/Grass_Delineation_py3.py"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "DbNsIhvb4128"
      },
      "source": [
        "\n",
        "\n",
        "---\n",
        "\n",
        "\n",
        "# **Calculation of hydraulic geometry properties**\n",
        "\n",
        "We can now use the segments, the delineated catchments, and the HAND raster to obtain the hydraulic geometry properties of each of these segments. We will start by running the GeoFlood script `River_Attribute_Estimation.py`, which is used to compute bed slope, river length, and catchment area. The inputs to this script are the DEM in the input folder, the river segments we have previously identified, and the segment catchment raster we just got from GRASS in our `OUTPUT` folder. This script will output a segment catchment polygon shapefile and a river attribute text file, which will be placed in the Hydraulics folder in your `OUTPUT/Hydraulics` folder.\n",
        "\n",
        "You will find new outputs in your Outputs folder. In particular, in the GIS folder you will find the _segmentCatchment.shp and in the Hydraulics folder you will find the river attribute text file."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "rxvCiB-d46rh"
      },
      "outputs": [],
      "source": [
        "!python geoflood_demo/GeoFlood/GeoFlood/River_Attribute_Estimation.py"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "4cw9oEk2498T"
      },
      "source": [
        "---\n",
        "# **Associate channel roughness to segments**\n",
        "We will now run the script `Network_Mapping.py` which takes the channel roughness information from the original river network and saves it in the segmented river network. The script takes as input the original catchment shapefile, which is the new input we just added as part of this lab (`Catchment.shp`) and the river segment shapefile we have created (`OC1mTest_channelSegment.shp`).\n",
        "\n",
        "This script will create a new output in your `OUTPUT` GIS folder which has csv format and contains the COMID of each of the segments in our flowline."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "sKcAnChz5Cr1"
      },
      "outputs": [],
      "source": [
        "!python geoflood_demo/GeoFlood/GeoFlood/Network_Mapping.py"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "1wqNVZg_5Kb8"
      },
      "source": [
        "---\n",
        "# **Derivation of channel hydraulic properties**\n",
        "The next script is part of TauDEM and is called CatchHydroGeo. This script is used to derive the channel hydraulic properties using the terrain analysis products. The inputs are the HAND raster we computed (`OC1mTest_HAND.tif`), the segment catchment raster (`OC1mTest_segmentCatchment.tif`), the Dinf slope raster (`OC1mTest_slp.tif`), the river attribute text file (`OC1mTest_River_Attribute.txt`) and the stage range text file we input.\n",
        "\n",
        "Once the execution is finished, the resulting channel characteristic table with hydraulic properties directly derived from the HAND raster (`hydroprop-basetable.csv`) will be found in the `OUTPUT/Hydraulics` folder and contains the hydraulic geometry properties for each stage value we input."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "13O0Fqx25PhW"
      },
      "outputs": [],
      "source": [
        "!mpiexec -n $numcores catchhydrogeo \\\n",
        "    -hand geoflood_demo/OUTPUT/GIS/OnionCreek/OC1mTest_hand_GeoFlood.tif \\\n",
        "    -catch geoflood_demo/OUTPUT/GIS/OnionCreek/OC1mTest_segmentCatchment.tif \\\n",
        "    -catchlist geoflood_demo/OUTPUT/Hydraulics/OnionCreek/OC1mTest_River_Attribute.txt \\\n",
        "    -slp geoflood_demo/OUTPUT/GIS/OnionCreek/OC1mTest_slp.tif \\\n",
        "    -h geoflood_demo/INPUT/Hydraulics/OnionCreek/stage.txt \\\n",
        "    -table geoflood_demo/OUTPUT/Hydraulics/OnionCreek/hydroprop-basetable.csv"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "AHbNVAg85TSy"
      },
      "source": [
        "---\n",
        "# **Derivation of synthetic rating curves**\n",
        "We have a final GeoFlood script to run, called `Hydraulic_Property_Postprocess.py`, which combines the channel roughness coefficient and the properties that we have obtained in the above steps to derive reach average geometric properties for the reach and synthetic rating curves. The inputs are the COMID-based Manning’s values (`COMID_Roughness.csv`), the network ID mapping csv table (`OC1mTest_networkMapping.csv`), and the base geometric properties (`hydroprop-basetable.csv`).\n",
        "\n",
        "You will find as output the full channel properties table with synthetic rating curves in the `OUTPUT` Hydraulics folder, called `hydroprop-fulltable.csv`."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "_mvUVk2b5YfA"
      },
      "outputs": [],
      "source": [
        "!python geoflood_demo/GeoFlood/GeoFlood/Hydraulic_Property_Postprocess.py"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "pipY123QHXxI"
      },
      "source": [
        "# **(13) Estimation of flood inundation extent and depth**\n",
        "\n",
        "And finally, what we were all waiting for! We will integrate all the results above to estimate flood inundation extent and depth based on the hydrological prediction provided as input which comes from the National Water Model. \n",
        "\n",
        "https://water.noaa.gov/about/nwm\n",
        "\n",
        "---\n",
        "# **Assign NWM input forecasted discharge to channel segments**\n",
        "\n",
        "This first operation takes the NWM forecast that we provided in input and will assign it to each channel segment based on COMID. The script performs the discharge-stage conversion based on the NWM netCDF output (general format: `nwm.txxz.product_type.channel_rt.fxxx.conus.nc`) from the NWM input folder. The script also uses the network id mapping table, which records how the COMIDs are linked to the ids of the extracted river segments. The script also uses the full channel property table, which contains stage-discharge data pairs for each segment used in stage height interpolation. \n",
        "\n",
        "The script produces the depth corresponding to the discharge for each segment in both netCDF and CSV formats. The depths files have the same file name as the NWM inputs and they will be stored in the corresponding OUPUT NWM folder. \n",
        "\n",
        "Since our stream does not flood very often we are going to overwrite the discharge forecast so that we can make this system flood. \n",
        "\n",
        "<mark>**Modify forecasted discharge in Foreast_Table.py**</mark>\n",
        "\n",
        "Open `Forecast_Table.py` from the left menu - the script is stored within the GeoFlood folder. Double click on the name of the script and uncomment line 63 by removing the symbol `#`. As you can see we are imposing a discharge of 500 to obtain a flood inundation map (do not worry about any of the symbols that are written on that same line after 500). \n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "CuB9GkGDJA9f"
      },
      "outputs": [],
      "source": [
        "!python geoflood_demo/GeoFlood/GeoFlood/Forecast_Table.py \\\n",
        "    geoflood_demo/INPUT/NWM/OnionCreek/nwm.t00z.analysis_assim.channel_rt.tm01.conus.nc"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "DmkxJ5C_dhdg"
      },
      "source": [
        "\n",
        "\n",
        "---\n",
        "# **Compute flood inundation map corresponding to the forecasted discharge**\n",
        "\n",
        "We will now quantify the flood extent and depth by using TauDEM's InunMap. \n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "0GJFXVo-f-bF"
      },
      "outputs": [],
      "source": [
        "!mpiexec -n $numcores inunmap \\\n",
        "         -hand geoflood_demo/OUTPUT/GIS/OnionCreek/OC1mTest_hand_GeoFlood.tif \\\n",
        "         -catch geoflood_demo/OUTPUT/GIS/OnionCreek/OC1mTest_segmentCatchment.tif \\\n",
        "         -forecast geoflood_demo/OUTPUT/NWM/OnionCreek/nwm.t00z.analysis_assim.channel_rt.tm01.conus.nc \\\n",
        "         -mapfile geoflood_demo/OUTPUT/Inundation/OnionCreek/OC1mTest_NWM_inunmap.tif"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "6Ki4fp2gTP4G"
      },
      "source": [
        "# **(14) Plot the inundation map!**\n",
        "\n",
        "We can now have a look at the results here and check our inundation map. "
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "3PiwJSlgT3Y9"
      },
      "outputs": [],
      "source": [
        "import matplotlib.pyplot as plt\n",
        "import numpy as np\n",
        "import rasterio\n",
        "import rasterio.plot\n",
        "\n",
        "# open the HAND raster\n",
        "with rasterio.open('geoflood_demo/OUTPUT/GIS/OnionCreek/OC1mTest_hand_GeoFlood.tif') as src:\n",
        "    oc_hand_image = src.read(1) # read band 1 into a numpy array\n",
        "\n",
        "# open the inundation raster\n",
        "with rasterio.open('geoflood_demo/OUTPUT/Inundation/OnionCreek/OC1mTest_NWM_inunmap.tif') as src:\n",
        "    oc_inun_image = src.read(1) # read band 1 into a numpy array\n",
        "\n",
        "fig, ax = plt.subplots(1, 2, figsize=(16,5))\n",
        "\n",
        "hand_im = ax[0].imshow(oc_hand_image, cmap='Blues_r', vmin=0, vmax=80)\n",
        "fig.colorbar(hand_im, ax=ax[0], shrink=0.4)\n",
        "ax[0].set_title(\"GeoFlood HAND\")\n",
        "ax[0].set_xlabel('x (m)')\n",
        "ax[0].set_ylabel('y (m)')\n",
        "\n",
        "inun_im = ax[1].imshow(oc_inun_image, cmap='Blues_r', vmin=0, vmax=8)\n",
        "fig.colorbar(inun_im, ax=ax[1], shrink=0.4)\n",
        "ax[1].set_title(\"Inundation map from NWM forecast\")\n",
        "ax[1].set_xlabel('x (m)')\n",
        "ax[1].set_ylabel('y (m)')\n",
        "\n",
        "fig.tight_layout()\n",
        "\n",
        "# save figure \n",
        "fig.savefig('demo.png', bbox_inches='tight')\n",
        "# fig.show()"
      ]
    }
  ],
  "metadata": {
    "colab": {
      "collapsed_sections": [],
      "include_colab_link": true,
      "name": "GeoFlood_demo.ipynb",
      "provenance": []
    },
    "kernelspec": {
      "display_name": "Python 3.8.6 ('geoflood-env')",
      "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.8.6"
    },
    "vscode": {
      "interpreter": {
        "hash": "d060aa1cd3282031c5d92abe05edcc06d33f6053b78819f748f49e53977cc352"
      }
    }
  },
  "nbformat": 4,
  "nbformat_minor": 0
}