{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# 07 - Making and using hostart files"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Create a hotstart file to resume a simulation from given hydrological conditions\n",
    "\n",
    "Hydrological models have state variables that describe the snow pack, soil moisture, underground reservoirs, etc. Typically, those cannot be measured empirically, so one way to estimate those values is to run the model for a period before the period we are actually interested in, and save the state variables at the end of this *warm-up* simulation.\n",
    "\n",
    "This notebook shows how to save those state variables and use them to configure another Raven simulation. These *states* are configured by the `:HRUStateVariableTable` and `:BasinStateVariables` commands, but `ravenpy` has a convenience function `set_solution` to update those directly from the `solution.rvc` simulation output.\n",
    "\n",
    "In the following, we run the model on two years then save the final states. Next, we use those final states to configure the initial state of a second simulation over the next two years. If everything is done correctly, these two series should be identical to a simulation over the full four years.\n",
    "\n",
    "## Model configuration\n",
    "\n",
    "At this point the following blocks of code should be quite familiar! If not, please go back to notebook \"04 - Emulating hydrological models\" to understand what is happening.\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Import packages\n",
    "import datetime as dt\n",
    "import warnings\n",
    "\n",
    "from matplotlib import pyplot as plt\n",
    "\n",
    "from ravenpy import Emulator, RavenWarning\n",
    "\n",
    "# Import the GR4JCN model\n",
    "from ravenpy.config import commands as rc\n",
    "from ravenpy.config import emulators\n",
    "\n",
    "# Utility that simplifies working with test data hosted on GitHub\n",
    "from ravenpy.testing.utils import get_file"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Start and end date for full simulation\n",
    "# Make sure the end date is before the end of the hydrometeorological data NetCDF file.\n",
    "start_date = dt.datetime(1986, 1, 1)\n",
    "end_date = dt.datetime(1988, 1, 1)\n",
    "\n",
    "# Define HRU\n",
    "hru = dict(\n",
    "    area=4250.6,\n",
    "    elevation=843.0,\n",
    "    latitude=54.4848,\n",
    "    longitude=-123.3659,\n",
    "    hru_type=\"land\",\n",
    ")\n",
    "\n",
    "# Get dataset:\n",
    "ERA5_full = get_file(\"notebook_inputs/ERA5_weather_data.nc\")\n",
    "\n",
    "# Set alternative names for netCDF variables\n",
    "alt_names = {\n",
    "    \"TEMP_MIN\": \"tmin\",\n",
    "    \"TEMP_MAX\": \"tmax\",\n",
    "    \"PRECIP\": \"pr\",\n",
    "}\n",
    "\n",
    "# Data types to extract from netCDF\n",
    "data_type = [\"TEMP_MAX\", \"TEMP_MIN\", \"PRECIP\"]\n",
    "data_kwds = {\n",
    "    \"ALL\": {\n",
    "        \"elevation\": hru[\"elevation\"],\n",
    "        \"latitude\": hru[\"latitude\"],\n",
    "        \"longitude\": hru[\"longitude\"],\n",
    "    }\n",
    "}\n",
    "\n",
    "# Model configuration\n",
    "config = emulators.GR4JCN(\n",
    "    params=[0.529, -3.396, 407.29, 1.072, 16.9, 0.947],\n",
    "    Gauge=[\n",
    "        rc.Gauge.from_nc(\n",
    "            ERA5_full,\n",
    "            data_type=data_type,\n",
    "            alt_names=alt_names,\n",
    "            data_kwds=data_kwds,\n",
    "        )\n",
    "    ],\n",
    "    HRUs=[hru],\n",
    "    StartDate=start_date,\n",
    "    EndDate=end_date,\n",
    "    RunName=\"full\",\n",
    ")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Silence the Raven warnings\n",
    "warnings.simplefilter(\"ignore\", category=RavenWarning)\n",
    "\n",
    "# Run the model and get the outputs.\n",
    "out1 = Emulator(config=config).run()\n",
    "\n",
    "# Plot the model output\n",
    "out1.hydrograph.q_sim.plot(label=\"Part 1\")\n",
    "plt.legend()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Now let's run the model for the next two years, setting the initial conditions to the final states of the first simulation.\n",
    "\n",
    "The path to the `solution.rvc` file can be found in `out1.files[\"solution\"]`.\n",
    "The content itself can be displayed with `out1.solution`"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# The path to the solution (final model states)\n",
    "hotstart = out1.files[\"solution\"]\n",
    "\n",
    "# Configure and run the model, this time with the next two years first 3 years (1988-1990).\n",
    "conf2 = config.set_solution(hotstart)\n",
    "conf2.start_date = dt.datetime(1988, 1, 1)\n",
    "conf2.end_date = dt.datetime(1990, 1, 1)\n",
    "conf2.run_name = \"part_2\"\n",
    "\n",
    "out2 = Emulator(config=conf2).run()\n",
    "\n",
    "# Plot the model output\n",
    "out2.hydrograph.q_sim.plot(label=\"Part 2\")\n",
    "plt.legend()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Compare with simulation over entire period\n",
    "\n",
    "Now in theory, those two simulations should be identical to one simulation over the whole period of four years, let's confirm this."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "full = config.copy()\n",
    "full.end_date = dt.datetime(1990, 1, 1)\n",
    "full.run_name = \"full\"\n",
    "\n",
    "out = Emulator(config=full).run(overwrite=True)\n",
    "\n",
    "out.hydrograph.q_sim.plot(label=\"Full\", color=\"gray\", lw=4)\n",
    "out1.hydrograph.q_sim.plot(\n",
    "    label=\"Part 1\",\n",
    "    color=\"blue\",\n",
    "    lw=0.5,\n",
    ")\n",
    "out2.hydrograph.q_sim.plot(label=\"Part 2\", color=\"orange\", lw=0.5)\n",
    "plt.legend()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "And now if we look at the difference between both hydrographs, we can see that there are differences in the second part at machine precision levels, due to rounding in the hotstart file (note that the y-axis is 1e-6, which is essentially 0!). But the rest is perfect!\n",
    "\n",
    "Therefore, we can provide forecasting abilities by saving simulation final states and using those to initialize model states for the forecasting runs. This will be used in other notebooks such as notebook #12 on hindcasting."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "\n",
    "delta1 = np.abs(out1.hydrograph.q_sim - out.hydrograph.q_sim)\n",
    "delta2 = np.abs(out2.hydrograph.q_sim - out.hydrograph.q_sim)\n",
    "\n",
    "delta1.plot(label=\"Part 1\")\n",
    "delta2.plot(label=\"Part 2\")\n",
    "plt.title(\"Difference between two parts and full simulation\")"
   ]
  }
 ],
 "metadata": {
  "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.10.9"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}