examples/Notebooks/flopy3_swi2package_ex4.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# FloPy\n",
    "\n",
    "### SWI2 Example 4. Upconing Below a Pumping Well in a Two-Aquifer Island System\n",
    "\n",
    "This example problem is the fourth example problem in the SWI2 documentation (http://pubs.usgs.gov/tm/6a46/) and simulates transient movement of the freshwater-seawater interface beneath an island in response to recharge and groundwater withdrawals. The island is 2,050$\\times$2,050 m and consists of two 20-m thick aquifers that extend below sea level. The aquifers are confined, storage changes are not considered (all MODFLOW stress periods are steady-state), and the top and bottom of each aquifer is horizontal. The top of the upper aquifer and the bottom of the lower aquifer are impermeable.\n",
    "\n",
    "The domain is discretized into 61 columns, 61 rows, and 2 layers, with respective cell dimensions of 50 m (`DELR`), 50 m (`DELC`), and 20 m. A total of 230 years is simulated using three stress periods with lengths of 200, 12, and 18 years, with constant time steps of 0.2, 0.1, and 0.1 years, respectively. \n",
    "\n",
    "The horizontal and vertical hydraulic conductivity of both aquifers are 10 m/d and 0.2 m/d, respectively. The effective porosity is 0.2 for both aquifers. The model is extended 500 m offshore along all sides and the ocean boundary is represented as a general head boundary condition (GHB) in model layer 1. A freshwater head of 0 m is specified at the ocean bottom in all general head boundaries. The GHB conductance that controls outflow from the aquifer into the ocean is 62.5 m$^{2}$/d and corresponds to a leakance of 0.025 d$^{-1}$ (or a resistance of 40 days).\n",
    "\n",
    "The groundwater is divided into a freshwater zone and a seawater zone, separated by an active ZETA surface between the zones (`NSRF=1`) that approximates the 50-percent seawater salinity contour. Fluid density is represented using the stratified density option (`ISTRAT=1`). The dimensionless density difference ($\\nu$) between freshwater and saltwater is 0.025. The tip and toe tracking parameters are a `TOESLOPE` and `TIPSLOPE` of 0.005, a default `ALPHA` of 0.1, and a default `BETA` of 0.1. Initially, the interface between freshwater and saltwater is 1 m below land surface on the island and at the top of the upper aquifer offshore. The SWI2 `ISOURCE` parameter is set to -2 in cells having GHBs so that water that infiltrates into the aquifer from the GHB cells is saltwater (zone 2), whereas water that flows out of the model at the GHB cells is identical to water at the top of the aquifer. `ISOURCE` in layer 2, row 31, column 36 is set to 2 so that a saltwater well may be simulated in the third stress period of simulation 2. In all other cells, the SWI2 `ISOURCE` parameter is set to 0, indicating boundary conditions have water that is identical to water at the top of the aquifer and can be either freshwater or saltwater, depending on the elevation of the active `ZETA` surface in the cell.\n",
    "\n",
    "A constant recharge rate of 0.4 millimeters per day (mm/d) is used in all three stress periods. The development of the freshwater lens is simulated for 200 years, after which a pumping well having a withdrawal rate of 250 m$^3$/d is started in layer 1, row 31, column 36. For the first simulation (simulation 1), the well pumps for 30 years, after which the interface almost reaches the top of the upper aquifer layer. In the second simulation (simulation 2), an additional well withdrawing\n",
    "saltwater at a rate of 25 m$^3$/d is simulated below the freshwater well in layer 2 , row 31, column 36, 12 years after the freshwater groundwater withdrawal begins in the well in layer 1. The saltwater well is intended to prevent the interface from\n",
    "upconing into the upper aquifer (model layer)."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Import `numpy` and `matplotlib`, set all figures to be inline, import `flopy.modflow` and `flopy.utils`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "flopy is installed in /Users/jdhughes/Documents/Development/flopy_git/flopy_fork/flopy\n",
      "3.7.3 | packaged by conda-forge | (default, Jul  1 2019, 14:38:56) \n",
      "[Clang 4.0.1 (tags/RELEASE_401/final)]\n",
      "numpy version: 1.17.3\n",
      "matplotlib version: 3.1.1\n",
      "flopy version: 3.3.1\n"
     ]
    }
   ],
   "source": [
    "import os\n",
    "import sys\n",
    "import platform\n",
    "import numpy as np\n",
    "import matplotlib as mpl\n",
    "import matplotlib.pyplot as plt\n",
    "\n",
    "# run installed version of flopy or add local path\n",
    "try:\n",
    "    import flopy\n",
    "except:\n",
    "    fpth = os.path.abspath(os.path.join('..', '..'))\n",
    "    sys.path.append(fpth)\n",
    "    import flopy\n",
    "\n",
    "print(sys.version)\n",
    "print('numpy version: {}'.format(np.__version__))\n",
    "print('matplotlib version: {}'.format(mpl.__version__))\n",
    "print('flopy version: {}'.format(flopy.__version__))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Define model name of your model and the location of MODFLOW executable. All MODFLOW files and output will be stored in the subdirectory defined by the workspace. Create a model named `ml` and specify that this is a MODFLOW-2005 model."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Set name of MODFLOW exe\n",
    "#  assumes executable is in users path statement\n",
    "exe_name = 'mf2005'\n",
    "if platform.system() == 'Windows':\n",
    "    exe_name = 'mf2005.exe'\n",
    "\n",
    "workspace = os.path.join('data')\n",
    "#make sure workspace directory exists\n",
    "if not os.path.exists(workspace):\n",
    "    os.makedirs(workspace)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Define the number of layers, rows and columns. The heads are computed quasi-steady state (hence a steady MODFLOW run) while the interface will move. There are three stress periods with a length of 200, 12, and 18 years and 1,000, 120, and 180 steps. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "ncol = 61\n",
    "nrow = 61\n",
    "nlay = 2\n",
    "\n",
    "nper = 3\n",
    "perlen = [365.25 * 200., 365.25 * 12., 365.25 * 18.]\n",
    "nstp = [1000, 120, 180]\n",
    "save_head = [200, 60, 60]\n",
    "steady = True"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Specify the cell size along the rows (`delr`) and along the columns (`delc`) and the top and bottom of the aquifer for the `DIS` package. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "# dis data\n",
    "delr, delc = 50.0, 50.0\n",
    "botm = np.array([-10., -30., -50.])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Define the `IBOUND` array and starting heads for the `BAS` package. The corners of the model are defined to be inactive."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [],
   "source": [
    "# bas data\n",
    "# ibound - active except for the corners\n",
    "ibound = np.ones((nlay, nrow, ncol), dtype= np.int)\n",
    "ibound[:, 0, 0] = 0\n",
    "ibound[:, 0, -1] = 0\n",
    "ibound[:, -1, 0] = 0\n",
    "ibound[:, -1, -1] = 0\n",
    "# initial head data\n",
    "ihead = np.zeros((nlay, nrow, ncol), dtype=np.float)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Define the layers to be confined and define the horizontal and vertical hydraulic conductivity of the aquifer for the `LPF` package."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [],
   "source": [
    "# lpf data\n",
    "laytyp = 0\n",
    "hk = 10.\n",
    "vka = 0.2"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Define the boundary condition data for the model"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [],
   "source": [
    "# boundary condition data\n",
    "# ghb data\n",
    "colcell, rowcell = np.meshgrid(np.arange(0, ncol), np.arange(0, nrow))\n",
    "index = np.zeros((nrow, ncol), dtype=np.int)\n",
    "index[:, :10] = 1\n",
    "index[:, -10:] = 1\n",
    "index[:10, :] = 1\n",
    "index[-10:, :] = 1\n",
    "nghb = np.sum(index)\n",
    "lrchc = np.zeros((nghb, 5))\n",
    "lrchc[:, 0] = 0\n",
    "lrchc[:, 1] = rowcell[index == 1]\n",
    "lrchc[:, 2] = colcell[index == 1]\n",
    "lrchc[:, 3] = 0.\n",
    "lrchc[:, 4] = 50.0 * 50.0 / 40.0\n",
    "# create ghb dictionary\n",
    "ghb_data = {0:lrchc}\n",
    "\n",
    "# recharge data\n",
    "rch = np.zeros((nrow, ncol), dtype=np.float)\n",
    "rch[index == 0] = 0.0004\n",
    "# create recharge dictionary\n",
    "rch_data = {0: rch}\n",
    "\n",
    "# well data\n",
    "nwells = 2\n",
    "lrcq = np.zeros((nwells, 4))\n",
    "lrcq[0, :] = np.array((0, 30, 35, 0))\n",
    "lrcq[1, :] = np.array([1, 30, 35, 0])\n",
    "lrcqw = lrcq.copy()\n",
    "lrcqw[0, 3] = -250\n",
    "lrcqsw = lrcq.copy()\n",
    "lrcqsw[0, 3] = -250.\n",
    "lrcqsw[1, 3] = -25.\n",
    "# create well dictionary\n",
    "base_well_data = {0:lrcq, 1:lrcqw}\n",
    "swwells_well_data = {0:lrcq, 1:lrcqw, 2:lrcqsw}"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [],
   "source": [
    "# swi2 data\n",
    "nadptmx = 10\n",
    "nadptmn = 1\n",
    "nu = [0, 0.025]\n",
    "numult = 5.0\n",
    "toeslope = nu[1] / numult  #0.005\n",
    "tipslope = nu[1] / numult  #0.005\n",
    "z1 = -10.0 * np.ones((nrow, ncol))\n",
    "z1[index == 0] = -11.0\n",
    "z = np.array([[z1, z1]])\n",
    "iso = np.zeros((nlay, nrow, ncol), dtype=np.int)\n",
    "iso[0, :, :][index == 0] = 1\n",
    "iso[0, :, :][index == 1] = -2\n",
    "iso[1, 30, 35] = 2\n",
    "ssz=0.2\n",
    "# swi2 observations\n",
    "obsnam = ['layer1_', 'layer2_']\n",
    "obslrc=[[0, 30, 35], [1, 30, 35]]\n",
    "nobs = len(obsnam)\n",
    "iswiobs = 1051"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Create output control (OC) data using words "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [],
   "source": [
    "# oc data\n",
    "spd = {(0,199): ['print budget', 'save head'],\n",
    "       (0,200): [],\n",
    "       (0,399): ['print budget', 'save head'],\n",
    "       (0,400): [],\n",
    "       (0,599): ['print budget', 'save head'],\n",
    "       (0,600): [],\n",
    "       (0,799): ['print budget', 'save head'],\n",
    "       (0,800): [],\n",
    "       (0,999): ['print budget', 'save head'],\n",
    "       (1,0): [],\n",
    "       (1,59): ['print budget', 'save head'],\n",
    "       (1,60): [],\n",
    "       (1,119): ['print budget', 'save head'],\n",
    "       (1,120): [],\n",
    "       (2,0): [],\n",
    "       (2,59): ['print budget', 'save head'],\n",
    "       (2,60): [],\n",
    "       (2,119): ['print budget', 'save head'],\n",
    "       (2,120): [],\n",
    "       (2,179): ['print budget', 'save head']}"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Create the model with the freshwater well (Simulation 1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "ModflowSwi2: specification of nobs is deprecated.\n"
     ]
    }
   ],
   "source": [
    "modelname = 'swiex4_s1'\n",
    "ml = flopy.modflow.Modflow(modelname, version='mf2005', exe_name=exe_name, model_ws=workspace)\n",
    "\n",
    "discret = flopy.modflow.ModflowDis(ml, nlay=nlay, nrow=nrow, ncol=ncol, laycbd=0,\n",
    "                                   delr=delr, delc=delc, top=botm[0], botm=botm[1:],\n",
    "                                   nper=nper, perlen=perlen, nstp=nstp)\n",
    "bas = flopy.modflow.ModflowBas(ml, ibound=ibound, strt=ihead)\n",
    "lpf = flopy.modflow.ModflowLpf(ml, laytyp=laytyp, hk=hk, vka=vka)\n",
    "wel = flopy.modflow.ModflowWel(ml, stress_period_data=base_well_data)\n",
    "ghb = flopy.modflow.ModflowGhb(ml, stress_period_data=ghb_data)\n",
    "rch = flopy.modflow.ModflowRch(ml, rech=rch_data)\n",
    "swi = flopy.modflow.ModflowSwi2(ml, nsrf=1, istrat=1, toeslope=toeslope, tipslope=tipslope, nu=nu,\n",
    "                                zeta=z, ssz=ssz, isource=iso, nsolver=1,\n",
    "                                nadptmx=nadptmx, nadptmn=nadptmn, \n",
    "                                nobs=nobs, iswiobs=iswiobs, obsnam=obsnam, obslrc=obslrc, iswizt=55)\n",
    "oc = flopy.modflow.ModflowOc(ml, stress_period_data=spd)\n",
    "pcg = flopy.modflow.ModflowPcg(ml, hclose=1.0e-6, rclose=3.0e-3, mxiter=100, iter1=50)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Write the simulation 1 MODFLOW input files and run the model"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(True, [])"
      ]
     },
     "execution_count": 11,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "ml.write_input()\n",
    "ml.run_model(silent=True)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Create the model with the saltwater well (Simulation 2)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "ModflowSwi2: specification of nobs is deprecated.\n"
     ]
    }
   ],
   "source": [
    "modelname2 = 'swiex4_s2'\n",
    "ml2 = flopy.modflow.Modflow(modelname2, version='mf2005', exe_name=exe_name, model_ws=workspace)\n",
    "\n",
    "discret = flopy.modflow.ModflowDis(ml2, nlay=nlay, nrow=nrow, ncol=ncol, laycbd=0,\n",
    "                                   delr=delr, delc=delc, top=botm[0], botm=botm[1:],\n",
    "                                   nper=nper, perlen=perlen, nstp=nstp)\n",
    "bas = flopy.modflow.ModflowBas(ml2, ibound=ibound, strt=ihead)\n",
    "lpf = flopy.modflow.ModflowLpf(ml2, laytyp=laytyp, hk=hk, vka=vka)\n",
    "wel = flopy.modflow.ModflowWel(ml2, stress_period_data=swwells_well_data)\n",
    "ghb = flopy.modflow.ModflowGhb(ml2, stress_period_data=ghb_data)\n",
    "rch = flopy.modflow.ModflowRch(ml2, rech=rch_data)\n",
    "swi = flopy.modflow.ModflowSwi2(ml2, nsrf=1, istrat=1, \n",
    "                                toeslope=toeslope, tipslope=tipslope, nu=nu,\n",
    "                                zeta=z, ssz=ssz, isource=iso, nsolver=1,\n",
    "                                nadptmx=nadptmx, nadptmn=nadptmn,\n",
    "                                nobs=nobs, iswiobs=iswiobs, obsnam=obsnam, obslrc=obslrc, iswizt=55)\n",
    "oc = flopy.modflow.ModflowOc(ml2, stress_period_data=spd)\n",
    "pcg = flopy.modflow.ModflowPcg(ml2, hclose=1.0e-6, rclose=3.0e-3, mxiter=100, iter1=50)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Write the simulation 2 MODFLOW input files and run the model"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(True, [])"
      ]
     },
     "execution_count": 13,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "ml2.write_input()\n",
    "ml2.run_model(silent=True)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Load the simulation 1 `ZETA` data and `ZETA` observations."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {},
   "outputs": [],
   "source": [
    "# read base model zeta\n",
    "zfile = flopy.utils.CellBudgetFile(os.path.join(ml.model_ws, modelname+'.zta'))\n",
    "kstpkper = zfile.get_kstpkper()\n",
    "zeta = []\n",
    "for kk in kstpkper:\n",
    "    zeta.append(zfile.get_data(kstpkper=kk, text='ZETASRF  1')[0])\n",
    "zeta = np.array(zeta)\n",
    "# read swi obs\n",
    "zobs = np.genfromtxt(os.path.join(ml.model_ws, modelname+'.zobs.out'), names=True)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Load the simulation 2 `ZETA` data and `ZETA` observations."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {},
   "outputs": [],
   "source": [
    "# read saltwater well model zeta\n",
    "zfile2 = flopy.utils.CellBudgetFile(os.path.join(ml2.model_ws, modelname2+'.zta'))\n",
    "kstpkper = zfile2.get_kstpkper()\n",
    "zeta2 = []\n",
    "for kk in kstpkper:\n",
    "    zeta2.append(zfile2.get_data(kstpkper=kk, text='ZETASRF  1')[0])\n",
    "zeta2 = np.array(zeta2)\n",
    "# read swi obs\n",
    "zobs2 = np.genfromtxt(os.path.join(ml2.model_ws, modelname2+'.zobs.out'), names=True)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Create arrays for the x-coordinates and the output years"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {},
   "outputs": [],
   "source": [
    "x = np.linspace(-1500, 1500, 61)\n",
    "xcell = np.linspace(-1500, 1500, 61) + delr / 2.\n",
    "xedge = np.linspace(-1525, 1525, 62)\n",
    "years = [40, 80, 120, 160, 200, 6, 12, 18, 24, 30]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Define figure dimensions and colors used for plotting `ZETA` surfaces"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {},
   "outputs": [],
   "source": [
    "# figure dimensions\n",
    "fwid, fhgt = 8.00, 5.50\n",
    "flft, frgt, fbot, ftop = 0.125, 0.95, 0.125, 0.925\n",
    "\n",
    "# line color definition\n",
    "icolor = 5\n",
    "colormap = plt.cm.jet  #winter\n",
    "cc = []\n",
    "cr = np.linspace(0.9, 0.0, icolor)\n",
    "for idx in cr:\n",
    "    cc.append(colormap(idx))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Recreate **Figure 9** from the SWI2 documentation (http://pubs.usgs.gov/tm/6a46/)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAiAAAAFuCAYAAAC4M6OdAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8QZhcZAAAgAElEQVR4nOzdeVyUZff48Q8wLOJamWlhKpIkyGKouYuGWyjot9z3Ussnl18L6qNmWCJm9Vi2mOZuJj6aKVpqaVEuj5oLLpCCCy6ohQtuyX7//pi4mRmWYZkNOO/Xy9eLWe8DM3M8c1/XdS47RVEUhBBCCCEsyN7aAQghhBCi8pECRAghhBAWJwWIEEIIISxOChAhhBBCWJwUIEIIIYSwOClAhBBCCGFxVilA1q9fj7e3N/b29hw6dEjvtsjISDw8PPD09GTHjh3WCE8IYUWSH4SoHDTWOGizZs3YuHEjr7zyit718fHxREVFERcXx5UrVwgKCiIhIQEHBwdrhCmEsALJD0JUDlY5A9K0aVM8PT3zXb9582YGDhyIs7MzjRo1wsPDg4MHD1ohQiGEtUh+EKJysKk5IMnJydSvX1+97ObmRnJyshUjEkLYCskPQlQsZhuCCQoK4tq1a/muj4iIIDQ0tMDHFNQV3s7OrsD7Ll68mMWLFwNw6tQpnn76aRIvXicjowxBF8HOHpo+fCP3EjyopnejXd265jmwqNSUa9dAycm7oso9QPs5+ePmI3o3mUpVzT2uX79u+ifWUd7zg509eLnXNs2TCWEBSmoioE0YZckdpswPZitAdu7cWeLHuLm5cenSJfXy5cuXefzxxwu879ixYxk7diwALVq0yDdZzdSaj9rOg7QsAJwdb3N01RD1trTZ4bjMCDfr8UXlZPjeaj58DemZNQHwdNFwdHkPkx+zRYsWJn9OQ+U9PzTvt4T7OsWMixMcXT/apMcQoiy6Ba4iLSOvyrjhnYXiqv2iXJbcYcr8YJVJqIUJCQlh8ODBvPHGG1y5coXExERatWpl7bAA9F4sr/5fWzESUZmlpynE/7eXtcOwClvKD4bFRvN+S2gaugSQYkTYhuQaGcxc97d62QkH+laxrdxhlQLku+++Y8KECaSkpBAcHIy/vz87duzA29ub/v374+XlhUaj4fPPP5cZ7kJUMuUxP+gWHLmFiBCW1rzfEtL+OTPn4gwDqky0bkBGWKUA6du3L3379i3wtunTpzN9+nQLRySEsBWSH4QonbQM+GNz+Tn7ZlNDMOWGnR3cClcvOnjtsV4sokJz8Nqj917D7imrxSKEEKYkBUhp2NeEh8LVi3YOgVYLRVRsdg5Zeu817LdaLRYhhDAlKUBKQWMPL0b/qV52ufUmMi1VmMPIW2+SpvNe09hU5x4hhCg9KUBKYbBPLdKO/KFe3nAyBUZaLx5RccVuTeHFZnnvtRY+tawYjSgOFyf9iaiyKkaYS49Re0j7pz0EgJOTFYMpBSlASmHq1PZ6l9eFnrFSJKKiy87MITw80NphiBIwLDZkVYwwl0s3TjFjXrx62Smn4MZ8tkoKECGEEKIcyrGDQU//x9phlJqMKAshhBDC4qQAEUIIIYTFyRCMSdiRNjs875KzC85hU60WjSjf0j+Yi5Ke9s+l+kXeVwghyispQEzAztlZb8Mw3WJEiJJS0tPU95Od7DskhKigZAhGCCGEEBYnZ0CEEEJUGHPn5vXGcHHR5GubIGyHFCBCCCEqjLS0LLV3Tnh4jFVjEUWTIRghhBBCWJycATEF2R1XmJDeDriy+60QooKSAsQE7Oxrse6h/1Mva56+xQtWjEeUb9FPe5P1z/vJzv6SlaMRQtiS5qO2k547x8VJsXI0ZSMFiAnUd6nJwkE31cu3snx4Yb0VAxLl2rvv+PCQRvt+qu9S08rRiLKSzemEKaXdTObES6u1F+xcgDFWjacsjBYge/fuxd/fn6pVq/L1119z5MgRJk2aRIMGDSwRX7mwfbn+LOumsjmdKIOMDIhZH2jtMIpF8oNxsjmdMC0FTWiMtYMwCaOTUMeNG4erqyvHjh1j3rx5NGjQgOHDh1siNiGEjZP8IIQoLaMFiEajwc7Ojs2bNzNp0iQmTZrE3bt3LRGbEMLGSX4QtmZPIoRHa//tSbR2NKIoRodgqlevTmRkJF9//TW//fYb2dnZZGZmWiI2IYSNk/wgbE1WDoSHaH8OXGvdWETRjJ4BWbduHc7OzixdupS6deuSnJxMWFiYJWITQtg4yQ9CiNIq8gxIdnY2Q4cOZefOnep1Tz75pIzxCiEkPwghyqTIMyAODg64urpy+/ZtS8UjhCgnJD8IIcrC6BwQFxcXfHx86Nq1K1WrVlWvX7BggVkDE0LYPskPQojSMlqABAcHExwcbIlYhBDljOQHIURpGS1ARowYwYMHD7h48SKenp4mOWhYWBhbtmzBycmJxo0bs3z5cmrVqgVAZGQkS5cuxcHBgQULFtC9e3eTHFMIYXqSH4QQpWV0FcyWLVvw9/enR48eAMTGxhISElKmg3bt2pWTJ09y/PhxmjRpQmRkJADx8fFERUURFxfH9u3b+de//kV2dnaZjiWEMB/JD0KI0jJ6BiQ8PJyDBw8SGBgIgL+/P+fPny/TQbt166b+3Lp1azZs2ADA5s2bGThwIM7OzjRq1AgPDw8OHjxImzZtynQ8y7MjbXa49idnF5zDplo1GmHb0j+Yi5KepnNNfavFUlKSH4QQpWW0ANFoNNSsqb8hlp2dnckCWLZsGQMGDAAgOTmZ1q1bq7e5ubmRnJxssmNZip2zMy4zwgHUQkSIwijpaer7BcCu/9fWC6aEJD8IIUrLaAHSrFkzvvnmG7Kzs0lMTGTBggW0bdvW6BMHBQVx7dq1fNdHREQQGhqq/qzRaBgyZAgAipJ/a+HCktnixYtZvHgxACkpKUbjEUKYnuQHIURpGS1APv30UyIiInB2dmbw4MF0796dt99+2+gT6zYnKsjKlSvZunUru3btUpOIm5sbly5dUu9z+fJlHn/88QIfP3bsWMaOHQtAixYtjMYjhDA9yQ9CiNIyOgn1+++/JyIigt9//53ff/+d2bNnEx0dXaaDbt++nffff5/o6GhcXV3V60NCQoiKiiI9PZ3z58+TmJhIq1atynQsIYT5SH4QQpSW0TMgkZGR9OvXz+h1JTF+/HjS09Pp2rUroJ1o9uWXX+Lt7U3//v3x8vJCo9Hw+eef4+DgUOrjCCHMS/KDEOY1d+4e0tKy1MsOjkbPG5QbhRYg27Zt44cffiA5OZmJEyeq19+5cweNxmjdUqQzZ84Uetv06dOZPn16mZ7f6uzs4FY4AA5ee6wbi7B5Dl571PcLAHZPWS2W4pL8IIRlfHMilSYDWqqX/RvetWI0plVopnj88cdp0aIF0dHRBAQEqNdXr16d+fPnWyS4csu+JjwUDoCdQ6BVQxG2z84hS32/AGC/1WqxFJfkByEsIyvrNhs6LFQvZ/6yB+htvYBMqNACxM/PDz8/PwYPHkxWVpZJOx1WdHbYsY7jAGie8uMFK8cjbNvmp/zI+uf9Atr3j62T/FB6zo7g9c9Sa2cXe46uGmzliIQtU9LSSPs877Kdc5D1gjExo+dKt2/fzltvvUVGRgbnz58nNjaWmTNnlnmiWUVW36UmCwfdBOBWlg8vrLdyQMKmvfuODw9pbqqX67vULOLetkXyQ8kdePa62niu+e/lp+mcsBZFr09QRVKqTqhJSUlmDqt82768vfpz09Czes3IpDOqAP3upxkZ9YlZH2jdgEpJ8kPJ6X3+Q5dYLxAhrKxUnVBF8el2RQXpjCq0dLuflqfOp4YkPwghSstsnVCFEBWf5Adha1K4r87BS+G+laMRRTG6oPjTTz8lLi4OZ2dnBg0aRI0aNfj4448tEZsQwsZJfhC2RkFhAL4MwBeF/O37he0wegbE1dWViIgIIiIiLBGPEKIckfwghCgtowXIoUOHmDNnDklJSWRl5XVjO378eBGPEkJUBpIfhBClZbQAGTJkCB988AE+Pj7Y21ecFrDWYufsok5ElRUxlYfuqhfQvvYVgeQHIURpGS1AHn30UUJCQiwRS6WgW3DIipjKQ3fVS0Ui+UEI0+sxKm//FycnKwdjRkYLkFmzZjF69Giee+45nJ2d1ev/7//+z6yBVRg6+8JoL7tALTnrUemlzgXlnzMi5WDvl8JIfhDC9C6l3WbmWm2TOk38CStHYz5GC5Dly5dz6tQpMjMz1VOsdnZ2kmCKS2dfGEC/GBGVl5KW974oB3u/FEbygxCmp+SkMuDW7wBknooDLysHZCZGC5Bjx45x4kTFrcCEEKUn+aFsXJwUvAblFaDOLhqOLu9hxYiELdDd/6Ui7f1iyGgB0rp1a+Lj4/HyqqAlmJnpbkwH4OT8EH2tGI+wDd85P0TGP++L8rD5XGEkP5TN78O+zRuKA3yWDbNiNMJ2VNz9X3QZLUD27NnDypUradSoEc7OziiKgp2dnSyzKybdjekAUh1r03eV9mfdFTG5l2VVTMWhu/LFcNXLrFdrUytT+74oT5vPGZL8UDaakO36VyyTvWFE5VGs3XBF6eluTAd523AD+YoNWRVTsRS18iUjTSHmv4EWjcccJD8IIUrLaAHSoEEDS8QhhCiHJD8IIUrLaAEiLEealJVvFbXZmBDCvObOzev7AeDgWDma+kkBYkOkSVn5VlGbjQkhzGvniSzaDwhUL/s0ume9YCxIChBL021MJk3JKg/dxmNQrpuPCSFM62rOfZqG5E3c9on/CehlvYAsxGgBsnHjRqZMmcJff/2FoijqLPc7d+5YIr4yWbFiBbVr16ZXr7K/kDk5OWqjpVWrVrFkyRImT55c8ufWbUxWRFMyWSFTPhS10kWPbuMxKNfNx3RJftAyWX4QlZJu4zGo2M3HdBktQCZPnsyWLVto2rSpJeIxq71797JlyxauXbvGjBkz2LRpEx07dqRVq1b079+fNWvWMGPGDNLT08nOzubjjz+mY8eO9O7dmy5dutCqVSsAhg8fTk5OjlljlRUy5UNlH3aR/GDa/ODkBIGDYgBwcdHkW0UnKibdxmNQsZuP6TJagDz22GMVIrkAODk5kZGRgaurKxs3buSVV15h6tSpPHjwgE6dOrFz506SkpJo2rQpFy5cIDk5GVdXV6ZONd2ZB93GZNKUrPLQbTwG5bv5mC7JD6bNDzNn/UGGfTwAsyd7AVKAVA6Vo/GYIaMFSIsWLRgwYAB9+vQp95tNvf/++6xdu5Z9+/bxyy+/ULNmTapVq8b8+fNZvXo1v/32G+3atWPixInqY2rWNG2TqHr2Vfkj2heABO7St5gbicoKGdtQ2pUu63e2oAm+6uV69jeLuHf5IfnBtPnhBa+P1J/ftZOmZKJiM1qA3LlzB1dXV3788Uf1uvK02dQXX3zB1q1bcXd3p2PHjrzzzjvcv3+fhx56CIAXXniBtWvXUr16dbp3786rr75KWFgYqampfPrppwU+59atW1m5ciVVqlShWrVqBAYGFjueIB8NaUdiADhx5iqEtCvW42SFjG0o7ZDLif9epIlHpno5yKdizP+W/JBfWfKDqDx0l95WlmW3huwURVGsHURZtWjRgkOHDpX4cfHx8cyaNYt58+ZZpaGSV/+vif/v0BI/rqBv4XJGxHwMJ5qW5m9d2tfaGkr7ebJV5TU/NA1dwh+bR1v8uOWd16CtxK/tle9nWxM0JEZdenv61lbWjrDNOA2ZMj8U+jVs3rx5TJ48mQkTJmBnl3+8esGCBaU+6Ntvv83mzZuxt7enTp06rFixgscffxyAyMhIli5dioODAwsWLKB79+6lPo4xXl5erFu3zmzPb5Tuklwo9rJcw/8A0z+YK8MzJlRQgVeq8VndpbcVbNmt5AchykZ36W1lWXZrqNACJHdiWYsWLUx+0LCwMN577z1Am6jeffddvvzyS+Lj44mKiiIuLo4rV64QFBREQkICDg4OJo/BJuguyYUil+UWRbfg0C1GQAqS4jI8y2GSCWG6S28ryLLbXJIfhCgb3aW3lWXZraFCC5DevXsDMGLECJMftEaNGurP9+/fV79Bbd68mYEDB+Ls7EyjRo3w8PDg4MGDtGnTxuQxVFRydqR4THaWo5KS/CBE2eguva0sy27zUaxk2rRpipubm+Lt7a389ddfiqIoymuvvaasXr1avc9LL72krF+/vsDHL1q0SAkICFACAgIUQP7JP/lnon8BAQEWyQFFkfwg/+Sfbf4zZX4w21T8oKAgrl27lu/6iIgIQkNDiYiIICIigsjISD777DNmzZqFUsB82ILGlwHGjh3L2LFjgfI1aa7HqLyZz6mOV4j9JCHvRgu3Zjc8C6CrvJ4tsanfSWcOiP+kJtTKfFy9yZabTJljWMWQ5Aet5v2WkJaRd9nFCY6ul4mnZaE78fTF6D/ZEPKYelt4eAzh4YFWiqxwzUdtJ11nMzrHrGtMW/U3AE440LfKa9YKLR9T5gezFSA7d+4s1v0GDx5McHAws2bNws3NjUuXLqm3Xb58WZ18VlHo/qfjNWirSeaAlFZR/xkbziUpjCX+Uy+qqCgoHpsZStGZA5KRuZWYtYFWDceWSH7QMiw2mvdbQtNQbf8PKUbKzskxh/DovMt7Eq0XiyHdL6M52ckcH7I670Y7FzRVtlspMsspcQHyxRdf8Mgjj/DCCy+g0ZSufklMTOSpp7SrAqKjo3n66acBCAkJYfDgwbzxxhtcuXKFxMREtb1xRaTbFRVsqzNqcYuK4hYqZWFTRUUJ6HY/rSidT42R/FA2ugVHbiEiSi+0ZwoDqKdeDlxrvVh0+34AXMm8yvRlKQA4nD2PplmMlSKznhJnCEVR2LNnD2vWrCE6Otr4AwowdepUTp8+jb29PQ0aNODLL78EwNvbm/79++Pl5YVGo+Hzzz+v0DPcH6UqA3S6Y66z31Pudsotj8M0ZmOw421GlTrq67uQitH51BjJD6bj7KjtH6NedrHn6KrBVoxIlMXOE1lq3w8A78bfE/rRaQDsnB+DZlYKzIpKXIC89lrZx6K+/fbbQm+bPn0606dPL/MxwHy7XYaFhXHz5k1SU1NZsWIF1atXL9Vzurho1I2nAFJdGjFg+b+0Fyw8HCNMwGDH28hR21iYFgNoX+vKQPKD6fLDgWev6w09Nv+9fpnjrIgMzyw4ONrmZy0rB8J1tt5Ye2oXLsP+Y72AbIDRVyo9PZ1vv/2WpKQksrLyXuSZM2eaNTBzMNVulx988AEA8+fPJzY2lg4dOpQqHsNJiF6DKlaviMouIy27ws/7kPxgvvyQ7+yiDMkUKC0tS29i6cnoP60XTBFSuK8/5J6TbcVobIPRAiQ0NJSaNWsSEBCgt9lUeWTK3S6vXbvGoUOH9DamEqKykfwg+cHmZF2EWwu1P1epAy6+Rd/fjHTPztg73GXArY3qbZW1+ZguowXI5cuX2b69YszGNdVul8nJyUydOtXkY9C6k1KdXB6lbynatAsL05n38Z3Lo+qkU6gcE08lP1guP4hiUnLUoVCnvz/RO+tw0zGN8PAYQDssOnWqeZfC68778HlmLzw0UL0tOz4cx/KxZ6PZGC1A2rZty4kTJ/Dx8bFEPCZnjt0uQ0JCaNKkCVOnTmXcuHH4+fmZJFbdSanrqgBVdCp3mRNim3SX2nJcb1JxZZh4KvkhP3PlB91JqTIhtXj6pt8C17zP5PoXHyP8n74guYWIOenO+1ib8B1ps39Sb7NzdjH78W2d0QJkz549rFixgkaNGuHs7IyiKNjZ2XH8+HFjD7W6kSNHMnLkyEJvj4+PZ/78+cybNw8AjUbDkiX646wbNmzI97jDhw+bNM5cupNSU10eMGC59U4dipKLHJXMwrS8oqMyTDyV/GC5/KA7KVUmpObZk4herw8nx6xC76vbF8QcPUEMJ8TedEzLO6udlVku2wmYk9EMuW3bNkvEYRW2tttlviZluuxcyt0S3QrJYKktdnnfYirDpFNDkh8sR29SqkxIVRmuLln34DDQrsD76vYFMUdPEMOlts2e+Z8670PmfORntABp0KABx44dY/fu3QB06NDBZKcURQnoFhwyHGM9BkttKzvJD6K80jhq9IZhSjsnRO+sh4PGoBj6HaqEAzLnoyBGC5BPPvmEr776iv/7P+1fbujQoYwdO5YJEyaYPbjKLF+XVBzoi7cVIxIF+Y44MtAup6sMk04NSX6wDicn9HoI2fLeQlale+YY9FbFtH+xvV6xUNo5IbrLgAdE/6WXtzUXzpO2QXt8mfORn9ECZOnSpRw4cICqVasCMGXKFNq0aSMJxszq2Vflj+i8OSAJjlfp2/OfC4YfKhmSMS/dYRc7/SSyftvDNMnUntKtZ1/xJ50akvxgHW+/E0d247/Vy++OeASonAVIvv4aupsWGuRF3VUxZx0fIzw6b6O6/Rc0xVohYzjPY/8FjTqvxMExW28ietbpo2hmxJTm16oUjBYgiqLoLSVzcHAocFdKYVpBPhrSjsSol0+cS4Oe/+xpYFhsyJCMeRUx7HJi7TGauGvbKQf5VPxJp4YkP1hH722PoaTnFbzh6VWtGI1lGRYA9o4Oev/pk74RXAt+rN6qmJ7HGYDOTrnknRGZO3dPoWdEXFw0eo3PwqPz5qCsS/sObv2l3qZkV76cUBJG/zqjRo3i2WefpW9f7TZpmzZt4uWXXzZ7YJWdYfX9X+mSapOyM7NscntvS5H8YB2VuUuquTqfujjmraZx8WtPeM+C7zd3m/6qm7OOf7IObQxOOZl6X1Zk3kfRjBYgb7zxBoGBgezZswdFUVi+fDnNmze3RGxCh16TMsP5ILJCxrSKWOmiO+cDKue8D12SH2yDg6O9RRtsVQROOOgN3Xj2zMurhkWGLhdHg1U3/Kmegcn65Q3Sjt5Sb5N5H0UrtAC5c+cONWrU4ObNmzRs2JCGDRuqt928eZOHH37YEvGJf+g1KcOgx4KskDGtIoZcMsiudM3GCiL5wbb4BdfhZI0qACRsSGBqBZ4Pkq/vB+fyWq9DvnlaenS+rPU1+LKmm1enFnL2A7RfQtbpfAlxerAb0rRLbZV0pNdHCRRagAwePJitW7cSEBCAnV3et7zcRkPnzp2zSIBCq9hNymSCaukUMdFUV2VsNlYQyQ+2Zdn1kyjJ/zQp+7tiNynL3/cjb6krQPoHc1HStZftnF30h6uK+LJmeEakME4YzDlJ26h+YZEhl5IpNHtu3aqdc3D+/HmLBSMKp7vEznvQ94Uv0ZUJqqVTxFkP3WGXzLScStdsrCCSH2yLNCnLo6SnqWch0maHF/txxW5zkDoXlPXqxazd+8k6rj2ODLmUjNGvb8899xy7du0yep2wHN3hGChgSEaXzA8pWBHzPAzpDrtU1iGXwkh+sD0VvUdIkctuS6K0Z4sNvqxkHQ+XYZdSKrQASUtL4++//+b69evcunVLXVp3584drly5YrEARX4ae/0x0ATHRxlQ2Jil7gcqdW7lHp4xHGYpZkfTzdse5Y9M7c8ae/OEVt5IfrBdM2f9QYZ9vHp59mQvynuPEL1t7Uuw7NbO2UXvLIjekExJzhYXc4hWlEyhBciiRYv4+OOPuXLlCgEBAWqCqVGjBq+99prFAhT5GfYIOXkujXU9UwAjHVMr+/BMMduoG650Obn2SqXu9VEQyQ+26wWvj/Quv2tX/odkdJfelmTZreFy5SKHZAzPiBje9k/u0J1jAjLsUhaFZtNJkyYxadIkPv30U+lqaGMMl9jFDIopfIVMUSr68EwJhll05VvpknmzUvf6KIjkh/JDd0imvA7H6K58cXLMMc9Bipn/dOeYiLIx+nVuwoQJnDx5kvj4eNLS8pL58OHDzRqYKD7dIZkih2MMVcThmVIOs+jSHXIBGXYpiuQH26fbtr28tmzXXfli2G1UdxIoFH1GQndIJt8KmSJoz3qkGX1+UTJGC5BZs2YRExNDfHw8zz//PNu2baN9+/aSYGyI7pCM7nAMlGATO8NiQ7cgseVipKCzHKUoOnSHXXSHXECGXYoi+cH26bZtLy8t2w3brWscdT6DSlapJ4HqFhwlWSEjZz3Mw2hm3bBhA8eOHaN58+YsX76cP//8k9GjR1siNlFMukMyusMxUMIhGV1FnR3RZYnixLDIMDx+KQoOQ3orXWTIpdgkP9g+3f90nfotKRcrZAzbrevuMlvqVS8GDCeoGruvMD2jBUiVKlWwt7dHo9Fw584d6tSpU+maDH3//fcsWrSI6OhofvjhB44ePcrt27d5//339Zow2QLDFTJnHR+DnoW0cC+uogqMoooTUzFRkaHLcKLp1m2PyUqXUpD8oHXo0CF27NhBlSpVeOONN6wdTqHKyy66ht1O9XaZLWLVS0kUd/hFmI/RAqRFixakpqYyZswYAgICqFatGq1atbJEbDbh6NGjpKWl4e7uDsBPP/3E/PnzWbVqFceOHcPf35+///6biIgIIiIimDhxIuPHj+err77CycmJ4OBgsrOzOXLkCJcuXeLDDz+kXbt2TJgwgZs3b5KZmYmHhwfBwcEmiddwhczVcxoG9NQmmFKfDSmKrQ7NGGE40fSrtXto7P4HIEMuJVHZ80OuZcuW0bhxYxRFUbvB2lJeyGW4i+57OdVtZg8ZvWEXB43BfispQD2rxCXMx2im/eKLLwB49dVX6dGjB3fu3MHXt5A24DYmY/UKci4kFfv+9g0a4jRspN5127Ztw9XVlaNHj3L06FG923LPfri6ulKjRg3WrVtH586diY6Opnbt2tSuXZszZ87wyCOPULVqVeLi4gDw9fVl4MCBHDhwgN27d3P79u0y/Z66DBOI7ulWw1bDpT4jUk7pnvVwwkHvtqxKvqttaVX2/JDr6tWrzJs3jzVr1nDkyBECAgJsKi/kMvzW779iCycfbgpAwrrfsebXCd1hl3wbwaXtgQfa/VakD0fFYbQACQ0NZcCAAYSGhuptOFUeFJYsSmLatGkAJCUl0bx5c65evcqcOXO4ffs2w4YNU+/38ssv07NnTw4cOEBCQgKLFy+mevXqBAQEsGPHDurXr092tvY/P2dnZ0C7aVeVKlU4efJkmeMsjO6QjIujt94mS2Y5I2LDdMO/B2cAACAASURBVM96zN0G4bLSpcwqe37I9eqrr6p5YciQIer1tpoXcq14+HNQtH1D/NNHqGdDwPJnRHSHXXS3uAf9be6lD0fFYacoRc/o+fXXX1m3bh3ff/89rVq1YsCAAfTq1QsXl7K/6B9++CFhYWGkpKRQu3ZtACIjI1m6dCkODg4sWLCA7t27G32eFi1acOjQoTLHUxZffvkltWrVYuDAgVaNw5Duac095zTsXJWXUAzPCFS0syGG8zx0f8eg4Xto7543y97ap59tSUk+T5IfimareaEgQ1dsIe3hvOGzs98dJbSB9nU0x+fDcKWLbn5ax3H9bqe3wtUCJG22tD63JlN+noyeAenUqROdOnUiOzubn3/+ma+++oqXXnqJO3fulOnAly5d4qeffuLJJ59Ur4uPjycqKoq4uDiuXLlCUFAQCQkJODg4FPFMtuHVV1+1dggF0k0aQUNi9E5t6p4R0W4xXf6HZwyLKt0kpnfWI1uGXExB8kPRbDUvFET3bAjAsM6ToVYgADvXxZhkeEa36HBx0eh9BvMNu4gKr1iz7R48eMCWLVtYt24dR44cYcSIEWU+8Ouvv868efMIDQ1Vr9u8eTMDBw7E2dmZRo0a4eHhwcGDB2nTpk2ZjyfyT1DdeU7D1H8mqBoWG7oFiS0XIwWd5dD75qRj59q8sx4y0dR0JD9UDJqQ7XqXX4yfRJaXGwAbN6SZZHhGd57H3G2GK/byhl2cHuzWbnOfS+Z9VEhGs/CAAQM4cOAAPXr04LXXXiMwMBB7+7INmEdHR/PEE0/g5+end31ycjKtW7dWL7u5uZGcnFymY4k8hglD94yIiyN680N0Cw5bOztS1FkOXXO3QZrOPA8562F6kh8qrtAzpyFxIgCbgvLOhgDEfLeHNJ2CpLj2X9Do5Rz9lS5/5n2W0zbqLb3Xnfchcz4qDqMFyKhRo/jmm29KfJozKCiIa9eu5bs+IiKCOXPm8OOPP+a7raDpKIX12Vi8eDGLFy8GICUlpcD7mMKWLVv43//+R2pqKm+++SanT5+26T4gJaF7RiTmgoa0zLwCRbcgsfbZkZKc5dAtOvZ/t4fWDfLGmOWsh+lV9vyQkJDAnDlz6NOnD3369GHr1q2cPHmShIQEFi1ahKOjo9mObW66Z0ReOPk62Y3PqJfP0pLS9A8JfBaDYV/9z3VhpBNpxWQ0I3fs2JHIyEguXrzI4sWLSUxM5PTp0/Tq1avIx+3cubPA60+cOMH58+fVbzeXL1/mmWee4eDBg7i5uXHp0iX1vpcvX+bxxx8v8HnGjh3L2LFjAe2kGHNxcXHh6tWrpKenU6dOHT777DOb7gNSErpnRMLDY/S+jRQ1HlvU2ZGi6BYrhkWFsccVVnAYSsvM+1YVfkTOeJhbZc8PTZo0YeTIkaSmpgLQq1cvevXqxfjx48nIyMDR0bHc5gddhv1DNC/+l6ynC27rXtSXEt2iI9/nOnUuKOu1P8uQS6VQrDMgAQEB7Nu3D9Ce9uzXr5/RBFMYHx8f/vorbyOhhg0bcujQIWrXrk1ISAiDBw/mjTfe4MqVKyQmJpatqdHdFZCVVPz7axpC9ZF6V8XHx/PZZ5/x22+/8euvv+rdZot9QEpLOyEsRr28/4KG8H++4RgOz+gqydkPwzMnxS0qimI4zLL/uz2EH8mb5CbMq7Lnh4J88sknhISEULWq9j/oipAfDPuHhEb3gNPa4RnsXPTOlhT1pUTvc69bcPzzPLLUtnIxmqHPnj3LunXrWLt2LaBtvWxk5W6peXt7079/f7y8vNBoNHz++edlm+FejGRhTN26dQkPD+fOnTtMnz4de3v7ctUHpLgM54fMnbunWMMzJWGqoZqihlkCPTVMnRpokuMI4yp7frh27RobNmzgwYMHNG/enPXr13P48GEURaFly5Y89NBDQPnPD4Z0C47MjUF62zH0NfponbMchWyxIEMulYPRPiBt27Zl165dtGvXjiNHjnD27FkGDRrEwYMHLRWjUdIHxLwM1+vvv6Chdd+Cx39LW5zoHc9w8qgB3aJD+neYXkk+T5Ifiqci5wfdreoNlWXLe9mrxTaZ8vNktAD56aefmD17NvHx8XTr1o29e/eyYsUKAgMDTRKAKdhCgqlMDNfy6xYAxoqH4jAsYgwLICk6zKsknyfJD6IoRRUnhqToKB8sWoAA3Lhxg/3796MoCq1bt1a7EtoKSTDWY1gc6CpJoWCq5xFlV9LPk+QHISoPi3RCPXLkiN7levW0OxFevHiRixcv8swzz5gkAFG+FVUYzJ27R29ia1EMuyIK2yb5QQhRVoUWIG+++WahD7Kzs+Pnn382S0C2xtg6/2+++Ybr169z//59Zs6cae1wbYqctai4JD9oGeaHRYsWcfPmTc6cOcN//vMfatasae0QhbBZhRYgv/zyiyXjsFnG1vnHxsYyf/583n33XVJTU6lVqxbJyclERUUxYcIEpk2bRt++fdm2bRsPHjxg3LhxJCQkcPr0ae7fv8/QoUOZNGkSY8aM4dixY9SsWRNfX186duxo5d9ciMJJftAyzA9VqlTh0qVLZGZmUqNGDYAK0QdECHMotACZN28ekydPBmD9+vX069dPvW3atGnMmTPH/NGV0S+c5S/uF/v+dahKZxobvZ/hOn/Q78j4xBNPcP36dZYvX86QIUNYtmwZfn5+PHjwgDNnznD79m1q1qzJvn37GDp0KG3btqVXr144OTlx4sQJ7t69W7JfVAgLk/xQsMuXL/PFF1+wcOFC4uLiaNasWYXoAyKEORRagERFRakJJjIyUi/BbN++vVwkmOIUE8YYW+fv7+/PRx9pd5CsVauW+rghQ4YwdepUXnnlFfr378+2bdtwc3OjSZMmrFixAk9Pz3zr/lNTU3F2dubEiRPyjUfYNMkPWob5oUqVKoSHh5OcnMzQoUPV+1W0PiBCmEKhq2CaN2/O0aNH8/1c0GVrs7VZ7tnZ2cyePZvOnTvLUIood4rzeZL8UDIVuQ+IqFwssgpGd0jBcMOn8rwBmyU4ODjwzjvvWDsMIcxG8kPJvPrqq9YOQQibU2gBcuzYMWrUqIGiKDx48ECdUKUoCmlpxWssI4SomCQ/CCHKqtACJHc8UgghDEl+EEKUlb21A7B1CQkJjBw5kk2bNgFw6NAhIiIi+M9//gPAypUr+eijj3j33XetGWalEhMTQ4MGDejcuTNdu3blxo0bJXrsjBkzzBidaeW2NJ87dy7JycnExsaqTcBiY2NZunSpFaMTu3fv5v3332fEiBHcvKndrn7lypWMGjXKypHZjtzPa2BgIIGBgaVa1WOKz+2KFSvIyckp03MYoxtn+/bF74P0//7f/yM7O5uYmBjOnTsHaONdsmRJoY8JDw8nJiamTPEWxNhxTUkKECNy1/nnWrZsGS4uLiiKgqIoxMbGqk2ZcnsBJCcn89FHH5GRkcFbb73F3r17mTFjBm+++SZnzpzhhx9+YP78+cyePZukpCRCQ0PZunUrERERfPbZZ/z222/W+FXLlWHDhvHLL78wYsQIdSdWczF30iqOqVOn8sQTT+gVIP7+/rz88stWjqxy69ChA1OmTMHDw4PU1FR++uknGjZsqNeA7O+//2b69OkATJw4kYSEBMLCwpg+fTr79u1j9+7dfPLJJ7z11lsAtGvXjqioKL744gs++eQTvv/+e6v8bqY0bNgwYmJiiImJ0fvbWPKzVdwCxBqf948//hgHBwe9AqQ0bCFXlUSFLkBW7IXw6OL/W7HX+HNevXqVV155hWrVqum1oy6qD0hUVBQNGzakYcOGen1Ajh07BqD2AWnZsiXp6enSB6QEcos+gAMHDhAYGEi7du1Yvnw5AHv37qVdu3Z07tyZdevWAXDixAl69+5Nu3btuHfvHleuXKFz5860b9+ef/3rX4D2m0xISAi9e/dmx44dfPnll7Ru3ZopU6aoZyUKOl6unJwcRo8eTadOnejZU7uz3s6dO2ndujWtW7dm586dgPYMx7///W9atmypns3YunUrAQEBvPLKK2RlaffHGTlyJGfOnGHx4sV88MEHDBkyRO/b1vvvv0+7du3o0qULFy9eBMDPz4/hw4fj5+dHbGwsN2/eJDAwkM6dOzNx4kSTvxbljanywzfffIO7uzvu7u7s2rWLo0ePcvToUc6ePQtQaB8Qd3d3zpw5w507dwrsAxIQEEBmZmaF7AMSGBjI5MmTGT58OA8ePGDQoEF06dKFAQMGkJmZyXfffUerVq3o0qULP/zwA5D/cxsWFsbJkyf58ccfad68OQAjRozgr7/+IjIykk6dOvHss89y9OhRDh48SGxsLM899xyrV6/mzJkzdOvWjU6dOjF79mxA+xkbP348PXr0UOO8c+cOAwYMALSF4aZNmzh37hz/+te/UBSFcePG0aVLF4KDg7l161aRv3NR8QYGBpKRkcGKFSt488031S+1u3bt4vnnn+f5559HURRu3rxJ586def7554mNjQXy56pJkybRqVMnOnTowMWLF4mLi1OL29wvMT///DPz5s0jNjaWTp060bp1a6ssnS90DkhFMLJd2Z/DcJ3/q6++ypw5c7h9+zZDhgyRPiBWsnr1arZs2UJOTg579uwBYObMmURHR1O9enW6du2qvgabN2+mdu3a5OTkqGeXtmzZQkREBLt27aJnz5789NNPaDQahg4dSmJiIgAZGRls376drKws2rdvz969ezl06BAHDhwo9HhOTk4AbN68mTp16rBkyRL1W0l4eDg//vgjAD169CAoKAiA/v37Ex4eTteuXXn55ZeJjIzk119/5datW3Tu3Fnv9x47dixZWVmMHj1aPf167do1fv75Z/bu3cuePXuIjIxk4cKFXLt2jQMHDnD48GFWrlxJcHAwgYGBhIeHU4w9KCs8U+SH9evXs2rVKnr06MGFCxeYO3cuAElJSTRunNdnpLL3AVm9ejV79uyhUaNGarHet29f2rRpw6effkpISAiDBg1i4cKFbNiwge+//57//ve/NGzYEEVR+PXXXwH9z23btm3Zu3cvV69epV69ety9e5c///yTOnXqMGnSJP79739z5swZ3nnnHdasWYO/vz87d+5Eo9EwYMAAli5dSv369Rk0aBCXL18GtEXGZ599psZdo0YNbt++TXp6Oo888gj79u3j3r17tG3blq1bt/Lkk0+ycOFCtm3bxpdffkmbNm0K/RsUFS+Avb09I0eOpH379gQFBbFixQrq1KnD2rVrGTNmDMePH2fHjh2MHj2aIUOG0L17d/W5c3MVQKdOnXB1dWXnzp0sWrRI3a06KSmJZs2asXfvXm7dukVgYCCenp7ExMRgZ2dH586def311037whtRoQsQU6hbt67eG7JBgwZ6L/yIESPyPSY7O5tvv/1WbdTUoUMHOnTooN6eO19kyJAhgHb8D7T/EYniGTZsGLNmzWL06NFcvHgRHx8fjh07RkhICADXr18nJSUFQN2d1d5ee8KvWbNmgPbbQGpqKjdu3GDcuHGkpqaSlJTElStXANQN1a5fv86TTz6Jg4MD/v7+agwFHe+JJ54AtHOH2rZtq3dcOzs7dbWIg4OD+jzNmjXD0dFRvZ+9vT3VqlWjWrVqPProo0b/FklJSfj6+gLaNfqzZs0CwMPDAxcXF/X37NSpE7/99huDBw+mZ8+eDBs2rLh/blGIfv366TVhy/Xxxx/rXd6wYQNvvvkm9vb2PP300+ocMkB9Tw0ePFjvsblnziqCYcOGqWcacgUEBADwxx9/cPjwYRYtWkRaWhqDBg1ixowZzJ49m6ysLHX4yvBz27NnT8LCwsjJyWHIkCFs3ryZxx57DNAWPGvWrMHe3r7AZeGnT59W3/+pqakkJyfrxaTL3d2d9evX8/zzz7N//3727dtHWFgY69evJyoqih07dpCVlVVk8QHa4qaweAtj+DufO3eOXr16Aeht+Kj787x589i1axeZmZk0bdoUOzs7nJ2d+fnnnxk/fjzR0dFcv36dt956i3PnzvHmm2/y999/c/r0af76668i4zE1KUDMQPqAWIaDgwNTp04lPDyctWvX0rx5czZs2EDVqlXJzMzE0dEROzs7bty4wSOPPKKeidBNSIqi8M0339CnTx9GjhzJkCFD1LMDuQVB7dq1uXTpEjk5ORw/flx9bEHHy+Xp6cn+/fvp1asXOTk52Nvbk5OTw507dwD9VSSGCTInJ4f79+9z69YttYjK5ejoSHp6ut51DRs2VIfzDh06pH7zNvw9s7Oz1eLX399fChALkj4g+eV+vjw9PXnuued44YUXAMjMzCQrK4slS5awb98+/vOf/zBgwIB87+c6depw9epV3NzcaNeuHX369GHcuHEAfPHFF+ow2JgxYwDtZyc7OxuNRoOnpycff/wx9erVIzs7Gzs7OxYuXKjGpKtt27Z8+OGHrF69mqNHj3Lq1CkaNWqEp6cnw4cPV4dLMjMz2bu38HH8ouLNlRtjLsPfuVGjRhw7dgwvLy+OHj2qfhnOjfvGjRvExMSwe/dufvrpJ9asWQNoC6vPPvuMX375hQ0bNpCeno6LiwsLFy5Uh5Xbt29v8TOjUoCIcs3T05OUlBSuXbvGrFmzCAkJIScnh4cffphvv/2WyMhIevfujbOzM6+++mqB3zi6dOnC8OHD1ZVOhjQaDSNGjKBt27a0adNGLTQKOl6ukJAQtmzZQseOHalWrRo//PADM2fOpFu3biiKUuSqqSlTptCxY0eeeeYZ6tatq3db69atGTlyJCdPnlQTdt26dencuTNt27bFycmJlStXFvi8Bw8eZNq0aWRmZqrDP0JY29ixYxkzZgxffPEFiqIQGRnJxo0b2b9/P/fu3VOHuAtSr149fHx8aNiwISkpKepZx1atWtGxY0e9TtTBwcH06dOH0aNHExERwUsvvUR6ejqOjo56n11DuZsDNmvWjBYtWnD9+nVA+xmfOHEiXbp0AbRnsnPPcJY03lyBgYFMmzaNAwcO8OSTT+Z7/OjRo3nhhRdYtWqVOlSn66GHHqJatWp06dJFPSua+zusX7+emjVr4ubmpp4VDg4OZvz48Xh5eanDx5ZUaCv28sQWWi2Lii0rKwuNRsOBAwdYtmwZixYtsnZIZlPRPk8V7fcRwpos0opdCJHn008/ZdOmTWRkZBR6hkEIIUTxSQFixO7du9m3bx/x8fHMnz+f/fv3c/ToUW7fvs3777/PqlWruH79Ovfv32fmzJnWDleYyeuvv27xGeLC9sXFxbFjxw4SExN57733OHjwoF5+kH1xhChchS5AVqyIJSkp1fgd/9GwYS1GjvTXuy53Bct7772nNhqaP38+q1at4tixY8TGxjJ//nzeffddUlNTqVWrFsnJyURFRTFhwgSmTZtG37592bZtGw8ePGDcuHEkJCRw+vRp7t+/z9ChQ5k0aRJjxozh2LFj1KxZE19fX9lFVwgzM0V+8Pb2JjY2lt27d+Po6JgvP/j7+/P3338TERFBREQEEydOZPz48Xz11Vc4OTkRHBxMdnY2R44c4dKlS3z44YfqnIObN2+SmZmJh4eHLMsXFVKFLkAMk0Vp6TYa0lXYjqCGjciWLVuGn58fDx480GtEtm/fPoYOHao2InNycuLEiRPSiEwICzBVfhgyZAi1atVSG8Dlys0JhTUiq127NmfOnOGRRx4psBHZgQMH2L17d4VsRCYEVPBOqKaQ22goJSWFCxcu0LVrV+bMmUNcXBy+vr5FNiLbsmULzZs3p3///pw7dw4HBweaNGnCH3/8gbOzc5GNyIQQtm/79u3MmzeP6OhoHn300Xz5IdfLL7/Mhx9+SGhoKL169eLPP/8kPT0dX19fTp48iaura4VuRCZEQWQVjBlkZ2cze/ZsOnfuLEMpotyxtc9TWdnC7/Pll19Sq1YtBg4caNU4hCgrWQVj46QRmRBClzQiEyI/qwzBhIeH88QTT+Dv74+/v7+62RBAZGQkHh4eeHp6smPHDmuEJ4SwIskPQlQOVjsD8vrrr6s79OWKj48nKiqKuLg4rly5QlBQEAkJCXr7ZgghKj7JD0JUfDY1CXXz5s0MHDgQZ2dnGjVqhIeHBwcPHrR2WEIIGyD5QYiKxWoFyGeffYavry8vvfQSt27dAiA5OZn69eur93Fzc1N3KRRCVB6SH4So+Mw2BBMUFMS1a9fyXR8REcG4ceN4++23sbOz4+233+bNN99k2bJlBe7EV1gnwcWLF7N48WIATp48SYsWLUz7C5RBSkpKsbZRtwSJpWASS+FOnTpl9mNIfrCN11tiKZjEUjiT5gfFys6fP694e3sriqIoc+bMUebMmaPe1q1bN2Xfvn1GnyMgIMBs8ZWGLcUjsRRMYimcLcUj+cG8JJaCSSyFM2U8VhmCuXr1qvrzd999R7NmzQDt9sZRUVGkp6dz/vx5EhMTadWqlTVCFEJYieQHISoHq6yCmTx5MrGxsdjZ2dGwYUN1a3Nvb2/69++Pl5cXGo2Gzz//XGa4C1HJSH4QonKwSgGyevXqQm+bPn0606dPL9HzjR07tqwhmZQtxSOxFExiKZy145H8YDkSS8EklsKZMp4K0YpdCCGEEOWLTfUBEUIIIUTlYPMFyPr16/H29sbe3l5vA5ykpCSqVKmitmvW3Wvh8OHD+Pj44OHhwcSJE9Xle+np6QwYMAAPDw+effZZkpKSTBYPFN4m2pzx5CpN++rC4jK17du34+npiYeHB3PnzjXLMQw1bNgQHx8f/P391SWYN2/epGvXrjz11FN07dpV7S8Bpm/x/dJLL1GnTh11AmVpj2+K16igWGz5/VISkh+Kx5Zfb8kPlPr45T4/mGw9jZnEx8crp06dUjp16qT8/vvv6vW6y/MMtWzZUtm3b5+Sk5Oj9OjRQ/nhhx8URVGUzz//XHnllVcURVGUtWvXKv379zdZPHFxcYqvr6+SlpamnDt3TnF3d1eysrLMHk+ud955R/nggw/yXV+auEwpKytLcXd3V86ePaukp6crvr6+SlxcnMmPY6hBgwZKSkqK3nVhYWFKZGSkoiiKEhkZqUyePFlRlKL/RqX166+/KocPH9Z7j5bm+KZ4jQqKxVbfLyUl+aF4bPX1lvxQufODzZ8Badq0KZ6ensW+/9WrV7lz5w5t2rTBzs6O4cOHs2nTJkDbynnEiBEAvPjii+zatavEFWNh8RTWJtrc8RhTmrhM6eDBg3h4eODu7o6TkxMDBw5k8+bNJj9Ocej+vUeMGKH3Opi6xXfHjh15+OGHy3R8U71GBcVSGGu/X0pK8kPZWPv1lvxQ+uNXhPxg8wVIUc6fP0/z5s3p1KkTu3fvBrTtmt3c3NT76LZr1m3lrNFoqFmzJjdu3DBJLIW1ibZkPCVpX11UXKZkrfbZdnZ2dOvWjYCAALUj5p9//km9evUAqFevHn/99ZdFYyzp8c39Gtni+8WUJD/os8XXW/JDnsqYH6y2G66uotoyh4aGFviYevXqcfHiRR555BEOHz5Mnz59iIuLK7Jdc1G3lTWewp7bFPEUJ66Stq8u6bFLy1LHMbR3714ef/xx/vrrL7p27crTTz9d6H2tFaOx45szLlt9vxRE8kPR8RQnLlt9vSU/GFeR84NNFCA7d+4s8WOcnZ1xdnYGICAggMaNG5OQkICbmxuXL19W73f58mUef/xxQFuVXbp0CTc3N7Kysrh9+3aBp55KE0/ucxse1xTxlDSuMWPG0KtXr1LHZUqFHd/cco9Rp04d+vbty8GDB3nssce4evUq9erV4+rVq9SpU8eiMZb0+OZ8jR577DH1Z1t6vxRE8oPkB1OT/GA8llzmfL+U2yGYlJQUsrOzATh37hyJiYm4u7tTr149qlevzv79+1EUhVWrVqnfSkJCQli5ciUAGzZsoEuXLiarGAtrE22peEravrqouEypZcuWJCYmcv78eTIyMoiKiiIkJMTkx9F1//597t69q/78448/0qxZM72/98qVK/VeB0u0+C7p8c35Gtnq+8VUJD/os9XXW/JDnkqZH4zPkbWujRs3Kk888YTi5OSk1KlTR+nWrZuiKIqyYcMGxcvLS/H19VWaN2+uREdHq4/5/fffFW9vb8Xd3V157bXXlJycHEVRFOXBgwfKiy++qDRu3Fhp2bKlcvbsWZPFoyiKMnv2bMXd3V1p0qSJ3gxgc8aTa+jQoUqzZs0UHx8fpXfv3sqVK1dKHZepff/998pTTz2luLu7K7NnzzbLMXSdPXtW8fX1VXx9fRUvLy/1mNevX1e6dOmieHh4KF26dFFu3LihPqawv1FpDRw4UKlbt66i0WiUJ554QlmyZEmpjm+K16igWGz5/VISkh+Kx5Zfb8kPlTc/SCdUIYQQQlhcuR2CEUIIIUT5JQWIEEIIISxOChAhhBBCWJwUIEIIIYSwOClAhBBCCGFxUoCYQLVq1fQur1ixgvHjx5foOaKjo026E2RqaipffPFFse5rGH9BwsPD+fDDDwGYOXNmkU2PNm3aRHx8fPECNaNDhw4xceJEix0vKSmJb775xmLHE+WD5Ad9kh9ELilAbEBWVhYhISFMnTrVZM9ZkgRTUu+++y5BQUGF3m4rCaZFixYsWLDAYscrTYLJbZYlRGEkP5iH5AfrkwLEzC5cuMBzzz2Hr68vzz33HBcvXgRg5MiRvPHGG3Tu3JkpU6bofSvy9/dX/1WpUoVff/2Vmzdv0qdPH3x9fWndujXHjx8HtN88XnrpJQIDA3F3d1c/UFOnTuXs2bP4+/sTFhbGvXv3eO6553jmmWfw8fEp1o6TEREReHp6EhQUxOnTp9XrR44cyYYNG9TjeHl54evry1tvvcW+ffuIjo4mLCwMf39/zp49y1dffUXLli3x8/PjhRde4O+//1afZ+LEibRt2xZ3d3f1OQHmzZuHj48Pfn5+auI9e/YsPXr0ICAggA4dOnDq1Kki44+JI8JAygAAIABJREFUiVFbCBf2dzJUrVo1pkyZQkBAAEFBQRw8eFB9THR0NKBNCmFhYbRs2RJfX18WLVqk/i12796Nv78/8+fPL/R+MTExdO7cmcGDB+Pj48P9+/cJDg7Gz8+PZs2asW7dOqOvjagYJD9IfqjU+aHEbdNEPvb29oqfn5/6r379+sprr72mKIqi9OrVS1mxYoWiKIqydOlSJTQ0VFEURRkxYoQSHBysZGVlKYqiKMuXL1cfkys6Olpp3769kpGRoYwfP14JDw9XFEVRdu3apfj5+SmKoijvvPOO0qZNGyUtLU1JSUlRHn74YSUjI0M5f/684u3trT5XZmamcvv2bUVRFCUlJUVp3Lix2qmuatWq+X6nQ4cOKc2aNVPu37+v3L59W2ncuLHywQcfqLGvX79euXHjhtKkSRP1eW7duqV3e67r16+rP0+fPl1ZsGCBer8XX3xRyc7OVuLi4pTGjRsriqIoP/zwg9KmTRvl/v37iqIoakfALl26KAkJCYqiKMr+/fuVzp07F/m6/PLLL0pwcHCRfydDgNrdr0+fPkrXrl2VjIwMJTY2Vv2bL1q0SHnvvfcURVGUtLQ0JSAgQDl37pze8Yzdz9XVVTl37pyiKNqunaNHj1Yfl5qaWuTvJcoXyQ+SHyQ/FMwmNqMr76pUqUJsbKx6ecWKFRw6dAiA//3vf2zcuBGAYcOGMXnyZPV+/fr1w8HBocDnTExMJCwsjJ9//hlHR0f27NnDt99+C0CXLl24ceMGt2/fBiA4OFjdfKtOnTr8+eef+Z5PURSmTZvGb7/9hr29PcnJyfz555/UrVu3wOPv3r2bvn374urqClDg/gw1atTAxcWF0aNHExwcrH6bMHTy5ElmzJhBamoq9+7do3v37uptffr0wd7eHi8vLzXunTt3MmrUKPXYDz/8MPfu3WPfvn3069dPfWx6enqBxytMQX8n3S2kAZycnOjRowcAPj4+ODs74+joiI+PD0lJSQD8+OOPHD9+XP1Gdvv2bRITE3FyctJ7rqLu16pVKxo1aqQe56233mLKlCn06tWLDh06lOj3ErZN8oPkB8kPBZMCxMJ0N5OqWrVqgfe5f/8+/fv356uvvlJ3FFSK2O44d9dPAAcHB7KysvLdd82aNaSkpHD48GEcHR1p2LAhaWlpxY61IBqNhoMHD7Jr1y6ioqL47LPP+Pnnn/Pdb+TIkWzatAk/Pz9WrFhBTEyMeptu7Lm/o6Io+Y6dk5NDrVq19BJ5SRXn7+To6Kge297eXn2Mvb29en9FUfj000/1EiWg93sZu5/ua9+kSRMOHz7MDz/8wL///W+6devGzJkzS/17ivJL8oPkh8qUH2QOiJm1bduWqKgoQPshb9++vdHHjBo1ilGjRulVuh07dmTNmjWA9g1au3ZtatSoUehzVK9eXd3xEbTVdZ06dXB0dOSXX37hwoULRcbQsWNHvvvuOx48eMDdu3fZsmVLvvvcu3eP27dv8/zzz/Pxxx+rH37DY9+9e5d69eqRmZmp/g5F6datG8uWLVPHgm/evEmNGjVo1KgR69evB7Qf3mPHjgHa3Rr//e9/G31eU+nevTsLFy4kMzMTgISEBO7fv5/v9y7sfoauXLmCq6srQ4cO5a233uLIkSOW+UWE1Ul+kPxQmfODnAExswULFvDSSy/xwQcf8Oijj7J8+fIi73/hwgU2bNhAQkICy5YtA2DJkiWEh4czatQofH19cXV1VbdtLswjjzxCu3btaNasGT179mTKlCn07t2bFi1a4O/vz9NPP13k45955hkGDBiAv78/DRo0KPC03927dwkNDSUtLQ1FUZg/fz4AAwcOZMyYMSxYsIANGzbw3nvv8eyzz9KgQQN8fHz0PoQF6dGjB7GxsbRo0QInJyeef/555syZw5o1axg3bhyzZ88mMzOTgQMH4ufnx9mzZ4tMtqY2evRokpKSeOaZZ1AUhUcffZRNmzbh6+uLRqPBz8+PkSNHMmnSpALvZ+jEiROEhYVhb2+Po6MjCxcutNjvIqxL8oPkh8qcH2Q3XFHuDR06lPnz5/Poo49aOxQhhI2R/GC7pAARQgghhMVZZQ7I+vXr8fb2xt7eXp0NnisyMhIPDw88PT3ZsWOHNcITQliR5AchKgerzAFp1qwZGzdu5JVXXtG7Pj4+nqioKOLi4rhy5QpBQUEkJCQUuhRNCFHxSH4QonKwyhmQpk2b4unpme/6zZs3M3DgQJydnWnUqBEeHh4cPHjQChEKIaxF8oMQlYNNLcNNTk6mfv366mU3NzeSk5OtGJEQwlZIfhCiYjHbEExQUBDXrl3Ld31ERAShoaEFPqaoZjqGFi9ezOLFiwE4deoUTz/9NIkXr5ORUYagdY9rD17utU3zZEJYgJKaCOQA8MfNR1BySvc8VTX3uH79uukCK4A18oMQouySkpJMlh/MVoAUtR1zYdzc3Lh06ZJ6+fLly2qnP0Njx45l7NixgHZXQ8PJamXVvN8S7usUMy5OcHT9aJMeQ4iymP38bNLS8ro1Rrs9TlYV7efF00XD0eU9SvW8LVq0MEl8RSnv+UGIysqU+cGmGpGFhIQwePBg3njjDa5cuUJiYiKtWrWySiyGxUbzfktoGroEkGJE2IZ1jnXwfNlPveznmMPqkDZWjMi8bCk/CCHKzioFyHfffceECRNISUkhODgYf39/duzYgbe3N/3798fLywuNRsPnn39uMzPcdQuO3EJECEtr3m8Jaf+cmXNxsmfDC89aNyAzKI/5QQhRchWiEZmlT7E2DV3CH5vlDIiwPEu89yrakEVF+32EsCZTfp5sahWMEEIIISoHKUCEEEIIYXFSgAghhBDC4mxqFUx54eKkPxFVVsUIc+kxao/eUlsnJysGI4QQJiQFSCkYFhuyKkaYy6Ubp5gxL1697JRTcOMtIYQob6QAEcKG5djBoKf/Y+0whBDC5GQOiBBCCCEsTgoQIYQQQlicFCBCCCGEsDgpQIQQQghhcVKACCGEEMLiZBWMEEIIIbRycuD2bUhIgLNnITkZzp+HM2fg/n2THkoKECGEEKIyURS4cQP++ANOncorMK5d09720EPg6QkNG0KTJtCzJ7i7Q5Uq0LKlycKQAkQIIYSoiLKz4dIlOH4cTp/WFhznzkFWFtSuDd7e8NRT0L07TJgAdeuCneWaHUoBIoQQQpR39+7BkSNw7Ji20Dh9GjIzoUEDaNZMe0YjNBQaNwYHB2tHC0gBIoTN0d3/RfZ+EULk8+AB/P67tuA4elQ7V8PVFQICwN8fOnXSFhzOztaOtEhSgJiAbE4nTEl3/xfZ+0UIwd27sG8f/PYb7N+vnafRqpW24OjdGxo1Avvyt6jVaAGyd+9e/P39qVq1Kl9//TVHjhxh0qRJNGjQwBLxlQuyOZ0wpfK0/4vkByHM4NYt2LNHW3D8/rv2VGi7dhAUBDNmaCeDVgBGS6Zx48bh6urKsWPHmDdvHg0aNGD48OGWiE0IYeMkPwhhArdvw8aNMGkSdOgAQ4bAyZPQpw/s2AE//gjvvAOdO1eY4gOKcQZEo9FgZ2fH5s2bmTRpEi+//DIrV660RGxCCBsn+UGIUlAUiI+H7dvhp5+0vTe6doWhQ+Gjj0BTOWZHGP0tq1evTmRkJF9//TW//fYb2dnZZGZmWiI2IYSNk/wgRDEpinbC6Lffws8/Q9Om2v4aUVFQq5a1o7MKo0Mw69atw9nZmaVLl1K3bl2Sk5MJCwuzRGxCiP/f3p3Hx3D/Dxx/5SbOEuqsm0rkIOJISOIsDUERR6ijqtVq9RC0lChB8aujzqizVZQqoWhRaYOSuomvirjjaIokQu58fn+MjCSyNonNHsnn+Xjksdmd2Zn3zs5+9r0zn3l/jJxsHyTpOYSAo0fh00+hdWv49lvlCpU//4RVq6Bv32KbfICWIyDp6ekMGjSIffv2qY+98sor8hyvJEmyfZAkTS5dgk2b4JdfwNlZSTTmztVrkS9T8NwExMLCAltbW+Li4ihXrpy+YpIkyQTI9kGSsnj0CHbtUo5ylCoFAwbA+PHFpj9HQWjdMiVKlMDR0ZFOnTpRqlQp9fGFCxcWamCSJBk/2T5Ixd6NGzBvnlKnw9cX1q9XypxLWmlNQHx8fPDx8dFHLJIkmRjZPkjFkhDK1SvLlinjqnz4IcyeLY925JPWrTVkyBASExO5fv06jRo10slKAwIC2LFjB9bW1tSrV4/Vq1dT/klHnJkzZ7Jy5UosLCxYuHAhr732mk7WKUmS7sn2QSpWUlJgxw5YsACcnGDaNGVAN6lAtF4Fs2PHDlxcXOjSpQsAp06dwtfX94VW2qlTJ86dO8eZM2do2LAhM2fOBOD8+fNs3LiRiIgI9uzZw3vvvUd6evoLrUuSpMIj2wepWEhNfXoFy8mTShKyaJFMPl6Q1gQkMDCQ8PBw9ReIi4sLV65ceaGVdu7cGcsnh6patWrFzZs3Adi+fTv9+/fHxsaGOnXqUL9+fcLDw19oXZIkFR7ZPkhFWkYGrF2rlEGPi1MSj+nTQXa61gmtCYilpeUzPdzNdHgp0apVq+jatSsA0dHR1KxZU51Wo0YNoqOjdbYuSZJ0S7YPUpGUnq5cQuvurgxrHxam1PKQnUt1SmsfkCZNmvDDDz+Qnp5OZGQkCxcuxN3dXeuCO3bsyJ07d555PCgoiB49eqj/W1pa4u/vD4AQ4pn5NTVmwcHBBAcHAxATE6M1HkmSdE+2D1KRExYGM2Yop1e2boVq1QwdUZGlNQH55ptvCAoKwsbGhoEDB/Laa6/xxRdfaF1w1uJEuVm7di07d+5k//79aiNSo0YNbty4oc5z8+ZNqml480eOHMnIkSMBaN68udZ4JEnSPdk+SEXGvXsQEAAxMcp4LPb2ho6o6BNa/Pjjj3l6LD92794tGjduLP79999sj587d044OTmJpKQkcfnyZVGnTh2RlpamdXmurq4vFI+uveq7wtAhSCbM0PtPfj5Psn2QioQlS4RwcxNi3z4hMjIMHY1R0+XnSesRkJkzZ9K3b1+tj+XH6NGjSU5OplOnToDS0WzZsmU4ODjg5+eHvb09lpaWLF68GAsLiwKvR5JMwaxZB0lKSlPvW1hp7ZplNGT7IJm06GhlBNqmTeGvv0DuT3qlMQHZvXs3u3btIjo6mg8//FB9PD4+Xu2hXlCXLl3SOG3ixIlMnDjxhZYvSabkh7OxNOznpt53qf3QgNHkjWwfJJOWnq7U8ti2Tblt2tTQERVLGluKatWq0bx5c0JCQnB1dVUfL1OmDPPmzdNLcJJUHKSlxbGl7VL1fuqBg0B3wwWUB7J9kEzW3bvK4HCennDggDzqYUAaExBnZ2ecnZ0ZOHAgaWlpOq10WNTZWIG93/fK/yXMObluoIEjkoyZSEoiafHT+2Y2HQ0XTB7J9kEyOUIoV7XMnQtffaUkIJJBaT1WumfPHsaOHUtKSgpXrlzh1KlTTJ48mZCQEH3EZ5KOtvwPkZwEQNO/a2qZW5IEJSYFGjqIApHtg2QS0tNh8GAoWVIZw6V0aUNHJFHASqhXr14t7LhMmk3ABEpMCnzypfJs7QJJKipk+yAZvehoaNkS2rSBlStl8mFEtB4Bya3SoSRJEsj2QTJyW7bAnDkQHAzNmhk6GimHQquEKklS0SfbB8loffYZ/P03hIYqp14ko6P1FMw333xDREQENjY2DBgwgLJlyzJ//nx9xCZJkpGT7YNkdB4+hHfeUW737ZPJhxHTegTE1taWoKAggoKC9BGPJEkmRLYPktGIi1MGkNu+HXx8lAJjklHTmoAcO3aMGTNmcPXqVdLSnlZrPHPmTKEGJkmS8ZPtg2QUhIB58+DBAxg4EJ4MaCgZN60JiL+/P3PmzMHR0RFzc9MpES1JUuGT7YNkFPz8IDYW9uyRhcVMiNYEpFKlSvj6+uojFkkqNroMezr+i7W1gYN5AbJ9kAwqJQVmzAAzM6W+h2RStCYgU6dOZcSIEXTo0AEbGxv18TfeeKNQA5OkouxGUhyTNyhF6izPnzVwNAUn2wfJYISAH3+Ef/9VxnORTI7WBGT16tVcuHCB1NRU9RCrmZmZbGAk6QWIjFj6PfgbgNQLEWBv4IAKSLYPksFMmgS//qqcdrGzM3Q0UgFoTUBOnz7N2bOm+wvN0EpYC+wH7FTv25Sw5OTqLgaMSDIGWcd/MYWxXzSR7YNkEOvXw++/Q3g4yL5HJktrAtKqVSvOnz+Pvb2J/kQzsL8H/wQiSb3vuGqwAaORjIfpjv+SlWwfJL2LioJVq2DxYpl8mDitCcjBgwdZu3YtderUwcbGBiEEZmZm8jK7PLL03ZP9gVXfGiYQSSoEsn2Q9Oq33+Dtt5UERJZWN3l5Gg1XkiQpN7J9kPTm8mWls+kff0Dt2oaORtIBrQlIrVq19BGHJEkmSLYPkl4IAf7+0LWrTD6KEK0JiCRJL27WrKd1PwAsrOS5a0nKk8ePlYHlmjaFyZMNHY2kQzIBkSQ92Hc2jTb9vNX7jnUSDBeMJJmSbdsgLU257FYqUmQCIkl6cDvjEY19n3bMdDy/F+hmuIAkyRSsWAFBQfDXX1C1qqGjkXRM63HgrVu30qBBA8qVK0fZsmUpU6YMZcuW1UdsL2zNmjXs3LlT+4x5kJGRof6/bt06PD09dbZsqehTCo9tVf98L0QYOiSdkO2DQrYPhUAIWL0ajh2TyUcRpfUIyLhx49ixYweNGzfWRzyF6tChQ+zYsYM7d+4wadIktm3bhqenJy1atMDPz4/169czadIkkpOTSU9PZ/78+Xh6etK9e3fat29PixYtAHjzzTezNTj5YW0N3gNCAShRwpI9q9vo6uVJRixr4TEw7eJjWcn2Qbftg/RESorS56NzZ1nltAjTmoC8/PLLRaJxAbC2tiYlJQVbW1u2bt3KO++8w4QJE0hMTMTLy4t9+/Zx9epVGjduzLVr14iOjsbW1pYJEyboLIbJU/9Hivl5AKaPswdkAlI8FI3CYznJ9kG37YMEZGTAzp3Qty8EBho6GqkQaU1AmjdvTr9+/ejZs6fJDzb11VdfsWHDBg4fPsyBAwcoV64cpUuXZt68eXz33Xf8+eefeHh48OGHH6rPKVeunE5j6G3/f+r/X5rJomSSaZPtg27bBwnltMucObBjh6EjkQqZ1gQkPj4eW1tbfvvtN/UxUxpsasmSJezcuZO6devi6enJlClTePToES+99BIAvXv3ZsOGDZQpU4bXXnuNd999l4CAAGJjY/nmm29yXebOnTtZu3YtJUuWpHTp0nh7e+vxFUmmosuwp5feWlsbOJhCItuHZ8n24QXEx8OWLRAWBpUqGToaqZCZCSGEoYN4Uc2bN+fYsWP5ft758+eZOnUqs2fPNkhBpcY9vuV/20fofb2SfjgM+IXJG2oCYHl+Jb3tTWPI8IJ+noyVqbYPxVKbNuDiAosWGToSSQNdtg8aj4DMnj2bcePG8cEHH2BmZvbM9IULFxZ4pV988QXbt2/H3NycypUrs2bNGqpVqwbAzJkzWblyJRYWFixcuJDXXnutwOvRxt7enk2bNhXa8qXiTSDohxMAaZGnoQiN1ybbB0nnzp+HGjVk8lGMaExAMjuWNW/eXOcrDQgIYNq0aYDSUH355ZcsW7aM8+fPs3HjRiIiIrh16xYdO3bk4sWLWFhY6DwGSSp0GXHwIBAAkV60Su7I9kHSqfBw6NMHli83dCSSHmlsFbt37w7AkCFDdL7SrHUCHj16pP6C2r59O/3798fGxoY6depQv359wsPDad26tc5jMAYWVuYEBoaq90uUsGTCBHlVTFGR9dLbonLZbSbZPkg6tWkTfPcdeHkZOhJJjwz2s2zixImsW7eOcuXKceDAAQCio6Np1aqVOk+NGjWIjo7O9fnBwcEEBwcDcPz48VwPA5uCiJ+y3//sM8PEIRWOkiFZ7owzjTfX1dXV0CHI9qE4+vprQ0cg5YEu24dCS0A6duzInTt3nnk8KCiIHj16EBQURFBQEDNnzmTRokVMnTqV3PrDamo4Ro4cyciRIwHj7jTXtO+3JKU8vV/CGk5ulh1Pi5ucg9H9dPIi6Ta2ANiUMOfkuoGGCu0ZhXFaJSfZPhRTqanw66+wbh3cvw9+fspf+fKGjkzKI122D4WWgOzbty9P8w0cOBAfHx+mTp1KjRo1uHHjhjrt5s2bauczU5Uz2Wja91sa91Dqf8hkpGjLehnu/ZTLNPCprk5zrl2d74f6GCo0g5PtQzGSkgK//w6bN8P//gedOsGMGVC/vqEjkwws3wnIkiVLqFixIr1798bSsmD5S2RkJA0aNAAgJCSEV199FQBfX18GDhzIJ598wq1bt4iMjFTLGxcVWROOzEREKppuPLrG5LX3ALCIukKfJsMNHFHhk+2DBChJx969Sk2PixehQwcYMwYcHUGeDpOeyHcLIYTg4MGDrF+/npCQEO1PyMWECRP4559/MDc3p1atWixbtgwABwcH/Pz8sLe3x9LSksWLFxfpHu42VmDv9/3T+0Z2KN5YZT2dYUwdd3OeZjHLSKHH/91X/rd5GZoYKjL9ke1DMXb3Lvz2G+zZA9euKUc6xo4FBwdDRyYV1MOHEB0NlpZw/DjUravTxRfpQmRr1qzBzs6Obt1efNjzjIwMzM2VwYMDAgK4f/8+sbGxrFmzhjJlyhRomclzZiGSk9T7Tf+uyf+2v/XCsRZ1gYGhBAZ6P/O/oXX0D6VNP2/1/j8PdrJhyIvve/pW1PpMmGr7YPRSU+Gvv5SE4+BBeOkleO015a9ePUNHJ+VHWhpcvQqVK8PPPyuDAK5dC2+9BSdOKFcnJSdDuXL6KUSWKTk5mZ9++omrV6+Slvb0193kyZN1EoA+6Wq0yzlz5gAwb948Tp06Rdu2bQsUj01AjkGs5CkZk3Y74xGNfc+o9x3P7wVMLwHJD9k+FF77YJQuXYL9+2HJEihTBjw8oEsXZdC4ojreQFFz+jQ4OSnv4bvvQkAAfPUVHDgAw4eDpydUrQqZgyxmFvsrUULnoWhNQHr06EG5cuVwdXXNNtiUKdLlaJd37tzh2LFj2QamkvQk6SA8CH3yvyXgbbBQsp52Mbd4SL8HW9VpqRciilT109zI9qGItw+PH8Mff8DChRAXB/fuKV9MS5YoVUt374Z27QwdpZRVWppyyuTwYXB3V5KL8eOVPjgLFiinUpydlaMcFhZPL39++23ltk4dvYWqNQG5efMme/bs0UcshU5Xo11GR0czYcIEeQ7aUEQavBQIQAmbQIMWc9t3Nk097eLY7BC81F+dln4+ECvTGJOtwGT7UMTaByHgwgUlwfjf/+DyZaVC6YAByuWyOX8Fv/kmbNgArVrp9YtLQnlvatdWri7q00dJDGfNUm7nzIHMK8YGD1ZuFzwZi2r4k87wTzp6G5LWBMTd3Z2zZ8/i6Oioj3h0rjBGu/T19aVhw4ZMmDCBUaNG4ezsrJNYs3ZKlR1S82bCR8BL3ur9rMmIPqRlQKCv8v+Giz+TNH2vOs3MRveHLI2NbB+eVVjtQ6GJi4Pvv1dGoD1+HF59Fdq2hY8+0t6Xw9ZW6WxavjzMnQuffiqvctGV+Hgl4Tt1Cho3ho0boWdPJTn84APl0uYBA5T3SAiYOVM5ojF3rvL8fv2UWyO+VF1rJ1R7e3suXbpEnTp1sLGxQQiBmZkZZ86ced7T9KqojHaZtVOq7JCqWeD4QAK/ClTuPAhUj4ZA4XdKzXmly9bLSUxcp3zATWnE2+fJz+dJtg8mKCNDOa2ybRscOqR0Ju3cWTn3361bwROImzehVCn48Ud45x3dxlwUPXyoJAyXLkGlSsp4OK++CqGhynsRHg7e3nD9OjRpAunpSqJXwMvbdUWvnVB3796tkxUZI2Mb7TJbp1TZIVWjg1dfITDzCs9EbyiZZVpk4a476ykXgCbN/lL7fRSHPh85yfbBRERHK0c5wsPh77+VuhwNGyodD3V1lU6NGsov8R49lP4Ht28r/xv4C9MghIDYWKU/xv37yv27d5VtcfcuVKwICQlQ/UlxwowMaNNGSeAyk7cnAz5SlJLfHLTuGbVq1eL06dOEhYUB0LZtW+M/pCgVaWnp5uppDx6EZjsF472hkNed5ZQLwKbEv6FkIFA8+nzkJNsHI5WaqhQBO3tWSTycnJTLZGfNUg7ZP7lkWOfMzKBKFeUvNVXpn9CggXKFTN26ULp04axXn5KTITFRKbYWEwNWVkp/C1tb5WhF1arw339K4pCcrNwXApo2VRIQW1tDvwKjoTUBWbBgAStWrOCNN5SWddCgQYwcOZIPPvig0IMrzqytwXtAqHq/RAlL9qw2joJbRsWshDrkPYClWX21H4iuOqRmPe1y3yqJTTw9vWB57QpJW5T1F4c+HznJ9sGI/P03nDwJ8+Ypv7BtbODzz5VOiVlGGNYbKyuljwLArl1Kh8lJk5SYzpyBFi0KLxHKDyGUZALg1i2oUEFJ3OrVU/rFuLoq8XfqpBwx6tBBuZTV3V05ivHqq8pRpCpVQI7MnD9CC0dHR5GQkKDeT0hIEI6Ojtqepleurq6GDkHnNp/9SGx8vED9s+/7vaFDMpiZM8PElCkH1L8OA9aq07aKc2KjOK3++W35S502ZcoBnay/w8ADYsp2IaZsF2LArlvZpqVu89LJOoxJfj5Psn0woAcPhPjlFyE6dhSieXMhnJ2FWLZMiIgIQ0emWVqacrthgxAZGUKMGaPc//JL5XbbNuX2+HHl9vZt5fbBA+U2KUm5TUxUbuPilNv//lOWd+mScvv338q69u4VIiVFiE2bhHj8WIiFC4WIjxdi0iQh7t8XYsoUIWJihPjqK+X+hg3K9MOHlXXcuKEsLyOj0DaJqdHl50nrERAhRLZLySwsLHIdlVLSre67X0Yk31fvByaXMmA0hpWUlJatY2m/n46oRyGssaAfTuq07VaH1P4hR65ZFuhTOB3oAAAgAElEQVQS3ZwdTbGwVE+7bEr6GR78q04S6cXw/HYWsn3Qs6goWLRIuTLi1i3o2FGp79CxY6EUitK5zH2l/5PL1efPV27ffVe5zeyPcu4cNGsGq1fDZ58ptSq+/BK++AJmz4YpU5T6FnPmwLRpsHgxTJ4Mv/yiXCFy4YLy/NRUZZ116ihHhPr0UU4DffmlcrooMFBZ37hx2ePKPJJRo0ahbo7iTmvrOWzYMFq2bEmvXr0A2LZtG2+9Ja/OKGyySupTByN52ukUsLBKz5Z0ZNWj03H6lfQAIJA22fprzJp1UOPpmaxJx5FrlrTq5a1Oq9ry7tOEJyM121U3xbHfR1ayfShkSUlKP4pdu5TD/pUqQa9eyhd2o0aGjk53KlVSbtu3V27ffFO5/ewz5fbLL5Xb2bOV26++Um6nTVNuMyvvZtZoGTRIue3aVbl1c1Nuq1bVbdzSC9GagHzyySd4e3tz8OBBhBCsXr2apk2b6iM2KYusfUKKW3+QZzt+Hgc8cp3XWgg1WYiyepnAkJfVaSWc2xD4pD3KmoyAsk3V8WVCcqyPu2rCk3bgE5JOPlCnFcd+H1nJ9qEQXLgAK1cq/TmuXoXu3cHXF777rnheUSIVWRr35vj4eMqWLcv9+/epXbs2tWvXVqfdv3+fChUq6CM+6YkvpkSQXu8xAF8OqQgUnwSEtOvZOppSsrLGWXulJELyk3Lo7SrTr8R76rRZu7McSbHPvv2SeDotyuoum7irTrNODIMkZZkiGUpMCqS4k+2DDgmhlDQ/ckQ55eDmpnTYXL8eXn5Z69MlyVRpTEAGDhzIzp07cXV1xSxLYRrxpNDQ5cuX9RKgpMjaJ6S49QeJMTNn00tPz3OY/xFGUlggoByByHa6qvzT/60fL8h2xUqjrhb0IvehwX8mghTSAXDGgl5ZT/EkbVVPuxT3Uy6ZZPvwghITlVoZgYHKFRjm5kr/hnHjisalqpKUBxoTkJ07dwJw5coVvQUjaVaci5QJyNbnIylsq3oUIml6oMbn9Up+ALZPn/czEdkSkqyydWaNnQViszotLewIaWeU9RT3Uy6ZZPtQAMnJSmfKU6eU5MPPDz75RDm9YmpjxkiSDmg9odihQwf279+v9TFJf0pYC+wH7FTv25Sw5OTqLgaMSLdyXoViYV3AxjlHjZBeZiWyHSHRSCRl62iadiZQnnbRQLYPWiQkwJ498M038OCBMmjbJ59AUamwKkkvQGMCkpSUxOPHj/nvv/948OCBemldfHw8t27d0luA0rP+HvyT8iX5hOOqwQaMRvdyXnZ77qejGuc1symR7ShItlMyOZONrP1Icoqd9XSbmsmjHNrI9uE50tOV8VC+/x6uXVMuB12wAFxcDB2ZJBkVjQnI8uXLmT9/Prdu3cLV1VVtYMqWLcv777+vtwClZ1n65hj+fFXROiWT87Jba6s0jfPmvFz5eadkch4ReWbak6MeyqCAT+eTp12eJduHXBw/rnQc/fFH6N0bRoxQRi+Vo8NKUq40JiBjxoxhzJgxfPPNN7KsspErapfo3s54RGPfp301nB8fQ9Nlt/mSl9MvgEhOkqdctJDtwxPx8Urn0bNnlUHG5s5VhkW3sTF0ZJJk9LT2Afnggw84d+4c58+fJynp6WH/NzMLxUgGVxQu0c3a78Pc4qE6wiwonUCTzuSt9kbWUzLPXCHzHMpRjySty5eyK7btw/nzMHKkMtpp8+ZKsbCKFQ0dlSSZFK0JyNSpUwkNDeX8+fO8/vrr7N69mzZt2hT9BsaEZL1Ed1pGmQKVHze0rP0+zv10FF7qr07LTyfQrAnHc0/H5CCPehRMsWofkpJg61blCEf16sqgal26GMeAapJkgrQmIFu2bOH06dM0bdqU1atXc/fuXUaMGKGP2KQ8yvql67h2Z7bh6fdtCiVvxwAMK2u/j+f1+ciPnB1Utc0r5V+xaB/i45UxSzZtUjqS/vorVKtm6KgkyeRpTUBKliyJubk5lpaWxMfHU7ly5WJZZOjYsWP8+uuvlCxZkk8++cTQ4WjU23U/6fWevj+bvzeNw8JZy60/r9R6fuT19ItUcEW6fUhIgI8/VkqiOzgow93b2ho6KkkqMrQmIM2bNyc2Npa3334bV1dXSpcuTYsWLfQRm1FZtWoV9erVQwihVnt8/PgxQUFBBAUF8eGHHzJ69GhWrFiBtbU1Pj4+pKenc+LECW7cuMHcuXPx8PDggw8+4P79+6SmplK/fn18fHx0GmfOUXRnYGuUHVRz1vq4b5X0dMA3OZqqySiS7UNaGowapQz+5u4Ohw7JTqWSVAi0JiBLliwB4N1336VLly7Ex8fj5JT7SKTGJuW7NWRcu5rn+c1r1cZ68NBcp92+fZvZs2ezfv16Tpw4gaurK7a2tpQtW5ZNmzbRrl07QkJCsLOzw87OjkuXLlGxYkVKlSpFREQEAE5OTvTv35+jR48SFhZGXFycDl5ldjl/9U8+P4Y0e+UL4csBN3S+voJ6ptbH1r+fdjw1KwHyh6ZJMOX24RlCKJVK161Thrfftw/KljV0VJJUZGlNQHr06EG/fv3o0aNHtgGnTIGmZKIg3n33XWbMmEFcXBz+/v7q42+99RZdu3bl6NGjXLx4keDgYMqUKYOrqyu//vorNWvWJD1dGWPE5smvqPv371OyZEnOnTuns/g06XHpH4hUhqieljrEaDqoPlPrwzJF1uEwQabcPmRz7hwMHw6OjrBlCzRoYOiIJKnoE1qEhoaKUaNGiVdeeUX06dNHbN68WSQmJmp7Wp7MmTNHACImJkZ9bMaMGaJevXqiYcOGYs+ePXlajqurq07iKYilS5eKDRs2GGz9+eG/OkT03n5H/XMcuEOv6585M0xMmXJATJlyQDgO3i02itPq39ZH89X5EqdN0WtcUnb5+TyZfPuQkSHEqFFCuLsL8f33LxqyJBV5uvy+1XoExMvLCy8vL9LT0/n9999ZsWIFw4cPJz4+/oUSnxs3brB3715eeeUV9bHz58+zceNGIiIiuHXrFh07duTixYtYGPFATe+++66hQ8izNRUWg/g/9b5r6mC99g/ZdzaNNv28AWjS8jb9qPp0YvJWedrFBJl0+3DtGrz/PlSoAHv3yg6mkqRnWhMQgMTERHbs2MGmTZs4ceIEQ4YMeeEVf/zxx8yePZsePXqoj23fvp3+/ftjY2NDnTp1qF+/PuHh4bRu3fqF1yc9W8J9coMxpNm7ARD05i2dn57J2dEUC8unV7ok/QwP/n06TY6/YrJMsn04exYGDIB33oHiXMlVkgxIawLSr18/jh49SpcuXXj//ffx9vbG/AUL74SEhFC9enWcnZ2zPR4dHU2rVq3U+zVq1CA6OvqF1iVplrV/yLaO46C8tzot9OeDJGVJSLLKmpw8k2RkceSaJa16PV1m1ZZ3n17pkpGabcTZrP0+ZJ8P02GS7cPhw8pYLYcPQ506LxSrJEkFpzUBGTZsGD/88EO+D3N27NiRO3fuPPN4UFAQM2bM4Lfffntmmsjl8kszDQM5BQcHExwcDEBMTEy+YsuPixcvMmPGDHr27EnPnj3ZuXMn586d4+LFiyxfvhwrK6tCW3dhy3pEpPe5j0mvd0m9H4Ubmkq6Z01OSpSwzHY1y6zdkJSq/F+p5V0ad80ypgsW9OLJFRI5TrnISqSmyeTah6VLlb/Tp6Fy5XzFLEmSbmlNQDw9PZk5cybXr18nODiYyMhI/vnnH7p16/bc5+3bty/Xx8+ePcuVK1fUXzc3b96kWbNmhIeHU6NGDW7ceHqp6M2bN6mmoeLgyJEjGTlyJKDUIigsDRs2ZOjQocTGxgLQrVs3unXrxujRo0lJScHKysoo64HkV876IZZ9fiTt1VK5zhtl1QBSXwYgiexXs0RZ3aWb710gR8IBT4a836z8L0+5FAkm1T78+adyiW1IiEw+JMkI5OkIiKurK4cPHwaUw559+/bV2sBo4ujoyL//Pj33X7t2bY4dO4adnR2+vr4MHDiQTz75hFu3bhEZGfliRY0eroG0q3mf37I2lBmqdbYFCxbg6+tLqVLKF7Qx1gPJr5z1Q3qEdIF/lNMzmJXIdrTEumsEKdzNdTnZko6sCceT5chLbYsWk2kf7t+HwYNhzx4w5cuFJakI0ZqAREVFsWnTJjZs2AAopZdzOxSqCw4ODvj5+WFvb4+lpSWLFy9+sStg8pBMaHPnzh22bNlCYmIiTZs2ZfPmzRw/fhwhBG5ubrz00kuA8dcDya+sCUfq1o7wIFC930vrs7Mc5cjSzyMrecqlaDCZ9qFjR1i5Eho3LpTYJEnKPzOhpbVwd3dn//79eHh4cOLECaKiohgwYADh4eH6ilGr5s2bc+zYMYPGsGzZMsqXL0///v21z2xisg5Vn9OLDHkvx2oxTvn5PJlE+9C0KcfeeAO++MLQoUiSydPl963WIyBTp06lS5cu3LhxA39/fw4dOsSaNWt0svKixJTqgeTX8xKF5Dmz8jXirDzqUbSYRPtgYSGTD0kyQloTkE6dOtGsWTOOHDmCEIIFCxZgZ2enj9gkEyCPYhRvsn2QJKmgNCYgJ06cyHa/alWlauX169e5fv06zZo1K9zIJEkyWrJ9kCTpRWlMQD799FONTzIzM+P3338vlICMTc46IMuXL+f+/ftcunSJr7/+mnLlyhk6REnSO9k+SJL0ojQmIAcOHNBnHEYrZx2QkiVLcuPGDVJTUyn7ZKjuolAHRJLyQ7YPkiS9KI01k2fPnq3+v3nz5mzTPv/888KLSIcOEMUmzuT57wBRWpd58+ZNlixZQuvWrdW6HprqgNStW5dLly4RHx+fax0QV1dXUlNTjaIOiCTlR1FoHyRJMiyNCcjGjRvV/2fOnJlt2p49e3LObpTaUY9+OOX5rx31nllGZh2QHTt2cO3aNUqWLElgYCAnTpygVq1a6nxvvfUWc+fOpUePHnTr1o27d++SnJyMk5MT586dw9bW1qTqgOhSaGgotWrVwtvbG29v7wIlXKGhoUyaNOmF4lizZg0ZGRkvtIy88vb2LtDz2rRRyt8XNNbY2Fi2bt1aoHXnR1FoHyRJMiyNp2CylgfJWSqksAoNGaMqVaqwaNEi9f7HH3+c63xbtmzh008/xdzcnFdffZWvv/5anebi4gLAwIEDAZg/fz4AXbt2Laywjc7gwYOZPn36M49nZGS88OBlebVmzRoGDRqkdX36jEmTvMaaU2YC8sYbb2id90Vep2wfJEl6URpbn6yDPOUc8EnTAFDF2bvvvlski5AVFm9vb8aNG8ebb75JYmIiAwYMoH379vTr14/U1FR+/vlnWrRoQfv27dm1axegjBPSvXt3PDw8SEhIICAggHPnzvHbb7/RtGlTAIYMGcK///7LzJkz8fLyomXLlpw8eZLw8HBOnTpFhw4d+O6777h06RKdO3fGy8tLTYyGDh3K6NGj6dKlS7ZYcy4rM/7PPvsMNzc3Vq5cCcDOnTtxdXXlnXfeIS0t+wjBly5dwt3dnXbt2jFjxgwA/Pz88PLyonPnzsTHx6vz5ow1M1H19/dn/vz5JCUl0atXL1JTU+nQoQOenp707t2b9PR0goOD2bt3L97e3sTExPDtt9/Stm1b2rZtq1654uzszKBBg7KdRskv2T5IkvTChAbm5uaiTJkyonTp0sLCwkKUKVNGvW9paanpaQbh6upq6BCk5zhw4IB45ZVXhJeXlxg6dKgQQggvLy9x+PBhIYQQCxcuFD/88IMQQoglS5aIH374Qfj7+4srV64IIYTIyMgQBw4cEL6+vkIIIaZPny62bdsmtm7dKpYtWyamTJkiunbtKuLj48Vrr70mhBDi0aNHQgghIiMjxcCBA9V1pqamCiGE8PPzE9evXxdCCNG/f39x48YNMWTIEDWOrDQt68SJEyIpKUm0bdtWCCGEu7u7ePjwobh+/bqoV69etmWsWLFCrF69Wn09WZe7YsUKERwcLIQQwsPD45lYe/fuLRISEkT//v3F8OHDRVhYmJg2bZrIyMgQjx8/FkIIMXHiRPHbb7+JK1euCH9/fyGEEDExMaJ79+4iIyND3Lt3T/To0UMIIUSFChVEQkKCxvcrL58n2T5IUvGky8+TxlMwmf0VJEkXcjsF4+rqCsD//vc/jh8/zvLly0lKSmLAgAFMmjSJ6dOnk5aWxsSJEwFo0qQJANWrVyc2NpauXbsSEBBARkYG/v7+bN++nZdfVkbp/e6771i/fj3m5ua5/iL/559/GDx4MKCctoiOjs4WU1aaltWkSROsrKzU0xjm5uaULl2a0qVLU6lSpWzL8PPzIzAwEH9/fwYNGkTnzp0JCAjg7NmzxMfH06uX5hF2WrRowbZt26hVqxZ37tzh8OHDeHh48OjRI0aOHEl0dDR3796lQYMGNGjQQH3e5cuXOX36NO3atcu2vEaNGqkDKRaUbB8kSXpRWiuhFndhYWEcPnyY8+fPM2/ePCpUqMDatWsJDQ1l9erVhg7PpGV+cTdq1IgOHTrQu3dvAFJTU0lLS+Pbb7/l8OHDfP311/Tr1y/bl78QgsqVK3P79m1q1KiBh4cHPXv2ZNSoUQAsWbKEkydPEhUVxdtvvw2AlZUV6enpWFpa0qhRI+bPn0/VqlVJT0/HzMyMpUuX5tonIrdlwbOnGjIyMnj06BEPHjwgJiYm2zQrKyu+/vprUlJS8PDwoHLlyjx69Ig///yTFStWqAlQ1vkzY3V3d2f06NFMmzaN48eP8+OPP/L++++zZ88eGjZsyA8//MDEiRMRQqjPA6hTpw5ubm5s2bJF3a5Zt7skmZo1a9YQERFB9erVGTFiBKVLl9b6nI8++kjtd5fXdbi4uKh9917E1atX2bZtG+XLl9fZMosSmYBokXn+fNq0acTGxnL8+HFq166drQCZrAPyYkaOHMnbb7/NkiVLEEIwc+ZMtm7dypEjR0hISOD//u//ND63atWqODo6Urt2bWJiYnB3dweUowaenp54enqq8/r4+NCzZ09GjBhBUFAQw4cPJzk5GSsrK3766SeN68htWbkZP348np6eNGvWjCpVqmSbFhISwqJFi3j8+DGDBg2iUaNGXLp0iS5dulCzZk2qV6+ebf6ssfr4+HDhwgU8PDwoUaIEO3bsoFSpUrRs2ZKgoCCOHTtGuXLlaNCgAVWqVOH+/fv06dOH4OBgfHx88PT0xMLCgvbt2/OFHBNFKgJKlixJ6dKl1ba0Zs2a7N69m8TEREaNGsWKFSuoVasWrVu3Jjk5mYULF3L06FFWrVrFtGnTqF+/PtHR0VSuXJkOHTqwf/9+bt26hYODA5GRkVy/fh07OzuWLl1KhQoVSEtLo2XLlqxfv17tr5acnJzrslauXKk+p1+/fobeVEZN62i4pkDT6HxrDsHVe3lfTu2KMNTj2cd/+OEHhBD4+/szYcIEqlSpws8//8yqVauoV0+5dPerr76idu3aWFtbExUVRXp6OnZ2dlhZWVGxYkVu377NTz/9xO7duxk1ahRLly7l6NGjhIWFUa1aNfUKGUkyNGMYXVqXitrrMWn9+8OdO3mfv0oVyHLJd86jE5lt6QcffICzszOJiYk0aNCAhw8fcuXKFdq2bcv69ev55ptvmDdvHn369GHhwoXY2dmRnp7Ow4cPCQoKYtGiRZQrV47jx4/j6uqKi4sLNWrUYNiwYfTu3ZuLFy/SuXNn7t69my2pGDt2bLZljR07NttzRo4cWeSOgOh1NFxTllsykV+bN29m3bp1dOnShWvXrjFr1ixAObSWmXyAUgeka9euHD16lIsXLxIcHEyZMmVwdXXl119/pWbNmsW2DogkSRKQLZnQhcy21M/Pj927d1OjRg0aNmzI77//jrW1NWfPnsXS0hIzMzPMzc3JyMigfv363Lt3j1q1ahEREUFiYiI3btygQoUKZGRkUK9ePTZt2sSoUaNwcHAgNjYWJyenbOvLlHNZdnZ2zzxH0qxIHwHRp2XLllG+fHl5Ka5k8ozh86RLRe31SJIhySMgRujdd981dAiSJEmSZDJkd3hJeg5vb2/S0tIIDAxk3759eXrOqVOn1OJkq1atUh/PLLP+vHUVBm3rlSRJMgSZgEiSjrm4uPDWW28B2ROQ/NLXuDWSJEmGIBMQqVh5Xvn2o0eP4u3tjYeHh9YaL/Hx8WpveA8PD7Zt28bly5d577331IHzgoODOXv2LN7e3pw9e5b09HRGjx6Ni4uLOmBbcHAwrVq1YsKECeqys5apP3XqFF5eXrRq1Uot4d6vXz/i4uIIDg5WC5hlXu49ZswYvLy8aNu2LdevX88W8+eff46Hhwft2rXj1q1butmgkiRJBST7gGgRERHBr7/+SmRkJNOmTSM8PJyTJ08SFxfHV199Jce9MDHu7u4cOnSI27dvU7VqVR4+fMjdu3epXLkygwcPJiQkhDJlytCpUyf8/f01Lqds2bLExcWRnJxMxYoVOXz4MAkJCWodElDqm6xbt47Q0FAA7t27x+TJk0lNTWX06NF07NiRlStXcujQIcLDwzly5Ij63F69etG6dWsSExMJDQ3FzMyMdu3a8fHHH9OqVSuOHDnCsWPHsLKyIjU1FQsLCywsLJg5cya2trbs27eP5cuXExQUpC7z0KFDhIWFYW5uLgeMk0ySqRUiA2Uk79jYWHr27Kl13o8++oiePXvmeX5TL5pWpBOQNWtOcfVqbJ7nr127PEOHZt/oDg4OnDp1irCwMKysrNi7dy/z5s1j3bp1nD59GhcXF1mIzIR4eHhoLN9++vRpfH19Afjvv/+eqWaaU926ddm8eTOvv/46R44c4fDhwwQEBHDt2rVc569UqRKVK1cGlPLv//33H7Vq1cLS0vKZEvCZ969cucKnn37K48eP+eeff/j333/x8PDgl19+ITExEWdnZzZt2qQezZk9ezb79+8nNTWVxo0bZ1vmuHHjGDJkCBUrViQoKOiFy7FLUr6tWQMuLnD8OLi5wdGj0Lo1HDwIXl5w4AB06AB790KXLlCpEmQp+pjJGAqRAWzcuJEbN25ga2uLg4MDmzZt4uWXX6ZNmzYIITh06BAxMTF06tRJjX38+PFqbAkJCdm+GwAiIyM5e/YsderUYfPmzUW6aFqRTkByJhMF5e/vT/ny5Z85pJ159MPW1payZcuyadMm2rVrR0hICHZ2dtjZ2XHp0iUqVqxIqVKliIiIAMDJyYn+/furhcji4uJ0Eqek3fPKtzdt2pQtW7ZQqlQpUlNTsbKyeu6y3N3dmTt3Lt999x0nT57kwoUL1KlTJ1sComnUWCEEdnZ2XLt2jfT0dHWU3UyZ5dKXLl3K+PHj8fb2Vhu1pk2bMmbMGDw8PPDw8OCjjz7iq6++4t69e4SGhhIWFsbevXtZv359tmW2b98eHx8fZsyYwc6dO2WVRkn/hg5VbjN/XWfWyngyzhOZSXOjRhoX4e/vr/46z2xLsxYiu3TpEs2bN+fKlSskJSVhYWHBBx98QFpaGnfv3iUpKYk7T4qhRUVFMXz4cB4/fkytWrX4448/1EJkJUqU4MyZM2pRMYCOHTuqyQcop2Lt7Oz4888/cXBwoEOHDnTo0IFvvvmGmJgYvv76a0JCQrLFnzW2+Pj4bN8NAA0aNKBSpUo4Ozvz3Xff6SRWGxubZ153r169nnmOvsk+IFrs2bOH2bNnExISQqVKlejUqRMzZswgIiIiW6GZt956i7lz59KjRw+6devG3bt3SU5OxsnJiXPnzmFraysLkRkJTeXbp06diq+vL+3atctTPRcPDw+uXbtGkyZNaN68OS+99NIz89SsWZPevXtz4cKFZ6ZZWloybNgw3N3d2bp1a67r8PHxYfTo0fj5+WFtbQ0o48RYW1vj4eFBixYtuHDhAq1ateKll16idOnStG/fnl9++eWZZfXs2ZO2bduye/duvLy8tL4+STJ2WQuRXb58GQsLCxo2bEhcXNxzC5Glp6dTq1YtzMzM1EJkFhYW2QqRPX78WGshsn/++Qdra2u1bc9cV0ZGBt26dWPp0qUcPHgw23OyxpbzuwGUNmPfvn1EREToNNacyzKGommyEJmOyEJkUlFhDJ8nXSpqr0eSDEkWIjNCshCZJEmSJOWdQU7BBAYGUr16dbWn7a5du9RpM2fOpH79+jRq1Ihff/3VEOFJkmRAsn2QpOLBYEdAPv74Y8aOHZvtsfPnz7Nx40YiIiK4desWHTt25OLFi1hYWBgoSkmSDEG2D5JU9BlVJ9Tt27fTv39/bGxsqFOnDvXr1yc8PNzQYUmSZARk+yBJRYvBEpBFixbh5OTE8OHDefDgAQDR0dHUrFlTnadGjRpER0cbKkRJkgxEtg+SVPQV2imYjh07qtccZxUUFMSoUaP44osvMDMz44svvuDTTz9l1apVuVZn1FRpNDg4mODgYADOnTtH8+bNdfsCXkBMTAyVKlUydBiAjEUTGYtmuV0yrGuyfTCO91vGkjsZi2Y6bR+EgV25ckU4ODgIIYSYMWOGmDFjhjqtc+fO4vDhw1qX4erqWmjxFYQxxSNjyZ2MRTNjike2D4VLxpI7GYtmuozHIKdgbt++rf7/888/0+RJFTxfX182btxIcnIyV65cITIykhYtWhgiREmSDES2D5JUPBjkKphx48Zx6tQpzMzMqF27NsuXLweUcVf8/Pywt7fH0tKSxYsXyx7uklTMyPZBkooHgyQg3333ncZpEydOZOLEifla3siRI180JJ0ypnhkLLmTsWhm6Hhk+6A/MpbcyVg002U8RaIUuyRJkiRJpsWo6oBIkiRJklQ8GH0CsnnzZhwcHDA3N882AM7Vq1cpWbKkWq4561gsx48fx9HRkfr16/Phhx+ql+8lJyfTr18/6tevT8uWLbl69arO4gHNZaILM55MBSlfrSkuXduzZw+NGjWifv36zJo1q1DWkVPt2rVxdHTExcVFvQTz/v37dOrUiXWSUsUAABFwSURBVAYNGtCpUye1vgTovsT38OHDqVy5stqBsqDr18V7lFssxry/6JIh9j1Nctsn9SW/+6Mh4nnePllYbty4Qbt27WjcuDEODg4sWLAAMNy20RSPIbZNUlISLVq0wNnZGQcHB6ZMmQLoeNvo7HqaQnL+/Hlx4cIF4eXlJf7++2/18ayX5+Xk5uYmDh8+LDIyMkSXLl3Erl27hBBCLF68WLzzzjtCCCE2bNgg/Pz8dBZPRESEcHJyEklJSeLy5cuibt26Ii0trdDjyTRlyhQxZ86cZx4vSFy6lJaWJurWrSuioqJEcnKycHJyEhERETpfT061atUSMTEx2R4LCAgQM2fOFEIIMXPmTDFu3DghxPO3UUH98ccf4vjx49n20YKsXxfvUW6xGOv+okuG2vc0yW2f1Jf87I+GikfTPlmYbt26JY4fPy6EECI+Pl40aNBAREREGGzbaIrHENsmIyNDPHz4UAghREpKimjRooX466+/dLptjP4ISOPGjWnUqFGe5799+zbx8fG0bt0aMzMz3nzzTbZt2wYopZyHDBkCQJ8+fdi/f3++f8VpikdTmejCjkebgsSlS+Hh4dSvX5+6detibW1N//792b59u87XkxdZt/eQIUOyvQ+6LvHt6elJhQoVXmj9unqPcotFE0PvL7pkTPueoeVnfzRUPIZQtWpVmjVrBkCZMmVo3Lgx0dHRBts2muIxBDMzM0qXLg1AamoqqampmJmZ6XTbGH0C8jxXrlyhadOmeHl5ERYWBijlmmvUqKHOk7Vcc9ZSzpaWlpQrV4579+7pJBZNZaL1GU9+ylc/Ly5dMlT5bDMzMzp37oyrq6taEfPu3btUrVoVUD7o//77r15jzO/6C/s9Msb9RZeMrXR7bvukIWnaHw0pt31SX65evcrJkydp2bKlUWybrPGAYbZNeno6Li4uVK5cmU6dOul82xhFAtKxY0eaNGnyzN/zfq1UrVqV69evc/LkSb7++msGDhxIfHz8c8s1P2/ai8ajadm6iCcvcY0aNYqoqChOnTpF1apV+fTTTwscly7paz05HTp0iBMnTrB7924WL17Mn3/+qXFeQ8Wobf2FGZex7i+6ZGwx52efLI407ZP6kJCQQO/evZk/fz5ly5bV23rzGo+hto2FhQWnTp3i5s2bhIeHc+7cOZ0u3yB1QHLat29fvp9jY2ODjY0NAK6urtSrV4+LFy9So0YNbt68qc538+ZNqlWrBii/gG7cuEGNGjVIS0sjLi4u18OABYknc9k516uLePIb19tvv023bt0KHJcuaVp/YctcR+XKlenVqxfh4eG8/PLL3L59m6pVq3L79m0qV66s1xjzu/7CfI9efvll9X9j2l90yVD7nia57ZOenp4Gi0fT/mjIeDJl3ScLW2pqKr1798bf35833nhDjcVQ20ZTPJn0uW0ylS9fHm9vb/bs2aPTbWMUR0AKIiYmhvT0dAAuX75MZGQkdevWpWrVqpQpU4YjR44ghGDdunX06NEDUEo5r127FoAtW7bQvn17nf0i0lQmWl/x5Ld89fPi0iU3NzciIyO5cuUKKSkpbNy4EV9fX52vJ6tHjx7x8OFD9f/ffvuNJk2aZNvea9euzfY+6KPEd37XX5jvkbHuL7pkiH1PE037pCFp2h8NRdM+WZiEELz11ls0btyYTz75RH3cUNtGUzyG2DYxMTHExsYCkJiYyL59+3j11Vd1u20K3H1VT7Zu3SqqV68urK2tReXKlUXnzp2FEEJs2bJF2NvbCycnJ9G0aVMREhKiPufvv/8WDg4Oom7duuL9998XGRkZQgghEhMTRZ8+fUS9evWEm5ubiIqK0lk8Qggxffp0UbduXdGwYcNsVwgUZjyZBg0aJJo0aSIcHR1F9+7dxa1btwocl6798ssvokGDBqJu3bpi+vTphbKOrKKiooSTk5NwcnIS9vb26jr/++8/0b59e1G/fn3Rvn17ce/ePfU5mrZRQfXv319UqVJFWFpaiurVq4tvv/22QOvXxXuUWyzGvL/okr73PU007ZP6kt/90RDxPG+fLCxhYWECEI6OjsLZ2Vk4OzuLX375xWDbRlM8htg2p0+fFi4uLsLR0VE4ODiIqVOnCiGe347ml6yEKkmSJEmS3pnsKRhJkiRJkkyXTEAkSZIkSdI7mYBIkiRJkqR3MgGRJEmSJEnvZAIiSZIkSZLeyQREBzLr5Wdas2YNo0ePztcyQkJCdDpaZ2xsLEuWLMnTvDnjz01gYCBz584FYPLkyc8tirZt2zbOnz+ft0AL0bFjx/jwww/1tr6rV6/yww8/6G19kvQi7t27p46uWqVKlWyjrbq7uxfqum/dukWfPn0KdR2S8TOKSqjFXVpaGr6+vjotkpSZgLz33ns6W2amL7/88rnTt23bRrdu3bC3t9f5uvOjefPmeh32PDMBGThwYJ6fk56ejoWFRSFGJUm5q1ixIqdOnQKUHxilS5dm7Nixell3tWrV2LJli17W9TxpaWlYWsqvQUORR0AK2bVr1+jQoQNOTk506NCB69evAzB06FA++eQT2rVrx/jx47MdNcn8FeLi4kLJkiX5448/uH//Pj179sTJyYlWrVpx5swZQGk4hg8fjre3N3Xr1mXhwoUATJgwgaioKFxcXAgICCAhIYEOHTrQrFkzHB0d8zQqaFBQEI0aNaJjx478888/6uNDhw5VG48JEyZgb2+Pk5MTY8eO5fDhw4SEhBAQEICLiwtRUVGsWLECNzc3nJ2d6d27N48fP1aX8+GHH+Lu7k7dunWzNUizZ8/G0dERZ2dnJkyYAEBUVBRdunTB1dWVtm3bcuHChefGHxoaqpYs1rSdcipdujTjx4/H1dWVjh07Eh4erj4nJCQEUJKGgIAA3NzccHJyYvny5eq2CAsLw8XFhXnz5mmcLzQ0lHbt2jFw4EAcHR159OgRPj4+ODs706RJEzZt2qT1vZGkwpR5VDQ0NBQvLy/8/Pxo2LAhEyZMYP369bRo0QJHR0eioqIApWpm7969cXNzw83NjUOHDj13+VevXlWrea5Zs4Y33niDLl260KBBA8aNG/fM/Pv376dXr17q/b1796plyn/77Tdat25Ns2bN6Nu3LwkJCYDyQ8nNzY0mTZowcuRIdWwgb29vPv/8c7y8vFiwYAGbN2+mSZMmODs7G7Q0frGkk5JpxZy5ublatc7Z2VnUrFlTvP/++0IIIbp16ybWrFkjhBBi5cqVokePHkIIIYYMGSJ8fHxEWlqaEEKI1atXq8/JFBISItq0aSNSUlLE6NGjRWBgoBBCiP379wtnZ2chhBBTpkwRrVu3FklJSSImJkZUqFBBpKSkiCtXrggHBwd1WampqSIuLk4IIURMTIyoV6+eWsmyVKlSz7ymY8eOiSZNmohHjx6JuLg4Ua9ePTFnzhw19s2bN4t79+6Jhg0bqst58OBBtumZ/vvvP/X/iRMnioULF6rz9enTR6Snp4uIiAhRr149IYQQu3btEq1btxaPHj0SQgi10l779u3FxYsXhRBCHDlyRLRr1+6578uBAweEj4/Pc7dTToBa/bNnz56iU6dOIiUlRZw6dUrd5suXLxfTpk0TQgiRlJQkXF1dxeXLl7OtT9t8tra24vLly0IIparviBEj1OfFxsY+93VJkq5NmTJF/XwL8bRNOHDggChXrpy4deuWSEpKEtWqVROTJ08WQggxf/58MWbMGCGEEAMGDBBhYWFCCCGuXbsmXn311eeuL2v7tHr1alGnTh0RGxsrEhMTxSuvvCKuX7+ebf6MjAzRqFEj8e+//6rrCwkJETExMaJt27YiISFBCCHErFmz1IqdWSt0Dho0SK2W7eXlJUaNGqVOa9Kkibh586YQ4mkbJumHPPakAyVLllQPZYKS0R87dgyAv/76i61btwIwePDgbNl93759NR5+j4yMJCAggN9//x0rKysOHjzITz/9BED79u25d+8ecXFxAPj4+KiD81WuXJm7d+8+szwhBJ9//jl//vkn5ubmREdHc/fuXapUqZLr+sPCwujVqxe2trYAuZ4eKlu2LCVKlGDEiBH4+PhoHCDp3LlzTJo0idjYWBISEnjttdfUaT179sTc3Bx7e3s17n379jFs2DB13RUqVCAhIYHDhw/Tt29f9bnJycm5rk+T3LZT1iHmAaytrenSpQsAjo6O2NjYYGVlhaOjI1evXgWUX1xnzpxRj9jExcURGRmJtbV1tmU9b74WLVpQp04ddT1jx45l/PjxdOvWjbZt2+brdUlSYXJzc1OHX69Xrx6dO3cGlP32wIEDgPKZzdrvKz4+nocPH1KmTJk8raNDhw6UK1cOAHt7e65du0bNmjXV6WZmZgwePJjvv/+eYcOG8ddff7Fu3Tr27NnD+fPn8fDwACAlJYXWrVsDcODAAWbPns3jx4+5f/8+Dg4OdO/eHYB+/fqpy/bw8GDo0KH4+fmpR1Uk/ZAJiJ5lHWyuVKlSuc7z6NEj/Pz8WLFihTqCpnjO0OKZowKDMnxyWlraM/OuX7+emJgYjh8/jpWVFbVr1yYpKSnPsebG0tKS8PBw9u/fz8aNG1m0aBG///77M/MNHTqUbdu24ezszJo1awgNDVWnZY098zUKIZ5Zd0ZGBuXLl8+W6OVXXraTlZWVum5zc3P1Oebm5ur8Qgi++eabbIkUkO11aZsv63vfsGFDjh8/zq5du/jss8/o3LkzkydPLvDrlCRdyvq50fSZyMjI4K+//qJkyZIvvA5Nn81hw4bRvXt3SpQoQd++fbG0tEQIQadOndiwYUO2eZOSknjvvfc4duwYNWvWJDAwMFt7l/Xzt2zZMo4ePcovv/yCi4sLp06domLFigV6HVL+yD4ghczd3Z2NGzcCShLQpk0brc8ZNmwYw4YNy/ZL2NPTk/Xr1wPKF5idnR1ly5bVuIwyZcqoo2+C8uu7cuXKWFlZceDAAa5du/bcGDw9Pfn5559JTEzk4cOH7Nix45l5EhISiIuL4/XXX2f+/PlqcpBz3Q8fPqRq1aqkpqaqr+F5OnfuzKpVq9S+Ivfv36ds2bLUqVOHzZs3A8qX++nTpwFldMjPPvtM63J15bXXXmPp0qWkpqYCcPHiRR49evTM69Y0X063bt3C1taWQYMGMXbsWE6cOKGfFyJJOtK5c2cWLVqk3s9sC8LDw3nzzTd1so5q1apRrVo1pk+fztChQwFo1aoVhw4d4tKlSwA8fvyYixcvqsmGnZ0dCQkJz+3wGhUVRcuWLfnyyy+xs7Pjxo0bOolX0k4eASlkCxcuZPjw4cyZM4dKlSqxevXq585/7do1tmzZwsWLF1m1ahUA3377LYGBgQwbNgwnJydsbW3V4ZA1qVixIh4eHjRp0oSuXbsyfvx4unfvTvPmzXFxceHVV1997vObNWtGv379cHFxoVatWrmeFnj48CE9evQgKSkJIQTz5s0DoH///rz99tssXLiQLVu2MG3aNFq2bEmtWrVwdHTM9iWdmy5dunDq1CmaN2+OtbU1r7/+OjNmzGD9+vWMGjWK6dOnk5qaSv/+/XF2diYqKuq5yZiujRgxgqtXr9KsWTOEEFSqVIlt27bh5OSEpaUlzs7ODB06lDFjxuQ6X05nz54lICAAc3NzrKysWLp0qd5eiyTpwsKFC3n//fdxcnIiLS0NT09Pli1bxvXr1wt8VCQ3/v7+xMTEqFfYVapUiTVr1jBgwAD1lOz06dNp2LAhb7/9No6OjtSuXRs3NzeNywwICCAyMhIhBB06dMDZ2Vln8UrPJ0fDlUzeoEGDmDdvHpUqVTJ0KJIkZREQEMDgwYNxcnLSyfJGjx5N06ZNeeutt3SyPMmwZAIiSZIkGT1XV1dKlSrF3r17s/UZkUyXTEAkSZIkSdI72QlVkiRJkiS9kwmIJEmSJEl6JxMQSZIkSZL0TiYgkiRJkiTpnUxAJEmSJEnSO5mASJIkSZKkd/8PRGs6G5R4PM0AAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 576x396 with 4 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "plt.rcParams.update({'legend.fontsize': 6, 'legend.frameon' : False})\n",
    "fig = plt.figure(figsize=(fwid, fhgt), facecolor='w')\n",
    "fig.subplots_adjust(wspace=0.25, hspace=0.25, left=flft, right=frgt, bottom=fbot, top=ftop)\n",
    "# first plot\n",
    "ax = fig.add_subplot(2, 2, 1)\n",
    "# axes limits\n",
    "ax.set_xlim(-1500, 1500)\n",
    "ax.set_ylim(-50, -10)\n",
    "for idx in range(5):\n",
    "    # layer 1\n",
    "    ax.plot(xcell, zeta[idx, 0, 30, :], drawstyle='steps-mid', \n",
    "            linewidth=0.5, color=cc[idx], label='{:2d} years'.format(years[idx]))\n",
    "    # layer 2\n",
    "    ax.plot(xcell, zeta[idx, 1, 30, :], drawstyle='steps-mid',\n",
    "            linewidth=0.5, color=cc[idx], label='_None')\n",
    "ax.plot([-1500, 1500], [-30, -30], color='k', linewidth=1.0)\n",
    "# legend\n",
    "plt.legend(loc='lower left')\n",
    "# axes labels and text\n",
    "ax.set_xlabel('Horizontal distance, in meters')\n",
    "ax.set_ylabel('Elevation, in meters')\n",
    "ax.text(0.025, .55, 'Layer 1', transform=ax.transAxes, va='center', ha='left', size='7')\n",
    "ax.text(0.025, .45, 'Layer 2', transform=ax.transAxes, va='center', ha='left', size='7')\n",
    "ax.text(0.975, .1, 'Recharge conditions', transform=ax.transAxes, va='center', ha='right', size='8')\n",
    "\n",
    "# second plot\n",
    "ax = fig.add_subplot(2, 2, 2)\n",
    "# axes limits\n",
    "ax.set_xlim(-1500, 1500)\n",
    "ax.set_ylim(-50, -10)\n",
    "for idx in range(5, len(years)):\n",
    "    # layer 1\n",
    "    ax.plot(xcell, zeta[idx, 0, 30, :], drawstyle='steps-mid', \n",
    "            linewidth=0.5, color=cc[idx-5], label='{:2d} years'.format(years[idx]))\n",
    "    # layer 2\n",
    "    ax.plot(xcell, zeta[idx, 1, 30, :], drawstyle='steps-mid',\n",
    "            linewidth=0.5, color=cc[idx-5], label='_None')\n",
    "ax.plot([-1500, 1500], [-30, -30], color='k', linewidth=1.0)\n",
    "# legend\n",
    "plt.legend(loc='lower left')\n",
    "# axes labels and text\n",
    "ax.set_xlabel('Horizontal distance, in meters')\n",
    "ax.set_ylabel('Elevation, in meters')\n",
    "ax.text(0.025, .55, 'Layer 1', transform=ax.transAxes, va='center', ha='left', size='7')\n",
    "ax.text(0.025, .45, 'Layer 2', transform=ax.transAxes, va='center', ha='left', size='7')\n",
    "ax.text(0.975, .1, 'Freshwater well withdrawal', transform=ax.transAxes, va='center', ha='right', size='8')\n",
    "\n",
    "# third plot\n",
    "ax = fig.add_subplot(2, 2, 3)\n",
    "# axes limits\n",
    "ax.set_xlim(-1500, 1500)\n",
    "ax.set_ylim(-50, -10)\n",
    "for idx in range(5, len(years)):\n",
    "    # layer 1\n",
    "    ax.plot(xcell, zeta2[idx, 0, 30, :], drawstyle='steps-mid', \n",
    "            linewidth=0.5, color=cc[idx-5], label='{:2d} years'.format(years[idx]))\n",
    "    # layer 2\n",
    "    ax.plot(xcell, zeta2[idx, 1, 30, :], drawstyle='steps-mid',\n",
    "            linewidth=0.5, color=cc[idx-5], label='_None')\n",
    "ax.plot([-1500, 1500], [-30, -30], color='k', linewidth=1.0)\n",
    "# legend\n",
    "plt.legend(loc='lower left')\n",
    "# axes labels and text\n",
    "ax.set_xlabel('Horizontal distance, in meters')\n",
    "ax.set_ylabel('Elevation, in meters')\n",
    "ax.text(0.025, .55, 'Layer 1', transform=ax.transAxes, va='center', ha='left', size='7')\n",
    "ax.text(0.025, .45, 'Layer 2', transform=ax.transAxes, va='center', ha='left', size='7')\n",
    "ax.text(0.975, .1, 'Freshwater and saltwater\\nwell withdrawals', transform=ax.transAxes,\n",
    "        va='center', ha='right', size='8')\n",
    "\n",
    "# fourth plot\n",
    "ax = fig.add_subplot(2, 2, 4)\n",
    "# axes limits\n",
    "ax.set_xlim(0, 30)\n",
    "ax.set_ylim(-50, -10)\n",
    "t = zobs['TOTIM'][999:] / 365 - 200.\n",
    "tz2 = zobs['layer1_001'][999:]\n",
    "tz3 = zobs2['layer1_001'][999:]\n",
    "for i in range(len(t)):\n",
    "    if zobs['layer2_001'][i+999] < -30. - 0.1:\n",
    "        tz2[i] = zobs['layer2_001'][i+999]\n",
    "    if zobs2['layer2_001'][i+999] < 20. - 0.1:\n",
    "        tz3[i] = zobs2['layer2_001'][i+999]\n",
    "ax.plot(t, tz2, linestyle='solid', color='r', linewidth=0.75, label='Freshwater well')\n",
    "ax.plot(t, tz3, linestyle='dotted', color='r', linewidth=0.75, label='Freshwater and saltwater well')\n",
    "ax.plot([0, 30], [-30, -30], 'k', linewidth=1.0, label='_None')\n",
    "# legend\n",
    "leg = plt.legend(loc='lower right', numpoints=1)\n",
    "# axes labels and text\n",
    "ax.set_xlabel('Time, in years')\n",
    "ax.set_ylabel('Elevation, in meters')\n",
    "ax.text(0.025, .55, 'Layer 1', transform=ax.transAxes, va='center', ha='left', size='7')\n",
    "ax.text(0.025, .45, 'Layer 2', transform=ax.transAxes, va='center', ha='left', size='7');"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Use `ModelCrossSection` plotting class and `plot_fill_between()` method to fill between zeta surfaces."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAg8AAADNCAYAAAA/mpM0AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8QZhcZAAAcs0lEQVR4nO3debgkdX3v8ffHWQ7DMIDss4VBRBEIEkBMgjHJE5TlEUfD42XiAqhXNIFcMZgomRszMZm4RMNoFC4Q7gUU5YILIIrAeEUvXCbjoGzDokMYYFhkERg44gwz871//H4NdXr6LHV6q+7+vJ6nn1Nd1V317eqqb3/rV7+qo4jAzMzMbKJe1u0AzMzMrLe4eDAzM7NSXDyYmZlZKS4ezMzMrBQXD2ZmZlaKiwczMzMrxcVDj5B0vaT/2u04zGxwSXqXpGu7HYd1n4uHFpO0VtLzkp6T9KikCyRt1+24+pGkIUnnS7pf0rOSfibp6LrX/ImkuyX9WtIPJe1ZmCZJn5H0ZH58VpI6/0nMQNIbJP0/Sc9I+pWkGyW9rttxFUXExRHx5m7HUVbe9x+XtF7SrZIW1k1/Z84jw5Iul7RTt2LtFS4e2uPYiNgOOAj4HeCMLsczgqQp3Y6hrPxDX7+9TgUeBP4Q2AH4O+BSSQvye3YBvpXH7wSsAv534f0nA28DXgscCLwF+GDbPoTZKCRtD1wF/BtpW50L/AOwoZtxFUma2u0YxjNGjB8GZkfE9qT9/quSZuf37A+cA7wH2B34NXBWB8LtaS4e2igiHgWuIRURwItHy5+T9ICkX0r6H5JmFKYvlHRLrpDvlXRUYZZ75qORZyVdm38ca++7LLd0PCPpx3mHqE27QNLZkr4naRj4Y0k7S/pOXs5PJP2TpBsK79lX0nX5COgeSf+l0WeU9A5JN9eNO13S5eN9Xkkvl3RVPiJ4Kg/PK8zneklLJd1I2qFfUbd+hyNiSUSsjYgtEXEVcB9wSH7JnwKrI+KyiPgNsAR4raR98/QTgc9HxLqIeAj4PHBSo89p1mavAoiIr0fE5oh4PiKujYjbai+Q9AFJd+X9/05JB+fxcyR9M+9H90n6b4X3LJF0qaSL8vtWSzq0MP3jOc/U5vn2wrSTcr45U9KvgCV5XDFP/H7OH8/kv7/f6MNJ+mtJ36wb92+SluXhHZRaER+R9FDOR1PytL0l/R+l1sEnJF0sacfCfNZK+pik24DhRgVERNwWEZtqT4FpwPz8/F3AdyLixxHxHOlg408lzRrj+7KI8KOFD2AtcEQengfcDnyhMH0ZcCXp6GIW8B3gU3naYcAzwJtIhd1cYN887XrgXlKSmZGff7ow3/fl+Q3lZdxSmHZBnu/heb7bAJfkx7bAfqQj+Bvy62fm5+8lHd0fDDwB7N/g8w4BvwJeUxj3M+C4CXzenYHjcgyzgMuAywvzuR54ANg/xzFtnHW/O/Cbwjr7AnB23WvuKMT2DPD6wrRDgWe7vQ35MXgPYHvgSeBC4Gjg5XXT3wE8BLwOEPBKYM+8P98MfAKYTiqw/xM4Mr9vSd4njgGmAJ8CVtTNd06ez/HAMOkIHVIhvQn4y7z/zcjjanliJ+Ap0hH7VODP8vOdG3y+2XneO+bnU4HHgEPy88tJR/8zgd2AlcAH87RXknLiELAr8GNgWWHea4FbSMXAjDHW8VV5XQTwfeBlefwVwMfqXvtcLTY/Rlmf3Q6g3x55Q34OeDZvpD8o7DDKO9Dehdf/HnBfHj4HOHOU+V4P/PfC878Avj/Ka3fMy94hP78AuKgwfQrwAvDqwrh/KiSF44H/WzfPc4C/H2V5ZwNL8/D+OYEMjfd5G8znIOCpus/8yQmu92nAcuCcwrjzKRRYedyNwEl5eDO50MjP98nrTd3ejvwYvAfwmryvriP9aF8J7J6nXQN8uMF7Xg88UDfuDOB/5eElwPLCtP2A58eI4RZgYR4+qcG8TyrkifcAK+um31TbvxrM+2rgA3n4LcCdeXh30umZGYXX/hnww1Hm8zbgZ4Xna4H3TXAdTyMVZx8pjPsB8KG61z0E/FG3t4kqP3zaoj3eFhGzgD8C9gVqpxd2JR1l3yzpaUlPkyrgXfP0+aTWhdE8Whj+NbAdpD4Mkj6dmx/Xk3YmCsuF1JJQsysv9RdoNH1P4PW1GHOc7wL2GCWuC4F3ShIpoVwaERvG+7yStpV0jlJHpfWkI4odNbJPxoOMQ6kvxFeAjcCphUnPkY7oirYnFXaNpm8PPBc5e5h1UkTcFREnRcQ84ABSi8CyPHm03LAnMKduX/1b0g9yTX3e2KbWtC/phHyatPbeAxg9b9SbA9xfN+5+UotpIxcC787D7ybts7XPMA14pBDHOaQWCCTtJumSfDpjPfDVuhjHi/NFEfFCRFwNHCnprXn0eHnCGnDx0EYR8SPSkcTn8qgngOdJzf875scOkTpXQtoB9p7Eot4JLASOIHUcXJDHF68cKP4gPk46splXGDe/MPwg8KNCjDtGxHYR8eeNFh4RK0g/3H+QY6klhfE+7+nAq0mnDrYH3jhO3FvJBcv5pGR5XES8UJi8mtQZsvbamaT1u7rR9Dy8GrMui4i7SbnjgDxqtNzwIKklr7ivzoqIY8ZbhtKVR+eRCu6dI2JH0mm9ie5/D5N++It+i3TU3sjlwIGSDiC1PFxc+AwbgF0Kn2H7iKj12/pUjuPAnCfeXRfjeHE2MpWX1md9nngFqeX05yXnOVBcPLTfMuBNkg6KiC2knfVMSbWqeq6kI/Nrzwfeq3R54cvytH1HmW/RLNLO9yTpSP+fx3pxRGwmXYWwJB/97wucUHjJVcCrJL1H0rT8eJ2k14wx24uALwGbIuKGvJzxPu8sUnHxtNKlUX8/gc9a72xSc++xEfF83bRvAwdIOk7SNqTzwrflxFyL+a9yTHNIxcwFk4jBrClKHZRPr3UYljSf1HS/Ir/k34GPSjpEySvzj/9KYH3uMDgjt0IeoIld4jmT9KP7eF7me3mpWJmI75HyxDslTZV0POm0yFWNXhyp0/I3gK+RTnc8kMc/AlwLfF7S9jn37S3pD/NbZ5FaB56WNBf46xIx1tbt0Xn9TJP0btKByo/ySy4GjpX0B/kA45PAtyLCLQ9jcPHQZhHxOOlH6u/yqI8Ba4AVuQluOenom4hYSeqkeCapM9+P2Lqyb+QiUnPhQ8CdvJRwxnIqqZXiUVJLwdfJl4XlnebNwCLS0cWjwGdI1fhovkJKPF+pGz/q5yUVVjNILRQrSKc0Jiwnzw+S+ko8qnRvjeckvSt/jsdJHTKXkvphvD5/pppzSB04bycdcX03jzPrtGdJ2+d/KF0RtYK0TZ4OEBGXkbbjr+XXXg7slA8EjiXtA/eR9qV/J+3bY4qIO0lXGN0E/BL4bVKfoAmJiCdJLQinkw5c/gZ4S0Q8McbbLszLqc8TJ5A6fN5J2le/QepkCemS1YNJOfG7pAOfMkTq+/EYqVD6MHB8RPw0f47VwIdIRcRjpGLlL0ouY+DIp3cNQNJngD0i4sRJvn8Gacc7OCJ+0dLgzKwvSPot4G5Srlnf7Xhs8tzyMKByU96BuQn0MOD9pGb+yfpz4CcuHMyskdyx+a+AS1w49L6migelGwStlrRFhRuP5GlnSFqjdIOhI0ebh3XNLFLz3zBwKan58orJzEjSWlJT4OmtCs56n/OD1eS+BOtJ92uYTN8mq5hmbzd6B+kufiPOE0vaj3RueX/S5TzLJb0qn5+zCoiIn5BuvtKKeS1oxXys7zg/GJDuBku+tNz6Q1MtD/m65HsaTFpIapraEBH3kTrMHdbMssystzg/mPWvdv2jk7mM7PG/jlFuHCLpZNI/KmHGjBmH7LZDup34C/k+QdNi84jhVkzbNHU6Uzal2wFsnjptxDDQE9OqGtcgfZ5m5iGiZdtzo2mbNIXnn3qCJ5/fWMX/Elqp/LBp6nTA22aVp1U1rm7E3Irc0Yr8MG7xIGk5je8suDgiRjtH3iighpd1RMS5wLkAC2bvHncfm+4Lsmz24QCc9siNI4ZbMe2adyzmDVekG7fdsPC0EcNAT0yralyD9HmamcfMjcMt254bTfvaLq/ly+edT7v1en4I4COFaQFc6/xQuWlVjatTMW+7cZiPjLLNQnfyw7jFQ0QcMYn5rmPkHQvnke4X0HUBXP3FY178Esw6ZduNw3237fV6fqivYqrYTGODp5YrgK3yRVW20XZdqnklsEjp3zHvRfqHQyvbtKxSNMqwWbt523tRZfND7eBi5sZhtt043O1wbIDUCoaZDba7KuaLZi/VfLukdaT/lPhdSdfAi3fsupR0t7DvA6e4J7XZYOnF/OACz7ql17a9pjpMRsS3GeXGQhGxlHQ7VTMbQM4PZv3Ld5g0MzOzUlw8mJmZWSkuHszMzKyUgS4e3LPa2qnWe7q2jfXL5ZmDwvnB2qnX88NAFw+91rvVeovvIdDbnB+snXo9Pwx08WBmZmbluXgwMzOzUlw8mJmZWSkuHszMzKwUFw9mZmZWiosHMzMzK8XFg5mZmZXi4sHMzMxKcfFgZmaVUbvzou/sWW0uHszMrDJ8Z8/e4OLBzMzMSnHxYGZm1gHFUzK99o+w6rl4MDMz64B+OiUztdsBVEXt3+8CLJt9eM9XhdYdtSMLb0P9xfnBbCS3PGS9/u9RrRr66cjCXuL8YDaSiwczMzMrxcWDmZmZleLiwczMzEpx8WBmZmaluHgwMzOzUlw8mJmZWSkuHszMzKwUFw9mZmZWiosHMzMzK8XFg5mZmZXi4sHMzMxKaap4kPQvku6WdJukb0vasTDtDElrJN0j6cjmQzWzXuL8YNa/mm15uA44ICIOBH4OnAEgaT9gEbA/cBRwlqQpTS7LzHqL84NZn2qqeIiIayNiU366ApiXhxcCl0TEhoi4D1gDHNbMssystzg/mPWvVvZ5eB9wdR6eCzxYmLYuj9uKpJMlrZK06oUt0cJwzKxCnB/M+sjU8V4gaTmwR4NJiyPiivyaxcAm4OLa2xq8vuGeHxHnAucCLJi9u7ODWQ9xfjAbTOMWDxFxxFjTJZ0IvAX4k4io7dzrgPmFl80DHp5skGZWTc4PZoOp2astjgI+Brw1In5dmHQlsEjSkKS9gH2Alc0sy8x6i/ODWf8at+VhHF8ChoDrJAGsiIgPRcRqSZcCd5KaK0+JiM1NLsvMeovzgw20bTcOc/UXjwFg2ezDG5+b61FNFQ8R8coxpi0FljYzfzPrXc4PNujqO/c06uzTq3yHSTOzkgK4+ovHMHPjMNtuHO52OGYd5+LBzKwkjTJsNihcPIzCRxY2EbVzmrVtpZ/OaZqZjcbFwyh8ZGET0c/nNM3MRuPiwczMzEpx8WBmZmaluHgwMzOzUlw8mJmZWSkuHszMzKwUFw9mZmZWiosHMzMzK8XFg5mZWYvUbhzX7zeNc/FgZmbWIoNyg0EXD2ZmZlZKU/+S28xs0NX+Dw7AstmH93VTtVmNWx7MzJrg/29ig8jFwwTUjiz8XzYNBqdDlJnZaFw8TICPLKxoUDpEmZmNxsWDmZmZleLiwczMzEpx8WBmZtaEQewH5eLBzMysCYPYD8rFg5mZmZXi4sHMzMxKcfFgZtZCtfvC+J4w1s9cPJiZtdAgnv+2wePiwczMzEpx8TAJbpYcPIN4KZaZNVbLB4OcE1w8TIKbJQePv3Mzq/G/LHDxYGZmZiU1VTxI+kdJt0m6RdK1kuYUpp0haY2keyQd2XyoZtZLnB/M+lezLQ//EhEHRsRBwFXAJwAk7QcsAvYHjgLOkjSlyWWZWW9xfjDrU00VDxGxvvB0JrzYb2QhcElEbIiI+4A1wGHNLMvMeovzg1n/UkRz/UQlLQVOAJ4B/jgiHpf0JWBFRHw1v+Z84OqI+MZY85o7Z0785QfeD8CGfCAyFJtHDLdjWrPz2DRtGwCmbNoAwOapQyOG2zGt3fP35xk5TREt327aOW2jpnDFl/+Vmx9f39W+XM4PW0/bNG2bgd+fxppWv691Yn1NJq5ObzdVyw/jtjxIWi7pjgaPhQARsTgi5gMXA6fW3tZgVg2rFEknS1olaVUvXe4yPH0mIREaxH62gyekF7/zXrG5A7E6PzRWzA/D02c6XzRh5sZhFIGaPNBtVjEHNPpeBy0/jFs8RMQREXFAg8cVdS/9GnBcHl4HzC9Mmwc8PMr8z42IQyPi0OmxmdMeuZHTHrmRodjMUH5eHG7HtLLzmB6bWXz/tWyzZRNDWzZxw8LTuGHhaWyeOsTmqUNbDbdjWrvn788zctrQlk0vfufTO7RdNjNt3sZnmNKBZOv84PzQ8n1tjO+gKjHX8kHte95mAPNDs1db7FN4+lbg7jx8JbBI0pCkvYB9gJXNLKtKfM3/4PF3Xp7zg7eVfjXW9zoo3/nUJt//aUmvBrYA9wMfAoiI1ZIuBe4ENgGnROSTLWY2KJwfrG/U7iq5bPbhA3lHyXpNFQ8RcdwY05YCS5uZv5n1LueHl25lD/hHp8e5RWmkZlseBp6TQ3+qHWWAv1ebPN/GeHze13qTb0/dJCeH/uTv1awzvK/1Jrc8mJmZ1XGLyNhcPJiZWeV1+hSxW0TG5tMWZmYdUvsBnLlxmG03Dnc7nJ7iH/NqcfFgZtYh7rFfbbVTFTM3Dvs0xThcPLSYjyx6lxOH2WBzcTdxLh5azBtf7/J3Z9Y7fKDWXe4waWZmPacVxb6vqJg8Fw9mZl3gG8x1nzthTp5PW5iZdYF/uFqnVohN5DSG+za1houHNiqzQVvn1ZKIE4lZZ7X6B7xR4VWb/8y6/Xy899nEuHhoIx9ZVJu/H7PuaHfnZP/L7PZz8WBmVgG+esB6iYsHM7MK8KXC1ktcPJiZmVkpLh46yM2S3eee1mZmzfN9HjrIzZLd5+/AeoHvAWFV5+LBzKxi+u1KIN/Jsf/4tIWZmbVVvxVD5paHrnGzZGf4iMfMrPXc8tAlrsQ7w+vZ+oE7W1vVuHgwM6s4d/S1qnHxYGZmZqW4eKgIN0u2ju/lYP2sV/7hnvfD/ubioSLcLNk6XpfWz3qlH4/3w/7m4sHMzMxK8aWaFeTLOMvx5Zg2yGr5otvbvvfDweKWhwrqlWbJqvD6skFWldMD3g8Hi4sHM7M+4Y7X1ikuHnqAE8LW3JPbbGudboXwfji4WlI8SPqopJC0S2HcGZLWSLpH0pGtWM6gqkqzZJV4nfQO54fu6MQlnd4PB1fTxYOk+cCbgAcK4/YDFgH7A0cBZ0ma0uyybHBbIWpHOD7K6S3OD93T6Me82dzh/dBqWtHycCbwNzBiO1oIXBIRGyLiPmANcFgLljXwBrXSd2esnuX8UBH1+8xEWyWKpybGm6cNjqaKB0lvBR6KiFvrJs0FHiw8X5fHNZrHyZJWSVr1whbXsWb9wvmhusZqlZhZ17ow3vtsMI17nwdJy4E9GkxaDPwt8OZGb2swruGeHxHnAucCLJi9u7NDCf1+P4jaEU8/frZ+4fzQH8YqClwwWCPjFg8RcUSj8ZJ+G9gLuFUSwDzgp5IOIx1JzC+8fB7wcNPR2gijHT306o/tWDeZcQKrJucHs8E06dMWEXF7ROwWEQsiYgEpIRwcEY8CVwKLJA1J2gvYB1jZkohtVL3eH8L9GvqH84NZf2vL7akjYrWkS4E7gU3AKRGxuR3LssZ65ZSGT00MHucHs97XsptE5SOMJwrPl0bE3hHx6oi4ulXLsYlpx2VarTJa7223NPQv5wez/uJ/jDUgGl2mBZ1plXBfBjOz/uLbUw+gTrRK+NpwM7P+5ZYHG7NVova8ODyRaW5dMDPrX255sBFa9UPvgsHMrH+5eDAzM7NSXDyYmZlZKS4ezMzMrBQXD2ZmZlaKiwczMzMrxcWDmZmZleLiwczMzEpx8WBmZmaluHgwMzOzUlw8mJmZWSkuHszMzKwUFw9mZmZWiosHMzMzK8XFg5mZmZXi4sHMzMxKcfFgZmZmpbh4MDMzs1JcPJiZmVkpLh7MzMysFBcPZmZmVoqLBzMzMyvFxYOZmZmV4uLBzMzMSnHxYGZmZqW4eDAzM7NSXDyYmZlZKS4ezMzMrJSmigdJSyQ9JOmW/DimMO0MSWsk3SPpyOZDNbNe4vxg1r+mtmAeZ0bE54ojJO0HLAL2B+YAyyW9KiI2t2B5ZtY7nB/M+pAiYvJvlpYAzzVIDmcARMSn8vNrgCURcdNY85s3Z0588qTjAVg3fYc0buMzI4bbMa3d83fM/fl5qhoXwGPTtuOys77AzY+vF13i/OBtc5A/T1Xjgtbkh1a0PJwq6QRgFXB6RDwFzAVWFF6zLo/biqSTgZMBpr9MfPm881sQUvMe/81Gdt1merfDAKoTS1XigOrEUpU4YOtY1j77/OSPDFrH+aHNqhJLVeKA6sRSlTig9flh3JYHScuBPRpMWkxKAE8AAfwjMDsi3ifpy8BNEfHVPI/zge9FxDfHWdaqiDi0/MdoPcdS3TigOrFUJQ7oTizOD91XlViqEgdUJ5aqxAGtj2XcloeIOGIiM5J0HnBVfroOmF+YPA94uHR0ZlZpzg9mg6nZqy1mF56+HbgjD18JLJI0JGkvYB9gZTPLMrPe4vxg1r+a7fPwWUkHkZol1wIfBIiI1ZIuBe4ENgGnTLAn9blNxtNKjmVrVYkDqhNLVeKAasUCzg+dUpVYqhIHVCeWqsQBLY6lqastzMzMbPD4DpNmZmZWiosHMzMzK6UyxYOko/KtatdI+ngHlrdW0u35trmr8ridJF0n6Rf578sLr2/Z7XQl/U9Jj0m6ozCu9LIlHZI/wxpJX5RU+oYfo8RS+rbCzcYiab6kH0q6S9JqSR/uxnoZI45urJNtJK2UdGuO5R+6sU6qwPnB+cH5YatYupsfIqLrD2AKcC/wCmA6cCuwX5uXuRbYpW7cZ4GP5+GPA5/Jw/vlmIaAvXKsU5pY9huBg4E7mlk2qYf67wECrgaOblEsS4CPNnht22IBZgMH5+FZwM/z8jq6XsaIoxvrRMB2eXga8B/A73ZrW+nWA+eHSS27Fd/5KLF0Y19wfth63l3ND1VpeTgMWBMR/xkRG4FLgIVdiGMhcGEevhB4W2H8JRGxISLuA9aQYp6UiPgx8Ktmlq10Gdz2EXFTpG//osJ7mo1lNG2LJSIeiYif5uFngbtIdx3s6HoZI47RtHOdREQ8l59Oy4+gS9tKFzk/OD84P2wdS1fzQ1WKh7nAg4Xno96utoUCuFbSzUq3wAXYPSIegbSRALt1ML6yy56bh9sV06mSbsvNlrVmr47EImkB8DukSrpr66UuDujCOpE0RdItwGPAdRHR1XXSJc4P1fvOnR8GPD9UpXhodH6l3deQHh4RBwNHA6dIeuMYr+1GfOMtu50xnQ3sDRwEPAJ8vlOxSNoO+CZwWkSsH+ul7YylQRxdWScRsTkiDiLdhfEwSQeMFXY7Y+ki54fyy3Z+aGMszg/VKR46frvaiHg4/30M+DapmfGXuQmndne8xzoYX9llr8vDLY8pIn6ZN8otwHm81ATb1lgkTSPtkBdHxLfy6I6vl0ZxdGud1ETE08D1wFFUaFvpEOeHCn3nzg/OD1Cd4uEnwD6S9pI0HVhEuoVtW0iaKWlWbRh4M+nWuVcCJ+aXnQhckYc7cTvdUsvOzVHPSvrd3DP2hMJ7mqKStxVuRSz5fecDd0XEvxYmdXS9jBZHl9bJrpJ2zMMzgCOAu6nQttIhzg8V+s6dH5wfgGpcbZH6aHAMqefqvcDiNi/rFaRep7cCq2vLA3YGfgD8Iv/dqfCexTm2e2iypzrwdVLT1gukqu/9k1k2cChpI70X+BL5jqEtiOUrwO3AbXmDm93uWIA3kJrKbgNuyY9jOr1exoijG+vkQOBneZl3AJ+Y7Hbaim2lmw+cH5wfnB/qY+lqfvDtqc3MzKyUqpy2MDMzsx7h4sHMzMxKcfFgZmZmpbh4MDMzs1JcPJiZmVkpLh7MzMysFBcPZmZmVsr/B0E+HV+TiGlAAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 576x198 with 2 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "fig = plt.figure(figsize=(fwid, fhgt/2))\n",
    "fig.subplots_adjust(wspace=0.25, hspace=0.25, left=flft, right=frgt, bottom=fbot, top=ftop)\n",
    "\n",
    "colors = ['#40d3f7', '#F76541']\n",
    "\n",
    "ax = fig.add_subplot(1, 2, 1)\n",
    "modelxsect = flopy.plot.PlotCrossSection(model=ml, line={'Row': 30}, \n",
    "                                         extent=(0, 3050, -50, -10))\n",
    "modelxsect.plot_fill_between(zeta[4, :, :, :], colors=colors, ax=ax,\n",
    "                              edgecolors='none')\n",
    "linecollection = modelxsect.plot_grid(ax=ax)\n",
    "ax.set_title('Recharge year {}'.format(years[4]));\n",
    "\n",
    "ax = fig.add_subplot(1, 2, 2)\n",
    "ax.set_xlim(0, 3050)\n",
    "ax.set_ylim(-50, -10)\n",
    "modelxsect.plot_fill_between(zeta[-1, :, :, :], colors=colors, ax=ax)\n",
    "linecollection = modelxsect.plot_grid(ax=ax)\n",
    "ax.set_title('Scenario year {}'.format(years[-1]));"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "anaconda-cloud": {},
  "kernelspec": {
   "display_name": "Python 3",
   "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.7.3"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 1
}