10 Sensitivity analysis of Wofost.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<img style=\"float: right;\" src=\"data:image/png;base64,iVBORw0KGgoAAAANSUhEUgAAAOIAAAAjCAYAAACJpNbGAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAABR0RVh0Q3JlYXRpb24gVGltZQAzLzcvMTNND4u/AAAAHHRFWHRTb2Z0d2FyZQBBZG9iZSBGaXJld29ya3MgQ1M26LyyjAAACMFJREFUeJztnD1y20gWgD+6nJtzAsPhRqKL3AwqwQdYDpXDZfoEppNNTaWbmD7BUEXmI3EPMFCR2YI1UDQpdAPqBNzgvRZA/BGUZEnk9FeFIgj0z2ugX7/XP+jGer2mLv/8b6d+4Efgf/8KG0+Zn8XyXLx+bgEslqegcfzxSY3Irrx6bgEsFssBWsRGowGufwHAYtq7u+H6fUCOxTTWax4wBAbr+SRqNDKesOv3gN/133sW0yh927j1mucIaFWINl7PJ+OcvMcfW8Bol3iN44+mLIOsTCp3UJFfAETr+WRQcG8EOJpunEnTyDlYzycbeWr5xxq3jOF6PglK8ix9buv5xCsrAzBkMV1l5OwD/aJ4BXzV3+8F9z4gz/hTSbz8cxc84FuNvDc4VIsYA7+qohmGwAnycA194G22YqUYlZxv4vpN4AuwBv4oON5m8k3TVLnK4sYFcRyN86dWvCwnlCvFCeUVvwX8CkSZZ5eWs5mLJWE/VZThBMgpfirPk5J4f1SU4QsQ6LNP4+j9OkSUKdRiGlD87CWe3PcyR5PFdAhc1cz/joOziMoIeVF95GX1EGVY6bWhvsAeZQrm+kON80PDneD6PRbTi4LQpmJfsZieFaR1qXlXURh3y2BaBPyG63sspv0t6e+CKJTrf2YxHe8Qr6z8AXBdGbMoHgCTshgr4AiItfxljenPJGv5roCi+rGVw1TExTTWl99ThRsglfYHUnF7SMv+Bhjn4idxbhFLGiAu6gjXD3LuUBF5VzWi3CoAfMP1kxe7mNYZMT5DLFgf13eAXi3ZtvMOsUb3V3J5/mmqy+/66RbnTC1LFdfIu/kd8Qx2bTQeg2GBTPfiUF1TgHNE0QaIq/JDX9RKr/WBy/V8EhfEHWncWMO2EKV8S7UypYnYdE2r+o8gyj5MHXVYsZh+JnG7A+3LPQxR5g9II/UJ148ockmrybqm2+Qapo6gppwB8J7EM6jqaz8u0lhfkXgB58BKPam6rvEdh2kRARbTMa7/HXEfVqnW8hxxWwE+5+JJRTYd9CM90gxw/XFuMKMo/yTNDzUkLnbr6rCYnuH6N8igQ3CvNPJproDPuH6MKMd4Z5kMUjnrh98tn1if72/Ie729Vzq708L0YV3/HGmgB4iHsjOProhhd1lrEr4zaz/FvM4lolTnqWum/6jKmeuDmFb1jHylNg96hPQbhcU0wPVBXESvQI4W5aNshsK4jeOPhSOcOaThMVb48dhU8m2UlR+29ZHzrqyhLL0EaTROteGt67EYIsT6F1HXC/ikcvS00dl51PRwLaIwQtzCxGWRFnRMkT8v/SyAy8I+iliHJtDUsHHq7imipE42GtJanxdcB6mgQcm9MmKNs1m5F9MI13+n+cXZSEpAeV8mQgZqNkmU/HsuT7kf4PrGhXcK0h1SXv7iPKsJKCrDYvoV17+meMqhiDFlll7GEb4U3iseAf+k7mqksmU9qUoaj73E7TEtol3iZnks7Moai8WylUN3TS0WANbzyYv2rqxFtFheANYi7iGNRoPOrO2QGTQIu8vhU8vSmbWNDAHQD7vLYWfWbgFx2F3ee3FBZ9ZuIgMpTWAQdpeRXm9pPoPOrD3UMCtkQM4BRmF3ubG6ZZdxkOfCWsT9pU96CuX56KfOjeIFVC8Ar8NI0xuyOQJsVkWl8xzptQGPNY/6xFiLuL+0gIu0FVTrNESmbK7C7tLrzNpmPW0EeGF32UyFN19UnCAT4ZHGWWnYqDNrB4jViZBK/kbD9sLuMiBZSD8AVp1Z+0LD/NmZta+BIzOS3pm1xwBhd9kvkeEGUbQeqSmIdHhkXnGs5fIQRUxPV1x0Zm2zMuoq7C69rU/yBWAt4v7iAd86s/ZaDweZP+wBvwBOZ9b2SCrrmPzk+AWizA09j1QxMK4gZumcWKUWMvkdA56mfxN2l7GmHWk6V2F32Qi7yxaIsmnYHvkJ9zEQqAwBotQXwK2m0c+EN/Kk8zPTZiOkIWrp/xNTnpeOtYh7iFauN+k5W+0vXab6UsbyecAw229SxWiG3aVZ7NBCKrGHuneazy2iyBeIuxkjk9UDE1bzOtJ4IzbdwysNN0D6dnf9Rk3/iKSBWOnhUbASSWW+DbvLWM+HKreZ3O/r77gza5u842w6LxFrEfcTj+Jv3mK4q7Co63hE+fI6E94hUaT0cry+XushSuvoNZO2CdsCrlXJHDYVMUIUJso2BmhfL+wuV6rMvVR6AXnS1428XupaE7Hwnrqkg4cMGD0lr3NfpVegrUw1m2sN0+crNirEX1uTqiPbPoyI/QSKKmqA9I9aer+fcR2zxIj7GiMV+EYVIkZc3r5eH2rYI+0vnpBYIE/vGwUCdYM7s3agbqXJu58VIOwug86sfd2ZtSPNKwi7S9PHy4UnscCmXKuUZQRdsqbPwCHp2754pKYnW0akcZBO/x2df29XnvA//6iV8T3TSluBmOQlR+v5JNvaHixlDZRalRZifbZaAg3vIIrkmP6YVu6owI1M9x2r0vVIFCBGXNLS96Ph45IGY2ey6e1DY20UMaLGItUXoIhVvCv5tvDg2MWLqYNaoKBKWe6Z7gBR8OwAzZOyD4poBmtidlwt/gIxw/QHz0+oWKIoj19fRz8p3YOjoV8195F5l31ltZ5PfnluISyW+/IK6SPstRIiH/FaLHvLa2R+6F6f978AVsD7v0vf0HK4vNK9VfbVojSBceP4o/PcglgsD8GMmjaRbRCc1PEQIrbv45nlIfleIrs778XkrcWSZXMcXPZyqbvfxy7ckuyqHJPslJzH9c3We2ZRbx1O/07ziJbDI1FE2Qwp4n4DNzHJhkZF16+3bnwrCmi40U2eWoj7KZvobn7+YtKO1vPJVyyWPSZrER1kNU0TqfienpvlaWZR7oX+3tba6lxcX7MK3tNfo2RlpNc8tthsIFbAKYtpsA+TtRbLNp5/H4/EFXX0MOfbOGUxvbCKaDkEnl8Rq0jc1ayFjhFFjKwiWg6B/wNk+JCXXNBIXQAAAABJRU5ErkJggg==\">\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<h1>Applying sensitivity analysis to WOFOST</h1>\n",
    "\n",
    "Sensitivity analysis is the study of how the uncertainty in the output of a mathematical model or system (numerical or otherwise) can be apportioned to different sources of uncertainty in its inputs (<a href=\"https://en.wikipedia.org/wiki/Sensitivity_analysis\">Wikipedia</a>). The sensitivity of each input is often represented by a numeric value, called the sensitivity index.\n",
    "\n",
    "In this notebook we demonstrate how to run a sensitivity analysis for the WOFOST model for a number of important parameters. The problem with WOFOST is that the sensitivity of the model can differ for different crops and regions depending on whether you are close to certain thresholds in the model. Therefore it is more useful to provide a general guideline for sensitivity analysis through this notebook rather than predefined set of sensitive parameters.\n",
    "\n",
    "A limitation of this notebook is that is only applies for the scalar parameters in WOFOST. WOFOST also has a number of tabular parameters which are defined as as a function of development stage or temperature. The sensitivity of these parameters cannot be properly analyzed with the approach used in this notebook.\n",
    "\n",
    "Below the WOFOST model is run for a sugar beet crop in the Netherlands and then sensitivity analysis is applied. The target variable for sensitivity is the total weight of the storage organs (TWSO) which represents the crop yield. However, other variables could be choosen as well such as total crop biomass or maximum leaf area index.\n",
    "\n",
    "*Steven Hoek and Allard de Wit, Wageningen Environmental Sciences, June 2020*\n",
    "\n",
    "## Dependencies for running this notebook\n",
    "This notebook require the following package to be installed:\n",
    "- [PCSE](https://pcse.readthedocs.io/en/stable/) and its dependencies (numpy, pandas, PyYAML, SQLAlchemy,...)\n",
    "- [SAlib](https://salib.readthedocs.io/en/latest/)\n",
    "- [IPython](http://ipython.org/)\n",
    "- [tqdm](https://tqdm.github.io/) for showing a progressbar\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "This notebook was built with:\n",
      "python version: 3.10.13 | packaged by conda-forge | (main, Oct 26 2023, 18:01:37) [MSC v.1935 64 bit (AMD64)]\n",
      "PCSE version: 6.0.6\n"
     ]
    }
   ],
   "source": [
    "import sys, os.path\n",
    "import yaml\n",
    "import numpy as np\n",
    "import pandas as pd\n",
    "from tqdm import tqdm\n",
    "from SALib.sample import saltelli\n",
    "from SALib.analyze import sobol\n",
    "import matplotlib.pyplot as plt\n",
    "from IPython.display import display\n",
    "pd.set_option(\"display.max_rows\", None)\n",
    "pd.set_option(\"display.max_colwidth\", 250)\n",
    "\n",
    "import pcse\n",
    "from pcse.models import Wofost72_PP\n",
    "from pcse.base import ParameterProvider\n",
    "from pcse.input import YAMLCropDataProvider, NASAPowerWeatherDataProvider\n",
    "from pcse.util import WOFOST72SiteDataProvider, DummySoilDataProvider\n",
    "from progressbar import printProgressBar\n",
    "\n",
    "print(\"This notebook was built with:\")\n",
    "print(f\"python version: {sys.version}\")\n",
    "print(f\"PCSE version: {pcse.__version__}\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## WOFOST parameters\n",
    "The table below gives an overview of important scalar WOFOST parameters related to the potential production situation."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "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>Name</th>\n",
       "      <th>Typical Value</th>\n",
       "      <th>Category</th>\n",
       "      <th>Description</th>\n",
       "      <th>Unit</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>CFET</td>\n",
       "      <td>1.000</td>\n",
       "      <td>Water use</td>\n",
       "      <td>correction factor for evapotranspiration</td>\n",
       "      <td>-</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>CVL</td>\n",
       "      <td>0.720</td>\n",
       "      <td>Conversion of assimilates into biomass</td>\n",
       "      <td>efficiency conversion of assimilates into leaf dry matter</td>\n",
       "      <td>kg / kg</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>CVO</td>\n",
       "      <td>0.820</td>\n",
       "      <td>Conversion of assimilates into biomass</td>\n",
       "      <td>efficiency conversion of assimilates into storage organ dry matter</td>\n",
       "      <td>kg / kg</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>CVR</td>\n",
       "      <td>0.720</td>\n",
       "      <td>Conversion of assimilates into biomass</td>\n",
       "      <td>efficiency conversion of assimilates into root dry matter</td>\n",
       "      <td>kg / kg</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <td>CVS</td>\n",
       "      <td>0.690</td>\n",
       "      <td>Conversion of assimilates into biomass</td>\n",
       "      <td>efficiency conversion of assimilates into stem dry matter</td>\n",
       "      <td>kg / kg</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>5</th>\n",
       "      <td>DLC</td>\n",
       "      <td>16.000</td>\n",
       "      <td>Phenology</td>\n",
       "      <td>critical day length for development (lower threshold)</td>\n",
       "      <td>hours</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>6</th>\n",
       "      <td>DLO</td>\n",
       "      <td>8.000</td>\n",
       "      <td>Phenology</td>\n",
       "      <td>optimum day length for development</td>\n",
       "      <td>hours</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>7</th>\n",
       "      <td>PERDL</td>\n",
       "      <td>0.030</td>\n",
       "      <td>Death rates</td>\n",
       "      <td>max. relative death rate of leaves due to water stress</td>\n",
       "      <td>-</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>8</th>\n",
       "      <td>Q10</td>\n",
       "      <td>2.000</td>\n",
       "      <td>Maintenance respiration</td>\n",
       "      <td>increase of the respiration rate per 10 °C temperature increase</td>\n",
       "      <td>-</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>9</th>\n",
       "      <td>RDI</td>\n",
       "      <td>10.000</td>\n",
       "      <td>Water use</td>\n",
       "      <td>initial rooting depth</td>\n",
       "      <td>cm</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>10</th>\n",
       "      <td>RDMCR</td>\n",
       "      <td>120.000</td>\n",
       "      <td>Water use</td>\n",
       "      <td>crop-dependent maximum rooting depth</td>\n",
       "      <td>cm</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>11</th>\n",
       "      <td>RGRLAI</td>\n",
       "      <td>0.016</td>\n",
       "      <td>Initial state</td>\n",
       "      <td>max. relative increase in leaf area index</td>\n",
       "      <td>day-1</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>12</th>\n",
       "      <td>RML</td>\n",
       "      <td>0.030</td>\n",
       "      <td>Maintenance respiration</td>\n",
       "      <td>maintenance respiration rate coefficient of leaves per day</td>\n",
       "      <td>day-1</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>13</th>\n",
       "      <td>RMO</td>\n",
       "      <td>0.003</td>\n",
       "      <td>Maintenance respiration</td>\n",
       "      <td>maintenance respiration rate coefficient of storage per day</td>\n",
       "      <td>day-1</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>14</th>\n",
       "      <td>RMR</td>\n",
       "      <td>0.015</td>\n",
       "      <td>Maintenance respiration</td>\n",
       "      <td>maintenance respiration rate coefficient of roots per day</td>\n",
       "      <td>day-1</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>15</th>\n",
       "      <td>RMS</td>\n",
       "      <td>0.015</td>\n",
       "      <td>Maintenance respiration</td>\n",
       "      <td>maintenance respiration rate coefficient of stems per day</td>\n",
       "      <td>day-1</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>16</th>\n",
       "      <td>SPA</td>\n",
       "      <td>0.000</td>\n",
       "      <td>Crop green area</td>\n",
       "      <td>specific pod area</td>\n",
       "      <td>ha / kg</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>17</th>\n",
       "      <td>SPAN</td>\n",
       "      <td>37.000</td>\n",
       "      <td>Crop green area</td>\n",
       "      <td>life span of leaves growing at 35 Celsius</td>\n",
       "      <td>day</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>18</th>\n",
       "      <td>TBASE</td>\n",
       "      <td>3.000</td>\n",
       "      <td>Crop green area</td>\n",
       "      <td>lower threshold temperature for physiological ageing of leaves</td>\n",
       "      <td>˚C</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>19</th>\n",
       "      <td>TBASEM</td>\n",
       "      <td>3.000</td>\n",
       "      <td>Crop green area</td>\n",
       "      <td>lower threshold temperature for emergence</td>\n",
       "      <td>˚C</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>20</th>\n",
       "      <td>TDWI</td>\n",
       "      <td>0.510</td>\n",
       "      <td>Initial state</td>\n",
       "      <td>initial total crop dry weight</td>\n",
       "      <td>kg / ha</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>21</th>\n",
       "      <td>TEFFMX</td>\n",
       "      <td>20.000</td>\n",
       "      <td>Emergence</td>\n",
       "      <td>max. effective temperature for emergence</td>\n",
       "      <td>˚C</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>22</th>\n",
       "      <td>TSUM1</td>\n",
       "      <td>800.000</td>\n",
       "      <td>Phenological development</td>\n",
       "      <td>temperature sum from emergence to anthesis</td>\n",
       "      <td>˚C day</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>23</th>\n",
       "      <td>TSUM2</td>\n",
       "      <td>1400.000</td>\n",
       "      <td>Phenological development</td>\n",
       "      <td>temperature sum from anthesis to maturity</td>\n",
       "      <td>˚C day</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>24</th>\n",
       "      <td>TSUMEM</td>\n",
       "      <td>90.000</td>\n",
       "      <td>Emergence</td>\n",
       "      <td>temperature sum from sowing to emergence</td>\n",
       "      <td>˚C day</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>25</th>\n",
       "      <td>VERNBASE</td>\n",
       "      <td>14.000</td>\n",
       "      <td>Phenological development</td>\n",
       "      <td>base vernalization requirement</td>\n",
       "      <td>day</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>26</th>\n",
       "      <td>VERNDVS</td>\n",
       "      <td>0.300</td>\n",
       "      <td>Phenological development</td>\n",
       "      <td>critical DVS for vernalization to switch off</td>\n",
       "      <td>-</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>27</th>\n",
       "      <td>VERNSAT</td>\n",
       "      <td>70.000</td>\n",
       "      <td>Phenological development</td>\n",
       "      <td>saturated vernalization requirement</td>\n",
       "      <td>day</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "        Name  Typical Value                                Category  \\\n",
       "0       CFET          1.000                               Water use   \n",
       "1        CVL          0.720  Conversion of assimilates into biomass   \n",
       "2        CVO          0.820  Conversion of assimilates into biomass   \n",
       "3        CVR          0.720  Conversion of assimilates into biomass   \n",
       "4        CVS          0.690  Conversion of assimilates into biomass   \n",
       "5        DLC         16.000                               Phenology   \n",
       "6        DLO          8.000                               Phenology   \n",
       "7      PERDL          0.030                             Death rates   \n",
       "8        Q10          2.000                 Maintenance respiration   \n",
       "9        RDI         10.000                               Water use   \n",
       "10     RDMCR        120.000                               Water use   \n",
       "11    RGRLAI          0.016                           Initial state   \n",
       "12       RML          0.030                 Maintenance respiration   \n",
       "13       RMO          0.003                 Maintenance respiration   \n",
       "14       RMR          0.015                 Maintenance respiration   \n",
       "15       RMS          0.015                 Maintenance respiration   \n",
       "16       SPA          0.000                         Crop green area   \n",
       "17      SPAN         37.000                         Crop green area   \n",
       "18     TBASE          3.000                         Crop green area   \n",
       "19    TBASEM          3.000                         Crop green area   \n",
       "20      TDWI          0.510                           Initial state   \n",
       "21    TEFFMX         20.000                               Emergence   \n",
       "22     TSUM1        800.000                Phenological development   \n",
       "23     TSUM2       1400.000                Phenological development   \n",
       "24    TSUMEM         90.000                               Emergence   \n",
       "25  VERNBASE         14.000                Phenological development   \n",
       "26   VERNDVS          0.300                Phenological development   \n",
       "27   VERNSAT         70.000                Phenological development   \n",
       "\n",
       "                                                           Description  \\\n",
       "0                             correction factor for evapotranspiration   \n",
       "1            efficiency conversion of assimilates into leaf dry matter   \n",
       "2   efficiency conversion of assimilates into storage organ dry matter   \n",
       "3            efficiency conversion of assimilates into root dry matter   \n",
       "4            efficiency conversion of assimilates into stem dry matter   \n",
       "5                critical day length for development (lower threshold)   \n",
       "6                                   optimum day length for development   \n",
       "7               max. relative death rate of leaves due to water stress   \n",
       "8      increase of the respiration rate per 10 °C temperature increase   \n",
       "9                                                initial rooting depth   \n",
       "10                                crop-dependent maximum rooting depth   \n",
       "11                           max. relative increase in leaf area index   \n",
       "12          maintenance respiration rate coefficient of leaves per day   \n",
       "13         maintenance respiration rate coefficient of storage per day   \n",
       "14           maintenance respiration rate coefficient of roots per day   \n",
       "15           maintenance respiration rate coefficient of stems per day   \n",
       "16                                                   specific pod area   \n",
       "17                           life span of leaves growing at 35 Celsius   \n",
       "18      lower threshold temperature for physiological ageing of leaves   \n",
       "19                           lower threshold temperature for emergence   \n",
       "20                                       initial total crop dry weight   \n",
       "21                            max. effective temperature for emergence   \n",
       "22                          temperature sum from emergence to anthesis   \n",
       "23                           temperature sum from anthesis to maturity   \n",
       "24                            temperature sum from sowing to emergence   \n",
       "25                                      base vernalization requirement   \n",
       "26                        critical DVS for vernalization to switch off   \n",
       "27                                 saturated vernalization requirement   \n",
       "\n",
       "       Unit  \n",
       "0         -  \n",
       "1   kg / kg  \n",
       "2   kg / kg  \n",
       "3   kg / kg  \n",
       "4   kg / kg  \n",
       "5     hours  \n",
       "6     hours  \n",
       "7         -  \n",
       "8         -  \n",
       "9        cm  \n",
       "10       cm  \n",
       "11    day-1  \n",
       "12    day-1  \n",
       "13    day-1  \n",
       "14    day-1  \n",
       "15    day-1  \n",
       "16  ha / kg  \n",
       "17      day  \n",
       "18       ˚C  \n",
       "19       ˚C  \n",
       "20  kg / ha  \n",
       "21       ˚C  \n",
       "22   ˚C day  \n",
       "23   ˚C day  \n",
       "24   ˚C day  \n",
       "25      day  \n",
       "26        -  \n",
       "27      day  "
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "df = pd.read_excel(\"ScalarParametersOfWofost-Potential.xlsx\")\n",
    "display(df)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Settings for crop and location\n",
    "In the following, standard data are used to run the WOFOST model for sugarbeet in the Netherlands. Feel free to carry out this sensitivity analysis for a different crop and / or a different location."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Define location, crop and season\n",
    "latitude = 52.0\n",
    "longitude = 5.0\n",
    "crop_name = 'sugarbeet'\n",
    "variety_name = 'Sugarbeet_601'\n",
    "campaign_start_date = '2006-01-01'\n",
    "emergence_date = \"2006-03-31\"\n",
    "harvest_date = \"2006-10-20\"\n",
    "max_duration = 300"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Constructing the cropping calendar"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\n",
      "- 2006-01-01:\n",
      "    CropCalendar:\n",
      "        crop_name: sugarbeet\n",
      "        variety_name: Sugarbeet_601\n",
      "        crop_start_date: 2006-03-31\n",
      "        crop_start_type: emergence\n",
      "        crop_end_date: 2006-10-20\n",
      "        crop_end_type: harvest\n",
      "        max_duration: 300\n",
      "    TimedEvents: null\n",
      "    StateEvents: null\n",
      "\n"
     ]
    }
   ],
   "source": [
    "# Here we define the agromanagement for sugar beet\n",
    "agro_yaml = \"\"\"\n",
    "- {start}:\n",
    "    CropCalendar:\n",
    "        crop_name: {cname}\n",
    "        variety_name: {vname}\n",
    "        crop_start_date: {startdate}\n",
    "        crop_start_type: emergence\n",
    "        crop_end_date: {enddate}\n",
    "        crop_end_type: harvest\n",
    "        max_duration: {maxdur}\n",
    "    TimedEvents: null\n",
    "    StateEvents: null\n",
    "\"\"\".format(cname=crop_name, vname=variety_name, \n",
    "           start=campaign_start_date, startdate=emergence_date, \n",
    "           enddate=harvest_date, maxdur=max_duration)\n",
    "agro = yaml.safe_load(agro_yaml)\n",
    "print(agro_yaml)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Retrieving weather input and model parameters"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Unknown parameter values provided to WOFOSTSiteDataProvider: {'CO2': 360.0}\n"
     ]
    }
   ],
   "source": [
    "# Weather data for Netherlands\n",
    "wdp = NASAPowerWeatherDataProvider(latitude=latitude, longitude=longitude)\n",
    "\n",
    "# Parameter sets for crop, soil and site\n",
    "# Standard crop parameter library\n",
    "cropd = YAMLCropDataProvider()\n",
    "# We don't need soil for potential production, so we use dummy values\n",
    "soild = DummySoilDataProvider()\n",
    "# Some site parameters\n",
    "sited = WOFOST72SiteDataProvider(WAV=50, CO2=360.)\n",
    "\n",
    "# Retrieve all parameters in the form of a single object. \n",
    "# In order to see all parameters for the selected crop already, we\n",
    "# synchronise data provider cropd with the crop/variety: \n",
    "firstkey = list(agro[0])[0]\n",
    "cropcalendar = agro[0][firstkey]['CropCalendar'] \n",
    "cropd.set_active_crop(cropcalendar['crop_name'], cropcalendar['variety_name'])\n",
    "params = ParameterProvider(cropdata=cropd, sitedata=sited, soildata=soild)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Defining the analysis\n",
    "Sensitivity analysis involves that we vary the model parameters a bit, to see how sensitive the model is to inaccuracies / variations in those values. We focus on variations in the model parameters and not on inaccuracies in driving variables usch as rainfall or temperature. \n",
    "\n",
    "In the following, intervals are defined for parameters. It is necessary to acquire knowledge of what are realistic values for the parameters which is often crop specific and region specific."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Define the target variable\n",
    "target_variable = \"TWSO\"\n",
    "\n",
    "# For each scalar parameter, determine a sensible interval \n",
    "problem_yaml = \"\"\"\n",
    "    num_vars: 5\n",
    "    names: \n",
    "    - TSUM1\n",
    "    - TSUM2\n",
    "    - SPAN\n",
    "    - Q10\n",
    "    - TDWI\n",
    "    bounds:\n",
    "    - [500, 800]\n",
    "    - [1200, 1600]\n",
    "    - [28, 37]\n",
    "    - [1.8, 2.2]\n",
    "    - [0.4, 0.6]\n",
    "\"\"\"\n",
    "problem = yaml.safe_load(problem_yaml)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## First and second order sensitivities\n",
    "When running a sensitivity analysis, we may compute first and second order sensitivities. Calculating second order sensitivities may require a large number of models runs if the number of parameters is large.\n",
    "\n",
    "In this examples we use the Saltelli sampler to sample a number of combinations of parameter values and we will derive first and second order sensitivity. The variable `nsamples` indicates the number of samples to generate, higher values will provide better estimates of sensitivity but increase computation time."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "We are going to do 600 simulations\n"
     ]
    }
   ],
   "source": [
    "calc_second_order = True\n",
    "nsamples = 50\n",
    "paramsets = saltelli.sample(problem, nsamples, calc_second_order=calc_second_order)\n",
    "print(\"We are going to do %s simulations\" % len(paramsets))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Running the model\n",
    "Below we are going to loop over the parameter sets sampled by the Saltelli sampler."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "100%|████████████████████████████████████████████████████████████████████████████████| 600/600 [03:22<00:00,  2.96it/s]\n"
     ]
    }
   ],
   "source": [
    "# Loop over the samples of the parameter values and run WOFOST for each of the samples\n",
    "target_results = []\n",
    "with tqdm(total=len(paramsets)) as pbar:\n",
    "    for i, paramset in enumerate(paramsets):\n",
    "        params.clear_override()\n",
    "        for name, value in zip(problem[\"names\"], paramset):\n",
    "            params.set_override(name, value)\n",
    "        wofost = Wofost72_PP(params, wdp, agro)\n",
    "        wofost.run_till_terminate()\n",
    "        r = wofost.get_summary_output()\n",
    "        target_result = r[0][target_variable]\n",
    "        if target_result is None:\n",
    "            print(\"Target variable is not available in summary output!\")\n",
    "        target_results.append(target_result)\n",
    "        pbar.update()\n",
    "target_results = np.array(target_results)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Computing sensitivity indices\n",
    "We will use the [Sobol](https://en.wikipedia.org/wiki/Variance-based_sensitivity_analysis) method to assess the response of the model to the variations in the parameters. Besides Sobol, other known methods are Morris and FAST. See also: https://salib.readthedocs.io/en/latest/."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Use the Sobol method to analyse the results of the simulations\n",
    "Si = sobol.analyze(problem, target_results, calc_second_order=calc_second_order)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Plotting results\n",
    "After running the cell below, a chart with first order and total order sensitivity is drawn. The first order  sensitivity index for a parameter represents the fraction of the variance in the target variable caused by the variation in the input parameter. Total order indices represent the total variance in the target variable for a given parameter, including all variance caused by its interactions, of any order, with any other input variable. \n",
    "\n",
    "In many models there are statistical dependencies between parameters or physical constraints on the model. In such cases part of the variance in the output is also caused by interactions between variables which is not represented by the first order indices. In such cases, second order indices need to be calculated. If the total-order indices are substantially larger than the first-order indices, then it is likely that there are higher-order interactions occurring. It means that if both parameters are changed at the same time, it will have a large effect on the target variable.\n",
    "\n",
    "Note that indices do not necessarily add up to unity and often indices can become slightly negative due to sampling error. Such negative index values would disappear if the number of samples would be increased."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<Axes: xlabel='Parameter'>"
      ]
     },
     "execution_count": 16,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAzoAAAIvCAYAAAC83a4oAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjguMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/SrBM8AAAACXBIWXMAAA9hAAAPYQGoP6dpAAA7J0lEQVR4nO3df5xVBb3v//cMwgw/B5UANQR/cETKAEGN/BaQGGja1ewR/jiCyMHMuOIhLbkVoKZoKqFmailYpjeyo540D6UklgapKP4kbpmElgOaygBeGGXm+0fXqTkyygDOZpbP5+OxHw/22mvt/RlmifOatdfaZfX19fUBAAAokPJSDwAAALC9CR0AAKBwhA4AAFA4QgcAACgcoQMAABSO0AEAAApH6AAAAIWzU6kH2BJ1dXX561//ms6dO6esrKzU4wAAACVSX1+ftWvXZvfdd095edPHbVpF6Pz1r39Nr169Sj0GAACwg3j++efzwQ9+sMnHW0XodO7cOcnfv5guXbqUeBoAAKBUampq0qtXr4ZGaEqrCJ233q7WpUsXoQMAALzrKS0uRgAAABSO0AEAAApH6AAAAIXTKs7RAQDg/W3Tpk154403Sj0GLaBt27Zp06bNNj+P0AEAYIdVX1+f6urqvPbaa6UehRbUtWvX9OzZc5s+Q1PoAACww3orcrp3754OHTr48PiCq6+vz+uvv57Vq1cnSXbbbbetfi6hAwDADmnTpk0NkbPrrruWehxaSPv27ZMkq1evTvfu3bf6bWwuRgAAwA7prXNyOnToUOJJaGlvfc+35bysrQqdq6++On369EllZWUOOeSQPPTQQ02ue+ONN6asrKzRrbKycqsHBgDg/cXb1d5/tsf3vNmhM2/evEyZMiXTp0/Po48+mgEDBmTUqFEN76PbnC5duuTFF19suP35z3/epqEBAADeSbNDZ9asWZk4cWLGjx+f/v3759prr02HDh0yZ86cJrcpKytLz549G249evTYpqEBAADeSbMuRlBbW5slS5Zk6tSpDcvKy8szcuTILFq0qMnt1q1bl969e6euri4HHnhgLrroonzoQx9qcv2NGzdm48aNDfdramqaMyYAAAXX59yft9hrrbj409vtuYYPH56BAwdm9uzZ2+05W9KMGTNyxx13ZOnSpaUe5V0164jOyy+/nE2bNr3tiEyPHj1SXV292W3222+/zJkzJ//5n/+ZH/3oR6mrq8vHPvaxvPDCC02+zsyZM1NVVdVw69WrV3PGBACAkjrllFPedp56WVlZvvWtb+WCCy7YpucuKyvLHXfcsX0GLbD3/KprQ4cOzdixYzNw4MAMGzYst912Wz7wgQ/kuuuua3KbqVOnZs2aNQ23559//r0eEwAAtqvRo0c3Ok/9xRdfzODBg9O5c+cmt6mtrW3BCZu2LVc7K+Vz/7NmhU63bt3Spk2brFq1qtHyVatWpWfPnlv0HG3bts2gQYPyxz/+scl1Kioq0qVLl0Y3AABoTSoqKhqdp96zZ88cdthhOeussxrW6dOnTy644IKMHTs2Xbp0yWmnnZba2tpMmjQpu+22WyorK9O7d+/MnDmzYf0kOfbYY1NWVtZw/91cc8012WeffdKuXbvst99+uemmmxo9XlZWlmuuuSaf+cxn0rFjx1x44YVJkosvvjg9evRI586dM2HChGzYsOFtz3399ddn//33T2VlZfr165fvfve7DY+tWLEiZWVlmTdvXoYNG5bKysrcfPPNzfhb3HrNCp127dpl8ODBWbBgQcOyurq6LFiwIEOHDt2i59i0aVOefPLJbfqUUwAAKIrLLrssAwYMyGOPPZZvfOMbufLKK/Ozn/0sP/nJT7J8+fLcfPPNDUHz8MMPJ0nmzp2bF198seH+O7n99tszefLkfPnLX85TTz2VL3zhCxk/fnzuu+++RuvNmDEjxx57bJ588smceuqp+clPfpIZM2bkoosuyiOPPJLddtutUcQkyc0335xp06blwgsvzLJly3LRRRflG9/4Rn7wgx80Wu/cc8/N5MmTs2zZsowaNWob/ra2XLMuRpAkU6ZMybhx4zJkyJAcfPDBmT17dtavX5/x48cnScaOHZs99tijoTrPP//8fPSjH82+++6b1157LZdeemn+/Oc/59/+7d+271cCAAA7kLvuuiudOnVquH/EEUdsdr1PfvKT+fKXv9xwf+XKlenbt2/+v//v/0tZWVl69+7d8NgHPvCBJEnXrl23+B1Vl112WU455ZScccYZSf7+8/zixYtz2WWXZcSIEQ3rnXjiiQ0/0yfJ8ccfnwkTJmTChAlJkm9+85u59957Gx3VmT59ei6//PJ89rOfTZLstddeeeaZZ3Lddddl3LhxDeudddZZDeu0lGaHzpgxY/LSSy9l2rRpqa6uzsCBAzN//vyGCxSsXLky5eX/OFD06quvZuLEiamurs7OO++cwYMH57e//W369++//b4KAADYwYwYMSLXXHNNw/2OHTvmhBNOeNt6Q4YMaXT/lFNOyeGHH5799tsvo0ePzlFHHZVPfepTWz3HsmXLctpppzVaduihh+aKK654xzmWLVuW008/vdGyoUOHNhwJWr9+fZ599tlMmDAhEydObFjnzTffTFVV1Ts+d0todugkyaRJkzJp0qTNPrZw4cJG97/97W/n29/+9ta8DAAAtFodO3bMvvvuu0Xr/bMDDzwwzz33XP7rv/4r9957bz7/+c9n5MiR+elPf/pejbrZOd7NunXrkiTf//73c8ghhzR6rE2bNtv03NvDVoUOrcSMqndfp1RmrCn1BAAAO6wuXbpkzJgxGTNmTD73uc9l9OjReeWVV7LLLrukbdu22bRp0xY/1/77758HH3yw0VvJHnzwwXd9h9X++++f3/3udxk7dmzDssWLFzf8uUePHtl9993zpz/9KSeddFIzvrqWIXQAAGAHMmvWrOy2224ZNGhQysvLc+utt6Znz57p2rVrkr9feW3BggU59NBDU1FRkZ133vkdn++cc87J5z//+QwaNCgjR47MnXfemdtuuy333nvvO243efLknHLKKRkyZEgOPfTQ3HzzzXn66aez9957N6xz3nnn5cwzz0xVVVVGjx6djRs35pFHHsmrr76aKVOmbPPfxbYQOgAAtDorLv50qUd4z3Tu3Dnf+ta38oc//CFt2rTJQQcdlLvvvrvhPPjLL788U6ZMyfe///3sscceWbFixTs+3zHHHJMrrrgil112WSZPnpy99torc+fOzfDhw99xuzFjxuTZZ5/NV77ylWzYsCHHHXdcvvjFL+YXv/hFwzr/9m//lg4dOuTSSy/NOeeck44dO+aAAw5odAntUimrr6+vL/UQ76ampiZVVVVZs2aNz9RpDm9dAwBasQ0bNuS5557LXnvtlcrKylKPQwt6p+/9lrZBsz5HBwAAoDUQOgAA0Ip96EMfSqdOnTZ7u/nmm0s9Xsk4RwcAAFqxu+++O2+88cZmH3vrsy7fj4QOAAC0Yr179y71CDskb10DAAAKR+gAAACFI3QAAIDCEToAAEDhCB0AAKBwXHUNAIDWZ0ZVC77WmpZ7ra1UVlaW22+/Pcccc0yLveaKFSuy11575bHHHsvAgQNb7HW3lCM6AACwnZSVlb3jbcaMGU1uu2LFipSVlWXp0qUtNm+ROaIDAADbyYsvvtjw53nz5mXatGlZvnx5w7JOnTqVYqxtUltbm3bt2rW653ZEBwAAtpOePXs23KqqqlJWVtZwv3v37pk1a1Y++MEPpqKiIgMHDsz8+fMbtt1rr72SJIMGDUpZWVmGDx+eJHn44Ydz+OGHp1u3bqmqqsqwYcPy6KOPbvWMTz75ZD75yU+mffv22XXXXXPaaadl3bp1DY+fcsopOeaYY3LhhRdm9913z3777ZckeeihhzJo0KBUVlZmyJAheeyxx9723E899VSOOOKIdOrUKT169MjJJ5+cl19+ueHx4cOHZ9KkSTnrrLPSrVu3jBo1aqu/jncjdAAAoAVcccUVufzyy3PZZZfliSeeyKhRo/KZz3wmf/jDH5L8PSSS5N57782LL76Y2267LUmydu3ajBs3Lg888EAWL16cvn375sgjj8zatWubPcP69eszatSo7Lzzznn44Ydz66235t57782kSZMarbdgwYIsX74899xzT+66666sW7cuRx11VPr3758lS5ZkxowZOfvssxtt89prr+WTn/xkBg0alEceeSTz58/PqlWr8vnPf77Rej/4wQ/Srl27PPjgg7n22mub/TVsKW9dAwCAFnDZZZflq1/9ao4//vgkySWXXJL77rsvs2fPztVXX50PfOADSZJdd901PXv2bNjuk5/8ZKPn+d73vpeuXbvm/vvvz1FHHdWsGW655ZZs2LAhP/zhD9OxY8ckyXe+850cffTRueSSS9KjR48kSceOHXP99dc3vK3se9/7Xurq6nLDDTeksrIyH/rQh/LCCy/ki1/8YsNzf+c738mgQYNy0UUXNSybM2dOevXqlf/zf/5P/uVf/iVJ0rdv33zrW99q1txbwxEdAAB4j9XU1OSvf/1rDj300EbLDz300Cxbtuwdt121alUmTpyYvn37pqqqKl26dMm6deuycuXKZs+xbNmyDBgwoCFy3pqhrq6u0blEBxxwQKNzZ5YtW5aPfOQjqaysbFg2dOjQRs/9+OOP57777kunTp0abv369UuSPPvssw3rDR48uNlzbw1HdAAAYAc2bty4/O1vf8sVV1yR3r17p6KiIkOHDk1tbe179pr/HEJbat26dQ1Hhv673XbbbZuee2s4ogMAAO+xLl26ZPfdd8+DDz7YaPmDDz6Y/v37J0nDEZRNmza9bZ0zzzwzRx55ZD70oQ+loqKi0Qn+zbH//vvn8ccfz/r16xs9f3l5ecNFB5ra7oknnsiGDRsali1evLjROgceeGCefvrp9OnTJ/vuu2+jW0vFzT8TOgAA0ALOOeecXHLJJZk3b16WL1+ec889N0uXLs3kyZOTJN27d0/79u0bTuJfs+bvH1Tat2/f3HTTTVm2bFl+97vf5aSTTkr79u23aoaTTjoplZWVGTduXJ566qncd999+Z//83/m5JNPbjg/Z3NOPPHElJWVZeLEiXnmmWdy991357LLLmu0zpe+9KW88sorOeGEE/Lwww/n2WefzS9+8YuMHz/+bfHWErx1DQCA1mfGmlJP0Gxnnnlm1qxZky9/+ctZvXp1+vfvn5/97Gfp27dvkmSnnXbKlVdemfPPPz/Tpk3Lxz/+8SxcuDA33HBDTjvttBx44IHp1atXLrroordd8WxLdejQIb/4xS8yefLkHHTQQenQoUOOO+64zJo16x2369SpU+68886cfvrpGTRoUPr3759LLrkkxx13XMM6bx2x+upXv5pPfepT2bhxY3r37p3Ro0envLzlj6+U1dfX17f4qzZTTU1NqqqqsmbNmnTp0qXU47QeM6pKPUHTWuE/TgBAy9qwYUOee+657LXXXo1Ogqf43ul7v6Vt4K1rAABA4QgdAAAoiIsuuqjR5Z3/+XbEEUeUerwW5RwdAAAoiNNPPz2f//znN/vY1l7AoLUSOgAAUBC77LJLdtlll1KPsUPw1jUAAHZoreDaWWxn2+N7LnQAANghtW3bNkny+uuvl3gSWtpb3/O39oGt4a1rAADskNq0aZOuXbtm9erVSf7+GTBlZWUlnor3Un19fV5//fWsXr06Xbt2TZs2bbb6uYQOAAA7rJ49eyZJQ+zw/tC1a9eG7/3WEjoAAOywysrKsttuu6V79+554403Sj0OLaBt27bbdCTnLUIHAIAdXps2bbbLD7+8f7gYAQAAUDiO6GyjPuf+vNQjNGlFZaknAACA0nBEBwAAKByhAwAAFI7QAQAACkfoAAAAhSN0AACAwhE6AABA4QgdAACgcIQOAABQOEIHAAAoHKEDAAAUjtABAAAKR+gAAACFI3QAAIDCEToAAEDhCB0AAKBwhA4AAFA4QgcAACgcoQMAABSO0AEAAApH6AAAAIUjdAAAgMIROgAAQOEIHQAAoHCEDgAAUDhCBwAAKByhAwAAFI7QAQAACkfoAAAAhSN0AACAwhE6AABA4QgdAACgcIQOAABQOEIHAAAoHKEDAAAUjtABAAAKR+gAAACFI3QAAIDCEToAAEDhCB0AAKBwhA4AAFA4QgcAACgcoQMAABTOVoXO1VdfnT59+qSysjKHHHJIHnrooS3a7sc//nHKyspyzDHHbM3LAgAAbJFmh868efMyZcqUTJ8+PY8++mgGDBiQUaNGZfXq1e+43YoVK3L22Wfn4x//+FYPCwAAsCWaHTqzZs3KxIkTM378+PTv3z/XXnttOnTokDlz5jS5zaZNm3LSSSflvPPOy957771NAwMAALybZoVObW1tlixZkpEjR/7jCcrLM3LkyCxatKjJ7c4///x07949EyZM2KLX2bhxY2pqahrdAAAAtlSzQufll1/Opk2b0qNHj0bLe/Tokerq6s1u88ADD+SGG27I97///S1+nZkzZ6aqqqrh1qtXr+aMCQAAvM+9p1ddW7t2bU4++eR8//vfT7du3bZ4u6lTp2bNmjUNt+eff/49nBIAACianZqzcrdu3dKmTZusWrWq0fJVq1alZ8+eb1v/2WefzYoVK3L00Uc3LKurq/v7C++0U5YvX5599tnnbdtVVFSkoqKiOaMBAAA0aNYRnXbt2mXw4MFZsGBBw7K6urosWLAgQ4cOfdv6/fr1y5NPPpmlS5c23D7zmc9kxIgRWbp0qbekAQAA74lmHdFJkilTpmTcuHEZMmRIDj744MyePTvr16/P+PHjkyRjx47NHnvskZkzZ6aysjIf/vCHG23ftWvXJHnbcgAAgO2l2aEzZsyYvPTSS5k2bVqqq6szcODAzJ8/v+ECBStXrkx5+Xt66g8AAMA7Kquvr68v9RDvpqamJlVVVVmzZk26dOlS6nEa6XPuz0s9QpNWVJ5Y6hGaNmNNqScAAKAV2tI2cOgFAAAoHKEDAAAUjtABAAAKR+gAAACFI3QAAIDCEToAAEDhCB0AAKBwhA4AAFA4QgcAACgcoQMAABSO0AEAAApH6AAAAIUjdAAAgMIROgAAQOEIHQAAoHCEDgAAUDhCBwAAKByhAwAAFI7QAQAACkfoAAAAhSN0AACAwhE6AABA4QgdAACgcIQOAABQOEIHAAAoHKEDAAAUjtABAAAKR+gAAACFI3QAAIDCEToAAEDhCB0AAKBwhA4AAFA4QgcAACgcoQMAABSO0AEAAApH6AAAAIUjdAAAgMIROgAAQOEIHQAAoHCEDgAAUDhCBwAAKByhAwAAFI7QAQAACkfoAAAAhSN0AACAwhE6AABA4QgdAACgcIQOAABQOEIHAAAoHKEDAAAUjtABAAAKR+gAAACFI3QAAIDCEToAAEDhCB0AAKBwhA4AAFA4QgcAACgcoQMAABSO0AEAAApH6AAAAIUjdAAAgMIROgAAQOEIHQAAoHCEDgAAUDhCBwAAKByhAwAAFI7QAQAACkfoAAAAhSN0AACAwhE6AABA4QgdAACgcIQOAABQOEIHAAAoHKEDAAAUjtABAAAKR+gAAACFI3QAAIDCEToAAEDhbFXoXH311enTp08qKytzyCGH5KGHHmpy3dtuuy1DhgxJ165d07FjxwwcODA33XTTVg8MAADwbpodOvPmzcuUKVMyffr0PProoxkwYEBGjRqV1atXb3b9XXbZJV/72teyaNGiPPHEExk/fnzGjx+fX/ziF9s8PAAAwOY0O3RmzZqViRMnZvz48enfv3+uvfbadOjQIXPmzNns+sOHD8+xxx6b/fffP/vss08mT56cj3zkI3nggQe2eXgAAIDNaVbo1NbWZsmSJRk5cuQ/nqC8PCNHjsyiRYvedfv6+vosWLAgy5cvzyc+8Ykm19u4cWNqamoa3QAAALZUs0Ln5ZdfzqZNm9KjR49Gy3v06JHq6uomt1uzZk06deqUdu3a5dOf/nSuuuqqHH744U2uP3PmzFRVVTXcevXq1ZwxAQCA97mdWuJFOnfunKVLl2bdunVZsGBBpkyZkr333jvDhw/f7PpTp07NlClTGu7X1NSIHYAW0ufcn5d6hCatuPjTpR4BgFaiWaHTrVu3tGnTJqtWrWq0fNWqVenZs2eT25WXl2ffffdNkgwcODDLli3LzJkzmwydioqKVFRUNGc0AACABs1661q7du0yePDgLFiwoGFZXV1dFixYkKFDh27x89TV1WXjxo3NeWkAAIAt1uy3rk2ZMiXjxo3LkCFDcvDBB2f27NlZv359xo8fnyQZO3Zs9thjj8ycOTPJ38+3GTJkSPbZZ59s3Lgxd999d2666aZcc8012/crAQAA+H+aHTpjxozJSy+9lGnTpqW6ujoDBw7M/PnzGy5QsHLlypSX/+NA0fr163PGGWfkhRdeSPv27dOvX7/86Ec/ypgxY7bfVwEAAPBPyurr6+tLPcS7qampSVVVVdasWZMuXbqUepxGduiTditPLPUITZuxptQTAE3Yof9dczECgPe9LW2DZn9gKAAAwI5O6AAAAIUjdAAAgMIROgAAQOEIHQAAoHCEDgAAUDhCBwAAKByhAwAAFI7QAQAACkfoAAAAhSN0AACAwhE6AABA4QgdAACgcIQOAABQOEIHAAAoHKEDAAAUjtABAAAKR+gAAACFI3QAAIDCEToAAEDhCB0AAKBwhA4AAFA4QgcAACgcoQMAABSO0AEAAApH6AAAAIUjdAAAgMIROgAAQOEIHQAAoHCEDgAAUDhCBwAAKByhAwAAFI7QAQAACkfoAAAAhSN0AACAwhE6AABA4QgdAACgcIQOAABQODuVegAA2GIzqko9QdNmrCn1BAD8E0d0AACAwhE6AABA4QgdAACgcIQOAABQOEIHAAAoHKEDAAAUjtABAAAKR+gAAACFI3QAAIDCEToAAEDhCB0AAKBwhA4AAFA4QgcAACgcoQMAABSO0AEAAApH6AAAAIUjdAAAgMIROgAAQOEIHQAAoHCEDgAAUDhCBwAAKByhAwAAFI7QAQAACkfoAAAAhSN0AACAwhE6AABA4QgdAACgcIQOAABQOEIHAAAoHKEDAAAUjtABAAAKR+gAAACFI3QAAIDCEToAAEDhCB0AAKBwhA4AAFA4QgcAACgcoQMAABSO0AEAAApH6AAAAIWzVaFz9dVXp0+fPqmsrMwhhxyShx56qMl1v//97+fjH/94dt555+y8884ZOXLkO64PAACwrZodOvPmzcuUKVMyffr0PProoxkwYEBGjRqV1atXb3b9hQsX5oQTTsh9992XRYsWpVevXvnUpz6Vv/zlL9s8PAAAwOY0O3RmzZqViRMnZvz48enfv3+uvfbadOjQIXPmzNns+jfffHPOOOOMDBw4MP369cv111+furq6LFiwYJuHBwAA2JxmhU5tbW2WLFmSkSNH/uMJysszcuTILFq0aIue4/XXX88bb7yRXXbZpcl1Nm7cmJqamkY3AACALdWs0Hn55ZezadOm9OjRo9HyHj16pLq6eoue46tf/Wp23333RrH0382cOTNVVVUNt169ejVnTAAA4H2uRa+6dvHFF+fHP/5xbr/99lRWVja53tSpU7NmzZqG2/PPP9+CUwIAAK3dTs1ZuVu3bmnTpk1WrVrVaPmqVavSs2fPd9z2sssuy8UXX5x77703H/nIR95x3YqKilRUVDRnNAAAgAbNOqLTrl27DB48uNGFBN66sMDQoUOb3O5b3/pWLrjggsyfPz9DhgzZ+mkBAAC2QLOO6CTJlClTMm7cuAwZMiQHH3xwZs+enfXr12f8+PFJkrFjx2aPPfbIzJkzkySXXHJJpk2blltuuSV9+vRpOJenU6dO6dSp03b8UgAAAP6u2aEzZsyYvPTSS5k2bVqqq6szcODAzJ8/v+ECBStXrkx5+T8OFF1zzTWpra3N5z73uUbPM3369MyYMWPbpgcAANiMZodOkkyaNCmTJk3a7GMLFy5sdH/FihVb8xIAAABbrUWvugYAANAShA4AAFA4QgcAACgcoQMAABSO0AEAAApH6AAAAIUjdAAAgMIROgAAQOEIHQAAoHCEDgAAUDhCBwAAKByhAwAAFI7QAQAACkfoAAAAhSN0AACAwhE6AABA4QgdAACgcIQOAABQOEIHAAAoHKEDAAAUjtABAAAKR+gAAACFI3QAAIDCEToAAEDhCB0AAKBwhA4AAFA4QgcAACgcoQMAABSO0AEAAApH6AAAAIUjdAAAgMIROgAAQOEIHQAAoHCEDgAAUDhCBwAAKByhAwAAFI7QAQAACkfoAAAAhSN0AACAwhE6AABA4QgdAACgcIQOAABQOEIHAAAoHKEDAAAUjtABAAAKR+gAAACFI3QAAIDCEToAAEDhCB0AAKBwhA4AAFA4QgcAACgcoQMAABSO0AEAAApH6AAAAIUjdAAAgMIROgAAQOEIHQAAoHCEDgAAUDhCBwAAKByhAwAAFI7QAQAACkfoAAAAhSN0AACAwhE6AABA4QgdAACgcIQOAABQOEIHAAAoHKEDAAAUjtABAAAKR+gAAACFI3QAAIDCEToAAEDhCB0AAKBwhA4AAFA4QgcAACgcoQMAABSO0AEAAApH6AAAAIWzVaFz9dVXp0+fPqmsrMwhhxyShx56qMl1n3766Rx33HHp06dPysrKMnv27K2dFQAAYIs0O3TmzZuXKVOmZPr06Xn00UczYMCAjBo1KqtXr97s+q+//nr23nvvXHzxxenZs+c2DwwAAPBumh06s2bNysSJEzN+/Pj0798/1157bTp06JA5c+Zsdv2DDjool156aY4//vhUVFRs88AAAADvplmhU1tbmyVLlmTkyJH/eILy8owcOTKLFi3abkNt3LgxNTU1jW4AAABbqlmh8/LLL2fTpk3p0aNHo+U9evRIdXX1dhtq5syZqaqqarj16tVruz03AABQfDvkVdemTp2aNWvWNNyef/75Uo8EAAC0Ijs1Z+Vu3bqlTZs2WbVqVaPlq1at2q4XGqioqHA+DwAAsNWadUSnXbt2GTx4cBYsWNCwrK6uLgsWLMjQoUO3+3AAAABbo1lHdJJkypQpGTduXIYMGZKDDz44s2fPzvr16zN+/PgkydixY7PHHntk5syZSf5+AYNnnnmm4c9/+ctfsnTp0nTq1Cn77rvvdvxSAAAA/q7ZoTNmzJi89NJLmTZtWqqrqzNw4MDMnz+/4QIFK1euTHn5Pw4U/fWvf82gQYMa7l922WW57LLLMmzYsCxcuHDbvwIAAID/ptmhkySTJk3KpEmTNvvYf4+XPn36pL6+fmteBgAAYKvskFddAwAA2BZCBwAAKByhAwAAFI7QAQAACkfoAAAAhSN0AACAwhE6AABA4QgdAACgcIQOAABQOEIHAAAoHKEDAAAUjtABAAAKR+gAAACFI3QAAIDCEToAAEDh7FTqAYACmFFV6gmaNmNNqScAAErAER0AAKBwhA4AAFA4QgcAACgc5+hAK9Hn3J+XeoQmrags9QQAAI05ogMAABSO0AEAAApH6AAAAIUjdAAAgMIROgAAQOEIHQAAoHCEDgAAUDhCBwAAKByhAwAAFI7QAQAACkfoAAAAhSN0AACAwhE6AABA4QgdAACgcIQOAABQOEIHAAAoHKEDAAAUjtABAAAKR+gAAACFs1OpBwAA3p/6nPvzUo/QpBUXf7rUIwDbyBEdAACgcIQOAABQOEIHAAAoHKEDAAAUjtABAAAKR+gAAACFI3QAAIDCEToAAEDhCB0AAKBwhA4AAFA4QgcAACgcoQMAABSO0AEAAApH6AAAAIUjdAAAgMIROgAAQOEIHQAAoHCEDgAAUDhCBwAAKByhAwAAFI7QAQAACkfoAAAAhSN0AACAwtmp1AMAAOxwZlSVeoKmzVhT6gmgVXBEBwAAKByhAwAAFI7QAQAACkfoAAAAhSN0AACAwhE6AABA4QgdAACgcIQOAABQOEIHAAAoHKEDAAAUjtABAAAKR+gAAACFI3QAAIDCEToAAEDhbFXoXH311enTp08qKytzyCGH5KGHHnrH9W+99db069cvlZWVOeCAA3L33Xdv1bAAAABbotmhM2/evEyZMiXTp0/Po48+mgEDBmTUqFFZvXr1Ztf/7W9/mxNOOCETJkzIY489lmOOOSbHHHNMnnrqqW0eHgAAYHN2au4Gs2bNysSJEzN+/PgkybXXXpuf//znmTNnTs4999y3rX/FFVdk9OjROeecc5IkF1xwQe6555585zvfybXXXruN4wMAwDvrc+7PSz1Ck1Zc/OlSj1BYzQqd2traLFmyJFOnTm1YVl5enpEjR2bRokWb3WbRokWZMmVKo2WjRo3KHXfc0eTrbNy4MRs3bmy4v2bNmiRJTU1Nc8ZtEXUbXy/1CE2qKasv9QhN2wG/lzs6+9pWsq81m31tK9nXms2+tpXsa822Q+9rvp/N9tbfWX39O/932qzQefnll7Np06b06NGj0fIePXrk97///Wa3qa6u3uz61dXVTb7OzJkzc955571tea9evZoz7vteVakHeCcX79DT0Uw79HfTvlYoO/R3075WKDv0d9O+VihVs0s9Qeu1du3aVFU1/d9Ds9+61hKmTp3a6ChQXV1dXnnlley6664pKysr4WStR01NTXr16pXnn38+Xbp0KfU4FJh9jZZiX6Ol2NdoKfa1rVNfX5+1a9dm9913f8f1mhU63bp1S5s2bbJq1apGy1etWpWePXtudpuePXs2a/0kqaioSEVFRaNlXbt2bc6o/D9dunTxHw4twr5GS7Gv0VLsa7QU+1rzvdORnLc066pr7dq1y+DBg7NgwYKGZXV1dVmwYEGGDh262W2GDh3aaP0kueeee5pcHwAAYFs1+61rU6ZMybhx4zJkyJAcfPDBmT17dtavX99wFbaxY8dmjz32yMyZM5MkkydPzrBhw3L55Zfn05/+dH784x/nkUceyfe+973t+5UAAAD8P80OnTFjxuSll17KtGnTUl1dnYEDB2b+/PkNFxxYuXJlysv/caDoYx/7WG655ZZ8/etfz//6X/8rffv2zR133JEPf/jD2++r4G0qKioyffr0t70FELY3+xotxb5GS7Gv0VLsa++tsvp3uy4bAABAK9Osc3QAAABaA6EDAAAUjtABAAAKR+gAAACFI3QAAIDCEToAAEDhNPtzdGhdHn/88Rx44IHZtGlTqUehAO6+++7cdttt2WWXXXLqqaemX79+DY+9+uqrOe644/KrX/2qhBMCbJ2HHnooixYtSnV1dZKkZ8+eGTp0aA4++OAST0YRfPazn92i9W677bb3eJL3F6HzPuCjktgebrnllowdOzajR4/O8uXLc9VVV+X666/PSSedlCSpra3N/fffX+IpKYoRI0akrKzsHdcpKyvLggULWmgiimr16tU57rjj8uCDD2bPPfds+AD0VatW5d///d9z6KGH5j/+4z/SvXv3Ek9Ka1ZVVVXqEd6XhE4r926/IVizZs27/rAAW+LSSy/NrFmzcuaZZyZJfvKTn+TUU0/Nhg0bMmHChBJPR9EMHDiwycfWrl2bW265JRs3bmy5gSisM844I5s2bcqyZcuy3377NXps+fLlOfXUU/OlL30pt956a4kmpAjmzp1b6hHel8rq/bq/VWvbtm0OP/zwht9A/XevvPJK7rrrLm9dY5t16tQpTz75ZPbaa6+GZffdd18+85nP5NJLL82xxx6b3Xff3b7Ge+bNN9/M1VdfnQsvvDBVVVW54IILcvzxx5d6LFq5zp0759e//nUGDRq02ceXLFmS4cOHZ+3atS08GUUyd+7cHHbYYdlzzz1LPcr7iiM6rdz++++f4447rsnfqC9dujR33XVXC09FEXXp0iWrVq1qFDojRozIXXfdlaOOOiovvPBCCaej6G6++eZMmzYt//f//t/MmDEjp512Wnbayf/C2HYVFRWpqalp8vG1a9emoqKiBSeiiM4444zU1tamd+/eGTFiRMNtjz32KPVoheaqa63c4MGD8+ijjzb5eEVFhd8esF0cfPDB+a//+q+3LR82bFjuvPPOzJ49u+WHovDmz5+fgQMH5owzzsgpp5ySP/zhDznjjDNEDtvNmDFjMm7cuNx+++2Ngqempia33357xo8fnxNOOKGEE1IEr732Wu69996MHTs2f/zjHzNx4sTsueee2W+//XL66adn3rx5WbVqVanHLBxvXWvlNm7cmE2bNqVDhw6lHoWCu//++/Pb3/42U6dO3ezj9913X374wx96HzLbxUMPPZSvfvWrWbx4cU4//fR87WtfS7du3Uo9FgW0cePGnHXWWZkzZ07efPPNtGvXLsnfL7Cy0047ZcKECfn2t7/tqA7b1YYNG7Jo0aLcd999WbhwYR5++OG88cYbefPNN0s9WqEIHQB2OOXl5Wnfvn1OO+20Rm+X/O/eujgGbKuamposWbKk0eWlBw8enC5dupR4MoqotrY2ixYtyq9+9assXLgwv/vd77L77rvnT3/6U6lHKxShA8AOp0+fPlt0eWk/FACtQW1tbRYvXpyFCxfmV7/6VX73u9+ld+/e+cQnPpFPfOITGTZsWHr16lXqMQtH6LRybdq02aL1XAmLbWVfA96PVq1aleuuuy7Tpk0r9Si0Yu3bt0/37t1z9NFHZ9iwYfn4xz+enj17lnqswnM2ZytXX1+f3r17Z9y4cU1eGhO2B/saO5LXXnstP/rRjzJp0qRSj0LBVVdX57zzzhM6bJMBAwbksccey69//euUl5envLw8w4cPz6677lrq0QrNEZ1W7pFHHskNN9yQH//4x9lrr71y6qmn5qSTTsrOO+9c6tEoGPsaO4IFCxbkhhtuyO23354OHTrkb3/7W6lHopV74okn3vHx3//+9znhhBMcrWabrVu3Lg888EDDBQgee+yx/Mu//EuGDx+eYcOGZdiwYenevXupxywUoVMQGzZsyE9/+tPMnTs3ixcvztFHH50JEybk8MMPL/VoFIx9jZb2/PPPZ+7cuZk7d25WrlyZ448/PieffHIOO+ywtG3bttTj0cqVl5enrKwsm/tx6K3lZWVlQoftbu3atfnNb36Te+65J3Pnzs26detcdW07EzoF9Nxzz2XChAm5//7789JLL2WXXXYp9UgUlH2N98obb7yRO+64I9dff31+85vfZPTo0TnxxBNzwgkn5PHHH0///v1LPSIF0a1bt3zrW9/KYYcdttnHn3766Rx99NFCh+2mrq4uDz/8cBYuXJj77rsvDz74YNavX5/evXvnueeeK/V4heIcnQJ54YUXcuONN+bGG2/M66+/nnPOOcdlMXlP2Nd4r+2xxx7p169f/vVf/zU//vGPG94i6YMb2d4GDx6cv/71r+ndu/dmH3/ttdc2e7QHmuOhhx7KwoULs3DhwjzwwANZt25dPvjBD2b48OG58sorM2LEiPTp06fUYxaO0Gnlamtrc/vtt+eGG27Ib37zmxxxxBGZPXt2jjjiiC2+ShZsCfsaLenNN99MWVlZysrK7F+8p04//fSsX7++ycf33HNPH4TMNvvoRz+anj17ZsSIEZk1a1ZGjBiRffbZp9RjFZ63rrVyu+66azp37pxx48bl5JNPbvIkNr9tZ1vZ12hJGzZsyH/8x3/khhtuyOLFi3PEEUfkX//1XzNmzJgsXbrUW9d4z7z88stp166df8vYrpYvX5799tuv1GO87widVq68vLzhz5v7cD0nUbK92NcolWeffTZz5szJD3/4w/zlL3/JCSeckFNOOSWf/OQnHe1hu3jttdfyta99LfPmzcurr76aJPnABz6Q8ePH5xvf+EY6dOhQ4gkpivr6+ixZsiQrVqxIWVlZ9tprrwwaNOhdPyCZrSN0Wrn7779/i9YbNmzYezwJRWdfo9Tq6uryi1/8IjfccEPuvPPOdO7cOS+//HKpx6KVe+WVVzJ06ND85S9/yUknnZT9998/SfLMM8/klltuSb9+/fLAAw/kiSeeyOLFi3PmmWeWeGJaq/vuuy8TJkzIn//854bzvt6KnTlz5uQTn/hEiScsHqEDwA5p8eLFufPOO1NbW5vDDjsso0ePbnjspZdeyk033ZQpU6aUcEKK4KyzzsqCBQty7733pkePHo0eq66uzqc+9anst99++eUvf5krr7wy48aNK9GktGZ//OMfM2DAgBxyyCGZPHly+vXrl/r6+jzzzDO58sor88gjj+SJJ57I3nvvXepRC0XotHJNfdBZVVVV9txzT4dC2W7sa7Skn/70pxkzZkzat2+ftm3bpqamJpdccknOPvvsUo9GwfTp0yfXXXddRo0atdnH58+fnyOPPDLTp0/P9OnTW3g6imLSpElZtmxZFixY8LbH6uvrM3LkyPTv3z9XXXVVCaYrLqHTyjX1QWdlZWWprKzMWWedlfPPP9/72Nlm9jVa0uDBg3PQQQfl6quvTps2bTJz5sxceumleeWVV0o9GgVTUVGRZ599Nh/84Ac3+/gLL7yQPn36+CBHtsmHP/zhzJw5M0cfffRmH7/zzjszderUPPXUUy08WbG5vHQr19QHS7322mtZsmRJvvGNb2TnnXf2W1C2mX2NlrR8+fLMmzevIZy//OUvZ9q0aVm9enWTV/yDrdGtW7esWLGiydB57rnn7HNss5UrV+aAAw5o8vEPf/jD+fOf/9yCE70/OKJTcD/96U9z3nnn5cknnyz1KBScfY3tqby8PNXV1Y1+wOzcuXMef/xx72Fnuzr11FPz7LPP5p577km7du0aPbZx48aMGjUqe++9d+bMmVOiCSmCzf2b9s9WrVqV3Xff3ZVLtzNHdApu8ODBTf4mHrYn+xrb2/XXX59OnTo13H/zzTdz4403plu3bg3LXAGLbXX++ednyJAh6du3b770pS81nCS+bNmyfPe7383GjRvzwx/+sNRjUgDPPPNMqqurN/uYK0i+NxzRKbhFixblxBNP9AMo7zn7GttTnz593vUCF2VlZfnTn/7UQhNRZM8991zOOOOM/PKXv2x02d/DDz883/nOd7LvvvuWeEJau6bOc03SsNxn0W1/QqfAXnrppZxwwgnZc889HXLnPWVfA4rg1VdfzR/+8Ickyb777ptddtmlxBNRFFt6/k3v3r3f40neX4ROK9fUp+muWbMmL7zwQsO1/3v27FmC6SgS+xotadGiRfnb3/6Wo446qmHZD3/4w0yfPj3r16/PMccck6uuuioVFRUlnBJgy5x//vk5++yz06FDh1KP8r4idFq58847b7PLu3Tpkv322y+jRo1yuV+2C/saLWn06NEZMWJEvvrVryZJnnzyyRx44IE55ZRTsv/+++fSSy/NF77whcyYMaO0gwJsgTZt2uTFF190Bb8WJnQA2OHstttuufPOOzNkyJAkyde+9rXcf//9eeCBB5Ikt956a6ZPn55nnnmmlGMCbJF3u+oa7w1XXSugDRs2ZN68eVm/fn0OP/zw9O3bt9QjUVD2Nd4rr776anr06NFw//77788RRxzRcP+ggw7K888/X4rRALbKu11ghe3PEZ1WbsqUKXnjjTdy1VVXJUlqa2tz8MEH55lnnkmHDh3y5ptv5p577snQoUNLPCmtnX2NltS7d+/cdNNN+cQnPpHa2tp07do1d955Zw477LAkf38r27Bhw/LKK6+UeFKAd1deXp6qqqp3jR3/pm1fjui0cr/85S9z0UUXNdy/+eabs3LlyvzhD3/InnvumVNPPTXf/OY38/Of/7yEU1IE9jVa0pFHHplzzz03l1xySe6444506NAhH//4xxsef+KJJ7LPPvuUcEKA5jnvvPNSVVVV6jHeV4ROK7dy5cr079+/4f4vf/nLfO5zn2u4POHkyZNz5JFHlmo8CsS+Rku64IIL8tnPfjbDhg1Lp06d8oMf/KDRp9bPmTMnn/rUp0o4IUDzHH/88c7RaWFCp5UrLy9v9OFTixcvzje+8Y2G+127ds2rr75aitEoGPsaLalbt2759a9/nTVr1qRTp05vu6Lfrbfemk6dOpVoOoDmcX5OaZSXegC2zf77758777wzSfL0009n5cqVGTFiRMPjf/7znxud0Atby75GKVRVVW32suW77LJLoyM8ADsyp8SXhiM6rdxXvvKVHH/88fn5z3+ep59+OkceeWT22muvhsfvvvvuHHzwwSWckKKwrwHA1qmrqyv1CO9Ljui0cscee2zuvvvufOQjH8m///u/Z968eY0e79ChQ84444wSTUeR2NcAgNbE5aVbufPPPz9nn312OnToUOpRKDj7GgDQmgidVq5NmzZ58cUXXcWD95x9DQBoTbx1rZXTqbQU+xoA0JoInQJwyUJain0NAGgtvHWtlSsvL09VVdW7/gD6yiuvtNBEFJV9DQBoTVxeugDOO++8VFVVlXoM3gfsawBAa+GITitXXl6e6upqJ4jznrOvAQCtiXN0WjnnTNBS7GsAQGsidFo5B+RoKfY1AKA18dY1AACgcBzRAQAACkfoAAAAhSN0AACAwhE6AABA4QgdAACgcIQOAO/olFNOSVlZWcrKytKuXbvsu+++Of/88/Pmm2+WerStUlZWljvuuKPUYwDwHtup1AMAsOMbPXp05s6dm40bN+buu+/Ol770pbRt2zZTp05t1vNs2rQpZWVlKS9v/b9ne+ONN9K2bdtSjwFAE1r//2kAeM9VVFSkZ8+e6d27d774xS9m5MiR+dnPfpZZs2blgAMOSMeOHdOrV6+cccYZWbduXcN2N954Y7p27Zqf/exn6d+/fyoqKrJy5co8/PDDOfzww9OtW7dUVVVl2LBhefTRRxu9ZllZWa677rocddRR6dChQ/bff/8sWrQof/zjHzN8+PB07NgxH/vYx/Lss8822u4///M/c+CBB6aysjJ77713zjvvvIajT3369EmSHHvssSkrK2u4/27bvTXPNddck8985jPp2LFjLrzwwu38twzA9iR0AGi29u3bp7a2NuXl5bnyyivz9NNP5wc/+EF+9atf5Stf+UqjdV9//fVccskluf766/P000+ne/fuWbt2bcaNG5cHHnggixcvTt++fXPkkUdm7dq1jba94IILMnbs2CxdujT9+vXLiSeemC984QuZOnVqHnnkkdTX12fSpEkN6//mN7/J2LFjM3ny5DzzzDO57rrrcuONNzZEycMPP5wkmTt3bl588cWG+++23VtmzJiRY489Nk8++WROPfXU7f73CsB2VA8A72DcuHH1/+N//I/6+vr6+rq6uvp77rmnvqKiov7ss89+27q33npr/a677tpwf+7cufVJ6pcuXfqOr7Fp06b6zp071995550Ny5LUf/3rX2+4v2jRovok9TfccEPDsv/9v/93fWVlZcP9ww47rP6iiy5q9Nw33XRT/W677dboeW+//fZG62zpdmedddY7fh0A7DicowPAu7rrrrvSqVOnvPHGG6mrq8uJJ56YGTNm5N57783MmTPz+9//PjU1NXnzzTezYcOGvP766+nQoUOSpF27dvnIRz7S6PlWrVqVr3/961m4cGFWr16dTZs25fXXX8/KlSsbrffP2/Xo0SNJcsABBzRatmHDhtTU1KRLly55/PHH8+CDDzY6ErNp06a3zfTfbel2Q4YM2Zq/PgBKQOgA8K5GjBiRa665Ju3atcvuu++enXbaKStWrMhRRx2VL37xi7nwwguzyy675IEHHsiECRNSW1vbEAft27dPWVlZo+cbN25c/va3v+WKK65I7969U1FRkaFDh6a2trbRev98sv9bz7G5ZXV1dUmSdevW5bzzzstnP/vZt30NlZWVTX59W7pdx44dm3wOAHYsQgeAd9WxY8fsu+++jZYtWbIkdXV1ufzyyxuuovaTn/xki57vwQcfzHe/+90ceeSRSZLnn38+L7/88jbPeeCBB2b58uVvm/WftW3bNps2bWr2dgC0LkIHgK2y77775o033shVV12Vo48+Og8++GCuvfbaLdq2b9++uemmmzJkyJDU1NTknHPOSfv27bd5pmnTpuWoo47Knnvumc997nMpLy/P448/nqeeeirf/OY3k/z9ymsLFizIoYcemoqKiuy8885btB0ArYurrgGwVQYMGJBZs2blkksuyYc//OHcfPPNmTlz5hZte8MNN+TVV1/NgQcemJNPPjlnnnlmunfvvs0zjRo1KnfddVd++ctf5qCDDspHP/rRfPvb307v3r0b1rn88stzzz33pFevXhk0aNAWbwdA61JWX19fX+ohAAAAtidHdAAAgMIROgAAQOEIHQAAoHCEDgAAUDhCBwAAKByhAwAAFI7QAQAACkfoAAAAhSN0AACAwhE6AABA4QgdAACgcP5/k7Pmja17WhYAAAAASUVORK5CYII=",
      "text/plain": [
       "<Figure size 1000x600 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "df_sensitivity = pd.DataFrame({\n",
    "    \"Parameter\": problem[\"names\"],\n",
    "    \"First_order\": Si[\"S1\"],\n",
    "    \"Total_order\": Si[\"ST\"]}).set_index(\"Parameter\")\n",
    "fig, axes = plt.subplots(figsize=(10,6))\n",
    "df_sensitivity.plot(kind=\"bar\", ax=axes)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Second order sensitivities\n",
    "Second order sensitivities deal with the combined sensitivity of two parameters. The results can therefore most easily be visualized through a matrix "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "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>TSUM1</th>\n",
       "      <th>TSUM2</th>\n",
       "      <th>SPAN</th>\n",
       "      <th>Q10</th>\n",
       "      <th>TDWI</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>TSUM1</th>\n",
       "      <td>NaN</td>\n",
       "      <td>0.097104</td>\n",
       "      <td>0.110361</td>\n",
       "      <td>0.112299</td>\n",
       "      <td>0.118296</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>TSUM2</th>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>0.018845</td>\n",
       "      <td>0.038701</td>\n",
       "      <td>0.033748</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>SPAN</th>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>-0.088541</td>\n",
       "      <td>-0.064614</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>Q10</th>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>-0.032382</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>TDWI</th>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "       TSUM1     TSUM2      SPAN       Q10      TDWI\n",
       "TSUM1    NaN  0.097104  0.110361  0.112299  0.118296\n",
       "TSUM2    NaN       NaN  0.018845  0.038701  0.033748\n",
       "SPAN     NaN       NaN       NaN -0.088541 -0.064614\n",
       "Q10      NaN       NaN       NaN       NaN -0.032382\n",
       "TDWI     NaN       NaN       NaN       NaN       NaN"
      ]
     },
     "execution_count": 17,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Show second-order sensitivities\n",
    "second_order = np.array(Si['S2'])\n",
    "pd.DataFrame(second_order, index=problem[\"names\"], columns=problem[\"names\"])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "or through an image showing the values in the matrix values graphically. From the results, it is clear that the second order interactions are small as all indices have a value < 0.1."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAArcAAAMjCAYAAAC/HTiWAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjguMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/SrBM8AAAACXBIWXMAAA9hAAAPYQGoP6dpAABfbElEQVR4nO3de3zO9f/H8efnunZy2pw3K8wphyIia/KNZZlDB30lK8KcvilFRE1C5Luor1C++SaibyolfPlWDsk6aCFaR0Q5FRvRNoxtdn1+f+y36+vKoXHZdX189rjfbp9b9rne1+d6fXx07bXXXtfrY5imaQoAAACwAYe/AwAAAAAuFZJbAAAA2AbJLQAAAGyD5BYAAAC2QXILAAAA2yC5BQAAgG2Q3AIAAMA2SG4BAABgGyS3AAAAsA2SWwAAANgGyS0AAAD0ySef6LbbblNkZKQMw9CyZcv+9DkpKSm67rrrFBwcrPr162v+/PlnrJk1a5aioqIUEhKi6Ohobdy48dIHfxqSWwAAAOj48eO69tprNWvWrGKt37Vrl7p27arY2FilpaVp+PDhGjhwoFatWuVes2jRIo0YMULjx4/Xli1bdO211yo+Pl4HDx4sqdOQYZqmWWJHBwAAwGXHMAwtXbpU3bp1O+eaxx57TO+9956+++47976EhARlZmZq5cqVkqTo6Ghdf/31evHFFyVJLpdLNWvW1EMPPaTHH3+8RGIPKJGjAgAAlBInT55UXl6ev8M4g2maMgzDY19wcLCCg4MvyfFTU1MVFxfnsS8+Pl7Dhw+XJOXl5Wnz5s1KSkpyP+5wOBQXF6fU1NRLEsPZkNwCAABcpJMnT6pCaLhO5Wf7O5QzlC9fXseOHfPYN378eE2YMOGSHD89PV3h4eEe+8LDw5Wdna0TJ07o999/V0FBwVnXbNu27ZLEcDYktwAAABcpLy9Pp/KzdXXzp+V0hvg7HLeCgpP6Pm2s9u3bp9DQUPf+S1W1tTKSWwAAAC85nSFyOsv4O4wzhIaGeiS3l1JERIQyMjI89mVkZCg0NFRlypSR0+mU0+k865qIiIgSiUliWgIAAIDXDMOQ4bDQ9ode25IQExOjtWvXeuxbs2aNYmJiJElBQUFq2bKlxxqXy6W1a9e615QEklsAAADo2LFjSktLU1pamqTCUV9paWnau3evJCkpKUl9+vRxr7///vv1888/a/To0dq2bZv++c9/6u2339YjjzziXjNixAjNmTNHCxYs0NatWzVkyBAdP35ciYmJJXYetCUAAABAX375pWJjY91fjxgxQpLUt29fzZ8/XwcOHHAnupJUp04dvffee3rkkUc0Y8YMXXnllXrllVcUHx/vXtOzZ08dOnRI48aNU3p6upo3b66VK1ee8SGzS4k5twAAABcpOztbYWFhah49Tc4A6/TcFpw6obQNI5SVlVViPbdWRVsCAAAAbIPkFgAAALZBzy0AAIC3HI7CzSqsFIuPld4zBwAAgO2Q3AIAAMA2aEsAAADwkuE0ZDhL/sYJxWWY1onF16jcAgAAwDZIbgEAAGAbtCUAAAB4yXAYMhzWaQWwUiy+RuUWAAAAtkFyCwAAANugLQEAAMBLTEuwDiq3AAAAsA2SWwAAANgGbQkAAADecjgKN6uwUiw+VnrPHAAAALZDcgsAAADboC0BAADAWw5rTUuQy0Kx+BiVWwAAANgGyS0AAABsg7YEAAAALxkOQ4bDOq0AVorF16jcAgAAwDZIbgEAAGAbtCUAAAB4yXBaa1qCwbQEAAAA4PJHcgsAAADboC0BAADAS0xLsA4qtwAAALANklsAAADYBm0JAAAA3nIaktNCNUOmJQAAAACXP5JbAAAA2AZtCQAAAF5iWoJ1ULkFAACAbZDcAgAAwDZoSwAAAPCS4TRkOK3TCmAwLQEAAAC4/JHcAgAAwDZoSwAAAPAS0xKsg8otAAAAbIPkFgAAALZBWwIAAICXDIfFpiUUWCcWX6NyCwAAANsguQUAAIBt0JYAAADgLYejcLMKK8XiY6X3zAEAAGA7JLcAAACwDdoSAAAAvMRNHKyDyi0AAABsg+QWAAAAtkFbAgAAgJcMp8Vu4mChWHyNyi0AAABsg+QWAAAAtkFbAgAAgJeYlmAdVG4BAABgGyS3AAAAsA3aEgAAALzldBRuVmGlWHys9J45AAAAbIfkFgAAALZBWwIAAICXmJZgHVRuAQAAYBsktwAAALAN2hIAAAC8ZDgMGU7rtALQlgAAAADYAMktAAAAbIO2BAAAAC8xLcE6qNwCAADANkhuAQAAYBu0JQAAAHjJcFpsWoKFYvE1KrcAAACwDZJbAAAA2AZtCQAAAF4yHA4ZDuvUDK0Ui6+V3jMHAACA7ZDcAgAAwDZoSwAAAPCSw2HIYaEbJ1gpFl+jcgsAAADbILkFAACAbdCWAAAA4CXDYa1WAKMUly9L8akDAADAbkhuSznTNP0dAgAAwCVDW0IpdeDAARUUFOjKK6+UaZoyDOv8KqU02rVrl1atWqW9e/cqJiZGsbGxKl++vL/DAgAUk+E0ZDit873USrH4GpXbUmjbtm2Kjo7W3/72N/38888yDIMKrh99++23atOmjd577z29+eabmjRpkp599lkVFBT4OzQAAC47JLelzK+//qqBAweqSpUqysnJ0ZgxY0hw/WjPnj3q3r27EhMTtWzZMu3cuVM333yzVq5cqfz8fH+HBwDAZYfktpT5+uuvFRAQoNmzZ+u+++7TgQMHPBJcl8vl7xBLjYKCAv3nP//RNddco+HDh0uSnE6nHnzwQf3888/aunWrfwMEABRb0U0crLSVVvTcljJdunRR2bJlFR0drejoaLlcLv373//WmDFjNHnyZNWrV8+jB7egoEBOp9PPUdtXWFiYOnXqpOrVq7v3BQUFqaCgQMePH/djZAAAXJ6o3JYC27dv18cff+z+ul27du4/Dxw4UH369NH+/fv1xBNPuCu4EyZM0KFDh0hsS8D27duVkpIip9Opnj17avDgwZLkrppXrVpVNWrUUHBwsPs5S5cuVXZ2tl/iBQCUHrNmzVJUVJRCQkIUHR2tjRs3nnNt+/btZRjGGVvXrl3da/r163fG4506dSrRc6Bya3Pp6elq3LixJOm9995T586d3Y+5XC45HA4NGDBAhmFowYIFeuKJJxQYGKjXX39dPXr0ULVq1fwVui2dfj1WrFihrl27yjRNuVwu9w8SBQUFys3NVV5eniRp7Nixeumll7RlyxaFhob6LXYApU9OTo7Kli3r7zAuCw7DWq0AjouYgrRo0SKNGDFCs2fPVnR0tKZPn674+Hht377d4zeMRZYsWeL+XiVJhw8f1rXXXqsePXp4rOvUqZNeffVV99enF29KApVbmzNNU9dff70aNGigrl27avHixe4Pj53eY9u/f3/17t1bK1as0IoVK7RlyxZdffXVfo7efk6/Hrfddpv7ejgc//tf8dixY/r999/ldDr1zDPP6LnnntPq1atVu3ZtP0YOoLT59ddf1atXL61YscLfocBHpk2bpkGDBikxMVFNmjTR7NmzVbZsWc2bN++s6ytXrqyIiAj3tmbNGpUtW/aM5DY4ONhjXaVKlUr0PEhuba5SpUqqX7++GjdurEmTJunuu+/WkiVL3D21pye43333nZxOpz777DM1b97cj1Hb1/muR9G0irJly6pu3bpKSkrSxIkT9dlnn6lly5Z+jtze9uzZozlz5vg7DMBSDhw4oCNHjujFF1/UqlWr/B0OLlJ2drbHlpube9Z1eXl52rx5s+Li4tz7HA6H4uLilJqaWqzXmjt3rhISElSuXDmP/SkpKapevboaNmyoIUOG6PDhwxd/QsVAcmtzISEhGjt2rL755htFRUVp1KhRuuuuu7R06VJ3gutwOPTRRx/p1Vdf1bp166jYlqA/ux6maSovL08HDx7U5s2blZqaqlatWvk7bFs7deqUXnrpJf3973/XrFmz/B0O4Hc7d+7UokWL1KpVK02YMEEhISGaOnUqCe6f8PdkhHNNS6hZs6bCwsLcW3Jy8lnj/+2331RQUKDw8HCP/eHh4UpPT//T89+4caO+++47DRw40GN/p06d9Nprr2nt2rWaMmWKPv74Y3Xu3LlEZ7nTc2szP/74o7Zt26bY2FhVqFBBUuE/7I4dOyo9PV3jx49XTk6OunfvrnfffVd33nmnJOnmm2/WTz/9RI/tJXYx1yM4OFiDBw/W3XffrQYNGvj5DOwvICBA/fv318mTJzV37ly5XC499NBD/g4L8JvJkydr8eLF6tatm2JjY+VyuTRz5kxNnTpVkhQfH+/nCHEh9u3b5/F5jZLqd507d66aNm2q1q1be+xPSEhw/7lp06Zq1qyZ6tWrp5SUFHXo0KFEYiG5tZH09HQ1atRIkjRo0CCVL19ekydPVvny5dWxY0cNHDhQvXr10pQpU+RwOJSQkKAFCxa4/+GR2F5aF3M9Xn31Vd17770aM2YMt0T2oauuukoPPvigCgoK3B96IMFFaTVlyhR9/vnnGjBggF5//XV16NBBDodD06dPJ8G9DIWGhhbrw8hVq1aV0+lURkaGx/6MjAxFRESc97nHjx/XW2+9pYkTJ/7p69StW1dVq1bVzp07Syy5pS3BRnJzc3X77bdLkqpUqaIff/xRrVq1UlJSkho3bqzbb79dr7/+usqWLauxY8eqT58+euCBB3Ts2DHuTlYCLuZ6DB06VEePHuV6lLCdO3dq2bJl2rJli6TCD/o1aNBADz30kNq0aaNXXnlFM2fO9HOUgH9UrFhRf/vb3/Ttt9/qP//5jyQpNjZWjzzyiMqXL6+pU6dq9erVfo7SehwOh+W2CxEUFKSWLVtq7dq17n0ul0tr165VTEzMeZ/7zjvvKDc3V7179/7T1/nll190+PBh1ahR44LiuxCGyXfRy97hw4dVuXJlGYahXbt2adiwYUpLS9MXX3yh1NRUrVy5Uu+8847y8vLUsmVLffTRRwoMDDxnfw28w/WwtoMHDyoyMlIul0uhoaG68cYbddVVV6lfv36qX7++TNNUUlKSNm/erLvvvtt99ziUvKLxhPCtrVu36oUXXtCTTz6psLAwlS1bVnv37tVtt92mBg0aaPHixe61H3/8sWbMmKG9e/dqypQpJVZ5u5xkZ2crLCxMtz+5XIEh5f78CT6Sf/K4lk+6XVlZWcUeI7lo0SL17dtX//rXv9S6dWtNnz5db7/9trZt26bw8HD16dNHV1xxxRl9u3/5y190xRVX6K233vLYf+zYMT311FPq3r27IiIi9NNPP2n06NE6evSovv322xJrkaAt4TJ39OhRNWnSRL169dI//vEP1alTRy+88IL69euntm3bat26derevbvuu+8+LV68WDExMQoMDJRU+CsIXFpcD+urXr26evfurffee0/33nuvMjMz9eOPPyo2NlbVqlVTQkKCKleurGuuuUYLFixQUFCQHnjgAX+HDZSIU6dO6YknntCyZcuUl5en+vXr66677lL9+vU1c+ZM3XzzzXr55ZfdN5tp166d8vLy9Oqrr6p+/fp+jh6XWs+ePXXo0CGNGzdO6enpat68uVauXOkuuuzdu/eMH0C3b9+uzz777KzVfKfTqW+++UYLFixQZmamIiMj1bFjR02aNKlEZ91SubWBOXPm6OGHH9ajjz6qiRMnyjAM7d27V4mJidq6das++eQT1a9fXydOnFCZMmX8Ha7tcT2s6/TbSffq1Us7duzQkCFDlJiYqE2bNumbb77R/PnzlZOTo6+++kqS1KBBA23cuFFhYWH+DN3WZs6cqS+//FIZGRm6/fbb3VUe+MamTZvUqVMntWvXTo0bN9bixYs1a9YstW3bVk899ZQ++OADvfLKKx6TW3j/+p+iym23CSssV7ldNuG2C6rc2gWVWxsYNGiQQkJClJiYKEmaOHGiatWqpVdffVUDBgzQTTfdpM8++0x169bl134+wPWwJtM05XQ6lZ+fr8DAQC1cuFC9evVScnKyHA6H7rrrLl1//fW6++67ZZqmlixZoh07dqh3794ktiXoiSee0Msvv6yEhAQFBQVp1KhRWrt2rQYNGuRxR0WUnOuvv15jx47Vl19+qR49eqhKlSrq37+/Bg0aJMMwVLlyZa1Zs0bNmzeXYRhyOp0ktrA0ktvL2KlTpxQQUHgJ77vvPklSYmKiTNPU008/rVq1amnOnDkaMmSImjRpom3btikqKsqPEdsb18Oadu7cqd9//13NmjVTcHCwAgMDlZeXp6CgIC1cuFD33XefkpOTVVBQoB49erhHtvXr1899Jz+UjJ07d2rp0qVauHChOnbsKEm6//77NXToUM2YMUOmaapLly5+jtJ+tm/frlWrVumOO+5w3/mwRYsW+vjjj7Vr1y6NGDFCjRo10urVq7Vr1y6lpKQoLS1NPXv2VN26df0cPfDnKBldZvbs2aNHHnlEUuF8ztOHIN9333169dVXNWXKFE2fPl2SFBUVpVmzZqlLly7Kz8/3R8i2xvWwtkOHDumqq65Sly5d9MQTT7jHGBX1OUvSv//9b7Vs2VLPPvuslixZomPHjrkfI7EtWWXKlNHRo0d18uRJSYVtI61atdKcOXOUnZ2tl19+WXv27PFzlPbz4Ycfavjw4RoyZIj7xiXt27fXddddp4ceekjHjh1Tly5d9Oijj2rUqFFq0qSJMjMz/Rv0ZcDfN2w4100cSiOS28uIaZpavXq13n33Xd1///2SCpu1/5hQTZkyRY899pg2b94sqXCm3Ntvv80NAS4xrof1FX1ALC4uTldddZXmz5+v9u3ba/Lkyfrll1/c6xYuXKhWrVpp1KhR+s9//sMoNh9yOp3aunWr++uCggI1bdpUL774oj766CMtWrTIj9HZ04MPPqj169erUqVKmjRpkm699Vbt2LFD48aNU7t27TR48GAVFBToyiuvVNu2bbVp0ybt27ePqi0uGyS3lwmXyyXDMHTXXXfpscce0xdffOG+xd0fE6oePXroyiuv9Kh4FP26HJcG18P6Tp06JUm64YYbFBoaqsGDB2vLli3q1q2bdu7cqebNm2vKlCl6//33JRVWcBMSEhQdHU3F1keuuOIKjR49WmPGjNHKlSvdH/Y7deqUrrvuOg0dOlRvvvmmTpw4IZfL5edo7eH48ePKz89XTEyMpk+frvnz5+unn35S586d9cQTTyg+Pl6VKlXSypUrJUn5+fkqU6aMrrjiCj9HDhQfye1lYPfu3UpOTtaOHTtUqVIl9enTRwMGDNCXX37pkVAVfTMPDg5WtWrVVLZsWX+GbVtcj8tD0Q8Qt99+uxYvXqwZM2YoJCREw4cP19GjR3Xq1Cl9+umnGjRokDp27KgNGzZo5syZjDcqQb///rv7z0XV8d69e2vAgAH661//qlWrVsnpdLqvXbly5VSjRg2VKVOGD156Yffu3Zo4caL+8pe/KDo6Wh06dNAnn3yiSpUqqVOnTvruu+/UtWtXffLJJxo+fLiWL1/unm17egsPzs9wGJbbSiveLSzu22+/VadOnbRp0yalpKTINE1VqFBBiYmJ7oSqb9++kv73zXzGjBnKzs7Wtdde68/QbYnrYW27du3S888/r5EjR+qDDz7QwYMHFRUVpfHjx2vTpk1yuVxKTExUamqqNm3apLlz5+rVV1/V0aNHmTNcwsaMGaPY2Fj9+OOPkv7XzxwaGqqkpCT16tVLt956q2bMmKH169frhx9+0BtvvEH7jpeK3rO+/vprtWjRQh07dtSxY8fUuXNnzZo1S4cPH5bT6dSMGTM0depU9e/fX7/++qsWLFigw4cP+zt84KIw59bCtm7dqptuukn9+/fX6NGjVaVKFY/Hs7Ky9NZbb2natGmqUKGCYmNjlZGRoTVr1uiDDz5Q8+bN/RO4TXE9rO3rr79Wp06dVKtWLe3cuVN5eXl64IEHNG7cOG3ZskWDBw9WhQoVdPDgQb3zzju6/vrr3c9lJFvJevnll/Xkk0/K6XSqevXqWrx48RkV8iNHjmjBggVKTk5WQECAKlSooGuuuUbvvvuuJDG54iJ8/fXXatOmjYYNG6ZRo0apUqVKkqTffvtNjz/+uF5//XXNmTPHPd2lyH//+181aNBADRs29EfYl52iObfdn37PcnNu3x3btVTOuSW5tajc3Fz17dtXlSpV0ksvveTef/LkSWVlZeno0aPuW4Vu3LhRzz//vHJychQVFaUHH3yQN6VLjOthbd99952io6OVlJSkYcOGqUKFCho8eLAWLVqkdevWufs3582bpw8//FBt2rTxeD6JU8nZv3+/nnrqKdWoUUP9+vXT3XffrZycHC1btuysLSA//fSTTpw4oby8PF133XWS+OHjYmzdulXNmjVTUlKSJk6cKMnz33nRe9pHH32k77//XtWqVeP/g4tUlNz2+Pv7lktu3xnTpVQmt7xbWJDL5VJgYKC2bdumpk2buvevWrVKI0eOVOPGjdW0aVP17dtXu3btUnR0tN566y0tX75cM2bMIJG6xLge1paRkaE2bdqobdu2Gjt2rHtO7YQJExQQEOD+JH737t3VokUL96iv0z+gxDf0klOtWjV17txZ3bt3V1RUlN5//32VKVPG/cG+07lcLtWrV0/XXHONO7E1TZPE9gKdOnVKn3/+uQoKCtSyZUv3/tP/nQcGBmrIkCHKycnRxx9/fMbjwOWMdwyL+f777zV58mT9/vvvioqK0ocffqgff/xRycnJGjZsmA4fPqypU6fqtdde05IlS7Rs2TJ/h2xrXA/rCw8PV1xcnA4ePKh///vfysrKkiTt27dPOTk5ioyMlCS1bdtWAQEB+sc//iFJJEw+UHQ3uG7duqlp06Y6deqUqlatqg8++OCMBHf37t1avHix+4OYRUi4Lsz333+v5ORkJSQk6LHHHlOPHj30xhtveKwp+oHhxhtv1KlTp3TkyBE/RQuUDN7dLeTrr79W06ZNFRgYqCpVqujOO+/UL7/8ohtvvFHPPfecRo4cqcmTJ2vgwIHq0aOH2rdvr08//dRjJiffCC4droe17d69W7NmzdKOHTu0ZMkSNWjQQM8884w++eQTbd26Vd27d9ff/vY3xcbGyjRNBQYGauTIkdq+fbt+++03ZtmWoDfffFMjR45U27ZtNWDAAP3zn/+UVPghS5fLpapVq7oruN27d9fKlSvVpk0bffLJJ4zJ80LRe1ZAQIDKlSun5ORkjRgxQv369dObb77pXmcYhgoKCvTFF1+ocePGuvHGG/0YtX0YFrhpw+lbaZ6WwLuIRfzwww+KiYnRuHHj9Pjjj0uS+vbtq3bt2unQoUOqVauWwsPDJRX+1J2bmyuHw+G+1zcuLa6HtX377be66667dPXVV+vKK69UgwYN9Pbbb+uuu+7SyJEjdeTIEd15553uO8MV9RK2bt1aGzZsYDJCCRo1apTeeecd3XDDDWratKk+/fRTLVy4UGvWrNHbb7+twMBA5efnq1q1alq9erXatWunLl26qEuXLnrxxRf9Hf5l6/T3rKSkJPf+Z555RlLh+5dpmrr33nslFY4r/M9//qNq1aopIiLCLzEDJYXk1gK+++47xcbGKioqShMmTJBU2OwfHBysqKgoRUVFeaw3TVOTJ0/WV199peeee873Adsc18Patm3bpnbt2ulvf/ubHnroIXfbgSQtXrxYAwcO1KJFi/SXv/xFJ06cUJkyZdyP8028ZE2bNk2vv/66VqxYoebNmysgIED79u3TO++8o3Hjxqlnz55asmSJe3bqb7/9piNHjqhXr17697//LYkPj12Ms71nFbWESP9LcPv16ydJuvfeezV+/HjNnTtXn3zyyRmTX4DLHcmtnxWNamndurV+/PFHDRs2TDNmzFBwcLAKCgrcd+wpsmTJEq1atUpLly7VqlWrmAF5iXE9rO3kyZMaN26c7r33XiUnJ7v35+fn69dff1VQUJBeeeUVOZ1OTZ48WU6nU926dVO5ctb5BLMdmaapnJwcrVy5UklJSWrVqpVM05RpmqpZs6YGDBggh8OhMWPG6IUXXtBDDz2k48eP66mnnlL9+vVJbL1wrveswMBAj/esogR38ODBWrBggdavX6+PP/5Y11xzjT/DtxWHw5DDaZ3f3DlKcVsC7yJ+9OWXX+r666/X6NGj9eGHH2r8+PF64403NGzYMEln3sZ17dq1WrFihTIzM/Xxxx+rRYsW/grdlrge1hcQEKD09HQ1atTIvW/VqlUaPXq0rr32WkVHR6t79+7617/+pXbt2mnkyJF67733/Bhx6WAYhjIzM7Vx40aPH/CKWnTCwsJ09913q1GjRtq4caMkKSQkRI8++qhSUlIkkdhejAt9z3rmmWc0ZMgQpaSk6NNPP/WYpADYCZVbP8rJydGQIUM0fvx4SVLPnj0lSU888YSkwjtbFb05OZ1OdejQQQ0aNFDFihVL3cw6X+B6WF9OTo4OHTqkb775Rtu3b9eSJUu0YMECXXPNNZo0aZLKly+viRMnauLEiXr55ZeVl5enVq1a+TvsUiE0NFRBQUH66quv1LlzZ4/ec9M0FRkZqa5du2rp0qXKy8tTUFCQ+8YmjPu6OBf6niVJzz77rJKSklS5cmX/BA34AMmtH91000266aabJBW+uYeFhSkhIUHSmW9ORd8MatWq5bd47Y7rYX2hoaGaNWuW4uPjtXr1ah05ckTPPvusOnTooPr16ys/P1+LFi3S9u3bJUnz58/3b8CliGEYql27tt577z317NlT9erVkySPqRS///67YmJiFBQUpFOnTrknI/AhzItzIe9Zp/99k9iWjKIpBVZhpVh8jeTWIk6/z/rpb05Op1PTpk1TUFCQP8Mrdbge1nXzzTfr559/1sGDB1W7dm2PyQdOp1NhYWGqW7euO6kicfKN8uXLa+rUqerYsaMmTZqkcePGqW7duu6//4MHD2rt2rU6cOCANm3apPvuu08PPPCAxwf+cPH+7D2LEWsoTfjXbkFFb04Oh0ODBw9WcHCwx4dn4FtcD+upWbOmatas6bEvLy9PkyZN0vr16zV58mSSWj+IjY3VjBkzNHz4cO3du1d33nmnYmNjtW3bNk2aNEkREREaP368nE6nbrrpJhLbEsJ7Fko7kluLCg0NVY8ePRQYGKiYmBh/h1PqcT2s7fXXX9emTZu0aNEiffDBB0yt8KMhQ4aoXr16Gj58uMaMGaMTJ06oVatWiomJ0ezZs/0dXqnBe5bvORwOS/WOWykWXzNMbtNjaUXD52ENXA/r2b59u+6//35VqlRJkydPVuPGjf0dElTYX5uTk6ODBw/qiiuuUPXq1SXprCP1UHJ4zyp52dnZCgsLU+9paxRUxjpjB/NOHNfrI25RVlZWqfvQM5Vbi+NNyVq4HtbTsGFDLVq0SMHBwQoLC/N3OPh/lSpVUqVKlXTFFVe495mmSWLrY7xnoTQiuQVw2SuqCsLaSLRgZ4bDkGGhCQVWisXXSm9DBgAAAGyH5BYAAAC2QVsCAACAl7iJg3VQuQUAAIBtkNzaUG5uriZMmKDc3Fx/h1LqcS2sg2thHVwL6+BawI6Yc2tDRTP3SuNsO6vhWlgH18I6uBbWwbXwXtHfYeKstQoqU97f4bjlnTimVx/sUCqvLZVbAAAA2AbJLQAAAGyDaQmncblc2r9/vypUqHBZDxvPzs72+C/8h2thHVwL6+BaWIddroVpmjp69KgiIyPlcPinbse0BOsguT3N/v37VbNmTX+HccnY6Vwud1wL6+BaWAfXwjrsci327dunK6+80t9hwM9Ibk9ToUIFSYX/c5S25msrysz7yd8h4DQuk7cL4Awm3X1WcPToMTWrf4P7+zhKN75bnaaoFSE0NJTk1gJcedb51Ckklxno7xAAyzFJbi3Fny2FtCVYB/9XAgAAwDZIbgEAAGAbtCUAAAB4yWFYrC3hMp765C0qtwAAALANklsAAADYBm0JAAAA3nI4ZPjpBhJnZaVYfKz0njkAAABsh+QWAAAAtkFbAgAAgJccTkMOp3UmFFgpFl+jcgsAAADbILkFAACAbdCWAAAA4CWHw2I3cbBQLL5G5RYAAAC2QXILAAAA26AtAQAAwEu0JVgHlVsAAADYBsktAAAAbIO2BAAAAC/RlmAdVG4BAABgGyS3AAAAsA3aEgAAALzkMAw5DOu0AlgpFl+jcgsAAADbILkFAACAbdCWAAAA4CWmJVgHlVsAAADYBsktAAAAbIO2BAAAAC85HYacFmoFsFIsvkblFgAAALZBcgsAAADboC0BAADAS0xLsA4qtwAAALANklsAAADYBm0JAAAAXjJkyGFYpxXAkHVi8TUqtwAAALANklsAAADYBm0JAAAAXmJagnVQuQUAAIBtkNwCAADANmhLAAAA8BJtCdZB5RYAAAC2QXILAAAA26AtAQAAwEsOR+FmFVaKxddK8akDAADAbkhuAQAAYBvFTm4NwzjvNmHCBEnS0qVLdcMNNygsLEwVKlTQ1VdfreHDh7uPM2HCBDVv3vyM4+/evVuGYSgtLU2SlJKSIsMwVKlSJZ08edJj7aZNm9yvW+TkyZPq16+fmjZtqoCAAHXr1q3YfwkAAADecBqG5bbSqtjJ7YEDB9zb9OnTFRoa6rHv0Ucf1dq1a9WzZ091795dGzdu1ObNmzV58mTl5+dfdIAVKlTQ0qVLPfbNnTtXtWrV8thXUFCgMmXK6OGHH1ZcXNxFvx4AAAAuX8X+QFlERIT7z2FhYTIMw2OfJK1YsUI33nijRo0a5d531VVXeVVF7du3r+bNm6d77rlHknTixAm99dZbevjhhzVp0iT3unLlyumll16SJK1fv16ZmZl/euzc3Fzl5ua6v87Ozr7oOAEAAOB/l7TnNiIiQt9//72+++67S3bM++67T59++qn27t0rSXr33XcVFRWl6667zutjJycnKywszL3VrFnT62MCAIDSp3BagmGhzd9/I/5zSU/9oYce0vXXX6+mTZsqKipKCQkJmjdvnkd19EJVr15dnTt31vz58yVJ8+bNU//+/S9JvElJScrKynJv+/btuyTHBQAAgH9c0uS2XLlyeu+997Rz506NHTtW5cuX18iRI9W6dWvl5ORc9HH79++v+fPn6+eff1Zqaqp69ep1SeINDg5WaGioxwYAAFBazZo1S1FRUQoJCVF0dLQ2btx4zrXz588/Y8BASEiIxxrTNDVu3DjVqFFDZcqUUVxcnHbs2FGi51AiRet69epp4MCBeuWVV7Rlyxb98MMPWrRokSQpNDRUWVlZZzynqEc2LCzsjMc6d+6sEydOaMCAAbrttttUpUqVkggbAADgovi/DeHM7UItWrRII0aM0Pjx47VlyxZde+21io+P18GDB8/5nD8OGNizZ4/H41OnTtXMmTM1e/ZsbdiwQeXKlVN8fPwZk7AupRLvyIiKilLZsmV1/PhxSVLDhg31yy+/KCMjw2Pdli1bFBIScsYUBEkKCAhQnz59lJKScslaEgAAAPA/06ZN06BBg5SYmKgmTZpo9uzZKlu2rObNm3fO5xQNGCjawsPD3Y+Zpqnp06dr7NixuuOOO9SsWTO99tpr2r9/v5YtW1Zi53FJk9sJEyZo9OjRSklJ0a5du/TVV1+pf//+ys/P1y233CJJio+PV8OGDXXPPffo888/188//6zFixdr7NixGjZsmJxO51mPPWnSJB06dEjx8fHnfP0ffvhBaWlpOnLkiLKyspSWluaemwsAAFDaZGdne2zn+hxUXl6eNm/e7DFO1eFwKC4uTqmpqec8/rFjx1S7dm3VrFlTd9xxh77//nv3Y7t27VJ6errHMcPCwhQdHX3eY3qr2KPAiqNdu3aaNWuW+vTpo4yMDFWqVEktWrTQ6tWr1bBhw8IXDAjQ6tWrNWbMGN1zzz06dOiQ6tSpo2HDhmnEiBHnPHZQUJCqVq163tfv0qWLRzm8RYsWkgp/cgAAACgpDsOQw0I3TiiK5Y+ToMaPH+++8dbpfvvtNxUUFHhUXiUpPDxc27ZtO+trNGzYUPPmzVOzZs2UlZWl5557Tm3atNH333+vK6+8Uunp6e5j/PGYRY+VhItKbvv166d+/fqdsT82NlaxsbF/+vzIyEj39INzad++/XmT0m7dup3x+O7du//0tQEAAEqLffv2eXxgPjg4+JIdOyYmRjExMe6v27Rpo8aNG+tf//qXx70IfK0UT0EDAACwtz9OhTpXclu1alU5nc4zPhOVkZFxxk27ziUwMFAtWrTQzp07Jf3vBmDeHPNikNwCAAB4ybDAdITTN+MCpyUEBQWpZcuWWrt2rXufy+XS2rVrPaqz51NQUKBvv/1WNWrUkCTVqVNHERERHsfMzs7Whg0bin3Mi3FJe24BAABweRoxYoT69u2rVq1aqXXr1po+fbqOHz+uxMRESVKfPn10xRVXKDk5WZI0ceJE3XDDDapfv74yMzP17LPPas+ePRo4cKCkwkkKw4cP19NPP60GDRqoTp06evLJJxUZGalu3bqV2HmQ3AIAAEA9e/bUoUOHNG7cOKWnp6t58+ZauXKl+wNhe/fuleO0+/r+/vvvGjRokNLT01WpUiW1bNlSn3/+uZo0aeJeM3r0aB0/flyDBw9WZmam2rZtq5UrV55xs4dLyTAZJeCWnZ2tsLAwZWVlcbcyC8jMK9k7mODCuMxAf4cAWI5p0t1nBUezj6pO+DV++f5dlDs8++E3KlOugk9f+3xOHD+qUXHNSmVOw/+VAAAAsA2SWwAAANgGPbcAAABechqGnBc4oaAkOS10Qwlfo3ILAAAA2yC5BQAAgG3QlgAAAOAlh2HIYaFWACvF4mtUbgEAAGAbJLcAAACwDdoSAAAAvORwGHJYaFqClWLxNSq3AAAAsA2SWwAAANgGbQkAAABeoi3BOqjcAgAAwDZIbgEAAGAbtCUAAAB4ibYE66ByCwAAANsguQUAAIBt0JYAAADgJYdhyGFYpxXASrH4GpVbAAAA2AbJLQAAAGyDtgQAAAAvMS3BOqjcAgAAwDZIbgEAAGAbtCUAAAB4yeEw5LRQKwBtCQAAAIANkNwCAADANmhLAAAA8BI3cbAOKrcAAACwDZJbAAAA2AZtCQAAAF7iJg7WQeUWAAAAtkFyCwAAANugLQEAAMBLDkfhZhVWisXXSvGpAwAAwG5IbgEAAGAbtCUAAAB4ySGL3cRB1onF10huYVkVgxr4OwSc5pRrvb9DwP9zmYH+DgGwFDPwmL9DgIXQlgAAAADboHILAADgJYdRuFmFlWLxNSq3AAAAsA2SWwAAANgGbQkAAABechgWm5ZgoVh8jcotAAAAbIPkFgAAALZBWwIAAICXaEuwDiq3AAAAsA2SWwAAANgGbQkAAABe4iYO1kHlFgAAALZBcgsAAADboC0BAADASw5DclpoQgFtCQAAAIANkNwCAADANmhLAAAA8BLTEqyDyi0AAABsg+QWAAAAtkFbAgAAgJcchiGHpaYlWCcWX6NyCwAAANsguQUAAIBt0JYAAADgJdoSrIPKLQAAAGyD5BYAAAC2QVsCAACAl7iJg3VQuQUAAIBtkNwCAADANmhLAAAA8BLTEqyDyi0AAABsg+QWAAAAtkFbAgAAgJeYlmAdVG4BAABgGyS3AAAAsA3aEgAAALzEtATroHILAAAA2yC5BQAAgG3QlgAAAOAlw2JtCYaFYvE1KrcAAACwDZJbAAAA2AZtCQAAAF7iJg7WQeUWAAAAtkFyCwAAANugLQEAAMBL3MTBOqjcAgAAwDZIbgEAAGAbxU5uDcM47zZhwgRJ0tKlS3XDDTcoLCxMFSpU0NVXX63hw4e7jzNhwgQ1b978jOPv3r1bhmEoLS1NkpSSkiLDMFSpUiWdPHnSY+2mTZvcr1skJSVFd9xxh2rUqKFy5cqpefPmWrhwYfH/JgAAAC5S0bQEK22lVbGT2wMHDri36dOnKzQ01GPfo48+qrVr16pnz57q3r27Nm7cqM2bN2vy5MnKz8+/6AArVKigpUuXeuybO3euatWq5bHv888/V7NmzfTuu+/qm2++UWJiovr06aP//ve/F/3aAAAAuLwU+wNlERER7j+HhYXJMAyPfZK0YsUK3XjjjRo1apR731VXXaVu3bpddIB9+/bVvHnzdM8990iSTpw4obfeeksPP/ywJk2a5F43ZswYj+cNGzZMq1ev1pIlS3Trrbee9di5ubnKzc11f52dnX3RcQIAAMD/LmnPbUREhL7//nt99913l+yY9913nz799FPt3btXkvTuu+8qKipK11133Z8+NysrS5UrVz7n48nJyQoLC3NvNWvWvGRxAwCA0qNoWoKVtosxa9YsRUVFKSQkRNHR0dq4ceM5186ZM0d/+ctfVKlSJVWqVElxcXFnrO/Xr98ZraydOnW6qNiK65Imtw899JCuv/56NW3aVFFRUUpISNC8efM8qqMXqnr16urcubPmz58vSZo3b5769+//p897++23tWnTJiUmJp5zTVJSkrKystzbvn37LjpOAACAy9miRYs0YsQIjR8/Xlu2bNG1116r+Ph4HTx48KzrU1JSdM8992jdunVKTU1VzZo11bFjR/36668e6zp16uTRyvrmm2+W6Hlc0uS2XLlyeu+997Rz506NHTtW5cuX18iRI9W6dWvl5ORc9HH79++v+fPn6+eff1Zqaqp69ep13vXr1q1TYmKi5syZo6uvvvqc64KDgxUaGuqxAQAA2EV2drbHdr6C47Rp0zRo0CAlJiaqSZMmmj17tsqWLat58+addf3ChQv1wAMPqHnz5mrUqJFeeeUVuVwurV271mNdcHCwIiIi3FulSpUu6Tn+UYmMAqtXr54GDhyoV155RVu2bNEPP/ygRYsWSZJCQ0OVlZV1xnMyMzMlFfbz/lHnzp114sQJDRgwQLfddpuqVKlyztf++OOPddttt+n5559Xnz59Ls0JAQAAnIdTktOw0Pb/cdWsWdOjBTM5Ofms8efl5Wnz5s2Ki4tz73M4HIqLi1Nqamqx/g5ycnKUn59/RktoSkqKqlevroYNG2rIkCE6fPjwxfwVF1uJ36EsKipKZcuW1fHjxyVJDRs21C+//KKMjAyFh4e7123ZskUhISFnTEGQpICAAPXp00dTp07VBx98cM7XSklJ0a233qopU6Zo8ODBl/5kAAAALiP79u3z+M10cHDwWdf99ttvKigo8MjNJCk8PFzbtm0r1ms99thjioyM9EiQO3XqpL/+9a+qU6eOfvrpJ40ZM0adO3dWamqqnE7neY528S5pcjthwgTl5OSoS5cuql27tjIzMzVz5kzl5+frlltukSTFx8erYcOGuueee/T0008rIiJCW7Zs0dixYzVs2LBznuikSZM0atSoc1Zt161bp1tvvVXDhg1T9+7dlZ6eLkkKCgo674fKAAAA7MpXbZfPPPOM3nrrLaWkpCgkJMS9PyEhwf3npk2bqlmzZqpXr55SUlLUoUOHEonlkrYltGvXTj///LP69OmjRo0aqXPnzkpPT9fq1avVsGFDSYVV2NWrV6tWrVq65557dM0112j8+PEaNmyYx2ivPwoKClLVqlU9btxwugULFignJ0fJycmqUaOGe/vrX/96KU8RAADgDP6ejODttISqVavK6XQqIyPDY39GRsYZo1//6LnnntMzzzyj1atXq1mzZuddW7duXVWtWlU7d+68oPguhGGaplliR7/MZGdnKywsTFlZWXy4DPiDU671/g4B/89lBvo7BMBSsrOPqVrlDn75/l2UO2zcfUDlLZQ7HMvOVuuoGhf0dxIdHa3WrVvrhRdekCS5XC7VqlVLQ4cO1eOPP37W50ydOlWTJ0/WqlWrdMMNN/zpa/zyyy+qVauWli1bpttvv734J3QBSuQDZQAAALi8jBgxQnPmzNGCBQu0detWDRkyRMePH3ePVe3Tp4+SkpLc66dMmaInn3xS8+bNU1RUlNLT05Wenq5jx45Jko4dO6ZRo0bpiy++0O7du7V27Vrdcccdql+/vuLj40vsPEr8A2UAAAB2582NE0rCxcTSs2dPHTp0SOPGjVN6erqaN2+ulStXuj9ktnfvXjkc/6uLvvTSS8rLy9Ndd93lcZzx48drwoQJcjqd+uabb7RgwQJlZmYqMjJSHTt21KRJk875wbZLgbaE09CWAJwbbQnWQVsC4MkKbQlf7km3XFtCq9oRpTKnoS0BAAAAtkFbAgAAgJccRuFmFVaKxdeo3AIAAMA2SG4BAABgG7QlAAAAeMmw2LSEc930qjSgcgsAAADbILkFAACAbdCWAAAA4CWmJVgHlVsAAADYBsktAAAAbIO2BAAAAC85LDYtwUqx+BqVWwAAANgGyS0AAABsg7YEAAAALxlyyZDL32G4WSkWX6NyCwAAANsguQUAAIBt0JYAAADgJcNwyTCs0wpgpVh8jcotAAAAbIPkFgAAALZBWwIAAICXHCqQQwX+DsPNSrH4GpVbAAAA2AbJLQAAAGyDtgQAAAAvMS3BOqjcAgAAwDZIbgEAAGAbtCUAAAB4yZBLhqzTCmClWHyNyi0AAABsg+QWAAAAtkFbAgAAgJcMw7TUhALDMP0dgt9QuQUAAIBtkNwCAADANmhLAAAA8JJDBXKowN9huFkpFl+jcgsAAADbILkFAACAbdCWAAAA4CXDcFlsWoJ1YvE1KrcAAACwDZJbAAAA2AZtCQAAAF4y5JIh67QCWCkWX6NyCwAAANugcgugWAIcN/o7BPw/M/ttf4eAIs4gf0cASYEnc/wdAiyE5BYAAMBLTEuwDtoSAAAAYBsktwAAALAN2hIAAAC8ViBDBf4O4jRWisW3qNwCAADANkhuAQAAYBu0JQAAAHjJYbjksNCEAivF4mtUbgEAAGAbJLcAAACwDdoSAAAAvGTIJUPWaQWwUiy+RuUWAAAAtkFyCwAAANugLQEAAMBLhmHKsNCEAsMw/R2C31C5BQAAgG2Q3AIAAMA2aEsAAADwkqECGSrwdxhuVorF16jcAgAAwDZIbgEAAGAbtCUAAAB4yTBcFpuWYJ1YfI3KLQAAAGyD5BYAAAC2QVsCAACAlxxyySHrtAJYKRZfo3ILAAAA2yC5BQAAgG3QlgAAAOAlQy4ZFmoFsFIsvkblFgAAALZBcgsAAADboC0BAADAS9zEwTqo3AIAAMA2SG4BAABgG7QlAAAAeMlQgQwV+DsMNyvF4mtUbgEAAGAbJLcAAACwDdoSAAAAvOQwXHJYaEKBlWLxNSq3AAAAsA2SWwAAANgGbQkAAABeMmTKkHVaAQyZ/g7Bb6jcAgAAwDZIbgEAAGAbtCUAAAB4yTBcMiw0ocBKsfgalVsAAADYBsktAAAAbIO2BAAAAC8ZKpChAn+H4WalWHyNyi0AAABsg+QWAAAAtlHiye2hQ4c0ZMgQ1apVS8HBwYqIiFB8fLzWr18vSYqKipJhGDIMQ+XKldN1112nd955x+MYJ06cUOXKlVW1alXl5uae8RpFx/jiiy889g8fPlzt27cvsXMDAACQ/jctwUrbxZg1a5aioqIUEhKi6Ohobdy48bzr33nnHTVq1EghISFq2rSp3n//fY/HTdPUuHHjVKNGDZUpU0ZxcXHasWPHRcVWXCWe3Hbv3l1fffWVFixYoB9//FHLly9X+/btdfjwYfeaiRMn6sCBA/rqq690/fXXq2fPnvr888/dj7/77ru6+uqr1ahRIy1btuysrxMSEqLHHnuspE8HAADAlhYtWqQRI0Zo/Pjx2rJli6699lrFx8fr4MGDZ13/+eef65577tGAAQP01VdfqVu3burWrZu+++4795qpU6dq5syZmj17tjZs2KBy5copPj5eJ0+eLLHzKNHkNjMzU59++qmmTJmi2NhY1a5dW61bt1ZSUpJuv/1297oKFSooIiJCV111lWbNmqUyZcpoxYoV7sfnzp2r3r17q3fv3po7d+5ZX2vw4MH64osvzviJAQAAAH9u2rRpGjRokBITE9WkSRPNnj1bZcuW1bx58866fsaMGerUqZNGjRqlxo0ba9KkSbruuuv04osvSiqs2k6fPl1jx47VHXfcoWbNmum1117T/v37z1msvBRKNLktX768ypcvr2XLlp21neBsAgICFBgYqLy8PEnSTz/9pNTUVN199926++679emnn2rPnj1nPK9OnTq6//77lZSUJJereKX43NxcZWdne2wAAAB28cc851z5WF5enjZv3qy4uDj3PofDobi4OKWmpp71OampqR7rJSk+Pt69fteuXUpPT/dYExYWpujo6HMe81Io0eQ2ICBA8+fP14IFC1SxYkXdeOONGjNmjL755puzrs/Ly1NycrKysrJ08803S5LmzZunzp07q1KlSqpcubLi4+P16quvnvX5Y8eO1a5du7Rw4cJixZecnKywsDD3VrNmzYs7UQAAUKoZpstymyTVrFnTI9dJTk4+a/y//fabCgoKFB4e7rE/PDxc6enpZ31Oenr6edcX/fdCjnkp+KTndv/+/Vq+fLk6deqklJQUXXfddZo/f757zWOPPaby5curbNmymjJlip555hl17dpVBQUFWrBggXr37u1e27t3b82fP/+s1dlq1arp0Ucf1bhx49yV3/NJSkpSVlaWe9u3b98lOWcAAAAr2Ldvn0euk5SU5O+QSpxPRoGFhITolltu0ZNPPqnPP/9c/fr10/jx492Pjxo1Smlpafrll1/0+++/uz8YtmrVKv3666/q2bOnAgICFBAQoISEBO3Zs0dr164962uNGDFCJ06c0D//+c8/jSs4OFihoaEeGwAAgF38Mc8JDg4+67qqVavK6XQqIyPDY39GRoYiIiLO+pyIiIjzri/674Uc81Lwy5zbJk2a6Pjx4+6vq1atqvr16ysiIkKGYbj3z507VwkJCUpLS/PYEhISzvnBsvLly+vJJ5/U5MmTdfTo0RI/FwAAAJku620XICgoSC1btvQoHrpcLq1du1YxMTFnfU5MTMwZxcY1a9a419epU0cREREea7Kzs7Vhw4ZzHvNSKNHb7x4+fFg9evRQ//791axZM1WoUEFffvmlpk6dqjvuuOO8zz106JBWrFih5cuX65prrvF4rE+fPrrzzjt15MgRVa5c+YznDh48WM8//7zeeOMNRUdHX9JzAgAAsKMRI0aob9++atWqlVq3bq3p06fr+PHjSkxMlFSYf11xxRXuvt1hw4apXbt2+sc//qGuXbvqrbfe0pdffqmXX35ZkmQYhoYPH66nn35aDRo0UJ06dfTkk08qMjJS3bp1K7HzKNHktnz58oqOjtbzzz+vn376Sfn5+apZs6YGDRqkMWPGnPe5r732msqVK6cOHTqc8ViHDh1UpkwZvf7663r44YfPeDwwMFCTJk3Svffee8nOBQAAwM569uypQ4cOady4cUpPT1fz5s21cuVK9wfC9u7dK4fjf7/0b9Omjd544w2NHTtWY8aMUYMGDbRs2TKPouTo0aN1/PhxDR48WJmZmWrbtq1WrlypkJCQEjsPwzRNs8SOfpnJzs5WWFiYsrKy6L8FYFlm9tv+DgFFnEH+jgCSsrNzVDGyl1++f7tzhyNLFRpazqevfT7Z2ccVVvnOUpnT+KXnFgAAACgJJLcAAACwjRLtuQUAACgVTLNwsworxeJjVG4BAABgGyS3AAAAsA3aEgAAALx1ETdOKFFWisXHqNwCAADANkhuAQAAYBu0JQAAAHjLZUouC7UCuJiWAAAAAFz2SG4BAABgG7QlAAAAeItpCZZB5RYAAAC2QXILAAAA26AtAQAAwFu0JVgGlVsAAADYBsktAAAAbIO2BAAAAG/RlmAZVG4BAABgGyS3AAAAsA3aEgAAALzlchVuVmGlWHyMyi0AAABsg+QWAAAAtkFbAgAAgLdMs3CzCivF4mNUbgEAAGAbJLcAAACwDdoSAAAAvMVNHCyDyi0AAABsg+QWAAAAtkFbAgAAgLdoS7AMKrcAAACwDZJbAAAA2AZtCQAAAF4yTZdMC7UCWCkWX6NyCwAAANsguQUAAIBt0JYAAADgLZercLMKK8XiY1RuAQAAYBsktwAAALAN2hIAAAC8ZlrsxgmmvwPwGyq3AAAAsA2SWwAAANgGbQkAAADeMs3CzSqsFIuPkdwCwGXGCL3b3yHg/5l5H/g7BEhSYIG/I4CF0JYAAAAA26ByCwAA4C3TZa1pCVaKxceo3AIAAMA2SG4BAABgG7QlAAAAeIu2BMugcgsAAADbILkFAACAbdCWAAAA4C2Xq3CzCivF4mNUbgEAAGAbJLcAAACwDdoSAAAAvMW0BMugcgsAAADbILkFAACAbdCWAAAA4C3aEiyDyi0AAABsg+QWAAAAtkFbAgAAgLdMs3CzCivF4mNUbgEAAGAbJLcAAACwDdoSAAAAvOVyFW5WYaVYfIzKLQAAAGyD5BYAAAC2QVsCAACAt7iJg2VQuQUAAIBtkNwCAADANmhLAAAA8JZpWqsVgJs4AAAAAJc/klsAAADYBm0JAAAA3uImDpZB5RYAAAC2QXILAAAA26AtAQAAwFu0JVgGlVsAAADYBsktAAAAbIO2BAAAAG+ZFmtLsNINJXyMyi0AAABsg+QWAAAAtkFbAgAAgLdcZuFmFVaKxceo3AIAAMA2SG4BAABgG7QlAAAAeIubOFgGlVsAAADYBsktAAAALsiRI0fUq1cvhYaGqmLFihowYICOHTt23vUPPfSQGjZsqDJlyqhWrVp6+OGHlZWV5bHOMIwztrfeeuuCYqMtAQAAwFulrC2hV69eOnDggNasWaP8/HwlJiZq8ODBeuONN866fv/+/dq/f7+ee+45NWnSRHv27NH999+v/fv3a/HixR5rX331VXXq1Mn9dcWKFS8oNpJbAAAAFNvWrVu1cuVKbdq0Sa1atZIkvfDCC+rSpYuee+45RUZGnvGca665Ru+++67763r16mny5Mnq3bu3Tp06pYCA/6WkFStWVERExEXHR1sCAACATWVnZ3tsubm5Xh8zNTVVFStWdCe2khQXFyeHw6ENGzYU+zhZWVkKDQ31SGwl6cEHH1TVqlXVunVrzZs3T6Z5YTN7LZHc7tu3T/3791dkZKSCgoJUu3ZtDRs2TIcPH3avWbJkiTp27KgqVarIMAylpaWdcZyTJ0/qwQcfVJUqVVS+fHl1795dGRkZPjwTAABQKhW1JVhpk1SzZk2FhYW5t+TkZK9PNT09XdWrV/fYFxAQoMqVKys9Pb1Yx/jtt980adIkDR482GP/xIkT9fbbb2vNmjXq3r27HnjgAb3wwgsXFJ/fk9uff/5ZrVq10o4dO/Tmm29q586dmj17ttauXauYmBgdOXJEknT8+HG1bdtWU6ZMOeexHnnkEa1YsULvvPOOPv74Y+3fv19//etffXUqAAAAlrJv3z5lZWW5t6SkpHOuffzxx8/6ga7Tt23btnkdU3Z2trp27aomTZpowoQJHo89+eSTuvHGG9WiRQs99thjGj16tJ599tkLOr7fe24ffPBBBQUFafXq1SpTpowkqVatWmrRooXq1aunJ554Qi+99JLuu+8+SdLu3bvPepysrCzNnTtXb7zxhm6++WZJhQ3JjRs31hdffKEbbrjhjOfk5uZ6lOezs7Mv8dkBAAD4T2hoqEJDQ4u1duTIkerXr99519StW1cRERE6ePCgx/5Tp07pyJEjf9ore/ToUXXq1EkVKlTQ0qVLFRgYeN710dHRmjRpknJzcxUcHFys8/BrcnvkyBGtWrVKkydPdie2RSIiItSrVy8tWrRI//znP2UYxnmPtXnzZuXn5ysuLs69r1GjRqpVq5ZSU1PPmtwmJyfrqaeeujQnAwAASi8bTEuoVq2aqlWr9qfrYmJilJmZqc2bN6tly5aSpI8++kgul0vR0dHnfF52drbi4+MVHBys5cuXKyQk5E9fKy0tTZUqVSp2Yiv5uS1hx44dMk1TjRs3PuvjjRs31u+//65Dhw796bHS09MVFBR0xriI8PDwc/Z/JCUleZTq9+3bd8HnAAAAUJo0btxYnTp10qBBg7Rx40atX79eQ4cOVUJCgntSwq+//qpGjRpp48aNkgoT244dO+r48eOaO3eusrOzlZ6ervT0dBUUFEiSVqxYoVdeeUXfffeddu7cqZdeekl///vf9dBDD11QfH5vS5D0p5+CCwoKKpHXDQ4OvqCfBAAAACAtXLhQQ4cOVYcOHeRwONS9e3fNnDnT/Xh+fr62b9+unJwcSdKWLVvckxTq16/vcaxdu3YpKipKgYGBmjVrlh555BGZpqn69etr2rRpGjRo0AXF5tfktn79+jIMQ1u3btWdd955xuNbt25VtWrVijW8NyIiQnl5ecrMzPRYn5GR4dWsNAAAgD/lMi3WlnBh47MuVOXKlc95wwZJioqK8ihetm/f/k+LmZ06dfK4ecPF8mtbQpUqVXTLLbfon//8p06cOOHxWHp6uhYuXPinjc1FWrZsqcDAQK1du9a9b/v27dq7d69iYmIuZdgAAACwKL+PAnvxxReVm5ur+Ph4ffLJJ9q3b59WrlypW265RVdddZXGjRsnqfDDZ2lpafrhhx8kFSauaWlp7n7asLAwDRgwQCNGjNC6deu0efNmJSYmKiYm5qwfJgMAAID9+D25bdCggTZt2qS6devq7rvvVu3atdW5c2ddddVVWr9+vcqXLy9JWr58uVq0aKGuXbtKkhISEtSiRQvNnj3bfaznn39et956q7p3766bbrpJERERWrJkiV/OCwAAlCIu03pbKWWYF3pPMx8YP368pk2bpjVr1vi06pqdna2wsDD37eAAADgfM+8Df4cASdnZOapY7S6/fP8uyh0yUx5TaHnrfEg9+1iuKrafUipzGktMS/ijp556SlFRUfriiy/UunVrORx+LzADAADgMmDJ5FaSEhMT/R0CAABA8djgJg52QUkUAAAAtkFyCwAAANuwbFsCAADAZcO0WFuCaaFYfIzKLQAAAGyD5BYAAAC2QVsCAACAl0zTJdNCrQBWisXXqNwCAADANkhuAQAAYBu0JQAAAHiLmzhYBpVbAAAA2AbJLQAAAGyDtgQAAABv0ZZgGVRuAQAAYBsktwAAALAN2hIAAAC85TILN6uwUiw+RuUWAAAAtkFyCwAAANugLQEAAMBbTEuwDCq3AAAAsA2SWwAAANgGbQkAAADeoi3BMqjcAgAAwDZIbgEAAGAbtCUAAAB4y2VaqxWAmzgAAAAAlz+SWwAAANgGbQkAAADeYlqCZVC5BQAAgG2Q3AIAAMA2aEsAAADwFm0JlkHlFgAAALZBcgsAAADboC0BAADAW6bF2hJMC8XiY1RuAQAAYBsktwAAALAN2hIAAAC85TILN6uwUiw+RuUWAAAAtkFyCwAAANugLQEAAMBb3MTBMkhuAQC4SEZQZ3+HAElGULa/Q4CF0JYAAAAA26ByCwAA4C3aEiyDyi0AAABsg+QWAAAAtkFbAgAAgLdoS7AMKrcAAACwDZJbAAAA2AZtCQAAAF4yC0yZBaa/w3CzUiy+RuUWAAAAtkFyCwAAANugLQEAAMBbLrNwsworxeJjVG4BAABgGyS3AAAAsA3aEgAAALxVYBZuVmGlWHyMyi0AAABsg+QWAAAAtkFbAgAAgJdM05RpoQkFpmmdWHyNyi0AAABsg+QWAAAAtkFbAgAAgLcKZK0JBQX+DsB/qNwCAADANkhuAQAAYBu0JQAAAHirwFW4WYWVYvExKrcAAACwDZJbAAAA2AZtCQAAAF4yXRa7iYOFYvE1KrcAAACwDZJbAAAA2AZtCQAAAN4qMC12EwcLxeJjVG4BAABgGyS3AAAAsA3aEgAAALzlMgs3q7BSLD5G5RYAAAC2QXILAAAA26AtAQAAwEtmgSnTQhMKrBSLr1G5BQAAgG2Q3AIAAMA2aEsAAADwlumSXC5/R/E/poVi8TEqtwAAALANklsAAABckCNHjqhXr14KDQ1VxYoVNWDAAB07duy8z2nfvr0Mw/DY7r//fo81e/fuVdeuXVW2bFlVr15do0aN0qlTpy4oNtoSAAAAvFVgFm5WUcKx9OrVSwcOHNCaNWuUn5+vxMREDR48WG+88cZ5nzdo0CBNnDjR/XXZsmXdfy4oKFDXrl0VERGhzz//XAcOHFCfPn0UGBiov//978WOjeQWAAAAxbZ161atXLlSmzZtUqtWrSRJL7zwgrp06aLnnntOkZGR53xu2bJlFRERcdbHVq9erR9++EEffvihwsPD1bx5c02aNEmPPfaYJkyYoKCgoGLFR1sCAACATWVnZ3tsubm5Xh8zNTVVFStWdCe2khQXFyeHw6ENGzac97kLFy5U1apVdc011ygpKUk5OTkex23atKnCw8Pd++Lj45Wdna3vv/++2PFRuQUAAPCS6TJluqzTllAUS82aNT32jx8/XhMmTPDq2Onp6apevbrHvoCAAFWuXFnp6ennfN69996r2rVrKzIyUt98840ee+wxbd++XUuWLHEf9/TEVpL76/Md949IbgEAAGxq3759Cg0NdX8dHBx8zrWPP/64pkyZct7jbd269aJjGTx4sPvPTZs2VY0aNdShQwf99NNPqlev3kUf949KpC3hj5+E++M2YcIE7d6922NfhQoVdPXVV+vBBx/Ujh073Mfatm2bDMPQF1984fEaN9xwg0JCQnTy5En3vpMnTyokJERz586VJPXr10/dunUriVMEAACwvNDQUI/tfMntyJEjtXXr1vNudevWVUREhA4ePOjx3FOnTunIkSPn7Kc9m+joaEnSzp07JUkRERHKyMjwWFP09YUct0QqtwcOHHD/edGiRRo3bpy2b9/u3le+fHn99ttvkqQPP/xQV199tXJycvTtt99qxowZuvbaa7VixQp16NBBjRo1UkREhFJSUnTDDTdIko4ePaotW7YoPDxcX3zxhdq3by+psFcjNzdXN998c0mcFgAAwNnZYFpCtWrVVK1atT9dFxMTo8zMTG3evFktW7aUJH300UdyuVzuhLU40tLSJEk1atRwH3fy5Mk6ePCgu+1hzZo1Cg0NVZMmTYp93BKp3EZERLi3sLAwGYbhsa98+fLutVWqVFFERITq1q2rO+64Qx9++KGio6M1YMAAFRQUSJJiY2OVkpLifs5nn32mq666SrfddpvH/pSUFNWuXVt16tQpidMCAAAo9Ro3bqxOnTpp0KBB2rhxo9avX6+hQ4cqISHBPSnh119/VaNGjbRx40ZJ0k8//aRJkyZp8+bN2r17t5YvX64+ffropptuUrNmzSRJHTt2VJMmTXTffffp66+/1qpVqzR27Fg9+OCD5604/5HlpiU4HA4NGzZMe/bs0ebNmyUVJrefffaZe4jvunXr1L59e7Vr107r1q1zP3fdunWKjY0t9mvl5uae8SlCAAAAnN/ChQvVqFEjdejQQV26dFHbtm318ssvux/Pz8/X9u3b3dMQgoKC9OGHH6pjx45q1KiRRo4cqe7du2vFihXu5zidTv33v/+V0+lUTEyMevfurT59+njMxS0OS36grFGjRpKk3bt3q3Xr1oqNjdXx48e1adMmxcTEKCUlRaNGjVLbtm3Vt29fnTx5UqZpauPGjRo4cGCxXyc5OVlPPfVUSZ0GAAAoLWzQlnAhKleufN4bNkRFRck0/xdDzZo19fHHH//pcWvXrq3333/fq9gsV7mV5P7LMAxDklS/fn1deeWVSklJUXZ2tr766iu1a9dONWrUUK1atZSamurut72Qym1SUpKysrLc2759+0rkfAAAAOAblqzcFo2ZOL13tn379lq3bp2aNWumBg0auBuNi1oTTNNU/fr1z5jndj7BwcEX1MMBAAAAa7NccutyuTRz5kzVqVNHLVq0cO+PjY3Vww8/rCZNmrinI0jSTTfdpDlz5sg0zQuq2gIAAFwqVr2JQ2nk97aEw4cPKz09XT///LOWL1+uuLg4bdy4UXPnzpXT6XSvK+q7nTdvntq1a+fe365dO23YsEEbN24kuQUAACjl/F65jYuLkySVLVtWtWvXVmxsrF5++WXVr1/fY12dOnVUu3Zt7dmzxyO5rVWrliIjI7V7926Pii4AAABKH8M8/aNspVx2drbCwsKUlZXlcas6AABgXf78/l302hlJ8QoNCfTpa59P9sl8hSevKpU5jd/bEgAAAIBLheQWAAAAtuH3nlsAAIDLnWlabFpCKe46pXILAAAA2yC5BQAAgG3QlgAAAOCtArNwsworxeJjVG4BAABgGyS3AAAAsA3aEgAAALzlMgs3q7BSLD5G5RYAAAC2QXILAAAA26AtAQAAwEtmgWRaaEKBWeDvCPyHyi0AAABsg+QWAAAAtkFbAgAAgLeYlmAZVG4BAABgGyS3AAAAsA3aEgAAALxV4CrcrMJKsfgYlVsAAADYBsktAAAAbIO2BAAAAC+ZLlOmhSYUWCkWX6NyCwAAANsguQUAAIBt0JYAAADgLZcpFVioFYC2BAAAAODyR3ILAAAA26AtAQAAwEtMS7AOKrcAAACwDZJbAAAA2AZtCQAAAF4yC0yZFpqWYKVYfI3KLQAAAGyD5BYAAAC2QVsCAACAl5iWYB1UbgEAAGAbJLcAAACwDdoSAAAAvOQqMOWy0IQCK8Xia1RuAQAAYBsktwAAALAN2hIAAAC8xLQE66ByCwAAANsguQUAAIBt0JYAAADgJdPlkuly+TsMNyvF4mskt6cxzcL+lOzsbD9HAgAAiqvo+3bR93GUbiS3pzl69KgkqWbNmn6OBAAAXKijR48qLCzM32HAz0huTxMZGal9+/apQoUKMgzD3+FctOzsbNWsWVP79u1TaGiov8Mp1bgW1sG1sA6uhXXY5VqYpqmjR48qMjLSf0EUmDKtdOMEK8XiYyS3p3E4HLryyiv9HcYlExoaelm/WdkJ18I6uBbWwbWwDjtcCyq2KMK0BAAAANgGlVsAAAAvmabFbuJQij9cR+XWhoKDgzV+/HgFBwf7O5RSj2thHVwL6+BaWAfXAnZkmKU5tQcAAPBCdna2wsLC9HOP1qoQaJ1fiB/NP6W672xUVlbWZd9PfaGscxUAAAAuU2aBKdNhnXqhpSY3+BhtCQAAALANklsAAADYBm0JAAAAXjJdFpuWYKFYfI3KLQAAAGyD5BYAAAC2QVsCAACAl1wuUy4LtQJYKRZfo3ILAAAA2yC5BQAAgG3QlgAAAOAls0AWu4mDvyPwHyq3AAAAsA2SWwAAANgGbQkAAABe4iYO1kHlFgAAALZBcgsAAADboC0BAADAS7QlWAeVWwAAANgGyS0AAABsg7YEAAAAL5kFpsVu4mCdWHyNyi0AAABsg+QWAAAAtkFbAgAAgJdM0yXTZfg7DDfTdPk7BL+hcgsAAADbILkFAACAbdCWAAAA4CWzwJRpWGdCAdMSAAAAABsguQUAAIBt0JYAAADgJdNlynRZpxXASrH4GpVbAAAA2AbJLQAAAGyD5BYAAMBLLpdpua0kHTlyRL169VJoaKgqVqyoAQMG6NixY+dcv3v3bhmGcdbtnXfeca872+NvvfXWBcVGzy0AAAAuSK9evXTgwAGtWbNG+fn5SkxM1ODBg/XGG2+cdX3NmjV14MABj30vv/yynn32WXXu3Nlj/6uvvqpOnTq5v65YseIFxUZyCwAAYFPZ2dkeXwcHBys4ONirY27dulUrV67Upk2b1KpVK0nSCy+8oC5duui5555TZGTkGc9xOp2KiIjw2Ld06VLdfffdKl++vMf+ihUrnrH2QtCWAAAA4CWzwLTcJhVWTMPCwtxbcnKy1+eampqqihUruhNbSYqLi5PD4dCGDRuKdYzNmzcrLS1NAwYMOOOxBx98UFWrVlXr1q01b948meaFtVhQuQUAALCpffv2KTQ01P21t1VbSUpPT1f16tU99gUEBKhy5cpKT08v1jHmzp2rxo0bq02bNh77J06cqJtvvllly5bV6tWr9cADD+jYsWN6+OGHix0fyS0AAIBNhYaGeiS35/P4449rypQp512zdetWr2M6ceKE3njjDT355JNnPHb6vhYtWuj48eN69tlnSW4BAAB8yQ43cRg5cqT69et33jV169ZVRESEDh486LH/1KlTOnLkSLF6ZRcvXqycnBz16dPnT9dGR0dr0qRJys3NLXbVmeQWAAAAqlatmqpVq/an62JiYpSZmanNmzerZcuWkqSPPvpILpdL0dHRf/r8uXPn6vbbby/Wa6WlpalSpUoX1E5BcgsAAIBia9y4sTp16qRBgwZp9uzZys/P19ChQ5WQkOCelPDrr7+qQ4cOeu2119S6dWv3c3fu3KlPPvlE77///hnHXbFihTIyMnTDDTcoJCREa9as0d///nc9+uijFxQfyS0AAICXzAJTpmGhtoSCko1l4cKFGjp0qDp06CCHw6Hu3btr5syZ7sfz8/O1fft25eTkeDxv3rx5uvLKK9WxY8czjhkYGKhZs2bpkUcekWmaql+/vqZNm6ZBgwZdUGyGeaHzFQAAACCpcI5sWFiYtrRqpAoBTn+H43b0VIGu+3KbsrKyiv2BMrtgzi0AAABsg7YEAAAAb5nWmpagUvyLeSq3AAAAsA2SWwAAANgGbQkAAABeMl0Wm5ZgpRYJH6NyCwAAANsguQUAAIBt0JYAAADgJbPAlCnrtAKU9E0crIzKLQAAAGyD5BYAAAC2QVsCAACAl1wuUy4LTUtwMS0BAAAAuPyR3AIAAMA2aEsAAADwkssluQx/R/E/Lpe/I/AfKrcAAACwDZJbAAAA2AZtCQAAAF6iLcE6qNwCAADANkhuAQAAYBu0JQAAAHiJtgTroHILAAAA2yC5BQAAgG3QlgAAAOAll1m4WYWVYvE1KrcAAACwDZJbAAAA2AZtCQAAAF5iWoJ1ULkFAACAbZDcAgAAwDZoSwAAAPCSy7RWKwDTEgAAAAAbILkFAACAbdCWAAAA4CXTJVmoK0GmlYLxMSq3AAAAsA2SWwAAANgGbQkAAABeclmsLcFKkxt8jcotAAAAbIPkFgAAALZBWwIAAICXaEuwDiq3AAAAsA2SWwAAANgGbQkAAABeoi3BOqjcAgAAwDZIbgEAAGAbtCUAAAB4ibYE66ByCwAAANsguQUAAIBt0JYAAADgJdoSrIPKLQAAAGyD5BYAAAC2QVsCAACAl0zTlGma/g7DzUqx+BqVWwAAANgGyS0AAABsg7YEAAAALzEtwTqo3AIAAMA2SG4BAABgG7QlAAAAeIm2BOugcgsAAADbILkFAACAbdCWAAAA4CWXaa1WAFfpvYcDlVsAAADYB8ktAAAAbIO2BAAAAC+ZLsll+DuK/zFpSwAAAAAufyS3AAAAsA3aEgAAALzkslhbAtMSAAAAABsguQUAAIBt0JYAAADgJdoSrIPKLQAAAGyD5BYAAAC2QVsCAACAl2hLsA4qtwAAALANklsAAADYBm0JAAAAXnKZksvfQZyGtgQAAADABkhuAQAAYBu0JQAAAHiJaQnWQeUWAAAAtkFyCwAAANugLQEAAMBLtCVYB5VbAAAA2AbJLQAAAGyDtgQAAAAv0ZZgHVRuAQAAYBsktwAAALANklsAAAAvucz/b02wylbCbQmTJ09WmzZtVLZsWVWsWLFYzzFNU+PGjVONGjVUpkwZxcXFaceOHR5rjhw5ol69eik0NFQVK1bUgAEDdOzYsQuKjeQWAAAAFyQvL089evTQkCFDiv2cqVOnaubMmZo9e7Y2bNigcuXKKT4+XidPnnSv6dWrl77//nutWbNG//3vf/XJJ59o8ODBFxSbYZpmKW45BgAAuHjZ2dkKCwvTHNVTWTn9HY5bjgo0SD8pKytLoaGhJfY68+fP1/Dhw5WZmXnedaZpKjIyUiNHjtSjjz4qScrKylJ4eLjmz5+vhIQEbd26VU2aNNGmTZvUqlUrSdLKlSvVpUsX/fLLL4qMjCxWTExLAAAA8NIJufwdgoeieLKzsz32BwcHKzg42Ofx7Nq1S+np6YqLi3PvCwsLU3R0tFJTU5WQkKDU1FRVrFjRndhKUlxcnBwOhzZs2KA777yzWK9FcgsAAHCRgoKCFBERoYfTd/k7lDOUL19eNWvW9Ng3fvx4TZgwweexpKenS5LCw8M99oeHh7sfS09PV/Xq1T0eDwgIUOXKld1rioPkFgAA4CKFhIRo165dysvL83coZzBNU4bhOXz3fFXbxx9/XFOmTDnvMbdu3apGjRpdkvhKCsktAACAF0JCQhQSEuLvMLw2cuRI9evX77xr6tate1HHjoiIkCRlZGSoRo0a7v0ZGRlq3ry5e83Bgwc9nnfq1CkdOXLE/fziILkFAACAqlWrpmrVqpXIsevUqaOIiAitXbvWncxmZ2drw4YN7okLMTExyszM1ObNm9WyZUtJ0kcffSSXy6Xo6OhivxajwAAAAHBB9u7dq7S0NO3du1cFBQVKS0tTWlqax0zaRo0aaenSpZIkwzA0fPhwPf3001q+fLm+/fZb9enTR5GRkerWrZskqXHjxurUqZMGDRqkjRs3av369Ro6dKgSEhKKPSlBonILAACACzRu3DgtWLDA/XWLFi0kSevWrVP79u0lSdu3b1dWVpZ7zejRo3X8+HENHjxYmZmZatu2rVauXOnR0rFw4UINHTpUHTp0kMPhUPfu3TVz5swLio05twAAALAN2hIAAABgGyS3AAAAsA2SWwAAANgGyS0AAABsg+QWAAAAtkFyCwAAANsguQUAAIBtkNwCAADANkhuAQAAYBsktwAAALANklsAAADYxv8BqecWDPHaJJYAAAAASUVORK5CYII=",
      "text/plain": [
       "<Figure size 800x1000 with 2 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "# Visualise second-order sensitivities\n",
    "fig, axes = plt.subplots(figsize=(8,10))\n",
    "aximage = axes.matshow(second_order, vmin=-1.0, vmax=1.0, cmap=\"RdYlBu\")\n",
    "cbar = fig.colorbar(aximage)\n",
    "aximage.axes.set_xticks(range(len(problem[\"names\"])))\n",
    "aximage.axes.set_xticklabels(problem[\"names\"], rotation=45)\n",
    "r = aximage.axes.set_yticklabels([\"\"] + problem[\"names\"])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Explanation\n",
    "The background of the sensitivity of sugarbeet yield to TSUM1 is related to the fact that a larger TSUM1 means that the crop will remain in its vegetative stage for a longer period. Given that there is a fixed harvest date this means that the time available for accumulating sugar in the beet root is shorter and thus simulated yields are lower.  Apparently, the other parameters that are part of the analysis do not have large impact within the selected ranges and thus TSUM1 dominates in terms of sensitivity over the other parameters."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAA1AAAAHWCAYAAAB5bWjdAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjguMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/SrBM8AAAACXBIWXMAAA9hAAAPYQGoP6dpAACAJUlEQVR4nOzdeXhM1x8G8Hcmy2Sd7Psuloh9F0sihCCo0qJUKarUlihVbSndKFW6KD9tlba2alE7sUQsQSyxZZFESEIWWzLZt7m/P9LcmiaYkJgs7+d55qm558yd7x036bzuPedIBEEQQERERERERE8l1XQBREREREREtQUDFBERERERkZoYoIiIiIiIiNTEAEVERERERKQmBigiIiIiIiI1MUARERERERGpiQGKiIiIiIhITQxQREREREREamKAIiIiIiIiUhMDFBFRHfHbb7/Bw8MDOjo6MDU11XQ5VIcsWLAAEokE9+7d03QpVaZHjx7o0aPHM73W1dUVY8eOfWq/devWQSKR4ObNm8/0PkRUMzFAEVGNV/Yl5Ny5c5ou5bnt3bsXCxYsqPL9RkdHY+zYsXB3d8ePP/6INWvWlOtz8+ZNSCQStR5lX/hu3ryJN998E+7u7tDT04OtrS28vb3x8ccfq+zb1dUVAwYMqLC2c+fOQSKRYN26deK2si/kUqkUSUlJ5V6jUCigr68PiUSCqVOnqrStWrUKr776KpydnSGRSNT6Iku116lTp7BgwQJkZGRouhQiIgCAtqYLICKqT/bu3YuVK1dWeYgKCQmBUqnEN998g4YNG1bYx8rKCr/99pvKtmXLliE5ORnLly8v1zcuLg4dOnSAvr4+xo0bB1dXV6SkpODChQv48ssvsXDhwueuWyaTYdOmTXjvvfdUtm/btu2xr/nyyy+RlZWFjh07IiUl5blroJrt1KlTWLhwIcaOHVulV1YPHjxYZfsiovqFAYqIqA5IT08HgCd+wTQ0NMTrr7+usm3z5s14+PBhue0AsHz5cmRnZyMiIgIuLi4Vvt/z6t+/f4UBauPGjQgICMBff/1V7jXHjh0Trz4ZGRlVSR0E5OTkwNDQUNNlvDC6urqaLoGIainewkdEtdLYsWNhZGSExMREDBgwAEZGRnBwcMDKlSsBAFeuXEHPnj1haGgIFxcXbNy4UeX1ZbcFhoaG4u2334aFhQXkcjneeOMNPHz4UKXv33//jYCAANjb20Mmk8Hd3R2ffvopSkpKytV15swZ9O/fH2ZmZjA0NETLli3xzTffiDWX1ffo7XJP88MPP6BZs2aQyWSwt7fHlClTVG5ncnV1FW+ps7KygkQiqZIrXPHx8XB0dCwXngDA2tr6ufcPACNHjkRERASio6PFbampqThy5AhGjhxZ4WtcXFzU+twe58aNG3j11Vdhbm4OAwMDdO7cGXv27FHpExISAolEgj/++AOff/45HB0doaenh169eiEuLu6p75GVlYXAwEC4urpCJpPB2toavXv3xoULF8Q+jxtHU9HYnFu3bmHQoEEwNDSEtbU1goKCcODAAUgkEoSEhIj9jh8/Lt7eKJPJ4OTkhKCgIOTl5ansr+znJz4+Hv3794exsTFGjRr11OO6d+8ehg0bBrlcDgsLC8yYMQP5+fnl+v3+++9o164d9PX1YW5ujhEjRlR4q+aZM2fQt29fmJiYwMDAAD4+Pjh58qTYvmDBAsyePRsA4ObmVu4W0//6+OOPoaOjg7t375ZrmzhxIkxNTcV6K/qcCwoK8PHHH6Nhw4bi5/fee++hoKDgqZ/NtWvX0LNnT+jr68PR0RGfffYZlErlU19HRLUPr0ARUa1VUlKCfv36wdvbG0uWLMGGDRswdepUGBoa4sMPP8SoUaMwZMgQrF69Gm+88Qa8vLzg5uamso+pU6fC1NQUCxYsQExMDFatWoVbt26JX6CB0rBlZGSEmTNnwsjICEeOHMH8+fOhUCiwdOlScV/BwcEYMGAA7OzsMGPGDNja2iIqKgq7d+/GjBkz8Pbbb+POnTsIDg4udyvd4yxYsAALFy6En58fJk+eLNYYHh6OkydPQkdHBytWrMCvv/6K7du3Y9WqVTAyMkLLli2f+/N1cXHBoUOHcOTIEfTs2fO591cRb29vODo6YuPGjfjkk08AAFu2bIGRkRECAgKq/P3S0tLQpUsX5ObmYvr06bCwsMD69esxaNAg/Pnnn3j55ZdV+i9evBhSqRSzZs1CZmYmlixZglGjRuHMmTNPfJ9Jkybhzz//xNSpU+Hp6Yn79+/jxIkTiIqKQtu2bStVc05ODnr27ImUlBTxvNq4cSOOHj1aru/WrVuRm5uLyZMnw8LCAmfPnsV3332H5ORkbN26VaVvcXEx/P390a1bN3z11VcwMDB4ai3Dhg2Dq6srFi1ahNOnT+Pbb7/Fw4cP8euvv4p9Pv/8c8ybNw/Dhg3DhAkTcPfuXXz33Xfw9vbGxYsXxaukR44cQb9+/dCuXTt8/PHHkEql+OWXX9CzZ08cP34cHTt2xJAhQ3D9+nVs2rQJy5cvh6WlJYDSfyioyOjRo/HJJ59gy5YtKmPnCgsL8eeff2Lo0KHQ09Or8LVKpRKDBg3CiRMnMHHiRDRt2hRXrlzB8uXLcf36dezYseOxn0tqaip8fX1RXFyM999/H4aGhlizZg309fWf+pkSUS0kEBHVcL/88osAQAgPDxe3jRkzRgAgfPHFF+K2hw8fCvr6+oJEIhE2b94sbo+OjhYACB9//HG5fbZr104oLCwUty9ZskQAIPz999/ittzc3HI1vf3224KBgYGQn58vCIIgFBcXC25uboKLi4vw8OFDlb5KpVL885QpUwR1f/Wmp6cLurq6Qp8+fYSSkhJx+/fffy8AENauXStu+/jjjwUAwt27d9Xad5mAgADBxcWlwrarV68K+vr6AgChdevWwowZM4QdO3YIOTk55fq6uLgIAQEBFe4nPDxcACD88ssvFdY7a9YsoWHDhmJbhw4dhDfffFMQBEEAIEyZMuWx9RsaGgpjxox5+oH+IzAwUAAgHD9+XNyWlZUluLm5Ca6uruLnfPToUQGA0LRpU6GgoEDs+8033wgAhCtXrjzxfUxMTJ5YtyCUfmYV1e7j4yP4+PiIz5ctWyYAEHbs2CFuy8vLEzw8PAQAwtGjR8XtFZ2rixYtEiQSiXDr1i1xW9nPz/vvv//EGsuU/X0NGjRIZfs777wjABAuXbokCIIg3Lx5U9DS0hI+//xzlX5XrlwRtLW1xe1KpVJo1KiR4O/vr/LzkZubK7i5uQm9e/cWty1dulQAICQkJKhVq5eXl9CpUyeVbdu2bSv3Wf33c/7tt98EqVSqcm4IgiCsXr1aACCcPHlS3Pbfv7uy8+rMmTPitvT0dMHExKRStRNR7cBb+IioVpswYYL4Z1NTUzRp0gSGhoYYNmyYuL1JkyYwNTXFjRs3yr1+4sSJ0NHREZ9PnjwZ2tra2Lt3r7jt0X9FzsrKwr1799C9e3fk5uaKt55dvHgRCQkJCAwMLDcO6VlvNzt06BAKCwsRGBgIqfTfX9dvvfUW5HJ5udvOqlqzZs0QERGB119/HTdv3sQ333yDwYMHw8bGBj/++GOVvc/IkSMRFxeH8PBw8b+Pu33vee3duxcdO3ZEt27dxG1GRkaYOHEibt68icjISJX+b775pspYme7duwNAhefSo0xNTXHmzBncuXPnuWvev38/HBwcMGjQIHGbnp4e3nrrrXJ9Hz1Xc3JycO/ePXTp0gWCIODixYvl+k+ePLlStUyZMkXl+bRp0wBA/HnZtm0blEolhg0bhnv37okPW1tbNGrUSLxqFhERgdjYWIwcORL3798X++Xk5KBXr14IDQ195tvf3njjDZw5cwbx8fHitg0bNsDJyQk+Pj6Pfd3WrVvRtGlTeHh4qNRedvW1oit+Zfbu3YvOnTujY8eO4jYrKyu1boskotqHt/ARUa2lp6dX7lYeExMTODo6lgstJiYm5cY2AUCjRo1UnhsZGcHOzk5ljMW1a9fw0Ucf4ciRI1AoFCr9MzMzAUD8sta8efNnPp7/unXrFoDSAPgoXV1dNGjQQGyvTo0bN8Zvv/2GkpISREZGYvfu3ViyZAkmTpwINzc3+Pn5qb2vxwXJNm3awMPDAxs3boSpqSlsbW2r7ZbBW7duoVOnTuW2N23aVGx/9O/Q2dlZpZ+ZmRkAVHguPWrJkiUYM2YMnJyc0K5dO/Tv3x9vvPEGGjRo8Ew1u7u7l/v8KpptMTExEfPnz8fOnTvL1Vh2rpbR1taGo6NjpWr578+Lu7s7pFKp+PMSGxsLQRDK9StT9o8VsbGxAIAxY8Y89r0yMzPFz7syhg8fjsDAQGzYsAHz589HZmYmdu/ejaCgoCf+Y0ZsbCyioqIee3vgkyZOedx59d+fXSKqGxigiKjW0tLSqtR2QRAq/R4ZGRnw8fGBXC7HJ598Iq6HdOHCBcyZM6feDBLX0tJCixYt0KJFC3h5ecHX1xcbNmwQA5Senl65iQrK5Obmin0eZ+TIkVi1ahWMjY0xfPhwlStumvSs59KwYcPQvXt3bN++HQcPHsTSpUvx5ZdfYtu2bejXrx+AxwfKkpKSx77vk5SUlKB379548OAB5syZAw8PDxgaGuL27dsYO3ZsuXNVJpM99+f832NQKpWQSCTYt29fhcdQNmtiWS1Lly5F69atK9z3s86waGZmhgEDBogB6s8//0RBQUGFM03+t/YWLVrg66+/rrDdycnpmeohorqHAYqI6rXY2Fj4+vqKz7Ozs5GSkoL+/fsDKJ2N7f79+9i2bRu8vb3FfgkJCSr7cXd3BwBcvXr1iVdlKnM7X9nsdzExMSpXLgoLC5GQkFCpqz9VqX379gCgsgaTi4tLudvfysTExIh9HmfkyJGYP38+UlJS1J5g41m4uLiI9Tyq7FbMJ9VYWXZ2dnjnnXfwzjvvID09HW3btsXnn38uBigzM7MKF4e9deuWyt932WcrCILK+fPf2QCvXLmC69evY/369XjjjTfE7cHBwVV2TLGxsSoTscTFxUGpVMLV1RVA6c+BIAhwc3ND48aNH7ufsp8XuVz+1PP4WW6BfeONN/DSSy8hPDwcGzZsQJs2bdCsWbMnvsbd3R2XLl1Cr169Kv2eLi4u4lW1R1V0rhFR7Vcz/omPiEhD1qxZg6KiIvH5qlWrUFxcLH7JLftX9EevOBQWFuKHH35Q2U/btm3h5uaGFStWlPtS/Ohry9bZqeiL83/5+flBV1cX3377rco+fv75Z2RmZlbLLHWPOn78uMpnU6ZsvMujtyf1798fycnJ5WYqKygowE8//QRra+snzj7n7u6OFStWYNGiRSrjSKpa//79cfbsWYSFhYnbcnJysGbNGri6usLT0/O536OkpKTc7XLW1tawt7dXmQ7b3d0dp0+fRmFhobht9+7d5ab79vf3x+3bt7Fz505xW35+frlxaBWdq4IgiNPoV4WyafjLfPfddwAg/rwMGTIEWlpaWLhwYbmrdIIg4P79+wCAdu3awd3dHV999RWys7PLvc+j05BX5memTL9+/WBpaYkvv/wSx44de+rVJ6D0quHt27crHN+Xl5eHnJycx762f//+OH36NM6ePatyDBs2bFC7ZiKqPXgFiojqtcLCQvTq1QvDhg1DTEwMfvjhB3Tr1k0csN+lSxeYmZlhzJgxmD59OiQSCX777bdyXw6lUilWrVqFgQMHonXr1njzzTdhZ2eH6OhoXLt2DQcOHABQ+sURAKZPnw5/f39oaWlhxIgRFdZmZWWFuXPnYuHChejbty8GDRok1tihQwe1vhQ+jy+//BLnz5/HkCFDxGnRL1y4gF9//RXm5uYIDAwU+06cOBFr167Fq6++inHjxqFNmza4f/8+tmzZgqtXr+LXX3996sKlM2bMUKuuXbt24dKlSwCAoqIiXL58GZ999hkAYNCgQU+cwv3999/Hpk2b0K9fP0yfPh3m5uZYv349EhIS8Ndff1XJrYNZWVlwdHTEK6+8glatWsHIyAiHDh1CeHg4li1bJvabMGEC/vzzT/Tt2xfDhg1DfHw8fv/9d/HqTJm3334b33//PV577TXMmDEDdnZ22LBhg3hLZNnVEg8PD7i7u2PWrFm4ffs25HI5/vrrr6eO16qMhIQEDBo0CH379kVYWBh+//13jBw5Eq1atQJQGgo/++wzzJ07Fzdv3sTgwYNhbGyMhIQEbN++HRMnTsSsWbMglUrx008/oV+/fmjWrBnefPNNODg44Pbt2zh69Cjkcjl27doF4N+fmQ8//BAjRoyAjo4OBg4c+MRFf3V0dDBixAh8//330NLSwmuvvfbUYxs9ejT++OMPTJo0CUePHkXXrl1RUlKC6Oho/PHHHzhw4IB49fW/3nvvPfz222/o27cvZsyYIU5j7uLigsuXL1f2Yyaimk4TU/8REVXG46YxNzQ0LNfXx8dHaNasWbnt/51mu2yfx44dEyZOnCiYmZkJRkZGwqhRo4T79++rvPbkyZNC586dBX19fcHe3l547733hAMHDpSbFlkQBOHEiRNC7969BWNjY8HQ0FBo2bKl8N1334ntxcXFwrRp0wQrKytBIpGoNaX5999/L3h4eAg6OjqCjY2NMHny5HJTpVfHNOYnT54UpkyZIjRv3lwwMTERdHR0BGdnZ2Hs2LFCfHx8uf4PHz4UgoKCBDc3N0FHR0eQy+WCr6+vsG/fvnJ91a0XFUxjXjYFd0WPR6dKf5z4+HjhlVdeEUxNTQU9PT2hY8eOwu7du1X6lE1jvnXrVpXtCQkJT32fgoICYfbs2UKrVq3E86BVq1bCDz/8UK7vsmXLBAcHB0Emkwldu3YVzp07V256bUEQhBs3bggBAQGCvr6+YGVlJbz77rvCX3/9JQAQTp8+LfaLjIwU/Pz8BCMjI8HS0lJ46623hEuXLpWr+XE/P49T9vcVGRkpvPLKK4KxsbFgZmYmTJ06VcjLyyvX/6+//hK6desmGBoaCoaGhoKHh4cwZcoUISYmRqXfxYsXhSFDhggWFhaCTCYTXFxchGHDhgmHDx9W6ffpp58KDg4OglQqVXta8LNnzwoAhD59+lTYXtHnXFhYKHz55ZdCs2bNBJlMJpiZmQnt2rUTFi5cKGRmZor9KpqC/vLly4KPj4+gp6cnODg4CJ9++qnw888/cxpzojpIIgjPMKqaiKiWW7duHd58802Eh4c/9l+ViWqyFStWICgoCMnJyXBwcNB0OTXOpUuX0Lp1a/z6668YPXq0psshojqEY6CIiIhquP/OcJifn4///e9/aNSoEcPTY/z4448wMjLCkCFDNF0KEdUxHANFRERUww0ZMgTOzs5o3bo1MjMz8fvvvyM6OpqTFFRg165diIyMxJo1azB16tQnjpUiInoWDFBEREQ1nL+/P3766Sds2LABJSUl8PT0xObNmzF8+HBNl1bjTJs2DWlpaejfvz8WLlyo6XKIqA7iGCgiIiIiIiI1cQwUERERERGRmhigiIiIiIiI1FSvx0AplUrcuXMHxsbG4kKERERERERU/wiCgKysLNjb2z9xYfV6HaDu3LkDJycnTZdBREREREQ1RFJSEhwdHR/bXq8DlLGxMYDSD0kul2u4GiIiIiIi0hSFQgEnJycxIzxOvQ5QZbftyeVyBigiIiIiInrq0B5OIkFERERERKQmBigiIiIiIiI1MUARERERERGpqV6PgSIiIiIiqgkEQUBxcTFKSko0XUqdpaWlBW1t7edevogBioiIiIhIgwoLC5GSkoLc3FxNl1LnGRgYwM7ODrq6us+8DwYoIiIiIiINUSqVSEhIgJaWFuzt7aGrq/vcV0ioPEEQUFhYiLt37yIhIQGNGjV64mK5T8IARURERESkIYWFhVAqlXBycoKBgYGmy6nT9PX1oaOjg1u3bqGwsBB6enrPtB9OIkFEREREpGHPejWEKqcqPmf+TREREREREamJAYqIiIiIiEhNDFBERERERERqYoAiIiIiIiK1SSSSJz4WLFgAOzs7LF68WOV177//PiQSCUJCQlS29+jRA6NHjwYAlJSUYPHixfDw8IC+vj7Mzc3RqVMn/PTTTyqvSUpKwrhx48SZC11cXDBjxgzcv3+/Wo8dYICqMUqUgqZLICIiIiJ6qpSUFPGxYsUKyOVylW2zZs1Cjx49ygWlo0ePwsnJSWV7fn4+Tp8+jZ49ewIAFi5ciOXLl+PTTz9FZGQkjh49iokTJyIjI0N8zY0bN9C+fXvExsZi06ZNiIuLw+rVq3H48GF4eXnhwYMH1Xr8nMa8BihRCvBechQetsbw87RBr6bWsDZ+tmkViYiIiKj2EgQBeUUlGnlvfR0ttdagsrW1Ff9sYmICiUSisg0AfH198e6776K4uBja2trIysrCxYsXsXz5cmzdulXsFxYWhoKCAvj6+gIAdu7ciXfeeQevvvqq2KdVq1Yq+54yZQp0dXVx8OBB6OvrAwCcnZ3Rpk0buLu748MPP8SqVasq/wGoiQGqBricnIHbGXm4nZGHw9HpAIDWTqbo7WmD3p42aGRtxAXViIiIiOqBvKISeM4/oJH3jvzEHwa6VRMPfH19kZ2djfDwcHh5eeH48eNo3Lgxhg4dilmzZiE/Px96eno4evQoXF1d4erqCqA0nB05cgTvvPMOrKysyu33wYMHOHDgAD7//HMxPJWxtbXFqFGjsGXLFvzwww/V9v2Zt/DVAK2dTLE/sDtm9WmMVk6mAICIpAwsPRCDPstD4bM0BJ/sisSp+HsoLlFqtlgiIiIioqdo1KgRHBwcxNv1QkJC4OPjA1tbWzg7OyMsLEzcXnb1CQC+/vpr3L17F7a2tmjZsiUmTZqEffv2ie2xsbEQBAFNmzat8H2bNm2Khw8f4u7du9V2bLwCVQNIJBJ42MrhYSvH1J6NkKbIx+GodByKSsOJuHtIfJCLtScTsPZkAkz0deDbxAp+njbwaWwFYz0dTZdPRERERFVEX0cLkZ/4a+y9q1LZOKi5c+ciJCQEs2fPBgD4+PggJCQEnTt3xpkzZ/DWW2+Jr/H09MTVq1dx/vx5nDx5EqGhoRg4cCDGjh2rMpGEIGhu/gAGqBrIRq6HkZ2cMbKTM3ILixF6/R4ORaXhSHQ6HuQUYkfEHeyIuAMdLQk6N7BAH08b9GpqA3tT/afvnIiIiIhqLIlEUmW30Wmar6+vODPexYsX4ePjA6A0QP3vf/+Dt7c3CgsLxQkkykilUnTo0AEdOnRAYGAgfv/9d4wePRoffvghGjZsCIlEgqioKLz88svl3jMqKgpmZmYV3v5XVerG304dZqCrjb7NbdG3uS1KlAIuJD7Eocg0BEem4ca9HByPvYfjsfcw7+9raGYvh1/T0nFTzezlHDdFRERERBrj6+uLnJwcfP3112jUqBGsra0BAN7e3hg/fjz27dsn3ur3JJ6engCAnJwcuLm5oXfv3vjhhx8QFBSkMg4qNTUVGzZswBtvvFGt34MZoGoRLakEHVzN0cHVHHP7N0X83WwcikzDoag0nLv1ENfuKHDtjgLfHI6FnYmeGKY6N7CArjaHuxERERHRi9OgQQM4Ozvju+++w6hRo8TtTk5OsLe3x5o1a/Daa6+pvOaVV15B165d0aVLF9ja2iIhIQFz585F48aN4eHhAQD4/vvv0aVLF/j7++Ozzz6Dm5sbrl27htmzZ8PBwQGff/55tR4Xv1XXYu5WRnjbxx1bJ3XBuQ/9sPSVlujjaQN9HS2kZObjt9O38Mbas2j7aTCmbLyAHRdvIyO3UNNlExEREVE94evri6ysLPTo0UNlu4+PD7KyslQmkAAAf39/7Nq1CwMHDkTjxo0xZswYeHh44ODBg9DWLr3206hRI5w7dw4NGjTAsGHD4O7ujokTJ8LX1xdhYWEwNzev1mOSCJocgaVhCoUCJiYmyMzMhFwu13Q5VSa/qASn4u8hOLJ0Ioq7WQViW+lVLDP09rRF76Y2cLYw0GClRERERPVbfn4+EhIS4ObmBj09rgNa3Z70eaubDRig6mCAepRSKeDy7UwER6biUGQ6YtKyVNob2xiht6cN/JraoJWjKaRSjpsiIiIielEYoF6sqghQlb6Fr2wqQXt7e0gkEuzYsUOlPS0tDWPHjoW9vT0MDAzQt29fxMbGlit8ypQpsLCwgJGREYYOHYq0tDSVPomJiQgICICBgQGsra0xe/ZsFBcXq/QJCQlB27ZtIZPJ0LBhQ6xbt66yh1PnSaUStHYyxWx/DxwI8kbobF/MG+AJrwYW0JJKcD0tGyuPxuPlH06h06LDmLvtMg5HpSFfQytgExERERHVZJUOUDk5OWjVqhVWrlxZrk0QBAwePBg3btzA33//jYsXL8LFxQV+fn7IyckR+wUFBWHXrl3YunUrjh07hjt37mDIkCFie0lJCQICAlBYWIhTp05h/fr1WLduHebPny/2SUhIQEBAAHx9fREREYHAwEBMmDABBw5oZuXm2sLZwgDju7lh08TOuPBRb3wzojUCWtrBSKaNu1kF2HQ2CePXn0ObT4Ix8ddz+ONcEu5nFzx9x0RERERE9cBz3cInkUiwfft2DB48GABw/fp1NGnSBFevXkWzZs0AAEqlEra2tvjiiy8wYcIEZGZmwsrKChs3bsQrr7wCAIiOjkbTpk0RFhaGzp07Y9++fRgwYADu3LkDGxsbAMDq1asxZ84c3L17F7q6upgzZw727NmDq1evivWMGDECGRkZ2L9/v1r114db+NRVWKzE6Rv3cSgqDYci03AnM19sk0iAds5m8PMsndXP3cpIg5USERER1R28he/F0sgtfE9SUFB6peLRYqRSKWQyGU6cOAEAOH/+PIqKiuDn5yf28fDwgLOzM8LCwgAAYWFhaNGihRiegNIZORQKBa5duyb2eXQfZX3K9vG4+hQKhcqDSulqS+Hd2AqfvNQcJ9/vid3TuiHQrxGaO8ghCMC5Ww+xeF80ei07hp5fheCLvVE4m/AAJcp6O4SOiIiIqMrU42kJXqiq+JyrdB2osiA0d+5c/O9//4OhoSGWL1+O5ORkpKSkAChd4EpXVxempqYqr7WxsUFqaqrY59HwVNZe1vakPgqFAnl5eSqLapVZtGgRFi5cWCXHWpdJJBI0dzBBcwcTBPo1xp2MPByOSkNwVDrC4u/hxr0crAm9gTWhN2BmoIOeHjbo7WmN7o2sYCjj0mJERERE6tLR0QEA5ObmVvj9lapWbm4ugH8/92dRpd92dXR0sG3bNowfPx7m5ubQ0tKCn58f+vXrVyNS9dy5czFz5kzxuUKhgJOTkwYrqh3sTfUx2ssVo71ckZVfhNDr93AoKg1HotPxMLcIf11Ixl8XkqGrLUVXdwv09rRFr6bWsJHzMjQRERHRk2hpacHU1BTp6ekAAAMDA0gknBW5qgmCgNzcXKSnp8PU1BRaWlrPvK8qv1zQrl07REREIDMzE4WFhbCyskKnTp3Qvn17AICtrS0KCwuRkZGhchUqLS0Ntra2Yp+zZ8+q7Ldslr5H+/x35r60tDTI5fLHpneZTAaZTFYlx1lfGevpIKClHQJa2qGoRIlzNx/iUFQagiPTkPggF0dj7uJozF1gO9DK0QR+TW3Qu5kNmtgY85cBERERUQXKvt+WhSiqPqampuLn/ayq7X4rExMTAEBsbCzOnTuHTz/9FEBpwNLR0cHhw4cxdOhQAEBMTAwSExPh5eUFAPDy8sLnn3+O9PR0WFtbAwCCg4Mhl8vh6ekp9tm7d6/KewYHB4v7oOqnoyWFl7sFvNwt8FFAU8SmZyM4sjRMRSRl4FJyJi4lZ2JZ8HU4munDr6kN+njaoIObOXS0qnT4HREREVGtJZFIYGdnB2traxQVFWm6nDpLR0fnua48lan0LHzZ2dmIi4sDALRp0wZff/01fH19YW5uDmdnZ2zduhVWVlZwdnbGlStXMGPGDLRr1w5//fWXuI/Jkydj7969WLduHeRyOaZNmwYAOHXqFIDSacxbt24Ne3t7LFmyBKmpqRg9ejQmTJiAL774AkDpNObNmzfHlClTMG7cOBw5cgTTp0/Hnj174O/vr9axcBa+6pOelY8jUekIjkzDibh7KChWim3GetrwbWKN3p428GliBbnes9+DSkRERERUFdTNBpUOUCEhIfD19S23fcyYMVi3bh2+/fZbLF26FGlpabCzs8Mbb7yBefPmQVdXV+ybn5+Pd999F5s2bUJBQQH8/f3xww8/qFxOu3XrFiZPnoyQkBAYGhpizJgxWLx4MbS1/71oFhISgqCgIERGRsLR0RHz5s3D2LFj1T4WBqgXI7ewGCdiS8dNHY5Kx/2cQrFNWypB5wYW6O1pg15NreFoZqDBSomIiIiovqq2AFWXMEC9eCVKARFJD3EwsnS9qfi7OSrtnnby0vWmmtqguYOc46aIiIiI6IVggFIDA5Tm3bib/c/ivek4d+sBHl1WylauBz9Pa/g1tYGXuwVk2s9/zyoRERERUUUYoNTAAFWzPMgpxNHo0nFTobF3kVtYIrYZ6mrBp4kV/JraoKeHNUwNdJ+wJyIiIiKiymGAUgMDVM2VX1SCsBv3EfzPrX7pWQVim5ZUgvYuZujtaQO/pjZwtTTUYKVEREREVBcwQKmBAap2UCoFXLmdKa43FZ2apdLeyNoIfv+EqTZOppBKOW6KiIiIiCqHAUoNDFC1U9KD3NJxU1FpOHPjAYofGThlaaSLXh428PO0QbeGltDX5bgpIiIiIno6Big1MEDVfpm5RQi5no5DUekIiU5HVkGx2KanI0W3hlbo7WmNnh42sDKWabBSIiIiIqrJGKDUwABVtxQWK3E24YF4q9/tjDyxTSIB2jiZws/TBn08beBuZcQp0omIiIhIxAClBgaouksQBESlZIlh6srtTJV2VwsDcRKKdi5m0NaSaqhSIiIiIqoJGKDUwABVf6Rk5uFwVOkU6WHx91FYohTbTA100LOJNXp72qB7YysYybQ1WCkRERERaQIDlBoYoOqn7IJiHL9+F8FRaTgSnY6M3CKxTVdLCi93C/HqlK2JngYrJSIiIqIXhQFKDQxQVFyixPlbDxEcmYbgqDTcup+r0t7S0QR+TUvDVFM7Y46bIiIiIqqjGKDUwABFjxIEAfF3s3Hwn8V7LyZl4NGfDgdTffHKVEc3c+hqc9wUERERUV3BAKUGBih6krtZBTgSnYbgyHSciLuL/KJ/x00Z62mjRxNr+DW1Ro8m1jDR19FgpURERET0vBig1MAARerKKyzBibh7OBSZhsPRabiXXSi2aUsl6NTAXLzVz8ncQIOVEhEREdGzYIBSAwMUPQulUsDFpAwciiq91S82PVul3cPWWLzVr4WDCaRSjpsiIiIiqukYoNTAAEVV4ea9HHG9qfCbD6B85CfKRi5Dr6Y26N3UBl7uFtDT0dJcoURERET0WAxQamCAoqr2MKcQR2PScSgqDcdi7iKnsERsM9DVgncjK/h52qCnhzXMDXU1WCkRERERPYoBSg0MUFSdCopLEBZ//59b/dKRqsgX26QSoL2LOfw8rdHb0xZuloYarJSIiIiIGKDUwABFL4ogCLh6W4Hgf271i0pRqLS7Wxmit6ctentao7WTGbQ4boqIiIjohWKAUgMDFGlK8sNcHI5KR3BkGk7fuI/iRwZOWRjqoqeHNXp72qBbI0sY6GprsFIiIiKi+oEBSg0MUFQTKPKLcCzmLoIj03A0Jh1Z+cVim0xbim4NLdHb0wY9m1rD2lhPg5USERER1V0MUGpggKKapqhEifCEBzgYmYZDUWlIfpgntkkkQGsnU/g1tUFvTxs0sjaCRMJb/YiIiIiqAgOUGhigqCYTBAExaVkIvlYapi4lZ6q0u1gY4OU2DhjXzQ1yPR0NVUlERERUNzBAqYEBimqTNEW+uHjvyfj7KCxWAgBM9HXwtk8DjO3iyvFSRERERM+IAUoNDFBUW+UUFONQVBq+OxKHuPRsAIClkS7e6dEQIzs5c8FeIiIiokpigFIDAxTVdiVKAX9H3MaKQ7FIfJALALAz0cP0Xo3wSjtH6GhJNVwhERERUe3AAKUGBiiqK4pKlNh6LhnfHo4VF+x1sTBAoF8jDGrlwHWliIiIiJ6CAUoNDFBU1+QXlWDDmUT8cDQO93MKAQCNbYwws3dj+Dez5ax9RERERI/BAKUGBiiqq3IKirHu1E3871g8FP+sK9XcQY53+zRBj8ZWDFJERERE/8EApQYGKKrrMvOK8NPxG1h7IgE5hSUAgPYuZni3TxN4uVtouDoiIiKimoMBSg0MUFRf3M8uwOpj8fg17BYK/pn+vFtDS7zbpzHaOJtpuDoiIiIizWOAUgMDFNU3qZn5+P5oLLaEJ6GopPRH36+pNWb2bgJPe/4MEBERUf3FAKUGBiiqr5Ie5OKbw7HYdiEZyn9+AwxoaYeg3o3hbmWk2eKIiIiINIABSg0MUFTfxaVnY8Wh69h9OQUAIJUAQ9s6YnqvRnAyN9BwdUREREQvDgOUGhigiEpF3lHg6+AYHIpKBwDoaEkwooMzpvZsCBu5noarIyIiIqp+DFBqYIAiUnUx8SGWHbyOE3H3AAAybSne8HLB5B4NYW6oq+HqiIiIiKqPutlAWtkdh4aGYuDAgbC3t4dEIsGOHTtU2rOzszF16lQ4OjpCX18fnp6eWL16tUqfHj16QCKRqDwmTZqk0icxMREBAQEwMDCAtbU1Zs+ejeLiYpU+ISEhaNu2LWQyGRo2bIh169ZV9nCI6BFtnM3w+4RO2PRWZ7RzMUNBsRI/Hk9A9y+P4OuDMcjMK9J0iUREREQaVekAlZOTg1atWmHlypUVts+cORP79+/H77//jqioKAQGBmLq1KnYuXOnSr+33noLKSkp4mPJkiViW0lJCQICAlBYWIhTp05h/fr1WLduHebPny/2SUhIQEBAAHx9fREREYHAwEBMmDABBw4cqOwhEdF/eLlb4M9JXvjlzQ5o7iBHTmEJvj0SB+8lR7HyaBxyC4ufvhMiIiKiOui5buGTSCTYvn07Bg8eLG5r3rw5hg8fjnnz5onb2rVrh379+uGzzz4DUHoFqnXr1lixYkWF+923bx8GDBiAO3fuwMbGBgCwevVqzJkzB3fv3oWuri7mzJmDPXv24OrVq+LrRowYgYyMDOzfv1+t+nkLH9HTCYKA/VdT8XXwdcSmZwMALI10MblHQ4zq5Aw9HS0NV0hERET0/KrtFr6n6dKlC3bu3Inbt29DEAQcPXoU169fR58+fVT6bdiwAZaWlmjevDnmzp2L3NxcsS0sLAwtWrQQwxMA+Pv7Q6FQ4Nq1a2IfPz8/lX36+/sjLCzssbUVFBRAoVCoPIjoySQSCfq1sMP+QG8sH94KzuYGuJddiE93R8L3qxBsPJOIohKlpsskIiIieiG0q3qH3333HSZOnAhHR0doa2tDKpXixx9/hLe3t9hn5MiRcHFxgb29PS5fvow5c+YgJiYG27ZtAwCkpqaqhCcA4vPU1NQn9lEoFMjLy4O+vn652hYtWoSFCxdW6fES1RdaUglebuOIAS3t8ef5ZHx7OBYpmfn4YPsVrD4Wj0C/RniptQO0pBJNl0pERERUbaolQJ0+fRo7d+6Ei4sLQkNDMWXKFNjb24tXjCZOnCj2b9GiBezs7NCrVy/Ex8fD3d29qksSzZ07FzNnzhSfKxQKODk5Vdv7EdVFOlpSvNbRGS+3ccDGM4n4ISQOiQ9yMfOPS1gVEo+ZvRvDv5ktpAxSREREVAdVaYDKy8vDBx98gO3btyMgIAAA0LJlS0REROCrr74qd8tdmU6dOgEA4uLi4O7uDltbW5w9e1alT1paGgDA1tZW/G/Ztkf7yOXyCq8+AYBMJoNMJnv2AyQikZ6OFsZ1c8OIjk5Yd+om/nfsBmLTszF5wwU0s5djVp8m6NHEChIJgxQRERHVHVU6BqqoqAhFRUWQSlV3q6WlBaXy8WMkIiIiAAB2dnYAAC8vL1y5cgXp6elin+DgYMjlcnh6eop9Dh8+rLKf4OBgeHl5VcWhEJGaDHS18U6Phgh9zxfTezaEoa4Wrt1R4M114XhldRhOxd/TdIlEREREVabSV6Cys7MRFxcnPk9ISEBERATMzc3h7OwMHx8fzJ49G/r6+nBxccGxY8fw66+/4uuvvwYAxMfHY+PGjejfvz8sLCxw+fJlBAUFwdvbGy1btgQA9OnTB56enhg9ejSWLFmC1NRUfPTRR5gyZYp4BWnSpEn4/vvv8d5772HcuHE4cuQI/vjjD+zZs6cqPhciqiQTfR3M7NMEY7q44n+hN7D+1E2cv/UQI388g64NLfBunyZo62ym6TKJiIiInkulpzEPCQmBr69vue1jxozBunXrkJqairlz5+LgwYN48OABXFxcMHHiRAQFBUEikSApKQmvv/46rl69ipycHDg5OeHll1/GRx99pDJd4K1btzB58mSEhITA0NAQY8aMweLFi6Gtra1SS1BQECIjI+Ho6Ih58+Zh7Nixah8LpzEnqj5pinx8fyQOm8MTUVRS+muml4c1ZvZpjGb2JhqujoiIiEiVutngudaBqu0YoIiqX9KDXHx7OBZ/XUiG8p/fNgEt7RDk1xgNrY00WxwRERHRPxig1MAARfTixN/NxopDsdh16Q4AQCoBXm7jiEC/RnAyN9BwdURERFTfMUCpgQGK6MWLSlFg2cHrOBRVOoumjpYEwzs4YVrPRrCR62m4OiIiIqqvGKDUwABFpDkRSRlYdjAGx2NLZ+mTaUsxurMLJvdwh4URlxsgIiKiF4sBSg0MUESad/rGfXx1IAbnbj0EABjqlq4vNaF7A5jo62i4OiIiIqovGKDUwABFVDMIgoBj1+9i2cHruHI7EwAg19PG2z7uGNvFFYayKl3zm4iIiKgcBig1MEAR1SyCIODAtVQsO3gdsenZAAALQ11M7uGO1zu7QE9HS8MVEhERUV3FAKUGBiiimqlEKWDXpTtYfug6bt3PBQDYyvUwrVdDvNrOCbraUg1XSERERHUNA5QaGKCIaraiEiX+Op+Mbw/H4k5mPgDAyVwfgb0aY3AbB2hJJRqukIiIiOoKBig1MEAR1Q75RSXYdDYRK4/G4152AQCgobURZvZujL7NbCFlkCIiIqLnxAClBgYootolt7AY60/dwupj8cjMKwIANLOX490+jeHbxBoSCYMUERERPRsGKDUwQBHVTor8Ivx0PAE/H7+BnMISAEBbZ1PM6tMEXRpaarg6IiIiqo0YoNTAAEVUuz3IKcT/jsVjfdhN5BcpAQBd3C0wy78J2jqbabg6IiIiqk0YoNTAAEVUN6Qr8vH90ThsOpuIopLSX2k9Pazxbp/GaGZvouHqiIiIqDZggFIDAxRR3ZL8MBffHo7FXxduo0RZ+qstoIUdgno3QkNrYw1XR0RERDUZA5QaGKCI6qYbd7Ox4lAsdl2+A0EApBJgcBsHBPZqDGcLA02XR0RERDUQA5QaGKCI6rboVAWWHbyO4Mg0AIC2VILhHZwwrWcj2Jroabg6IiIiqkkYoNTAAEVUP1xKysBXB2NwPPYeAEBXW4rRnV0wuYc7LI1kGq6OiIiIagIGKDUwQBHVL2du3MdXB2MQfvMhAMBAVwvjurrhre4NYGKgo+HqiIiISJMYoNTAAEVU/wiCgNDYe/jqQAyu3M4EAMj1tDHRuwHe7OoGQ5m2hiskIiIiTWCAUgMDFFH9JQgCDlxLw9fBMbielg0AMDfUxTs93PF6Zxfo6WhpuEIiIiJ6kRig1MAARUQlSgG7L9/B8uDruHk/FwBgI5dhWs9GGNbeCbraUg1XSERERC8CA5QaGKCIqExRiRLbLiTjm0OxuJOZDwBwMtfHjF6N8XIbB2hJJRqukIiIiKoTA5QaGKCI6L8Kikuw6Uwivj8aj3vZBQAAdytDBPVujP7N7SBlkCIiIqqTGKDUwABFRI+TV1iC9WE3sfpYPDJyiwAAnnZyvNunMXp6WEMiYZAiIiKqSxig1MAARURPo8gvws/HE/DziQRkFxQDANo4m2JWnybo2tBSw9URERFRVWGAUgMDFBGp62FOIVaHxmP9qZvIL1ICALwaWGCWfxO0czHTcHVERET0vBig1MAARUSVla7Ixw8h8dh4JhGFJaVByreJFd7t0wTNHUw0XB0RERE9KwYoNTBAEdGzSn6Yi+8Ox+HPC8koUZb+Gu3fwhZBfo3RyMZYw9URERFRZTFAqYEBioieV8K9HKw4dB07L92BIABSCTC4tQMC/RrD2cJA0+URERGRmhig1MAARURVJTpVgeXB13HgWhoAQFsqwavtnTC9V0PYmehruDoiIiJ6GgYoNTBAEVFVu5ycga8OXkfo9bsAAF1tKV7v5IJ3fN1haSTTcHVERET0OAxQamCAIqLqcjbhAb46EIOzNx8AAAx0tfBmV1dM7O4OEwMdDVdHRERE/8UApQYGKCKqToIg4HjsPXx1MAaXkzMBAMZ62pjYvQHe7OYGI5m2hiskIiKiMgxQamCAIqIXQRAEHIxMw9cHryMmLQsAYG6oi8k+7hjt5QI9HS0NV0hEREQMUGpggCKiF0mpFLDr8h2sOBSLhHs5AABrYxmm9WyI4R2coast1XCFRERE9RcDlBoYoIhIE4pLlNh24Ta+ORyL2xl5AABHM33M6NUIL7dxgLYWgxQREdGLpm42qPT/pUNDQzFw4EDY29tDIpFgx44dKu3Z2dmYOnUqHB0doa+vD09PT6xevVqlT35+PqZMmQILCwsYGRlh6NChSEtLU+mTmJiIgIAAGBgYwNraGrNnz0ZxcbFKn5CQELRt2xYymQwNGzbEunXrKns4REQvnLaWFMM6OOHILB8sHNQMVsYyJD/Mw+w/L6PPilDsunQHSmW9/bctIiKiGq3SASonJwetWrXCypUrK2yfOXMm9u/fj99//x1RUVEIDAzE1KlTsXPnTrFPUFAQdu3aha1bt+LYsWO4c+cOhgwZIraXlJQgICAAhYWFOHXqFNavX49169Zh/vz5Yp+EhAQEBATA19cXERERCAwMxIQJE3DgwIHKHhIRkUbItLUwposrQmf7Ym4/D5ga6ODG3RxM23QRAd+dwKHINNTjmwSIiIhqpOe6hU8ikWD79u0YPHiwuK158+YYPnw45s2bJ25r164d+vXrh88++wyZmZmwsrLCxo0b8corrwAAoqOj0bRpU4SFhaFz587Yt28fBgwYgDt37sDGxgYAsHr1asyZMwd3796Frq4u5syZgz179uDq1avi+4wYMQIZGRnYv3+/WvXzFj4iqkmy8ovw84kE/HQ8AdkFpVfcWzuZYlafJuja0AISiUTDFRIREdVd1XYL39N06dIFO3fuxO3btyEIAo4ePYrr16+jT58+AIDz58+jqKgIfn5+4ms8PDzg7OyMsLAwAEBYWBhatGghhicA8Pf3h0KhwLVr18Q+j+6jrE/ZPipSUFAAhUKh8iAiqimM9XQQ6NcYx9/zxSQfd+jpSBGRlIHXfz6D1348jfO3Hmi6RCIionqvygPUd999B09PTzg6OkJXVxd9+/bFypUr4e3tDQBITU2Frq4uTE1NVV5nY2OD1NRUsc+j4amsvaztSX0UCgXy8vIqrG3RokUwMTERH05OTs99vEREVc3MUBfv9/NA6Hu+GNvFFbpaUpy+8QBDV4Vh7C9ncfV2pqZLJCIiqreqJUCdPn0aO3fuxPnz57Fs2TJMmTIFhw4dquq3qrS5c+ciMzNTfCQlJWm6JCKix7I21sOCQc1wdHYPjOjgBC2pBCExdzHguxOY/Pt5XP9nTSkiIiJ6cbSrcmd5eXn44IMPsH37dgQEBAAAWrZsiYiICHz11Vfw8/ODra0tCgsLkZGRoXIVKi0tDba2tgAAW1tbnD17VmXfZbP0PdrnvzP3paWlQS6XQ19fv8L6ZDIZZDJZlRwrEdGL4mCqj8VDW2KSjztWHLqOvy/dwb6rqdh/LRWDWzsg0K8RXCwMNV0mERFRvVClV6CKiopQVFQEqVR1t1paWlAqlQBKJ5TQ0dHB4cOHxfaYmBgkJibCy8sLAODl5YUrV64gPT1d7BMcHAy5XA5PT0+xz6P7KOtTtg8iorrG1dIQK0a0wf4Z3ujbzBaCAGy/eBs9lx3D3G2XcSej4tuXiYiIqOpU+gpUdnY24uLixOcJCQmIiIiAubk5nJ2d4ePjg9mzZ0NfXx8uLi44duwYfv31V3z99dcAABMTE4wfPx4zZ86Eubk55HI5pk2bBi8vL3Tu3BkA0KdPH3h6emL06NFYsmQJUlNT8dFHH2HKlCniFaRJkybh+++/x3vvvYdx48bhyJEj+OOPP7Bnz56q+FyIiGqsJrbGWD26Ha4kZ+KrgzE4dv0uNp1Nwl8XbmNUJ2e806MhrIx5tZ2IiKg6VHoa85CQEPj6+pbbPmbMGKxbtw6pqamYO3cuDh48iAcPHsDFxQUTJ05EUFCQOAVvfn4+3n33XWzatAkFBQXw9/fHDz/8IN6eBwC3bt3C5MmTERISAkNDQ4wZMwaLFy+Gtra2Si1BQUGIjIyEo6Mj5s2bh7Fjx6p9LJzGnIjqgvCbD7D0QAzOJpTO0qevo4WxXV3xtncDmBroarg6IiKi2kHdbPBc60DVdgxQRFRXCIKAE3H38NWBGFxKLp2lz1imjbe8G2BcNzcYyap0yCsREVGdwwClBgYoIqprBEFAcGQavg6+jujU0ln6zAx0MLmHO97wcoWejpaGKyQiIqqZGKDUwABFRHWVUilg95UUrAi+jhv3cgAA1sYyzPBrhBEdnKEllWi4QiIiopqFAUoNDFBEVNcVlyix7eJtfHMoFrf/maXPw9YY8wd4oktDSw1XR0REVHMwQKmBAYqI6ouC4hJsOJ2Ibw7HIjOvCADQ29MGH/ZvCldLriFFRETEAKUGBigiqm8e5hRixaHr+P1MIkqUAnS0JHizqxum9mwIuZ6OpssjIiLSGAYoNTBAEVF9FZuWhU/3RCH0+l0AgIWhLt7t0wTDOzhxfBQREdVLDFBqYIAiovpMEASExNzFp3siceNu6UQTHB9FRET1FQOUGhigiIiAohIlfgu7hRWHrkORXwwA6ONpgw84PoqIiOoRBig1MEAREf3rYU4hlh+6jg2PjI8a19UNUzg+ioiI6gEGKDUwQBERlXc9LQuf7o7E8dh7AABLo9LxUcPac3wUERHVXQxQamCAIiKqmCAIOBqTjs92R4kL8XrYGmP+QE90cef4KCIiqnsYoNTAAEVE9GQVjY/yb1Y6PsrFguOjiIio7mCAUgMDFBGReh78s35U2fgoXS0p3uzmiqm+DWHM8VFERFQHMECpgQGKiKhyKhofNatPE7zK8VFERFTLMUCpgQGKiKjyBEHAkeh0fL7n3/FRnnZyzBvgCS93Cw1XR0RE9GwYoNTAAEVE9OwKi5X47fQtfPPI+Ki+zWzxQf+mcLYw0HB1RERElcMApQYGKCKi5/cgpxDLg69jw5lbUArg+CgiIqqVGKDUwABFRFR1YlJLx0ediCsbHyXDbP/GeKUdx0cREVHNxwClBgYoIqKqJQgCDkel4/O9UUh4ZHzU/IGe6NyA46OIiKjmYoBSAwMUEVH1KCxW4tewm/jmcCyyOD6KiIhqAQYoNTBAERFVr/vZBVh+6Do2nkkUx0eN7+6Gd3q4c3wUERHVKAxQamCAIiJ6MTg+ioiIajoGKDUwQBERvTgVjY9qZi/H/AGe6MTxUUREpGEMUGpggCIievEqGh/Vv4Ut5vZrCidzjo8iIiLNYIBSAwMUEZHm3M8uwNfB17HprOr4qCm+DWEk09Z0eUREVM8wQKmBAYqISPOiUxX4dHckTsbdB1A6Puo9/yZ4pZ0jpBwfRURELwgDlBoYoIiIagZBEHAoKh2f74nEzfu5AErHR308sBk6uplruDoiIqoPGKDUwABFRFSzFBSX4NdTt/Dt4VhkFZSOjwpoYYf3+3lwfBQREVUrBig1MEAREdVM9/4ZH7W5bHyUthQTurnhHY6PIiKiasIApQYGKCKimi0qpXR81Kn40vFRVsYyzPZvglfacnwUERFVLQYoNTBAERHVfIIgIDgyDZ/vjcKtf8ZHNXeQY/4Ajo8iIqKqwwClBgYoIqLao6C4BOtP3cR3h+M4PoqIiKocA5QaGKCIiGqfe9kFWHbwOraE/zs+6q3ubninR0MYcnwUERE9IwYoNTBAERHVXpF3SsdHhd34d3zUe/5NMJTjo4iI6BkwQKmBAYqIqHYTBAEHI9PwxSPjo1o4mGD+QE90cOX4KCIiUp+62UBa2R2HhoZi4MCBsLe3h0QiwY4dO1TaJRJJhY+lS5eKfVxdXcu1L168WGU/ly9fRvfu3aGnpwcnJycsWbKkXC1bt26Fh4cH9PT00KJFC+zdu7eyh0NERLWYRCKBfzNbHAzyxtx+HjCSaePK7Uy8ujoMUzZeQPLDXE2XSEREdUylA1ROTg5atWqFlStXVtiekpKi8li7di0kEgmGDh2q0u+TTz5R6Tdt2jSxTaFQoE+fPnBxccH58+exdOlSLFiwAGvWrBH7nDp1Cq+99hrGjx+PixcvYvDgwRg8eDCuXr1a2UMiIqJaTqathbd93HF0Vg+81tEZEgmw53IKei47hq8OxCDnn0kniIiIntdz3cInkUiwfft2DB48+LF9Bg8ejKysLBw+fFjc5urqisDAQAQGBlb4mlWrVuHDDz9EamoqdHV1AQDvv/8+duzYgejoaADA8OHDkZOTg927d4uv69y5M1q3bo3Vq1erVT9v4SMiqpsi7yjwye5rOH3jAQDA2liG9/p6YEgbB46PIiKiClXbLXyVkZaWhj179mD8+PHl2hYvXgwLCwu0adMGS5cuRXHxv/86GBYWBm9vbzE8AYC/vz9iYmLw8OFDsY+fn5/KPv39/REWFvbYegoKCqBQKFQeRERU93jay7Hprc5Y/Xo7OJsbID2rALO2XsLgH07i3M0Hmi6PiIhqsWoNUOvXr4exsTGGDBmisn369OnYvHkzjh49irfffhtffPEF3nvvPbE9NTUVNjY2Kq8pe56amvrEPmXtFVm0aBFMTEzEh5OT03MdHxER1VwSiQR9m9sieKY33v9nfNTl5Ey8sjoMUzk+ioiInlG1Lpixdu1ajBo1Cnp6eirbZ86cKf65ZcuW0NXVxdtvv41FixZBJpNVWz1z585VeW+FQsEQRURUx8m0tTDJxx1D2zpi2cEYbDmXhN2XUxAcmYaJ3g0wyced60cREZHaqu0K1PHjxxETE4MJEyY8tW+nTp1QXFyMmzdvAgBsbW2Rlpam0qfsua2t7RP7lLVXRCaTQS6XqzyIiKh+sDKWYfHQltg9rRs6uZmjoFiJ747EoeeyEPx1PhlKZb1d1YOIiCqh2gLUzz//jHbt2qFVq1ZP7RsREQGpVApra2sAgJeXF0JDQ1FUVCT2CQ4ORpMmTWBmZib2eXRiirI+Xl5eVXgURERU1zSzN8HmiZ2x+vW2cDLXR5qiAO9yfBQREamp0gEqOzsbERERiIiIAAAkJCQgIiICiYmJYh+FQoGtW7dWePUpLCwMK1aswKVLl3Djxg1s2LABQUFBeP3118VwNHLkSOjq6mL8+PG4du0atmzZgm+++Ubl9rsZM2Zg//79WLZsGaKjo7FgwQKcO3cOU6dOrewhERFRPVM6PsoOwUE+mNNXdXzUtE0XcTsjT9MlEhFRDVXpacxDQkLg6+tbbvuYMWOwbt06AMCaNWsQGBiIlJQUmJiYqPS7cOEC3nnnHURHR6OgoABubm4YPXo0Zs6cqTL+6fLly5gyZQrCw8NhaWmJadOmYc6cOSr72rp1Kz766CPcvHkTjRo1wpIlS9C/f3+1j4XTmBMREQCkZ+Vj2YHr+ON8EgQBkGlL8bZ3A0zq4Q4DXY6PIiKqD9TNBs+1DlRtxwBFRESPuno7E5/sjsTZhNJb+WzkMszp64HBrbl+FBFRXccApQYGKCIi+i9BELD/aiq+2BeFpAelt/K1cjTB/IGeaOdiruHqiIioujBAqYEBioiIHie/qARrTyZg5ZE45BSWAAAGtbLHnH4ecDDV13B1RERU1Rig1MAARURET/Pf8VF6OlJM9HbHJJ8GHB9FRFSHMECpgQGKiIjUxfFRRER1GwOUGhigiIioMsrGR32+NwrJD/8ZH+VkivkDPNHOxUzD1RER0fNggFIDAxQRET2LisZHvdTaHnP6esCe46OIiGolBig1MEAREdHzSFfk46uDMdh6Ppnjo4iIajkGKDUwQBERUVW4kpyJT3dH4uzN0vFRtnI9vN/PA4Na2XN8FBFRLcEApQYGKCIiqiqCIGDf1VR88cj4qNZOppg/0BNtnTk+ioiopmOAUgMDFBERVbX8ohL8fCIBK4/GIZfjo4iIag0GKDUwQBERUXVJV+Rj6YEY/Hnh3/FRk3zc8ba3O/R1tTRdHhER/QcDlBoYoIiIqLpdSc7EJ7uvIfzmQwCAnYke5vT1wEut7SGRcHwUEVFNwQClBgYoIiJ6EQRBwN4rpeOjbmeUjo9q41y6flQbjo8iIqoRGKDUwABFREQvUkXjo15u44D3+jaBnQnHRxERaRIDlBoYoIiISBPSFflYciAGf55PBsDxUURENQEDlBoYoIiISJMuJ2fgk12ROHfr3/FRZetHcXwUEdGLxQClBgYoIiLSNEEQsOdKChbtjeb4KCIiDWKAUgMDFBER1RT5RSX46fgN/BASrzI+ak5fD9ia6Gm4OiKiuo8BSg0MUEREVNOkKfKxZH8M/rpQOj5KX0cLk3zcMdG7AcdHERFVIwYoNTBAERFRTfXf8VH2JnqYw/FRRETVhgFKDQxQRERUkwmCgN2XU7B437/jo9o6m2L+wGZo7WSq2eKIiOoYBig1MEAREVFtkF9Ugh9DS8dH5RWVjo8a0sYB73F8FBFRlWGAUgMDFBER1SZpinx8uT8a2y7cBlA6Pmpyj9LxUXo6HB9FRPQ8GKDUwABFRES10aWkDHyyOxLnHxkf9X7/phjY0o7jo4iInhEDlBoYoIiIqLYSBAG7Lqdg8d4o3MnMBwC0czHD/AGeaMXxUURElcYApQYGKCIiqu3yCkvw4/EbWPXo+Ki2petH2cg5PoqISF0MUGpggCIioroiNTMfSw6ojo96p4c73uL4KCIitTBAqYEBioiI6pqIpAx8susaLiRmAAAcTPXxfj8PDOD4KCKiJ2KAUgMDFBER1UWCIGDnpTv4cl+0OD6qvYsZ5g/0REtHU80WR0RUQzFAqYEBioiI6rKKxkcNbeuI9/o24fgoIqL/YIBSAwMUERHVBymZeVi6PwbbLpaOjzLQLR0fNaE7x0cREZVhgFIDAxQREdUnFxMf4pPdkbj4yPiouf09ENCC46OIiBig1MAARURE9U3Z+KjF+6KR8s/4qA6uZpg/oBlaOJpouDoiIs1hgFIDAxQREdVXeYUl+F9oPFYfi0d+kRISSen4qNn+HB9FRPUTA5QaGKCIiKi+S8nMw5L9Mdj+yPioKb4NMb6bG8dHEVG9wgClBgYoIiKiUhcSH+KTXZGISMoAwPFRRFT/qJsNpJXdcWhoKAYOHAh7e3tIJBLs2LFDpV0ikVT4WLp0qdjnwYMHGDVqFORyOUxNTTF+/HhkZ2er7Ofy5cvo3r079PT04OTkhCVLlpSrZevWrfDw8ICenh5atGiBvXv3VvZwiIiICEBbZzNsm9wFK4a3hq1cD7cz8jB140X4rwjFzycS8DCnUNMlEhHVCJUOUDk5OWjVqhVWrlxZYXtKSorKY+3atZBIJBg6dKjYZ9SoUbh27RqCg4Oxe/duhIaGYuLEiWK7QqFAnz594OLigvPnz2Pp0qVYsGAB1qxZI/Y5deoUXnvtNYwfPx4XL17E4MGDMXjwYFy9erWyh0REREQApFIJBrdxwJFZPgj0awR9HS1cT8vGp7sj0emLw5i26SJOxt2DUllvb14hInq+W/gkEgm2b9+OwYMHP7bP4MGDkZWVhcOHDwMAoqKi4OnpifDwcLRv3x4AsH//fvTv3x/Jycmwt7fHqlWr8OGHHyI1NRW6uroAgPfffx87duxAdHQ0AGD48OHIycnB7t27xffq3LkzWrdujdWrV6tVP2/hIyIiejxFfhH+jriDLeGJuHpbIW53MtfH8PZOeKWdE2xNOOEEEdUN1XYLX2WkpaVhz549GD9+vLgtLCwMpqamYngCAD8/P0ilUpw5c0bs4+3tLYYnAPD390dMTAwePnwo9vHz81N5P39/f4SFhT22noKCAigUCpUHERERVUyup4PRnV2we1p37J7WDaM7u8BYTxtJD/Lw1cHr6LL4MCasD0dwZBqKS5SaLpeI6IXQrs6dr1+/HsbGxhgyZIi4LTU1FdbW1qpFaGvD3NwcqampYh83NzeVPjY2NmKbmZkZUlNTxW2P9inbR0UWLVqEhQsXPtcxERER1UfNHUzQ3MEEH/Rvir1XUrAlPAlnbz7Aoah0HIpKh7WxDK+2d8Sw9k5wsTDUdLlERNWmWq9ArV27FqNGjYKeXs24vD937lxkZmaKj6SkJE2XREREVKvo62phaDtH/DHJC4dm+mCidwNYGOoiPasAK4/Gw2dpCEb+eBp/R9xGflGJpsslIqpy1XYF6vjx44iJicGWLVtUttva2iI9PV1lW3FxMR48eABbW1uxT1pamkqfsudP61PWXhGZTAaZTPZsB0REREQqGlob4YP+TTGrTxMcjkrD5vAkhMbexan4+zgVfx+mBjp4uY0DRnRwRhNbY02XS0RUJartCtTPP/+Mdu3aoVWrVirbvby8kJGRgfPnz4vbjhw5AqVSiU6dOol9QkNDUVRUJPYJDg5GkyZNYGZmJvYpm5ji0T5eXl7VdUhERERUAV1tKfq1sMP6cR1x/D1fzOjVCPYmesjILcIvJ2/Cf0UoXv7hJLaEJyKnoFjT5RIRPZdKz8KXnZ2NuLg4AECbNm3w9ddfw9fXF+bm5nB2dgZQOoOFnZ0dli1bhkmTJpXbR79+/ZCWlobVq1ejqKgIb775Jtq3b4+NGzcCADIzM9GkSRP06dMHc+bMwdWrVzFu3DgsX75cnO781KlT8PHxweLFixEQEIDNmzfjiy++wIULF9C8eXO1joWz8BEREVWPEqWA47F3sflsEg5FpaH4n6nPDXW1MLCVPYZ3cEJrJ1Mu0ktENYa62aDSASokJAS+vr7lto8ZMwbr1q0DAKxZswaBgYFISUmBiYlJub4PHjzA1KlTsWvXLkilUgwdOhTffvstjIyMxD6XL1/GlClTEB4eDktLS0ybNg1z5sxR2c/WrVvx0Ucf4ebNm2jUqBGWLFmC/v37q30sDFBERETV725WAbZdSMaW8CTcuJcjbvewNcbwDk54uY0DTA10n7AHIqLqV20Bqi5hgCIiInpxBEHA2YQH2BKehD1XUlBQXDr1ua62FP2a22J4Byd0drOAVMqrUkT04jFAqYEBioiISDMy84qwM+I2Np1NQmTKv+syulgYYFh7J7zazhHW8poxiy8R1Q8MUGpggCIiItIsQRBw9bYCm8MT8XfEHWT/M8mEllSCnh7WGNHBCT6NraCtVa0rrxARMUCpgwGKiIio5sgtLMaey6WL9J679VDcbivXExfpdTI30GCFRFSXMUCpgQGKiIioZopLz8Lms0nYdvE2HuQUitu7NbTE8A5O6NPMBjJtLQ1WSER1DQOUGhigiIiIaraC4hIcikzH5vBEnIi7h7JvLWYGOhjS1hEjOjihkQ0X6SWi58cApQYGKCIiotoj6UEutp5Lwh/nkpGqyBe3t3Mxw/AOThjQ0g4GutoarJCIajMGKDUwQBEREdU+JUoBodfvYtPZRByOTkfJP4v0Gsm0MbCVPUZ0cEJLRxMu0ktElcIApQYGKCIiototPSsff52/jS3hibh5P1fc3tROjhEdnDC4tQNMDHQ0WCER1RYMUGpggCIiIqobBEHA6RsPsCU8EXuvpqLwn0V6ZdpS9G9hh+EdnNDJzZxXpYjosRig1MAARUREVPdk5hZhR8RtbDqbiOjULHG7m6UhhrV3wtB2DrA25iK9RKSKAUoNDFBERER1lyAIuJycic3hSdgZcRs5hSUAAG2pBL2aWmNEB2d4N7aClpRXpYiIAUotDFBERET1Q05B6SK9m8MTcSExQ9xuZ6KHV9s74dV2jlykl6ieY4BSAwMUERFR/XM9LQtbwpOw7UIyHuYWAQAkktJFekd0cEZvTxvoaks1XCURvWgMUGpggCIiIqq/CopLcPBaGraEJ+FE3D1xu7mhLoa2dcDwDk5oaM1FeonqCwYoNTBAEREREQAk3s/FH+eSsPV8EtIUBeL29i5mGNHRGf1b2HKRXqI6jgFKDQxQRERE9KjiEiWOXb+LTWeTcDTm30V6jWXaGNTaHiM6OKOFo4mGqySi6sAApQYGKCIiInqcNEU+/jyfjC3hSUh88O8ivc3sSxfpHdTaASb6XKSXqK5ggFIDAxQRERE9jVIp4PSN+9gcnoT9V1NRWFK6SK+eTukivSM6OKODqxkX6SWq5Rig1MAARURERJXxMKcQOyJuY/PZJMSk/btIbwNLQwzv4IQhbR1hZSzTYIVE9KwYoNTAAEVERETPQhAERCRlYEt4EnZeuoPcRxbp7e1pg+EdnNC9ERfpJapNGKDUwABFREREzyu7oBi7L93B5vAkRCRliNsdTPXxantHvNreCQ6m+porkIjUwgClBgYoIiIiqkrRqYp/Fum9jcy8fxfp9W5khREdnNCrKRfpJaqpGKDUwABFRERE1SG/qAQHrqViS3gSTsXfF7dbGuliaFtHDOvgBHcrIw1WSET/xQClBgYoIiIiqm437+X8s0hvMu5m/btIb0dXc4zo6IR+ze2gr6ulwQqJCGCAUgsDFBEREb0oxSVKHI25iy3hiTgSnY5/1uiFsZ42Brd2wPAOTmjuwEV6iTSFAUoNDFBERESkCamZ+fjzfBK2nEtC0oM8cXtzBzlGdHDGoNb2kOtxkV6iF4kBSg0MUERERKRJSqWAU/H3sTk8EQevpaks0hvQwh6vdXRCOxcu0kv0IjBAqYEBioiIiGqKBzmF2H7xNjafTURsera43d3KECM6OGNIWwdYGHGRXqLqwgClBgYoIiIiqmkEQcCFxAxsCU/ErkspyCsqXaRXR0uCPp62GN7BCd0aWkLKRXqJqhQDlBoYoIiIiKgmy8ovwq5LKdgSnohLyZnidgdTfQzv4IRX2zvCzoSL9BJVBQYoNTBAERERUW0ReUeBP84lYduFZCjyiwEAUgng09gKwzs4o1dTa+hocZFeomfFAKUGBigiIiKqbfKLSrD/aio2hyfi9I0H4nZLIxleaeeI4R2c4GZpqMEKiWonBig1MEARERFRbZZwLwdbwpPw5/lk3Mv+d5HeTm7meK2jM/o2t4WeDhfpJVIHA5QaGKCIiIioLigqUeJIdDq2hCchJObfRXrletp4uY0Dhndwhqc9v+sQPQkDlBoYoIiIiKiuuZORhz/PJ2NLeBJuZ/y7SG9LRxOM6OCMga3sYMxFeonKYYBSAwMUERER1VVKpYATcfewJTwJByNTUVRS+pVPX0cLA1raYURHZ7R1NuUivUT/UDcbVHqqltDQUAwcOBD29vaQSCTYsWNHuT5RUVEYNGgQTExMYGhoiA4dOiAxMVFs79GjByQSicpj0qRJKvtITExEQEAADAwMYG1tjdmzZ6O4uFilT0hICNq2bQuZTIaGDRti3bp1lT0cIiIiojpJKpXAu7EVVo5qi9Nze+GjgKZwtzJEXlEJtp5PxtBVp9BneSh+On4DD3IKNV0uUa1R6QCVk5ODVq1aYeXKlRW2x8fHo1u3bvDw8EBISAguX76MefPmQU9PT6XfW2+9hZSUFPGxZMkSsa2kpAQBAQEoLCzEqVOnsH79eqxbtw7z588X+yQkJCAgIAC+vr6IiIhAYGAgJkyYgAMHDlT2kIiIiIjqNAsjGSZ0b4BDM33w5yQvvNLOEXo6UsSmZ+OzPVHo/MVhTN14ASdi70GprLc3JxGp5blu4ZNIJNi+fTsGDx4sbhsxYgR0dHTw22+/PfZ1PXr0QOvWrbFixYoK2/ft24cBAwbgzp07sLGxAQCsXr0ac+bMwd27d6Grq4s5c+Zgz549uHr1qsp7Z2RkYP/+/WrVz1v4iIiIqL5S5BdhZ8QdbAlPwpXb/y7S62imj+HtnfBqeyfYmug9YQ9EdUu13cL3JEqlEnv27EHjxo3h7+8Pa2trdOrUqcLb/DZs2ABLS0s0b94cc+fORW5urtgWFhaGFi1aiOEJAPz9/aFQKHDt2jWxj5+fn8o+/f39ERYW9tj6CgoKoFAoVB5ERERE9ZFcTwevd3bBrmndsHtaN4zu7AJjPW0kP8zDsuDr6LL4MMavC8fBa6koKlFqulyiGkO7KneWnp6O7OxsLF68GJ999hm+/PJL7N+/H0OGDMHRo0fh4+MDABg5ciRcXFxgb2+Py5cvY86cOYiJicG2bdsAAKmpqSrhCYD4PDU19Yl9FAoF8vLyoK+vX66+RYsWYeHChVV5yERERES1XnMHEzR3MMEH/Zti39UUbA5PwtmEBzgcnY7D0emwMpbh1X8W6XWx4CK9VL9VaYBSKkv/deKll15CUFAQAKB169Y4deoUVq9eLQaoiRMniq9p0aIF7Ozs0KtXL8THx8Pd3b0qS1Ixd+5czJw5U3yuUCjg5ORUbe9HREREVJvo62phSFtHDGnriPi72fjjn0V672YV4IeQePwQEg+vBhYY0dEJ/s24SC/VT1UaoCwtLaGtrQ1PT0+V7U2bNsWJEyce+7pOnToBAOLi4uDu7g5bW1ucPXtWpU9aWhoAwNbWVvxv2bZH+8jl8gqvPgGATCaDTCar3EERERER1UPuVkaY278p3u3TBEei07DpbBJCY+8i7MZ9hN24DxN9HbzcxgEjOjrBw5Zjyan+qNIxULq6uujQoQNiYmJUtl+/fh0uLi6PfV1ERAQAwM7ODgDg5eWFK1euID09XewTHBwMuVwuhjMvLy8cPnxYZT/BwcHw8vKqikMhIiIiIgC62lL0bW6H9eM64sScngj0awR7Ez1k5hVh3amb6LviOF5aeRKbzyYiu6D46TskquUqPQtfdnY24uLiAABt2rTB119/DV9fX5ibm8PZ2Rnbt2/H8OHDsXLlSvj6+mL//v0IDAxESEgIunXrhvj4eGzcuBH9+/eHhYUFLl++jKCgIDg6OuLYsWMASqcxb926Nezt7bFkyRKkpqZi9OjRmDBhAr744gsApdOYN2/eHFOmTMG4ceNw5MgRTJ8+HXv27IG/v79ax8JZ+IiIiIgqr0Qp4HjsXWwJT0JwZBqK/5n63EBXC+/0cMckH3doa1Xpv9MTVTt1s0GlA1RISAh8fX3LbR8zZoy4kO3atWuxaNEiJCcno0mTJli4cCFeeuklAEBSUhJef/11XL16FTk5OXBycsLLL7+Mjz76SKXQW7duYfLkyQgJCYGhoSHGjBmDxYsXQ1v737sOQ0JCEBQUhMjISDg6OmLevHkYO3as2sfCAEVERET0fO5lF+Cv88nYEp6EG/dyAABtnU2xfHhrTjhBtUq1Bai6hAGKiIiIqGoIgoBtF25jwc5ryCoohoGuFuYN8MSIDk6QSCSaLo/oqTSyDhQRERER1U8SiQRD2zliX2B3dHIzR25hCeZuu4IJ68/hblaBpssjqjIMUERERERUZRzNDLDxrc74oL8HdLWkOBydjr4rQhEcmfb0FxPVAgxQRERERFSltKQSTPR2x99Tu8LD1hj3cwrx1q/nMOfPy5ypj2o9BigiIiIiqhZN7eT4e2pXvO3dABIJsOVcEvp/cxznbj7QdGlEz4wBioiIiIiqjUxbC3P7N8WmtzrDwVQfiQ9yMex/YVh6IBqFxUpNl0dUaQxQRERERFTtOjewwL7A7hja1hFKAVh5NB4v/3ASsWlZmi6NqFIYoIiIiIjohZDr6WDZsFZYNaotTA10cO2OAgHfncDaEwlQKuvtyjpUyzBAEREREdEL1a+FHQ4GesOnsRUKi5X4ZHckRq89g5TMPE2XRvRUDFBERERE9MJZy/Ww7s0O+HRwc+jpSHEy7j78l4di56U7mi6N6IkYoIiIiIhIIyQSCUZ3dsGe6d3RytEEivxiTN90EdM3XURmbpGmyyOqEAMUEREREWmUu5UR/pzcBTN6NYKWVIKdl+7Af0UoTsTe03RpROUwQBERERGRxuloSRHUuzH+nOQFN0tDpCry8frPZ7Bw1zXkF5VoujwiEQMUEREREdUYbZzNsGd6N7ze2RkA8MvJmxjw3QlcvZ2p4cqISjFAEREREVGNYqCrjc8Gt8AvYzvA0kiGuPRsDF55EiuPxqGE052ThjFAEREREVGN5OthjYNB3vBvZoNipYClB2Iw/H9hSLyfq+nSqB5jgCIiIiKiGsvcUBerX2+Hr15tBSOZNs7deoh+34RiS3giBIFXo+jFY4AiIiIiohpNIpHglXaO2DejOzq6miOnsARz/rqCib+dx73sAk2XR/UMAxQRERER1QpO5gbYNLEz5vbzgI6WBMGRaei7IhSHItM0XRrVIwxQRERERFRraEkleNvHHX9P6YYmNsa4l12ICb+ew9xtl5FTUKzp8qgeYIAiIiIiolrH016Ov6d2xVvd3SCRAJvOJqHfN8dx/tZDTZdGdRwDFBERERHVSno6WvgwwBMbJnSCvYkeEh/k4tXVp7DsYAyKSpSaLo/qKAYoIiIiIqrVurhbYl+gN4a0cYBSAL47EochP5xCXHqWpkujOogBioiIiIhqPRN9HXw9vDVWjmwLUwMdXLmdiYBvT2DdyQQoufguVSEGKCIiIiKqMwJa2uFAoDe6N7JEQbESC3ZFYswvZ5Gama/p0qiOYIAiIiIiojrFRq6HX8d1xCcvNYNMW4rjsffgvyIUuy7d0XRpVAcwQBERERFRnSORSPCGlyv2TO+Olo4myMwrwrRNFzFj80Vk5hVpujyqxRigiIiIiKjOamhthL8md8H0ng0hlQB/R9xB3xWhOBV3T9OlUS3FAEVEREREdZqOlhQz+zTBn5O7wNXCACmZ+Rj50xl8ujsS+UUlmi6PahkGKCIiIiKqF9o6m2HP9O4Y2ckZAPDziQQM+v4Ert3J1HBlVJswQBERERFRvWEo08YXL7fA2rHtYWkkw/W0bAxeeRKrQuJRwunOSQ0MUERERERU7/T0sMGBwO7o42mDohIBX+6Pxog1YUh6kKvp0qiGY4AiIiIionrJwkiG/41uhyWvtIShrhbCbz5E3xWh+ONcEgSBV6OoYgxQRERERFRvSSQSDGvvhP2B3ujgaoacwhK89+dlTPr9PO5nF2i6PKqBGKCIiIiIqN5zMjfA5olemNPXAzpaEhy4lgb/FcdxJDpN06VRDcMARUREREQEQEsqweQe7tgxpSsa2xjhXnYBxq07hw+2X0FOQbGmy6MaggGKiIiIiOgRzexNsHNqN0zo5gYA2HgmEQHfHseFxIcaroxqgkoHqNDQUAwcOBD29vaQSCTYsWNHuT5RUVEYNGgQTExMYGhoiA4dOiAxMVFsz8/Px5QpU2BhYQEjIyMMHToUaWmql0cTExMREBAAAwMDWFtbY/bs2SguVk3+ISEhaNu2LWQyGRo2bIh169ZV9nCIiIiIiMrR09HCRwM8sXFCJ9iZ6OHm/Vy8suoUvj4Yg6ISpabLIw2qdIDKyclBq1atsHLlygrb4+Pj0a1bN3h4eCAkJASXL1/GvHnzoKenJ/YJCgrCrl27sHXrVhw7dgx37tzBkCFDxPaSkhIEBASgsLAQp06dwvr167Fu3TrMnz9f7JOQkICAgAD4+voiIiICgYGBmDBhAg4cOFDZQyIiIiIiqlCXhpbYH+iNwa3toRSAb4/EYeiqU4i/m63p0khDJMJzzNEokUiwfft2DB48WNw2YsQI6Ojo4LfffqvwNZmZmbCyssLGjRvxyiuvAACio6PRtGlThIWFoXPnzti3bx8GDBiAO3fuwMbGBgCwevVqzJkzB3fv3oWuri7mzJmDPXv24OrVqyrvnZGRgf3791f43gUFBSgo+Hc2FYVCAScnJ2RmZkIulz/rx0BERERE9cCuS3fw4fYrUOQXQ09Hig/6N8Xozi6QSCSaLo2qgEKhgImJyVOzQZWOgVIqldizZw8aN24Mf39/WFtbo1OnTiq3+Z0/fx5FRUXw8/MTt3l4eMDZ2RlhYWEAgLCwMLRo0UIMTwDg7+8PhUKBa9euiX0e3UdZn7J9VGTRokUwMTERH05OTlVx2ERERERUDwxsZY+DQT7o3sgS+UVKzP/7Gsb8Eo40Rb6mS6MXqEoDVHp6OrKzs7F48WL07dsXBw8exMsvv4whQ4bg2LFjAIDU1FTo6urC1NRU5bU2NjZITU0V+zwansray9qe1EehUCAvL6/C+ubOnYvMzEzxkZSU9NzHTERERET1h62JHta/2RELBnpCpi1F6PW78F8Rij2XUzRdGr0g2lW5M6WydEDdSy+9hKCgIABA69atcerUKaxevRo+Pj5V+XaVJpPJIJPJNFoDEREREdVuUqkEY7u6oVsjSwRuicDV2wpM2XgBh6IcsGBQM5jo62i6RKpGVXoFytLSEtra2vD09FTZ3rRpU3EWPltbWxQWFiIjI0OlT1paGmxtbcU+/52Vr+z50/rI5XLo6+tX2TEREREREVWkobUxtk3uiqm+DSGVANsv3ka/FaEIi7+v6dKoGlVpgNLV1UWHDh0QExOjsv369etwcXEBALRr1w46Ojo4fPiw2B4TE4PExER4eXkBALy8vHDlyhWkp6eLfYKDgyGXy8Vw5uXlpbKPsj5l+yAiIiIiqm662lLM8m+CrZO84GxugDuZ+Rj502l8vicS+UUlmi6PqkGlb+HLzs5GXFyc+DwhIQEREREwNzeHs7MzZs+ejeHDh8Pb2xu+vr7Yv38/du3ahZCQEACAiYkJxo8fj5kzZ8Lc3BxyuRzTpk2Dl5cXOnfuDADo06cPPD09MXr0aCxZsgSpqan46KOPMGXKFPEWvEmTJuH777/He++9h3HjxuHIkSP4448/sGfPnir4WIiIiIiI1NfOxRz7ZnTHZ3siselsEn48noDjsfewfHhrNLXjbM91SaWnMQ8JCYGvr2+57WPGjBEXsl27di0WLVqE5ORkNGnSBAsXLsRLL70k9s3Pz8e7776LTZs2oaCgAP7+/vjhhx/E2/MA4NatW5g8eTJCQkJgaGiIMWPGYPHixdDW/jfzhYSEICgoCJGRkXB0dMS8efMwduxYtY9F3akKiYiIiIjUdSgyDe9vu4x72YXQ1ZLi3T6NMaF7A2hJOd15TaZuNniudaBqOwYoIiIiIqoO97IL8P5fV3AoqnTMfkc3cyx7tRWczA00XBk9jkbWgSIiIiIiIsDSSIYf32iHJUNbwlBXC2cTHqDfN8fx5/lk1OPrF3UCAxQRERERUTWQSCQY1sEJ+2Z4o52LGbILijFr6yVM/v0CHuQUaro8ekYMUERERERE1cjZwgB/vO2F2f5NoC2VYP+1VPRZHoqj0elPfzHVOAxQRERERETVTEsqwRTfhtgxpSsaWhvhXnYB3lwXjg+3X0FuYbGmy6NKYIAiIiIiInpBmjuYYPe0bnizqysAYMOZRAR8ewIXEx9qtjBSGwMUEREREdELpKejhY8HNsPv4zvBVq6HhHs5eGV1GJYHX0dRiVLT5dFTMEAREREREWlAt0aWOBDojUGt7FGiFPDN4Vi8suoU4u9ma7o0egIGKCIiIiIiDTEx0MG3r7XBNyNaQ66njUvJmQj49jh+C7vJ6c5rKAYoIiIiIiINe6m1A/YHeqNrQwvkFykx7+9rGPtLONIV+Zoujf6DAYqIiIiIqAawN9XHb+M6Yf4AT+hqS3Hs+l30WRGKfVdSNF0aPYIBioiIiIiohpBKJRjXzQ17pnVDM3s5MnKLMHnDBcz8IwKK/CJNl0dggCIiIiIiqnEa2Rhj+ztdMcXXHVIJsO3CbfRbcRynb9zXdGn1HgMUEREREVENpKstxWx/D/zxtheczQ1wOyMPr/14Gov2RqGguETT5dVbDFBERERERDVYe1dz7J3RHcPbO0EQgP+F3sBL359EVIpC06XVSwxQREREREQ1nJFMG1++0hJrRreDhaEuolOz8NL3J7EmNB4lSk53/iIxQBERERER1RJ9mtlif6A3enlYo7BEiS/2RmPkj6eR/DBX06XVGwxQRERERES1iJWxDD+NaY/FQ1rAQFcLZxIeoN+K49h2IZmL774ADFBERERERLWMRCLBiI7O2Du9O9o6myKroBgz/7iEKRsv4GFOoabLq9MYoIiIiIiIailXS0P88bYXZvVpDG2pBHuvpMJ/RShCYtI1XVqdxQBFRERERFSLaWtJMbVnI2x/pyvcrQyRnlWAsb+EY96Oq8gr5HTnVY0BioiIiIioDmjhaILd07pjbBdXAMBvp28h4NvjiEjK0GhddQ0DFBERERFRHaGvq4UFg5rh13EdYSOX4ca9HAxddQrfHIpFcYlS0+XVCQxQRERERER1jHdjKxwI9EZASzuUKAUsP3Qdr6wOQ8K9HE2XVusxQBERERER1UGmBrr4/rU2+GZEaxjraSMiKQP9vzmODWducbrz58AARURERERUR0kkErzU2gEHAr3Rxd0CeUUl+HD7VYxbF470rHxNl1crMUAREREREdVx9qb6+H18J3wU0BS62lIcjbkL/+Wh2H81RdOl1ToMUERERERE9YBUKsGE7g2wa2o3NLWT42FuESb9fgGztl5CVn6RpsurNRigiIiIiIjqkSa2xtgxpQsm93CHRAL8eT4ZfVccx9mEB5ourVZggCIiIiIiqmdk2lqY09cDWyZ6wdFMH7cz8jB8TRgW74tGQTEX330SBigiIiIionqqo5s59s3ojmHtHSEIwOpj8Ri88hRiUrM0XVqNxQBFRERERFSPGevpYMkrrbD69XYwN9RFVIoCA787gZ+O34BSyenO/4sBioiIiIiI0Le5LfYHdkdPD2sUlijx2Z4ojPrpDG5n5Gm6tBqFAYqIiIiIiAAA1sZ6+HlMe3z+cnPo62gh7MZ99F0eiu0Xk7n47j8YoIiIiIiISCSRSDCqkwv2zuiO1k6myCooRtCWS5i68SIycgs1XZ7GMUAREREREVE5bpaG+HOSF2b2bgwtqQR7rqTAf0UoQq/f1XRpGlXpABUaGoqBAwfC3t4eEokEO3bsUGkfO3YsJBKJyqNv374qfVxdXcv1Wbx4sUqfy5cvo3v37tDT04OTkxOWLFlSrpatW7fCw8MDenp6aNGiBfbu3VvZwyEiIiIiosfQ1pJieq9G2Da5CxpYGSJNUYA31p7Fx39fRV5h/ZzuvNIBKicnB61atcLKlSsf26dv375ISUkRH5s2bSrX55NPPlHpM23aNLFNoVCgT58+cHFxwfnz57F06VIsWLAAa9asEfucOnUKr732GsaPH4+LFy9i8ODBGDx4MK5evVrZQyIiIiIioido5WSKPdO6Y4yXCwBgfdgtDPjuOC4nZ2i2MA2QCM8xGkwikWD79u0YPHiwuG3s2LHIyMgod2XqUa6urggMDERgYGCF7atWrcKHH36I1NRU6OrqAgDef/997NixA9HR0QCA4cOHIycnB7t37xZf17lzZ7Ru3RqrV69Wq36FQgETExNkZmZCLper9RoiIiIiovrs2PW7mL31EtKzCqAtlWBGr0aY3MMd2lq1e3SQutmgWo4yJCQE1tbWaNKkCSZPnoz79++X67N48WJYWFigTZs2WLp0KYqLi8W2sLAweHt7i+EJAPz9/RETE4OHDx+Kffz8/FT26e/vj7CwsMfWVVBQAIVCofIgIiIiIiL1+TS2woFAb/RvYYtipYBlwdfx6v/CcPNejqZLeyGqPED17dsXv/76Kw4fPowvv/wSx44dQ79+/VBS8u89ktOnT8fmzZtx9OhRvP322/jiiy/w3nvvie2pqamwsbFR2W/Z89TU1Cf2KWuvyKJFi2BiYiI+nJycnvt4iYiIiIjqGzNDXawc2RbLh7eCsUwbFxMz0O+b49h4JrHOT3euXdU7HDFihPjnFi1aoGXLlnB3d0dISAh69eoFAJg5c6bYp2XLltDV1cXbb7+NRYsWQSaTVXVJorlz56q8t0KhYIgiIiIiInoGEokEL7dxREc3C7z7RwRO33iAD7ZfweGoNCwe2hJWxtX3vV6Tqv1GxQYNGsDS0hJxcXGP7dOpUycUFxfj5s2bAABbW1ukpaWp9Cl7bmtr+8Q+Ze0VkclkkMvlKg8iIiIiInp2Dqb62DihMz4KaApdLSkOR6fDf0UoDlx7/J1htVm1B6jk5GTcv38fdnZ2j+0TEREBqVQKa2trAICXlxdCQ0NRVFQk9gkODkaTJk1gZmYm9jl8+LDKfoKDg+Hl5VUNR0FERERERI8jlUowoXsD7JzWFR62xniQU4i3fzuP9/68hOyC4qfvoBapdIDKzs5GREQEIiIiAAAJCQmIiIhAYmIisrOzMXv2bJw+fRo3b97E4cOH8dJLL6Fhw4bw9/cHUDr5w4oVK3Dp0iXcuHEDGzZsQFBQEF5//XUxHI0cORK6uroYP348rl27hi1btuCbb75Ruf1uxowZ2L9/P5YtW4bo6GgsWLAA586dw9SpU6vgYyEiIiIiosrysJXj76ld8bZPA0gkwB/nktHvm1CE33yg6dKqTKWnMQ8JCYGvr2+57WPGjMGqVaswePBgXLx4ERkZGbC3t0efPn3w6aefihM+XLhwAe+88w6io6NRUFAANzc3jB49GjNnzlQZ/3T58mVMmTIF4eHhsLS0xLRp0zBnzhyV99y6dSs++ugj3Lx5E40aNcKSJUvQv39/tY+F05gTEREREVWPMzfuY+Yfl3A7Iw8SCTDJxx1Bfo2hq10zpztXNxs81zpQtR0DFBERERFR9cnKL8KCnZH460IyAMDTTo4VI1qjsY2xhisrT6PrQBERERERERnr6WDZsFZYNaotzAx0EJmiwIDvTuDnEwlQKmvndRwGKCIiIiIiqlb9WtjhQKA3ejSxQmGxEp/ujsTrP5/BnYw8TZdWaQxQRERERERU7azlevhlbAd8Nrg59HW0cCr+PvxXhOLviNuaLq1SGKCIiIiIiOiFkEgkeL2zC/ZM74ZWTqbIyi/GxcQMTZdVKdqaLoCIiIiIiOqXBlZG+HOSF34Lu4WRnZw1XU6lMEAREREREdELp6Mlxbhubpouo9J4Cx8REREREZGaGKCIiIiIiIjUxABFRERERESkJgYoIiIiIiIiNTFAERERERERqYkBioiIiIiISE0MUERERERERGpigCIiIiIiIlITAxQREREREZGaGKCIiIiIiIjUxABFRERERESkJgYoIiIiIiIiNTFAERERERERqYkBioiIiIiISE3ami5AkwRBAAAoFAoNV0JERERERJpUlgnKMsLj1OsAlZWVBQBwcnLScCVERERERFQTZGVlwcTE5LHtEuFpEasOUyqVuHPnDoyNjSGRSDRai0KhgJOTE5KSkiCXyzVaC9UOPGeosnjOUGXxnKHK4jlDlVHTzhdBEJCVlQV7e3tIpY8f6VSvr0BJpVI4OjpqugwVcrm8RpxAVHvwnKHK4jlDlcVzhiqL5wxVRk06X5505akMJ5EgIiIiIiJSEwMUERERERGRmhigagiZTIaPP/4YMplM06VQLcFzhiqL5wxVFs8ZqiyeM1QZtfV8qdeTSBAREREREVUGr0ARERERERGpiQGKiIiIiIhITQxQREREREREamKAIiIiIiIiUhMDVDVasGABJBKJysPDw0Nsz8/Px5QpU2BhYQEjIyMMHToUaWlpKvtITExEQEAADAwMYG1tjdmzZ6O4uPhFHwq9IE87Z3r06FGufdKkSSr74DlT/9y+fRuvv/46LCwsoK+vjxYtWuDcuXNiuyAImD9/Puzs7KCvrw8/Pz/Exsaq7OPBgwcYNWoU5HI5TE1NMX78eGRnZ7/oQ6EX4Gnny9ixY8v9nunbt6/KPni+1C+urq7lzgmJRIIpU6YA4PcZUvW086UufJfR1nQBdV2zZs1w6NAh8bm29r8feVBQEPbs2YOtW7fCxMQEU6dOxZAhQ3Dy5EkAQElJCQICAmBra4tTp04hJSUFb7zxBnR0dPDFF1+88GOhF+NJ5wwAvPXWW/jkk0/E5wYGBuKfec7UPw8fPkTXrl3h6+uLffv2wcrKCrGxsTAzMxP7LFmyBN9++y3Wr18PNzc3zJs3D/7+/oiMjISenh4AYNSoUUhJSUFwcDCKiorw5ptvYuLEidi4caOmDo2qgTrnCwD07dsXv/zyi/j8v1MM83ypX8LDw1FSUiI+v3r1Knr37o1XX30VAL/PkKqnnS9AHfguI1C1+fjjj4VWrVpV2JaRkSHo6OgIW7duFbdFRUUJAISwsDBBEARh7969glQqFVJTU8U+q1atEuRyuVBQUFCttZNmPOmcEQRB8PHxEWbMmPHYdp4z9c+cOXOEbt26PbZdqVQKtra2wtKlS8VtGRkZgkwmEzZt2iQIgiBERkYKAITw8HCxz759+wSJRCLcvn27+oqnF+5p54sgCMKYMWOEl1566bHtPF9oxowZgru7u6BUKvl9hp7q0fNFEOrGdxnewlfNYmNjYW9vjwYNGmDUqFFITEwEAJw/fx5FRUXw8/MT+3p4eMDZ2RlhYWEAgLCwMLRo0QI2NjZiH39/fygUCly7du3FHgi9MI87Z8ps2LABlpaWaN68OebOnYvc3FyxjedM/bNz5060b98er776KqytrdGmTRv8+OOPYntCQgJSU1NVfteYmJigU6dOKr9rTE1N0b59e7GPn58fpFIpzpw58+IOhqrd086XMiEhIbC2tkaTJk0wefJk3L9/X2zj+VK/FRYW4vfff8e4ceMgkUj4fYae6L/nS5na/l2Gt/BVo06dOmHdunVo0qQJUlJSsHDhQnTv3h1Xr15FamoqdHV1YWpqqvIaGxsbpKamAgBSU1NVTp6y9rI2qnuedM4YGxtj5MiRcHFxgb29PS5fvow5c+YgJiYG27ZtA8Bzpj66ceMGVq1ahZkzZ+KDDz5AeHg4pk+fDl1dXYwZM0b8e6/ovHj0d421tbVKu7a2NszNzXne1DFPO1+A0tv3hgwZAjc3N8THx+ODDz5Av379EBYWBi0tLZ4v9dyOHTuQkZGBsWPHAgC/z9AT/fd8AVAnvsswQFWjfv36iX9u2bIlOnXqBBcXF/zxxx/Q19fXYGVUUz3pnBk/fjwmTpwotrdo0QJ2dnbo1asX4uPj4e7uromSScOUSiXat28v3hfepk0bXL16FatXrxa/EBOVUed8GTFihNi/RYsWaNmyJdzd3RESEoJevXpppG6qOX7++Wf069cP9vb2mi6FaoGKzpe68F2Gt/C9QKampmjcuDHi4uJga2uLwsJCZGRkqPRJS0uDra0tAMDW1rbcLDZlz8v6UN326DlTkU6dOgGA2M5zpv6xs7ODp6enyramTZuKt36W/b1XdF48+rsmPT1dpb24uBgPHjzgeVPHPO18qUiDBg1gaWmp8nuG50v9dOvWLRw6dAgTJkwQt/H7DD1ORedLRWrjdxkGqBcoOzsb8fHxsLOzQ7t27aCjo4PDhw+L7TExMUhMTISXlxcAwMvLC1euXFH5H1VwcDDkcnm5/wFS3fToOVORiIgIABDbec7UP127dkVMTIzKtuvXr8PFxQUA4ObmBltbW5XfNQqFAmfOnFH5XZORkYHz58+LfY4cOQKlUin+j43qhqedLxVJTk7G/fv3VX7P8Hypn3755RdYW1sjICBA3MbvM/Q4FZ0vFamV32U0PYtFXfbuu+8KISEhQkJCgnDy5EnBz89PsLS0FNLT0wVBEIRJkyYJzs7OwpEjR4Rz584JXl5egpeXl/j64uJioXnz5kKfPn2EiIgIYf/+/YKVlZUwd+5cTR0SVbMnnTNxcXHCJ598Ipw7d05ISEgQ/v77b6FBgwaCt7e3+HqeM/XP2bNnBW1tbeHzzz8XYmNjhQ0bNggGBgbC77//LvZZvHixYGpqKvz999/C5cuXhZdeeklwc3MT8vLyxD59+/YV2rRpI5w5c0Y4ceKE0KhRI+G1117TxCFRNXra+ZKVlSXMmjVLCAsLExISEoRDhw4Jbdu2FRo1aiTk5+eL++H5Uv+UlJQIzs7Owpw5c8q18fsM/dfjzpe68l2GAaoaDR8+XLCzsxN0dXUFBwcHYfjw4UJcXJzYnpeXJ7zzzjuCmZmZYGBgILz88stCSkqKyj5u3rwp9OvXT9DX1xcsLS2Fd999VygqKnrRh0IvyJPOmcTERMHb21swNzcXZDKZ0LBhQ2H27NlCZmamyj54ztQ/u3btEpo3by7IZDLBw8NDWLNmjUq7UqkU5s2bJ9jY2AgymUzo1auXEBMTo9Ln/v37wmuvvSYYGRkJcrlcePPNN4WsrKwXeRj0gjzpfMnNzRX69OkjWFlZCTo6OoKLi4vw1ltvqUwnLAg8X+qjAwcOCADK/e4QBH6fofIed77Ule8yEkEQBE1fBSMiIiIiIqoNOAaKiIiIiIhITQxQREREREREamKAIiIiIiIiUhMDFBERERERkZoYoIiIiIiIiNTEAEVERERERKQmBigiIiIiIiI1MUARERERERGpiQGKiIiIiIhITQxQRERU40kkkic+FixYAADYvn07OnfuDBMTExgbG6NZs2YIDAwU97NgwQK0bt263P5v3rwJiUSCiIgIAEBISAgkEgnMzMyQn5+v0jc8PFx83zL5+fkYO3YsWrRoAW1tbQwePLiKPwEiIqopGKCIiKjGS0lJER8rVqyAXC5X2TZr1iwcPnwYw4cPx9ChQ3H27FmcP38en3/+OYqKip75fY2NjbF9+3aVbT///DOcnZ1VtpWUlEBfXx/Tp0+Hn5/fM78fERHVfNqaLoCIiOhpbG1txT+bmJhAIpGobAOAXbt2oWvXrpg9e7a4rXHjxs91NWjMmDFYu3YtXnvtNQBAXl4eNm/ejOnTp+PTTz8V+xkaGmLVqlUAgJMnTyIjI+OZ35OIiGo2XoEiIqI6wdbWFteuXcPVq1erbJ+jR4/G8ePHkZiYCAD466+/4OrqirZt21bZexARUe3CAEVERHXCtGnT0KFDB7Ro0QKurq4YMWIE1q5di4KCgmfep7W1Nfr164d169YBANauXYtx48ZVUcVERFQbMUAREVGdYGhoiD179iAuLg4fffQRjIyM8O6776Jjx47Izc195v2OGzcO69atw40bNxAWFoZRo0ZVYdVERFTbMEAREVGd4u7ujgkTJuCnn37ChQsXEBkZiS1btgAA5HI5MjMzy72mbMySiYlJubZ+/fohLy8P48ePx8CBA2FhYVGt9RMRUc3GAEVERHWWq6srDAwMkJOTAwBo0qQJkpOTkZaWptLvwoUL0NPTKze7HgBoa2vjjTfeQEhICG/fIyIizsJHRER1w4IFC5Cbm4v+/fvDxcUFGRkZ+Pbbb1FUVITevXsDAP7f3h2iKBQEcBz+L4JdMFhNL78TGAVv8IpXeCBYxKZX8DLvEEa7UQSLzfI2LJhcmbBhd/m+NmWGKQM/hmHm83mqqkrTNNnv95lMJjkej9lut2nbNoPB4OXcu90u6/X67e3T6XTK4/HI7XbL/X5//in16t8pAP4uAQXAvzCbzXI4HLJcLnO5XDIajVLXdbquS1VVSb5uk7quy2azSdM0uV6vmU6nads2q9Xq27mHw2HG4/Hb9ReLRc7n83Nc13WSpO/7H9gdAL/FR+9kBwAAKOINFAAAQCEBBQAAUEhAAQAAFBJQAAAAhQQUAABAIQEFAABQSEABAAAUElAAAACFBBQAAEAhAQUAAFBIQAEAABT6BJj68E4PEUF1AAAAAElFTkSuQmCC",
      "text/plain": [
       "<Figure size 1000x500 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "TSUM1_values = list(range(500, 800, 50))\n",
    "results = []\n",
    "for tsum1 in TSUM1_values:\n",
    "    params.clear_override()\n",
    "    params.set_override(\"TSUM1\", tsum1)\n",
    "    wofost = Wofost72_PP(params, wdp, agro)\n",
    "    wofost.run_till_terminate()\n",
    "    r = wofost.get_summary_output()\n",
    "    results.append(r[0][target_variable])\n",
    "\n",
    "df = pd.DataFrame({\"TSUM1\": TSUM1_values,\n",
    "                   target_variable: results}\n",
    "                 ).set_index(\"TSUM1\")\n",
    "r = df.plot(title=\"Impact of TSUM1 on sugar beet yield\", figsize=(10,5))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Speeding up simulating the samples\n",
    "Running over the samples will quickly become a lengthy process when many variables are being assessed for sensitivity. However, because the simulation for each sample is independent we can easily speed up the simulations by running on multiple CPUs. Python provides the `multiprocessing` module for starting a pool of worker processes that can carry out the simulations.\n",
    "\n",
    "Demonstrating the use of `multiprocessing` within a jupyter notebook is complicated because the notebook infrastructure does not support it. Therefore a separate python script is provided that does exactly the same calculations as this notebook both using sequential and parallel processing. It can be downloaded [here](SA_mp_example.py)."
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3 (ipykernel)",
   "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.10.13"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}