{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "52f09a4b",
   "metadata": {},
   "source": [
    "# Spatially allocating a `ws3` harvest schedule to raster grids\n",
    "\n",
    "This tutorial demonstrates how to take an aspatial harvest schedule produced by `ws3`, disaggregate it to annual steps, and allocate it across a rasterized forest inventory using the [`ws3.spatial.ForestRaster`](https://ws3.readthedocs.io/en/latest/api.html#ws3.spatial.ForestRaster) class. The workflow mirrors the spatial allocation utilities used by `spades_ws3`, but keeps things self-contained inside this repository.\n",
    "\n",
    "We will:\n",
    "\n",
    "- build a `ForestModel` from the sample Woodstock-format files\n",
    "- generate a quick heuristic harvest schedule\n",
    "- rasterize the sample stand inventory shapefile to a three-band GeoTIFF (hash, age, block id)\n",
    "- use `ForestRaster` to convert the 10-year schedule into annual GeoTIFF disturbance maps\n",
    "- inspect the resulting rasters and summarise disturbed area per year\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "5071bdac",
   "metadata": {},
   "source": [
    "## Recommended environment\n",
    "\n",
    "As with the other notebooks in this folder, running inside a sandboxed virtual environment keeps dependencies tidy (see `000_venv_python_kernel_setup.ipynb`). The commands below assume this notebook lives inside the WS3 repository checkout."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "id": "88bf14a0",
   "metadata": {},
   "outputs": [],
   "source": [
    "%load_ext autoreload\n",
    "%autoreload 2"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "id": "1a51f590",
   "metadata": {},
   "outputs": [],
   "source": [
    "clobber_ws3 = False\n",
    "if clobber_ws3:\n",
    "    %pip uninstall -y ws3\n",
    "    %pip install -e ..\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "id": "3c3bead8",
   "metadata": {},
   "outputs": [],
   "source": [
    "import sys\n",
    "from pathlib import Path\n",
    "import numpy as np\n",
    "import pandas as pd\n",
    "import rasterio\n",
    "import matplotlib.pyplot as plt\n",
    "\n",
    "# allow notebooks to reuse helper utilities in this folder\n",
    "EXAMPLES_DIR = Path.cwd() / 'examples'\n",
    "if str(EXAMPLES_DIR) not in sys.path:\n",
    "    sys.path.insert(0, str(EXAMPLES_DIR))\n",
    "\n",
    "import ws3.forest\n",
    "from ws3.common import rasterize_stands, hash_dt\n",
    "from ws3.spatial import ForestRaster\n",
    "from util import schedule_harvest_areacontrol, compile_scenario, plot_scenario\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "5936bb22",
   "metadata": {},
   "source": [
    "## Data paths\n",
    "\n",
    "We re-use the sample TSA 24 dataset bundled with the repository. The Woodstock formatted files feed the forest model, while the companion shapefile provides a spatial inventory for rasterization."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "id": "dacaded0",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(PosixPath('data/raster_spatial_demo/tsa24_inventory.tif'),\n",
       " PosixPath('data/spatial_outputs/tsa24_demo'))"
      ]
     },
     "execution_count": 4,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "#DATA_DIR = Path('examples/data')\n",
    "DATA_DIR = Path('data')\n",
    "WOODSTOCK_DIR = DATA_DIR / 'woodstock_model_files_tsa24_clipped'\n",
    "SHAPE_PATH = DATA_DIR / 'shp' / 'tsa24_clipped.shp' / 'stands.shp'\n",
    "\n",
    "# output directories (created on demand)\n",
    "RASTER_DIR = DATA_DIR / 'raster_spatial_demo'\n",
    "RASTER_DIR.mkdir(parents=True, exist_ok=True)\n",
    "INVENTORY_TIF = RASTER_DIR / 'tsa24_inventory.tif'\n",
    "\n",
    "SPATIAL_DIR = DATA_DIR / 'spatial_outputs' / 'tsa24_demo'\n",
    "SPATIAL_DIR.mkdir(parents=True, exist_ok=True)\n",
    "\n",
    "INVENTORY_TIF, SPATIAL_DIR\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "23a00d32",
   "metadata": {},
   "source": [
    "## Build a `ForestModel`\n",
    "\n",
    "We load the Woodstock-format inputs, initialise state, and keep the horizon short (5 periods of 10 years) so that the spatial allocation step produces a manageable number of rasters."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "id": "afd381d0",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<ws3.forest.ForestModel at 0x7a02b8306de0>"
      ]
     },
     "execution_count": 5,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "BASE_YEAR = 2020\n",
    "HORIZON = 5\n",
    "PERIOD_LENGTH = 10\n",
    "\n",
    "fm = ws3.forest.ForestModel(\n",
    "    model_name='tsa24_clipped',\n",
    "    model_path=str(WOODSTOCK_DIR),\n",
    "    base_year=BASE_YEAR,\n",
    "    horizon=HORIZON,\n",
    "    period_length=PERIOD_LENGTH,\n",
    "    max_age=1000,\n",
    ")\n",
    "\n",
    "fm.import_landscape_section()\n",
    "fm.import_areas_section(convert_periods_to_years=period_length)\n",
    "fm.import_yields_section(convert_periods_to_years=period_length)\n",
    "fm.import_actions_section(convert_periods_to_years=period_length)\n",
    "fm.import_transitions_section(convert_periods_to_years=period_length)\n",
    "fm.initialize_areas()\n",
    "fm.add_null_action()\n",
    "fm.reset_actions()\n",
    "\n",
    "fm\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "9cea2bcd",
   "metadata": {},
   "source": [
    "## Heuristic harvest schedule\n",
    "\n",
    "We call the `schedule_harvest_areacontrol` helper (defined in `examples/util.py`) to generate an aspatial schedule. The function compiles the internal `applied_actions` structure, which `ForestRaster` will later consume."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "id": "1d5d9360",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>period</th>\n",
       "      <th>oha</th>\n",
       "      <th>ohv</th>\n",
       "      <th>ogs</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>1</td>\n",
       "      <td>103.411088</td>\n",
       "      <td>14069.080670</td>\n",
       "      <td>134541.524980</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>2</td>\n",
       "      <td>105.655079</td>\n",
       "      <td>14690.387934</td>\n",
       "      <td>130879.310887</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>3</td>\n",
       "      <td>103.411088</td>\n",
       "      <td>14352.583697</td>\n",
       "      <td>126127.588834</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>4</td>\n",
       "      <td>103.411088</td>\n",
       "      <td>14403.805319</td>\n",
       "      <td>121730.561271</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <td>5</td>\n",
       "      <td>103.411088</td>\n",
       "      <td>14167.526587</td>\n",
       "      <td>119035.437396</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "   period         oha           ohv            ogs\n",
       "0       1  103.411088  14069.080670  134541.524980\n",
       "1       2  105.655079  14690.387934  130879.310887\n",
       "2       3  103.411088  14352.583697  126127.588834\n",
       "3       4  103.411088  14403.805319  121730.561271\n",
       "4       5  103.411088  14167.526587  119035.437396"
      ]
     },
     "execution_count": 6,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "sch = schedule_harvest_areacontrol(fm)\n",
    "scenario = compile_scenario(fm)\n",
    "scenario\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "id": "577074f3",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAA9oAAAGHCAYAAABPmCpHAAAAOnRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjEwLjYsIGh0dHBzOi8vbWF0cGxvdGxpYi5vcmcvq6yFwwAAAAlwSFlzAAAPYQAAD2EBqD+naQAAi7dJREFUeJzs3XdYFNf7NvAbhAWkLKICoojEELFgCRjEXghgR7GgJDYsUbDGmihiC7GLJaImliQQlUSJGoNiSUApCootippgiX6BGAQUFRDO+4fvzs+RpaiLaLw/1zXX5Z7z7MyZWTjuw5k5R0sIIUBEREREREREGqFd2Q0gIiIiIiIi+i9hok1ERERERESkQUy0iYiIiIiIiDSIiTYRERERERGRBjHRJiIiIiIiItIgJtpEREREREREGsREm4iIiIiIiEiDmGgTERERERERaRATbSIiIiIiIiINYqJNRG+F3377DVpaWvjtt980sr+tW7dCS0sL165de+l9BQYGQktLC3fu3Hn5hv1HXLlyBW5ublAqldDS0kJERIRGr/mrcPPmTejr6+P48eOV3ZTXlpaWFgIDAyu7GW8VXnOgY8eOaNKkSakxBQUFsLa2xldfffWKWkVE/zVMtIlII1RJUGJiotr68nyxeR198cUXiIiIqOxmvNEePHiAwMDA5/ojx9ChQ3Hu3DksWrQI3333HZycnCqugRVk/vz5cHZ2Rps2bSq7KURvhNjYWAQGBiIrK6uymwJdXV1MmTIFixYtwqNHjyq7OUT0BmKiTURvhfbt2+Phw4do3779c72vpET7448/xsOHD2FjY6OhFv53PXjwAPPmzSt3ov3w4UPExcXB19cX/v7++Oijj1CnTp2KbaSG/fPPP9i2bRs++eSTym4KkczDhw8xe/bsym6GWrGxsZg3b95rkWgDwPDhw3Hnzh2EhYVVdlOI6A3ERJuI/tMePXqEoqIiaGtrQ19fH9ramun2qlSpAn19fWhpaWlkfxVNCIGHDx9WdjPK5Z9//gEAmJqaVm5DXsL3338PHR0d9OzZU2P7zM3N1di+6PXzqj5ffX196OjovJJjvelMTU3h5uaGrVu3VnZTiOgNxESbiCrV999/D0dHRxgYGMDMzAze3t64efOmLKZevXoYNmxYsfd27NgRHTt2lF6rnsPevn07Zs+ejdq1a6Nq1arIyclR+4z2lStX4OXlBUtLS+jr66NOnTrw9vZGdnY2gCfPMubm5mLbtm3Q0tKClpaW1I6Snhf+9ddf0aFDBxgbG8PExAQtW7Ys92hIVlYWhg0bBlNTUyiVSgwfPhwPHjyQxWzZsgWdO3eGubk59PT00KhRI6xfv77YvurVq4cePXrgwIEDcHJygoGBATZs2IAmTZqgU6dOxeKLiopQu3Zt9OvXT1a2atUqNG7cGPr6+rCwsMCYMWNw9+5d2XsTExPh7u6OGjVqwMDAALa2thgxYgQA4Nq1a6hZsyYAYN68edJ1LOkZ0cDAQOkugWnTpkFLSwv16tUr9bp99dVXaNy4MfT09GBlZQU/Pz/ZiNjq1atRpUoVWdny5cuhpaWFKVOmSGWFhYUwNjbGjBkzpLLt27fD0dFR+jwdHBwQHBxcansAICIiAs7OzjAyMpKVx8TEoH///qhbty709PRgbW2NyZMnF/sjyLBhw2BkZIQ///wT3bp1g7GxMXx8fACU/3P5+eef0b17d1hZWUFPTw/169fHggULUFhYWGb7n1VQUAAzMzMMHz68WF1OTg709fUxdepUqSwjIwO+vr6wsLCAvr4+mjVrhm3btpV5nGHDhqn9vFXzGDxNS0sL/v7+CA8PR6NGjWBgYAAXFxecO3cOALBhwwa8++670NfXR8eOHdU+25+QkAAPDw8olUpUrVoVHTp0eKln6ouKihAYGAgrKytUrVoVnTp1wh9//FGsD1P1H7///jvGjRsHc3Nz2V0bFfkz/ezvn+raXr16tcz+5+HDh5gwYQJq1KgBY2Nj9OrVC7du3Sr3c99r1qxB48aNUbVqVVSrVg1OTk5S/xgYGIhp06YBAGxtbaW+QvW5PX78GAsWLED9+vWhp6eHevXq4bPPPkNeXl6x47xIP3zw4EFUrVoVgwYNwuPHj6XyDz/8EMeOHUNmZmaZ50dE9DT+SZOINCo7O1vtpF4FBQXFyhYtWoQ5c+ZgwIABGDlyJP755x+sWbMG7du3x+nTp194RHPBggVQKBSYOnUq8vLyoFAoisXk5+fD3d0deXl5GD9+PCwtLXHr1i3s27cPWVlZUCqV+O677zBy5Eh88MEHGD16NACgfv36JR5369atGDFiBBo3boxZs2bB1NQUp0+fRmRkJAYPHlxmuwcMGABbW1sEBQXh1KlT+Prrr2Fubo7FixdLMevXr0fjxo3Rq1cv6OjoYO/evRg3bhyKiorg5+cn219KSgoGDRqEMWPGYNSoUWjQoAEGDhyIwMBApKWlwdLSUoo9duwYbt++DW9vb6lszJgx2Lp1K4YPH44JEyYgNTUVa9euxenTp3H8+HHo6uoiIyMDbm5uqFmzJmbOnAlTU1Ncu3YNu3btAgDUrFkT69evx9ixY9GnTx/07dsXANC0aVO116Bv374wNTXF5MmTMWjQIHTr1q1Ysvq0wMBAzJs3D66urhg7dixSUlKwfv16nDx5Umpju3btUFRUhGPHjqFHjx4AniS82traiImJkfZ1+vRp3L9/X3q8ICoqCoMGDUKXLl2kz+DixYs4fvw4Jk6cWGKbCgoKcPLkSYwdO7ZYXXh4OB48eICxY8eievXqOHHiBNasWYO///4b4eHhstjHjx/D3d0dbdu2xbJly1C1atVyfy7Ak59HIyMjTJkyBUZGRjhy5AgCAgKQk5ODpUuXlth+dXR1ddGnTx/s2rULGzZskP1ORUREIC8vT/rZefjwITp27IirV6/C398ftra2CA8Px7Bhw5CVlVXqtXteMTEx2LNnj/SzHxQUhB49emD69On46quvMG7cONy9exdLlizBiBEjcOTIEem9R44cQdeuXeHo6Ii5c+dCW1tb+kNWTEwMPvjgg+duz6xZs7BkyRL07NkT7u7uOHPmDNzd3Ut8xnfcuHGoWbMmAgICpBHtivyZLk15+p9hw4Zh586d+Pjjj9GqVSv8/vvv6N69e7muzaZNmzBhwgT069cPEydOxKNHj3D27FkkJCRg8ODB6Nu3Ly5fvowffvgBK1euRI0aNQBA+kPdyJEjsW3bNvTr1w+ffvopEhISEBQUhIsXL2L37t3ScV6kH963bx/69euHgQMHYvPmzahSpYpU5+joCCEEYmNjpWtNRFQugohIA7Zs2SIAlLo1btxYir927ZqoUqWKWLRokWw/586dEzo6OrJyGxsbMXTo0GLH7NChg+jQoYP0+ujRowKAeOedd8SDBw9ksaq6o0ePCiGEOH36tAAgwsPDSz0vQ0NDtcdWnW9qaqoQQoisrCxhbGwsnJ2dxcOHD2WxRUVFpR5j7ty5AoAYMWKErLxPnz6ievXqsrJnz0sIIdzd3cU777wjK7OxsREARGRkpKw8JSVFABBr1qyRlY8bN04YGRlJ+4+JiREARGhoqCwuMjJSVr57924BQJw8ebLE8/vnn38EADF37twSY56WmpoqAIilS5fKyp+95hkZGUKhUAg3NzdRWFgoxa1du1YAEJs3bxZCCFFYWChMTEzE9OnThRBPPo/q1auL/v37iypVqoh79+4JIYRYsWKF0NbWFnfv3hVCCDFx4kRhYmIiHj9+XK52q1y9elXtNRZC/ecXFBQktLS0xPXr16WyoUOHCgBi5syZstjyfi4lHWvMmDGiatWq4tGjR891TkIIceDAAQFA7N27V1berVs32c/fqlWrBADx/fffS2X5+fnCxcVFGBkZiZycHKn82Z+LoUOHChsbm2LHVv2OPA2A0NPTk34ehBBiw4YNAoCwtLSUHWfWrFmyn52ioiJhZ2cn3N3dZb+fDx48ELa2tuLDDz8s1zV5WlpamtDR0RGenp6y8sDAQAFA1o+ofpbbtm0r+/mq6J9p1XV7+pqXt/9JSkoSAMSkSZNkccOGDSvX73fv3r1l/weos3TpUtnnpJKcnCwAiJEjR8rKp06dKgCII0eOCCHK3w936NBBastPP/0kdHV1xahRo2TXXOX27dsCgFi8eHGpbSciehZvHScijVq3bh2ioqKKbc+OYO7atQtFRUUYMGAA7ty5I22Wlpaws7PD0aNHX7gNQ4cOhYGBQakxSqUSAHDgwIFit0e+iKioKNy7dw8zZ86Evr6+rK68z3E/O3FWu3bt8O+//yInJ0cqe/q8VHcPdOjQAX/99Zd0y7uKra0t3N3dZWXvvfcemjdvjh07dkhlhYWF+PHHH9GzZ09p/+Hh4VAqlfjwww9ln4+joyOMjIykz0d118G+ffvU3rVQkQ4dOoT8/HxMmjRJ9uz9qFGjYGJigl9++QUAoK2tjdatWyM6OhrAk1Hpf//9FzNnzoQQAnFxcQCejAg2adJEOidTU1Pk5uYiKirqudr177//AgCqVatWrO7pzy83Nxd37txB69atIYTA6dOni8U/Oype3s/l2WPdu3cPd+7cQbt27fDgwQNcunTpuc4JADp37owaNWrIfnbu3r2LqKgoDBw4UCrbv38/LC0tMWjQIKlMV1cXEyZMwP379/H7778/97FL0qVLF9mt5s7OzgAALy8vGBsbFyv/66+/AADJycm4cuUKBg8ejH///Ve6jrm5uejSpQuio6NRVFT0XG05fPgwHj9+jHHjxsnKx48fX+J7Ro0aJRs9reif6dKU1f9ERkYCwHOd39NMTU3x999/4+TJk+WKf9r+/fsBQHZbPAB8+umnACBdl+fth3/44QcMHDgQY8aMwYYNG9TO4aH6Pebyi0T0vJhoE5FGffDBB3B1dS22PZt0XLlyBUII2NnZoWbNmrLt4sWLyMjIeOE22NralitmypQp+Prrr1GjRg24u7tj3bp1xZLV8vrzzz8B4KWWMKtbt67steqaPf3s7fHjx+Hq6gpDQ0OYmpqiZs2a+OyzzwBAbaKtzsCBA3H8+HHcunULwJNn2zMyMmTJ0pUrV5CdnQ1zc/Nin8/9+/elz6dDhw7w8vLCvHnzUKNGDfTu3RtbtmxR+9ykpl2/fh0A0KBBA1m5QqHAO++8I9UDT5KGpKQkPHz4EDExMahVqxbef/99NGvWTLrV9tixY2jXrp30nnHjxuG9995D165dUadOHYwYMUJKNspDCFGs7MaNGxg2bBjMzMxgZGSEmjVrokOHDgCKf346OjrFZlsv7+cCABcuXECfPn2gVCphYmKCmjVr4qOPPlJ7rPLQ0dGBl5cXfv75Z+nz3bVrFwoKCmQ/O9evX4ednV2xpKVhw4ZSvaY8+zuj+gOatbW12nLV79KVK1cAPPmj3LPX8euvv0ZeXt5zXyPVeb377ruycjMzM7V/dAGK/45W9M90acrqf65fvw5tbe1ibX72fEsyY8YMGBkZ4YMPPoCdnR38/PzK/Ty86tjPHsvS0hKmpqbSdXmefjg1NRUfffQRvLy8sGbNmhL/IKr6PX5TJr4kotcHn9EmokpRVFQELS0t/Prrr7IRHZWnn8st6QtOYWGh2veWNZqtsnz5cgwbNgw///wzDh48iAkTJiAoKAjx8fGVspyUunMB/u+L3p9//okuXbrA3t4eK1asgLW1NRQKBfbv34+VK1cWG4Er6ToMHDgQs2bNQnh4OCZNmoSdO3dCqVTCw8NDiikqKoK5uTlCQ0PV7kP13KSWlhZ+/PFHxMfHY+/evThw4ABGjBiB5cuXIz4+vtTnq1+ltm3boqCgAHFxcYiJiZGSj3bt2iEmJgaXLl3CP//8I0tKzM3NkZycjAMHDuDXX3/Fr7/+ii1btmDIkCGlTuxVvXp1ACg2OVlhYSE+/PBDZGZmYsaMGbC3t4ehoSFu3bqFYcOGFfv89PT0iiWr5f1csrKy0KFDB5iYmGD+/PmoX78+9PX1cerUKcyYMeO5R2tVvL29sWHDBvz666/w9PTEzp07YW9vj2bNmr3Q/p5V2u+6OiX9zpT1u6Q6/6VLl6J58+ZqY1/Fz255+yp1XuRnujRlXbOX1bBhQ6SkpGDfvn2IjIzETz/9hK+++goBAQGYN29eufahyWS3Vq1aqFWrFvbv34/ExEQ4OTmpjVP9HqueGSciKi8m2kRUKerXrw8hBGxtbfHee++VGlutWjW166pev34d77zzzku1w8HBAQ4ODpg9ezZiY2PRpk0bhISEYOHChQDK/8VONUna+fPnyz3C87z27t2LvLw87NmzRzb69Ly32dva2uKDDz7Ajh074O/vj127dsHT0xN6enpSTP369XHo0CG0adOmXMlAq1at0KpVKyxatAhhYWHw8fHB9u3bMXLkyAobCVLNTp6SkiL7OcjPz0dqaipcXV2lsg8++AAKhQIxMTGIiYmRZjdu3749Nm3ahMOHD0uvn6ZQKNCzZ0/07NkTRUVFGDduHDZs2IA5c+aU+DnXrVsXBgYGSE1NlZWfO3cOly9fxrZt2zBkyBCp/HluTS/v5/Lbb7/h33//xa5du2Tn9Gybnlf79u1Rq1Yt7NixA23btsWRI0fw+eefy2JsbGxw9uxZaVk9FdXt6qWtPV/a77omqX5fTUxMZD8nL0N1XlevXpWN+v7777/F/uhS1j4q8mf6RdnY2KCoqAipqamws7OTyq9evVrufRgaGmLgwIEYOHAg8vPz0bdvXyxatAizZs0qdblE1bGvXLki3RkBAOnp6cjKypKu2/P0w/r6+ti3bx86d+4MDw8P/P7772jcuHGxONXvzNPHJSIqD946TkSVom/fvqhSpQrmzZtXbMRECCE95wo8+fIUHx+P/Px8qWzfvn3FlgF7Hjk5ObIlXIAnSbe2trbstmdDQ0O1X/yf5ebmBmNjYwQFBRWbYVhTI0KqEaen95ednY0tW7Y8974GDhyI+Ph4bN68GXfu3JHd+gs8mYG4sLAQCxYsKPbex48fS9fk7t27xc5PNUKouo6q2bLLcx2fh6urKxQKBVavXi1rwzfffIPs7GzZbMj6+vpo2bIlfvjhB9y4cUM2+vfw4UOsXr0a9evXR61ataT3PP0zCDx5LlY110Bpt8br6urCyckJiYmJsnJ1n58QolzLhamU93NRd6z8/Hx89dVX5T6WOtra2ujXrx/27t2L7777Do8fPy72s9OtWzekpaXJnuV+/Pgx1qxZAyMjI+lWeXXq16+P7OxsnD17Vir73//+J5tVWhMcHR1Rv359LFu2DPfv3y9Wr1rL/Xl06dIFOjo6xZbbW7t2bbn3UdE/0y9DNd/Dsz9Da9asKdf7n/19UigUaNSoEYQQ0vwOhoaGAIr3Fd26dQMArFq1Sla+YsUKAJCuy/P2w0qlEgcOHIC5uTk+/PBD6dbzpyUlJUFLSwsuLi7lOk8iIhWOaBNRpahfvz4WLlyIWbNm4dq1a/D09ISxsTFSU1Oxe/dujB49WlqXd+TIkfjxxx/h4eGBAQMG4M8//8T3339f6lJbZTly5Aj8/f3Rv39/vPfee3j8+DG+++47VKlSBV5eXlKco6MjDh06hBUrVsDKygq2trbSxEpPMzExwcqVKzFy5Ei0bNkSgwcPRrVq1XDmzBk8ePCgXGsIl8XNzU0aYR0zZgzu37+PTZs2wdzcHP/73/+ea18DBgzA1KlTMXXqVJiZmRUb1evQoQPGjBmDoKAgJCcnw83NDbq6urhy5QrCw8MRHByMfv36Ydu2bfjqq6/Qp08f1K9fH/fu3cOmTZtgYmIifTk2MDBAo0aNsGPHDrz33nswMzNDkyZNXup5duDJbdKzZs3CvHnz4OHhgV69eiElJQVfffUVWrZsKT2PrNKuXTt8+eWXUCqVcHBwAPDk9vAGDRogJSWl2FrtI0eORGZmJjp37ow6derg+vXrWLNmDZo3b17m6Fbv3r3x+eefIycnByYmJgAAe3t71K9fH1OnTsWtW7dgYmKCn376qdyjnUD5P5fWrVujWrVqGDp0KCZMmAAtLS189913apON3377DZ06dcLcuXPLtRbywIEDsWbNGsydOxcODg7FrsXo0aOxYcMGDBs2DElJSahXrx5+/PFHHD9+HKtWrZJNUvYsb29vzJgxA3369MGECRPw4MEDrF+/Hu+99x5OnTpV7utUFm1tbXz99dfo2rUrGjdujOHDh6N27dq4desWjh49ChMTE+zdu1eK19LSQocOHfDbb7+VuE8LCwtMnDgRy5cvR69eveDh4YEzZ87g119/RY0aNcp1Z0dF/0y/DEdHR3h5eWHVqlX4999/peW9Ll++DKDsu3/c3NxgaWmJNm3awMLCAhcvXsTatWvRvXt36WfC0dERAPD555/D29sburq66NmzJ5o1a4ahQ4di48aN0mMRJ06cwLZt2+Dp6YlOnToBeLF+uEaNGoiKikLbtm3h6uqKY8eOoXbt2lJ9VFQU2rRpIz0SQkRUbq9whnMi+g9TLVdT0jJPTy+n8rSffvpJtG3bVhgaGgpDQ0Nhb28v/Pz8REpKiixu+fLlonbt2kJPT0+0adNGJCYmlri8l7olu55d3uuvv/4SI0aMEPXr1xf6+vrCzMxMdOrUSRw6dEj2vkuXLon27dsLAwMD2RI9zy41pbJnzx7RunVrYWBgIExMTMQHH3wgfvjhh1KvnWp5nX/++UdWru4Ye/bsEU2bNhX6+vqiXr16YvHixWLz5s3F4mxsbET37t1LPW6bNm3ULpnztI0bNwpHR0dhYGAgjI2NhYODg5g+fbq4ffu2EEKIU6dOiUGDBom6desKPT09YW5uLnr06CESExNl+4mNjRWOjo5CoVCUuRRQeZf3Ulm7dq2wt7cXurq6wsLCQowdO1a2nJHKL7/8IgCIrl27yspHjhwpAIhvvvlGVv7jjz8KNzc3YW5uLhQKhahbt64YM2aM+N///ldi21XS09OFjo6O+O6772Tlf/zxh3B1dRVGRkaiRo0aYtSoUeLMmTMCgNiyZYsUN3ToUGFoaFji/sv6XIQQ4vjx46JVq1bCwMBAWFlZienTp0tLdKl+D4QQYu/evQKACAkJKfO8hHiyTJK1tbUAIBYuXFji+Q8fPlzUqFFDKBQK4eDgIDs/FXU/CwcPHhRNmjQRCoVCNGjQQHz//fclLu/l5+cnKyvpZ6ekvuH06dOib9++onr16kJPT0/Y2NiIAQMGiMOHD0sx9+7dEwCEt7d3WZdGPH78WMyZM0dYWloKAwMD0blzZ3Hx4kVRvXp18cknn0hxZfWXFfUzLUTJy3uVp//Jzc0Vfn5+wszMTBgZGQlPT09pycAvv/yy1GuzYcMG0b59e+la169fX0ybNk1kZ2fL4hYsWCBq164ttLW1ZccvKCgQ8+bNE7a2tkJXV1dYW1uLWbNmqV2qrqx+WN3/R1evXhW1atUSDRs2lK5FVlaWUCgU4uuvvy713IiI1NESQkP3NBIREZHE19cXly9flmaAfl1Nnz4dP/zwA65evSp7Tp+e2L9/P3r06IEzZ85Io8bPIysrC9WqVcPChQuLPc/+X5CcnIwWLVrg+++/h4+PT2U3R6NWrVqFJUuW4M8//3ypieuI6O3EZ7SJiIgqwNy5c3Hy5MlyL2FUWY4ePYo5c+YwyS7B0aNH4e3tXa4k++HDh8XKVM8Vd+zYUcMte/VKOj9tbW2NTbr2uigoKMCKFSswe/ZsJtlE9EI4ok1ERESkAVu3bsXWrVvRrVs3GBkZ4dixY/jhhx/g5uaGAwcOVHbzXtq8efOQlJSETp06QUdHR1r2TvVcPhER/R8m2kREREQacOrUKUyfPh3JycnIycmBhYUFvLy8sHDhwtdmTfmXERUVhXnz5uGPP/7A/fv3UbduXXz88cf4/PPPoaPD+XWJiJ7GW8fpP2nr1q3Q0tLCtWvXKrspGrVz506YmZnJlqPR0tKCv7+/xo4RGRkJIyOjF1rehogqXmX1b4GBgRW2Jvqr0q1bN4waNarC9v/+++/j0KFDuHPnDvLz83Hz5k2sWrUKRkZGCAkJQd26dUtdGu519+GHH+LYsWPIzMxEfn4+rl69irlz5zLJphfy22+/QUtLq9TZ/N8Wqmvx448/vvA+Tpw4AYVCgevXr2uwZeVTUFAAa2vrl15C8r+GifZrSPUl6tk1WFU6duz40svivA4ePHiAwMBAdrDlVFhYiLlz52L8+PEVOjLi4eGBd999F0FBQRV2DHp7sX+jynL8+HEcPHgQM2bM0Oh+J0+ejPfffx9mZmaoWrUqGjZsiMDAwGLrcw8bNgz5+fm8xZpeidTUVPj7++O9995D1apVUbVqVTRq1Ah+fn6yderfZufOnUO/fv1gY2MDfX191K5dGx9++GGxteG/+OILREREVE4jn8Pnn3+OQYMGwcbGRmP7jI6ORq9evWBtbQ19fX1YWlrCw8Oj2Nwjurq6mDJlChYtWlRsDfu3GRNtqjQPHjzAvHnz+EW0nPbu3YuUlBSMHj26wo81ZswYbNiwAffu3avwYxH9F7F/e/0sXboUXbp0wbvvvqvR/Z48eRLt2rXDvHnzEBwcjE6dOuHLL7+Eh4cHioqKpDh9fX0MHToUK1asULueOZGm7Nu3D02aNMF3330HV1dXrFy5EsHBwejatSv279+P5s2bV8qo57Pat2+Phw8fVspEerGxsXBycsKZM2cwatQorF27FiNHjoS2tjaCg4NlsW9Cop2cnIxDhw7hk08+0eh+L1++DG1tbXzyySdYt24dpk6dirS0NLRv3x6RkZGy2OHDh+POnTsICwvTaBveZLzX5y0nhMCjR484o+YLyM3NhaGh4Ss73pYtW9CmTRvUrl27wo/l5eWF8ePHIzw8HCNGjKjw4xFVBPZvpJKRkYFffvkFISEhGt/3sWPHipXVr18fU6dOxYkTJ9CqVSupfMCAAViyZAmOHj2Kzp07a7wtRH/++Se8vb1hY2ODw4cPo1atWrL6xYsX46uvvoK2duljba/iO462tjb09fUr9BglWbRoEZRKJU6ePAlTU1NZXUZGRqW06WVs2bIFdevWlfU3mjBy5EiMHDlSVjZu3Di88847WLVqFTw8PKRyU1NTuLm5YevWrfzu+P9xRPs/YsuWLejcuTPMzc2hp6eHRo0aYf369cXi6tWrhx49euDAgQNwcnKCgYEBNmzYgCZNmqBTp07F4ouKilC7dm3069dPVrZq1So0btwY+vr6sLCwwJgxY3D37l3ZexMTE+Hu7o4aNWrAwMAAtra20i/etWvXULNmTQBPZjHV0tKClpYWAgMDpfdfunQJ/fr1g5mZGfT19eHk5IQ9e/YUa+OFCxfQuXNnGBgYoE6dOli4cKFsFKE0Z8+exbBhw/DOO+9It8SMGDEC//77ryxO9WziH3/8gcGDB6NatWpo27atVP/999/D0dERBgYGMDMzg7e3N27evCnbR0xMDPr374+6detCT08P1tbWmDx5strlUp716NEjREZGwtXVtcSYiIgINGnSBHp6emjcuHGxvzRev34d48aNQ4MGDWBgYIDq1aujf//+ap/zNDc3R9OmTfHzzz+X2Taiisb+7fn7t2XLlkFLS0vtqNWsWbOgUChk5xQeHi71YTVq1MBHH32EW7dulXqMa9euQUtLC1u3bi1W9+z5qvrQy5cv46OPPoJSqUTNmjUxZ84cCCFw8+ZN9O7dGyYmJrC0tMTy5cuL7TMvLw9z587Fu+++K/Wh06dPL9czz7/88gseP35crA9VPcpw7NgxTJgwATVr1oSpqSnGjBmD/Px8ZGVlYciQIahWrRqqVauG6dOnl2s0ul69egCerKH9NEdHR5iZmbFvpQqzZMkS5ObmYsuWLcWSbADQ0dHBhAkTYG1tLZUNGzYMRkZG+PPPP9GtWzcYGxtLa6Ln5ubi008/hbW1NfT09NCgQQMsW7ZM9nvQt29fvP/++7Lj9OzZE1paWrJ+LSEhAVpaWvj1118BqH9GW/X40B9//IFOnTqhatWqqF27NpYsWVLsXK5fv45evXrB0NAQ5ubmmDx5Mg4cOFCu577//PNPNG7cuFiSDTz5DqSipaWF3NxcbNu2TerLhw0bJtWfPn0aXbt2hYmJCYyMjNClSxfEx8cX22dWVhYmT56MevXqQU9PD3Xq1MGQIUNw586dEtuYl5eHHj16QKlUIjY2ttTziYiIQOfOnYvNo6H6f/G3336T/l90cHCQrs+uXbvg4OAAfX19ODo64vTp06UeBwCqVq2KmjVrFuvfAPk8DsQR7ddadna22l/AgoKCYmXr169H48aN0atXL+jo6GDv3r0YN24cioqK4OfnJ4tNSUnBoEGDMGbMGIwaNQoNGjTAwIEDERgYiLS0NFhaWkqxx44dw+3bt+Ht7S2VjRkzBlu3bsXw4cMxYcIEpKamYu3atTh9+jSOHz8OXV1dZGRkwM3NDTVr1sTMmTNhamqKa9euYdeuXQCAmjVrYv369Rg7diz69OmDvn37AgCaNm0K4MmXS9Xo7cyZM2FoaIidO3fC09MTP/30E/r06QMASEtLQ6dOnfD48WMpbuPGjeUewYqKisJff/2F4cOHw9LSEhcuXMDGjRtx4cIFxMfHF+uw+vfvDzs7O3zxxRfSfzKLFi3CnDlzMGDAAIwcORL//PMP1qxZg/bt2+P06dNSJx4eHo4HDx5g7NixqF69Ok6cOIE1a9bg77//Rnh4eKntTEpKQn5+frH/yJ7+nHbt2oVx48bB2NgYq1evhpeXF27cuIHq1asDeHJ7Y2xsLLy9vVGnTh1cu3YN69evR8eOHfHHH3+gatWqsn06Ojq+9rdK0ZuL/VvF9m8DBgzA9OnTsXPnTkybNk1Wt3PnTri5uaFatWoAIJ1vy5YtERQUhPT0dAQHB+P48eOyPkwTBg4ciIYNG+LLL7/EL7/8goULF8LMzAwbNmxA586dsXjxYoSGhmLq1Klo2bKldEtpUVERevXqhWPHjmH06NFo2LAhzp07h5UrV+Ly5ctl9lWxsbGoXr16ic8ujh8/HpaWlpg3bx7i4+OxceNGmJqaIjY2FnXr1sUXX3yB/fv3Y+nSpWjSpAmGDBkie//jx4+RlZWF/Px8nD9/HrNnz4axsTE++OCDYsd6//33X/u11enNtW/fPrz77rtwdnZ+rvc9fvwY7u7uaNu2LZYtW4aqVatCCIFevXrh6NGj8PX1RfPmzXHgwAFMmzYNt27dwsqVKwEA7dq1w88//4ycnByYmJhACIHjx49DW1sbMTEx6NWrF4AnAw7a2tpo06ZNqW25e/cuPDw80LdvXwwYMAA//vgjZsyYAQcHB3Tt2hXAkz8AdO7cGf/73/8wceJEWFpaIiwsDEePHi3X+drY2CAuLg7nz58vdV6Q7777DiNHjsQHH3wgPbpXv359AE/68nbt2sHExATTp0+Hrq4uNmzYgI4dO+L333+XPoP79++jXbt2uHjxIkaMGIH3338fd+7cwZ49e/D333+jRo0axY778OFD9O7dG4mJiTh06BBatmxZYhtv3bqFGzdulPgd8erVqxg8eDDGjBmDjz76CMuWLUPPnj0REhKCzz77DOPGjQMABAUFYcCAAUhJSSl2x0NOTg7y8/Nx584dfPvttzh//jw+++yzYsdydHSEEAKxsbHo0aNHiW1+awh67WzZskUAKHVr3Lix7D0PHjwoth93d3fxzjvvyMpsbGwEABEZGSkrT0lJEQDEmjVrZOXjxo0TRkZG0v5jYmIEABEaGiqLi4yMlJXv3r1bABAnT54s8Tz/+ecfAUDMnTu3WF2XLl2Eg4ODePTokVRWVFQkWrduLezs7KSySZMmCQAiISFBKsvIyBBKpVIAEKmpqSUeXwj11+2HH34QAER0dLRUNnfuXAFADBo0SBZ77do1UaVKFbFo0SJZ+blz54SOjo6sXN2xgoKChJaWlrh+/Xqp7fz6668FAHHu3LlidQCEQqEQV69elcrOnDlT7PNUd/y4uDgBQHz77bfF6r744gsBQKSnp5faNqLnwf7t1fVvLi4uwtHRUVZ24sQJ2e98fn6+MDc3F02aNBEPHz6U4vbt2ycAiICAAKlM1Q+qpKamCgBiy5YtxY797Lmr3jt69Gip7PHjx6JOnTpCS0tLfPnll1L53bt3hYGBgRg6dKhU9t133wltbW0RExMjO05ISIgAII4fP17qtWjbtm2xayHE//08uru7i6KiIqncxcVFaGlpiU8++aRYezt06FBsP6q+VLU1aNBAHD16VG1bRo8eLQwMDEptL9GLyM7OFgCEp6dnsbq7d++Kf/75R9qe7leHDh0qAIiZM2fK3hMRESEAiIULF8rK+/XrJ7S0tKTvHSdPnhQAxP79+4UQQpw9e1YAEP379xfOzs7S+3r16iVatGghvT569KgAIPtd6dChQ7HvJXl5ecLS0lJ4eXlJZcuXLxcAREREhFT28OFDYW9vX2yf6hw8eFBUqVJFVKlSRbi4uIjp06eLAwcOiPz8/GKxhoaGsv5IxdPTUygUCvHnn39KZbdv3xbGxsaiffv2UllAQIAAIHbt2lVsH6p+R3UtwsPDxb1790SHDh1EjRo1xOnTp0s9DyGEOHTokAAg9u7dW6xO9f9ibGysVHbgwAEBQBgYGMi+f27YsKHEa+fu7i71bwqFQowZM0b2f8bT5w9ALF68uMx2vw146/hrbN26dYiKiiq2qUZFnvb0CIdqpKhDhw7466+/kJ2dLYu1tbWFu7u7rOy9995D8+bNsWPHDqmssLAQP/74I3r27CntPzw8HEqlEh9++CHu3LkjbY6OjjAyMpL+kqgaAdm3b5/aEarSZGZm4siRIxgwYADu3bsnHePff/+Fu7s7rly5It3SuH//frRq1Uo2alCzZk3plqeyPH3dHj16hDt37kjPt5w6dapY/LOTTOzatQtFRUUYMGCA7HpYWlrCzs5O9pfVp4+Vm5uLO3fuoHXr1hBClHmrjupWdtUI1LNcXV2lv7ACT0bOTExM8Ndff6k9fkFBAf7991+8++67MDU1VXuuqmOVdlsT0Yti/1bx/dvAgQORlJSEP//8UyrbsWMH9PT00Lt3bwBPboHPyMjAuHHjZM9Kdu/eHfb29vjll1+e6/zK8vSzflWqVIGTkxOEEPD19ZXKTU1N0aBBA1n/FR4ejoYNG8Le3l722aiecy5rFOvff/8tsf8EAF9fX9kdTM7OzsXapWrv0+1SadSoEaKiohAREYHp06fD0NCw2KzjKtWqVcPDhw/x4MGDUttM9LxycnIAQO3KJB07dkTNmjWlbd26dcVixo4dK3u9f/9+VKlSBRMmTJCVf/rppxBCSLeAt2jRAkZGRoiOjgbwZORadWv0qVOn8ODBAwghcOzYMbRr167M8zAyMsJHH30kvVYoFPjggw9kv3uRkZGoXbu2NFoOPJlwsLzL93344YeIi4tDr169cObMGSxZsgTu7u6oXbu22sd4nlVYWIiDBw/C09MT77zzjlReq1YtDB48GMeOHZM+j59++gnNmjWT7lZ62rN3TmZnZ8PNzQ2XLl3Cb7/9hubNm5fZlrK+IzZq1AguLi7Sa9VIe+fOnVG3bt1i5er6uC+//BIHDx7EN998g1atWiE/Px+PHz8uFsfvjnK8dfw19sEHH8DJyalYebVq1Yr9AB8/fhxz585FXFxcsf+8s7OzoVQqpde2trZqjzdw4EB89tlnuHXrFmrXro3ffvsNGRkZGDhwoBRz5coVZGdny55feZpqAokOHTrAy8sL8+bNw8qVK9GxY0d4enpi8ODB0NPTK/W8r169CiEE5syZgzlz5pR4nNq1a+P69etqb49q0KBBqcdQyczMxLx587B9+/Zik188+wUeKH7trly5AiEE7Ozs1O5fV1dX+veNGzcQEBCAPXv2FHveU92x1BElPBv4dEepUq1aNdlxHj58iKCgIGzZsgW3bt2S7Uvd8VX1b/q6ufR6Yv9W8f1b//79MWXKFOzYsQOfffYZhBAIDw+XnicEID3DrW6f9vb2aif6ehnP9lVKpRL6+vrFbp1UKpWyuTKuXLmCixcvSs++P6s8kxeV1H+W1C4AsudYVeXP9t8AYGJiIj3/3bt3b4SFhaF37944deoUmjVrprYd7FtJ04yNjQFA7R95VCuJpKeny5JYFR0dHdSpU0dWdv36dVhZWUn7VWnYsKFUDzz5I5SLiwtiYmIAPEm027Vrh7Zt26KwsBDx8fGwsLBAZmZmuRLtOnXqFPv9qFatmmxZsuvXr6N+/frF4p5nVYGWLVti165dyM/Px5kzZ7B7926sXLkS/fr1Q3JyMho1alTie//55x88ePBAbd/ZsGFDFBUV4ebNm2jcuDH+/PNPeHl5latNkyZNwqNHj3D69Gk0bty43OcClP87Ymn9GwC1fdzTCf9HH32E999/H8OGDSu27jf7Nzkm2v8Bf/75J7p06QJ7e3usWLEC1tbWUCgU2L9/P1auXFls4pySnu8bOHAgZs2ahfDwcEyaNAk7d+6EUqmUzShYVFQEc3NzhIaGqt2H6kuQlpYWfvzxR8THx2Pv3r04cOAARowYgeXLlyM+Pr7UdaBV7Z06dWqxkSkVTS3PMmDAAMTGxmLatGlo3rw5jIyMUFRUVGxZFpVnr11RUZE0sUeVKlWKxavOs7CwEB9++CEyMzMxY8YM2Nvbw9DQELdu3cKwYcPKnNxI9Zz13bt3i/1HCEDtsQF5pzt+/Hhs2bIFkyZNgouLC5RKJbS0tODt7a32+KqOVt2zQ0SvCvu3F2dlZYV27dph586d+OyzzxAfH48bN25g8eLFGtl/SV+kCgsLS3yPur6qPP1XUVERHBwcsGLFCrWxz35hfFb16tXVfnksqw3qyktL2FX69u2Ljz/+GNu3by+WaN+9exdVq1blbPikcUqlErVq1cL58+eL1an+aKduAlQA0NPTK3Mm8tK0bdtWWkM5JiYGn3/+OUxNTdGkSRPExMTAwsICAMqVaJenT9AkhUKBli1bomXLlnjvvfcwfPhwhIeHY+7cuRVyvNL07t0b27dvx5dffolvv/22XJ/J098R1Xme/g0o+zorFAr06tULX375JR4+fCjry/jdUY6J9n/A3r17kZeXhz179sj+alXeCSFUbG1t8cEHH2DHjh3w9/fHrl274OnpKRuhqV+/Pg4dOoQ2bdqU60tCq1at0KpVKyxatAhhYWHw8fHB9u3bMXLkyBK/pKluwdHV1S11lm3gyWQWV65cKVaekpJSZtvu3r2Lw4cPY968eQgICJDK1e2vJPXr14cQAra2tnjvvfdKjDt37hwuX76Mbdu2ySbRiYqKKtdx7O3tAQCpqalwcHAod/ue9uOPP2Lo0KGy2XwfPXqkdtZI1bFq1KhR4ggS0avA/u3F+jeVgQMHYty4cUhJScGOHTtQtWpV9OzZU3YM1T6fXW4qJSWlxMnDgP+7RfDZPqQi1uetX78+zpw5gy5durzQSIm9vT1++uknjberJHl5eSgqKlJ7t1Bqaqo0Ikikad27d8fXX3+NEydOqJ2M73nY2Njg0KFDuHfvnmxU+9KlS1K9Srt27ZCfn48ffvgBt27dkhLq9u3bS4n2e++9JyXcL8vGxgZ//PEHhBCyPuHq1asvtV/VnVb/+9//pDJ1fU7NmjVRtWpVtf3xpUuXoK2tLf0BsH79+mr/+KGOp6cn3NzcMGzYMBgbG6tdYeNZT39HfFUePnwIIQTu3bsn+/9S1Qb2cU/wGe3/ANVfpJ69FXjLli3Pva+BAwciPj4emzdvxp07d2S3VQJPRoALCwuxYMGCYu9VzboKPElin/2LmOq2E9VSLKpZrp/9kmZubo6OHTtiw4YNso5O5Z9//pH+3a1bN8THx+PEiROy+pJGpJ6m7roBwKpVq8p8r0rfvn1RpUoVzJs3r9h+hBDSrY/qjiWEQHBwcLmO4+joCIVCgcTExHK37VlVqlQp1sY1a9aUOPqUlJQke6aHqDKwf3ux/k3Fy8sLVapUwQ8//IDw8HD06NFDtjauk5MTzM3NERISIlsm69dff8XFixfRvXv3EvdtYmKCGjVqSM9lqnz11Vflbl95DRgwALdu3cKmTZuK1T18+BC5ubmlvt/FxQV3795V++zhy8jKylL7nP7XX38NAGofjzh16hRat26t0XYQqUyfPh1Vq1bFiBEjkJ6eXqz+eUaFu3XrhsLCQqxdu1ZWvnLlSmhpaUkzgANPRsx1dXWxePFimJmZSbc9t2vXDvHx8fj999/LNZpdXu7u7rh165bseepHjx6p7SPUOXr0qNprsX//fgDyx2kMDQ2L9eVVqlSBm5sbfv75Z9ldAunp6QgLC0Pbtm2lR3S8vLykW9Ofpa4NQ4YMwerVqxESEoIZM2aUeS61a9eGtbX1S31HLIm6x3KysrLw008/wdrautijVklJSdDS0uL3x/+PI9r/AW5ublAoFOjZsyfGjBmD+/fvY9OmTTA3N1f7Ra40AwYMwNSpUzF16lSYmZkVG3Hp0KEDxowZg6CgICQnJ8PNzQ26urq4cuUKwsPDERwcjH79+mHbtm346quv0KdPH9SvXx/37t3Dpk2bYGJigm7dugF4cotno0aNsGPHDrz33nswMzNDkyZN0KRJE6xbtw5t27aFg4MDRo0ahXfeeQfp6emIi4vD33//jTNnzgB48h/Kd999Bw8PD0ycOFFa/sbGxkb2LI86JiYmaN++PZYsWYKCggLUrl0bBw8efK6/CNavXx8LFy7ErFmzcO3aNXh6esLY2BipqanYvXs3Ro8ejalTp8Le3h7169fH1KlTcevWLZiYmOCnn34q9VbGp+nr68PNzQ2HDh3C/Pnzy92+p/Xo0QPfffcdlEolGjVqhLi4OBw6dEi65ehpGRkZOHv2bLGlk4heNfZvL9a/qZibm6NTp05YsWIF7t27V+yPC6ovxsOHD0eHDh0waNAgaXmvevXqYfLkyaXuf+TIkfjyyy8xcuRIODk5ITo6GpcvXy7vR1JuH3/8MXbu3IlPPvkER48eRZs2bVBYWIhLly5h586d0trpJenevTt0dHRw6NAhaYkeTfjtt98wYcIE9OvXD3Z2dsjPz0dMTAx27doFJyenYs/CJiUlITMzU5qMjkjT7OzsEBYWhkGDBqFBgwbw8fFBs2bNIIRAamoqwsLCoK2trfYxtGf17NkTnTp1wueff45r166hWbNmOHjwIH7++WdMmjRJNglr1apV4ejoiPj4eGkNbeDJiHZubi5yc3M1mmiPGTMGa9euxaBBgzBx4kTUqlULoaGh0qSOZd35Mn78eDx48AB9+vSBvb098vPzERsbix07dqBevXoYPny4FOvo6IhDhw5hxYoVsLKygq2tLZydnbFw4UJERUWhbdu2GDduHHR0dLBhwwbk5eXJ1v2eNm0afvzxR/Tv3x8jRoyAo6MjMjMzsWfPHoSEhBR7vAQA/P39kZOTg88//xxKpVLtUlpP6927N3bv3l1shP9lde3aFXXq1IGzszPMzc1x48YNbNmyBbdv35ZNMKoSFRWFNm3aqP1u+VaqwBnN6QWplhspaemYDh06FFv+Zs+ePaJp06ZCX19f1KtXTyxevFhs3ry52BIwNjY2onv37qUev02bNgKAGDlyZIkxGzduFI6OjsLAwEAYGxsLBwcHMX36dHH79m0hhBCnTp0SgwYNEnXr1hV6enrC3Nxc9OjRQyQmJsr2ExsbKxwdHYVCoSi2HMyff/4phgwZIiwtLYWurq6oXbu26NGjh/jxxx9l+zh79qzo0KGD0NfXF7Vr1xYLFiwQ33zzTbmWv/n7779Fnz59hKmpqVAqlaJ///7S0gTqlqb5559/1O7np59+Em3bthWGhobC0NBQ2NvbCz8/P5GSkiLF/PHHH8LV1VUYGRmJGjVqiFGjRknLcKlbHudZu3btElpaWuLGjRuycgDCz8+vWLyNjY1sOYq7d++K4cOHixo1aggjIyPh7u4uLl26VCxOCCHWr18vqlatKnJycspsF9HzYP/2xKvo31Q2bdokAAhjY2O1y7EIIcSOHTtEixYthJ6enjAzMxM+Pj7i77//lsU8u7yXEE+WXvP19RVKpVIYGxuLAQMGiIyMjHL3oUOHDhWGhobF2qPu5yA/P18sXrxYNG7cWOjp6Ylq1aoJR0dHMW/ePJGdnV3mdejVq5fo0qWLrKykn8fytvfq1atiyJAh4p133hEGBgZCX19fNG7cWMydO1fcv3+/WBtmzJgh6tatK1tKjKgiXL16VYwdO1a8++67Ql9fXxgYGAh7e3vxySefiOTkZFlsSb+HQghx7949MXnyZGFlZSV0dXWFnZ2dWLp0qdqf4WnTpqld2undd98VAGTLYAlR8vJez/7uq9poY2MjK/vrr79E9+7dhYGBgahZs6b49NNPxU8//SQAiPj4+NIuj/j111/FiBEjhL29vTAyMhIKhUK8++67Yvz48cWWNb106ZJo3769MDAwEABk35lOnTol3N3dhZGRkahataro1KmTbCktlX///Vf4+/uL2rVrC4VCIerUqSOGDh0q7ty5I7sW4eHhsvdNnz5dABBr164t9XxOnTolABRbArGk/xfVfXdULdm4dOlSqWzt2rWibdu2okaNGkJHR0fUrFlT9OzZU7YErkpWVpZQKBTi66+/LrWtbxMtISpoZgEi0qjCwkI0atQIAwYMUHtrqya1aNECHTt2xMqVKyv0OEREr0pMTAw6duyIS5culbhSREXKy8tDvXr1MHPmTEycOPGVH5/obbBq1SpMnjwZf//9N2rXrl3ZzXmlunTpAisrK3z33XeVcvxVq1ZhyZIl+PPPPznZ4//HRJvoDbJjxw6MHTsWN27cKHVm45cRGRmJfv364a+//ipxmSMiojeR6jbI8j7HqUkhISH44osvcOXKlTKXgSOisj074/WjR4/QokULFBYWVsgjLK+7hIQEtGvXDleuXCl1IsuKUFBQgPr162PmzJkYN27cKz3264yJNhERERERvVG6du2KunXronnz5sjOzsb333+PCxcuIDQ0FIMHD67s5hFxMjQiIiIiInqzuLu74+uvv0ZoaKj0eN327duLTfpIVFk4ok1ERERERESkQVxHm4iIiIiIiEiD3shbx4uKinD79m0YGxtrdK04Inp7CCFw7949WFlZQVv7v/M3R/aPRPQ8jh8/jtWrVyM5ORlpaWkIDQ1F9+7d1faPn3zyCTZs2ICVK1di0qRJUnlmZibGjx+PvXv3QltbG15eXggODpZN2nn27Fn4+fnh5MmTqFmzJsaPH4/p06fL2hIeHo45c+bg2rVrsLOzw+LFi6W16YEn/fbcuXOxadMmZGVloU2bNli/fv1zzSLPPpKIXsZzfX+spGXFXsrNmzcFAG7cuHF76e3mzZuV3aVpFPtHbty4aWp7un/ctWuXaNasmbCyshIrV66U9TseHh6iWbNmIj4+XsTExIh3331XDBo0SKrPzs4WFhYWwsfHR5w/f1788MMPwsDAQGzYsEGKOX78uKhSpYpYsmSJ+OOPP8Ts2bOFrq6uOHfunBTz5ZdfCqVSKSIiIsSZM2dEr169hK2tbYnrw7OP5MaNW0Vt5fn++EY+o52dnQ1TU1PcvHkTJiYmld0cInoD5eTkwNraGllZWVAqlZXdHI1h/0hEL0qpVCI0NBTt27eX9Y+3bt2Cs7MzDhw4gO7du2PSpEnSiPbFixfRqFEjnDx5Ek5OTgCeLBPZrVs3/P3337CyssL69evx+eefIy0tDQqFAgAwc+ZMRERE4NKlSwCAgQMHIjc3F/v27ZPa06pVKzRv3hwhISEQQsDKygqffvoppk6dCuBJf2dhYYGtW7fC29u7XOfIPpKIXsbzfH98I28dV93qY2Jiwk6SiF7Kf+3WQfaPRPQyqlatKvUdWlpaKCoqwscff4xp06ahcePGxeLj4uJgamoqJdkA4OrqCm1tbSQkJKBPnz6Ii4tD+/btpSQbeDJj9OLFi3H37l1Uq1YNcXFxmDJlimzf7u7uiIiIAACkpqYiLS0Nrq6uUr1SqYSzszPi4uJKTLTz8vKQl5cnvb537x4A9pFE9HLK8/3xv/NgIhERERFp1OLFi6Gjo4MJEyaorU9LS4O5ubmsTEdHB2ZmZkhLS5NiLCwsZDGq12XFPF3/9PvUxagTFBQEpVIpbdbW1qWeLxGRpjDRJiIiIqJiTp8+jeDgYGzduvWNvftn1qxZyM7OlrabN29WdpOI6C3BRJuIiIiIiomLi0NGRgbq1q0LHR0d6Ojo4Pr16/j0009Rr149AIClpSUyMjJk73v8+DEyMzNhaWkpxaSnp8tiVK/Linm6/un3qYtRR09PT7pNnLeLE9GrxESbiIiIiIrx9vbG2bNnkZycLG1WVlaYNm0aDhw4AABwcXFBVlYWkpKSpPcdOXIERUVFcHZ2lmKio6NRUFAgxURFRaFBgwaoVq2aFHP48GHZ8aOiouDi4gIAsLW1haWlpSwmJycHCQkJUgwR0evkjZwMjYiIiIhe3v3793H16lXpdWpqKs6ePQsAMDMzk0auVXR1dWFpaYkGDRoAABo2bAgPDw+MGjUKISEhKCgogL+/P7y9vWFlZQUAGDx4MObNmwdfX1/MmDED58+fR3BwMFauXCntd+LEiejQoQOWL1+O7t27Y/v27UhMTMTGjRsBPJl4aNKkSVi4cCHs7Oxga2uLOXPmwMrKCp6enhV4hYiIXgwTbSIiIqK3VGJiIjp16iS9fnbm7/IIDQ2Fv78/unTpAm1tbXh5eWH16tVSvVKpxMGDB+Hn5wdHR0fUqFEDAQEBGD16tBTTunVrhIWFYfbs2fjss89gZ2eHiIgINGnSRIqZPn06cnNzMXr0aGRlZaFt27aIjIyEvr7+C549EVHFeSPX0c7JyYFSqUR2djaftSGiF/Jf7Uf+q+dFRK/Of7kf+S+fGxFVvOfpQ/iMNhEREREREZEGMdEmIiIiIiIi0iAm2kREREREREQaxESbiIiIiIiISIOYaBMRERERERFpEJf3ekPVm/lLZTehXK592b2ym0BEb6E3oY9k/0j0ZmB/QkQvgiPaRERERERERBrERJuIiIiIiIhIg5hoExEREREREWkQE20iIiIiIiIiDWKiTURERERERKRBb9Ws45w18vX1Jnw2QPk/nzfhfJ7nZ+2/dj5EVHnehP4EYJ9CREQv561KtImIiN40TEyJiIjePEy0iYiIiIjeEm/CH+/4hzv6L+Az2kREREREREQaxBFtInrrRUdHY+nSpUhKSsL//vc/7N69G56enmpjP/nkE2zYsAErV67EpEmTpPLMzEyMHz8ee/fuhba2Nry8vBAcHAwjIyMp5uzZs/Dz88PJkydRs2ZNjB8/HtOnT5ftPzw8HHPmzMG1a9dgZ2eHxYsXo1u3bhVx2kSkARwdJCIidTiiTURvvdzcXDRr1gzr1q0rNW737t2Ij4+HlZVVsTofHx9cuHABUVFR2LdvH6KjozF69GipPicnB25ubrCxsUFSUhKWLl2KwMBAbNy4UYqJjY3FoEGD4Ovri9OnT8PT0xOenp44f/685k6WiIiIiCrcc49olzXyI4TA3LlzsWnTJmRlZaFNmzZYv3497OzspJjyjPwQEb0qXbt2RdeuXUuNuXXrFsaPH48DBw6ge3f56NDFixcRGRmJkydPwsnJCQCwZs0adOvWDcuWLYOVlRVCQ0ORn5+PzZs3Q6FQoHHjxkhOTsaKFSukhDw4OBgeHh6YNm0aAGDBggWIiorC2rVrERISUgFn/n84KkdERESkOc89ol3WyM+SJUuwevVqhISEICEhAYaGhnB3d8ejR4+kmLJGfoiIXidFRUX4+OOPMW3aNDRu3LhYfVxcHExNTaUkGwBcXV2hra2NhIQEKaZ9+/ZQKBRSjLu7O1JSUnD37l0pxtXVVbZvd3d3xMXFldi2vLw85OTkyDYiIiIiqlzPPaJd2siPEAKrVq3C7Nmz0bt3bwDAt99+CwsLC0RERMDb27tcIz9ERK+TxYsXQ0dHBxMmTFBbn5aWBnNzc1mZjo4OzMzMkJaWJsXY2trKYiwsLKS6atWqIS0tTSp7Oka1D3WCgoIwb9685z4nIiIiIqo4Gn1GOzU1FWlpabIRGaVSCWdnZ2lEpjwjP8/iiA0RVZakpCQEBwdj69at0NLSquzmFDNr1ixkZ2dL282bNyu7SURERERvPY0m2qpRl9JGZMoz8vOsoKAgKJVKabO2ttZks4mIShQTE4OMjAzUrVsXOjo60NHRwfXr1/Hpp5+iXr16AABLS0tkZGTI3vf48WNkZmbC0tJSiklPT5fFqF6XFaOqV0dPTw8mJiayjYiIiIgq1xsx6zhHbIiosnz88cc4e/YskpOTpc3KygrTpk3DgQMHAAAuLi7IyspCUlKS9L4jR46gqKgIzs7OUkx0dDQKCgqkmKioKDRo0ADVqlWTYg4fPiw7flRUFFxcXCr6NImIiIhIgzS6jrZq1CU9PR21atWSytPT09G8eXMppqyRn2fp6elBT09Pk00lIpLcv38fV69elV6npqYiOTkZZmZmqFu3LqpXry6L19XVhaWlJRo0aAAAaNiwITw8PDBq1CiEhISgoKAA/v7+8Pb2luadGDx4MObNmwdfX1/MmDED58+fR3BwMFauXCntd+LEiejQoQOWL1+O7t27Y/v27UhMTJQtAUZERET/h6tm0OtKoyPatra2sLS0lI3I5OTkICEhQRqRKc/IDxHRq5SYmIgWLVqgRYsWAIApU6agRYsWCAgIKPc+QkNDYW9vjy5duqBbt25o27atLEFWKpU4ePAgUlNT4ejoiE8//RQBAQGyFRdat26NsLAwbNy4Ec2aNcOPP/6IiIgINGnSRHMnS0REREQV7rlHtMsa+Zk0aRIWLlwIOzs72NraYs6cObCyspLW2i7PyA8R0avUsWNHCCHKHX/t2rViZWZmZggLCyv1fU2bNkVMTEypMf3790f//v3L3RYiIiIiev08d6KdmJiITp06Sa+nTJkCABg6dCi2bt2K6dOnIzc3F6NHj0ZWVhbatm2LyMhI6OvrS+8JDQ2Fv78/unTpAm1tbXh5eWH16tUaOB0iIiIiIiKiyvXciXZZIz9aWlqYP38+5s+fX2JMeUZ+iIiIiIiIiN5Eb8Ss40RERERERERvCibaRERERERERBrERJuIiIiIiIhIg5hoExEREREREWkQE20iIiIiIiIiDXruWceJiIiI6L8hOjoaS5cuRVJSEv73v/9h9+7d6Ny5MwCgoKAAM2bMwP79+/HXX39BqVTC1dUVX375JaysrKR9ZGZmYvz48di7d6+0bGtwcDCMjIykmLNnz8LPzw8nT55EzZo1MX78eEyfPl3WlvDwcMyZMwfXrl2DnZ0dFi9ejG7dukn1QgjMnTsXmzZtQlZWFtq0aYP169fDzs6ugq8S0atRb+Yvld2Ecrn2ZffKbsIbgSPaRERERG+p3NxcNGvWDOvWrStW9+DBA5w6dQpz5szBqVOnsGvXLqSkpKBXr16yOB8fH1y4cAFRUVHYt28foqOjMXr0aKk+JycHbm5usLGxQVJSEpYuXYrAwEBs3LhRiomNjcWgQYPg6+uL06dPw9PTE56enjh//rwUs2TJEqxevRohISFISEiAoaEh3N3d8ejRowq4MkREL4cj2kRERERvqa5du6Jr165q65RKJaKiomRla9euxQcffIAbN26gbt26uHjxIiIjI3Hy5Ek4OTkBANasWYNu3bph2bJlsLKyQmhoKPLz87F582YoFAo0btwYycnJWLFihZSQBwcHw8PDA9OmTQMALFiwAFFRUVi7di1CQkIghMCqVaswe/Zs9O7dGwDw7bffwsLCAhEREfD29q6oS0RE9EI4ok1ERERE5ZKdnQ0tLS2YmpoCAOLi4mBqaiol2QDg6uoKbW1tJCQkSDHt27eHQqGQYtzd3ZGSkoK7d+9KMa6urrJjubu7Iy4uDgCQmpqKtLQ0WYxSqYSzs7MUo05eXh5ycnJkGxHRq8BEm4iIiIjK9OjRI8yYMQODBg2CiYkJACAtLQ3m5uayOB0dHZiZmSEtLU2KsbCwkMWoXpcV83T90+9TF6NOUFAQlEqltFlbWz/XORMRvSgm2kRERERUqoKCAgwYMABCCKxfv76ym1Nus2bNQnZ2trTdvHmzsptERG8JPqNNRERERCVSJdnXr1/HkSNHpNFsALC0tERGRoYs/vHjx8jMzISlpaUUk56eLotRvS4r5ul6VVmtWrVkMc2bNy+x7Xp6etDT03ue0yUiDXnbZ1HniDYRERERqaVKsq9cuYJDhw6hevXqsnoXFxdkZWUhKSlJKjty5AiKiorg7OwsxURHR6OgoECKiYqKQoMGDVCtWjUp5vDhw7J9R0VFwcXFBQBga2sLS0tLWUxOTg4SEhKkGCKi1wkTbSIiIqK31P3795GcnIzk5GQATyYdO3v2LIAnSXa/fv2QmJiI0NBQFBYWIi0tDWlpacjPzwcANGzYEB4eHhg1ahROnDiB48ePw9/fH97e3tJa24MHD4ZCoYCvry8uXLiAHTt2IDg4GFOmTJHaMXHiRERGRmL58uW4dOkSAgMDkZiYCH9/fwCAlpYWJk2ahIULF2LPnj04d+4chgwZAisrK3h6er66C0ZEVE68dZyIiIjoLZWYmIhOnTpJr59Ofm/fvo09e/YAQLHbs48ePYqOHTsCAEJDQ+Hv748uXbpAW1sbXl5eWL16tRSrVCpx8OBB+Pn5wdHRETVq1EBAQIBsre3WrVsjLCwMs2fPxmeffQY7OztERESgSZMmUsz06dORm5uL0aNHIysrC23btkVkZCT09fU1eUmIiDSCiTYRERHRW6pjx44QQsjKcnJyoFQqYWNjU6xOHTMzM4SFhZUa07RpU8TExJQa079/f/Tv37/Eei0tLcyfPx/z588vs01ERJWNt44TERERERERaRATbSIiIiIiIiINYqJNREREREREpEFMtImIiIiIiIg0iIk2ERERERERkQYx0SYiIiIiIiLSICbaRERERERERBrERJuIiIiIiIhIg5hoExEREREREWkQE20iIiIiIiIiDWKiTURERERERKRBTLSJ6K0XHR2Nnj17wsrKClpaWoiIiJDqCgoKMGPGDDg4OMDQ0BBWVlYYMmQIbt++LdtHZmYmfHx8YGJiAlNTU/j6+uL+/fuymLNnz6Jdu3bQ19eHtbU1lixZUqwt4eHhsLe3h76+PhwcHLB///4KOWciIiIiqjhMtInorZebm4tmzZph3bp1xeoePHiAU6dOYc6cOTh16hR27dqFlJQU9OrVSxbn4+ODCxcuICoqCvv27UN0dDRGjx4t1efk5MDNzQ02NjZISkrC0qVLERgYiI0bN0oxsbGxGDRoEHx9fXH69Gl4enrC09MT58+fr7iTJyIiIiKN06nsBhARVbauXbuia9euauuUSiWioqJkZWvXrsUHH3yAGzduoG7durh48SIiIyNx8uRJODk5AQDWrFmDbt26YdmyZbCyskJoaCjy8/OxefNmKBQKNG7cGMnJyVixYoWUkAcHB8PDwwPTpk0DACxYsABRUVFYu3YtQkJCKvAKEBEREZEmcUSbiOg5ZWdnQ0tLC6ampgCAuLg4mJqaSkk2ALi6ukJbWxsJCQlSTPv27aFQKKQYd3d3pKSk4O7du1KMq6ur7Fju7u6Ii4srsS15eXnIycmRbURERERUuZhoExE9h0ePHmHGjBkYNGgQTExMAABpaWkwNzeXxeno6MDMzAxpaWlSjIWFhSxG9bqsGFW9OkFBQVAqldJmbW39cidIRERERC+NiTYRUTkVFBRgwIABEEJg/fr1ld0cAMCsWbOQnZ0tbTdv3qzsJhERERG99fiMNhFROaiS7OvXr+PIkSPSaDYAWFpaIiMjQxb/+PFjZGZmwtLSUopJT0+XxahelxWjqldHT08Penp6L35iRERERKRxHNEmIiqDKsm+cuUKDh06hOrVq8vqXVxckJWVhaSkJKnsyJEjKCoqgrOzsxQTHR2NgoICKSYqKgoNGjRAtWrVpJjDhw/L9h0VFQUXF5eKOjUiIiIiqgBMtInorXf//n0kJycjOTkZAJCamork5GTcuHEDBQUF6NevHxITExEaGorCwkKkpaUhLS0N+fn5AICGDRvCw8MDo0aNwokTJ3D8+HH4+/vD29sbVlZWAIDBgwdDoVDA19cXFy5cwI4dOxAcHIwpU6ZI7Zg4cSIiIyOxfPlyXLp0CYGBgUhMTIS/v/8rvyZERERE9OKYaBPRWy8xMREtWrRAixYtAABTpkxBixYtEBAQgFu3bmHPnj34+++/0bx5c9SqVUvaYmNjpX2EhobC3t4eXbp0Qbdu3dC2bVvZGtlKpRIHDx5EamoqHB0d8emnnyIgIEC21nbr1q0RFhaGjRs3olmzZvjxxx8RERGBJk2avLqLQUREREQvjc9oE9Fbr2PHjhBClFhfWp2KmZkZwsLCSo1p2rQpYmJiSo3p378/+vfvX+bxiIiIiOj1xRFtIiIiIiIiIg1iok1ERERERESkQUy0iYiIiIiIiDSIiTYRERERERGRBjHRJiIiIiIiItIgJtpEREREREREGqTxRLuwsBBz5syBra0tDAwMUL9+fSxYsEC2PI4QAgEBAahVqxYMDAzg6uqKK1euaLopRERERERERK+cxhPtxYsXY/369Vi7di0uXryIxYsXY8mSJVizZo0Us2TJEqxevRohISFISEiAoaEh3N3d8ejRI003h4iIiIiIiOiV0tH0DmNjY9G7d290794dAFCvXj388MMPOHHiBIAno9mrVq3C7Nmz0bt3bwDAt99+CwsLC0RERMDb21vTTSIiIiIiIiJ6ZTQ+ot26dWscPnwYly9fBgCcOXMGx44dQ9euXQEAqampSEtLg6urq/QepVIJZ2dnxMXFqd1nXl4ecnJyZBsRERERERHR60jjifbMmTPh7e0Ne3t76OrqokWLFpg0aRJ8fHwAAGlpaQAACwsL2fssLCykumcFBQVBqVRKm7W1taabTURERPTWiY6ORs+ePWFlZQUtLS1ERETI6sszr05mZiZ8fHxgYmICU1NT+Pr64v79+7KYs2fPol27dtDX14e1tTWWLFlSrC3h4eGwt7eHvr4+HBwcsH///uduCxHR60LjifbOnTsRGhqKsLAwnDp1Ctu2bcOyZcuwbdu2F97nrFmzkJ2dLW03b97UYIuJiIiI3k65ublo1qwZ1q1bp7a+PPPq+Pj44MKFC4iKisK+ffsQHR2N0aNHS/U5OTlwc3ODjY0NkpKSsHTpUgQGBmLjxo1STGxsLAYNGgRfX1+cPn0anp6e8PT0xPnz55+rLURErwuNP6M9bdo0aVQbABwcHHD9+nUEBQVh6NChsLS0BACkp6ejVq1a0vvS09PRvHlztfvU09ODnp6epptKRERE9Fbr2rWr9Hjfs8ozr87FixcRGRmJkydPwsnJCQCwZs0adOvWDcuWLYOVlRVCQ0ORn5+PzZs3Q6FQoHHjxkhOTsaKFSukhDw4OBgeHh6YNm0aAGDBggWIiorC2rVrERISwjl+iOiNo/ER7QcPHkBbW77bKlWqoKioCABga2sLS0tLHD58WKrPyclBQkICXFxcNN0cIiIiInoB165dK3Nenbi4OJiamkpJNgC4urpCW1sbCQkJUkz79u2hUCikGHd3d6SkpODu3btSzNPHUcWojvMic/wAnOeHiCqPxhPtnj17YtGiRfjll19w7do17N69GytWrECfPn0AAFpaWpg0aRIWLlyIPXv24Ny5cxgyZAisrKzg6emp6eYQERER0QvIyMgAUPq8OmlpaTA3N5fV6+jowMzMTBajbh+qutJinq4vqy3qcJ4fIqosGr91fM2aNZgzZw7GjRuHjIwMWFlZYcyYMQgICJBipk+fjtzcXIwePRpZWVlo27YtIiMjoa+vr+nmEBEREdFbatasWZgyZYr0Oicnh8k2Eb0SGk+0jY2NsWrVKqxatarEGC0tLcyfPx/z58/X9OGJiIiISANUI9WlzatjaWkpjXyrPH78GJmZmdK8PJaWlkhPT5fFqF6XFfN0fVltUYfz/BBRZdH4reNERERE9OarV69emfPquLi4ICsrC0lJSVLMkSNHUFRUBGdnZykmOjoaBQUFUkxUVBQaNGiAatWqSTFPH0cVozoO5/ghojcNE20iIiKit9T9+/eRnJyM5ORkAE8mHTt79iyA8s2r07BhQ3h4eGDUqFE4ceIEjh8/Dn9/f3h7e8PKygoAMHjwYCgUCvj6+uLChQvYsWMHgoODZbd0T5w4EZGRkVi+fDkuXbqEwMBAJCYmwt/fv9xtISJ6nWj81nEiIiIiejMkJiaiU6dO0uunk1+gfPPqhIaGwt/fH126dIG2tja8vLywevVqqV6pVOLgwYPw8/ODo6MjatSogYCAANla261bt0ZYWBhmz56Nzz77DHZ2doiIiECTJk2eqy1ERK8LJtpEREREb6mOHTtCCCEry8nJgVKpBFC+eXXMzMwQFhZW6nGaNm2KmJiYUmP69++P/v37l1jPOX6I6E3CW8eJiIiIiIiINIiJNhEREREREZEGMdEmIiIiIiIi0iAm2kREREREREQaxESbiIiIiIiISIOYaBMRERERERFpEBNtIiIiIiIiIg1iok1ERERERESkQUy0ieitFx0djZ49e8LKygpaWlqIiIiQ1QshEBAQgFq1asHAwACurq64cuWKLCYzMxM+Pj4wMTGBqakpfH19cf/+fVnM2bNn0a5dO+jr68Pa2hpLliwp1pbw8HDY29tDX18fDg4O2L9/v8bPl4iIiIgqFhNtInrr5ebmolmzZli3bp3a+iVLlmD16tUICQlBQkICDA0N4e7ujkePHkkxPj4+uHDhAqKiorBv3z5ER0dj9OjRUn1OTg7c3NxgY2ODpKQkLF26FIGBgdi4caMUExsbi0GDBsHX1xenT5+Gp6cnPD09cf78+Yo7eSIiIiLSOJ3KbgARUWXr2rUrunbtqrZOCIFVq1Zh9uzZ6N27NwDg22+/hYWFBSIiIuDt7Y2LFy8iMjISJ0+ehJOTEwBgzZo16NatG5YtWwYrKyuEhoYiPz8fmzdvhkKhQOPGjZGcnIwVK1ZICXlwcDA8PDwwbdo0AMCCBQsQFRWFtWvXIiQk5BVcCSIiIiLSBI5oExGVIjU1FWlpaXB1dZXKlEolnJ2dERcXBwCIi4uDqamplGQDgKurK7S1tZGQkCDFtG/fHgqFQopxd3dHSkoK7t69K8U8fRxVjOo46uTl5SEnJ0e2EREREVHlYqJNRFSKtLQ0AICFhYWs3MLCQqpLS0uDubm5rF5HRwdmZmayGHX7ePoYJcWo6tUJCgqCUqmUNmtr6+c9RSIiIiLSMCbaRERvsFmzZiE7O1vabt68WdlNIiIiInrrMdEmIiqFpaUlACA9PV1Wnp6eLtVZWloiIyNDVv/48WNkZmbKYtTt4+ljlBSjqldHT08PJiYmso2IiIiIKhcTbSKiUtja2sLS0hKHDx+WynJycpCQkAAXFxcAgIuLC7KyspCUlCTFHDlyBEVFRXB2dpZioqOjUVBQIMVERUWhQYMGqFatmhTz9HFUMarjEBEREdGbgYk2Eb317t+/j+TkZCQnJwN4MgFacnIybty4AS0tLUyaNAkLFy7Enj17cO7cOQwZMgRWVlbw9PQEADRs2BAeHh4YNWoUTpw4gePHj8Pf3x/e3t6wsrICAAwePBgKhQK+vr64cOECduzYgeDgYEyZMkVqx8SJExEZGYnly5fj0qVLCAwMRGJiIvz9/V/1JSEiIiKil8DlvYjorZeYmIhOnTpJr1XJ79ChQ7F161ZMnz4dubm5GD16NLKystC2bVtERkZCX19fek9oaCj8/f3RpUsXaGtrw8vLC6tXr5bqlUolDh48CD8/Pzg6OqJGjRoICAiQrbXdunVrhIWFYfbs2fjss89gZ2eHiIgINGnS5BVcBSIiIiLSFCbaRPTW69ixI4QQJdZraWlh/vz5mD9/fokxZmZmCAsLK/U4TZs2RUxMTKkx/fv3R//+/UtvMBERERG91njrOBEREREREZEGMdEmIiIiIiIi0iAm2kREREREREQaxESbiIiIiIiISIOYaBMRERERERFpEBNtIiIiIiIiIg1iok1ERERERESkQUy0iYiIiIiIiDSIiTYRERERERGRBjHRJiIiIiIiItIgJtpEREREREREGsREm4iIiIiIiEiDmGgTERERERERaRATbSIiIiIiIiINYqJNRERERGoVFhZizpw5sLW1hYGBAerXr48FCxZACCHFCCEQEBCAWrVqwcDAAK6urrhy5YpsP5mZmfDx8YGJiQlMTU3h6+uL+/fvy2LOnj2Ldu3aQV9fH9bW1liyZEmx9oSHh8Pe3h76+vpwcHDA/v37K+bEiYheEhNtIiIiIlJr8eLFWL9+PdauXYuLFy9i8eLFWLJkCdasWSPFLFmyBKtXr0ZISAgSEhJgaGgId3d3PHr0SIrx8fHBhQsXEBUVhX379iE6OhqjR4+W6nNycuDm5gYbGxskJSVh6dKlCAwMxMaNG6WY2NhYDBo0CL6+vjh9+jQ8PT3h6emJ8+fPv5qLQUT0HJhoExEREZFasbGx6N27N7p374569eqhX79+cHNzw4kTJwA8Gc1etWoVZs+ejd69e6Np06b49ttvcfv2bURERAAALl68iMjISHz99ddwdnZG27ZtsWbNGmzfvh23b98GAISGhiI/Px+bN29G48aN4e3tjQkTJmDFihVSW4KDg+Hh4YFp06ahYcOGWLBgAd5//32sXbv2lV8XIqKyMNEmIiIiIrVat26Nw4cP4/LlywCAM2fO4NixY+jatSsAIDU1FWlpaXB1dZXeo1Qq4ezsjLi4OABAXFwcTE1N4eTkJMW4urpCW1sbCQkJUkz79u2hUCikGHd3d6SkpODu3btSzNPHUcWojqNOXl4ecnJyZBsR0augU9kNICIiIqLX08yZM5GTkwN7e3tUqVIFhYWFWLRoEXx8fAAAaWlpAAALCwvZ+ywsLKS6tLQ0mJuby+p1dHRgZmYmi7G1tS22D1VdtWrVkJaWVupx1AkKCsK8efOe97SJiF4aR7SJiIiISK2dO3ciNDQUYWFhOHXqFLZt24Zly5Zh27Ztld20cpk1axays7Ol7ebNm5XdJCJ6S3BEm4iIiIjUmjZtGmbOnAlvb28AgIODA65fv46goCAMHToUlpaWAID09HTUqlVLel96ejqaN28OALC0tERGRoZsv48fP0ZmZqb0fktLS6Snp8tiVK/LilHVq6Onpwc9Pb3nPW0iopfGEW0iIiIiUuvBgwfQ1pZ/XaxSpQqKiooAALa2trC0tMThw4el+pycHCQkJMDFxQUA4OLigqysLCQlJUkxR44cQVFREZydnaWY6OhoFBQUSDFRUVFo0KABqlWrJsU8fRxVjOo4RESvkwpJtG/duoWPPvoI1atXh4GBARwcHJCYmCjVl2e9RSIiIiKqXD179sSiRYvwyy+/4Nq1a9i9ezdWrFiBPn36AAC0tLQwadIkLFy4EHv27MG5c+cwZMgQWFlZwdPTEwDQsGFDeHh4YNSoUThx4gSOHz8Of39/eHt7w8rKCgAwePBgKBQK+Pr64sKFC9ixYweCg4MxZcoUqS0TJ05EZGQkli9fjkuXLiEwMBCJiYnw9/d/5deFiKgsGr91/O7du2jTpg06deqEX3/9FTVr1sSVK1ekv0YC/7fe4rZt22Bra4s5c+bA3d0df/zxB/T19TXdJCIiIiJ6AWvWrMGcOXMwbtw4ZGRkwMrKCmPGjEFAQIAUM336dOTm5mL06NHIyspC27ZtERkZKftOFxoaCn9/f3Tp0gXa2trw8vLC6tWrpXqlUomDBw/Cz88Pjo6OqFGjBgICAmRrbbdu3RphYWGYPXs2PvvsM9jZ2SEiIgJNmjR5NReDiOg5aDzRXrx4MaytrbFlyxap7OlZJJ9dbxEAvv32W1hYWCAiIkJ6BoiIiIiIKpexsTFWrVqFVatWlRijpaWF+fPnY/78+SXGmJmZISwsrNRjNW3aFDExMaXG9O/fH/379y81hojodaDxW8f37NkDJycn9O/fH+bm5mjRogU2bdok1ZdnvcVncQ1EIiIiIiIielNoPNH+66+/sH79etjZ2eHAgQMYO3YsJkyYIC0DUZ71Fp8VFBQEpVIpbdbW1ppuNhEREREREZFGaDzRLioqwvvvv48vvvgCLVq0wOjRozFq1CiEhIS88D65BiIRERERERG9KTSeaNeqVQuNGjWSlTVs2BA3btwAANl6i08rbR1EPT09mJiYyDYiIiIiIiKi15HGE+02bdogJSVFVnb58mXY2NgAKN96i0REr5PCwkLMmTMHtra2MDAwQP369bFgwQIIIaSY8ixbmJmZCR8fH5iYmMDU1BS+vr64f/++LObs2bNo164d9PX1YW1tjSVLlryScyQiIiIizdF4oj158mTEx8fjiy++wNWrVxEWFoaNGzfCz88PQPnWWyQiep0sXrwY69evx9q1a3Hx4kUsXrwYS5YswZo1a6QY1bKFISEhSEhIgKGhIdzd3fHo0SMpxsfHBxcuXEBUVBT27duH6Oho2dI1OTk5cHNzg42NDZKSkrB06VIEBgZi48aNr/R8iYiIiOjlaHx5r5YtW2L37t2YNWsW5s+fD1tbW6xatQo+Pj5STHnWWyQiel3Exsaid+/e6N69OwCgXr16+OGHH3DixAkA5Vu28OLFi4iMjMTJkyfh5OQE4Mn6tN26dcOyZctgZWWF0NBQ5OfnY/PmzVAoFGjcuDGSk5OxYsUKWUJORERERK83jY9oA0CPHj1w7tw5PHr0CBcvXsSoUaNk9ar1FtPS0vDo0SMcOnQI7733XkU0hYjopbVu3RqHDx/G5cuXAQBnzpzBsWPH0LVrVwDlW7YwLi4OpqamUpINAK6urtDW1kZCQoIU0759eygUCinG3d0dKSkpuHv3rtq2cflDIiIiotePxke0iYj+a2bOnImcnBzY29ujSpUqKCwsxKJFi6Q7dcqzbGFaWhrMzc1l9To6OjAzM5PF2NraFtuHqq5atWrF2hYUFIR58+Zp4CyJiIiISFMqZESbiOi/ZOfOnQgNDUVYWBhOnTqFbdu2YdmyZdi2bVtlN43LHxIRERG9hjiiTURUhmnTpmHmzJnw9vYGADg4OOD69esICgrC0KFDZcsW1qpVS3pfeno6mjdvDuDJ0oYZGRmy/T5+/BiZmZnS+y0tLdUufaiqU0dPTw96enovf5JEREREpDEc0SYiKsODBw+grS3vLqtUqYKioiIA5Vu20MXFBVlZWUhKSpJijhw5gqKiIjg7O0sx0dHRKCgokGKioqLQoEEDtbeNExEREdHriYk2EVEZevbsiUWLFuGXX37BtWvXsHv3bqxYsQJ9+vQBUL5lCxs2bAgPDw+MGjUKJ06cwPHjx+Hv7w9vb29YWVkBAAYPHgyFQgFfX19cuHABO3bsQHBwMKZMmVJZp05EREREL4C3jhMRlWHNmjWYM2cOxo0bh4yMDFhZWWHMmDEICAiQYsqzbGFoaCj8/f3RpUsXaGtrw8vLC6tXr5bqlUolDh48CD8/Pzg6OqJGjRoICAjg0l5EREREbxgm2kREZTA2NsaqVauwatWqEmNUyxbOnz+/xBgzMzOEhYWVeqymTZsiJibmRZtKRERERK8B3jpOREREREREpEFMtImIiIiIiIg0iIk2ERERERERkQYx0SYiIiIiIiLSICbaRERERERERBrERJuIiIiIiIhIg5hoExEREREREWkQE20iIiIiIiIiDWKiTURERERERKRBTLSJiIiIiIiINIiJNhEREREREZEGMdEmIiIiIiIi0iAm2kREREREREQaxESbiIiIiIiISIOYaBMRERERERFpEBNtIiIiIiIiIg1iok1ERERERESkQUy0iYiIiIiIiDSIiTYRERERlejWrVv46KOPUL16dRgYGMDBwQGJiYlSvRACAQEBqFWrFgwMDODq6oorV67I9pGZmQkfHx+YmJjA1NQUvr6+uH//vizm7NmzaNeuHfT19WFtbY0lS5YUa0t4eDjs7e2hr68PBwcH7N+/v2JOmojoJTHRJiIiIiK17t69izZt2kBXVxe//vor/vjjDyxfvhzVqlWTYpYsWYLVq1cjJCQECQkJMDQ0hLu7Ox49eiTF+Pj44MKFC4iKisK+ffsQHR2N0aNHS/U5OTlwc3ODjY0NkpKSsHTpUgQGBmLjxo1STGxsLAYNGgRfX1+cPn0anp6e8PT0xPnz51/NxSAieg46ld0AIiIiIno9LV68GNbW1tiyZYtUZmtrK/1bCIFVq1Zh9uzZ6N27NwDg22+/hYWFBSIiIuDt7Y2LFy8iMjISJ0+ehJOTEwBgzZo16NatG5YtWwYrKyuEhoYiPz8fmzdvhkKhQOPGjZGcnIwVK1ZICXlwcDA8PDwwbdo0AMCCBQsQFRWFtWvXIiQkRG378/LykJeXJ73OycnR7AUiIioBR7SJiIiISK09e/bAyckJ/fv3h7m5OVq0aIFNmzZJ9ampqUhLS4Orq6tUplQq4ezsjLi4OABAXFwcTE1NpSQbAFxdXaGtrY2EhAQppn379lAoFFKMu7s7UlJScPfuXSnm6eOoYlTHUScoKAhKpVLarK2tX+JqEBGVHxNtIiIiIlLrr7/+wvr162FnZ4cDBw5g7NixmDBhArZt2wYASEtLAwBYWFjI3mdhYSHVpaWlwdzcXFavo6MDMzMzWYy6fTx9jJJiVPXqzJo1C9nZ2dJ28+bN5zp/IqIXxVvHiYiIiEitoqIiODk54YsvvgAAtGjRAufPn0dISAiGDh1aya0rm56eHvT09Cq7GUT0FuKINhERERGpVatWLTRq1EhW1rBhQ9y4cQMAYGlpCQBIT0+XxaSnp0t1lpaWyMjIkNU/fvwYmZmZshh1+3j6GCXFqOqJiF4nTLSJiIiISK02bdogJSVFVnb58mXY2NgAeDIxmqWlJQ4fPizV5+TkICEhAS4uLgAAFxcXZGVlISkpSYo5cuQIioqK4OzsLMVER0ejoKBAiomKikKDBg2kGc5dXFxkx1HFqI5DRPQ6YaJNRERERGpNnjwZ8fHx+OKLL3D16lWEhYVh48aN8PPzAwBoaWlh0qRJWLhwIfbs2YNz585hyJAhsLKygqenJ4AnI+AeHh4YNWoUTpw4gePHj8Pf3x/e3t6wsrICAAwePBgKhQK+vr64cOECduzYgeDgYEyZMkVqy8SJExEZGYnly5fj0qVLCAwMRGJiIvz9/V/5dSEiKguf0SYiIiIitVq2bIndu3dj1qxZmD9/PmxtbbFq1Sr4+PhIMdOnT0dubi5Gjx6NrKwstG3bFpGRkdDX15diQkND4e/vjy5dukBbWxteXl5YvXq1VK9UKnHw4EH4+fnB0dERNWrUQEBAgGyt7datWyMsLAyzZ8/GZ599Bjs7O0RERKBJkyav5mIQET0HJtpEREREVKIePXqgR48eJdZraWlh/vz5mD9/fokxZmZmCAsLK/U4TZs2RUxMTKkx/fv3R//+/UtvMBHRa4C3jhMRlcOtW7fw0UcfoXr16jAwMICDgwMSExOleiEEAgICUKtWLRgYGMDV1RVXrlyR7SMzMxM+Pj4wMTGBqakpfH19cf/+fVnM2bNn0a5dO+jr68Pa2hpLlix5JedHRERERJrDRJuIqAx3795FmzZtoKuri19//RV//PEHli9fLk3QAwBLlizB6tWrERISgoSEBBgaGsLd3R2PHj2SYnx8fHDhwgVERUVh3759iI6Olt0WmZOTAzc3N9jY2CApKQlLly5FYGAgNm7c+ErPl4iIiIheDm8dJyIqw+LFi2FtbY0tW7ZIZba2ttK/hRBYtWoVZs+ejd69ewMAvv32W1hYWCAiIgLe3t64ePEiIiMjcfLkSTg5OQEA1qxZg27dumHZsmWwsrJCaGgo8vPzsXnzZigUCjRu3BjJyclYsWKFLCF/Wl5eHvLy8qTXOTk5FXEJiIiIiOg5cESbiKgMe/bsgZOTE/r37w9zc3O0aNECmzZtkupTU1ORlpYGV1dXqUypVMLZ2RlxcXEAgLi4OJiamkpJNgC4urpCW1sbCQkJUkz79u2hUCikGHd3d6SkpODu3btq2xYUFASlUilt1tbWGj13IiIiInp+TLSJiMrw119/Yf369bCzs8OBAwcwduxYTJgwAdu2bQMApKWlAQAsLCxk77OwsJDq0tLSYG5uLqvX0dGBmZmZLEbdPp4+xrNmzZqF7Oxsabt58+ZLni0RERERvSzeOk5EVIaioiI4OTnhiy++AAC0aNEC58+fR0hICIYOHVqpbdPT04Oenl6ltoGIiIiI5DiiTURUhlq1aqFRo0aysoYNG+LGjRsAAEtLSwBAenq6LCY9PV2qs7S0REZGhqz+8ePHyMzMlMWo28fTxyAiIiKi1x8TbSKiMrRp0wYpKSmyssuXL8PGxgbAk4nRLC0tcfjwYak+JycHCQkJcHFxAQC4uLggKysLSUlJUsyRI0dQVFQEZ2dnKSY6OhoFBQVSTFRUFBo0aCCb4ZyIiIiIXm8Vnmh/+eWX0NLSwqRJk6SyR48ewc/PD9WrV4eRkRG8vLyKjeIQEb0uJk+ejPj4eHzxxRe4evUqwsLCsHHjRvj5+QGA1MctXLgQe/bswblz5zBkyBBYWVnB09MTwJMRcA8PD4waNQonTpzA8ePH4e/vD29vb1hZWQEABg8eDIVCAV9fX1y4cAE7duxAcHAwpkyZUlmnTkREREQvoEIT7ZMnT2LDhg1o2rSprHzy5MnYu3cvwsPD8fvvv+P27dvo27dvRTaFiOiFtWzZErt378YPP/yAJk2aYMGCBVi1ahV8fHykmOnTp2P8+PEYPXo0WrZsifv37yMyMhL6+vpSTGhoKOzt7dGlSxd069YNbdu2la2RrVQqcfDgQaSmpsLR0RGffvopAgICSlzai4iIiIheTxU2Gdr9+/fh4+ODTZs2YeHChVJ5dnY2vvnmG4SFhaFz584AgC1btqBhw4aIj49Hq1atKqpJREQvrEePHujRo0eJ9VpaWpg/fz7mz59fYoyZmRnCwsJKPU7Tpk0RExPzwu0kIiIiospXYSPafn5+6N69u2xdWQBISkpCQUGBrNze3h5169aV1pt9Vl5eHnJycmQbERERERER0euoQka0t2/fjlOnTuHkyZPF6tLS0qBQKGBqaiorf3q92WcFBQVh3rx5FdFUIiIiIiIiIo3S+Ij2zZs3MXHiRISGhsqeTXwZs2bNQnZ2trTdvHlTI/slIiIiIiIi0jSNJ9pJSUnIyMjA+++/Dx0dHejo6OD333/H6tWroaOjAwsLC+Tn5yMrK0v2vqfXm32Wnp4eTExMZBsRERERERHR60jjt4536dIF586dk5UNHz4c9vb2mDFjBqytraGrq4vDhw/Dy8sLAJCSkoIbN25I680SERERERERvak0nmgbGxujSZMmsjJDQ0NUr15dKvf19cWUKVNgZmYGExMTjB8/Hi4uLpxxnIiIiIiIiN54Fba8V2lWrlwJbW1teHl5IS8vD+7u7vjqq68qoylEREREREREGvVKEu3ffvtN9lpfXx/r1q3DunXrXsXhiYiIiIiIiF6ZCltHm4iIiIiIiOhtxESbiIiIiIiISIOYaBMRERERERFpEBNtIiIiIiIiIg1iok1ERERERESkQUy0iYiIiIiIiDSIiTYRERERERGRBjHRJiIiIiIiItIgJtpEREREREREGsREm4iIiIiIiEiDmGgTERERERERaRATbSIiIiIiIiINYqJNREREREREpEFMtImIiIioXL788ktoaWlh0qRJUtmjR4/g5+eH6tWrw8jICF5eXkhPT5e978aNG+jevTuqVq0Kc3NzTJs2DY8fP5bF/Pbbb3j//fehp6eHd999F1u3bi12/HXr1qFevXrQ19eHs7MzTpw4URGnSUT00phoExEREVGZTp48iQ0bNqBp06ay8smTJ2Pv3r0IDw/H77//jtu3b6Nv375SfWFhIbp37478/HzExsZi27Zt2Lp1KwICAqSY1NRUdO/eHZ06dUJycjImTZqEkSNH4sCBA1LMjh07MGXKFMydOxenTp1Cs2bN4O7ujoyMjIo/eSKi58REm4iIiIhKdf/+ffj4+GDTpk2oVq2aVJ6dnY1vvvkGK1asQOfOneHo6IgtW7YgNjYW8fHxAICDBw/ijz/+wPfff4/mzZuja9euWLBgAdatW4f8/HwAQEhICGxtbbF8+XI0bNgQ/v7+6NevH1auXCkda8WKFRg1ahSGDx+ORo0aISQkBFWrVsXmzZtf7cUgIioHJtpEREREVCo/Pz90794drq6usvKkpCQUFBTIyu3t7VG3bl3ExcUBAOLi4uDg4AALCwspxt3dHTk5Obhw4YIU8+y+3d3dpX3k5+cjKSlJFqOtrQ1XV1cpRp28vDzk5OTINiKiV0GnshtARERERK+v7du349SpUzh58mSxurS0NCgUCpiamsrKLSwskJaWJsU8nWSr6lV1pcXk5OTg4cOHuHv3LgoLC9XGXLp0qcS2BwUFYd68eeU7USIiDeKINhERERGpdfPmTUycOBGhoaHQ19ev7OY8t1mzZiE7O1vabt68WdlNIqK3BBNtIiIiIlIrKSkJGRkZeP/996GjowMdHR38/vvvWL16NXR0dGBhYYH8/HxkZWXJ3peeng5LS0sAgKWlZbFZyFWvy4oxMTGBgYEBatSogSpVqqiNUe1DHT09PZiYmMg2IqJXgYk2EREREanVpUsXnDt3DsnJydLm5OQEHx8f6d+6uro4fPiw9J6UlBTcuHEDLi4uAAAXFxecO3dONjt4VFQUTExM0KhRIynm6X2oYlT7UCgUcHR0lMUUFRXh8OHDUgwR0euEiTYR0XOq7HVkiYheFWNjYzRp0kS2GRoaonr16mjSpAmUSiV8fX0xZcoUHD16FElJSRg+fDhcXFzQqlUrAICbmxsaNWqEjz/+GGfOnMGBAwcwe/Zs+Pn5QU9PDwDwySef4K+//sL06dNx6dIlfPXVV9i5cycmT54stWXKlCnYtGkTtm3bhosXL2Ls2LHIzc3F8OHDK+XaEBGVhpOhERE9h9LWkf3ll18QHh4OpVIJf39/9O3bF8ePHwfwf+vIWlpaIjY2Fv/73/8wZMgQ6Orq4osvvgDwf+vIfvLJJwgNDcXhw4cxcuRI1KpVC+7u7q/8XImIymPlypXQ1taGl5cX8vLy4O7ujq+++kqqr1KlCvbt24exY8fCxcUFhoaGGDp0KObPny/F2Nra4pdffsHkyZMRHByMOnXq4Ouvv5b1fQMHDsQ///yDgIAApKWloXnz5oiMjCw2QRoR0euAiTYRUTk9vY7swoULpXLVOrJhYWHo3LkzAGDLli1o2LAh4uPj0apVK2kd2UOHDsHCwgLNmzfHggULMGPGDAQGBkKhUMjWkQWAhg0b4tixY1i5cmWJiXZeXh7y8vKk11y6hogq2m+//SZ7ra+vj3Xr1mHdunUlvsfGxgb79+8vdb8dO3bE6dOnS43x9/eHv79/udtKRFRZeOs4EVE5VfY6suoEBQVBqVRKm7W19UufJxERERG9HCbaRETloFpHNigoqFjdq1pHVh0uXUNERET0+uGt40REZVCtIxsVFfXarSOrp6cnTSZERERERK8HjmgTEZXhdVlHloiIiIjeDEy0iYjK8LqsI0tEREREbwbeOk5EVAbVOrJPe3odWQDSOrJmZmYwMTHB+PHjS1xHdsmSJUhLS1O7juzatWsxffp0jBgxAkeOHMHOnTvxyy+/vNoTJiIiIqKXwkSbiEgDXtU6skRERET0+mOiTUT0AipzHVkiIiIier3xGW0iIiIiIiIiDWKiTURERERERKRBTLSJiIiIiIiINIiJNhEREREREZEGMdEmIiIiIiIi0iAm2kREREREREQaxESbiIiIiIiISIOYaBMRERERERFpEBNtIiIiIiIiIg1iok1ERERERESkQUy0iYiIiIiIiDSIiTYRERERERGRBjHRJiIiIiIiItIgjSfaQUFBaNmyJYyNjWFubg5PT0+kpKTIYh49egQ/Pz9Ur14dRkZG8PLyQnp6uqabQkRERERERPTKaTzR/v333+Hn54f4+HhERUWhoKAAbm5uyM3NlWImT56MvXv3Ijw8HL///jtu376Nvn37aropRERERERERK+cjqZ3GBkZKXu9detWmJubIykpCe3bt0d2dja++eYbhIWFoXPnzgCALVu2oGHDhoiPj0erVq003SQiIiIiIiKiV6bCn9HOzs4GAJiZmQEAkpKSUFBQAFdXVynG3t4edevWRVxcnNp95OXlIScnR7YRERERERERvY4qNNEuKirCpEmT0KZNGzRp0gQAkJaWBoVCAVNTU1mshYUF0tLS1O4nKCgISqVS2qytrSuy2UREREREREQvrEITbT8/P5w/fx7bt29/qf3MmjUL2dnZ0nbz5k0NtZCIiIiIiIhIszT+jLaKv78/9u3bh+joaNSpU0cqt7S0RH5+PrKysmSj2unp6bC0tFS7Lz09Pejp6VVUU4mIiIiIiIg0RuMj2kII+Pv7Y/fu3Thy5AhsbW1l9Y6OjtDV1cXhw4elspSUFNy4cQMuLi6abg4RERERERHRK6XxEW0/Pz+EhYXh559/hrGxsfTctVKphIGBAZRKJXx9fTFlyhSYmZnBxMQE48ePh4uLC2ccJyIiIiIiojeexhPt9evXAwA6duwoK9+yZQuGDRsGAFi5ciW0tbXh5eWFvLw8uLu746uvvtJ0U4iIiIiIiIheOY0n2kKIMmP09fWxbt06rFu3TtOHJyIiIiIiIqpUFb6ONhEREREREdHbhIk2EREREakVFBSEli1bwtjYGObm5vD09ERKSoos5tGjR/Dz80P16tVhZGQELy8vpKeny2Ju3LiB7t27o2rVqjA3N8e0adPw+PFjWcxvv/2G999/H3p6enj33XexdevWYu1Zt24d6tWrB319fTg7O+PEiRMaP2ciIk1gok1EREREav3+++/w8/NDfHw8oqKiUFBQADc3N+Tm5koxkydPxt69exEeHo7ff/8dt2/fRt++faX6wsJCdO/eHfn5+YiNjcW2bduwdetWBAQESDGpqano3r07OnXqhOTkZEyaNAkjR47EgQMHpJgdO3ZgypQpmDt3Lk6dOoVmzZrB3d0dGRkZr+ZiEBE9hwpbR5uIiIiI3myRkZGy11u3boW5uTmSkpLQvn17ZGdn45tvvkFYWBg6d+4M4MkEuA0bNkR8fDxatWqFgwcP4o8//sChQ4dgYWGB5s2bY8GCBZgxYwYCAwOhUCgQEhICW1tbLF++HADQsGFDHDt2DCtXroS7uzsAYMWKFRg1ahSGDx8OAAgJCcEvv/yCzZs3Y+bMma/wqhARlY0j2kRERERULtnZ2QAAMzMzAEBSUhIKCgrg6uoqxdjb26Nu3bqIi4sDAMTFxcHBwQEWFhZSjLu7O3JycnDhwgUp5ul9qGJU+8jPz0dSUpIsRltbG66urlKMOnl5ecjJyZFtRESvAhNtIqJyeN2eUyQietWKioowadIktGnTBk2aNAEApKWlQaFQwNTUVBZrYWGBtLQ0KebpJFtVr6orLSYnJwcPHz7EnTt3UFhYqDZGtQ91goKCoFQqpc3a2vr5T5yI6AUw0SYiKofX6TlFIqLK4Ofnh/Pnz2P79u2V3ZRymzVrFrKzs6Xt5s2bld0kInpL8BltIqJyeJ2eU3xaXl4e8vLypNe8LZKIKoK/vz/27duH6Oho1KlTRyq3tLREfn4+srKyZKPa6enpsLS0lGKenR1cdbfP0zHP3gGUnp4OExMTGBgYoEqVKqhSpYraGNU+1NHT04Oent7znzAR0UviiDYR0QuorOcUn8XbIomoIgkh4O/vj927d+PIkSOwtbWV1Ts6OkJXVxeHDx+WylJSUnDjxg24uLgAAFxcXHDu3DnZ7OBRUVEwMTFBo0aNpJin96GKUe1DoVDA0dFRFlNUVITDhw9LMURErxMm2kREz6kyn1N8Fm+LJKKK5Ofnh++//x5hYWEwNjZGWloa0tLSpP5IqVTC19cXU6ZMwdGjR5GUlIThw4fDxcUFrVq1AgC4ubmhUaNG+Pjjj3HmzBkcOHAAs2fPhp+fnzTa/Mknn+Cvv/7C9OnTcen/tXe/oVWWfxzHP8cT20x31NkfXctl+A+LTfLPOEQ1UhIpKfDBkNQ1owdyJupIY0886RNFUBOcS5KsJ0OHuARNbVhuiEVzJmiUaAQOc04I5jzRUTz378GPDo0t3Tm7tuu6b98vuCBPN+3zRfvA13Ou7ddftWfPHjU2NmrdunXpLDU1Nfrss8/05Zdf6pdfftGqVauUSCTS34UcAFzCR8cBIEP/3FM8c+aM7Sh8LBLAkKqvr5cklZeX93p9//79eu+99yRJO3fu1IgRI7RkyRIlk0ktXLhQe/bsST8bDod19OhRrVq1StFoVKNGjVJlZaU2b96cfmby5Mk6duyY1q1bp127dqmoqEj79u3rdWWmoqJCt27d0saNG9XZ2alZs2bpxIkTff5yEgBcwKINABmwfU8RAIaT53kPfSYvL091dXWqq6v7z2eKi4v19ddfP/C/U15erp9++umBz1RXV6u6uvqhmQDANj46DgAD4Mo9RQAAALiPd7QBYABisZgaGhp05MiR9D1F6f/3E0eOHNnrnmJBQYEikYhWr179n/cUt23bps7Ozn7vKe7evVsbNmzQypUr9e2336qxsVHHjh2zNjsAAAAywzvaADAA9fX16u7uVnl5uSZOnJg+Bw8eTD+zc+dOvfXWW1qyZIleffVVTZgwQYcPH07/+3/uKYbDYUWjUS1btkwrVqzo955ic3OzSktLtX379j73FAEAAOA23tEGgAFw7Z4iAAAA3MU72gAAAAAAGMSiDQAAAACAQSzaAAAAAAAYxKINAAAAAIBBLNoAAAAAABjEog0AAAAAgEEs2gAAAAAAGMSiDQAAAACAQSzaAAAAAAAYxKINAAAAAIBBLNoAAAAAABjEog0AAAAAgEEs2gAAAAAAGMSiDQAAAACAQSzaAAAAAAAYxKINAAAAAIBBLNoAAAAAABjEog0AAAAAgEEs2gAAAAAAGMSiDQAAAACAQSzaAAAAAAAYxKINAAAAAIBBLNoAAAAAABjEog0AAAAAgEEs2gAAAAAAGMSiDQAAAACAQSzaAAAAAAAYxKINAAAAAIBBVhfturo6Pffcc8rLy1NZWZl+/PFHm3EAwBn0IwD0j34E4AfWFu2DBw+qpqZG8Xhc58+fV2lpqRYuXKiuri5bkQDACfQjAPSPfgTgF9YW7R07duiDDz5QVVWVZs6cqU8//VSPP/64Pv/8c1uRAMAJ9CMA9I9+BOAXj9n4onfv3lV7e7tqa2vTr40YMUILFizQ999/3+f5ZDKpZDKZ/nV3d7ck6fbt2xl93VTyrywTD5+BzuSHWSTmcVkm//8EbZ5/P+953lDEyRr9+N8e1T+zfphFYh6XPar9KJnpyKD9HjPP8AtSn0iP7jz/fnZA/ehZcP36dU+Sd/bs2V6vr1+/3ps3b16f5+PxuCeJw+FwjJ+Ojo7hqr4BoR85HI4rx+/96Hl0JIfDGZozkH608o52pmpra1VTU5P+dSqV0p9//qnx48crFApZyXT79m09++yz6ujoUCQSsZLBJOZxW5DmcWUWz/PU09OjwsJCaxlMoB+HHvO4jXnMC0o/SnTkcAjSPEGaRWKeoZBJP1pZtJ944gmFw2HdvHmz1+s3b97UhAkT+jyfm5ur3NzcXq+NHTt2KCMOWCQSCcQf3H8wj9uCNI8Ls4wZM8bq1+8P/egu5nEb85gVhH6U6MjhFKR5gjSLxDymDbQfrXwztJycHM2ePVunTp1Kv5ZKpXTq1ClFo1EbkQDACfQjAPSPfgTgJ9Y+Ol5TU6PKykrNmTNH8+bN0yeffKJEIqGqqipbkQDACfQjAPSPfgTgF9YW7YqKCt26dUsbN25UZ2enZs2apRMnTujpp5+2FSkjubm5isfjfT6O5FfM47YgzROkWYYK/egW5nEb8zxa/N6PUvB+j4M0T5BmkZjHtpDnOfazGwAAAAAA8DErd7QBAAAAAAgqFm0AAAAAAAxi0QYAAAAAwCAWbQAAAAAADGLRBgAAAADAIBbtDLW2tmrx4sUqLCxUKBTSV199ZTvSoGzZskVz585Vfn6+nnrqKb3zzju6fPmy7VhZq6+vV0lJiSKRiCKRiKLRqI4fP247lhFbt25VKBTS2rVrbUfJyscff6xQKNTrzJgxw3YsGBakjqQf/YN+hB/Qj+6iH93m145k0c5QIpFQaWmp6urqbEcxoqWlRbFYTD/88IOam5t17949vfHGG0okErajZaWoqEhbt25Ve3u7zp07p9dff11vv/22fv75Z9vRBqWtrU179+5VSUmJ7SiD8sILL+jGjRvpc+bMGduRYFiQOpJ+9Af6EX5BP7qLfnSfLzvSQ9YkeU1NTbZjGNXV1eVJ8lpaWmxHMWbcuHHevn37bMfIWk9Pjzd16lSvubnZe+2117w1a9bYjpSVeDzulZaW2o6BYRS0jqQf3UM/wq/oR/fRj+7wa0fyjjZ66e7uliQVFBRYTjJ49+/f14EDB5RIJBSNRm3HyVosFtObb76pBQsW2I4yaFeuXFFhYaGef/55vfvuu7p27ZrtSMCA0Y/uoR8BN9CP7glSP0r+7MjHbAeAO1KplNauXauXX35ZL774ou04Wbt48aKi0aj+/vtvjR49Wk1NTZo5c6btWFk5cOCAzp8/r7a2NttRBq2srExffPGFpk+frhs3bmjTpk165ZVXdOnSJeXn59uOBzwQ/ege+hFwA/3oniD1o+TfjmTRRlosFtOlS5f8cefhAaZPn64LFy6ou7tbhw4dUmVlpVpaWnxXlh0dHVqzZo2am5uVl5dnO86gLVq0KP3PJSUlKisrU3FxsRobG/X+++9bTAY8HP3oFvoRcAf96Jag9aPk345k0YYkqbq6WkePHlVra6uKiopsxxmUnJwcTZkyRZI0e/ZstbW1adeuXdq7d6/lZJlpb29XV1eXXnrppfRr9+/fV2trq3bv3q1kMqlwOGwx4eCMHTtW06ZN09WrV21HAR6IfnQP/Qi4gX50T9D7UfJPR7JoP+I8z9Pq1avV1NSk06dPa/LkybYjGZdKpZRMJm3HyNj8+fN18eLFXq9VVVVpxowZ+uijj3xfknfu3NFvv/2m5cuX244C9It+dBf9CNhFP7or6P0o+acjWbQzdOfOnV5/e/L777/rwoULKigo0KRJkywmy04sFlNDQ4OOHDmi/Px8dXZ2SpLGjBmjkSNHWk6XudraWi1atEiTJk1ST0+PGhoadPr0aZ08edJ2tIzl5+f3ues0atQojR8/3pd3oD788EMtXrxYxcXF+uOPPxSPxxUOh7V06VLb0WBQkDqSfnQX/Qg/oh/dRT+6zbcdafm7nvvOd99950nqcyorK21Hy0p/s0jy9u/fbztaVlauXOkVFxd7OTk53pNPPunNnz/f++abb2zHMsbPP56hoqLCmzhxopeTk+M988wzXkVFhXf16lXbsWBYkDqSfvQX+hGuox/dRT+6za8dGfI8zxvybR4AAAAAgEcEP0cbAAAAAACDWLQBAAAAADCIRRsAAAAAAINYtAEAAAAAMIhFGwAAAAAAg1i0AQAAAAAwiEUbAAAAAACDWLQBAAAAADCIRRsAAAAAAINYtAEAAAAAMIhFGwAAAAAAg/4HbdqHKQ9Vls0AAAAASUVORK5CYII=",
      "text/plain": [
       "<Figure size 1200x400 with 3 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "fig, axes = plot_scenario(scenario)\n",
    "fig.suptitle('Heuristic harvest flows (area, volume, growing stock)')\n",
    "plt.show()\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "2120f892",
   "metadata": {},
   "source": [
    "### Ensure transition lookups cover scheduled ages\n",
    "\n",
    "Some Woodstock models only declare a generic transition for `age = -1`. The snippet below mirrors what `spades_ws3` does: it copies those generic transitions to every `(action, age)` pair actually referenced in the schedule, so that `ForestRaster` can resolve targets without raising a `KeyError`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "id": "39a182a9",
   "metadata": {},
   "outputs": [],
   "source": [
    "for period in fm.periods:\n",
    "    for acode, by_dtype in fm.applied_actions[period].items():\n",
    "        for dtk, by_age in by_dtype.items():\n",
    "            dt = fm.dtypes[dtk]\n",
    "            for age in list(by_age.keys()):\n",
    "                if (acode, age) not in dt.transitions:\n",
    "                    if (acode, -1) in dt.transitions:\n",
    "                        dt.transitions[acode, age] = dt.transitions[acode, -1]\n",
    "                    else:\n",
    "                        raise KeyError(f\"Missing transition for {(acode, age)} in {dtk}\")\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "aa658c71",
   "metadata": {},
   "source": [
    "## Rasterize the stand inventory\n",
    "\n",
    "`ForestRaster` expects a three-band raster: band 1 stores hashed development type theme values, band 2 stores age (years), and band 3 stores block IDs. The helper `ws3.common.rasterize_stands` creates this raster from the sample shapefile and returns the required hash map."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "id": "3ad04cf9",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Raster shape: 3 bands 41 x 41\n",
      "CRS: EPSG:3005\n",
      "Pixel size: 100.0 m\n",
      "Hash map entries: 12\n"
     ]
    }
   ],
   "source": [
    "theme_columns = ['theme0', 'theme1', 'theme2', 'theme3', 'curve1']\n",
    "\n",
    "hdt_map = rasterize_stands(\n",
    "    str(SHAPE_PATH),\n",
    "    str(INVENTORY_TIF),\n",
    "    theme_cols=theme_columns,\n",
    "    age_col='age',\n",
    "    blk_col='curve2',\n",
    "    age_divisor=1.0,\n",
    ")\n",
    "\n",
    "with rasterio.open(INVENTORY_TIF) as src:\n",
    "    print('Raster shape:', src.count, 'bands', src.width, 'x', src.height)\n",
    "    print('CRS:', src.crs)\n",
    "    print('Pixel size:', src.transform.a, 'm')\n",
    "    print('Hash map entries:', len(hdt_map))\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "3b806574",
   "metadata": {},
   "source": [
    "## Allocate the schedule to rasters\n",
    "\n",
    "We instantiate `ForestRaster` with the model, hash metadata, and destination directory. The planner horizon of five periods (10-year periods) becomes 50 annual GeoTIFFs. Each file encodes whether a cell is disturbed (`1`) or untouched (`nodata`/`0`) for the target year."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "id": "a3097890",
   "metadata": {},
   "outputs": [],
   "source": [
    "# clean older outputs from previous runs\n",
    "for tif in SPATIAL_DIR.glob('harvested_*.tif'):\n",
    "    tif.unlink()\n",
    "\n",
    "forest_raster = ForestRaster(\n",
    "    hdt_map=hdt_map,\n",
    "    hdt_func=hash_dt,\n",
    "    src_path=str(INVENTORY_TIF),\n",
    "    snk_path=str(SPATIAL_DIR),\n",
    "    acode_map={'harvest': 'harvested'},\n",
    "    forestmodel=fm,\n",
    "    base_year=BASE_YEAR,\n",
    "    horizon=HORIZON,\n",
    "    period_length=PERIOD_LENGTH,\n",
    "    time_step=1,\n",
    "    piggyback_acodes={},\n",
    ")\n",
    "\n",
    "forest_raster.allocate_schedule(verbose=False, sda_mode='randblk', nthresh=10)\n",
    "forest_raster.cleanup()\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "d8649bff",
   "metadata": {},
   "source": [
    "## Inspect the outputs\n",
    "\n",
    "The allocation step writes one GeoTIFF per year. Below we list the first few files, visualise the first raster, and summarise the disturbed area (in number of pixels) per year."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "id": "6470ac1b",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "[PosixPath('data/spatial_outputs/tsa24_demo/harvested_2020.tif'),\n",
       " PosixPath('data/spatial_outputs/tsa24_demo/harvested_2021.tif'),\n",
       " PosixPath('data/spatial_outputs/tsa24_demo/harvested_2022.tif'),\n",
       " PosixPath('data/spatial_outputs/tsa24_demo/harvested_2023.tif'),\n",
       " PosixPath('data/spatial_outputs/tsa24_demo/harvested_2024.tif')]"
      ]
     },
     "execution_count": 11,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "output_files = sorted(SPATIAL_DIR.glob('harvested_*.tif'))\n",
    "output_files[:5]\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "id": "3f0f764f",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAeEAAAH4CAYAAAB9k1VdAAAAOnRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjEwLjYsIGh0dHBzOi8vbWF0cGxvdGxpYi5vcmcvq6yFwwAAAAlwSFlzAAAPYQAAD2EBqD+naQAAPZpJREFUeJzt3Xlc1VX+x/E3qAgqmwLmhhpihrmUg1SWGrZNrlmmiEuamaUt44Brlstk2aJOqf3SUmdSSyVsFFMsp2xxiZrJZcSktNRMkbFwF5Lz+8PhjldA7nW4nou+no+Hj4ece+75nu/d3pzv93s/+BhjjAAAwCXna3sCAABcqQhhAAAsIYQBALCEEAYAwBJCGAAASwhhAAAsIYQBALCEEAYAwBJCGAAASwhhlIn58+fLx8dHP/zwg6Otffv2at++vbU5FefBBx9UtWrVbE/DYxo0aKBOnTrZnoZXKSgo0HXXXafnnnvO9lSuOP/+979VtWpVffDBB7an4rUI4RIUhspXX31V7O3t27fXddddd4lnVfZOnDih8ePH65NPPrE9FUnS+vXrNX78eP3666+2p4LLxDvvvKO9e/dq2LBhjraMjAwNGzZMTZs2VdWqVRUZGakHHnhAO3fuLHaMzMxM3X333apWrZqqV6+uvn376tChQ059duzYoREjRqhly5YKDAxUrVq11LFjxxI/Q3766Sc98MADCgkJUVBQkLp27apdu3a5tE+e2NbevXs1YcIEtW7dWqGhoQoLC1P79u310UcflTqfhx9+WD4+PkV+AaxRo4YGDRqkcePGubRfVySDYs2bN89IMhkZGcXe3q5dO9O0adNLPKuyd+jQISPJPPvss//TOIWP1+7dux1tp0+fNqdPn3ZrnJdeeqnIOGWpf//+pmrVqh4Z2xvUr1/fdOzY0fY0vEqLFi3M4MGDndruu+8+c9VVV5nHH3/czJkzx0yaNMnUrFnTVK1a1WzdutWp7969e01YWJiJiooyf/7zn81zzz1nQkNDTYsWLZxe33/84x9NSEiIeeihh8wbb7xhXnzxRRMVFWUqVKhgPvzwQ6cxjx49aqKjo01ERISZMmWKmTp1qqlXr56pW7euycnJKXWfPLGt1157zQQEBJiEhAQzY8YMM336dHPDDTcYSWbu3LklziUjI8NUrFjR+Pv7F/va2759u5Fk1q5dW+p+XYkI4RJcyhAuKCgwJ06cKJOx3OXJEL4YngrhY8eOGWMI4cvNmTNnzMmTJ0u8/R//+IeRZD766COn9i+++KLIL4g7d+40lStXNomJiU7tjz76qAkICDA//vijo+3DDz80kswbb7zhaPvqq6/M0aNHne6bk5NjwsPDTZs2bZzap0yZYiSZL7/80tGWmZlpKlSoYEaPHl3KXntmW9u2bTOHDh1yuu+pU6dMkyZNTN26dYudR0FBgbnpppvMwIEDL/jau+6660zfvn1L3a8rESFcgosJ4blz55rbbrvNhIeHGz8/P3PttdeaWbNmFblv4Yt19erVplWrVqZy5cpm2rRppmnTpqZ9+/ZF+p85c8bUrl3b3HfffU5t06ZNMzExMaZy5comIiLCDB482Bw+fNjpvhkZGebOO+80NWrUMP7+/qZBgwZmwIABxhhjdu/ebSQV+VdaIG/bts3cdtttxt/f39SpU8dMmjTJvPXWW0XCs127dqZdu3ZO93311VdNTEyMCQgIMCEhIaZVq1Zm4cKFxhhjnn322WLns3v3bsdc582bV2Q+58+5cJx//etfJiEhwYSEhJiWLVsaY/4bwt9//7258847TZUqVUytWrXMhAkTTEFBgdO4L730krnppptM9erVjb+/v7nhhhvM0qVLi93+0KFDzbJly0zTpk2Nn5+fiYmJMatWrSrSd9++fWbgwIGmVq1axs/PzzRo0MAMGTLEKRB++eUX8+STT5q6desaPz8/ExUVZV544QVz5syZCz4vxvz3tZWenm5atGhhKleubK699lrz3nvvOfp8//33RpKZOnVqkft/8cUXRpJZtGhRseMfPXrUVKlSxTzxxBNFbtu7d6/x9fU1kydPdntf3H2sFyxYYGJiYkzFihXNsmXLSnw8nnnmGePn52fy8vJK7HOuG264wdxwww1ObREREaZHjx5F+jZu3Nh06NCh1DG7d+9uqlev7tQWGxtrYmNji/S98847TVRUlEtzvVTbGj58uJFkjhw5UuS2v/zlLyYwMND8/PPPFwzhP/zhDyYkJKTIewzGVPT04e7yLjc3Vzk5OUXa8/Pzi7S9/vrratq0qbp06aKKFStqxYoVeuyxx1RQUKChQ4c69f3222+VkJCgRx55RA8//LCuueYa9ezZU+PHj9eBAwd01VVXOfp+/vnn2r9/v3r16uVoe+SRRzR//nwNGDBATzzxhHbv3q0ZM2bon//8p7744gtVqlRJ2dnZuvPOOxUeHq5Ro0YpJCREP/zwg1JTUyVJ4eHhev311/Xoo4/q3nvvVffu3SVJzZs3L/HxOHDggG677Tb99ttvGjVqlKpWrarZs2crICCg1Mdyzpw5euKJJ3T//ffrySef1KlTp7RlyxZt2rRJvXv3Vvfu3bVz50698847mjZtmsLCwhzzPP/8myt69Oih6OhoTZ48Weacv9h55swZ3X333brxxhv14osvavXq1Xr22Wf122+/aeLEiY5+f/7zn9WlSxclJiYqLy9P7777rnr06KG0tDR17NjRaVuff/65UlNT9dhjjykwMFCvvvqq7rvvPu3Zs0c1atSQJO3fv1+tW7fWr7/+qsGDB6tJkyb66aeflJKSohMnTsjPz08nTpxQu3bt9NNPP+mRRx5RZGSk1q9fr9GjR+vnn3/W9OnTS93vrKws9ezZU0OGDFH//v01b9489ejRQ6tXr9Ydd9yhq6++Wm3atNHChQv1hz/8wem+CxcuVGBgoLp27Vrs2NWqVdO9996rxYsXa+rUqapQoYLjtnfeeUfGGCUmJkqSW/vizmP997//XUuWLNGwYcMUFhamBg0alPhYrF+/Xtddd50qVapU6uNmjNHBgwfVtGlTR9tPP/2k7Oxs/e53vyvSv3Xr1i5dcHTgwAHHa1k6e6HYli1bNHDgwGLHXLNmjY4eParAwMBSx74U2zpw4ICqVKmiKlWqOLUfPXpUI0eO1JgxY5w+r4rTqlUrTZs2Tf/6178ui2tpypTlXwK8VuFK+EL/zl8JF3dI+a677jJXX321U1v9+vWNJLN69Wqn9m+//dZIMq+99ppT+2OPPWaqVavmGP+zzz4zkhwryEKrV692al+2bNkFV/PGuH84+qmnnjKSzKZNmxxt2dnZJjg4uNSVcNeuXUs9hF/S4eiLWQknJCQU6du/f38jyTz++OOOtoKCAtOxY0fj5+fndDju/OczLy/PXHfddSY+Pr7I9v38/Mx3333naNu8eXOR57Jfv37G19e32OejcIUwadIkU7VqVbNz506n20eNGmUqVKhg9uzZU+S+5yp8bZ278s3NzTW1atUy119/vaPtjTfeMJJMZmam0/6FhYWZ/v37X3Ab6enpRlKRlX7z5s2dnm939sWdx9rX19f861//uuAcC9WtW9fpCNKFvP3220aSeeuttxxtGRkZRpL561//WqR/cnKykWROnTpV4piffvqp8fHxMePGjXO0Fb7nJk6cWKT/zJkzjSSzY8cOl+bs6W1lZWUZf3//Yg8lJyUlmYYNGzr2/0Ir4fXr1xtJZvHixe7u1mWPq6NLMXPmTH344YdF/hW3Wjx3NVi4gm7Xrp127dql3Nxcp74NGzbUXXfd5dTWuHFjtWzZUosXL3a0nTlzRikpKercubNj/KVLlyo4OFh33HGHcnJyHP9atWqlatWq6eOPP5YkhYSESJLS0tKKXblfjA8++EA33nijWrdu7WgLDw93rH4uJCQkRPv27VNGRkaZzKU0Q4YMKfG2c6+U9fHx0bBhw5SXl+d0Jei5z+cvv/yi3Nxc3XrrrfrHP/5RZLzbb79dUVFRjp+bN2+uoKAgxxWoBQUFev/999W5c+diV1U+Pj6Szj63t956q0JDQ52e29tvv11nzpzRp59+Wup+165dW/fee6/j56CgIPXr10///Oc/deDAAUnSAw88IH9/fy1cuNDRLz09XTk5OerTp88Fx7/99ttVu3Ztp/tu27ZNW7ZscbqvO/vizmPdrl07xcTElPo4SGe/IhMaGlpqvx07dmjo0KG66aab1L9/f0f7yZMnJUmVK1cuch9/f3+nPufLzs5W79691bBhQ40YMaJMxiyJJ7Z14sQJ9ejRQwEBAXrhhRecbtu5c6f+/Oc/66WXXip27PMVPgfFHVW80nE4uhStW7cu9kOz8IPlXF988YWeffZZbdiwQSdOnHC6LTc3V8HBwY6fGzZsWOz2evbsqTFjxuinn35SnTp19Mknnyg7O1s9e/Z09MnKylJubq4iIiKKHSM7O1vS2Q+r++67TxMmTNC0adPUvn17devWTb1793bpjVOcH3/8UXFxcUXar7nmmlLvO3LkSH300Udq3bq1GjVqpDvvvFO9e/dWmzZtLmoupSnpMfb19dXVV1/t1Na4cWNJcvqec1pamv70pz/pm2++0enTpx3thYF5rsjIyCJtoaGh+uWXXyRJhw4d0pEjR0o9FJeVlaUtW7YoPDy82NsLn9sLadSoUZE5nrt/V111lUJCQtS5c2ctWrRIkyZNknT2UHSdOnUUHx9/wfF9fX2VmJio119/XSdOnFCVKlW0cOFC+fv7q0ePHhe1L+481iU9ryUx55yKKM6BAwfUsWNHBQcHKyUlxekQe+EvB+fOqdCpU6ec+pzr+PHj6tSpk44eParPP//c6bvp7oxZ+EtToeDg4CLbK6ttnevMmTPq1auXtm/frlWrVql27dpOtz/55JO6+eabdd999xW5b3EKn4Pins8rHSFcRr7//nt16NBBTZo00dSpU1WvXj35+fnpgw8+0LRp01RQUODUv6RzqD179tTo0aO1dOlSPfXUU1qyZImCg4N19913O/oUFBQoIiLCaSVyrsIPPR8fH6WkpGjjxo1asWKF0tPTNXDgQL3yyivauHHjJS9ace211+rbb79VWlqaVq9erffee0+zZs3SM888owkTJlzwviW9ec+cOVPifVw5T12Szz77TF26dFHbtm01a9Ys1apVS5UqVdK8efO0aNGiIv3P/eA+V2kBcL6CggLdcccdTquZcxWGaVno16+fli5dqvXr16tZs2Zavny5HnvsMfn6ln6ArF+/fnrppZf0/vvvKyEhQYsWLVKnTp2cftF0dV/cfazdeV5r1Kjh+EWoOLm5ufr973+vX3/9VZ999lmRsKlVq5Yk6eeffy5y359//lnVq1cv8gttXl6eunfvri1btig9Pb3IL16F9ylpTEmOeRRuv9C8efP04IMPemRb53r44YeVlpamhQsXFvml7O9//7tWr16t1NRUp19af/vtN508eVI//PCDqlevrqCgIMdthc/BueercRYhXEZWrFih06dPa/ny5U6rosJDw65q2LChWrdurcWLF2vYsGFKTU1Vt27dnN7oUVFR+uijj9SmTRuXPpBuvPFG3XjjjXruuee0aNEiJSYm6t1339WgQYPc/s20fv36ysrKKtL+7bffunT/qlWrqmfPnurZs6fjA+S5557T6NGj5e/vX+J8Cg9nnV/E48cff3Rr/tLZcNi1a5dToBUWaSi8yOe9996Tv7+/0tPTnR77efPmub096ewvRkFBQdq2bdsF+0VFRenYsWO6/fbbL2o7kvTdd9/JGOP0WJ6/f5J09913Kzw8XAsXLlRcXJxOnDihvn37urSN6667Ttdff70WLlyounXras+ePXrttdcual/K+rE+V5MmTbR79+5ibzt16pQ6d+6snTt36qOPPir2EHedOnUUHh5ebBGML7/8Ui1btnRqKygoUL9+/bR27VotWbJE7dq1K3I/X19fNWvWrNgxN23apKuvvtpxodSHH37odPu5F42V9bYKJScna968eZo+fboSEhKK3G/Pnj2S5LiQ81w//fSTGjZsqGnTpumpp55ytBc+B9dee22R+1zpOCdcRgpXQueufHJzcy/qg6Rnz57auHGj5s6dq5ycHKdD0dLZ83lnzpxxHEY812+//eYIql9++aXISqzwQ6Pw8FThFY+uVqi65557tHHjRn355ZeOtkOHDpW4Kj/Xv//9b6ef/fz8FBMTI2OM45x11apVi51PUFCQwsLCipwTnTVrlkvzPt+MGTMc/zfGaMaMGapUqZI6dOgg6ezz6ePj47TS/uGHH/T+++9f1PZ8fX3VrVs3rVixotgPxMLn6YEHHtCGDRuUnp5epM+vv/6q3377rdRt7d+/X8uWLXP8fOTIEf31r39Vy5Ytna5irVixohISErRkyRLNnz9fzZo1u+CV8efr27ev1qxZo+nTp6tGjRr6/e9/73S7q/tS1o/1uW666SZt27atyOHYM2fOqGfPntqwYYOWLl2qm266qcQx7rvvPqWlpWnv3r2OtrVr12rnzp1Oh98l6fHHH9fixYs1a9asYkOq0P3336+MjAyn18K3336rv//9705j3n777U7/zl0Zl/W2JOmll17Syy+/rDFjxujJJ58sdrz4+HgtW7asyL/w8HD97ne/07Jly9S5c2en+3z99dcKDg52+iUC/2HtkjAv5+73hHfs2GH8/PxMs2bNzIwZM8wLL7xgoqKiTIsWLYpc7VtaQYW9e/caHx8fExgYaKpXr17sdxwfeeQRI8n8/ve/N9OmTTMzZswwTz75pKldu7bj+5XTpk0z0dHRZsSIEeaNN94wL7/8srnmmmtMUFCQ2bVrl2OsmJgYc9VVV5mZM2ead955p0jFoHPt37/f1KhRw4SGhprx48ebl156yURHR5vmzZuXenX0DTfcYO655x7z3HPPmTfffNP88Y9/NJUrVzadO3d29Pnyyy+NJHPPPfeYv/71r+add95xFNoYNWqUkWQeeugh8/rrr5uEhATTqlWrEq+OPr/wgDFnr4729/c30dHRpl+/fmbmzJmmU6dORpIZM2aMo9/atWuNJHPrrbea119/3UyYMMFEREQ49vNc+s93V89Xv359pyuN9+3bZ6666ipTpUoV89RTT5k33njDjB8/3jRt2tT88ssvxhhjjh8/bm644QZTsWJFM2jQIPP666+bl19+2fH95uL26fxtNm7c2ISEhJhRo0aZadOmmWbNmhlfX98iV+Mbc7bog/5ztf+UKVMuOPb5Dhw4YCpWrGgkmUcffbTI7a7uS1k81iUp3L/09HSn9ieffNJIMp07dzZvv/12kX/n2rNnj6lRo4aJiooyr776qpk8ebIJDQ01zZo1c7oyetq0aUaSuemmm4ods/B1bIwxR44cMVFRUSYiIsK8+OKLZtq0aaZevXqmdu3aJjs7u9T98sS2UlNTjSQTHR1d7JgHDhy44JxKK9bRp0+fUvfrSkQIl+BiinUsX77cNG/e3FEUY8qUKWbu3Lluh7AxxrRp08ZIMoMGDSqxz+zZs02rVq1MQECACQwMNM2aNTMjRoww+/fvN8acrRaUkJBgIiMjHQU9OnXqZL766iuncdavX29atWpl/Pz8XPq60pYtW0y7du3cLtbxxhtvmLZt25oaNWqYypUrm6ioKJOcnGxyc3Odxp80aZKpU6eO8fX1dRrzxIkT5qGHHjLBwcEmMDDQPPDAAyY7O9vtED6/WEfNmjXNs88+W6SAxFtvvWWio6NN5cqVTZMmTcy8efMcY5/L1RA2xpgff/zR9OvXz4SHh5vKlSubq6++2gwdOtSpWMfRo0fN6NGjTaNGjYyfn58JCwszN998s3n55ZdLLTpxbrGO5s2bO+ZeXOGLQk2bNjW+vr5m3759Fxy7OPfcc4+RZNavX1/s7a7uy//6WF9I8+bNzUMPPeTU1q5duwt+/fB827Ztc7xeQkJCTGJiYpFQKvz6W0n/zv/a3d69e839999vgoKCTLVq1UynTp1MVlaWS/vkiW2VVCyn8N/HH398wTmV9LmWmZlppKJVy3CWjzFuXjkC4LJy/fXXq3r16lq7dq3b97333nu1detWfffddx6YWdl4++23NXToUO3Zs8fxtT1cOk899ZQ+/fRTff3111wdXQzOCQNXsK+++krffPON+vXr5/Z9f/75Z61cudLli7lsSUxMVGRkpGbOnGl7Klecf//733rzzTf1pz/9iQAuASth4Aq0bds2ff3113rllVeUk5OjXbt2OYo3lGb37t364osv9OabbyojI0Pff/99qWULARSPlTBwBUpJSdGAAQOUn5+vd955x+UAlqR169apb9++2r17t/7yl78QwMD/gJUwAACWsBIGAMASQhgAAEsIYQAALCGEAQCwhBAGAMASQhgAAEsIYQAALCGEAQCwhBAGAMASQhgAAEsIYQAALCGEAQCwhBAGAMASQhgAAEsqutrx119/9dgkQkJCPDZ2eePJxxkoC+68X3k9X1nc/Swvj6+Pss4rVsIAAFhCCAMAYAkhDACAJYQwAACWEMIAAFhCCAMAYAkhDACAJYQwAACWEMIAAFhCCAMAYAkhDACAJS7Xji6PymNdUsDb8b5CWfHU3w0oT69RVsIAAFhCCAMAYAkhDACAJYQwAACWEMIAAFhCCAMAYAkhDACAJYQwAACWEMIAAFhCCAMAYInLZSvdKS/mbsmw8lRiDLDFUyX+JN6DVxpPvpbgHlbCAABYQggDAGAJIQwAgCWEMAAAlhDCAABYQggDAGAJIQwAgCWEMAAAlhDCAABYQggDAGCJy2UrUf55svQoLg7PCWy43F9L5el9xUoYAABLCGEAACwhhAEAsIQQBgDAEkIYAABLCGEAACwhhAEAsIQQBgDAEkIYAABLCGEAACxxuWyl7dJeKMqd0mzwTryvLg1vea+483y7O2deSxfH9muDlTAAAJYQwgAAWEIIAwBgCSEMAIAlhDAAAJYQwgAAWEIIAwBgCSEMAIAlhDAAAJYQwgAAWEIIAwBgicu1o1H+UVsWKD94vzqzXePZU1gJAwBgCSEMAIAlhDAAAJYQwgAAWEIIAwBgCSEMAIAlhDAAAJYQwgAAWEIIAwBgCSEMAIAlLpet3HrtMJcHbZY546ImA6DsuFPm73IvkejO/l2u5RHLu8v1OWQlDACAJYQwAACWEMIAAFhCCAMAYAkhDACAJYQwAACWEMIAAFhCCAMAYAkhDACAJYQwAACW+BhjjCsdL/eydsDlxlOl+/gsgLfzZMnWsn5fsRIGAMASQhgAAEsIYQAALCGEAQCwhBAGAMASQhgAAEsIYQAALCGEAQCwhBAGAMASQhgAAEsqutox5vm9Lg+6fXS9i5oMJL9XQl3um/fHXzw4EwC4/HmqvKurWAkDAGAJIQwAgCWEMAAAlhDCAABYQggDAGAJIQwAgCWEMAAAlhDCAABYQggDAGAJIQwAgCWEMAAAlvgYY4wrHX/99VcPTwVXAk/Wab3cX6Pe8ti5M4/L/TnBlaes34eshAEAsIQQBgDAEkIYAABLCGEAACwhhAEAsIQQBgDAEkIYAABLCGEAACwhhAEAsIQQBgDAkoq2JwBcyOVe9tBbSlF6w7jAlYiVMAAAlhDCAABYQggDAGAJIQwAgCWEMAAAlhDCAABYQggDAGAJIQwAgCWEMAAAlhDCAABYQtlKXFKUPLx4PHbeyVOlR3m+rwyshAEAsIQQBgDAEkIYAABLCGEAACwhhAEAsIQQBgDAEkIYAABLCGEAACwhhAEAsIQQBgDAEspWAi5wpzQh5QavLDzf+F+wEgYAwBJCGAAASwhhAAAsIYQBALCEEAYAwBJCGAAASwhhAAAsIYQBALCEEAYAwBJCGAAASwhhAAAsoXY04ALqA/8XdbSvLDzfnsVKGAAASwhhAAAsIYQBALCEEAYAwBJCGAAASwhhAAAsIYQBALCEEAYAwBJCGAAASwhhAAAsoWwlALe4U5rQnZKH7o7tDkovXhre8nyXJ6yEAQCwhBAGAMASQhgAAEsIYQAALCGEAQCwhBAGAMASQhgAAEsIYQAALCGEAQCwhBAGAMASylYCZSzm+b0u990/JcTlvpQEvHg8FhePx86zWAkDAGAJIQwAgCWEMAAAlhDCAABYQggDAGAJIQwAgCWEMAAAlhDCAABYQggDAGAJIQwAgCWEMAAAlvgYY4wrHb2lfqg79XO9Zc64eH6vhLrcN++Pv3hwJvZRO7r8c/c5dAfP96VR1s8hK2EAACwhhAEAsIQQBgDAEkIYAABLCGEAACwhhAEAsIQQBgDAEkIYAABLCGEAACwhhAEAsKSi7QngyhLz/F63+m8ffXmXokT556lSlJShvDKwEgYAwBJCGAAASwhhAAAsIYQBALCEEAYAwBJCGAAASwhhAAAsIYQBALCEEAYAwBJCGAAASyhb+R/ulp4rjyXltl47zOW+zTJneGQO20fX88i4QFnxVBlKd10Jn0lgJQwAgDWEMAAAlhDCAABYQggDAGAJIQwAgCWEMAAAlhDCAABYQggDAGAJIQwAgCWEMAAAlvgYY4wrHSmJBm/3TN/uLved+HaqB2fiGZQxvDS8pWwlz593KuvXBythAAAsIYQBALCEEAYAwBJCGAAASwhhAAAsIYQBALCEEAYAwBJCGAAASwhhAAAsIYQBALCEEAYAwBJqR5dj1BIu/2Ke3+ty3+2j67k1tjuvj8v9teHJetCX+2MHZ9SOBgDgMkEIAwBgCSEMAIAlhDAAAJYQwgAAWEIIAwBgCSEMAIAlhDAAAJYQwgAAWEIIAwBgSUXbE/AWV0IJyPJYxrDd1A4u9103fK1HxpWkDeYfLvcd9c1tLvfd/naqW/OA53nytV8e34PwLFbCAABYQggDAGAJIQwAgCWEMAAAlhDCAABYQggDAGAJIQwAgCWEMAAAlhDCAABYQggDAGAJZSv/40ooEVce99GdUpSe5E4pyhdafuxy37yLmYyLao/c6nLf/VOaudy3PL6O3OHJEraX+2PnLcpTeVBWwgAAWEIIAwBgCSEMAIAlhDAAAJYQwgAAWEIIAwBgCSEMAIAlhDAAAJYQwgAAWEIIAwBgiY8xxrjS0XZpL1xaW68d5nLfZpkzPDaPZ/p2d7nvxLdTPTaP8sjd8ouu8pbPAm8pTegt88B/eeq17wmshAEAsIQQBgDAEkIYAABLCGEAACwhhAEAsIQQBgDAEkIYAABLCGEAACwhhAEAsIQQBgDAEkIYAABLqB0Nr+b3SqjLfW/yucHlvuuGr72Y6ZQr3lA/l5rNuNyU9fuKlTAAAJYQwgAAWEIIAwBgCSEMAIAlhDAAAJYQwgAAWEIIAwBgCSEMAIAlhDAAAJYQwgAAWELZSvzPKB948byhtKQneUvZSnfwGsWFULYSAIDLBCEMAIAlhDAAAJYQwgAAWEIIAwBgCSEMAIAlhDAAAJYQwgAAWEIIAwBgCSEMAIAlFW1PAAA8jVKU8FashAEAsIQQBgDAEkIYAABLCGEAACwhhAEAsIQQBgDAEkIYAABLCGEAACwhhAEAsIQQBgDAEspWolghISEu93WnJKA747o7Nq4snnpt8BrFpcRKGAAASwhhAAAsIYQBALCEEAYAwBJCGAAASwhhAAAsIYQBALCEEAYAwBJCGAAASwhhAAAsIYQBALCE2tHlmCdr3FIP99Lw5OPsqfrf7oh5fq9b/bePrueReXjDYwEUh5UwAACWEMIAAFhCCAMAYAkhDACAJR4N4ccee0yJiYlF2j///HOFhoYqNzfXk5v/n+3Zs0ehoaHaunWr7akAAC5D5XYlnJeXZ3sK1vEYAED55hUhfPjwYT300EOKiYlR7dq1dfPNNyslJcWpT6dOnZScnKzRo0crKipK9913nwYNGqSBAwc69cvPz1dUVJTeffddSVJBQYGmTp2qFi1aqFatWrrlllv0t7/9zdH/119/1cMPP6xGjRqpVq1aatWqlRYuXChJatGihSSpbdu2Cg0NVadOnYrM3RijG264Qa+99ppT+9atWxUaGqpdu3ZJknJzc/XEE0+oUaNGioyMVJcuXZxW2Lt371bv3r3VuHFj1a1bV/Hx8frkk0+cxmzevLleeuklDRkyRJGRkRo8eLA7DzMAwMt4xfeET506pZYtW+qpp55SYGCg1qxZoyFDhqhhw4Zq1aqVo9+7776rAQMGaPXq1ZKkXbt2acCAATp27JiqVasmSVq7dq1Onjypjh07SpKmTp2qpUuXaurUqYqKitL69ev1yCOPKCwsTG3atNFzzz2nb7/9VkuXLlWNGjW0a9cunTx50jFWhw4d9P7776tJkyby8/MrMncfHx8lJiZq0aJFevzxxx3tCxcu1M0336yrr75akvTggw/K399fS5cuVVBQkObPn69u3brpq6++UmhoqI4dO6Y77rhDTz/9tCpXrqx3331XCQkJ+vLLL1Wv3n+/O/naa69pxIgRGjlypIKCgsr4mQAAXEoeD+H09HTVrVvXqe3MmTNOP9euXdspwAYPHqy1a9fq/fffdwrhq6++WhMnTnT83LBhQ1WpUkVpaWnq1auXJCklJUV33323AgMDdfr0aU2bNk3Lli1T69atJUkNGjTQxo0bNW/ePLVp00b79u1T8+bNdf3110uSIiMjHeOHhYVJkqpXr66aNWuWuI+9e/fW888/r6+//lqtWrVSfn6+UlJSNGnSJEnShg0b9PXXXysrK0uVK1eWJE2aNEkrV67U3/72Nz344INq1qyZmjVr5hhz7NixSktL06pVq5xWvG3bttWwYcMkuV+sAwDgXTwewrfeeqteeeUVp7avvvpKjzzyiOPnM2fOaOrUqVq2bJl+/vln5efn6/Tp06pSpYrT/Vq2bOn0c8WKFdWtWzelpKSoV69eOn78uFatWqU333xT0tmV8okTJ9S9e3en++Xl5al58+aSpIEDB6p///7avHmzbrvtNnXs2FFxcXFu7WOtWrV05513asGCBWrVqpVWr16tvLw8de3aVZK0bds2HT9+XFFRUU73O3nypHbv3i1JOnbsmKZMmaI1a9bowIEDOnPmjE6ePKl9+/Zd8DEAAJRfHg/hKlWqOA7JFtq/f7/Tz6+++qr+7//+T5MnT1ZMTIyqVq2q0aNHF7nw6PxQlqQePXqoU6dOOnTokD7++GP5+/vr9ttvlyQdP35ckrR48WLVqlXL6X6Fh5bvuOMObdmyRR9++KE+/vhjdevWTYMGDXKsYl3Vt29fDRkyRJMnT9bChQt17733OuZ7/PhxXXXVVVqxYkWR+wUHB0uSxo0bp08++USTJk1Sw4YNFRAQoP79+ys/P9+pf9WqVR3/vxLK65XHcoPeMmdveDzcLUPpDUd3PFkOFjifV5wT3rRpk+655x717NlT0tmLqb7//ntdc801pd43Li5OderU0bJly/Thhx+qa9euqlSpkiTpmmuuUeXKlbV37161adOmxDHCwsKUkJCghIQEzZs3T88++6wmTZrkGOf8w+fFufPOO1W1alXNnTtXa9eu1cqVKx23tWjRQgcPHlTFihWdDnef/xj07t3bcfHXsWPHtGfPnlK3CwAov7wihKOiovS3v/1NmzZtUkhIiGbNmqXs7GyXQliS7r//fs2bN0/fffedli9f7mgPDAzUsGHDNHbsWBljdOONN+rIkSPatGmTAgMDlZCQoMmTJ6tly5Zq0qSJTp8+rfT0dDVu3FiSFB4eroCAAH300UeqXbu2Kleu7Fi5nq9ChQpKSEjQxIkTFRUV5TgHLUnt27dXbGysEhMTNWHCBDVq1Eg///yz1qxZo06dOun6669XVFSUVqxYobvvvluSNHnyZBljLvYhBQCUA17xFaWkpCS1aNFC999/vzp37qyIiAjH1c2u6NGjh3bs2KFatWrpxhtvdLpt7NixSk5O1rRp0xQXF6f7779fa9ascaxI/fz8NHHiRN1yyy3q2LGjKlSooLfeekvS2XPOL7zwgubPn69rr7222MIj5+rbt6/y8vLUu3dvp3YfHx8tWbJEN998s4YNG6bf/e53euihh7R3716Fh4dLkp577jmFhITorrvuUkJCguLj4x3nrQEAlycf4+Jyi/MepVu/fr26deumbdu2KSIiwvZ0vJInz/l5y2vUW84Jl0fecE7YXTyHV5ayfo16xeHo8u706dPKycnRlClT1LVrVwIYAOASrzgcXd699957at68uXJzczVhwgTb0wEAlBMcjsYlxeFoZ94yZ2/B4Wh4u7J+jbISBgDAEkIYAABLCGEAACzh6mj8z8rjOVBPlib01D4+07d76Z3OMfHtVJf7lsdzse5o8Uyr0jv9x7rhaz04E9fFPL/Xrf7ulgiFd2AlDACAJYQwAACWEMIAAFhi7Zxwu6kdLun23D3PM3XqVKWlpSkrK0v+/v5q3bq1xo8fr+joaEefU6dO6emnn1Zqaqry8vIUHx+vl19+2VExa+vWrZo+fbo2btyow4cPKzIyUgMGDNCQIUOctvX5559r7Nix2rFjh+rUqaOkpKQi9acBAJcfVsIlWL9+vQYNGqQ1a9YoNTVV+fn56t69u+NvFEvSmDFjtHr1as2fP19paWk6cOCA+vbt67h98+bNCg8P1+zZs7VhwwYNHz5cEydO1OzZsx19fvzxR/Xs2VO33nqrPv30Uw0ZMkRPPPGE1q71jotDAACew9XRJUhJSXH6edasWYqOjtY333yjNm3aKDc3VwsWLNCcOXPUtm1bSdKMGTMUFxenjIwMxcbGqk+fPk5jNGjQQBkZGUpLS9PgwYMlSXPnzlVkZKT+9Kc/STr7N5A3btyo119/XR06XNqjBQCAS4uVsIuOHDkiSQoNDZV0dpWbn5+v9u3bO/o0btxYdevWVUZGxgXHKRxDkjIyMpzGkKQOHTroyy+/LLvJAwC8EiHsgoKCAo0ePVpxcXGKiYmRJB08eFB+fn4KDg526hsREaGDBw8WO86mTZu0bNky9e/f39GWnZ3t+JvChcLDw3X06FGdPHmyjPcEAOBNOBztgqSkJGVmZmrVqlUXPcb27duVmJiokSNHKj4+vgxnBwAor1gJlyI5OVnp6elasWKF6tSp42ivWbOm8vLylJub69Q/OztbNWvWdGrbsWOHunXrpv79+yspKcnptoiICB06dMip7dChQwoMDFRAQEAZ7w0AwJsQwiUwxig5OVkrV67U8uXLVb9+fafbW7RooUqVKmndunWOtqysLO3bt0+xsbGOtszMTHXp0kW9evXSuHHjimwnNjbWaQxJ+vjjj9W6desy3iMAgLfhcHQJkpKSlJKSokWLFqlatWqO87xBQUEKCAhQcHCw+vTpo7Fjxyo0NFSBgYEaMWKEYmNjHSG8fft2de3aVfHx8Ro6dKhjjAoVKigsLEySNHDgQL355pt65pln1KdPH3366ad6//33tXjxYjs7/h/lsR60O57o7N4pAXfqMHvKuutzS+90kdx5Dj352vDU2N5SD9od4ZUfdKv/1muvdblvs8wZbs4GnkIIl2Du3LmSpE6dOjm1z5w501FIY/LkyfL19VW/fv2cinUUWr58uXJycrRkyRItWbLE0V6vXj1t2bJFklS/fn0tXrxYY8aM0RtvvKHatWvr1Vdf5etJAHAFsBbC3v6b6S+//FJqH39/f7388stOwXuuUaNGadSoUaWOc8stt+jTTz91e44AgPKNc8IAAFhCCAMAYAkhDACAJYQwAACWEMIAAFhCCAMAYAkhDACAJYQwAACWUDELxSqPpSjLI3fKNHpLgRtPlbiEM7ef7+GemQc8i5UwAACWWFsJxzy/95Jub/voem71nzp1qtLS0pSVlSV/f3+1bt1a48ePV3R0tKPPqVOn9PTTTys1NdWpdnRERIQkaevWrZo+fbo2btyow4cPKzIyUgMGDNCQIUMcYxw4cEBPP/20vvnmG+3atUuPPPKInn/++bLZaQCAV2MlXIL169dr0KBBWrNmjVJTU5Wfn6/u3bvr+PHjjj5jxozR6tWrNX/+fKWlpenAgQPq27ev4/bNmzcrPDxcs2fP1oYNGzR8+HBNnDhRs2fPdvTJy8tTWFiYkpKSdN11113SfQQA2MU54RKkpKQ4/Txr1ixFR0frm2++UZs2bZSbm6sFCxZozpw5atu2rSRpxowZiouLU0ZGhmJjY9WnTx+nMRo0aKCMjAylpaVp8ODBkqTIyEi98MILkqQFCxZcgj0DAHgLVsIuOnLkiCQpNDRU0tlVbn5+vtq3b+/o07hxY9WtW1cZGRkXHKdwDADAlY2VsAsKCgo0evRoxcXFKSYmRpJ08OBB+fn5KTg42KlvRESEDh48WOw4mzZt0rJly7R48WKPzxkA4P0IYRckJSUpMzNTq1atuugxtm/frsTERI0cOVLx8fFlODsAQHnF4ehSJCcnKz09XStWrFCdOnUc7TVr1lReXp5yc3Od+mdnZ6tmzZpObTt27FC3bt3Uv39/JSUlXZJ5AwC8HyFcAmOMkpOTtXLlSi1fvlz169d3ur1FixaqVKmS1q1b52jLysrSvn37FBsb62jLzMxUly5d1KtXL40bN+6SzR8A4P04HF2CpKQkpaSkaNGiRapWrZrjPG9QUJACAgIUHBysPn36aOzYsQoNDVVgYKBGjBih2NhYRwhv375dXbt2VXx8vIYOHeoYo0KFCgoLC3Nsa+vWrZKk48ePKycnR1u3blWlSpXUpEmTS7zXAIBLyccYY1zpWNZlDL29WEdJVzDPnDlTvXv3lvTfYh3vvfeeU7GOwsPRL7zwgqZMmVJkjHr16mnLli0X3Nb5fa5UW68d5nLfZpkzvGIePfuNdLnv/inNLmY65Ya7nxvulLmktCpsKOtSrNZCGHAFIVy+EcK43JR1CHNOGAAASwhhAAAsIYQBALCEEAYAwBJCGAAASwhhAAAsIYQBALCEEAYAwBJCGAAAS6gdDa827OFMl/uuK72Lg98rxZclLUmzzF9c7rvfQ1Wf3J1zlUkuFcPzqLKuLgRcbqyFsDtlAMuCuyUNp06dqrS0NGVlZcnf31+tW7fW+PHjFR0d7ehTWDs6NTXVqXZ0RESEpLN/mGH69OnauHGjDh8+rMjISA0YMEBDhgxxjLFixQrNnTtXW7duVV5enpo0aaKRI0eqQ4cOZbPjAACvxeHoEqxfv16DBg3SmjVrlJqaqvz8fHXv3l3Hjx939BkzZoxWr16t+fPnKy0tTQcOHFDfvn0dt2/evFnh4eGaPXu2NmzYoOHDh2vixImaPXu203bat2+vJUuW6OOPP9Ytt9yihIQE/ngDAFwBOBxdgpSUFKefZ82apejoaH3zzTdq06aNcnNztWDBAs2ZM0dt27aVJM2YMUNxcXHKyMhQbGys+vTp4zRGgwYNlJGRobS0NA0ePFiS9Pzzzzv1eeaZZ7Rq1SqtXr1azZs39+AeAgBsYyXsoiNHjkj6758d3Lx5s/Lz89W+fXtHn8aNG6tu3brKyMi44Dgl/ZlESSooKNDRo0c5lwYAVwBWwi4oKCjQ6NGjFRcXp5iYGEnSwYMH5efnp+DgYKe+EREROnjwYLHjbNq0ScuWLdPixYtL3NZrr72m48eP69577y27HQAAeCVC2AVJSUnKzMzUqlWrLnqM7du3KzExUSNHjlR8fHyxfZYuXaoXX3xRCxcuVHh4+EVvCwBQPnA4uhTJyclKT0/XihUrVKdOHUd7zZo1lZeXp9zcXKf+2dnZqlmzplPbjh071K1bN/Xv319JSUnFbue9997Tk08+qblz5zod4gYAXL4I4RIYY5ScnKyVK1dq+fLlql+/vtPtLVq0UKVKlbRu3X+/nZqVlaV9+/YpNjbW0ZaZmakuXbqoV69eGjduXLHbSklJ0bBhw/Tmm2/qrrvu8swOAQC8DoejS5CUlKSUlBQtWrRI1apVc5znDQoKUkBAgIKDg9WnTx+NHTtWoaGhCgwM1IgRIxQbG+sI4e3bt6tr166Kj4/X0KFDHWNUqFBBYWFhks4egn7sscf0/PPPq1WrVo4+/v7+Rc43AwAuL4RwCebOnStJ6tSpk1P7zJkz1bt3b0nS5MmT5evrq379+jkV6yi0fPly5eTkaMmSJVqyZImjvV69eo7vAf/lL3/Rb7/9puTkZCUnJzv6JCQkaNasWR7bPwCAfT7GGJdq27lTXg8oK+2mul45bN3wtS73dbcEZN4fXS9b6c7Xyy73spWexGcSbCjrr496ZCXsLd9x5U16abhTgrRnv5Fujf2d3z9c7pvnxrjuBlSeG6+lz2r1Kb3Tf7hTTvVyD1WJ9yyuPFyYBQCAJYQwAACWEMIAAFhCCAMAYAkhDACAJYQwAACWEMIAAFhCCAMAYAkhDACAJdZqR58Y5+Pxbfid8393yg5K0tSpU5WWlqasrCz5+/urdevWGj9+vKKjox19Tp06paefflqpqalOtaMjIiIkSVu3btX06dO1ceNGHT58WJGRkRowYICGDBniGGPDhg0aP368srKydPLkSdWrV08PPvigHnvssf9p3wEA3s/l2tFl7VKE8LncDeH7779f3bt31/XXX6/ffvtNkyZNUmZmpjZu3KiqVatKkoYPH641a9Zo1qxZCgoK0ogRI+Tj46P09HRJ0oIFC7Rt2zZ17txZderU0aZNm/SHP/xB48eP1+DBgyVJW7Zs0c6dO9W0aVNVrVpVS9s+qpePbNKwwFbqUiW6xPlJ7pU89BbltaRpzPN7Xe67f0ozN2dz+boSylB6qr45vFO5qB19OUhJSXH6edasWYqOjtY333yjNm3aKDc3VwsWLNCcOXPUtm1bSdKMGTMUFxenjIwMxcbGqk8f5xrCDRo0UEZGhtLS0hwh3Lx5czVv3tzR586Ahlp3aq+25GWXGsIAgPKNc8IuOnLkiCQpNPTsX7LZvHmz8vPz1b59e0efxo0bq27dusrIyLjgOIVjFGdn/mH9K/+QWvrVLJuJAwC8FithFxQUFGj06NGKi4tTTEyMJOngwYPy8/NTcHCwU9+IiAgdPHiw2HE2bdqkZcuWafHixUVua9q0qXJycvRbXr4GVGumTlUalf2OAAC8CiHsgqSkJGVmZmrVqlUXPcb27duVmJiokSNHKj4+vsjtH3zwgY4dO6a/3fNHvXH0n6pTIVC3BzT4H2YNAPB2hHApkpOTlZ6erg8++EB16tRxtNesWVN5eXnKzc11Wg1nZ2erZk3nQ8k7duxQt27d1L9/fyUlJRW7nfr160uSCqo00uGCk5p3bAshDACXOc4Jl8AYo+TkZK1cuVLLly93hGShFi1aqFKlSlq3bp2jLSsrS/v27VNsbKyjLTMzU126dFGvXr00btw417YtKd8UlMl+AAC8FyvhEiQlJSklJUWLFi1StWrVHOd5g4KCFBAQoODgYPXp00djx45VaGioAgMDNWLECMXGxjpCePv27eratavi4+M1dOhQxxgVKlRQWFiYJGnOnDmqW7euGjduLElKO/Gd3j2+XfdVucbCXgMALiVCuARz586VJHXq1MmpfebMmerdu7ckafLkyfL19VW/fv2cinUUWr58uXJycrRkyRItWbLE0V6vXj1t2bJF0tkV98SJE7Vnzx5VqFBBV52sqCGB16tLAF9PAoDLnbViHZdCeSwUsPXaYS73pVjHxaNYx6VRHt+D7qJYx5WlrD/DOCcMAIAlLq+EPfkbrTesjrzlN/Zn+nZ3ue/Et1M9OBPXeXL17s5rw53n0JOvOW+Zh6d4y3sFsIGVMAAAlwlCGAAASwhhAAAsIYQBALCEEAYAwBJCGAAASwhhAAAsIYQBALDEWu1od4pSlIVXV/zdrf5Tp05VWlqasrKy5O/vr9atW2v8+PGKjv5vTedTp07p6aefVmpqqlPt6IiICEnS1q1bNX36dG3cuFGHDx9WZGSkBgwYoCFDhhS7zY0bN2rGyk9UI7CqEtrGFtsHAHD5YCVcgvXr12vQoEFas2aNUlNTlZ+fr+7du+v48eOOPmPGjNHq1as1f/58paWl6cCBA+rbt6/j9s2bNys8PFyzZ8/Whg0bNHz4cE2cOFGzZ88usr3c3Fw9+uijqhcWcil2DwDgBayVrfS2lXBp+5eTk6Po6GilpaWpTZs2ys3NVXR0tObMmaOuXbtKknbu3Km4uDitWbPG6W8KnyspKUk7d+7U8uXLndoHDhyoqKgoVahQQStXrtRnn31W6j65Uzhecq94fHksp+gtKOt4ZSmPpWZx8ShbacmRI0ckSaGhoZLOrnLz8/PVvn17R5/GjRurbt26ysjIuOA4hWMUWrhwoX788UeNHDmy7CcOAPBa/D1hFxQUFGj06NGKi4tTTEyMJOngwYPy8/NTcHCwU9+IiAgdPHiw2HE2bdqkZcuWafHixY6277//XhMmTNAHH3ygihV5OgDgSsKnvguSkpKUmZmpVatWXfQY27dvV2JiokaOHKn4+HhJ0pkzZ/Twww9r1KhRatSoUVlNFwBQThDCpUhOTlZ6ero++OAD1alTx9Fes2ZN5eXlKTc312k1nJ2drZo1azqNsWPHDnXr1k39+/dXUlKSo/3YsWP65z//qS1btmjEiBGSzq66jTEKCwtTamqq2rZt6+E9BADYQgiXwBijESNGaOXKlVqxYoXq16/vdHuLFi1UqVIlrVu3Tl26dJEkZWVlad++fU4XZWVmZqpr167q1auXxo0b5zRGYGCgvvjiC6e2t956S5999pnmz59fZJsAgMsLIVyCpKQkpaSkaNGiRapWrZrjPG9QUJACAgIUHBysPn36aOzYsQoNDVVgYKBGjBih2NhYRwhv375dXbt2VXx8vIYOHeoYo0KFCgoLC5Ovr6/jHHOh8PBwVa5cuUg7AODyQwiXYO7cuZKkTp06ObXPnDlTvXv3liRNnjxZvr6+6tevn1OxjkLLly9XTk6OlixZoiVLljja69Wrpy1btlyCvQAAeDNr3xO+WJ76/qq37J87+J6wdyqPryVcPL4nfGXhe8IAAFwmCGEAACwhhAEAsIRzwv/h7v65M48T43xc7pv3x1/cmoencE7Y+3jyNerJeQCXE84JAwBwmSCEAQCwhBAGAMASQhgAAEsIYQAALCGEAQCwhBAGAMASQhgAAEsIYQAALCGEAQCwxOWylQAAoGyxEgYAwBJCGAAASwhhAAAsIYQBALCEEAYAwBJCGAAASwhhAAAsIYQBALCEEAYAwJL/Bwn4SVHVF9tMAAAAAElFTkSuQmCC",
      "text/plain": [
       "<Figure size 600x600 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "years_to_display = [BASE_YEAR + i for i in range(5)]\n",
    "files_by_year = {int(path.stem.split('_')[1]): path for path in output_files}\n",
    "available = [(year, files_by_year[year]) for year in years_to_display if year in files_by_year]\n",
    "\n",
    "if not available:\n",
    "    raise ValueError('No harvest rasters found for the requested years.')\n",
    "\n",
    "from matplotlib.colors import ListedColormap, BoundaryNorm\n",
    "from matplotlib.patches import Patch\n",
    "\n",
    "with rasterio.open(INVENTORY_TIF) as src:\n",
    "    inventory = src.read(1)\n",
    "    inv_nodata = src.nodata\n",
    "    template_shape = inventory.shape\n",
    "    if inv_nodata is None:\n",
    "        forest_mask = np.ma.array(np.ones(template_shape, dtype=int), mask=np.isnan(inventory))\n",
    "    else:\n",
    "        forest_mask = np.ma.array(np.ones(template_shape, dtype=int), mask=(inventory == inv_nodata))\n",
    "\n",
    "display = np.full(template_shape, fill_value=-1, dtype=int)\n",
    "\n",
    "for idx, (year, tif_path) in enumerate(available):\n",
    "    with rasterio.open(tif_path) as src:\n",
    "        arr = src.read(1)\n",
    "        nodata = src.nodata\n",
    "        if nodata is None:\n",
    "            harvested = ~np.isnan(arr)\n",
    "        else:\n",
    "            harvested = arr != nodata\n",
    "        display = np.where(harvested, idx, display)\n",
    "\n",
    "masked_display = np.ma.array(display, mask=display == -1)\n",
    "\n",
    "color_palette = ['#2E7D32', '#1976D2', '#C2185B', '#F57C00', '#6D4C41'][:len(available)]\n",
    "harvest_cmap = ListedColormap(color_palette)\n",
    "norm = BoundaryNorm(range(len(color_palette) + 1), harvest_cmap.N)\n",
    "background_cmap = ListedColormap(['#e0e0e0'])\n",
    "\n",
    "fig, ax = plt.subplots(figsize=(6, 6))\n",
    "ax.imshow(forest_mask, cmap=background_cmap, interpolation='nearest', alpha=0.5)\n",
    "ax.imshow(masked_display, cmap=harvest_cmap, norm=norm, interpolation='nearest')\n",
    "ax.set_title(f\"Harvest disturbance by year ({available[0][0]}-{available[-1][0]})\")\n",
    "ax.axis('off')\n",
    "\n",
    "handles = [Patch(facecolor=color_palette[i], label=str(year)) for i, (year, _) in enumerate(available)]\n",
    "ax.legend(handles=handles, title='Harvest year', loc='lower left', frameon=False)\n",
    "plt.show()\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "id": "594636bc",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>year</th>\n",
       "      <th>disturbed_pixels</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>2020</td>\n",
       "      <td>10</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>2021</td>\n",
       "      <td>10</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>2022</td>\n",
       "      <td>10</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>2023</td>\n",
       "      <td>9</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <td>2024</td>\n",
       "      <td>9</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "   year  disturbed_pixels\n",
       "0  2020                10\n",
       "1  2021                10\n",
       "2  2022                10\n",
       "3  2023                 9\n",
       "4  2024                 9"
      ]
     },
     "execution_count": 13,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "summary = []\n",
    "for tif in output_files:\n",
    "    year = int(tif.stem.split('_')[1])\n",
    "    with rasterio.open(tif) as src:\n",
    "        arr = src.read(1)\n",
    "        nodata = src.nodata\n",
    "        disturbed = np.count_nonzero(arr != nodata)\n",
    "    summary.append({'year': year, 'disturbed_pixels': disturbed})\n",
    "\n",
    "pd.DataFrame(summary).head()\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "acc99b9c",
   "metadata": {},
   "source": [
    "### Next steps\n",
    "\n",
    "- Supply your own management-unit mask by passing a tuple (e.g., `('tsa24_clipped', '1', '?', '?', '?')`) to `ForestRaster.allocate_schedule`.\n",
    "- Use finer time steps (e.g., quarterly) by setting `time_step=0.25` when building `ForestRaster`.\n",
    "- Combine multiple actions by expanding `acode_map` (each entry gets its own set of year-indexed rasters).\n",
    "\n",
    "Remember to tidy up raster outputs if the directory should remain clean for version control."
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": ".venv",
   "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.12.3"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}