{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "163f60be",
   "metadata": {},
   "source": [
    "# Generating model datasets\n",
    "\n",
    "This notebook contains two workflows to generate a model dataset and add layer\n",
    "information to it:\n",
    "\n",
    "- Get data (download or read from file)\n",
    "- Regrid (convert data to match your desired grid)\n",
    "\n",
    "or\n",
    "\n",
    "- Create grid (define your desired grid)\n",
    "- Add data to dataset (add data based on this grid)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "4ef9418d",
   "metadata": {},
   "outputs": [],
   "source": [
    "from shapely.geometry import LineString\n",
    "\n",
    "import nlmod"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "1b23f82f",
   "metadata": {},
   "outputs": [],
   "source": [
    "nlmod.util.get_color_logger(\"INFO\")\n",
    "nlmod.show_versions()"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "13fb7319",
   "metadata": {},
   "source": [
    "## Get REGIS data"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "67d3b65b",
   "metadata": {},
   "outputs": [],
   "source": [
    "extent = [100000, 101000, 400000, 401000]\n",
    "\n",
    "# downlaod regis\n",
    "regis = nlmod.read.regis.download_regis(extent)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "e96dbc81",
   "metadata": {},
   "outputs": [],
   "source": [
    "f, ax = nlmod.plot.get_map(extent, figsize=5)\n",
    "nlmod.plot.data_array(regis[\"top\"], edgecolor=\"k\");"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "49a18b1b",
   "metadata": {},
   "source": [
    "## Define some properties of the grid\n",
    "We choose a resolution of the calculation grid (delr and delc) larger than the resolution of regis (100 meter)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "73677bcd",
   "metadata": {},
   "outputs": [],
   "source": [
    "delr = 200.0\n",
    "delc = 200.0\n",
    "\n",
    "polygon = LineString([(extent[0], extent[2]), (extent[1], extent[3])]).buffer(10)\n",
    "refinement_features = [([polygon], \"polygon\", 3)]"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "d22c22ff",
   "metadata": {},
   "source": [
    "## Regrid data\n",
    "create a model dataset from regis and refine."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "efb9f87d",
   "metadata": {},
   "outputs": [],
   "source": [
    "ds = nlmod.to_model_ds(regis, delr=delr, delc=delc)\n",
    "ds = nlmod.grid.refine(ds, model_ws=\"01_layer_generation\", refinement_features=refinement_features)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "cf9fbfd6",
   "metadata": {},
   "source": [
    "When we plot the model top, we see that all information has a resolution of 200 meter, also in the smaller cells. The original resolution of the regis data is 100 meter however. So information is lost in the process."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "6f46ad4d",
   "metadata": {},
   "outputs": [],
   "source": [
    "f, ax = nlmod.plot.get_map(extent, figsize=5)\n",
    "nlmod.plot.data_array(ds[\"top\"], ds=ds, edgecolor=\"k\");"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "e976498d",
   "metadata": {},
   "source": [
    "## Add data to existing grid\n",
    "We can also first create a grid first, and then warp the information from regis to this grid. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "af0c984f",
   "metadata": {},
   "outputs": [],
   "source": [
    "ds = nlmod.get_ds(extent, delr=delr, delc=delc)\n",
    "ds = nlmod.grid.refine(ds, model_ws=\"01_layer_generation\", refinement_features=refinement_features)\n",
    "ds = nlmod.grid.update_ds_from_layer_ds(ds, regis, method=\"average\")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "98b6d820",
   "metadata": {},
   "source": [
    "When we plot the model top, we now see that in the cells with an equal resolution to that of Regis (or smaller) no information is lost."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "f65ea923",
   "metadata": {},
   "outputs": [],
   "source": [
    "f, ax = nlmod.plot.get_map(extent, figsize=5)\n",
    "nlmod.plot.data_array(ds[\"top\"], ds=ds, edgecolor=\"k\");"
   ]
  }
 ],
 "metadata": {
  "language_info": {
   "name": "python"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}