{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Exporting model data to GIS\n",
    "This notebook shows how to export model data so it can be viewed in GIS (QGIS or ArcMAP)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import os\n",
    "import xarray as xr\n",
    "from IPython.display import FileLink\n",
    "\n",
    "import nlmod"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "nlmod.util.get_color_logger(\"INFO\")\n",
    "nlmod.show_versions()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Model types"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### structured grid"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "source_ws_struc = os.path.join(\"..\", \"examples\", \"01_basic_model\")\n",
    "model_name = \"IJmuiden\"\n",
    "output_ws = \"08_gis\"\n",
    "os.makedirs(output_ws, exist_ok=True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "ds_struc = xr.load_dataset(\n",
    "    os.path.join(source_ws_struc, f\"{model_name}.nc\"), mask_and_scale=False\n",
    ")\n",
    "ds_struc"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# create gisdir\n",
    "gisdir_struc = os.path.join(output_ws, \"structured\")\n",
    "os.makedirs(gisdir_struc, exist_ok=True)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### vertex grid"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "source_ws_vert = os.path.join(\"..\", \"examples\", \"03_local_grid_refinement\")\n",
    "model_name = \"IJm_planeten\""
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "ds_vert = xr.load_dataset(\n",
    "    os.path.join(source_ws_vert, f\"{model_name}.nc\"), mask_and_scale=False\n",
    ")\n",
    "ds_vert"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# create gisdir\n",
    "gisdir_vert = os.path.join(output_ws, \"vertex\")\n",
    "os.makedirs(gisdir_vert, exist_ok=True)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Export vector data"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### structured"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# write model data to a geopackage\n",
    "fname_geopackage = nlmod.gis.ds_to_vector_file(ds_struc, gisdir=gisdir_struc)\n",
    "\n",
    "# get download link\n",
    "FileLink(fname_geopackage, result_html_prefix=\"klik hier om te downloaden -> \")"
   ]
  },
  {
   "cell_type": "raw",
   "metadata": {
    "raw_mimetype": "text/x-python",
    "scrolled": true,
    "tags": []
   },
   "source": [
    "# write model data to multiple shapefiles\n",
    "fnames = nlmod.gis.ds_to_vector_file(\n",
    "    ds_struc, driver=\"ESRI Shapefile\", gisdir=gisdir_struc\n",
    ")\n",
    "\n",
    "# get download link\n",
    "FileLinks(gisdir_struc, included_suffixes=\".shp\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### vertex grid"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# write model data to a geopackage\n",
    "fname_geopackage = nlmod.gis.ds_to_vector_file(ds_vert, gisdir=gisdir_vert)"
   ]
  },
  {
   "cell_type": "raw",
   "metadata": {
    "raw_mimetype": "text/x-python",
    "scrolled": true
   },
   "source": [
    "# write model data to multiple shapefiles\n",
    "nlmod.gis.ds_to_vector_file(\n",
    "    ds_vert, driver=\"ESRI Shapefile\", gisdir=gisdir_vert\n",
    ")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Export griddata\n",
    "\n",
    "The model data can be exported to a netcdf file that can be visualised in Qgis. For a structured model the standard model dataset (xarray.Dataset) can be exported to a netdf file. For a vertex model you have to convert the model dataset to a certain format before you can write it to a netcdf and read it with Qgis. With the code below we export the vertex model dataset to a netcdf file ('model_qgis.nc') that can be read using Qgis."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### structured grid\n",
    "\n",
    "In order to load a structured grid netcdf file in Qgis use the \"Add Raster Layer..\" option."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# write model data to a netcdf file\n",
    "fname = os.path.join(gisdir_struc, \"model_struc_qgis.nc\")\n",
    "ds_struc.to_netcdf(fname)\n",
    "\n",
    "# get download link\n",
    "FileLink(fname, result_html_prefix=\"klik hier om te downloaden -> \")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### vertex grid\n",
    "\n",
    "In order to load a vertex netcdf file in Qgis use the \"Add Mesh Layer..\" option."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# write model data to a netcdf file\n",
    "fname = os.path.join(gisdir_vert, \"model_vert_qgis.nc\")\n",
    "out = nlmod.gis.ds_to_ugrid_nc_file(ds_vert, fname)\n",
    "\n",
    "# get download link\n",
    "FileLink(fname, result_html_prefix=\"klik hier om te downloaden -> \")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Add symbology (QGIS)\n",
    "\n",
    "It is always nice to have automatic symbology for your vector data. Some thoughts:\n",
    "\n",
    "- QGIS can save symbology of a single shapefile in a .qml file\n",
    "- In QGIS you can add a .qml file to a geopackage thus saving the symbology to a single file.\n",
    "- You can create a .qml file in QGIS from existing symbology.\n",
    "- a .qml file is an xml file so theoretically it is possible to modify a .qml file with Python based on the properties of the data."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Some limitations of the current gis functions:\n",
    "\n",
    "- when exporting shapefiles to gis, attributes cannot have names longer\n",
    "than 10 characters. Now the automatic character shortening of fiona is\n",
    "used. This is not optimal.\n",
    "- when exporting data variables with dimension time only the mean values\n",
    "in time are exported in the shapefile to avoid extremely large files."
   ]
  }
 ],
 "metadata": {
  "language_info": {
   "name": "python"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}