Week3/original/FittingPeaks/FittingPeaks_original.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Fitting is an Art!\n",
    "\n",
    "## Description\n",
    "Python macro for fitting (possible) peaks (typically signal) on top of a distribution (typically background). It is also an exercise in making ChiSquare and likelihood fits work, and get a feel for which fitting proceedure is likely to give good results.\n",
    "\n",
    "This program contains most the necessary code for fitting with a ChiSquare. The case has been chosen carefully, and should illustrate several points. Please, _play around_ with it as much as you can.\n",
    "\n",
    "\n",
    "## Your Task\n",
    "\n",
    "This exercise is meant to improve your understanding of but also your experience with\n",
    "ChiSquare and Likelihood fits. The below are indicative questions, but you should\n",
    "also \"play around\" yourself to see, if you understand how the three behave, when you\n",
    "change conditions.\n",
    "\n",
    "\n",
    "### Authors\n",
    "- Troels Petersen ([email](mailto:petersen@nbi.dk))\n",
    "- Étienne Bourbeau (notebook conversion) ([email](mailto:etienne.bourbeau@icecube.wisc.edu))\n",
    "\n",
    "\n",
    "### Date\n",
    "1st of December 2019\n",
    "\n",
    "---"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let's load some libraries..."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "from iminuit import Minuit\n",
    "from scipy import stats\n",
    "import os, sys                                         # Modules to see files and folders in directories\n",
    "from os.path import dirname as parent_folder\n",
    "# plt.close('all')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "sys.path.append('../../../External_Functions')\n",
    "from ExternalFunctions import Chi2Regression, BinnedLH, UnbinnedLH\n",
    "from ExternalFunctions import nice_string_output, add_text_to_ax    # Useful functions to print fit results on figure"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Ignore LogWarning (Python may complain, if it takes the log of too small numbers, e.g. log-likelihood):"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "import warnings\n",
    "warnings.filterwarnings(\"ignore\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "r = np.random             # Random generator\n",
    "r.seed(42)                # Set a random (but fixed) seed\n",
    "\n",
    "font_size = 10"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [],
   "source": [
    "Verbose = False\n",
    "Saveplots = False"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "These are switches to decide, if signal should be included in the data, and also fitted for:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [],
   "source": [
    "IncludeSignalInData = False\n",
    "IncludeSignalInFit = False"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "One may want to assert, that one is not fitting for a signal, that is not there!\n",
    "\n",
    "However, in reality this is often what we actually do, to assert what the largest deviation from background is.\n",
    "\n",
    "This type of \"bump hunting\" can be challeging, and very common in many sciences."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {
    "lines_to_next_cell": 2
   },
   "outputs": [],
   "source": [
    "# assert not ((IncludeSignalInData == False) and (IncludeSignalInFit == True))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "---\n",
    "## CASE: Gaussian distribution on an exponential background:\n",
    " * $signal(x) = \\frac{N_{sig}}{\\sigma \\sqrt{2\\pi}} \\cdot \\exp{\\left [-0.5 \\cdot (\\frac{x-\\mu}{\\sigma})^{2}\\right]}$ for $x$ in [0,200]\n",
    " \n",
    " * $bkg(x) = N_{bkg} \\cdot  c \\cdot \\exp{[-cx]}$                   for $x$ in [0,200]\n",
    "---"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [],
   "source": [
    "if IncludeSignalInData:\n",
    "    # Produce signal (Gaussian peak):\n",
    "    Nsigx  =  50\n",
    "    mux    = 125.0\n",
    "    # mux    = r.uniform(10.0, 190.0)    # For signal peak of unknown position\n",
    "    mux_init = mux                       # Initial value to use for the peak position in fit\n",
    "    sigmax =   3.0                       # This is normally not really known either!\n",
    "    x_signal = r.normal(mux, sigmax, Nsigx)\n",
    "    if (Verbose) :\n",
    "        print(\"  Signal: \", x_signal)\n",
    "else:\n",
    "    x_signal = []"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Produce background (here exponential):\n",
    "Nbkgx =  25000\n",
    "taux  =   80.0    # \"Lifetime\" of background\n",
    "x_background = r.exponential(taux, Nbkgx)\n",
    "if (Verbose) :\n",
    "    print(\"  Background: \", x_background)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Histogram specifications:\n",
    "xmin =   0.0\n",
    "xmax = 200.0\n",
    "Nbins =  100       # Note that this gives the bins a natural width, and is much smaller than sigma of the signal!\n",
    "binwidth = (xmax-xmin) / Nbins"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {
    "lines_to_next_cell": 2
   },
   "outputs": [],
   "source": [
    "# Note: The x_signal and x_background needs to be in a parenthesis inside the other\n",
    "# parenthesis, because of the posibility of several other options.\n",
    "x_data = np.concatenate((x_signal, x_background))    # Combine the two numpy arrays\n",
    "if (Verbose) :\n",
    "    print(\"  Data: \", x_data)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "---\n",
    "Define PDFs:\n",
    "---"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Normalized Gaussian:\n",
    "def gauss_pdf(x, mu, sigma):\n",
    "    # return np.exp(-(x - mu)**2 / 2.0 / sigma**2)\n",
    "    return 1.0 / np.sqrt(2 * np.pi) / np.abs(sigma) * np.exp(-(x - mu)**2 / 2.0 / sigma**2)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Normalized Exponential:\n",
    "def exponential_pdf(x, tau):\n",
    "    # return np.exp(-x/tau)\n",
    "    return 1.0 / tau * np.exp(-x/tau)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Next we write up the full model PDF (one for ChiSquare [CS], one for Binned Likelihood [BL], and\n",
    "one for unbinned likelihood [UL]). They are identical, except for the normalisation, where the ChiSquare includes the bin width.\n",
    "\n",
    "We make two versions, one which is a background fit only, and the other which includes the signal PDF in the fit function."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {},
   "outputs": [],
   "source": [
    "if IncludeSignalInFit:\n",
    "    def PDFmodel_CS(x, Nsig, mu, sigma, Nbkg, tau) :\n",
    "        return Nbkg * binwidth * exponential_pdf(x, tau) + Nsig * binwidth * gauss_pdf(x, mu, sigma)\n",
    "        \n",
    "    #def PDFmodel_BL(x, Nsig, mu, sigma, Nbkg, tau) :\n",
    "    #    return ???\n",
    "    \n",
    "    #def PDFmodel_UL(x, Nsig, mu, sigma, Nbkg, tau) :\n",
    "    #    return ???\n",
    "    \n",
    "else:\n",
    "    def PDFmodel_CS(x, Nbkg, tau) :\n",
    "        return Nbkg * binwidth * exponential_pdf(x, tau) \n",
    "        \n",
    "    #def PDFmodel_BL(x, Nbkg, tau) :\n",
    "    #    return ???\n",
    "    \n",
    "    #def PDFmodel_UL(x, Nbkg, tau) :\n",
    "    #    return ???\n",
    "    "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "---\n",
    "Plot the data:\n",
    "---\n",
    "\n",
    "The data from the histogram is here plotted with error bars:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "[Text(0, 0.5, 'Frequency / 2 GeV'),\n",
       " (0.0, 200.0),\n",
       " Text(0.5, 0, 'Photon invariant mass [GeV]'),\n",
       " Text(0.5, 1.0, 'Distribution of diphoton invariant masses')]"
      ]
     },
     "execution_count": 16,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAA4QAAAGDCAYAAABz1N6qAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvnQurowAAIABJREFUeJzt3XmcZHV97//XGwZFEQFlCQIyosQEjSJOdKJxxR0jJtedJGg0aIIaE3PjGE1EE/ObbGr8eaNBMSxBcYsRhRj3JbkZdUBc0YBkkBECAw6L4sLyuX/UaahpqrtOd9fa9Xo+HvPoOqfOOfWpPl0z/Z7vlqpCkiRJkjR7dhp3AZIkSZKk8TAQSpIkSdKMMhBKkiRJ0owyEEqSJEnSjDIQSpIkSdKMMhBKkiRJ0owyEErSGCR5W5I/GdC17p7kB0l2brY/k+QFg7h2c71/TXLsoK63hNf98yRXJvmflsdXkns1j1t/f5OcnOTPV1LrMCR5WJJvj+m1B/bzKUmabHEdQkkarCRbgP2AG4GbgG8CpwInVtXNy7jWC6rqE0s45zPAP1XVO5byWs25JwD3qqpfX+q5g5TkIOC/gIOr6oqW5xRwaFVduMTXOhnYWlWvXkadn2GZ3+tZkGQt8N/ALlV143irkST1YguhJA3Hr1TV7sDBwEbgFcBJg36RJGsGfc0JcTBwVdswqMGZa2mWJM0GA6EkDVFVXVNVZwLPBI5Ncl/YsZtikr2TfCTJ1Um+n+TzSXZKchpwd+DDTZfQP0qytuka+fwk3wU+1bWvOxzeM8kXk1yT5ENJ7tK81iOTbO2uMcmWJI9J8gTgj4FnNq/3leb5W7qgNnW9OsnFSa5IcmqSPZrn5uo4Nsl3m+6er1roe5Nkj+b8bc31Xt1c/zHAx4G7NXWcvMD5/zvJZUkuTfJb857r/v4+MsnWJH/c1LQlyTHzLrdXkrOSXJfkC0nu2XWthyT5UvO9/FKShzT7Xw88DHhLU+dbFju+63v5Z0n+o3mtjyXZe4H3t8O9aur+wyRfba79niS7Ns+dn+TJXceuad7rEc32+5L8T3Pe55LcZ9736q1Jzk7yQ+BR875/ezU/n9uSbG8eH9jyPX2u+Xp18z36pR7v84Smvn9qzv9akp9N8srmZ+ySJI/rOv55zfu9LslFSV7Y9VzPz1Lz3CuSfK8579tJjmz275RkQ5LvJLkqyXtz6+dl16auq5prfinJfr3ulyRNKwOhJI1AVX0R2EonQMz38ua5feh0Nf3jzin1G8B36bQ23qmq/qrrnEcAPw88foGX/E3gt4C70em6+uYWNX4U+AvgPc3r3b/HYc9t/jwKOAS4E/CWecf8MnBv4EjgT5P8/AIv+f8DezTXeURT8/Oa7rFPBC5t6nju/BPTCa9/CDwWOBR4TJ+39zPA3sABwLHAiUnu3fX8s4HXAnsBFwKvb17nLsBZdL5/dwXeAJyV5K5V9Srg88CLmzpfvNjxXa/1HOB5wL7A7Zr30dYzgCcA9wDuR+deALy7eQ9zHg9cWVXnNtv/Suf7tC9wLnD6vOs+p3nPuwP/Pu+5nYB/pNNqe3fgR9z2ni/0nh7efN2z+R795wLv61eA0+h8/78M/FvzugcArwP+oevYK4AnA3duXvONc8GXBT5Lzb1+MfCLTcv944EtzTkvBZ5K52fwbsB24P80zx1L52f0IDr380XN+5ekVcNAKEmjcylwlx77bwD2pzNe7oaq+nz1H+B9QlX9sKoW+uX0tKr6elX9EPgT4BkZTFfAY4A3VNVFVfUD4JXAs7Jj6+Rrq+pHVfUV4CvAbYJlU8szgVdW1XVVtQX4W+A3WtbxDOAfu97jCS3O+ZOq+klVfZZOaHtG13P/XFVfbMa5nQ4c3uw/Crigqk6rqhur6t3At+gEmF7aHP+PVfVfzb17b9drtfHmqrq0qr4PfLjr3HcBT0lyx2b7Oc0+AKrqnc33+Sd0vlf3T9Oy2/hQVf1HVd1cVT/ufsGquqqqPlBV11fVdXSC4yPm1bWS9wTw+ar6t+b7/z46gW5jVd0AnAGsTbJnU89ZVfWd6vgs8DFu/Y+WhT5LNwG3Bw5LsktVbamq7zTnvBB4VVVt7fr+PK35mb6BThC8V1XdVFXnVNW1S3xvkjTRDISSNDoHAN/vsf+v6bRKfazpArehxbUuWcLzFwO70GkhW6m7NdfrvvYaOq0xc7pnBb2eTivifHvTaUmaf60DllDH/Pe4mO1NcOw+/m5d2wvVPP/99quzzfFtvj8L6XluM5HO+cCvNKHwKTSBMMnOSTY2XSKv5daWse6fhwV/npLcMck/pNOt91o63UD3nPcfDCt5TwCXdz3+EZ3WzZu6tpm7ZpInJtnUdAm9GnhS13vp+Vlqvj8voxP2rkhyRpK5+38w8MGmS+jVdL6PN9H5mT6NTmvlGel0Tf6rJLss8b1J0kQzEErSCCT5RTqhYH53PJqWm5dX1SF0WpL+YG58E7BQS2G/FsSDuh7fnU5Lx5XAD4G5VqS5lrp9lnDdS+n8At197RvZ8Rf6Nq5sapp/re+1PP8ybvseF7NXkt3mHX9pi9eZ/37nzp2rc/73q9/xwzTXbfRo4Jtds60+p9n3GDrdH9c2+9N17mL3/eV0ugA/uKruzK3dQLPwKa2uu2RJbg98APgbYL+q2hM4e66WxT5LVfWuqvplOvengL9sLnsJ8MSq2rPrz65V9b2mlfG1VXUY8BA6XVV/c5DvSZLGzUAoSUOU5M7NZB9n0Fme4Gs9jnlyknslCXAtndaJudaRy+mMsVuqX09yWNNa9Drg/U2Ly38BuyY5qmnpeDWdrnRzLqfTPW+hfx/eDfx+knskuRO3jjlc0pICTS3vBV6fZPckBwN/APxTy0u8F3hu13t8TYtzXpvkdkkeRucX+/e1OOds4GeTPCediVqeCRwGfKR5fv796Xf8MJ0BPA74Hbq6i9IZF/gT4Co6/xnwF0u87u50WumubsZItvlez9kG3MzyfoZ7uR2dn9dtwI1JnkjnPQMLf5aS3DvJo5tA+WM672fuM/Y2Oj+HBzfX2CfJ0c3jRyX5heY/Tq6l858Yc+dJ0qpgIJSk4fhwkuvotD68is7kIs9b4NhDgU8APwD+E/j7qvpM89z/B7y66c62lMlHTgNOptOVb1c6E2dQVdcAvwu8g06r1Q/pTMIxZy4kXZXkXG7rnc21P0dnfbkfAy9ZQl3dXtK8/kV0Wk7f1Vy/r6r6V+BNwKfodBH8VJ9T/ofOZCGX0hkj+KKq+laL17mKTnh8OZ1A9UfAk6vqyuaQv6Mz3mx7kje3OH5oquoyOj8/DwHe0/XUqXS6rX6PzpqYm5Z46TcBd6DTqrsJ+OgSarqezpjD/2h+htcv8bXnX+86Oj/L76VzP58DnNl1yEKfpdvTWf7lSjo/C/vSmXAGOvfwTDrdTK+j8x4f3Dz3M8D76YTB84HP0v4/LSRpKrgwvSRpVUvySDqtswf2O1aSpFljC6EkSZIkzSgDoSRJkiTNKLuMSpIkSdKMsoVQkiRJkmaUgVCSJEmSZtSacRewEnvvvXetXbt23GVIkiRJ0licc845V1bVPss9f6oD4dq1a9m8efO4y5AkSZKksUhy8UrOt8uoJEmSJM0oA6EkSZIkzSgDoSRJkiTNKAOhJEmSJM0oA6EkSZIkzSgDoSRJkiTNKAOhJEmSJM0oA6EkSZIkzSgDoSRJkiTNKAOhJEmSJM0oA6EkSZIkzSgDoSRJkiTNqFUTCNduOIu1G84a+zUkSZIkaVqsmkAoSZIkSVoaA6EkSZIkzSgDoSRJkiTNKAOhJEmSJM0oA6EkSZIkzSgDoSRJkiTNKAOhJEmSJM0oA6EkSZIkzSgDoSRJkiTNKAPhEqzdcBZrN5w17jIkSZIkaSBWXSA85+LtE3ENSZIkSZp0qyIQdge4Y96xaVmBbhDXkCRJkqRpsioC4aaLrrrl8Q033rzD9iivIUmSJEnTZFUEwvWH3PWWx7us2WmH7VFeQ5IkSZKmyaoIhA88eK9bHp/+gvU7bI/yGpIkSZI0TVZFIOw2iCBnGJQkSZI0C1ZdIJQkSZIktWMglCRJkqQZNdRAmGTPJO9P8q0k5yf5pSR3SfLxJBc0X/dqjk2SNye5MMlXkxwxzNokSZIkadYNu4Xw74CPVtXPAfcHzgc2AJ+sqkOBTzbbAE8EDm3+HAe8dci1LZtrFEqSJElaDYYWCJPcGXg4cBJAVf20qq4GjgZOaQ47BXhq8/ho4NTq2ATsmWT/YdW3VC5cL0mSJGm1GWYL4SHANuAfk3w5yTuS7AbsV1WXATRf922OPwC4pOv8rc2+HSQ5LsnmJJu3bds2xPJ35ML1kiRJklabYQbCNcARwFur6gHAD7m1e2gv6bGvbrOj6sSqWldV6/bZZ5/BVNqCC9dLkiRJWm2GGQi3Alur6gvN9vvpBMTL57qCNl+v6Dr+oK7zDwQuHWRBazecxdoNZy3rXBeulyRJkrTarBnWhavqf5JckuTeVfVt4Ejgm82fY4GNzdcPNaecCbw4yRnAg4Fr5rqWtrFl41ErrrntNQyDkiRJklaDoQXCxkuA05PcDrgIeB6dVsn3Jnk+8F3g6c2xZwNPAi4Erm+OlSRJkiQNyVADYVWdB6zr8dSRPY4t4Phh1iNJkiRJutWw1yGUJEmSJE0oA6EkSZIkzSgDoSRJkiTNKAOhJEmSJM0oA6EkSZIkzaiZDITnXLx93CVIkiRJ0tjNTCDsDoHHvGPTskLhlo1HtV68XpIkSZIm3cwEwk0XXXXL4xtuvHmHbUmSJEmaRTMTCNcfctdbHu+yZqcdtiVJkiRpFs1MIHzgwXvd8vj0F6zfYVuSJEmSZtHMBMJuhkFJkiRJmtFAKEmSJEkyEEqSJEnSzDIQSpIkSdKMMhAO2NoNZ7F2w1njLkOSJEmS+jIQSpIkSdKMMhBKkiRJ0owyEEqSJEnSjDIQSpIkSdKMWjPuAkZpy8ajxl2CJEmSJE0MWwinkDOZSpIkSRoEA6EkSZIkzSgDoSRJkiTNKAPhkJxz8fae++3uKUmSJGlSGAgHqDsEHvOOTQuGQkmSJEmaBAbCAdp00VW3PL7hxpt32JYkSZKkSWMgHKD1h9z1lse7rNlph21JkiRJmjQGwgF64MF73fL49Bes32FbkiRJkiaNgXBIxhkGnbhGkiRJUhsGwinmpDWSJEmSVsJAOGWcyVSSJEnSoBgIp4wzmUqSJEkaFAPhmCy3Zc+ZTCVJkiQNioFwhAbR3dOZTCVJkiQNioFwhNp091zKDKGGQUmSJEkrYSAcIbt7SpIkSZokBsIRsrunJEmSpEmyZtwFrDZbNh7V6rhRhMFzLt5u6JQkSZK0IFsIVxnXKZQkSZLUloFwlXGdQkmSJEltDTUQJtmS5GtJzkuyudl3lyQfT3JB83WvZn+SvDnJhUm+muSIYda2Wg1q4pqlzHYqSZIkaTqNooXwUVV1eFWta7Y3AJ+sqkOBTzbbAE8EDm3+HAe8dQS1rTptJ64x8EmSJEkax6QyRwOPbB6fAnwGeEWz/9SqKmBTkj2T7F9Vl42hxok2SRPXSJIkSZpew24hLOBjSc5Jclyzb7+5kNd83bfZfwBwSde5W5t9O0hyXJLNSTZv27ZtiKWP16RPBmMLoyRJkjT9hh0IH1pVR9DpDnp8kocvcmx67Kvb7Kg6sarWVdW6ffbZZ1B1TgRnCJUkSZI0SkPtMlpVlzZfr0jyQeBBwOVzXUGT7A9c0Ry+FTio6/QDgUuHWd84LNbds9cMoXb7lCRJkjQsQ2shTLJbkt3nHgOPA74OnAkc2xx2LPCh5vGZwG82s42uB66ZtfGDg5ohVJIkSZLaGGYL4X7AB5PMvc67quqjSb4EvDfJ84HvAk9vjj8beBJwIXA98Lwh1jaR2s4QKkmSJEmDMLRAWFUXAffvsf8q4Mge+ws4flj1TBvDoCRJkqRhG8U6hJpiTmwjSZIkrV4GwlVoy8ajWq9V2CvwOdupJEmSNBsMhDOoX+DrNdtpm2tJkiRJmi4GwhnUL/D1m+3UFkRJkiRpdTAQzqB+ga/fbKdLaUGUJEmSNLkMhDNoKctb9HrO9RIlSZKk1cFAOOOWs7yF6yVKkiRJq8MwF6bXMrSdHXRSGAYlSZKk6WULoSRJkiTNKAOhJEmSJM0oA6EkSZIkzSgDoSRJkiTNKCeVUU/TNrmNJEmSpKUzEM4oA58kSZIku4xqaNZuOIu1G84adxmSJEmSFmALoZbFFkZJkiRp+tlCKEmSJEkzykAoSZIkSTPKQKixGcQYQ8cpSpIkSctnIJQkSZKkGWUglCRJkqQZZSDUxLJLqSRJkjRcBkJJkiRJmlEGQg3dORdvH3cJkiRJknowEGooukPgMe/YZCiUJEmSJpCBUEOx6aKrbnl8w40377AtSZIkaTIYCDUU6w+56y2Pd1mz0w7bw2ALpCRJkrR0CwbCJOtGWYhWlwcevNctj09/wfodtudbbpizW6okSZK0Mou1EL49yQVJXpfksJFVpFWnVxgcRJizW6okSZK0MgsGwqp6APBk4Cbg/UnOS/KKJAePrDqtWksJcwuFxaV0S7X1UJIkSbqtRccQVtW3q+q1VXUYcCywJ/CpJP8xkuq0avULc21aEPt1S7VLqSRJkrS4VpPKJNkJ2BfYD9gN2DbMorT69QtzS+0O2qtbattrrN1wFms3nNWqbkmSJGk1WbPYk0keBjwbeCrwdeAM4Per6poR1KYZ0SvMDWKW0lHPdCpJkiRNmwUDYZJLgO/SCYGvrarLR1aVZt5SZikd5jUkSZKk1WyxFsJfrqqL5zaS7FZVPxxBTdIOBhHkDIOSJEnSbS0YCOfCYJJfAk4C7gTcPcn9gRdW1e+OpkRNqy0bjxp3CZIkSZIW0WZSmTcBjweuAqiqrwAPH2ZRkiRJkqThazXLaFVdMm/XTUOoRZIkSZI0QovOMtq4JMlDgEpyO+ClwPnDLUtqZ1q6pc4tazEt9UqSJGk2tGkhfBFwPHAAsBU4vNluJcnOSb6c5CPN9j2SfCHJBUne04RMkty+2b6weX7tUt+MpsuWjUdNVEBy4XpJkiTNmr6BsKqurKpjqmq/qtq3qn69qhZfJXxHv8eOLYp/Cbyxqg4FtgPPb/Y/H9heVfcC3tgcJw1Vdwg85h2beoZCF66XJEnSarXYOoT3Ae5ZVWc2228E9miefktVndvv4kkOBI4CXg/8QZIAjwae0xxyCnAC8Fbg6OYxwPuBtyRJVdUS35NWiUG0Hva7xqaLbv2/jRtuvJlNF13lEhWSJEmaGYu1EG4EruzafjxwFvBp4E9bXv9NwB8BNzfbdwWurqobm+2tdLqi0ny9BKB5/prmeGlo1h9y64/YLmt22mF7kthKKUmSpGFYLBDuX1X/t2v72qr6QFWdBuzd78JJngxcUVXndO/ucWi1eK77uscl2Zxk87Zt2/qVIS2quzXw9Best3VQkiRJM2WxQLh790ZVre/a3LfFtR8KPCXJFuAMOl1F3wTsmWSuq+qBwKXN463AQQDN83sA359/0ao6sarWVdW6ffbZp0UZUjuGQUmSJM2axQLhpUkePH9nkvXcGuIWVFWvrKoDq2ot8CzgU1V1DJ0up09rDjsW+FDz+Mxmm+b5Tzl+ULPCLqGSJEkah8XWIXwF8J4kJwNzE8g8kE5oe+YKXvMVwBlJ/hz4MnBSs/8k4LQkF9JpGXzWCl5DmkjnXLzdlkhJkiRNjAUDYVV9sWkhfDHw3Gb3N4D1VXX5Ul6kqj4DfKZ5fBHwoB7H/Bh4+lKuK02D+UtbOFZRkiRJk2KxFkKq6grazygqrWrLbd1zaQtJkiRNqr4L00uzbBAL1w9yaYtery9JkiQtl4FQM2/LxqMWXMC+V+veUq10aYs2oVSSJElaDgOhtIhBL1w/qC6nkiRJ0iAsGAiT7JzkhUn+LMlD5z336uGXJo3fJCxcP+hQKkmSJM1ZrIXwH4BHAFcBb07yhq7nfm2oVUkTaBRhsFd30EkIpZIkSVqdFguED6qq51TVm4AHA3dK8s9Jbg9kNOVJq99SxggaBiVJkjRIiwXC2809qKobq+o44DzgU8Cdhl2YNCtGNUaw32yokiRJmj2LBcLNSZ7QvaOqXgf8I7B2mEVJs8QxgpIkSRqXBQNhVf16VX20x/53VNUuwy1Lmj7LXQ7CMYKSJEkaF5edkFag7fi/xdY67GYYlCRJ0iitGXcB0qRbLMj1Gv9nqJMkSdK0WGwdQsOi1Ifj/yRJkjTNFgt9m5JsBT4KfLSqtoymJGl6jGr8X5vuppIkSdJSLRgIq2pdkoOBJwJvSnIA8O/AvwKfraqfjKhGaSrYVVSSJEnTZtFJZarq4qp6W1U9FXgI8GHgMcDnk7igmTRB2q4zuNzZUAfF9RAlSZImR+txglV1A51F6T8F0LQYShqAYXcJnT8bqstbSJIkCVaw7ERVfW+QhUganl6zoUqSJEmuQyjNAGdDlSRJUi99A2GS+46iEEnDM6rZUB0fKEmSNF3ajCF8W5LbAScD76qqq4dbkjRdpm1JCMcOSpIkaU7fFsKq+mXgGOAgYHOSdyV57NArkzRytvBJkiTNllZjCKvqAuDVwCuARwBvTvKtJL82zOIkLd24l5Voa1rqlCRJWs3ajCG8X5I3AucDjwZ+pap+vnn8xiHXJ6mF+ctKTGrYmpY6JUmSZkWbFsK3AOcC96+q46vqXICqupROq6GkMZuWZSWmpU5JkqRZ0WZSmScBP6qqmwCS7ATsWlXXV9VpQ61OUivTsqzEtNQpSZI0K9q0EH4CuEPX9h2bfZImRJtlJbZsPGpkM6Iu1BV0VMtfSJIkqZ02gXDXqvrB3Ebz+I7DK0nSSowrZC11fKBhUJIkafzaBMIfJjlibiPJA4EfDa8kSeO2nMleHB8oSZI0fdoEwpcB70vy+SSfB94DvHi4ZUkatZXOAOr4QEmSpOnTd1KZqvpSkp8D7g0E+FZV3TD0yiSNVK8WvqV063R8oCRJ0vRpM8sowC8Ca5vjH5CEqjp1aFVJGrlBtvAZBiVJkqZD30CY5DTgnsB5wE3N7gIMhNIqsppa+NZuOAtgZLOqSpIkTas2LYTrgMOqqoZdjKTlG2T4GWYYNKRJkiRNjjaTynwd+JlhFyJJkiRJGq02LYR7A99M8kXgJ3M7q+opQ6tK0kSyK6YkSdLq0iYQnjDsIiRJkiRJo9dm2YnPJjkYOLSqPpHkjsDOwy9N0rSx5VCSJGm6tJll9LeB44C70Jlt9ADgbcCRwy1N0qgZ6CRJkmZLm0lljgceClwLUFUXAPsOsyhJs2vthrNuGau4UudcvH0g15EkSVqt2gTCn1TVT+c2kqyhsw7hopLsmuSLSb6S5BtJXtvsv0eSLyS5IMl7ktyu2X/7ZvvC5vm1y3tLkmZZdwg85h2bDIWSJEmLaBMIP5vkj4E7JHks8D7gwy3O+wnw6Kq6P3A48IQk64G/BN5YVYcC24HnN8c/H9heVfcC3tgcJ2kCTXLI2nTRVbc8vuHGm3fYliRJ0o7aBMINwDbga8ALgbOBV/c7qTp+0Gzu0vwp4NHA+5v9pwBPbR4f3WzTPH9kkrSoT9IITEvL2/pD7nrL413W7LTDtiRJknbUNxBW1c1V9faqenpVPa153LfLKECSnZOcB1wBfBz4DnB1Vd3YHLKVziQ1NF8vaV7zRuAa4Da/ySU5LsnmJJu3bdvWpgxJAzApLW/9xhg+8OC9bnl8+gvW77AtSZKkHbWZZfS/6TFmsKoO6XduVd0EHJ5kT+CDwM/3OmzupRZ5rvuaJwInAqxbt65VMJW0ctPY8mYYlCRJWlybhenXdT3eFXg6nSUoWquqq5N8BlgP7JlkTdMKeCBwaXPYVuAgYGszcc0ewPeX8jqShmfULW/nXLzdQCdJkjRkbbqMXtX153tV9SY64wAXlWSfpmWQJHcAHgOcD3waeFpz2LHAh5rHZzbbNM9/qm3XVEmjNaygNi3jFCVJklaLNl1Gj+ja3IlOi+HuLa69P3BKkp2b895bVR9J8k3gjCR/DnwZOKk5/iTgtCQX0mkZfFb7tyFpNeg1TtFWQkmSpOFp02X0b7se3whsAZ7R76Sq+irwgB77LwIe1GP/j+l0R5U0o6ZpnOLcxDZbNh411a8hSZJmW99AWFWPGkUhkrSUcYqLjTE0QEmSJLXTpsvoHyz2fFW9YXDlSFJHr7A3f4yhy0pIkiStTJuF6dcBv0NnncADgBcBh9EZR9hmLKEkDcSg1kLst5ahJEnSrGgzhnBv4Iiqug4gyQnA+6rqBcMsTNLkGXdXzGkaYyhJkjQN2rQQ3h34adf2T4G1Q6lGkhYx6rUQJUmSVrs2gfA04ItJTkjyGuALwKnDLUuSFjcpYXAUayW6HqMkSRqWNgvTvx54HrAduBp4XlX9xbALkzSbtmw8auxdU/uZP7nNMALbKF5DkiSpTQshwB2Ba6vq74CtSe4xxJokaaINanKbcb+GJElS30DYdBN9BfDKZtcuwD8NsyhJGoWFWt36zUI6isltnEBHkiSNQpsWwl8FngL8EKCqLsXlJiRNqUF0xRzF5DZOoCNJkkahzbITP62qSlIASXYbck2StKCVji/s1RVzJWFroXPnWhgHMR7SMChJkoalTQvhe5P8A7Bnkt8GPgG8fbhlSdJwTFJXzH5dUyVJkoatbwthVf1NkscC1wL3Bv60qj4+9MokaQjsiilJknSrRQNhkp2Bf6uqxwCGQEmrimFQkiTNukW7jFbVTcD1SfYYUT2SNDFc+0+SJK12bSaV+THwtSQfp5lpFKCqXjq0qiRpTObPQrpQt9JBTBbTzyBeY5CT20iSpNWnTSA8q/kjSaveoGchPefi7XZNlSRJE2vBQJjk7lX13ao6ZZQFSdI4DWIW0ratjN3HGxolSdI4LDaG8F/mHiT5wAhqkaSR2LLxqAW7UA5iFtJerYzzzQ+NjleUJEnjsFggTNfjQ4ZdiCRNmuW22rVpZWwTGiVJkoZtsUBYCzyWJC2iTSvjILqmSpIkrdRigfD+Sa5Nch1wv+bxtUmuS3LtqAqUpGm2UCvjILqmrt1w1i2ziPZjl1RJktTLgoGwqnaJXGsFAAAX8ElEQVSuqjtX1e5VtaZ5PLd951EWKUmr2bAmlBnEOMWlhE5JkjR92iw7IUkzZbWs2TfoJTQkSdLqs1iXUUnSFHOcoiRJ6scWQkkagkloZRzEOEVJkrS62UIoSVOuzdhAw6AkSerFQChJY7Jl41HLbklcbQvbO3mNJEnjYSCUpCk06oXtpz1wSpKk3gyEkjSFRjFhzGprhZQkSbdlIJSkKTSKCWPatkLa3VOSpOnlLKOSNOUWC4Mrme10UK2Qc2FxEmZelSRJO7KFUJLU0ziWrVioW6qtkJIkDYeBUJLU1zDDoGMVJUkaHwOhJGmsRj1jqiRJupWBUJJm2CC7Yi63ZW8UM6bCaLqd2rVVkjRtDISSNKVWsrD9oCylu+dCz41jrKIkSeowEEqSFtQvdPbr7rnU8YH9wqDjCyVJGiwDoSRp2fp19xzE+EAnnZEkaXgMhJKkZYesft09BzE+cNCTzhgoJUm61dACYZKDknw6yflJvpHk95r9d0ny8SQXNF/3avYnyZuTXJjkq0mOGFZtkqTBt7z16u45iPGBgwiVo25lNHRKkqbFMFsIbwReXlU/D6wHjk9yGLAB+GRVHQp8stkGeCJwaPPnOOCtQ6xNkmbeqJd7WCwMLjZWcRChchTv1a6tkqRpNLRAWFWXVdW5zePrgPOBA4CjgVOaw04Bnto8Pho4tTo2AXsm2X9Y9UnSrBvVcg+DtNwZSEfxXl1PUZI0jUYyhjDJWuABwBeA/arqMuiERmDf5rADgEu6Ttva7Jt/reOSbE6yedu2bcMsW5JWtdW03EO/9f9G8V6nMWBLkrRm2C+Q5E7AB4CXVdW1SRY8tMe+us2OqhOBEwHWrVt3m+clSUu3koA07rUQl2pYwXc1BWxJ0uwYaiBMsgudMHh6Vf1zs/vyJPtX1WVNl9Armv1bgYO6Tj8QuHSY9UmShm/aAuMgGAYlSdNiaIEwnabAk4Dzq+oNXU+dCRwLbGy+fqhr/4uTnAE8GLhmrmupJGm2TUKonOuSOgm1SJI0KMMcQ/hQ4DeARyc5r/nzJDpB8LFJLgAe22wDnA1cBFwIvB343SHWJklahZzZU5KkpRlaC2FV/Tu9xwUCHNnj+AKOH1Y9kqTbWg2tXfOXe1ho/N4kvFdbGSVJk2bok8pIkjRMvZZ7GNcYPoOeJGnajGTZCUmShmXUyz3YLVWStJoYCCVJU20Uyz3M75Y6zFDYb01FSZIGyUAoSVo1htVVtFe31JWwlVGSNCkMhJIksXjL3CC6pY6ylVGSpLYMhJIk9TGIbqmDbmWUJGkQDISSJC3BcruljnryG0mS2nDZCUnS1JuG5R6W2sp4zsXbx7Z8Rj+upyhJq4cthJIkdRnF2L6Fgp7jDCVJo2YglCTNvEkJYqMaZzgJS1tMQg2SJLuMSpLUM4jNb8UbRfdIxxlKkkbNFkJJ0syblCDWZpzhJLWs2aVVkqafgVCSNPMGsaxEG1s2HtW6pXFSJ5SZlO61kqTBMBBKktRlUoPYMCwnzA16nKOBUpLGy0AoSdIq0q9LaZsWvsWuMYjutbYyStLkMBBKkjRB2nYrXW6IWmkL3yC6145qNlVJUn8GQkmSpsQgWtYGOYHOcrvXTsokPpIkl52QJAkYzbISK9VmeYx+RjWBzqTXIEnqsIVQkqQpsZSWtTath/2C2CjG9i1UwyQtryFJq5mBUJKkKdGvZW0QXUqd8EWSZotdRiVJmkK9WtYG0aW0zTWmoXutJKkdA6EkSatE2y6liwW6UU34MohQOdel1IAqSctnl1FJklaJQUzWMmkTvthlVZKGy0AoSdIUabtO4SCC3EqusZJJYQY1jtGJaSSpPwOhJEmaKC5cL0mjYyCUJEkTZanjGFfSrdRWREmzzkllJElaRQYxwcogJ2k55+LtS+562mYc4/xupZMw3lGSppEthJIkaaAGuZbhQiFvKd1KnZhGkhZmIJQkSQM1ijGA/bqVLjWUOnGNpFllIJQkSQM1irUM+3UrbRNKB9mSOUyDCJ0GV0kLcQyhJEkaqEGsZbiUcYy9rt8mlPYKjfOvNReiVjKuchDXkKRhsYVQkiQNzbgmemkTSgfZkulMp5KmlYFQkiStaguF0pW2ZI66y+kgrj+p3WIljY+BUJIkTaUtG48aWDfMfmGwV5Ba6uQ5ywljgwid0zJWUtJ4GAglSdLADTKsjauGfkGqTZfTpYSxQYTOXkYx66uk6WUglCRJ6qFfkGrT5bTfNQYROvsZxayvkqaXs4xKkqSZtVgL4lKC1EJdTvtdo99Mp4OYsXUQ15C0etlCKEmS1MMowtggQudy65nP2U6l2WQLoSRJUh8LBamVrpc4Ta13rqcorU5DC4RJ3gk8Gbiiqu7b7LsL8B5gLbAFeEZVbU8S4O+AJwHXA8+tqnOHVZskSdKkWUno7BfWJiXEGSqlyTPMLqMnA0+Yt28D8MmqOhT4ZLMN8ETg0ObPccBbh1iXJElSK5MwW+qoLXdZilF0ObVbqzR4QwuEVfU54Pvzdh8NnNI8PgV4atf+U6tjE7Bnkv2HVZskSdKo9AuVkxA6V7o8hqTpNepJZfarqssAmq/7NvsPAC7pOm5rs0+SJEktLTesrXR5jKUyVEqTY1JmGU2PfdXzwOS4JJuTbN62bduQy5IkSZpsgwhry1keo009S62zbZfQlQTKfq8xiG6pdm3VNBl1ILx8rito8/WKZv9W4KCu4w4ELu11gao6sarWVdW6ffbZZ6jFSpIkTbqlhLWFrHR5jDZhb6V1DrqVctwMjZoUow6EZwLHNo+PBT7Utf8307EeuGaua6kkSZIWtpS1DNtYzvIYbcLeSutsGyhH0co4CgZGjcrQAmGSdwP/Cdw7ydYkzwc2Ao9NcgHw2GYb4GzgIuBC4O3A7w6rLkmSpNVk1GsZ9rp+m7C3lDp7hbVBBN9RT54z6aFTguHOMvrsqtq/qnapqgOr6qSquqqqjqyqQ5uv32+Orao6vqruWVW/UFWbh1WXJEnSarWSMLiS2VCXGkp7Pd8vrC31NZbTbXUQ3VKXeo2VPi+t1KRMKiNJkqRVYLmhdCljDBd6jX5hbJCT5yykzTX61bnaxkv2Y/fY8Voz7gIkSZK0MuNex7BtDYsd06ZLaL/X6BXGusPjSifPaaPNNfrV2e/5STIX5CbhZ1DLYwuhJEmSxm4QYyGXEuiWM3nOnMVatNpco1+dg54oSFqMgVCSJEkTZbmtYYOcYGcQLXILXaNfnaOeKGgQnIRnehkIJUmStOosFKL6TZ6zFIMIMMuZgGdU+o3tG8ckPBo8xxBKkiRpIkzCOLR+NcwPML1a8FY6nrKNNmP3hj2+bxBjHadpvORqZSCUJEnSqjHsUDmqADMJ4XjOORdvX/b6k/0M4hpObLMydhmVJEmSWpq0CV/adLEcVlfONmMd+3U7ncbxkquNgVCSJElqaRICTJuwttKxeUtdk3GYk/C05fjD5TEQSpIkScswrtasNmGtzTGLtd4NuiV0sbC2kol+nJRm5QyEkiRJ0hRpE9ZWGujatoQuFuZGEdaW2pI56fp1sR0GJ5WRJEmSlmDck5e0CWtL6dq60KQxva61FKOYgGdQLZmzPDGNLYSSJEnSlGoTsHod06b1bqVrNg4yrC3UajaqMZ1tWu4G2bo3yq6vthBKkiRJU2alax2OovVu1BPwDOL6/VpLh6nNGpfDYAuhJEmSNGNGvXzGoMJaL21aMhdrvVvKWMdhLfMB4xsPaSCUJEmSZsyoWu9W2u100BPT9Dq/XxAbxTIfML41Lg2EkiRJ0gyb5MXgB9Fq1i+s9Qtio1jmA8a3xqVjCCVJkqQZNA0zag6i1azfeMl+QWzQy3y0Gac4ypBuC6EkSZKkiTSIVrOlhLXlLuHR75i2XUpX2sV2OWwhlCRJkjSxVhqQ2gS6tq+xWCAd96yuy2UglCRJkjQTlhvCVhpKB7km4yDq6WYglCRJkrSqjXu85KAnjBnkeomOIZQkSZKkEVlukBv0EhxzbCGUJEmSpCFbaSvlsMYh2kIoSZIkSRNuWAvXGwglSZIkacINa+F6u4xKkiRJ0hQYxuQ4thBKkiRJ0owyEEqSJEnSjDIQSpIkSdKMMhBKkiRJ0owyEEqSJEnSjDIQSpIkSdKMMhBKkiRJ0owyEEqSJEnSjDIQSpIkSdKMMhBKkiRJ0owyEEqSJEnSjDIQSpIkSdKMMhBKkiRJ0owyEEqSJEnSjEpVjbuGZUtyHfDtcdehRe0NXDnuIrQg78/k8x5NPu/RZPP+TD7v0WTz/ky+e1fV7ss9ec0gKxmDb1fVunEXoYUl2ew9mlzen8nnPZp83qPJ5v2ZfN6jyeb9mXxJNq/kfLuMSpIkSdKMMhBKkiRJ0oya9kB44rgLUF/eo8nm/Zl83qPJ5z2abN6fyec9mmzen8m3ons01ZPKSJIkSZKWb9pbCCVJkiRJyzS1gTDJE5J8O8mFSTaMu55Zl+SgJJ9Ocn6SbyT5vWb/CUm+l+S85s+Txl3rLEuyJcnXmnuxudl3lyQfT3JB83Wvcdc5i5Lcu+tzcl6Sa5O8zM/QeCV5Z5Irkny9a1/Pz0w63tz8u/TVJEeMr/LZscA9+usk32ruwweT7NnsX5vkR12fp7eNr/LZsMD9WfDvtSSvbD5D307y+PFUPVsWuEfv6bo/W5Kc1+z3MzRii/yOPbB/i6ayy2iSnYH/Ah4LbAW+BDy7qr451sJmWJL9gf2r6twkuwPnAE8FngH8oKr+ZqwFCugEQmBdVV3Zte+vgO9X1cbmP1f2qqpXjKtG3fJ33PeABwPPw8/Q2CR5OPAD4NSqum+zr+dnpvml9iXAk+jcu7+rqgePq/ZZscA9ehzwqaq6MclfAjT3aC3wkbnjNHwL3J8T6PH3WpLDgHcDDwLuBnwC+NmqummkRc+YXvdo3vN/C1xTVa/zMzR6i/yO/VwG9G/RtLYQPgi4sKouqqqfAmcAR4+5pplWVZdV1bnN4+uA84EDxluVWjoaOKV5fAqdv2Q0XkcC36mqi8ddyKyrqs8B35+3e6HPzNF0fqGqqtoE7Nn8Q64h6nWPqupjVXVjs7kJOHDkhQlY8DO0kKOBM6rqJ1X138CFdH7n0xAtdo+ShM5/7r97pEXpFov8jj2wf4umNRAeAFzStb0Vw8fEaP736AHAF5pdL26arN9pd8SxK+BjSc5Jclyzb7+qugw6f+kA+46tOs15Fjv+4+tnaLIs9Jnx36bJ9FvAv3Zt3yPJl5N8NsnDxlWUev695mdo8jwMuLyqLuja52doTOb9jj2wf4umNRCmx77p6/u6CiW5E/AB4GVVdS3wVuCewOHAZcDfjrE8wUOr6gjgicDxTTcRTZAktwOeAryv2eVnaHr4b9OESfIq4Ebg9GbXZcDdq+oBwB8A70py53HVN8MW+nvNz9DkeTY7/geln6Ex6fE79oKH9ti36OdoWgPhVuCgru0DgUvHVIsaSXah84N6elX9M0BVXV5VN1XVzcDbsevHWFXVpc3XK4AP0rkfl891JWi+XjG+CkUnrJ9bVZeDn6EJtdBnxn+bJkiSY4EnA8dUM2FC0xXxqubxOcB3gJ8dX5WzaZG/1/wMTZAka4BfA94zt8/P0Hj0+h2bAf5bNK2B8EvAoUnu0fxv+rOAM8dc00xr+pifBJxfVW/o2t/dZ/lXga/PP1ejkWS3ZjAySXYDHkfnfpwJHNscdizwofFUqMYO/xvrZ2giLfSZORP4zWaGt/V0JmG4bBwFzrokTwBeATylqq7v2r9PM2kTSQ4BDgUuGk+Vs2uRv9fOBJ6V5PZJ7kHn/nxx1PXpFo8BvlVVW+d2+BkavYV+x2aA/xatGXDNI9HMGvZi4N+AnYF3VtU3xlzWrHso8BvA1+amJgb+GHh2ksPpNFVvAV44nvIE7Ad8sPP3CmuAd1XVR5N8CXhvkucD3wWePsYaZ1qSO9KZPbn7c/JXfobGJ8m7gUcCeyfZCrwG2Ejvz8zZdGZ1uxC4ns4MsRqyBe7RK4HbAx9v/s7bVFUvAh4OvC7JjcBNwIuqqu2EJ1qGBe7PI3v9vVZV30jyXuCbdLr6Hu8Mo8PX6x5V1Uncdjw7+Bkah4V+xx7Yv0VTueyEJEmSJGnlprXLqCRJkiRphQyEkiRJkjSjDISSJEmSNKMMhJIkSZI0owyEkiRJkjSjDISSJEmSNKMMhJKkvpLclOS8JF9P8r4kd0yyNsnX+5+9w3Wem+RuA6ppXZI3D+JafV7n7CR7LvPcpyY5bNA1LVeSRya5JsnZXfsOTfKRJN9Jck6STyd5+CLX2C3JVUn2mLf/X5I8I8kzk1yY5CPDfC+SpMEwEEqS2vhRVR1eVfcFfgq8aJnXeS4wkEBYVZur6qWDuFYv6dipqp5UVVcv8zJPBSYmEDY+X1VPAkiyK3AWcGJV3bOqHgi8BDhkoZOr6ofAx+i8N5rr7AH8MvCRqnoP8IIh1i9JGiADoSRpqT4P3Kt5vHOStyf5RpKPJbkDQJLDk2xK8tUkH0yyV5KnAeuA05vWxjskOTLJl5N8Lck7k9y+OX9LktcmObd57ufmF9G0dn2keXxCc/5nklyU5KXN/r9M8rtd55yQ5OVJ7pTkk13XP7p5fm2S85P8PXAucFBTy97N8//StKJ9I8lxXdf9QZLXJ/lK8773S/IQ4CnAXzfv957z6j85yVubFrmLkjyieQ/nJzm567i3JtncvOZru/ZvTPLN5nv8N82+pzetuF9J8rkW9/IY4D+r6sy5HVX19ao6ubnebk1NX2ru09HNYe8GntV1nV8FPlpV17d4TUnSBDEQSpJaS7IGeCLwtWbXocD/qar7AFcD/6vZfyrwiqq6X3Psa6rq/cBm4JiqOhwo4GTgmVX1C8Aa4He6Xu7KqjoCeCvwhy3K+zng8cCDgNck2QU4A3hm1zHPAN4H/Bj41eb6jwL+NkmaY+4NnFpVD6iqi+e9xm81rWjrgJcmuWuzfzdgU1XdH/gc8NtV9X+BM4H/3bSufqdHzXsBjwZ+H/gw8EbgPsAvJDm8OeZVVbUOuB/wiCT3S3IXOiHsPs33+M+bY/8UeHxTx1NafM/uQyf4LuRVwKeq6hfpfJ/+OsluwEeBB3a9/2fRCYmSpCljIJQktXGHJOfRCXTfBU5q9v93VZ3XPD4HWNt0H9yzqj7b7D8F6DUm7d7N+f+1wHH/3H3dFjWeVVU/qaorgSuA/arqy8C+Se6W5P7A9qr6LhDgL5J8FfgEcACwX3Odi6tq0wKv8dIkXwE2AQfRCcTQ6UY7N2aubb0AH66qohOaL6+qr1XVzcA3uq7xjCTnAl+mE+AOA66lE2rfkeTXgLmWuf8ATk7y28DOLWu4RdOa+/Ukc9/7xwEbmnv/GWBX4O5V9VM6YfdpTevp4XS6kUqSpsyacRcgSZoKP2pa9W7RNKj9pGvXTcAdlnDN9Hl+7to30e7fq/m1zJ3zfuBpwM/QaTGETlfJfYAHVtUNSbbQCTsAP+xZbPJI4DHAL1XV9Uk+03XODU2wW0q93TXfPK/+m4E1Se5Bp3X0F6tqe9OVdNequjHJg4Aj6bTOvRh4dFW9KMmDgaOA85IcXlVXLfL636ArhFfVryZZB/zN3NsG/ldVfbvHue8GXt0c86GquqHle5YkTRBbCCVJA1VV1wDbkzys2fUbwFxr4XXA7s3jb9FpUbxXj+MG6Qw6oelpdMIhwB7AFU0YfBRwcIvr7EGnhfH6Zkzj+hbndL/f5bgznYB6TZL96HTXJcmdgD2q6mzgZXRa6Ehyz6r6QlX9KXAlnVbMxbwLeGiS7u6ld+x6/G/AS+a60yZ5QNdzn6bTQno8dheVpKllC6EkaRiOBd6W5I7ARcDzmv0nN/t/BPxSs/99zdjELwFvG3QhVfWNJLsD36uqy5rdpwMfTrIZOI9OOO3no8CLmm6m36bTbbSfM4C3N5PcPG2BcYSL1f6VJF+m05J3EZ0uodAJmR9KZ5bQ0BmDCJ0xfoc2+z4JfKXP9X+U5MnAG5K8CbicToidG5P4Z8CbgK82oXAL8OTm3JuTfAB4Op1xk5KkKZRbe7hIkqTVrOn2+odV9eTV8DqSpJWzy6gkSbPjp8B907Uw/aAleSbw98D2Yb2GJGlwbCGUJEmSpBllC6EkSZIkzSgDoSRJkiTNKAOhJEmSJM0oA6EkSZIkzSgDoSRJkiTNqP8HEFFQWr7iYhkAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 1080x432 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "fig, ax = plt.subplots(figsize=(15, 6))\n",
    "# For a normal histogram (without error bars) one would use:\n",
    "# hist_data = ax.hist(x_data, bins=Nbins, range=(xmin, xmax), histtype='step', linewidth=2, label='Data', color='blue')\n",
    "y, bin_edges = np.histogram(x_data, bins=Nbins, range=(xmin, xmax), normed=False)\n",
    "x = 0.5*(bin_edges[1:] + bin_edges[:-1])\n",
    "sy = np.sqrt(y)      # This is the standard for histograms - each bin is Poisson distributed!\n",
    "hist_data = ax.errorbar(x, y, sy, fmt='.', linewidth=2, label=\"Data\")\n",
    "ax.set(xlabel=\"Photon invariant mass [GeV]\", ylabel = \"Frequency / {0:1.0f} GeV\".format(binwidth),\n",
    "       title = \"Distribution of diphoton invariant masses\", xlim=(xmin, xmax))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "---\n",
    "Fit the data:\n",
    "---"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "NOTE, that for the starting values for the fit, we here (to begin with) use the values that the data was generated with!\n",
    "That is of course a great advantage, and one of the main challenges in general fitting is to find these good starting values."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## $\\chi^{2}$ fit:\n",
    "Do the fit with a Chi2 minimisation (only using bins with entries):"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {},
   "outputs": [],
   "source": [
    "chi2reg = Chi2Regression(PDFmodel_CS, x[y>0], y[y>0], sy[y>0])\n",
    "if IncludeSignalInFit:\n",
    "    minuit_cs = Minuit(chi2reg, pedantic=False, print_level=0, Nsig=Nsigx, mu=mux_init, sigma=sigmax, Nbkg=Nbkgx, tau=taux)\n",
    "else:\n",
    "    minuit_cs = Minuit(chi2reg, pedantic=False, print_level=0, Nbkg=Nbkgx, tau=taux)\n",
    "minuit_cs.migrad()                                                         # Perform the actual fit\n",
    "if (not minuit_cs.get_fmin().is_valid) :                                   # Check if the fit converged\n",
    "    print(\"  WARNING: The ChiSquare fit DID NOT converge!!!\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "  ChiSquare Fit result: Nbkg = 24887.23188 +/- 174.68929\n",
      "  ChiSquare Fit result: tau = 79.90851 +/- 0.84052\n"
     ]
    }
   ],
   "source": [
    "minuit_output = [minuit_cs.get_fmin(), minuit_cs.get_param_states()]       # Save the output parameters in case needed\n",
    "if IncludeSignalInFit:\n",
    "    csfit_Nsig, csfit_mu, csfit_sigma, csfit_Nbkg, csfit_tau = minuit_cs.args  # The fitted values of the parameters\n",
    "else:\n",
    "    csfit_Nbkg, csfit_tau = minuit_cs.args\n",
    "Chi2_value_cs = minuit_cs.fval                                             # The Chi2 value\n",
    "NvarModel_cs = len(minuit_cs.args)\n",
    "Ndof_cs = len(y[y>0]) - NvarModel_cs\n",
    "ProbChi2_cs = stats.chi2.sf(Chi2_value_cs, Ndof_cs)\n",
    "for name in minuit_cs.parameters:\n",
    "    print(\"  ChiSquare Fit result: {0} = {1:.5f} +/- {2:.5f}\".format(name, minuit_cs.values[name], minuit_cs.errors[name]))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Plotting the result on top of the previously plotted data:\n",
    "---"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAA4QAAAGDCAYAAABz1N6qAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvnQurowAAIABJREFUeJzs3Xd4VNXWx/HvIglNFBQRBYSIBUGpAUTECoiIvV0VFSvqtaNX8OpVvHoVrwVfbIhiB7tXUQR7RUFAsYIKGBAQpaqA9P3+sU7IJMwkk5BJCPP7PM88mbPnnDN7zgw6a9bee1kIAREREREREUk/VSq6AyIiIiIiIlIxFBCKiIiIiIikKQWEIiIiIiIiaUoBoYiIiIiISJpSQCgiIiIiIpKmFBCKiIiIiIikKQWEIiIVwMyGmtm/yuhcjc1smZllRNvvm9m5ZXHu6HxjzKxPWZ2vBM97s5ktNLP5Se4fzGy36H7S19fMHjOzmzelr6lgZvub2fcV9Nxl9vkUEZHNm6kOoYhI2TKzXKA+sBZYB3wHPAEMCyGsL8W5zg0hvF2CY94HngohPFyS54qOHQjsFkI4raTHliUz2xn4AWgSQvgtyWMCsHsIYXoJn+sxYE4I4bpS9PN9Snmt04GZZQM/AVkhhLUV2xsREYlHGUIRkdQ4MoSwNdAEGAT0B4aX9ZOYWWZZn3Mz0QRYlGwwKGUnL9MsIiLpQQGhiEgKhRB+DyGMAv4G9DGzvaHgMEUz297MXjOzpWa22Mw+MrMqZvYk0Bh4NRoSerWZZUdDI88xs9nAuzFtscHhrmb2mZn9bmavmNl20XMdZGZzYvtoZrlm1s3MDgP+Cfwter4vo8c3DEGN+nWdmc0ys9/M7Akzqx09ltePPmY2OxrueW2ia2NmtaPjF0Tnuy46fzfgLaBB1I/HEhz/DzP7xczmmdnZhR6Lvb4HmdkcM/tn1KdcM+td6HTbmtloM/vTzCaY2a4x5+psZhOjaznRzDpH7f8B9gfujfp5b1H7x1zLm8xsXPRcb5rZ9gleX4H3Kur3VWb2VXTuZ82sevTYVDM7ImbfzOi1tou2nzez+dFxH5rZXoWu1QNm9rqZLQcOLnT9to0+nwvMbEl0v1GSr+nD6O/S6BrtG+d1Doz691R0/NdmtoeZXRN9xn42s0Nj9j8rer1/mtlMMzs/5rG4/5aix/qb2dzouO/NrGvUXsXMBpjZDDNbZGbPWf6/l+pRvxZF55xoZvXjvV8iIpWVAkIRkXIQQvgMmIMHEIVdGT1WDx9q+k8/JJwOzMazjbVCCP+NOeZAoDnQI8FTngGcDTTAh64OSaKPY4FbgGej52sdZ7czo9vBQFOgFnBvoX26AM2ArsD1ZtY8wVPeA9SOznNg1OezouGxPYF5UT/OLHygefB6FdAd2B3oVszL2xHYHmgI9AGGmVmzmMdPAW4EtgWmA/+Jnmc7YDR+/eoCdwGjzaxuCOFa4CPg4qifFxe1f8xznQqcBewAVI1eR7JOAg4DdgFa4e8FwNPRa8jTA1gYQvg82h6DX6cdgM+BEYXOe2r0mrcGPi70WBXgUTxr2xj4i43f80Sv6YDob53oGn2a4HUdCTyJX/8vgDei520I/Bt4MGbf34AjgG2i5xycF/iS4N9S9F5fDHSIMvc9gNzomEuBY/DPYANgCXBf9Fgf/DO6M/5+XhC9fhGRLYYCQhGR8jMP2C5O+xpgJ3y+3JoQwkeh+AneA0MIy0MIib6cPhlC+CaEsBz4F3CSlc1QwN7AXSGEmSGEZcA1wMlWMDt5YwjhrxDCl8CXwEaBZdSXvwHXhBD+DCHkAncCpyfZj5OAR2Ne48AkjvlXCGFVCOEDPGg7Keaxl0IIn0Xz3EYAbaL2XsCPIYQnQwhrQwhPA9PwACaeZPZ/NITwQ/TePRfzXMkYEkKYF0JYDLwac+xI4Cgzqxltnxq1ARBCeCS6zqvwa9Xaosxu5JUQwrgQwvoQwsrYJwwhLAohvBhCWBFC+BMPHA8s1K9NeU0AH4UQ3oiu//N4QDcohLAGeAbINrM6UX9GhxBmBPcB8Cb5P7Qk+re0DqgGtDCzrBBCbghhRnTM+cC1IYQ5MdfnhOgzvQYPBHcLIawLIUwOIfxRwtcmIrJZU0AoIlJ+GgKL47Tfjmel3oyGwA1I4lw/l+DxWUAWniHbVA2i88WeOxPPxuSJXRV0BZ5FLGx7PJNU+FwNS9CPwq+xKEuiwDF2/wYx24n6XPj1FtfPZPZP5vokEvfYaCGdqcCRUVB4FFFAaGYZZjYoGhL5B/mZsdjPQ8LPk5nVNLMHzYf1/oEPA61T6AeGTXlNAL/G3P8Lz26ui9km75xm1tPMxkdDQpcCh8e8lrj/lqLrczke7P1mZs+YWd773wT4XzQkdCl+Hdfhn+kn8WzlM+ZDk/9rZlklfG0iIps1BYQiIuXAzDrgQUHh4XhEmZsrQwhN8UxSv7z5TUCiTGFxGcSdY+43xjMdC4HlQF4WKS9TV68E552Hf4GOPfdaCn6hT8bCqE+FzzU3yeN/YePXWJRtzWyrQvvPS+J5Cr/evGPz+ln4ehW3fyrlDRs9GvguZrXVU6O2bvjwx+yo3WKOLep9vxIfArxPCGEb8oeBWuJDkjpviZlZNeBF4A6gfgihDvB6Xl+K+rcUQhgZQuiCvz8BuC067c9AzxBCnZhb9RDC3CjLeGMIoQXQGR+qekZZviYRkYqmgFBEJIXMbJtosY9n8PIEX8fZ5wgz283MDPgDz07kZUd+xefYldRpZtYiyhb9G3ghyrj8AFQ3s15RpuM6fChdnl/x4XmJ/v/wNHCFme1iZrXIn3NYopICUV+eA/5jZlubWROgH/BUkqd4Djgz5jXekMQxN5pZVTPbH/9i/3wSx7wO7GFmp5ov1PI3oAXwWvR44fenuP1T6RngUOBCYoaL4vMCVwGL8B8DbinhebfGs3RLozmSyVzrPAuA9ZTuMxxPVfzzugBYa2Y98dcMJP63ZGbNzOyQKKBcib+evH9jQ/HPYZPoHPXM7Ojo/sFm1jL64eQP/EeMvONERLYICghFRFLjVTP7E88+XIsvLnJWgn13B94GlgGfAveHEN6PHrsVuC4azlaSxUeeBB7Dh/JVxxfOIITwO/B34GE8a7UcX4QjT16QtMjMPmdjj0Tn/hCvL7cSuKQE/Yp1SfT8M/HM6cjo/MUKIYwB7gbexYcIvlvMIfPxxULm4XMELwghTEvieRbhweOVeEB1NXBECGFhtMv/4fPNlpjZkCT2T5kQwi/456cz8GzMQ0/gw1bn4jUxx5fw1HcDNfCs7nhgbAn6tAKfczgu+gx3KuFzFz7fn/hn+Tn8/TwVGBWzS6J/S9Xw8i8L8c/CDviCM+Dv4Sh8mOmf+GvcJ3psR+AFPBicCnxA8j9aiIhUCipMLyIiWzQzOwjPzjYqbl8REZF0owyhiIiIiIhImlJAKCIiIiIikqY0ZFRERERERCRNKUMoIiIiIiKSphQQioiIiIiIpKnMiu7Apth+++1DdnZ2RXdDRERERESkQkyePHlhCKFeaY+v1AFhdnY2kyZNquhuiIiIiIiIVAgzm7Upx2vIqIiIiIiISJpSQCgiIiIiIpKmFBCKiIiIiIikKQWEIiIiIiIiaUoBoYiIiIiISJpSQCgiIiIiIpKmFBCKiIiIiIikKQWEIiIiIiIiaUoBoYiIiIiISJpSQCgiIiIiIpKmFBCKiIiIiIikKQWEIiIiIiIiaWqLCQizB4wme8DoCj+HiIiIiIhIZbHFBIQiIiIiIiJSMgoIRURERERE0pQCQpFUyciANm3ybx99VPwxd98NK1Yk/xyjRsGgQaXvo2ye5s6FLl2gRQto1w7eeiv/sdjP1aWX5rc/9xzsvjvssQe8+mrp20VERCStWAihovtQau3btw+TJk0C2DD3L3dQr1KfryzOIbJBrVqwbFnJjsnOhkmTYPvtU9IlqSR+/RXmz4fWrWHWLOjc2YNEiP+5Wr0amjWD8eNh1So4+GD48UdYu7Zk7VX0G6GIiEhlY2aTQwjtS3u8/u8vUt5q1YJrroG99oLzzvO2t9/2jM+8ef7lPO8+wPvvQ/fucMIJsPfecPnl3n7GGdC4MVx8ccHzv/oqdOzo5+jXL7/9//4PmjeHVq2gf/+Uv0zZBPXrezAI0KQJrFnjgVsiEyb456l+ff9M7LwzfPVVydtFREQk7WRWdAdEtlh//eVBWZ5XX/Uv3suXw5FHwr//Dbvu6oFft24wZYpnCN97b+MM4bhx8NlnHhAuXeptTzwBjz3mGcU8v/0GN93kQWTNmnDiifDuu3DIIXDDDfDzz7D11r6fVA5vvAFt20K1ar69cqUPI61RA269FQ44wLOJO+0EQ4fCdtvBjjvCL794JrEk7bGfVxEREUkLCghFUqVGDQ/yCqta1YcAAuyyi3+Zb9Cg6HPl5HgwCFCnTuL9xo+HmTPzz79smW8fcgh06ABnnw29esGxx5b89Uj5mz8frrwSXnklv23OHA/gJk709/HHHyFv6P8FF/jfl14Cs5K3i4iISNpRQChS3rKy8u+bwfr1xR9Tu3by5+/RA0aM2Lh97Fj4+GN45hm45x6YPDn5c0r5W7nShwnfcYdnkvPsuKP/7dDBs3yzZvkPCr/8kr/P/Pm+X61aJWsXERGRtFO5A8JKvCCOSFzbbAOLF5d+UZlOnXxO4Zw50KiRBwvVqvlcsdmz4cADoWVLX11SNl8hQJ8+cOqpcNhh+e2LF3vmuUYNyM314caNG0PTpvDttz4UeNUqf/9btfLFY0rSLiIiImmncgeEP/wAS5bAttuWy9NpFVIpkcJzCG+4ofihmpdcAscc4/O6XnghcdYmN9f3W7zYn+fjj+GWW+Dww+GBB3yO4tq1sNVWni0MAU4/Hf74A9atgzvvLLOXKSkwbpwP4/z+exg2zNtef93f97PO8iA/IwMeftjnioKXH9lvP78/eLCvGFq1asnaRUREJO1U7rITZmFSixYwZgzZ938NwIsXdianSekCxLyAL9E5FBCKiIiIiMjmJL3LTlSvDt99x+SjTsOCz8Pq/fB4Js9aUuJTxR5T2nOIiIiIiIhUJpU7INxzTzjgAMbX2BGLMp1r1q5n/MxFJT5V7DGlPYeIiIiIiEhlUrkDwowMeOMNOu26PVXXrSVj3VqyCHRqWrfEp4o9JiuzSqnOIbLBY4/5PK2FC337iCO8NiD4Co/xZGfn71/R5s6FLl2gRQuveffWWwUf//NPX9nyjjvy226/3UtjtGgBN95YdPsbb/j8yrxb1arxS3QU159Fi6B9ey/i3qYNvPxycv0vjfXrvZh7ad+jq67yxX3yyodA8dch3nVOZMIEXximeXM46aT89kTvy7XXenu7dvnXTURERNJO5V5UBqB6dXIev4dTD72Auit+p9Psr8mpfz5cfnmJThM7Z3DEuZ1KPQ9RZIOaNeHpp32hmMomMxPuu88DrVmzvK7h3Ln5j998s9dGzDNnDjz4IEyd6gvY7LknnHGGl9iI196jh9/Ayx8ccEDRRdET9WebbTzQrlXLA7WWLeGoo4rvfyKPPeYLtwwcuPFjr7/upR5KuwLsccfBySfDmWfmtxV3HQpf50TyFg165BEPhBcs8PZE78uiRfDmmx58Ll3qhe+7doWtty7daxMREZFKq3JnCPNkZPBI+6P5vXotcuZNgyuu8F/jk6nvFoeCQSkTRx4Jzz+/cXsIcN55np257LKNH58xw7M2P/7o27ff7l/kjzrK69Hl5qa024Bnslq39vtNmsCaNV6eAHx13wULNg5U1q71fVat8kxXXu3ERO15nnnG6+2Vpj9ZWfkZ199/97a1a4vuf2ndc8+mBfedOxcdTBa+DomuczyTJkG9eh4Mgt/PE+/6z5jhgWdmpvepYUMvdC8iIiJpZ8sICAHMGLbP8fDkk/4l58474bTTNv1LoEhp1a7tw/2mTi3YvmKFf/H/6isvDv/RR/mPzZnjw/0ef9xrBc6e7RmeyZPhv/+Fn34q39cAPqyxbVsvdQDQv7+X0IjVqJEHt40b+/0rr/TSGYnaY40Y4Zmz0vbnzz89M9iqFdx/vwc9Re1fGtOmwbJlyQVnpVX4OsS7zonMng077AA9e/prvf9+b090/Vu08CGmK1b4sVOnenF6ERERSTuVf8hoYaed5tmB447z4Xq//ur1vApnJUTKwxlneHAXKyMDunf3um/du3t2Z//9/bGjj4ZevTzAAQ8E99/f6wnuuadnu8rT/PkeRLzyim+/+irsscfG/Vi6FMaM8YB1zRqvb3fEEV5APV77Tjv5cT/8AMuX52fzStof8GGOX3/tQU3fvnD88Z45TLR/PIMH+/u0eDGsXu1z6rp39+wswL33wsUXJ9fH0ih8HRJd50RWroT33vMfGerU8XmVhx3mwV+869+ypQ9d7dzZs4MHH+yrNouIiEjaSWlAaGZ1gIeBvYEAnA18DzwLZAO5wEkhhCVmZsD/AYcDK4AzQwifl+qJu3eHDz/0It3vvuvzcsaM8WyNSHnq0cOzPLFD+Aozy78/bBj84x/w2WfQsaMPL60oK1d6JvOOO3yoKnhW6cUXPcBauNCD2p128uxb48b5P7y0bQtffOEZqHjteQFhSbKD8foTq3lzDwS//NIDouL2j3XFFX6LN4fwjz98vt3gwQWPmTABzj/f799/vwdXpVX4OiS6zr17xz++fn1//Y0b+3ZOjmc1i7r+/fr5DaBTp/xjRUREJK2kesjo/wFjQwh7Aq2BqcAA4J0Qwu7AO9E2QE9g9+jWF3hgk565bVv49FNo1sx/Ne/UybMIZUA1CiVpGRme4fvkk/y2devgnXd8juvbb/tCJXlycnyI6Nln+3Dn9u3h44/9i/20ab5ASnkIAfr0gVNP9UxTnptvhunTvS8XXwxXX+1Byo47+hy0Vavgr7/g889hl10St+cZOTJ+QHjNNX4rrj9z5/oCKeDZwO++84xXov1L49FH/TXmZR3z7LOPL8oyZcqmBYOw8XVIdJ3zFL4+HTr40M+8DOfXX0PTpkVf/7zr9sEHflwqh8OKiIjIZitlGUIz2wY4ADgTIISwGlhtZkcDB0W7PQ68D/QHjgaeCCEEYLyZ1TGznUIIv5S6E9nZMG6cL8bxySe+4MILL3gGsYQKF67XSqSStD59CmaXatb0BUSuuMJXdtxvv4L777uvf0b/9S+fN9i3ry8y06aNZ7o2ZS5cssaN86HW33/vWUvwVTYTZdm7dPFsaOvWns0691zPWEHi9gkTfEGYZs02Pt8vhf7ZJ+rP7Nl+fcAD7dtu8+zXxx+XrP95YlcABQ8sH3rIA/dNddFF3qeFC30+3/33+3+biroOiRS+PrVrw913wyGH+NDQU0/1IcaQ+PqfdZYHnVlZ8NRTBTPVIiIikjYspGhImpm1AYYB3+HZwcnAZcDcEEKdmP2WhBC2NbPXgEEhhI+j9neA/iGESYXO2xfPINK4ceOcWVHGJHvAaAByB/XauDN//eVfyp9/3hecGToUzjlno92KOsd9703n9je+ByDDoN+hzbjo4N1KcEVESunPP32e3MKFnsXJzdWX9/IyZowvVDVyZEX3RERERCQuM5scQmhf2uNTOWQ0E2gHPBBCaAssJ394aDzxvuFuFK2GEIaFENqHENrXK2peVqwaNTwj07+/L8F+7rnwz3+WqCyFCtdLhbniCl9B85BDYMgQBYPlqWdPBYMiIiKyRUvlojJzgDkhhAnR9gt4QPhr3lBQM9sJ+C1m/51jjm8EzCuz3lSpQjb7c3KP5Qx6+wG49VaYOdMXkUhidT0VrpcK8/DDFd0DEREREdlCpSxDGEKYD/xsZnkTY7riw0dHAX2itj5A3nrwo4AzzHUCfi/J/MHcQb3iDxct5Jk2h8Ho0T4E79lnfQ7XwoUlOoeCQRERERER2RKkepXRS4ARZvYV0Aa4BRgEdDezH4Hu0TbA68BMYDrwEPD3lPWqRw9fdKJRI19splMnrwMmUtbuvttXCJXkvPGGL56Td6ta1VfxBLj2Wth7b19g5+WXiz9Xov2few52393r/L36avH7l9b69bDXXht+cCqxRP1M5Pbbvf8tWsCNNxZ87M8/fUGdO+7Ib7vqKi9XsffepeufiIiIbBFSWocwhDAFiDfBsWucfQNwUSr7U0CrVr663xFHeF2ufff1L4F5BcJFysLdd8Npp/nKolK8Hj38Br6S5gEHeGA4aZLXApwyBZYu9bIyXbt6pj+eRPtXq+ZziceP91IMBx8MvXp5OYaSnD9PvLqFeV5/3ctBbL99ya/D6tXx+1klwW94c+Z4uZKpU31l1D33hDPOyC8xcfPNG5eVOO44L3VReGVVERERSSupzhBu3ho08AL2Rxzhdbi6ddMCElI23n7bA5l58/zLfN598FIDrVv7F/R7780/plat/PsHHeRBTTp75hkvLA8wY4Zfw8xMD7AaNvT6eokk2n/CBM/a1a/vhdh33tnrlJb0/Mm45x645JLSHZuon0VZu9aDx1WrPLOaV4z+hx9gwYKNA8LOnUsXrIqIiMgWJb0DQvAv4S+/7MWfV6/24s833eS/souUVrdunm1q0ADeey//Pnj9uS+/9OzPkCFeUF02NmJEfrH2Fi08SFqxwmsPTp1a9HVLtP/8+V6ncOhQH5K5446eiSzp+YszbRosW1b6Yu+J+plIo0Zw2WUePDZqBFdeCdtt54/17w833FC6foiIiMgWL6VDRiuNjAz/Yr7rrtCvH1x/vX8hHD7cS1aIlKXhw/1HiBA8azhvnn/hl3w//ADLl3smFaBlSx/a2LmzZ+8OPrjo1YET7b92rT9+wQX+96WXvIxHSc8/eDA8/riPLFi92t/P7t19Hh945vfii0v/+vN+kCrcz0SWLvWaiT/95IXp99vPRz5MmuRzEJs0KX1fREREZIumDGEeM7j8cv9it9VW8PTT/qVQ2RspS++/7wunjBvnWcM998yvhxn7hT8vcElXsdnBPP36+TUbPdr/XTZuXPQ54u3foEHBTNv8+fnBeEnOf8UVvu+//+1B25Qp+cHgH3/4fMS84a55JkzIXyznk0+K7ntR/Yzn7be9v7Vr+zDQtm19bvSECfDii/45u/de+O9//dqKiIiIRJQhLOyoo/zL2pFH+pepjh1h1Cj/EidSUtts41mkvLlaf/zh92vWhG+/9aGjeWrX9n1r1PAhh+ls5Eh47bWCbYsWQd268MEHfp1ih2Nec43/vfXWovdfs8av+2+/+Vy7OXN8ganizl8Sjz7qQ8+zsgq277NP/oqpxenYMXE/473eHXf0OY+rVvkPDJ9/7sNEDz/cF5QBX/imVi3vm4iIiEgkLQPCybOWFF1LsFUr+OwzOPZY+PRTH341YgSQlfgYkXguuQSOOcbnc73wAhx2GAwb5kMh99zTSxzk6d/fH2/f3hcRSVcTJnjg0qxZwfazzoLp0z3QeuqpghnVePPr4u1ftSoMGuT/psGHfuat3FnU+RMpvEJnCPDQQ56x2xRF9RM2fr1duvjqrK1b+37nngvNmxf9HBdd5ENRFy70eYf33+8/iImIiEhasVCJF09p3759mJTkSoyTZy3h+Ad8mFb1rCqMOLdT8QXmV66Evn3hySd9+9Zb/Ut7Ml8URST9jBnj/73QasUiIiJSTsxscgghXqm/pKTNHMLxMxdtuL9m7foC2wlVr+4LR+QNy7rmGujTxwNFEZHCevZUMCgiIiKVStoEhJ2a1t1wPyuzSoHtIpnBgAE+tKpmTf/1/5BD4NdfU9RTERERERGR8pE2AWHs8NCkhosWduyxvjJko0Y+r7Bjx+ILRUt6++EH2H9/r3HXtq3Xt4OCBegLGzXK547lmTvX54e1aOHzDd96a+Nj1q/3IuYLFybft733htzc/O2PPvLnaNMGvvsu+fNsSSZM8PnDzZvDSScV3x7PokU+B7R1a7+WL79c8PE///QVRO+4w7eTeX8TSfZ9v/Zaf7/btSvYn9tv9/YWLeDGG5N/XhEREdmipOWiMiUOBvO0aeMr+R1zjH9J7NzZh4dpIQaJ5/TT4c47/Qv/vHlQrVrxxxx1VMHPU2Ym3HefBxizZvlnbu7cgse8/jp06JC/kmlpjBgB//iHL6ySjkLw9+uRR/z9WrCg6PZEttnGS4vUquWBWsuW/n7mLQhz880FVy9N5v1NJJn3fdIkL4ExZYrXKmzbFrp2hd9/hwcf9HqrIfgCR2ecAbvsktxzi4iIyBYjbTKEZWbHHf0L36mneuHso4+Gm27KryUnAl5Oolo1DyLAs0J1Y4YpX3ONZ3fOOy+/7YwzvJZcbEHz+vXzi7M3aeJlE1atKvhc99zjq5nmqVUr/vn/+1/Pcp14Yv482KVL/YeO557zmnrpmiGcNAnq1ct/v+rVK7o9kays/Azw77/7e5VXU/KHHzygjA0Ik3l/Eyn8vsczY4a/p5mZHjg2bOg/aoH3a9Uqv1Wt6mVPREREJO0oICyN6tV9WfpbbvE5htdf71+yly2r6J7J5mLmTNh11/iPLV/udS6nTPEi9fPmefsTT3hQlsgbb3iGJzbTOG2af+5ig4x458/N9XIIn3/u9ehmzvR969Tx/Y46yocQTpniQwjTzezZsMMOvihM27ZegqGo9qL8+adnBlu18v2rVvX2/v29NmAi8d7fROK97/G0aOGjGVas8NcydaoXuW/UCC67zH+AaNQIrrzSS6OIiIhI2knLIaNlwsyzMK1aebbwpZc8A/DKK9C0aUX3TipaCInLk1St6kMDwYfozZ/vGcSizJ/vX9pfeaVg+733FswoJjp/bq5nuWrU8MxhkyYlfklbtJUr4b33fF5wnTo+D/CwwxK3F/VvfOut4euvPfjq2xeOPx7GjoU99kh83RO9v4nEe9/jadnSayV27uzZwYMP9h+0li71Ehk//eRZyf32gyOOgJ12Su75RUREZIuhDOGm6tXLi9jvuSd88w1L92pN75P/U9G9korWtKkXOY8nKyv/vlnxw41XroQTTvCFSGKzjn/84fPDTjih+PObFQxQK3H90ZSoX9+H0zZu7PMAc3I8C5eoPRnNm/t78eWXnqV78UX/78S99/rw3REjfL/5Au6eAAAgAElEQVRE728iid73RPr188zv6NEeeDZuDG+/7X9r1/ahpG3bwhdfJHc+ERER2aIoICwLzZrB+PFwxBHUWbmMJ567HgYP1pfudNamjQ/d/Phj3/7lF1+BsqRC8NqXp57qmalYjz4KvXsXDAATadfOV8ldudLnCOateCquQwe/JosXw+rVnuFr2jRxe55rrvFbnrlz89/n+fP9Wjds6IvJTJ/uweTFF8PVV/t7V9T7m0hJ3nfI788HH/jryMnxudATJ/r8wb/+8qHEWlBGREQkLWnIaFmpXRteeYV7upzKJZ8+67/Kf/GFr+RXo0ZF904qwpNP+qIuixf7l/fXXiu4sEys3FxfvXbxYv+C/vHHPkd1m218OPL338OwYb7v66/70L6HHvJMTzKaNIFzzvFMUIsW+vJfWO3acPfdXmN0zRoP0Pbc0x9L1A4e6MeaPduHiQKsWwe33Vb0MMxx4+K/v4mGEIdQsvcdfOXY6dP9M/jUU54p7tIFevTwBW2qVIFzz/WMpoiIiKQdC5U4i9W+ffswadKkpPfPHjAagNxBvVLVJbIHjKbntI954K0hvpBD+/bwv//5wg0iZWXMGA84R46s6J5IedL7LiIiIoWY2eQQQvvSHq8hoykwZs8uXrx+l1182fr27T0TIFJWevZUUJCO9L6LiIhIGUurgDB3UK+UZgcLaNXK5+gccgj8+quv7jd0qOYVioiIiIjIZiOtAsJyV7eu1xa74gqff3ThhXD22T5HbBNkDxi9YfiriIiIiIhIaSkgTLXMTLjrLp/3U6MGPPaY1/z66aeK7pmkWkaGrza6115w8sml/yGgVq2y7ZdUDs89B7vv7vULX321+P2vvRb23ttXlH355eLbAf780xewueOOsu27iIiIVBoKCMvLaad5aYpdd/XVR3NyvFi1bLlq1PD6b99849tDh1Zsf6TyWL0a+vf31Wbffhsuv7zoepWTJnltwilT/O8ll3iwl6g9z803+3+LREREJG0pIEyRybOWbNzYqhWtjvgPb+3WEZYsgcMPh5tuKr4wuVRuZnDggfmF6t9/H7p398Lie+/tX/YBhgzx7b339lpzeULw8hV77w2XXVbu3ZcKMGGCZ5br1/cC8jvvDF99lXj/GTM8G52Z6YXmGzb0OcyJ2gF++AEWLFBAKCIikuYUEJah2CCw98Pj4waFf1SvRd/jrvNf5gGuvx6OPhqWLi2vbkp5W7PGa8u1bJnfNm4cDBzo2cOBA70O4ZAhnkX+6CO44Qb47Tffd8UKDx6/+gomT/bHZcs2f77XLxw61IeO7rjjxjUPY7Vo4UHkihVeC3HqVD9HonbwDOQNN5TP6xEREZHNlgLCMjR+5qIN99esXV9gO1awKj6vZ8wY2G47L1jevn3RGQCpfP76y7MzHTr4UOFzzsl/LCfHM34Ader4MOL99/f5gttuCx075n8eMjI8o1iliv8tQe1NqaTyViO+4AI46SS/b5Z4/5Yt4cwzoXNnX7zq4IOhevXE7a++6nMTmzRJ9SsRERGRzVxmRXdgS9Kpad0N97MyqxTYjqtHD8/4HH88fP45dOoEw4b5fEOp/PLmEMZTu3bB7ZKUIykqMJAtQ4MGBTOC8+d7lrAo/fr5Dfy/JY0bJ25/+WV48UV45RVYuNB/bNhpJ+jdu+xfi4iIiGzWlCEsQzlNtt1wf8S5nQpsJ5Sd7QtHnHWWZ5ROP90Xfli9OnUdlc1Pu3b+OVi+3IcPT5zotSwB1q2Dd97xuaZvv+0ZR9mydewI337rw4Z//hnmzMn/PABcc43fYi2KRiR88AEsXpw/NzBe+803+5zWadPg4ovh6qsVDIqIiKQpZQhTJKlgME+NGjB8OOy7r385u/deDwieey7/V/4SyKtRmDuoV4mPlQqSne0/BOyzj2/fcAPssIPfr1kTnnnG61l27eplS2TLVrUqDBqU/14PHuxZvDzx5hOedZYHeVlZ8NRT+ZnkRO0iIiIiKCDcfJj5SpKtW8OJJ/pCEG3bwhNPQK/4gd3kWUtKFnhK+Vq2LH77QQf5rbBLL/VbYcuXl2WvpLL429/8Fs9jj23cNmpU/H0TtecZOLAkvRIREZEtjIaMbm46dvT5hL16+fCuI47woWFr1wLJrWQqIiIiIiKSDAWEm6O6df1X/dtu8xUmBw2CQw6BefOSXslURERERESkOAoIK0ixmb0qVXyhh3ff9dX/PvoI2rSh06IZG3ZJaiVTERERERGRBBQQlqNSDfc84AAvXdCtGyxYQM5JPTnpyzew9euTX8lUREREREQkDgWE5SiZ4Z7ZA0ZvWCV0gx12gLFjNyz+8N+x9/Dkc/8ip7pKU4iIiIiISOkpICxHJS5cHysjw0sRvPkmC2vWpsusL30V0g8/TEFPRUREREQkHSggLEelKlxfWLduHH7mECY02strkR1yCPznP168XEREREREpAQUEJax3EG9kioIvylz/37bui6nnnIL9O/vgeB118Ghh25UrFolKUREREREpCgKCCupdVWichRjx0K9er4aaevWTH7m9Q37qE6hiIiIiIgURQFhZdejB3z55YZVSMf/32PY+vWA6hSKiIiIiEjRUhoQmlmumX1tZlPMbFLUtp2ZvWVmP0Z/t43azcyGmNl0M/vKzNqlsm9blJ12gjfegFtvpdOcb6m2bg1V1q8jq4qVuk5h3NVORURERERki5JZDs9xcAhhYcz2AOCdEMIgMxsQbfcHegK7R7d9gAeiv5KMKlVgwAByDjyQu/tcyYy6jei0cAY5Lf8BTU7eaPe8YC+Z+Y4iIiIiIrJlqogho0cDj0f3HweOiWl/IrjxQB0z26kC+rfZK3Lhmn335erDL2OXxXPJmf45nHIKnHsuLF9evp0UEREREZHNXqoDwgC8aWaTzaxv1FY/hPALQPR3h6i9IfBzzLFzorYCzKyvmU0ys0kLFixIYdcr1qYsBvNH9Vr8/ZhrYOhQqF4dhg+HDh3g66/LrH8aUioiIiIiUvmlOiDcL4TQDh8OepGZHVDEvhanLWzUEMKwEEL7EEL7evXqlVU/NwuxQeAmrxBqBuefD599Bs2bw9SpHhTecw+EjS6riIiIiIikoZQGhCGEedHf34D/AR2BX/OGgkZ/f4t2nwPsHHN4I2BeKvtXEYoa7hm7ImiZrRDasiVMmgTnnQerVsGll8Lhh1NvmcpRiIiIiIiku5QFhGa2lZltnXcfOBT4BhgF9Il26wO8Et0fBZwRrTbaCfg9b2hpuohdETQrs0qpVwjdSM2aMGwYvPgibLcdjB3L2Ecuouv0CWVzfhERERERqZRSmSGsD3xsZl8CnwGjQwhjgUFAdzP7EegebQO8DswEpgMPAX9PYd82SzlNtt1wf8S5nQpsl4njjoOvvoJu3aj71x8Mf/EmuPBCWLGibJ9HREREREQqhZSVnQghzARax2lfBHSN0x6Ai1LVn8qmzIPBPA0bwhtvcFP387n6g8eoNnQovPcejBgBOTmpeU4REREREdksVUTZCaloVaowvMMxHHPGXbDXXvD999CpE9x2G6xbV2DXTVrYRkRERERENmsKCLdARdYpjDF1h6ZMfultuOQSWLsWBgyArl2ZPHHahn02ebVTERERERHZbCkgTEMFyls8+QWTr7wRXn8d6teHDz5g/NW3bChNUdxqpwoWRUREREQqLwWEaShueYuePX3BmSOPpNMPE6m+djW2fh1ZGRuvdlqm9RJFRERERKTCpGxRGdl8JSxvscMO8Mor5Dz4II/cMpAvGu5Jp1W/knPw9dAkfx2geAFlyhbBERERERGRlFGGMA0VWd7CDC64gOt6XMR+s74k54sPoVs3n2cYladIWb1EEREREREpVwoI01yizN7Muo04/rTb4aabIDMT7r0X2raF8eNTXy9RRERERETKhYaMbmaSWR20vKyrkgHXXQe9esEZZ8A338B++8GAAWSt68CajCwFgyIiIiIilZgyhFK8tm1h0iS4+mpfffSWW3jliX40W5Bb0T0TEREREZFNoIBQklOtmheu//BDaNqUFr/9xKjHL49bzF5ERERERCoHBYRSMl26wJdf8lSbnlRbFxWzP/BAmD69onsmIiIiIiIlpIBQSq5WLa7rcRF9TrwRGjSAceOgdWu4/35Yv76ieyciIiIiIklSQChx5Q7qVewCNx80zYGvv4ZTTvGSFBddBF27wowZ5dRLERERERHZFAoI01QyAV9SttsORo6EF17wwvbvvw+tWsGQIcoWioiIiIhs5hQQStk4/nj49ls49VTPFl52GZ81aclBfYdVdM9ERERERCQBCyFUdB9KrX379mHSpEkV3Q0p7JVX4IILYP58VmZWpfqgW+DyyyEjo6J7JiIiIiKyRTGzySGE9qU9XhlCKXtHHw3ffsuLex9C9bWr4aqrfHXSqVMrumciIiIiIhJDAaGkxnbbcWWvfpx1wg2+Eun48V7g/rbbYO3aiu6diIiIiIiggFBS7L1dO/jcwrPPhlWrvG5h587wzTdkDxhN9oDRm3T+sjiHiIiIiEi6UkAoqVenDgwfDmPGQKNGMHEitGvHpeOeJmvdmorunYiIiIhI2lJAKOXnsMM8W9i3L6xZQ7+PRzD60cvgk08qumciIiIiImlJAaGUr222gQcfhPfeY+a2Ddhj0WxfcObii+GPPwrsqiGlIiIiIiKppYBQKsZBB9HzrHu4r9OJXo7ivvugRQsYNaqieyYiIiIikjYUEErKTZ61JG77qqxq3H5gH5g8GTp0gLlzvWTFSSfB/Pnl3EsRERERkfSjgFBSIjYI7P3w+IRBIQCtWsGnn8LgwbDVVvD889C8OX/78g0IoRx6KyIiIiKSnhQQSkqMn7low/01a9cX2I4rIwMuv9wXnenZE5Yu5bax9/D0M/+EH39McW9FRERERNKTAkJJiU5N6264n5VZpcB2kZo0gdGjYcQIFtaszb6zv4aWLeGWW2D16oSHFZmBFBERERGRuBIGhGbWvjw7IluWnCbbbrg/4txOBbYL2yiYM4NTT6XbuQ/wwt5dvaD9tddC27bw4Ydxjyt2WKqIiIiIiGykqAzhQ2b2o5n928xalFuPZIsTLxhMJphbWmMbrup1Bbz1Fuy2G3z3HRx4IJx9NixcWPJhqSIiIiIiUkDCgDCE0BY4AlgHvGBmU8ysv5k1KbfeyRarJMHc5N1z4Ouv4YYboGpVePRRaNaMTt9+smHRmeKGpSp7KCIiIiKysSLnEIYQvg8h3BhCaAH0AeoA75rZuHLpnWyxiptjuFEG8de/YOBADwy7dYPFi8m5pA83v3Ef2y9bEndYqoaUioiIiIgULalFZcysCrADUB/YCliQyk7Jlq+4OYYJM4h77AFvvgkjR0L9+pz25Vg+feBMcu4fBMuXJ3eOQrIHjCZ7wOhNfUkiIiIiIpVOkQGhme1vZvcDc4B/AB8DzUIIx5RH5yQ9xJtjWGQG0QxOOQWmTeOJtr3IWL8e/vtfaNECRo1K7hwiIiIiIlLkKqM/A4OAqUDbEMKhIYRHQgi/l1vvJG0ltUppnTpcf+iFHHv6Hb4C6ezZcPTRcMwxMHt2iVY6FRERERFJR0VlCLuEEPYLIdwTQvjVzLYqt16JxCgukPuyQTP47DO4+27Yemt45RVo3hxuuYWqa9ckdQ4RERERkXRU1CqjswDMbF8z+w7PFGJmraNhpCJFyh3Ui9xBvcrnyTIz4bLLYOpUOOkkWLECrr2WNx75OwfNmFQ+fRARERERqWSSWVTmbqAHsAgghPAlcEAqOyVSag0bwrPPwttvQ/Pm7LLkFx57YaAPJZ05s6J7JyIiIiKyWUlqldEQws+FmtaloC8iZadrV/jyS24++GyWVa3hi820aOGlK/76q6J7JyIiIiKyWUgmIPzZzDoDwcyqmtlVRMNHRSpakcNSs7J4uONxHHzeg3DaabBqFdx4oweGL7+8oah9eVBpCxERERHZHCUTEF4AXAQ0xMtPtIm2k2JmGWb2hZm9Fm3vYmYTzOxHM3vWzKpG7dWi7enR49klfTFSuZTXHMMFtbaDJ5+EDz+EVq0gNxeOPRYOPxx++GHDfipcLyIiIiLpptiAMISwMITQO4RQP4SwQwjhtBBC/Arf8V1GwYzibcDgEMLuwBLgnKj9HGBJCGE3YHC0n0jZ2X9/mDwZ7rkHateGsWOZ3P14MtethRDo/fD4uEGhsnsiIiIisqUqqg7hXmZ2VMz2YDN7JLq1S+bkZtYI6AU8HG0bcAjwQrTL40Bekfujo22ix7tG+0uaKosM4kbnyMyEiy/2zOA55zC+QXPWm4EZa9asY/yMhZvYaxERERGRyqOoDOEgIPbbcQ9gNPAecH2S578buBpYH23XBZaGENZG23PwoahEf38GiB7/PdpfpOztsAM8/DCd/nUpmevXkbFuLVlrVtPp1gHw6acV3buNKEspIiIiIqlQVEC4Uwjhk5jtP0IIL4YQngS2L+7EZnYE8FsIYXJsc5xdQxKPxZ63r5lNMrNJCxYsKK4bIkXKObwLqzMy2XNBLiPeuouct1+Czp2hd2/4ufDiuiIiIiIiW5aiAsKtYzdCCJ1iNndI4tz7AUeZWS7wDD5U9G6gjpllRvs0AuZF9+cAOwNEj9cGFhc+aQhhWAihfQihfb169ZLohkgxrArf7rgbOZ+MhWuvhWrVYORIaNYMbryR6mtWVnQPRURERERSoqiAcJ6Z7VO40cw6kR/EJRRCuCaE0CiEkA2cDLwbQuiNDzk9IdqtD/BKdH9UtE30+LshlGNdAJFateDmm2HaNDjxRK9XOHAg7zx0IUd990FKy1RoSKiIiIiIVISiAsL+wDNmdoOZHRndBgJP4/MCS6s/0M/MpuNzBIdH7cOBulF7P2DAJjyHSOllZ8Nzz8EHH0DbtjT8cwFDXr0dunSBiRM36dQqbSEiIiIim5OEAWEI4TNgHyADODO6VQE6RY8lLYTwfgjhiOj+zBBCxxDCbiGEE0MIq6L2ldH2btHjM0v3kkTKyAEHwMSJXH3YpSyoWQc++QQ6doQ+fWBesUnyDWKDwESlLUREREREKkKRdQhDCL+FEK4PIRwf3a4PIfxaXp0TqXAZGTzX+lD2u/ARJve/GapWhSeegD32gBtvhOXLiz3F+Jn5ZTvXrF1fYFtEREREpCIVW5heJJ3lZfNWZ1ald1Y7Jr83GY491gPBgQNh993p3/NSml49KuE5OjXNr56SlVmlwHZp+yMiIiIiUhYUEEra26h4fYyNsnurqsNLL8GHH0KHDvDLL9w29h7GPHoJjBkTd+GZnCbbbrg/4txOBbaToSGnIiIiIpIqCghFipAwu7f//jB+PDz9ND/Xrk+zhbPh8MOhe3eYMiXh+UoaDIKGnIqIiIhI6iQMCM0sw8zON7ObzGy/Qo9dl/quiVS8IrN7VarAySfT9dyh3Hzw2VCnDrzzDrRr5wvPlFFh+7IccioiIiIiEquoDOGDwIHAImCImd0V89hxKe2VyGYoUXZvdWYWD3c8DmbMgCuugMzM/IVn/vlP+OOPpJ8j3nDQTR1yKiIiIiKSSFEBYccQwqkhhLvx8hO1zOwlM6sGWPl0T6QS2W47uOsumDoVTjoJVq6EW2+F3Xbj9M9fI3Pd2riHlWSOoIJBERERESlLRQWEVfPuhBDWhhD6AlOAd4Faqe6YSKW1667w7LM+x7BLF1iwgJveGsqbw//uBe/Xry+we3nNEcweMJrsAaNTcm4RERERqZyKCggnmdlhsQ0hhH8DjwLZqeyUyBZhn318NdKXXmLmtg1oumQe/O1vXtz+rbc2rEiqOYIiIiIiUlESBoQhhNNCCGPjtD8cQshKbbdEKp+4Qz3N4NhjOfSc+7mmx8XQoAFMngyHHgrdusFnn2mOoIiIiIhUGJWdENkEyc7/m3770dw69h748Ue47TZfkfTddz2LePzx7LrIVyRVMCgiIiIi5UkBoUgxSlS4vrj5fzVrwtVXw8yZMGAA1KgBL73Em8MvYtCYIWVWqkJEREREJBlF1SHMLM+OiFRGpZ7/t+22vgLp9OlwwQUE4OSv3oTdd4erroJFKj4vIiIiIqlXVIZwvJm9bGYXmFl2OfVHpFLZ5Pl/DRrAAw/Q7dwHeHXP/WHVKrjzTmjaFP7zH1i2DCg6SykiIiIiUlpFLSrTHrgs2rzbzCaa2WAzOzSqRSgiMTZl/l/udg255Oj+vuBMjx5ezP666zwwvPNOWLGiDHsqIiIiIuKKnEMYQpgVQhgaQjgG6Ay8CnQDPjIzFTQTKWvt2sHYsb7gzL77woIFPoR0113hnns8g5hAsnUGiyp8Xx5UD1FERERk85H0ojIhhDUhhHdDCFeHEDoCfVPYL5G0stGQ0IMPhnHj4PXXIScH5s+HSy/1OYbDhsGaNSU6f7KroYqIiIhIein1KqMhhLll2RERKcQMevaEiRPhf/+Dli19FdLzz4dmzeCxx2Dt2qROVeLVUEVEREQkLajshMjmzgyOOQamTIFnn4U994SffoKzzoK99oKRI2HduiJPUerVUEVERERki1ZsQGhme5dHR0SkGFWqwEknwTffwBNP+LzCH36A3r2hdWsO+34cFtbHPXSTV0NNkuYHioiIiFQuydQaHGpmVYHHgJEhhKWp7ZJI5VLu5SAyMuD00+Hkkz0w/Pe/4dtvGfrtt3y3wy6wzxrPKFaJ/3tPqoJBEREREal8is0QhhC6AL2BnYFJZjbSzLqnvGciUrSsLDjnHM8S3ncf82ttR4vffoLjj4c2beD552F9/IxhIsrwiYiIiKSXpOYQhhB+BK4D+gMHAkPMbJqZHZfKzolIEqpVg7//nf3PH84xp93O5Db7w9df+/DSli3hmWeKnWNYEbTSqYiIiEjFS2YOYSszGwxMBQ4BjgwhNI/uD05x/0QkCZNnLWFNZhZTGjand69rmDx4ODRuDN99B6ecAnvvzdHfvkfG+ooNDFX+QkRERGTzkkyG8F7gc6B1COGiEMLnACGEeXjWUEQqWIGyEuvWM771AfDjj16zMDsbpk3j/167k7cevhAefzzpchUp7afKX4iIiIhUuGQCwsPxxWT+AjCzKmZWEyCE8GQqOyciyYlbVqJqVTjvPJ9jOHw4s+rsSNMl8+DMM72O4SOPlLjAfUr6KSIiIiIVJpmA8G2gRsx2zahNRDYTRZaVyMqCs8+myYKfPTu4++4wc6YvSLPHHjB0KKxcWab9STQUtLzKX4iIiIhIcpIJCKuHEJblbUT3a6auSyKyKRIGWZmZcMYZPq/wqae8wH1uLlx4IeyyC9x+O1utWlHq5y3p/EAFgyIiIiIVL5mAcLmZtcvbMLMc4K/UdUlEUioz04vZf/MNPPusl6iYPx+uvpoHXr6F/X76gslf/VTi02p+oIiIiEjlk0xh+suB581sXrS9E/C31HVJRMpFRoaXpjjxRBg7lsn3PE7f5iexOiOT3o9/wYjwODn9zoVGjZI6neYHioiIiFQ+xQaEIYSJZrYn0AwwYFoIoXxXohCR1DGDnj0ZX313Vo6dBmasCYHxH39HTtOmPsy0f3+fe1gEzQ8UERERqXySKkwPdABaAW2BU8zsjNR1SUQqQqemdT04BLKqZtJp93penmL4cJ9vePLJMGVKUudSMCgiIiJSOSRTmP5J4A6gCx4YdgDap7hfIlLOCmT4+nYm54n7YNo0X400I8PnG7ZtC4cfDh99BCFUYG+Llj1gNNkDRld0N0REREQ2e8nMIWwPtAhhM/72JyLkDupVZufaEBzusQc8/DAMHAh33gkPPghjxvhtn33gH/+AY47xgLEC+ikiIiIimyaZIaPfADumuiMishlr1AgGD4ZZs+C662C77WDCBDjhBB9O+sAD8JcWHxYRERGpbJIJCLcHvjOzN8xsVN4t1R0Tkc1P9p2fkb2mE8yeDUOGQHY2TJ8Of/87NG7MZR+PZNsVv1d0N0VEREQkSckMGR2Y6k6ISCWz1VZwySVe1P6ll+D222HSJK5YOJILJrwIaz+Gfv1gt90quqciIiIiUoRiM4QhhA+AXCAruj8R+DzF/RKRyiAz02sZfvYZvPce9OpFjbWrfAjpHnv4kNIJEyq6lyIiIiKSQDKrjJ4HvAA8GDU1BF5OZadEpGLkDupVukVfzOCgg+C11+Cbb+CsszxYfPFF6NQJDjgARo2C9evLvM8iIiIiUnrJzCG8CNgP+AMghPAjsEMqOyUildhee8Ejj0BuLgwYALVre5mKo4/2BWjuvReWLUt4eFmWjJg8a0mZnEdERERkS5VMQLgqhLA6b8PMMoFiS1CYWXUz+8zMvjSzb83sxqh9FzObYGY/mtmzZlY1aq8WbU+PHs8u3UsSkc1CgwZw663w889w113QpAn8+KPPPWzUyEtWzJpV5k8bGwT2fni8gkIRERGRIiQTEH5gZv8EaphZd+B54NUkjlsFHBJCaA20AQ4zs07AbcDgEMLuwBLgnGj/c4AlIYTdgMHRfiKyGSpRkLX11nDFFb4a6QsvwH77we+/wx13QNOmcOKJMG5cmRW6Hz9z0Yb7a9auL7AtIiIiIgUlExAOABYAXwPnA68D1xV3UHB548KyolsADsHnJAI8DhwT3T862iZ6vKuZWRL9E5FysMmZt8xMOP54+PhjX4Smd2+oUsWDxC5doGNHGDmSrHVrNqmfnZrW3XA/K7NKgW0RERERKSiZVUbXhxAeCiGcGEI4Ibqf1E/5ZpZhZlOA34C3gBnA0hDC2miXOfgiNUR/f46ecy3wO7DRNzkz62tmk8xs0oIFC5LphoiUgTLNvHXoAE895UNGr70W6taFSZOgd28+GnoOf//0OVi4MO6hxc0xzGmy7Yb7I87tVGBbRERERApKZpXRn8xsZuFbMicPIawLIbQBGgEdgebxdst7qiIeiz3nsBBC+xBC+3r16iXTDREpAynJvDVoADff7IXuh/qz+/MAACAASURBVA2DFi3Ycdlirv7wCdh5Z+jb11ctLSUFgyIiIiJFS2bIaHugQ3TbHxgCPFWSJwkhLAXeBzoBdaKFacADxXnR/TnAzrBh4ZrawOKSPI+IpE5KM281a8J558E333DU6Xdy5eGXM3m7bHjoIWjZ0ktavPACrNm04aQiIiIiUlAyQ0YXxdzmhhDuxucBFsnM6plZneh+DaAbMBV4Dzgh2q0P8Ep0f1S0TfT4u8kOTRWR8pWqzNvk2Uv5qkEzXmzZjd5n3snki6+BrbaCDz7wxWd22YVLxj1NvWVaOVRERESkLCQzZLRdzK29mV0AbJ3EuXcC3jOzr4CJwFshhNeA/kA/M5uOzxEcHu0/HKgbtffDF7MRkTRSYJ7i+sD4486GuXNhyBBo1gzmzuXKj0cw7oGzfFGaTz8ts9VJRURERNJRZvG7cGfM/bVALnBScQeFEL4C2sZpn4nPJyzcvhI4MYn+iMgWKu48xdq1vXbhxRfDO+/w5t//Rdfpn8HIkX5r29YfO+UUqFGj3Pqat7BN7qBelfo5REREJL0VGxCGEA4uj46IiBQ5T9EMunWj73Gr2GblMh6tPoOch++CL76Ac86Bq67yvxdeqABKREREJEnFBoRm1q+ox0MId5Vdd0REXLx5inm1D/+oXoveWW0Z8em35Hz6Btx7L0yc6MXu77wTDj8cLroIDj0UMjLKu+siIiIilUayq4xeiNcJbAhcALTA5xEmM5dQRKRMbFQLcc6fcMYZXuh+wgS/n5UFo0d7ULjrrvCf/8AvvxQ4T3G1DEVERETSRTIB4fZAuxDClSGEK4EcoFEI4cYQwo2p7Z6IbE5yB/Wq0OGYRdZC7NgRHn8c5syBW26B7GwvfH/dddC4MZxwArz1FqxfX/4dFxEREdlMJRMQNgZWx2yvBrJT0hsRkSIkVQuxXj245hqYMQPGjIFjjvGVSF980YeQ7rEHF4x/gbrLl5Zjz0VEREQ2T8kEhE8Cn5nZQDO7AZgAPJHabomIFK3YWohVqsBhh8H//gezZ8O//w077wwzZjDgg8f49P4z4eST4b33Nql0Rd68xlQqj+cQEfn/9u47Turq3v/46zO7S+9IW3ABKSJNYBBWjb0gakQTO2o0RqM/NZqruSHxxmiChntjbNfEmlwsaw32rogVRmGRKkW6dJCyIm3L+f1xvttgZ+vMzizzfj4e5zGzZ7/znfP16+zum9NEJDVVZ2P6O4ErgK3ANuAK59xd8W6YiKSmuAxLzcyEP/wBli+HN97g/d4jSHNF8PzzcOKJ0K+fX4xm8+Zqna5sQBv7eCQuga0+3kNERESkOj2EAM2APOfc/cBqM+sZxzaJiMRHWhqccQZX/fQ2fnTNP+G223xYXLzYb1vRtavf8H7KlErnGu63uE2Zr2OlPt5DREREpMpAGAwT/S3wu6AqA3g6no0SEYm3da06kPvzm/zCM6+8AqNHQ36+3+z+xBNZ2b4rjB/vF6nZR6WL28RIfbyHiIiISHV6CM8BzgJ+AHDOrUXbTYhIA7XfUMw138OYMfDWW7BsGdx2G2tadqD7tvV+mGn37n4Li0mTYK9fX6tai9vUUX28h4iIiEiVG9MDe51zzswcgJk1j3ObRESiquv8woqGYpaErR494I47OGZXmB+tmMWTzPO9h2+/7UuHDnDppXDllSXniBbUivc5jMV8SIVBERERiZfq9BC+YGaPAG3M7CrgA+Cx+DZLRCQ+qjMUsyiUxieHhOGFF2DtWrj3Xhg4EDZtgnvugQEDeOmpm7lw1jvw/fe1bkuPcW+WBEcRERGRRKjOKqN3A/8GJgGHArc55/433g0TEYmHGg/FPOgguOkmmDMHvvgCrr4aWrZk2NpFTHj3QejcGa64Aj77rE7bV4iIiIgkQqWB0MzSzOwD59z7zrnfOOducc69X1+NExGJpxoNxTSDESPgkUdg3TpuPv3XfHHwQNi5EyZOhGOOgb594c9/hhUr4tVkERERkZiqNBA65wqBnWbWup7aIyKSNKLu/de8OZMGncQFF0/wW1aMGwddusCSJX4ri549efbZ33HenPfrNKRUREREJN6qs6jMbmCumb1PsNIogHPuV3FrlYhIguy7Cmm0YaXlFov5y198z+DkyfDEE/Dyyxy5ai5HrpoLnR+Dn/wELrsMTjzR74VYTbFYkCaWi9uIiIjIgac6gfDNoIiIHPAqXYW0MunpMGqUL9u3c+XP76bDD1s5b+5kwk8/DU8/Dd26wSWXwM9+Bv36xfEqRERERKonaiA0syzn3Crn3BP12SARkUSKxYbwuduKmNwnG4BXwqeTkzaf8NMPwfLlMGGCLyNGcFrmiXx0yHByV27V1hIiIiKSEJXNIXyl+ImZTaqHtoiI1IsVE86IOoQyFhvCl+tlLHJETjjHzy/85BO/h2HLluSuzuOjQ45gb1oGY//+CbkTX4Ldu2t+MSIiIiJ1UFkgtDLPD4l3Q0REkk1te+0q7GUMhfxKpI8/DuvXE/nNePampVMUSiPfQeSR56BTJ7+FxQcfQGFhrC5DREREJKrKAqGL8lxERCpRZS9js2Zkn3MiRaE0cI6MtBDZzfZCXp7fwuKUU/x8w1//GmbM0P6GIiIiEjeVBcLDzSzPzL4HBgfP88zsezPLq68Giog0ZNF6GUvqzci59keEJ78CCxf6bSt69YL16+G+++CII/wCNHfc4YedltFj3Jslq4hWJeoWGiIiIpLSogZC51yac66Vc66lcy49eF78dav6bKSIyIGsJBweeqgPft98A198Ab/6FXTs6Pc6vP126NMHRo6EBx6ADRuqPO++W2jUJhTWJHSKiIhIw1PpxvQiIqmoskVn6oUZjBgB998Pa9bAu+/6fQxbtIAvv4Qbb4TMTJ58/g+cO/cD2LatwtNUtIWGiIiISFkKhCIiySw9HU491W94v2EDPP88nHUWpKVx7IqvuPut+3wv4o9/7Pc6zCsd0R+LLTRERETkwKZAKCISB3HpZWzWDM4/H159Fdav53ejrufz7oP9iqRvvAGXXurD4TnnwHPPEW6fUfLS2m6hISIiIge2qBvTi4hIEmvXjmeHnMazQ05j0gX9CE971/cefvopvPKKL02b8mBWmDf7/YhwhxMS3WIRERFJQuohFBFJkLr0IpZbMOalxeSecRF8/DGsXu3nHh51FOzaxZmLPuOhVyf4nsOLL/ZBcffuWF1CzGjxGhERkcRQIBQRaYCiLhiTmelXJ/38c1i5kvEn/JxZXfrCDz/As8/64aSdOvlFal5/vdrhUNtWiIiIHJgUCEVEGqBqLRiTlcXjI37C2ZfdA8uWwYQJMGyYX3jmqaf84jQdOsBFF8GLL8KOHeVeHottK0RERCS5KRCKiDRAZReIqdaCMT17wm9/C7m5fl/D8eNh6FAfAp97zi9W06GD70F8+mnYtq3a21ZouKeIiEjDpUVlREQauMrCYIVzFPv0gVtv9WXZMnj5ZZg0CaZNK12QJiOD7B+PJa3XORRaiIyMtFpvW1EcFhO6t6OIiIhUSIFQRCSVHXII3HyzL2vW+HD40kvw8ceEX5rIC5kRpmUN5siD0gi/vtT3IHbtGrfm5K7cWmHAVagUERGJDw0ZFRERr2tXuP56+PBDWL8eHnuM7xs35+ovXyL8xrNwww3QrRsceSTcfbfvXYwBzVUUERFJHPUQiojI/jp0gF/8gsuXdKHV7h3MGZ7vew7ffhsiEV9+8xsYMIBbWg3g/d4joWg0hGr+74wVzVWsck6kiIiIxIQCoYhICqvOUMy8Ji3gkjPgkkv89hVvv+3nHL71Fsyfz5GZhZhz5A45lvBRA2HMGDjxRGjcuFptqNaKqTFQH8NONbRVREQaGgVCEZEGKiGho3lzOPdcX/buJfe1KVzwxS6KLESjwgJynruV8COPQIsWMGoUjBlD8z2N+aFx86jzA2u8YqqIiIjEjOYQiohIVCsmnBE9eDZqRKR9LwrSMigKpZHfqDGRy2+Eww/321lMmkTuuLsoDKUTKipk7EOfkTttfqXvV1UY1PxCERGR2FIgFBGRWttvuOfVF8CsWbB8OTzwAJETzmZvWroPjIVFRP7jDhg8GP7wB/jySygqqvI9tOiMiIhI/CgQiohIrUNW1OGePXrADTeQfed/UmQhcEVkGGRvWgJz58L48TByJHTpApdfzuiFn9Fyzw8VvkdFi87UhQKliIhIqbgFQjM72MymmNkCM5tvZjcG9e3M7H0z+yZ4bBvUm5k9YGZLzGyOmQ2LV9tERCT2PW9R5weagYXIue5YwvOnwTvvwLXXQlYWbNwITzzBQ69OYOYDF/vFaO65BxYtAueA2Cw6U9+9jAqdIiLSUMSzh7AAuNk5dxiQDVxnZv2BccBk51wfYHLwNcBooE9QrgYeimPbRERSXqx73qoS7t7Wrzw6ahT84x+wYgXMmQN/+Qv86Edk4GDKFLj5ZujXD/r0gRtvJLx4BukFBUDtF52pj2vV0FYREWmI4hYInXPrnHMzg+ffAwuArsAY4IngsCeAs4PnY4AnnRcB2phZl3i1T0Qk1dXXdg9RmcGgQTBuHHz6KWzaBM88A2PHQrt2sHQpPPAAnHoqsx+4kEdeGk/4/Umwdm2N36o+rrW+A7aIiEgs1Mu2E2bWAxgKfAF0cs6tAx8azaxjcFhX4NsyL1sd1K3b51xX43sQycrKimu7RUQOZEm33UO7dnDRRb4UFkIkAm++CW++SfM5cxj1TQSuivhjhw2D0aPhtNMgO5se//UuEH0rjvq41oQHbBERkVqIeyA0sxbAJOAm51yemUU9tII6t1+Fc48CjwIMHz58v++LiEjN1SUgxWU/xLQ0OPpoX+66i6P+3/9xwtIZ3NloFUyeDDNn+nLnndCqFQ91HsjHPYfBqkF+bmIl4hV8ky5gi4iIVENcA6GZZeDDYI5z7qWgeoOZdQl6B7sAG4P61cDBZV7eDaj5uCAREUkqsQiMa1t1JGfo6dw54QzYtQs++gjefdcvULNoEaPzpjJ68VTo/iAcdpjvORw1Co49Fpo2rftF1JDCoIiINBRxC4TmuwL/CSxwzt1T5luvAT8DJgSPr5apv97MngNGAtuLh5aKiEhqKxcqmzb1w0VHjw6+uYLfX3M3xy3LZdT6+bBggS/33gtNmsDxx8OoUfT6rglL23WrdRt6jHtz/7aIiIg0cPHsITwauBSYa2azgrrf44PgC2Z2JbAKOC/43lvA6cASYCdwRRzbJiIiB4oePXhmyGieGTKaSVeNILxmge85fOcd+OqrkueTgdWtOsK2s30P4kknQatWiW69iIhIQsUtEDrnPqPieYEAJ1VwvAOui1d7RERkfwdCb1e57R4mzvDz9+46Fu66CzZsgPfe88NL332Xbps3wqOP+pKWBiNHwimnwMkn++cZGXFtq3oZRUQk2dTLKqMiIiLxUtF2DyVz+Dp1gksv9aWoyC9EUzz3cNo0mDrVlzvugBYt4LjjfDg85RTo399vjVEDCnoiItLQKBCKiEiDVu3tHkIhGD7cl1tvhbw8vzjNBx/4smBByTYXAHTu7MPhySf74aWB3JVbtWiMiIgcMOK2Mb2IiEh9qPV2D61awVlnwQMPwNdfw+rV8MQTvjexSxdYvx6efhouv5zckafQuGAvVlTE2Ec+J3f+t1Wfv5Z6jHuzZGipiIhIvKmHUEREDhh16rnr2hUuu8wX53xIDHoPI3mtyQ+l4UIh8vMLiFx1C2G+9T2HJ5wARx5Zo+0t1MsoIiLJQj2EIiIi7NMzZwYDBsCNN8Lrr5P9r3txGFZURAaO7G/n+TmI48f7UNimjd/e4vbb4eOPYc+e/c5fbvGbxyPlvhYREUkUBUIREZEqhHt1wIVCuFCInOuPIzx/Grz2Gvz61zBkCOTn+yB4xx0+GLZp44Pi+PHw2Wewd2+Fi9+IiIgkmoaMioiI1EDJUM8f/9gXgC1bfCD86COYMgXmzoUPP/QFoFkzsk/9Kem9f0pBKI2M9LToi9+IiIjUIwVCERFp8BK+3UO7dnDOOb4AbNpUPiB+/TXhV57i+czpRLIGkb1pCeHcTn7+4QknwNChkF7+V3IyzzPUfooiIgcOBUIREZEyYhLEOnSAc8/1BWDDBvjoIxZMeILTFk+j15bVsBS/HyL4PRCPOorco08jVNSbIjPGPh6p2aqpIiIitaA5hCIikvLivuBLp05wwQX816jrOOmqh2HNGsjJgauugt69YccOeO89Im9PAxyYkb8nn8jDz/jQmJcX0+Ykw9YWydAGERFRD6GIiEiFC77s2zMX0+GRmZlw8cW+AKxbB59+Svbnc8koKqTAOTKKCsnOeRAmLIRQyC9ec8wxpaVjx9i1R0REUpZ6CEVEJOWVXeAlIz1U/wu+dOkC559P+P4/sye9EelFheSEGxG+5Cy/x2EoBDNnwv33+2GonTpBv35w9dXw1FOwYoXfO7GeaesMEZGGTz2EIiKS8sr2BsZz3l51exn3ZDQmfPEZcPGZvmLnTvjiC/jkEz77v5cZtnYhzRYtgkWL4LHH/DHdusHRR8NRR/nHwYMhIyPm17Dv8FrNcxQRadgUCEVERMpIynDTrFnJiqSX7BpOemEBS87tDJ9+Cp984vc6XL0ann/el+LXjBhRGhCzs/1qqGXUZgGd6gyvrYlkXk1VRCQVaMioiIhIA1OQlg4jR8Itt8Brr8HmzX7vw0cf5cWBJ7OsbabvVfzoI7jrLjjjDGjfHvr3J/ea/yStqBCci7qATmULvsRieG3cF/EREZFqUw+hiIhIEqnusNJyPWuhEAwcCAMH8pulmf48N4+AadNg6lT4/HOYPh0WLCDSeiDumB+VrmR6+72E+zb1vYjDh/uexUrEYnhtrHsZRUSk9hQIRUREGogazd/r0AHOOssXgL174auvyP5oJmmbijAgo7CA7PdehIkL/THp6TB0KH8s6sysLn1hyaHQqxeYVfgWtQ1xCV/ER0RESigQioiIEONtJeKkTj1rjRrByJGER44kf9yb4BwvnNePcPgPpb2Ic+bA9OlcUfyaN/7mh5qOHFlaRoyo83XU1yI+IiJSNQVCERGRBqImPWtVLtZiRji7P2T3L90P8fvv4csvuflvb1AUMi5ZESE8fxq89ZYvgcntuvoexINX+JA4eLAPnLUQrY3FcxgbQlAXEWnIFAhFREQaiKp61uq8JUTLluT2HsakwbsBeHvoKeQ83p3wyrl+24svvoCZM+m1ZQ29tqyB66f41zVuDMOG+XCYne0fu3ePOtRURESShwKhiIhIA1RR0IvFYi37nWNXI8IXXAAXXBBU5vuhpZFIaUhcvNgvYDNtWumJOnb0wXD48NLSsWPNLlJEROJOgVBEROQAUd0hpZUNw6zyHBkZEA77ct11vm7rVvjyy9KA+MUXsHEjvP66L8UOPrgkHK44cbg/Rx1oWKmISN0pEIqIiBwgYrFYS63O0bYtjBrlC4BzsGyZD4kzZvgycyZ8+60vL79c+tqePcv3Ig4bBm3alHxbG9eLiMSXAqGIiEgDUt3esFiEqFqfw4wejy0EWrHib3/zdYWFfmhpcUCcMQO++gqWL/flxRdLXp478mQyjrmOglAaYx+ZSs6lQwgf1rXGzVAPoohI1RQIRUREJP7S0uCww3y59FJfV1AACxaUD4mzZhGxthRaCGch8vMLiPz814S3z/O9h0OHlhYREakzBUIRERFJjPR0GDTIlyuC3Q/37iX746+w9zdgRUVkUET2mq/h2wU+PObklLz84T5H8k7fI8n940LCw/v6wJiZWaPVTdWLKCKpToFQRETkABKLYBPLcFTjOYCNGhE+ZSSFk31Qy7n2eMJ/nuFXNv3qq5KSu3E3N/34FvampfPOjgJyrvlPwmsXQocOJT2IZy5wzO/UC4qKIBSK2TWJiBxIFAhFREQkpuq8H2IZJa8bMcKXQGTyYna/txjMyDcjcvwYwm+th02b4L334L33uCKzH5GsQeT2fZhw52Y+KBYPO+3fHxo1qtN1iogcCBQIRUREJKZisR9iVbJ7d4D3vwEgo1E62Xf9Fp7+C6xc6XsQc7/hwt19KQyFaFRYQM5ztxL+/MHSE2RkQP/+3NBiAJHug8n9dzrhY4ZAp041aoeGnIpIQ6dAKCIiIjFV3f0Q6yLq9hg9ekCPHkTaLCH/3UUA5KelEfnjPYS3zC4ddvrNN+Ru2sNjJ/2UvWnpjJ32Azk3Hk+4YAsMHgyHH+4fBw/2C+E0bhzza6iuWIROBVcRiUaBUERERGIqFvsh1iS4VHT+/ULpqGzoPrr0gB07iLw0nd3zd4CFyAcifYYT/vhp+OADXwL5oTQy+h9WGhCLS2ZmtdqnMCYiyUyBUEREROImUZvKVxlKW7Qg+7gh8PVUADIaZ5A98X6w8TB7tl/EZs4clk6eRs8ta2DePF+eeab0HO3bc3tWNp8cEib371sIj+jnexNbtKhRWxUYRSSRFAhFRETkgBYtlFYcGttB9+5w1lkAnDTuTZrk72bhhVklIZE5c2D2bHIbd2DCCT9nb1o6U5cVkHPXZX6l0549YcAAGDgQBgyg/4atLG3frc7XUeMVW+N0DhE5sCgQioiISIMUyx61ykLS7owm5HbsTfjKI0ornSPy6kx2T1tXutLpiJMJv7kUli/35Y03APhzZj+mZQ0mN/sfhDNblguL9O1babtisWJrLFd9FZEDjwKhiIiIxFwyDH+saxsqDVJmZA89BCLrgWCl0/v+BC/cA0uW+OGl8+eTu2QDF3c+lfxQWulqpy+/XPom6enc1+doJvceQe4fFxI+POhd7NUL0tNjsmJrfaz6KiINlwKhiIiISAWqClJR5ykedpgv551HZMoS9pRd7fQ34wlvmgnz58O8eeTuzmDc6F+xNy2d93cUkHPDrX7YaePG0K8f2YcfQ3rHkykIhchISyP74JY1vo76WPVVRBouBUIRERFJWZX1ItYkSEXrcdvvHOecCN1/WlIXeW8BuycvLR12esq5hKc8CatWwezZhGfP5vnMD4hkDSJ71VzCExb73sN+/UqD52GH+a9bt66ybRouKiL7UiAUERERqUAsglRV58g+tDN8uAwIhp3e8R8w8c+Qlwdffw0LFjD9kTfpu3kl4WYF/kXffOPL66+Xf7MuXfYPiYcd5usraM++tNqpSGpSIBQRERGpQrQgVdf9EqMGxlatIDsbsrOZsKhj6Xvt3u3D4MKFsGBBaVm0CNat8+XDD8u/SatWvNy8i1/ptO08HxIPPdSvhtqoUbXbr8AocmCKWyA0s38BZwIbnXMDg7p2wPNAD2AFcL5zbquZGXA/cDqwE7jcOTczXm0TERERSTbVCp1NmsCgQb6Ucch/vkZm3iY+O6tLaUgsDo1btjA0L4+h6xbBuMmlL0pL86Gwb1/o25dLvt7L0nZdYfXhkJkJoVDMr1GhUiT5xLOHcCLwIPBkmbpxwGTn3AQzGxd8/VtgNNAnKCOBh4JHERERkYRpKMGlKJTG6jad4fTTfSnmHGzatH9I/OYbWLHCr4i6ZAm89RbnZPYjkjWI3JFPEd62Cvr08WHx0EOhb1+GrN3EsnZdo7ahPsKeAqVI7MUtEDrnPjGzHvtUjwGOD54/AXyED4RjgCedcw6ImFkbM+vinFsXr/aJiIiI1Ieqwktcw40ZdOzoy3HHlf/e7t2wdCksXkzu199y4bYsCi1Eo6ICcp69lfDs2TB7dsnhfygOjAPuJdwmrVxYpG9fGufvYU9G4/hdi4jERX3PIexUHPKcc+vMrGNQ3xX4tsxxq4M6BUIRERGRaspdubX6i980aeL3PBwwgEibJeSX3R7jsRcIH7QLFi/2gXHFd1zU8WQKyu6nOHViudMtwFjb6iCYPgh69/aroRY/9uoFLVrUrp0iElfJsqiMVVDnKjzQ7GrgaoCsrKx4tklEREQk6eWu3FryfOzjkVqtiLrf9hgDukH3tjBiBACRKUvYWzYw/ul+wvnLSgIjixczY2c60w8eSPbCuYT3XdgGyB10FI1OvYX8UDpjH/qMnN67CB/e04fFdu2A6g8JrUugrOo9YjEsVUNbpSGp70C4oXgoqJl1ATYG9auBg8sc1w1YW9EJnHOPAo8CDB8+vMLQKCIiIpIqIsu+K3meX1BEZNl3NQ5LVW6PsW9gPPkI6H5qSV3uyq2c/4/PAWiSZuS0X0v42/l+SOqSJbBsGZHmXSkIpeFCIfILC4j889+EIy/6E7RtC71788COpqxo0wUmbvK9i717Q6dOYBaT4JtMFBolWdR3IHwN+BkwIXh8tUz99Wb2HH4xme2aPygiIiJStf3CWpmva6NG22MEIsu+8/MVgXwHkWEnEL75qtIDCgvJzl0M/16CuSIyQkb2Ie1g1+E+MG7dCtOnc1bx8dOeL31t8+bQqxeR4edg7cM+UOYXEfliAeGOQ6Fp03JtqY9exvqgwCj1JZ7bTjyLX0DmIDNbDfwRHwRfMLMrgVXAecHhb+G3nFiC33biini1S0RERORAUlVYi+f7FasylKalER5xGEUvLfPtvPYowt2D+OccbNgAS5dy+X3vkx9K5z+abiT8dcT3MG7ZAnPmkL15L40vHEy+SyOjqJDsX18JFyyELl389hlBOW9OHqtbd4aVA6FbN7+9RqAmvYyxCIzJHjpFIL6rjF4U5VsnVXCsA66LV1tEREREUkFdwkddVkOtaSgt930z6NyZ3D2N+ajXNgByM0Lk3HGHP27LFli6lPDSpZz+1GfsaNyUXxauItxkL6Snw7p1vkydCsCFxauhHvUk4Y1LICurJCxGMkeC6wxmwfDazeXaEothqTU9R1WhUaFS4i1ZFpURERERkQNAbcNL1LmQ7dr5csQRvDSrJQCPTrjTH1hQAGvWwPLlsHw5uUs2cuHeQ4PtMwrJefb3hJcthGW+ZzI78zOaXHgn+aGgl/Gq86FpQWlgt4WQCwAAFF5JREFU7DAU3EFRA2OdrqOMqkLjgTZfsioaHptYCoQiIiIiDVwy/CFdnTZUdkx15kLu9/r0dOje3ZfjjycyZZ/tM/7vJcJZ5gPh8uWEly/nhs+fY0OLdoz5Npfw0ln+PPPm+TZk9isfGH9xHjTa48+flVX6XsUlM9O3oYbXUVVojMVCQfVFYa7hUyAUERERkYSLxVzI/cLYoZ399hn9+pXU/zUIMH96/2G/mE3Qu8jy5YRXruSmz3JY3/IgfvxtLuFls/2LFi6s+A3T0qBr13IhMZyVxZA137OqTWceu/7EWs25jPVCQSKVUSAUERERkaRS296wGofKtm19GTaspGpCEBhv/+AR2L4dVq2ClSt9CZ7P/HQWXfM20WnHFl+3ahV8+mnJOV4pfvIg0KFD+d7FrCzC3brRZ9NG1rY6iCd/uX9orO+FgmJBi/A0XAqEIiIiInLAiRYsajS0sXVrGDTIlzJ+EoTGSb8YTth2lAbGIDR+/mEuXbdvoscPm2HTJl9yc8ud4/3iJw+m+6Gn3brBwQeXPJ62aAPrWh5EOGMoFLYqt1pqfalqOGgiFuGR2FMgFBEREZGkkAzz0Kqzf2GxsU/M9AHmxF7ljjm6+ElREaxfX76X8dtvfVm92j9u3Fjay1jGw8VPnrrZz1Ps2rVcaLz9qzzWtTyIR35/tq/v1AlCoXLniPf8vljMdWxI8yUPVAqEIiIiInLAiHeorFGACYV8719mJmRnV3zMnj1+pdTigFg2LBY/btpUGigDtxc/eeUu/1gcGot7Gbt148qZ21jfoj181tp/r0sXaNKkxtccbShnLOY6xuIcWtimbhQIRURERESqKeYLvjRuDIcc4ks0u3f70FgmJI5/bwnfNW/NJT8sJTxvKmzevF9oPL14T8YLriK8NlgYp107Hw6Lg2rx87J1nTqRuzqv5DzRhnJWZ65jVWGtIc6XPNAoEIqIiIiIVFNCAkyTJtCrly/4HrvHt08F4O2MU8iZOJFwxyalofHbb8ldtY0LtmdRaGk0doXkTL6X8NypsGWLL3PnRn+/UIjIyZcTOnwMRaE08vcWEPnnvwn3TC8fHtuWXnss/jtoUZrEUCAUEREREamFRIWPioet9obeQQEiU5ZQULwno4WI/O2fhI87xPckrlkDa9cy7n/fofOO77ipfwtYu7akno0byZ43lUYDzyDfOb8n49P3wNp9tt9o0oSPm7RhY/N2sOT//JDUzp1LS5cutP1hO9uatqw0rNVlqKcWpak7BUIRERERkQakOsNWKzwmFIKOHX0ZOpTngp0ybto3kO3dS3j9erLueB0HTDgsg/BlZ5cGxuLwmJdH993r6b5tPUz6er825Gb2Y9eFd2I4xv7vFHKm/4twk/zywXHfENmqFZhV+7/FgbYoTSLmQyoQioiIiIjUQKIXL6nOsNWaDG3dr/euUSPIymJxhx7+XLdEud4dO3w4XL/el3Xryj2PNOvL3rR0P+zUOSJ7mxH+5MXKL65p04rDYtnnxaG2SZOYzelM5YVpFAhFRERERBqo6vSGVXRMdYZaVhmOWrSAvn19qUD2yq0UPTQVnCOjcQbZf7oZ8i+PGiBZtw527oTly32pSqtWhDt2ZHzz3sw4eACXZmUQ/lekNDB26lT6vG3bGvU87qs6gTGWobI+50MqEIqIiIiINDDVCR2VHVMfQy1LzmdGzlVHVu/8O3bsFxT//uxndPhhK+d3Tfd1Gzf6kpcHeXlcwhIumf1O5edNTy8Nh/uGxU6dGLp6E2tbdyR3+iLCg3v41V/rWaLmQyoQioiIiIikmJhvn1GFagebFi2gTx9fAn9d61dX7XXtUaXncQ62bi0Nhxs2VP48L690/uM+cjP7seDCO9mbls7Y5+aTc/bZhH9YVy44/vK7Vszv1Ivcv64mnNkCOnSAgw4qfWzatPR8tezdS9R8SAVCEREREZEUU1/bZ9R1+GTUXjMzv6diu3bQr1+V5+lzyys0KdjDxNEHEy7cVi40Rna2L53rCES6DyY8bSFs3w7ffENuZj+eCALjjPUF5Nx3a+m+jsXtPORwGp99G/lp6Yz9+8fkfPcR4dYhHxbLBsfix3bt/CI/ZdR3SC+mQCgiIiIiksKSeVXOWPSa5a7cSn56BvnpGYz9bDs5vziS8Oml5yiZ6wh+ruMz/4CW/yjpYYzM3sLutRlgIfLNiJw5lvCSKX4Lj82bYdMmIh37kl8cKgsdka/XEo5UsoBOKORDYZmQGO7QgYtXpLOoQw9+f2w3wgu+gE0HQfv2/rgWLeo0DzIaBUIRERERkRTUEFbUjEWvWVWhMmpvafv2cNhhZPfcCsWBsVE62b+/Drr/V+kbOEf2gjW4J77yC+hkpJN92Vnw0xGwaVNJaHxm5V6Wte/K6LVzCC+aURooF5b2Nt5V0pAKLiQjozQctm9f+ryOFAhFRERERCQpxWJoa01CZa228DAj3L8bzmb7Y3551H7H5K7cyu+DUPl0xrnkXD6ccPPCcr2M+z1u2gTffVdadu4sXWwnhhQIRUREREQkadW1J7M6obK671FZIK3Rqq6rthM+obffV7G6du0qHxCLy7XXVv8cFVAgFBERERGRlFDb+ZJ1DaUxWTCmaVN6PDh7//YoEIqIiIiIiESX6PmSsV7VNZYb14eqPkRERERERERiobZBbt8tOMp+XRfqIRQREREREYmzuvZSxmvjevUQioiIiIiIJLl4bVyvQCgiIiIiIpLkYj0PsZiGjIqIiIiIiDQA8VgcRz2EIiIiIiIiKUqBUEREREREJEUpEIqIiIiIiKQoBUIREREREZEUpUAoIiIiIiKSohQIRUREREREUpQCoYiIiIiISIpSIBQREREREUlRCoQiIiIiIiIpSoFQREREREQkRSkQioiIiIiIpCgFQhERERERkRSlQCgiIiIiIpKiFAhFRERERERSlDnnEt2GWjOz74FFiW6HVOogYHOiGyFR6f4kP92j5Kd7lNx0f5Kf7lFy0/1Jfoc651rW9sXpsWxJAixyzg1PdCMkOjOboXuUvHR/kp/uUfLTPUpuuj/JT/couen+JD8zm1GX12vIqIiIiIiISIpSIBQREREREUlRDT0QPproBkiVdI+Sm+5P8tM9Sn66R8lN9yf56R4lN92f5Fene9SgF5URERERERGR2mvoPYQiIiIiIiJSSw02EJrZaWa2yMyWmNm4RLcn1ZnZwWY2xcwWmNl8M7sxqL/dzNaY2aygnJ7otqYyM1thZnODezEjqGtnZu+b2TfBY9tEtzMVmdmhZT4ns8wsz8xu0mcosczsX2a20czmlamr8DNj3gPB76U5ZjYscS1PHVHu0V/NbGFwH142szZBfQ8z21Xm8/Rw4lqeGqLcn6g/18zsd8FnaJGZjUpMq1NLlHv0fJn7s8LMZgX1+gzVs0r+xo7Z76IGOWTUzNKAxcApwGpgOnCRc+7rhDYshZlZF6CLc26mmbUEcoGzgfOBHc65uxPaQAF8IASGO+c2l6n7H2CLc25C8I8rbZ1zv01UG6XkZ9waYCRwBfoMJYyZHQvsAJ50zg0M6ir8zAR/1N4AnI6/d/c750Ymqu2pIso9OhX40DlXYGb/DRDcox7AG8XHSfxFuT+3U8HPNTPrDzwLjAAygQ+Avs65wnptdIqp6B7t8/2/Adudc3/SZ6j+VfI39uXE6HdRQ+0hHAEscc4tc87tBZ4DxiS4TSnNObfOOTczeP49sADomthWSTWNAZ4Inj+B/yEjiXUSsNQ5tzLRDUl1zrlPgC37VEf7zIzB/0HlnHMRoE3wi1ziqKJ75Jx7zzlXEHwZAbrVe8MEiPoZimYM8Jxzbo9zbjmwBP83n8RRZffIzAz/j/vP1mujpEQlf2PH7HdRQw2EXYFvy3y9GoWPpBH869FQ4Iug6vqgy/pfGo6YcA54z8xyzezqoK6Tc24d+B86QMeEtU6KXUj5X776DCWXaJ8Z/W5KTj8H3i7zdU8z+8rMPjazYxLVKKnw55o+Q8nnGGCDc+6bMnX6DCXIPn9jx+x3UUMNhFZBXcMb+3oAMrMWwCTgJudcHvAQ0AsYAqwD/pbA5gkc7ZwbBowGrguGiUgSMbNGwFnAi0GVPkMNh343JRkzuxUoAHKCqnVAlnNuKPAfwDNm1ipR7Uth0X6u6TOUfC6i/D9Q6jOUIBX8jR310ArqKv0cNdRAuBo4uMzX3YC1CWqLBMwsA/8/ao5z7iUA59wG51yhc64IeAwN/Ugo59za4HEj8DL+fmwoHkoQPG5MXAsFH9ZnOuc2gD5DSSraZ0a/m5KImf0MOBMY64IFE4KhiN8Fz3OBpUDfxLUyNVXyc02foSRiZunAT4Dni+v0GUqMiv7GJoa/ixpqIJwO9DGznsG/pl8IvJbgNqW0YIz5P4EFzrl7ytSXHbN8DjBv39dK/TCz5sFkZMysOXAq/n68BvwsOOxnwKuJaaEEyv1rrD5DSSnaZ+Y14LJghbds/CIM6xLRwFRnZqcBvwXOcs7tLFPfIVi0CTM7BOgDLEtMK1NXJT/XXgMuNLPGZtYTf3++rO/2SYmTgYXOudXFFfoM1b9of2MTw99F6TFuc70IVg27HngXSAP+5Zybn+BmpbqjgUuBucVLEwO/By4ysyH4ruoVwC8T0zwBOgEv+58rpAPPOOfeMbPpwAtmdiWwCjgvgW1MaWbWDL96ctnPyf/oM5Q4ZvYscDxwkJmtBv4ITKDiz8xb+FXdlgA78SvESpxFuUe/AxoD7wc/8yLOuWuAY4E/mVkBUAhc45yr7oInUgtR7s/xFf1cc87NN7MXgK/xQ32v0wqj8VfRPXLO/ZP957ODPkOJEO1v7Jj9LmqQ206IiIiIiIhI3TXUIaMiIiIiIiJSRwqEIiIiIiIiKUqBUEREREREJEUpEIqIiIiIiKQoBUIREREREZEUpUAoIiIiIiKSohQIRUSkSmZWaGazzGyemb1oZs3MrIeZzav61eXOc7mZZcaoTcPN7IFYnKuK93nLzNrU8rVnm1n/WLeptszseDPbbmZvlanrY2ZvmNlSM8s1sylmdmwl52huZt+ZWet96l8xs/PN7AIzW2Jmb8TzWkREJDYUCEVEpDp2OeeGOOcGAnuBa2p5nsuBmARC59wM59yvYnGuipgXcs6d7pzbVsvTnA0kTSAMfOqcOx3AzJoAbwKPOud6OefCwA3AIdFe7Jz7AXgPf20E52kN/Ah4wzn3PPCLOLZfRERiSIFQRERq6lOgd/A8zcweM7P5ZvaemTUFMLMhZhYxszlm9rKZtTWzc4HhQE7Q29jUzE4ys6/MbK6Z/cvMGgevX2Fmd5jZzOB7/fZtRNDb9Ubw/Pbg9R+Z2TIz+1VQ/99m9v/KvOZ2M7vZzFqY2eQy5x8TfL+HmS0ws38AM4GDg7YcFHz/laAXbb6ZXV3mvDvM7E4zmx1cdyczOwo4C/hrcL299mn/RDN7KOiRW2ZmxwXXsMDMJpY57iEzmxG85x1l6ieY2dfBf+O7g7rzgl7c2Wb2STXu5VhgmnPuteIK59w859zE4HzNgzZND+7TmOCwZ4ELy5znHOAd59zOaryniIgkEQVCERGpNjNLB0YDc4OqPsDfnXMDgG3AT4P6J4HfOucGB8f+0Tn3b2AGMNY5NwRwwETgAufcICAduLbM2212zg0DHgJuqUbz+gGjgBHAH80sA3gOuKDMMecDLwK7gXOC858A/M3MLDjmUOBJ59xQ59zKfd7j50Ev2nDgV2bWPqhvDkScc4cDnwBXOeemAq8Bvwl6V5dW0Oa2wInAr4HXgXuBAcAgMxsSHHOrc244MBg4zswGm1k7fAgbEPw3Hh8cexswKmjHWdX4bzYAH3yjuRX40Dl3BP6/01/NrDnwDhAuc/0X4kOiiIg0MAqEIiJSHU3NbBY+0K0C/hnUL3fOzQqe5wI9guGDbZxzHwf1TwAVzUk7NHj94ijHvVT2vNVo45vOuT3Ouc3ARqCTc+4roKOZZZrZ4cBW59wqwIC7zGwO8AHQFegUnGelcy4S5T1+ZWazgQhwMD4Qgx9GWzxnrrrtBXjdOefwoXmDc26uc64ImF/mHOeb2UzgK3yA6w/k4UPt42b2E6C4Z+5zYKKZXQWkVbMNJYLe3HlmVvzf/lRgXHDvPwKaAFnOub34sHtu0Hs6BD+MVEREGpj0RDdAREQahF1Br16JoENtT5mqQqBpDc5pVXy/+NyFVO/31b5tKX7Nv4Fzgc74HkPwQyU7AGHnXL6ZrcCHHYAfKmys2fHAycCRzrmdZvZRmdfkB8GuJu0t2+aifdpfBKSbWU987+gRzrmtwVDSJs65AjMbAZyE7527HjjROXeNmY0EzgBmmdkQ59x3lbz/fMqEcOfcOWY2HLi7+LKBnzrnFlXw2meB/wqOedU5l1/NaxYRkSSiHkIREYkp59x2YKuZHRNUXQoU9xZ+D7QMni/E9yj2ruC4WHoOH5rOxYdDgNbAxiAMngB0r8Z5WuN7GHcGcxqzq/GastdbG63wAXW7mXXCD9fFzFoArZ1zbwE34XvoMLNezrkvnHO3AZvxvZiVeQY42szKDi9tVub5u8ANxcNpzWxome9NwfeQXoeGi4qINFjqIRQRkXj4GfCwmTUDlgFXBPUTg/pdwJFB/YvB3MTpwMOxbohzbr6ZtQTWOOfWBdU5wOtmNgOYhQ+nVXkHuCYYZroIP2y0Ks8BjwWL3JwbZR5hZW2fbWZf4XvyluGHhIIPma+aXyXU8HMQwc/x6xPUTQZmV3H+XWZ2JnCPmd0HbMCH2OI5iX8G7gPmBKFwBXBm8NoiM5sEnIefNykiIg2QlY5wERERkQNZMOz1FufcmQfC+4iISN1pyKiIiEjq2AsMtDIb08eamV0A/APYGq/3EBGR2FEPoYiIiIiISIpSD6GIiIiIiEiKUiAUERERERFJUQqEIiIiIiIiKUqBUEREREREJEUpEIqIiIiIiKSo/w+qWCChdI3nBgAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 1080x432 with 1 Axes>"
      ]
     },
     "execution_count": 19,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "x_csfit = np.linspace(xmin, xmax, 1000)      # Create the x-axis for the plot of the fitted function\n",
    "\n",
    "if IncludeSignalInFit:\n",
    "    y_csfit = PDFmodel_CS(x_csfit, csfit_Nsig, csfit_mu, csfit_sigma, csfit_Nbkg, csfit_tau)\n",
    "else:\n",
    "    y_csfit = PDFmodel_CS(x_csfit, csfit_Nbkg, csfit_tau)\n",
    "\n",
    "ax.plot(x_csfit, y_csfit, '-', color='red', linewidth=2, label='ChiSquare fit')\n",
    "\n",
    "\n",
    "# Adding in some text (useful histogram and fit information) around the plot before showing it:\n",
    "if IncludeSignalInFit:\n",
    "    names = ['Entries', 'Nsig', 'mu', 'sigma', 'Nbkg', 'tau', 'Chi2/ndf', 'Prob']\n",
    "    values = [\"{:d}\".format(len(x_data)),\n",
    "              \"{:.3f} +/- {:.3f}\".format(minuit_cs.values['Nsig'], minuit_cs.errors['Nsig']),\n",
    "              \"{:.3f} +/- {:.3f}\".format(minuit_cs.values['mu'], minuit_cs.errors['mu']),\n",
    "              \"{:.3f} +/- {:.3f}\".format(minuit_cs.values['sigma'], minuit_cs.errors['sigma']),\n",
    "              \"{:.3f} +/- {:.3f}\".format(minuit_cs.values['Nbkg'], minuit_cs.errors['Nbkg']),\n",
    "              \"{:.3f} +/- {:.3f}\".format(minuit_cs.values['tau'], minuit_cs.errors['tau']),\n",
    "              \"{:.3f} / {:d}\".format(Chi2_value_cs, Ndof_cs),\n",
    "              \"{:.3f}\".format(ProbChi2_cs),\n",
    "              ]\n",
    "else:\n",
    "    names = ['Entries', 'Nbkg', 'tau', 'Chi2/ndf', 'Prob']\n",
    "    values = [\"{:d}\".format(len(x_data)),\n",
    "              \"{:.3f} +/- {:.3f}\".format(minuit_cs.values['Nbkg'], minuit_cs.errors['Nbkg']),\n",
    "              \"{:.3f} +/- {:.3f}\".format(minuit_cs.values['tau'], minuit_cs.errors['tau']),\n",
    "              \"{:.3f} / {:d}\".format(Chi2_value_cs, Ndof_cs),\n",
    "              \"{:.3f}\".format(ProbChi2_cs),\n",
    "              ]\n",
    "    \n",
    "d={}\n",
    "for n,v in zip(names,values):\n",
    "    d[n]=v\n",
    "\n",
    "ax.text(0.2, 0.95, nice_string_output(d, 0), family='monospace', \n",
    "        transform=ax.transAxes, fontsize=font_size, color='red', verticalalignment='top')\n",
    "\n",
    "fig       # This prints plots the figure"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Plot everything together: \n",
    "\n",
    "\n",
    "---"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAABDAAAAGoCAYAAACwmRWfAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvnQurowAAIABJREFUeJzs3Xd4VVXWx/HvJglNmg1EEAKIUqUkINXOWMAyimVABRQRFRsWQJ0BHQuOzIDdURQbKm1eG1ZsIziUBAELiIhBEZEmVUog+/1jnSQ34d7kJuSm/j7Pkyf39H3OuUTPOmuv7bz3iIiIiIiIiIiUZpVKugEiIiIiIiIiIvlRAENERERERERESj0FMERERERERESk1FMAQ0RERERERERKPQUwRERERERERKTUUwBDREREREREREo9BTBERKTccM495Zz7axHtq5FzbrtzLi6Y/tQ5N7go9h3s713n3ICi2l8Bjnuvc26Dc25tlOt759zRweeor69z7nnn3L0H0tZYcM71dM59V0LHLrLvp4iISEXkvPcl3QYREZF8OefSgHrAXmAf8C3wIvC09z6jEPsa7L2fVYBtPgVe9t5PLMixgm3HAEd77y8t6LZFyTl3FLAcaOy9XxflNh5o7r1fUcBjPQ+s9t7fVYh2fkohr3VF4JxLBH4EErz3e0u2NSIiIsVHGRgiIlKWnO29rwk0BsYCI4Bni/ogzrn4ot5nKdEY2Bht8EKKTmYmj4iIiBSeAhgiIlLmeO+3eO/fBC4GBjjn2kDObgvOucOcc2875zY75zY55z53zlVyzr0ENALeCrqI3O6cSwy6SlzpnPsJ+DhkXmgwo5lzbr5zbotz7g3n3CHBsU5yzq0ObaNzLs05d5pz7gzgDuDi4HiLg+VZXVKCdt3lnFvlnFvnnHvROVc7WJbZjgHOuZ+C7h93Rro2zrnawfbrg/3dFez/NOBD4MigHc9H2P4259yvzrk1zrkrci0Lvb4nOedWO+fuCNqU5pzrn2t3BzvnZjrntjnn5jnnmoXsq5tzbkFwLRc457oF8+8DegKPBe18LK/1Q67l351zc4JjfeCcOyzC+eW4V0G7b3XOLQn2PcU5VzVYttQ51ydk3fjgXDsG09Occ2uD7f7rnGud61o96Zx7xzm3Azg51/U7OPh+rnfO/R58bhjlOf03+L05uEZdw5znmKB9Lwfbf+WcO8Y5Nyr4jv3snPtTyPqDgvPd5pxb6Zy7OmRZ2H9LwbIRzrlfgu2+c86dGsyv5Jwb6Zz7wTm30Tk31WX/e6katGtjsM8Fzrl64e6XiIhIKAUwRESkzPLezwdWYw+8ud0SLDsc63pyh23iLwN+wrI5anjv/xGyzYlAS+D0CIe8HLgCOBLryvJIFG18D7gfmBIcr12Y1QYGPycDTYEawGO51ukBHAucCvzNOdcywiEfBWoH+zkxaPOgoLvMmcCaoB0Dc2/oLNhyK9ALaA6cls/pHQEcBjQABgBPO+eODVn+F+Bu4GBgBXBfcJxDgJnY9TsU+Bcw0zl3qPf+TuBzYFjQzmF5rR9yrH7AIKAuUDk4j2hdBJwBNAGOw+4FwKvBOWQ6HdjgvV8YTL+LXae6wEJgcq799gvOuSYwO9eySsAkLCumEbCT/e95pHM6IfhdJ7hG/4twXmcDL2HX/0vg/eC4DYB7gH+HrLsO6APUCo45PjNQQ4R/S8G9HgZ0CjKjTgfSgm1uAM7DvoNHAr8DjwfLBmDf0aOw+zk0OH8REZE8KYAhIiJl3RrgkDDz04H6WL2HdO/95z7/wk9jvPc7vPeRHqZe8t5/7b3fAfwVuMgVTdeA/sC/vPcrvffbgVHAJS5n9sfd3vud3vvFwGJgv0BI0JaLgVHe+23e+zTgn8BlUbbjImBSyDmOiWKbv3rvd3vvP8OCDBeFLPuP935+UKdhMtA+mN8b+N57/5L3fq/3/lVgGfbAHU4060/y3i8P7t3UkGNF4xHv/Rrv/SbgrZBtXwHOcc5VD6b7BfMA8N4/F1zn3di1aueCzJnAG977Od77DO/9rtADeu83eu9neO//8N5vwwIdJ+Zq14GcE8Dn3vv3g+s/DQtAjPXepwOvAYnOuTpBe2Z673/w5jPgA7IDg5H+Le0DqgCtnHMJ3vs07/0PwTZXA3d671eHXJ++wXc6HQtcHO293+e9T/Xeby3guYmISAWkAIaIiJR1DYBNYeY/hL31/yBIiR8Zxb5+LsDyVUACloFwoI4M9he673jsbXem0FFD/sCyNHI7DHtTn3tfDQrQjtznmJffg0BH6PpHhkxHanPu882vndGsH831iSTstkHh0qXA2UEQ4xyCAIZzLs45NzboIrGV7MyD0O9DxO+Tc666c+7fzrr5bMW6hdTJFRA7kHMC+C3k804se2RfyDSZ+3TOnemcmxt0EdkMnBVyLmH/LQXX5yYsOLHOOfeacy7z/jcG/i/oIrIZu477sO/0S1g2yGvOuir9wzmXUMBzExGRCkgBDBERKbOcc52wh9jc6fkEb8Zv8d43xd7UD8/snw9EysTIL0PjqJDPjbA3yRuAHUDmW/rMTIjDC7DfNdgDX+i+95LzATQaG4I25d7XL1Fu/yv7n2NeDnbOHZRr/TVRHCf3+WZum9nO3Ncrv/VjKbMbybnAtyGjsfQL5p2GdYdIDOa7kG3zuu+3YF2Cjvfe1yK7W4iLvElU+y0w51wVYAYwDqjnva8DvJPZlrz+LXnvX/He98DujwceDHb7M3Cm975OyE9V7/0vQRbH3d77VkA3rOvK5UV5TiIiUj4pgCEiImWOc65WUFzxNWy4za/CrNPHOXe0c84BW7G3v5lvn3/DakQU1KXOuVbB2/h7gOnBG+3lQFXnXO/gTfJdWGp9pt+wdP1I/919FbjZOdfEOVeD7JoZBRoiM2jLVOA+51xN51xjYDjwcpS7mAoMDDnH0VFsc7dzrrJzrif2IDotim3eAY5xzvVzVhjzYqAV8HawPPf9yW/9WHoN+BNwDSHdR7C6FruBjVjw6v4C7rcmlgWxOajxEc21zrQeyKBw3+FwKmPf1/XAXufcmdg5A5H/LTnnjnXOnRIEQHZh55P5b+wp7HvYONjH4c65c4PPJzvn2gaBvq1Y0C1zOxERkYgUwBARkbLkLefcNuzt7p1YMcdBEdZtDswCtgP/A57w3n8aLHsAuCtIby9IsceXgOex1P6qWKFCvPdbgGuBiVhWwA6s6GGmzIf6jc65hezvuWDf/wV+xB4Gry9Au0JdHxx/JZaZ8kqw/3x5798FJgAfY10GPs5nk7VYccY1WI2Lod77ZVEcZyMW7LgFCwDcDvTx3m8IVnkYq5fwu3PukSjWjxnv/a/Y96cbMCVk0YtYN5ZfgG+BuQXc9QSgGpY1Mxd4rwBt+gOrmTEn+A53KeCxc+9vG/Zdnordz37AmyGrRPq3VAUbzngD9l2oixX4BLuHb2LdTrZh53h8sOwIYDoWvFgKfEb0QTYREanAXP71zERERERycs6dhGW/NMxvXREREZGioAwMERERERERESn1FMAQERERERERkVJPXUhEREREREREpNRTBoaIiIiIiIiIlHrxJd2AA3HYYYf5xMTEkm6GiIiIiIiIiOQjNTV1g/f+8MJuX6YDGImJiaSkpJR0M0REREREREQkH865VQeyvbqQiIiIiIiIiEippwCGiIiIiIiIiJR6CmCIiIiIiIiISKlXpmtgiIiIiIiIiKSnp7N69Wp27dpV0k0RoGrVqjRs2JCEhIQi3a8CGCIiIiIiIlKmrV69mpo1a5KYmIhzrqSbU6F579m4cSOrV6+mSZMmRbpvdSERERERERGRMm3Xrl0ceuihCl6UAs45Dj300JhkwyiAISIiIiIiImWeghelR6zuhQIYIiIiIiIiIlLqKYAhIiIiIiIicoDWrl3LJZdcQrNmzWjVqhVnnXUWy5cv59NPP6VPnz5htxk8eDDffvstAM899xxt27bluOOOo02bNrzxxhvF2fyIdu/ezWmnnUb79u2ZMmVKjjbff//9xdoWFfEUEREREREROQDee/785z8zYMAAXnvtNQAWLVrEb7/9lud2EydOBKwI6X333cfChQupXbs227dvZ/369TFt8759+4iLi8t3vS+//JL09HQWLVoEwMUXX5y17P777+eOO+6IWRtzUwaGiIiIiIiIyAH45JNPSEhIYOjQoVnz2rdvT8+ePQHYvn07ffv2pUWLFvTv3x/vPQAnnXQSKSkprFu3jpo1a1KjRg0AatSokTWCR2pqKu3ataNr167cdttttGnTBoDnn3+eYcOGZR2vT58+fPrppwBcc801JCcn07p1a0aPHp21TmJiIvfccw89evRg2rRp/PDDD5xxxhkkJSXRs2dPli1bluO81q1bx6WXXsqiRYto3749P/zwQ1abR44cyc6dO2nfvj39+/cv4isangIYIiIiIiIiUn44F5ufPHz99dckJSVFXP7ll18yYcIEvv32W1auXMmcOXNyLG/Xrh316tWjSZMmDBo0iLfeeitr2aBBg3jkkUf43//+F/UluO+++0hJSWHJkiV89tlnLFmyJGtZ1apVmT17NpdccglDhgzh0UcfJTU1lXHjxnHttdfm2E/dunWZOHEiPXv2ZNGiRTRr1ixr2dixY6lWrRqLFi1i8uTJUbftQCiAISIiIiIiIhJDnTt3pmHDhlSqVIn27duTlpaWY3lcXBzvvfce06dP55hjjuHmm29mzJgxbNmyhc2bN3PiiScCcNlll0V1vKlTp9KxY0c6dOjAN998k1WzArK7gGzfvp0vvviCCy+8kPbt23P11Vfz66+/Fs0Jx4hqYIiIiIiIiEj5EXTPKE6tW7dm+vTpEZdXqVIl63NcXBx79+7dbx3nHJ07d6Zz58706tWLQYMGcdNNN0UckjQ+Pp6MjIys6V27dgHw448/Mm7cOBYsWMDBBx/MwIEDs5YBHHTQQQBkZGRQp06drNoWZYEyMEREREREREQOwCmnnMLu3bt55plnsuYtWLCAzz77LKrt16xZw8KFC7OmFy1aROPGjalTpw61a9dm9uzZADm6aiQmJrJo0SIyMjL4+eefmT9/PgBbt27loIMOonbt2vz222+8++67YY9Zq1YtmjRpwrRp0wArRLp48eICnXdCQgLp6ekF2uZAlKsARuLImSSOnFlq9yciIiIiIiLlj3OO//u//+PDDz+kWbNmtG7dmjFjxnDkkUdGtX16ejq33norLVq0yBqu9OGHHwZg0qRJXHfddXTt2pVq1aplbdO9e3eaNGlC27ZtufXWW+nYsSNg9TQ6dOhA69atueKKK+jevXvE406ePJlnn32Wdu3a0bp16wIP3TpkyBCOO+64Yivi6XwJpNcUleTkZJ+SkpI1nRlsSBvbu0j2X9T7ExERERERkaK3dOlSWrZsWdLNiLm0tDT69OnD119/XdJNyVe4e+KcS/XeJxd2n+UqA0NEIoiLg/bts38+/zz/bSZMgD/+iP4Yb74JY8cWvo1SOv3yC/ToAa1aQceO8OGH2ctCv1c33JA9f+pUaN4cjjkGQipoF3i+iIiIiEgIFfEUqQiqVYOCFueZMAEuvRSqV49u/XPOsR8pX+Lj4fHHoV07WLUKunWzoAaE/17t2QMjRsDcubB7N5x8MvTuDXv3Fmx+JcXXRURERHJLTEwsE9kXsaL/QxSpyGrUgFGjoHVruOoqmzdrlr1RX7PGHiYzPwN8+in06gV9+0KbNnDTTTb/8suhUSMYNizn/t96Czp3tn0MH549/+GHoWVLOO44e3iV0qtePQteADRuDOnpFmiIZN48+z7Vq2ffiaOOgiVLCj5fRERERCQXZWCIVAQ7d1oQIdNbb9mD4o4dcPbZcM890KyZBSpOO83eqicmwiefwGGH5dzXnDkwf74FMDZvtnkvvgjPPw8hNWlYtw7+/ncLelSvDhdeCB9/DKecAqNHw88/Q82atp6UDe+/Dx06QOYwYLt2WbeSatXggQfghBNg7VqoXx+eegoOOQSOOAJ+/RW2by/Y/NDvq4iIiIgICmCIVAyRupBUrmxdAgCaNLGHz/wqJSclWfACoE6dyOvNnQsrV2bvf/t2mz7lFOjUCa64wroK/PnPBT8fKX5r18Itt0BoZerVqy3gsGCB3cfvv88ed33oUPv9n/+AcwWfLyIiIiKSiwIYIhVZQkL2Z+cgIyP/bWrXjn7/p58OIWNVZ3nvPZg9G157DR59FFJTo9+nFL9du6zb0LhxlqmT6Ygj7HenTpZFsWqVBcB+/TV7nbVrbb0aNQo2X0RERCTGNOpk2aMAhoiEV6sWbNq0fxeSaHXpYjUxVq+Ghg3t4bZKFat18NNPcOKJ0LatjT4hpZf3MGAA9OsHZ5yRPX/TJsvsqVYN0tKs+1GjRtC0KXzzjXUN2r3b7v9xx1mxzoLMFxERESlj4uLiaNu2Lenp6cTHxzNgwABuuukmKuVRnDwtLY0vvviCfv36FWNLyy4FMEQqgtw1MEaPzr/rxvXXw3nnWV2C6dMjvxVPS7P1Nm2y48yeDfffD2edBU8+aTU29u6Fgw6ybAzv4bLLYOtW2LcP/vnPIjtNiYE5c6xbx3ffwdNP27x33rH7PmiQBaXi4mDixOwRa8aOhe7d7fP48TaiSOXKBZsvIiIiUsZUq1aNRUG37XXr1tGvXz+2bNnC3XffHXGbtLQ0XnnlFQUwouR8Zv/jMii5TRufEjKETFGnAEW7P6UeiYiIiIiIlJylS5fSsmXLAm2T+Rw345puJDU++IDbUKNGDbZv3541vXLlSjp16sSGDRtYtWoVl112GTt27ADgscceo1u3bnTp0oWlS5fSpEkTBgwYwJ///Oew65VF4e6Jcy7Ve59c2H2W7ddcS5fCu+/uNzt11e9Fepii3p+IiIiIiIiUnNBnvP4T58bkma9p06ZkZGSwbt066taty4cffsjChQuZMmUKN9xwAwBjx46lZ8+eLFq0iJtvvjniemLKdgAjIwP69IF//YvUtE1Zs4viC1gcX2gREREREREpfnNXbsz6nL43I8d0Ucrs8ZCens5VV11F27ZtufDCC/n222/Drh/tehVV2a6BUb++Va+/5Rbm/rAHarQB57K+gAeSBhTuC10UaUUiIiIiIiJSsro0PTTrc0J8pRzTRWXlypXExcVRt25d7r77burVq8fixYvJyMigatWqYbcZP358VOtVVGU7A+PII2HaNKhenS6vv0DlfengfZF8AYvjCy1SLJ5/3golbthg0336wKef2ucaNcJvk5iYvX5J++UX6NEDWrWCjh3hww9zLt+2zf4WjBuXPe+hh6BNG9smtGhSuPnvv28FTjN/KleGoPhSgdqzcSMkJ0O7draf11+Prv2FkZEBrVsX/h7dequNBtOmTfa8/K5DuOscybx5NpJIy5Zw0UXZ8yPdlzvvtPkdO2ZfNxEREZEYCn05PXlwlyJ/Wb1+/XqGDh3KsGHDcM6xZcsW6tevT6VKlXjppZfYt28fADVr1mTbtm1Z20VaT0zZzsAA6NsXjj6apHPO4dVX7+CD5l34051DD/gLGOsvtEixql4dXn3VRhYpa+Lj4fHHLTCwahV062ZBgUz33gtJSdnTq1fDv/9tNXK8hxYt4PLLISEh/PzTT7cfsIyuE07IOWJLtO2pVcsCQzVqWGChbVs455z82x/J88/bSB9jxuy/7J13oFOnwg9xe/75cMklMHBg9rz8rkPu6xxJ5igzzz1ngZv1621+pPuycSN88IEFSzZvhg4d4NRToWbNwp2biIiISAEV1bPezp07ad++fdYwqpdddhnDhw8H4Nprr+WCCy5g2rRpnHzyyRx00EEAHHfcccTHx9OuXTsGDhwYcT0xZT+AAfY/2QsWQMeTGPXZ87BgKrzwggU3ioCCF1LmnX22ZSvlDmB4D1ddBf/7nz00PvxwzuU//AAXXghTpkDz5vYG/dln4Zhj4Jtv4KOPLFsjlurVsx+Axo0hPR1277bhO5cvtwfk3A/We/faOt5bJkHt2vDHH+Hnh3rttfz/buTVnoQEm79li83buzfv9Qvr0UdtqNrC6tbNgiOR5L4Oka5zOCkpcPjhFrwA+5wp3PWfP9/+hsfHW0CmQQP7e37KKYU6NREREZGSkle2RPPmzVmyZEnW9AMPPABAQkICH330UY51w60npmx3IQlVrx5/+csDTGtzmj2oXHihpShnZJR0y0RKXu3alv6/dGnO+X/8YQ+qS5ZAaip8/nn2stWrLf3/hRcsePHTT/YGPTUV/vEP+PHH4j0HsG4OHTpkP/yPGAGjR+dcp2FDuPFGaNTIPt9yCxxySOT5oSZPtsyEwrZn2zbLvDjuOHjiCXtIz2v9wli2DLZvjy6YUFi5r0O46xzJTz9B3bpw5pl2rk88YfMjXf9WrazLyR9/2LZLl8LatUV/TiIiIiK5pI3tTdrY3iXdDCmA8pGBEdgTn8BtZ93IhYPOhNtus9Trr7+2VGyl3khFd/nlFowIFRcHvXpBpUr2OyUFeva0ZeeeC7172wM5WOCiZ0/7t9SihWUTFKe1a+2h9403bPqttywTJHc7Nm+24ZV//NGyHbp3t7of1aqFn1+/vm23fDns2GFdPQrTHrBuD199ZQ/hQ4bABRdkZ2WEWz+c8ePtPm3aBHv2WE2IXr0s+wXgscdg2LDo2lgYua9DpOscya5d8MknFhSrU8fqgpxxhgUrwl3/tm2tK0u3bpZ9cfLJoGJVIiIiIhJGuQpgAOAcDB9uxeMuuQSmT4cVK+yhoVGjkm6dSMk5/XR7ix6a0p+bc9mfn37aAoHz50Pnzpb2X1J27bJMkXHjoFkzmzdvHsyYYf+2N2ywIEz9+pbd0KhRdveQDh3gyy/tDX+4+ZkBjIJkX4RrT6iWLS1wsXixPcDnt36om2+2n3A1MLZutXoR48fn3GbePLj6avv8xBMWDCis3Nch0nXu3z/89vXq2fln/r1NSrKskbyu//Dh9gPQpYv+VouIiEiheO9xof8/KyXGx+jZofx0IcntzDNh7lw4+mgrDtepE3zxRUm3SqTkxMVZBkXov4N9+6yORUYGzJpl/04yJSVZl5ErrrC6BcnJMHu2PYguW2YFKYuD9zBgAPTrZ2/yM917rwUnly2zjITbb7eH6iOOsBoKu3fDzp2wcCE0aRJ5fqZXXgkfwBg1yn7ya88vv1hBSrBsi2+/tYyCSOsXxqRJdo6ZWR2Zjj/e/s4tWnRgwQvY/zpEus6Zcl+fTp2sK0hmBslXX0HTpnlf/8zr9tlntl0su8eIiIhIuVS1alU2btwYswdniZ73no0bN8ZkCNiYZmA45+oAE4E2gAeuAL4DpgCJQBpwkff+d2ehsoeBs4A/gIHe+4UH1ICWLe3t4cUX28PZySfbaACDBx/QbkXKrAEDcr69r17dCjbefLMV8ezePef6Xbta94W//tXqXgwZYkNdtm9vmQQHUsshWnPmwH/+A999Z1khYKNwHHlk+PV79LBsk3btLFtg8GD7WwCR58+bZ6OHHHvs/vv79dfo2vPTT3Z9wAJDDz5o2QWzZxes/ZlCRwgBC4Q884z9LTtQ111nbdqwwepRPPGEjZiS13WIJPf1qV0bJkywIpzp6Ra4adHClkW6/oMGWZAkIQFefjlnJpCIiIhIFBo2bMjq1atZnzkCmpSoqlWr0rBhwyLfr4tlhMo59wLwufd+onOuMlAduAPY5L0f65wbCRzsvR/hnDsLuB4LYBwPPOy9Pz6v/ScnJ/uUlJSs6cSRMwH2L8Syd6+lJz/6qE0PHWqjLeQusJdLxP1FWG/GNd00YomUb9u2WZ2HDRvsLXlamh42i8u778JLL1mGhIiIiIhIGeScS/XeJxd2+5h1IXHO1QJOAJ4F8N7v8d5vBs4FMisJvgCcF3w+F3jRm7lAHedc/SJpTHw8PPKIpV9XqQJPPWVvB4ug0n3qqt+zPvefODfHtEi5c/PNNsLGKafYvykFL4rPmWcqeCEiIiIiFVosa2A0BdYDk5xzXzrnJjrnDgLqee9/BQh+1w3WbwD8HLL96mBeDs65Ic65FOdcSoHTgwYOtGEiGzSwNPCkJEuZPgBzV27M+py+NyPHtEi5M3GijS6xZImNUiIiIiIiIlJMYhnAiAc6Ak967zsAO4CReawf7lXufv1bvPdPe++TvffJh+c1mkIknTplDwe5Zg2ccAI8+2zB9xPo0vTQrM8J8ZVyTIuIiIiIiIhI0YhlAGM1sNp7n5niMB0LaPyW2TUk+L0uZP2jQrZvCKyJRcMSx6dw9PG3WEX9PXusmNy119rnAgqteTF5cBfVwJDSa8IEG0FEovP++1asNPOncmUb5QPgzjuhTRsraPr66/nvK9L6U6dC8+ZwzDHw1lv5r19YGRnQurXVLimMSO2M5KGHrP2tWsHdd+dctm2bFTAdNy573q232vCrbdoUrn0iIiIiUiHEbBQS7/1a59zPzrljvfffAacC3wY/A4Cxwe83gk3eBIY5517DinhuyexqEq38im2G2hsXb0U9O3a0op5PPmnD/U2bZsP9FXB/gIIXUrpNmACXXmojj0j+Tj/dfsBG2jjhBAtkpKTABx9YMGPzZujQwUZwqVkz/H4irV+lCowYYcM9795toyT17m3DixZk/5mef96Kqo4Zs/+yd96x7LPDDiv4ddizJ3w7K0WIf69ebcPvLl1qI6e0aAGXX549ZOq99+4/TOr559vQrblHXhERERERCRHLDAywUUUmO+eWAO2B+7HARS/n3PdAr2Aa4B1gJbACeAa4NsZtM4MGZdfFmD0bkpMPuC6GSKkya5Y9eK9ZYw+fmZ/Bhs5s184eKB97LHubGjWyP590kj2EV2SvvQZ9+9rnH36waxgfbwGBBg1gwYLI20Zaf948y4qoVw8aNYKjjrLaIgXdfzQefRSuv75w20ZqZ1727rVgx+7dlrlSu7bNX74c1q/fP4DRrVvhgisiIiIiUqHELAMDwHu/CAg3RMqpYdb1wHWxbE9EnTtbXYy+fS2IccIJlpFxxRUl0hyRInXaafY2PzERPvkk54PiE09Aw4aQnm4PqX37ZmUgSYjJk7Nr5bRqBffdZ91xNmywTIN0fZgEAAAgAElEQVS8RjSKtH5CAtSvb6MiHXKIXfdffy34/vOzbBls375/0CBaa9eGb2f79uHXb9gQbrzRgh379llXkUMOsWUjRlgm0KRJhWuLiIiIiFRoMQ1glCn16sFHH9kwkU88AVdeaUGN8ePtDaJIefTss1ZjwXvLylizRgGM3JYvhx07LFMFoG1b6+rQrZtlR5x8MlStGnn7SOvv3WvLhw613//5jw1LW9D9jx8PL7wAmzZZd4/XX4devawOBVhmzbBhhT9/78O3M5LNm+Hdd+HHHy0w1r079OljWTzHHAONGxe+LSIiIiJSocW6C0nZUrkyPP64PdRVrmyBjJNOsj7dIuXNp59aoco5cyxDo0ULK/YIOR9QMx+0K6rJk60+Q6jhw+2azZxpGQqNGuW9j3DrH3mkZTJkWrs2O3hUkP3ffLOte889FmRYtCg7eLF1q9XTyOz+kmnevOzipF98kXfb82pnOLNmWXtr17Zsnw4d4Msv7ZgzZtj37LHH4B//sGsrIiIiIhIlZWCEc8UVVg2/b1/43/+s0OeUKfYmVKSsqlXL3tJndiHZutU+V68O33wDixdnr1u7tq1brZp1QajIXnkF3n4757yNG+HQQ+Gzz+w6hXbPGDXKfj/wQN7rp6fbdV+3zmpFrF4Nxx2X//4LYtIk6N/fuquEOv747BFV8tO5c+R2hjvfI46wmh27d1tAbOFCGD0azjrLCniCFRqtUcPaJiIiIiISJQUwIsmsi9Gvn71RPO00uP9+uP32vNOnRUqr66+H886zegTTp8MZZ8DTT1vXiBYtLFCXacQIW56cbEUbK6p58+xB+9hjc84fNAhWrLDAwMsv5/yb8GuYwZPCrV+5Mowda10swLqCZI7skdf+I8k9gof38Mwz9vfrQOTVTtj/fHv0sNFb2rWz9QYPhpYt8z7GdddZ15QNG6yGxhNPWIFZEREREZEQzmf2by6DkpOTfUohRkdIHDkTiHKY1H377O3hfffZ9Hnn2XCFmVX1C7o/EZHi8O678NJLlkEiIiIiIlIKOOdSvffhBvqIimpg5CcuztKe33zTghavvw6dOsFXX5V0y0REIjvzTAUvRERERKRcUQAjWmefbV1K2rWD77+3PuQqQCciIiIiIiJSLCp0ACN11e8F26BZM6vYP2AA7NwJl14Kw4aRdk8vdR+R0m35cujZE1q1slEhfvrJ5teoEXmbN9+02geZfvnF6hu0amX1Mj78cP9tMjKgdWurZRCtNm0gLS17+vPP7Rjt28O330a/n/Jk3jwrlNmyJVx0Uf7zw9m40WqYtGtn1/L113Mu37bNRhgZN86mo7m/kUR73++80+53x4452/PQQza/VSu4++7ojysiIiIiFUqFK+IZGrToP3Eukwd3IanxwdHvoHp1q+zftSvccIMNu5qaCtOmWfE5kdLossvgn/+0B9Q1a6BKlfy3OeecnIUU4+Pt+96uHaxaBd262UNvqHfesS5WmSOdFMbkyXDbbVbIsiLy3u7Xc8/Z/Vq/Pu/5kdSqZUPl1qhhgYW2be1+ZhbgvPfenKObRHN/I4nmvqek2JCuixbB5s0WSDv1VNiyBf79b1i61M6xRQu4/HJo0iS6Y4uIiIhIhVHhMjDmrtyY9Tl9b0aO6ag5B1dfDbNnQ6NGMHeuvVH8+OMibKlIEVm82AIWPXrY9JFH2hCdmUaNsrfnV12VPe/yy+27PWxY9rx69ezhFqBxYxsGdPfunMd69FEb7SRTjRrh9/+Pf1gWwYUXwq5dNm/zZssUmDoV7rmn4mZgpKTA4Ydn36/DD897fiQJCdkZNlu22L3au9emly+3AEhoACOa+xtJ7vsezg8/2D2Nj7dAR4MGNtwqWLt277afypVzFEkWEREREclU4QIYXZpmP7glxFfKMV1gnTpZ9sWf/mQPA7162WglGRlF0FKRIrJypXV/CmfHDqvvsmgRvP++ZWcAvPiiBREief99e4MemsmxbBls357zoTjc/tPSbHjPhQthzBhrH0CdOrbeOedYl4JFi6xLQUXz009Qt64V4ezQwYYUzWt+XrZts8yL446z9StXtvkjRtjoSpGEu7+RhLvv4bRqZV1g/vjDzmXpUli71jLXbrzRAmYNG8Itt9hQvyIiIiIiuVS4LiSh3UUK3H0knMMOs/Tpu++Gv/8d7rrL+vC/9FL+b0hFioP3ljUUTuXK1lUALGV/7VrL0MjL2rX2kPnGGznnP/ZYzoyNSPtPS7MsgmrVLDOjceMCn1K5tmsXfPIJLFliQZ3kZDjjjMjzmzaNvK+aNW3EpKVLYcgQuOACeO89OOaYyNc90v2NJNx9D6dtWxg40L4PDRrAySdD1aqWefPuu/Djj5b10b079OkD9etHd3wRERERqTAqXAZGqAMOXmSKi7O31e+9ZwGNzLeXs2cXzf5FDkTTprBiRfhlCQnZn53LP3to1y7o29cKP4ZmdWzdavUN+vbNf//O5QyoeB/deVQU9epZ95pGjayORVKSZTlEmh+Nli3tXixebFkQM2ZYrYnHHrPuPJkjKkW6v5FEuu+RDB9umTUzZ1qgpFEjmDXLfteubX8/O3SAL7+Mbn8iIiIiUqFU6ABGkTv9dPsf7x49rPjdSSfBgw+qS4mUrPbtrStHZkDt119thIqC8t5G4OnXz978h5o0Cfr3zxmwiKRjR5gzxx6Wv/02e0QUMZ062TXZtAn27LEMiqZNI8/PNGqU/WT65Zfs+7x2rV3rBg2seOeKFRb8GDYMbr/d7l1e9zeSgtx3yG7PZ5/ZeSQlwRFHWC2M3bttdKeFC1XAU0RERETCqnBdSGKuYUMr5vnXv1rwYuRI+O9/rabAoQdQb0PkQLz0khXR3LTJHjbffjvy9zEtDc47z9bdudMCH/ffb2/9//Mf+O47ePppW/eddyzV/5ln7E16NBo3hiuvtDftrVrpYTW32rVhwgQ45RTrUtGvn2VLQOT5YIGpUD/9ZN1GAPbts79HeXXLmDMn/P2N1KXI+4Ldd7CRZVassO/gyy9bJk6PHhb8bdfORkgZPNgyRkREREREcnG+DKdvJycn+5SUlAJvlzhyJgBpY3sXdZNymjnTRnPYtAmOOgqmTLHhV0XKk3fftQDJK6+UdEukOOm+i4iIiEgBOedSvffJhd1eXUhiqXdv61LStSv8/DPpPXpy7ylXqs+/lC9nnqmH2IpI911EREREipkCGLHWqJH1977lFhIy9nHXJ89Zev7vv5d0y0RERERERETKDAUwikNCAowbx1Xn38WWKgfBm29a///580u6ZSIiIiIiIiJlggIYxejD5l3oPegRG01g1SorXjd+vLqUSOzFxdloJK1bwyWXWHHOwqhRo2jbJWXD1KnQvDkccwy89Vb+6995J7RpYyPOvP56/vMBtm2zgqHjxhVt20VERESk3FAAo5itrl3PRnW44QYbSWD4cOjTB9avL+mmSXlWrRosWgRff23TTz1Vsu2RsmPPHhgxwv5uzZoFN92U99DQKSnwwQf2ffvgA7j+egtORJqf6d57bVhVEREREZEIFMAoCZUrw8MP2xvIQw6xoQrbtYNPPinplkl55xyceKINZQnw6afQqxf07Wtvxm+6yeY/8ohNt2kDkyZlb++9Dcfapg3ceGOxN19KwLx5lrlTr57V9DnqKFiyJPL6P/xg2T7x8XDYYdCgASxYEHk+wPLlFsRVAENERERE8lAhAxhpY3vHfgjVaJx7rr2N7NkTfv0VTj0V7roL9u4t8kMljpyZNXysVGDp6RYwa9s2e96cOTBmjGVnjBkDaWkWwJg7Fz7/HEaPhnXrbN0//rBgx5IlkJpqy6V8W7sW6te3rJ2pU+GII+zvVSStWlnQ448/4KefYOlS20ek+WAZHqNHF8/5iIiIiEiZVSEDGKXKUUfBxx/b/7w7B/fdZ2/IV60q6ZZJebJzp7397tQJmjWDK6/MXpaUZBkVAHXq2NC/PXtavYuDD4bOnbPfuMfFWcZGpUr2OyWl+M9FildmjZ6hQ+Gii+yzc5HXb9sWBg6Ebt3gmmvg5JOhatXI8996y2prNG4c6zMRERERkTIuvqQbIFhK9Zgx9j/0/fvDF1/Yw+azz8L555d066Q8yKyBEU7t2jmnC1JUNq8HWSkfjjwyZ8bF2rWWhZGX4cPtB6BLF+t6Emn+66/DjBnwxhuwYYMFx+rXt7+FIiIiIiIhlIFRmpx4oj1knn02bN4MF1xgbyoLO2KESGF07GgFG3fssO/hggVw3HG2bN8++OgjK+I4a5ZldEj51rkzfPONdSP6+WdYvTr7+wAwapT9hNq40X5/9hls2pRd2yLc/HvvtZosy5bBsGFw++0KXoiIiIhIWMrAKG0OOwzeeIMxvYYy6tPnqPLUU/Yw+dprVkhPJNYSE22EiOOPt+nRo6FuXftcvbp9F2++2Wq2dO9eYs2UYlK5Mowdm32vx4+3LIlM4ephDBpkQYmEBHj55exMnUjzRURERESi4HxB0sVLmeTkZJ9ShvrgZxbRnHFNN5IaH5zvuq1/+4GZXzxuFfqrVbORSwYPLtT/9Gceu1QULxUREREREZEKxzmX6r1PLuz26kJSTFJX/Z71uf/EuTmmI/mmXjMb6WHgQOtGMmSIFdHbtCmGLRUREREREREpfRTAKCZzV27M+py+NyPHdJ5q1IBJkyzdumZNmD7d+p9/8kmMWioiIiIiIiJS+iiAUUy6ND0063NCfKUc01Hp398KfHbtCr/8YvUHbr8d9uwp4paKiIiIiIiIlD4KYBST0JoXkwd3ybcGRlhNm8J//2tDrjoHDz1kQxEuW1Zk7UwcOTOrXoaIiIiIiIhIaaEARgkoVPAiU3y8jQoxezY0aQJffmnDXj71FJThgqwiIiIiIiIieVEAo6zq2tW6lAwYYAU+r7kGzj0X1q/Pc7NoioeKiIiIiIiIlDYKYJRltWrB88/Da69BnTrw1lvQti28916O1QozAoqIiIiIiIhIaaIARnlw8cWweDGceCL89huceSbceCPs2gUcwAgoIiIiIiIiIqWEAhjlRaNG8NFH8MADVifjkUegUydYsuTAR0ARERERERERKWEKYJRyBeruERcHI0fC3LlwzDHw9dfQqRNJU57BZewDDmAEFBEREREREZESpABGKXTANSuSkmDhQrj6atizB0aMYMqrozhq81oFL0RERERERKRMUgCjFCpIzYrEkTNJHDlz/wUHHWRDq77zDtSvT+fV3/Lec8PgmWc03KqIiIiIiIiUOQpglEJFWrPizDPhq694u0VPDkrfBUOGwNlnw9q1RdBSERERERERkeKhAEYpFNrNo0hqVhx6KMPOHcENZ99mw63OnAlt2sD06QfYUhEREREREZHioQBGMUob25u0sb0LtE1R1qx4s9WJVtjzT3+CjRvhwgvhsstg8+b91i1w3Q0RERERERGRGFIAo6Jp0ADeew8efxyqVYOXX4a2bWHWrAMvHioiIiIiIiISIwpgVETOwbXXwuLF0KULrF4NvXox9+Hnswp85lc8VERERERERKQ4xTSA4ZxLc8595Zxb5JxLCeYd4pz70Dn3ffD74GC+c8494pxb4Zxb4pzrGMu2CdC8OXz+Odx7L8TH02XK01TZlw7eH3jx0EDEUVJERERERERECqA4MjBO9t63994nB9MjgY+8982Bj4JpgDOB5sHPEODJYmibxMfDnXfC/PkkHRzHK6/ewa3/fZHJOxeQVLdqxM0UmBAREREREZHiVBJdSM4FXgg+vwCcFzL/RW/mAnWcc/VLoH3lUr4FRDt0gNRU5jZqyzXzZpD0z9HQvj3873/F10gRERERERGRCGIdwPDAB865VOfckGBePe/9rwDB77rB/AbAzyHbrg7m5eCcG+KcS3HOpaxfvz6GTa+AqlThoRMH8OfLxkGrVrB8OfToAbfdBjt3xvTQyugQERERERGRvMQ6gNHde98R6x5ynXPuhDzWdWHm+f1meP+09z7Ze598+OGHF1U7y7yiHDFkSf1jIDUVRga9e8aNswwNZWOIiIiIiIhICYlpAMN7vyb4vQ74P6Az8Ftm15Dg97pg9dXAUSGbNwTWxLJ9ZV1Mhz2tWhUeeMCCFi1bwnffFVs2hoiIiIiIiEhuMQtgOOcOcs7VzPwM/An4GngTGBCsNgB4I/j8JnB5MBpJF2BLZleTiijfmhWQY5jTmA172rkzLFwII0bYdJCN0eGXZUV/LBEREREREZEI4mO473rA/znnMo/zivf+PefcAmCqc+5K4CfgwmD9d4CzgBXAH8CgGLatXAgd5rSohj0Nq2pVGDsW/vxnGDgQli1j+vLbeabTebDrVFsuIiIiIiIiEkMxC2B471cC7cLM3wicGma+B66LVXvKo6TGB2d9njy4S47pmDj+ePjySxg9Gh4ax9D5/7HaGJMmQZcusT22iIiIiIiIVGglMYyqxEDMgxeZqlaFBx/kgksfYsUhDWHZMujWDYYPhx07iqcNIiIiIiIiUuEogCGFsujIY+k96BG4/XZwDsaPh7ZtYdassOsXaYFRERERERERqXAUwJAcoikemml3fGVSrx0J8+dDu3bw44/QqxdccQX8/ntsR0kRERERERGRCkUBDCmw/QIThzWFBQvgvvugShWridGqFXPf/CxrvZiNkiIiIiIiIiIVggIYUmBhh29NSIA77oBFi6B7d1i7li5jR5GwNx28j3qUFGVpiIiIiIiISDgKYEiB5Tl8a4sW8N//wuOPk7R1Na+9OoobvniNyUdsIKlRnbD7U1cTERERERERyY8CGFJg+Q7fWqkSXHstfPMNW6sexPDZk0kadrnVx1i5cr/9hc3oEBEREREREQmhAIYckDyHb23UiEF9x3Bjn1vg0EPho4+gTRv4179g376s1fLM6BARERERERFBAQyJNed4o/XJsHQp/OUvsHMn3HILdO0KixcDUWR0iIiIiIiISIWnAEYZV5BhT0vU4YfDK6/A229Dw4Y2aklSEtx2G+zYkbWaghciIiIiIiISjgIYUrx694ZvvoEbboCMDBg3Dlq35qQfFpR0y0RERERERKQUUwBDil+tWvDwwzBvHnToAKtW8fz0u3ns9bGwZk1Jt05ERERERERKIQUwpOR06gTz58M//8kfCVXo891saNkSnngiR5FPEREREREREQUwpGTFx8Pw4fS68klmNesEW7fCdddB9+6wZElJt05ERERERERKCQUwpFT4pXZdBl/wN5gxA4480rqXdOwII0bkKPIpIiIiIiIiFZMCGBJTBRolxTk4/3wbcvX6663I5z/+AW3awLvvxrahIiIiIiIiUqopgCGFEtPhW2vVgkcegblzoX17SEuDs86Ciy9WkU8REREREZEKSgEMKb06d4YFC2yo1erVYepUtjc5mr+fMhjS00u6dSIiIiIiIlKMnPe+pNtQaMnJyT4lJaWkmyHFYdUquPFGeOMNm27bFh5/HHr2LNl2iYiIiIiISFScc6ne++TCbq8MDCkbGjeG119nUN/R/FS7Hnz1FZxwAgwYAL/9VtKtExERERERkRhTAEPKlE+adaLXlU/A3/4GVarAiy/CscdaNsa+fSXdPBEREREREYkRBTCkzNmdUAXuvhu+/hrOOAO2bIFhw6xmxrx5MTlm4siZJI6cGZN9i4iIiIiISP4UwJCy6+ij4Z13YMYMOOooWLgQunaFIUNof+OrCjiIiIiIiIiUIwpgSNnmHJx/PixdCiNGQFwcPPMMHz8zlIsXvw8ZGSXdQhERERERESkCCmBI+XDQQTB2LCxZAqecwiE7t/Lge49C9+6QmlrSrRMREREREZEDpACGlC8tW8KsWVx/9m38VuMQmDsXOnWCIUNg/fr9Vo9FbQvVyxARERERESl6CmBI+eMcb7U6kVMHPwW33prVrYTmzeGRRyA9vaRbKCIiIiIiIgWkAIaUSamrfs93ne1VqsNDD8FXX8Hpp9toJTfeCB06wMcfF0MrRUREREREpKgogCFlRmjQov/EuVEFMQBo0QLefRfefBOaNoVvvoFTT4W+fWm45bcYtVZERERERESKkgIYUmbMXbkx63P63owc0/lyDs4+24IX998P1avDjBnMmngNN85+BXbujEGLRUREREREpKgogCFlRpemh2Z9ToivlGM6alWrwqhR8N138Je/UHXvHm6e84oV/5wxA7zPc/Oosz5ERERERESkSCmAIWVGUuODsz5PHtwlx3QkEQMODRvCK69wYb+xfFu3CaxaBX37wmmnwddfR9xHgbquiIiIiIiISJFRAEPKpLyCFwUJOCw4qg19BkyAJ5+EQw6x4p7t28P118NG66JyQF1XREREREREpEgogCHlTkEDDhmV4mDoUPj+e7juOutG8thjcPTRMH48XY6qlbVuobuuiIiIiIiIyAFRAEPKncLUykhd9btlYDz2GCxeDL16webNMHw4SWf1oNPPX4P3UXddydqniIiIiIiIFAkFMKTcibZWRsSuJm3awPvvw9tvw7HHwvffM+2VkUyecidJm3/K89iqlyEiIiIiIhIbEQMYzrnk4myISCzklS2RZ1cT56B3b/jqK3jkETZXrUH3VUugQwe46ipYu7bg+8wlceRMEkfOLMDZiIiIiIiIVFx5ZWA845z73jl3j3OuVbG1SKSYRNXVJCEBrr+eE4c8w3NJ50BcHEycCM2bwwMPwK5dBd+niIiIiIiIFFjEAIb3vgPQB9gHTHfOLXLOjXDONS621onEUEGGZd1SrSb3nDbEhlg9+2zYvh3uuANatIApU6zwZwH3KSIiIiIiItHLswaG9/477/3d3vtWwACgDvCxc25OsbROpJhEHWg49lh480348ENo2xZWrYJLLoGePWH+/MLtU0RERERERPIVVRFP51wloC5QDzgIWB/LRolEkja2N2lje5d0M+C00+DLL+Hf/4bDD4c5c+D44+Hii2n0+68l3ToREREREZFyJ88AhnOup3PuCWA1cBswGzjWe39ecTROpFSLi4MhQ2DFChg5EqpWhalTmTXxGv4262nYGLmAp4iIiIiIiBRMXqOQ/AyMBZYCHbz3f/LeP+e931JsrRMpC2rVsoKey5fDgAHEZ+zjitQ3oVkzePBB2LmzpFsoIiIiIiJS5uWVgdHDe9/de/+o9/4359xBxdYqkbLoqKPg+efpPehhPmvSEbZsscyMY4+FF1+EfftKuoUiIiIiIiJlVl6jkKwCcM51dc59i2Vi4JxrF3QrESm1irpWRkH2t7RuUwZcdA988AG0awc//wwDBkBSkhX/LGGJI2eSOHJmSTdDRERERESkQKIp4jkBOB3YCOC9XwycEO0BnHNxzrkvnXNvB9NNnHPznHPfO+emOOcqB/OrBNMrguWJBT0ZkVKlVy9ITYUXXrDsjMWL4U9/gtNPt8+B1FW/l2AjRUREREREyoaoRiHx3v+ca1ZBcuFvJMjeCDwIjPfeNwd+B64M5l8J/O69PxoYH6wnUrbFxcHll8N338HYsVYv44MPSO19CQn70sF7+k+cm2cQQxkTIiIiIiIi0QUwfnbOdQO8c66yc+5WcgYkInLONQR6AxODaQecAkwPVnkByBzR5NxgmmD5qcH6IjETi2FZw+6zWjUYMQJ++AFuvJG5jduxz1UC50jfs5e536wu0jaIiIiIiIiUN9EEMIYC1wENsOFU2wfT0ZgA3A5kBNOHApu993uD6dXBfgl+/wwQLN8SrC9Sfhx2GEyYQJeH7iIuI4O4fXtJ2JtOl5sH2UgmO3aUdAsLTBkiIiIiIiJSHPINYHjvN3jv+3vv63nv63rvL/Xeb8xvO+dcH2Cd9z41dHa4Q0SxLHS/Q5xzKc65lPXr1+fXDJFSKalbG9LjEzh8x+9MXvF/JC1PhTvugKOPhieegD17SrqJIiIiIiIipUrEAIZzrrVz7pyQ6fHOueeCn45R7Ls7cI5zLg14Des6MgGo45yLD9ZpCKwJPq8GjgqOFQ/UBjbl3qn3/mnvfbL3Pvnwww+PohkipdfaWoeT9MZLNjpJcjKsXQvXXQctW8LkyZCRkf9OREREREREKoC8MjDGAhtCpk8HZgKfAH/Lb8fe+1He+4be+0TgEuBj733/YPu+wWoDgDeCz28G0wTLP/be75eBIVIunXYazJ8P06fDscfCypVw6aXQoQOnrJgP+qcgIiIiIiIVXF4BjPre+y9Cprd672d4718CDjuAY44AhjvnVmA1Lp4N5j8LHBrMHw6MPIBjiJQ9zsEFF8DXX8Ozz0LDhrBkCc/NuIdpk0fA558Xa3NU20JEREREREqTvAIYNUMnvPddQibrFuQg3vtPvfd9gs8rvfedvfdHe+8v9N7vDubvCqaPDpavLMgxRMqN+Hi44gr4/nv417/YVK0WnX75Fk44AXr3hkWLiuQweQ3dKiIiIiIiUtrkFcBY45w7PvdM51wXsutWiEisVK0KN9/MCVdPZEL3v0CNGvDOO9ChA/TrBytWFHiXoUGL/hPnKoghIiIiIiJlRl4BjBHAa8650c65s4OfMcCr2NCoIlIMtlepzoQe/Un9fDHceCNUrgyvvgotWsDgwbBqVdT7mrsyewCh9L0ZOaZFRERERERKs4gBDO/9fOB4IA4YGPxUAroEy0QkxnJkTMz4jtSbR8Py5TBwoBX2fPZZaN4crr2W4697Id+aFV2aHpr1OSG+Uo7pomyriIiIiIhIUcsrAwPv/Trv/d+89xcEP3/z3v9WXI0TKc/SxvYmbWzvPNcJmzHRuDFMmgRLl1pXkr174ckn+e+/r+Jvs562oVgjSGp8cNbnyYO75JguDHVJERERERGR4pJnAENESlaeGRPHHAOTJ8NXX0HfvlTZl84VqW9C06Zw++2wYUOYPWY70OAFqEuKiIiIiIgUHwUwREqxqDImWreGadM4c9AjfNC8C+zcCQ89BE2awF13we+xy4qIZU31Km4AACAASURBVJcUERERERGRUApgiJQR+WVMLK3blCHn3wULFsCZZ8L27XDffRbIuOce2Lq1UMfNq1tIUXdJERERERERiSRiAMM5F+ecu9o593fnXPdcy+6KfdNEpFCSk2241S++gNNOgy1bYPRoC2SMHUv1PTvz3UVhalsoeCEiIiIiIrGUVwbGv4ETgY3AI865f4UsOz+mrRKRA9e1K3z4IXz6KZxwAmzaBKNG8flTVzJ07nTL0IhAtS1ERERERKS0ySuA0dl73897PwEbTrWGc+4/zrkqgCue5onIATvxRAtifPghdOnCoTu3MvKz5yExER54IGzXkpKsbZE4cma+w8GKiIiIiEjFk1cAo3LmB+/9Xu/9EGAR8DFQI9YNE5Ei5Jx1J/niCy676B5SGrSEjRvhjjsskPH3v1tXk4BqW4iIiIiISGmTVwAjxTl3RugM7/09wCQgMZaNEpHCy7NehXO8NOWvJP/8DXz0kXUt+f13+NvfoHFjGDNmv1FLFLwQEREREZHSIGIAw3t/qff+vTDzJ3rvE2LbLBEpiAIX3XQOTjkFPvvMupecfLJlYNx9t2Vk/PWv1NlZuFFLREREREREYiG+pBsgInlLG9s733XCFd2MOnPixBPh44/h889tuNVZs+Dee5lduRovdOwDG46Hww4rbPNFRERERESKRF5dSESkjCiSops9e1qhzzlz4IwzqLFnJ9fNnWYZGbffDuvWFV2DRURERET+n707D4+6Ovs//r4nCXsA2fd9E1GRAQwirnVfcNcWqVqXqtjWR9vHpbbV9qfFWrVWq31UtGrjWrSiUBVFKSoBGRUBAQmYsC+yCIhAlvP740w2sk2Smcwk83ld17ky852Z75xhnOB8OOe+RWqo0gDDzLQ6Q6SBiGrRzaOOgv/8h3Mm3s97/UfBd9/Bfff5IOPGG2HNmjJ3z5l8RkSrREREREREROqiqhUYWWb2bzO71sz61NN8RKSOolV08/Nug7nygt/BggVw9tnw/ffw0EPQvz9cdRWsWBGV5xEREREREYlEVUU8RwK/CF/9i5l9YmYPmtnJZta0fqYnInEXDMLrr8PChfDDH0JBAUyZAkOGwCWX+OMxUG0hUhERERERSSpV1sBwzuU65/7unDsHOAp4A/gBMMfMptfHBEUkQRx2GDz/PCxfDldfDSkp8NJLMHw4nHkmfPxxlQ/vc+t0+txa9a+NGndTiZNIXouIiIiIiERXxEU8nXN5zrlZzrn/dc6NBq6J4bxEJM4qrW0xYAA8/jisWuVrYrRoAdOnw9ixcNxx8M474FytnrOibioiIiIiIiJQhy4kzrl10ZyIiDQwPXrAgw9CTg7ccQe0aQOzZ8Mpp8CoUfDqq1BYWKNTRqWbioiIiIiINEpqoyoiddOxI/zhD7B6NUyeDJ06QSgE558Pw4bBM89AXl5Ep4pqNxUREREREWlUqg0wzGxYfUxERBq41q3hllv8ioxHHoFevWDpUrj8cujfnys/+Tct9+2J+HTxCC9U20JEREREJHFFsgLj72Y238yuN7O2MZ+RiNRKpTUr6lvz5jBpEmRnwz/+AQcfDGvW8JtZT/LxY1fA7bfDxo3xnqWIiIiIiDQw1QYYzrmjgQlAT2CBmT1vZifFfGYi0rClpcFll8HixfDGG8zrcQht9n0Hf/wj9O4N11zjO5rUklZLiIiIiIgkl4hqYDjnVgB3ALcAxwJ/NbNlZnZeLCcnIo1AIABnnsnFE+7lvEvvg3PP9TUxnnjCr8447zzIyor3LEVEREREJMFFUgPjMDN7EFgKnACc5Zw7OHz5wRjPT0QakU+7H0zowSm+NsbVV0OTJvDaazBmDBxzDLzxBuZq1rkknkK52+M9BRERERGRpBHJCoxHgE+Bw51zk5xznwI459bjV2WIiFSp9Bf9CU9mEWrWCR5/3Bf8vP12aNsW5syBs8/mnSmTuPCLmbBvX/wmXIVyr0UhhoiIiIhIvYgkwDgdeN459z2AmQXMrAWAc+65WE5ORBqHrFVbiy/n5ReWXO/SBe6+27dgfeAB6NGDgVvXcN9/HoJ+/Xxb1m3b4jTrilX6WkREREREJKYiCTDeBZqXut4ifExEJCIZ/doXX05LDZS5DkB6OvzP/8CqVfDss3DoobB+Pdx2G/TsCTfcACtW1POsK1btaxERERERkZiIJMBo5pzbXXQlfLlF7KYkIo1NsPdBxZczr8ooc72MtDSYOBEWLoS33oJTToE9e+Bvf4PBg2H8eJg9G5yL6Xyr2hYS8WsREREREZGoiiTA+M7MRhRdMbMg8H3spiQijVlEX/jNfHjx1luwaBFceaUv+DltGhx3HASDnLPkfdIK8qI2r9rUtlB4ISIiIiJSfyIJMG4EXjGzOWY2B3gJuCG20xIRCRs2DJ580tfJuPNO6NgRPvuMiZ9O56r5rxH6w0NRqZOh2hYiIiIiIokttbo7OOc+MbMhwGDAgGXOuej9s6eISCQ6dYLf/Q5uuYXQlFf40detyAuk8PSOfDJHnUDw1LHwi1/AoEG1Or1qW4iIiIiIJLZIVmAAjAIOA44AfmhmP47dlEREqtCsGVlDx7AvJY3CQAp5qWlkdRoEjz4KQ4bA2WfD++/XuE6GaluIiIiIiCS2agMMM3sO+DNwND7IGAWMjPG8REQqldGvva+TAaQ1SSVj8q1w9dW+TsYbb8AJJ8Dw4X7ryZ49NT6/wgsRERERkcQTyQqMkcBY59z1zrmfhcfPYz0xEZHKlFstcewIePxxXyfjrrugc2f44gsfavTsCbfcArm5cZyxiIiIiIjUVSQBxmKgS6wnIiJSG2VWS3TqBL/9rQ8rnnsORo3yBT7/9Cfo1w/OOw8++CDmbVhjpc+t0+lz6/R4T0NEREREJC6qLeIJdAC+NLP5wL6ig865s2M2KxFpdHImn1F/T9a0KVx6qR/z5sHDD8PLL8Nrr/kxbBj87GcwYQK0bFnj09fraxERERERESCyFRh3AucA9wD3lxoiIonvyCPhn/8sacPapQssXgw//Sn06AG/+hV8/XW8ZykiIiIiItWoNsBwzs0GcoC08OVPgE9jPC8Rkago3nbRpYtvw5qbC5mZkJEBO3bAn/8M/fvDOedwVM7nDXZ7iYiIiIhIYxdJF5KrgX8B/xc+1B34dywnJSISM02awI9+BHPnwvz5MHEipKXB66/z/Et38M6USfDYY7BrV7xnKiIiIiIipUSyhWQSMBbYCeCcWwF0iuWkRETqxahR8OyzfnvJ738PXbsyaOtquP566NbN/1y0KN6zFBERERERIgsw9jnn9hddMbNUQGusRaTx6NwZfvMbv73kxRfhmGNg926/EuOww2DcOHj+edi3r/pziYiIiIhITEQSYMw2s9uB5mZ2EvAK8EZspyUiUrWcyWdEvxtIWhpcfDHMnu0LfU6aBOnp8OGHvmNJz55w222QkxPR6WLV9jSUuz3q5xQRERERSXSRBBi3AluARcBPgRnAHbGclIhI3B1yCDzyCKxfD3//Oxx+OGzZApMnQ79+cOaZMH06FBTUy3RKhxYTnsxSiCEiIiIiSSeSLiSFzrknnHMXOucuCF+udguJmTUzs/lmttDMlpjZXeHjfc1snpmtMLOXzKxJ+HjT8PXs8O196vriRETqrFUr33L1s8/g44/h0kv9So3p032IMWCADzU2b47pNLJWbS2+nJdfWOa6iIiIiEgyiKQLyddmturAEcG59wEnOOcOB4YDp5pZBnAv8KBzbiCwHbgyfP8rge3OuQHAg+H7iYhERZ1XLJjBmDHw3HOwdi3cey/07eu3k9x2m99eMmGC324Sg1asGf3aF19OSw2UuS4iIiIikgwi2UIyEhgVHuOAvwL/rO5BztsdvpoWHg44Ad+WFeAZ4Jzw5fHh64RvP9HMLIL5iYhUKGbbLjp2hP/9X8jOhhkz4KyzIC/PF/ocN85vP3nwQdp+vzOi00VSKyPY+6Diy5lXZZS5LiIiIiKSDCLZQrK11FjnnPsLPoSolpmlmNnnwGZgJrAS2OGcyw/fZS3QPXy5O7Am/Jz5wLdAuX9iNLNrzGyBmS3YsmVLJNMQkSQV820XgQCcdhpMmwZffw233w5dusDSpXDTTcz72495aNp98P77UV2VofBCRERERJJRJFtIRpQaI83sWiA9kpM75wqcc8OBHsBo4OCK7lb0VFXcVvqcjzvnRjrnRnbs2DGSaYhIkqrXbRe9e8Pdd8Pq1fDqq3DaaaQVFDB+6Ww44QQYPBjuuy/mtTJERERERBqrSLaQ3F9q/BEIAhfV5EmcczuAD4AMoK2ZpYZv6gGsD19eC/QECN/eBthWk+cRESktLtsu0tLg3HNhxgxGT3qGCRf9P0KHj4MVK/y2kx494KKLYOZMKCyM/XxERERERBqJSLaQHF9qnOScu9o5t7y6x5lZRzNrG77cHPgBsBR4H7ggfLfLgNfDl6eFrxO+fVYk3U5ERCJR39suQrnb+aZVOz7qO5wJZ91GKPMNXyujoABeeQVOPtl3MLnnHjruVlYrIiIiIlKd1OruYGY3VXW7c+6BSm7qCjxjZin4oORl59ybZvYl8KKZ/T/gM2BK+P5TgOfMLBu/8uKSCF+DiEjCKVd/o+sQgtOm+Q4mTz8NTz7p62b8+tfMtQDvDRgNxxiccgqkpMRx5iIiIiIiianaAIOSLiTTwtfPAv5LuOBmZZxzXwBHVHB8Fb4exoHH9wIXRjAfEZGEV2n9jR494De/8QU/Z86EJ57A/ft1TlmRBWecAd27w49/DFdcAQMHxmn2iaGoM0vO5DPiPBMRERERSQSR1MDoAIxwzt3snLsZXwOjh3PuLufcXbGdnohIw1Rt/Y2UFDj1VJg6laOu+wd/OubHfkvJunXwxz/CoEFwzDF+tcbu3cRTJG1eRURERERiLZIAoxewv9T1/UCfmMxGRKQRqq7+xpZWB/HomIsIzZwH//2vX33RsiXMmQM/+YlvzXrlleSc2YacP55eT7MWEREREUkskQQYzwHzzexOM/sdMA94NrbTEhFJDqHc7cWXJ0yZR6jXMHjqKdiwAaZMgaOPhu++88fGjfPtWCdPhvXrqziriIiIiEjjE0kXkruBK4DtwA7gCufcPbGemIhINORMPiOhayiUK/ZZdD093a++mDMHli+H226Dbt18O9bbboOePX3NjKlTYf/+MufUlg8RERERaYwiWYEB0ALY6Zx7CFhrZn1jOCcRkaRRabHP0gYNgnvugdxcmD4dLrjA19CYMcNf7tYNbrwRFi6sx5mLiIiIiNSvagOM8LaRW4DbwofSgH/GclIiIsmi2mKfpaWmwumnwyuv+C0kf/kLHHYYbN0KDz0Ew4fD8OFcOf81Ou7eXvl5REREREQaoEhWYJwLnA18B+CcWw+kx3JSIiLJqLpin2V06AC/+AV8/jmEQjBpErRrBwsX8pv3p5D16GVw2mnwwguwZ0/sJl0PStcJEREREZHkFUmAsd855wAHYGYtYzslEZHGoV7qb5jBiBHwyCN+Vcarr/L2wAwKLABvvQU/+lFxFxNmz4bCwlo/VX0GCWWKmz6ZpRBDRERERCIKMF42s/8D2prZ1cC7wBOxnZaIiNRY06Zw7rn89Lw7GH3Ds/C3v8GRR8KuXb6LyXHHQb9+cMcdvjBoBOIVJFRa3FREREREklYkXUj+DPwLmAoMBn7rnHs41hMTEZHa29G8NVx/PWRlwbJlPrTo3dsXAr37bhgyhM+6DeGOk6/3NTQqEa8gIaLipiIiIiKSVKoMMMwsxczedc7NdM79yjn3S+fczPqanIiI1F7xaonBg+EPf4BVq+CDD3x71vR0jtiwnP838zHo2hXOPRdeew327StzjngFCTUqbioiIiIiSaHKAMM5VwDsMbM29TQfERGpgyq3fAQCcOyxMGUKbNzIz8/6Fe/3C0JBAfz733Deeb5exlVXwfvvQ0FBjYKEPrdOp8+t06P+mhReiIiIiAhAagT32QssMrOZhDuRADjnfh6zWYmIJJFoFvqsaMtHhQFAixZMG3os04YeS84vRvhuJc8957uaTJniR7ducMklDNvai8Wd+ytIEBEREZG4iiTAmB4eIiKS4Gq15aNrV7jpJj+WLoXnn/dj1Sp44AHeBFa26w4trvZdTQYOrNMci1ZpxLxDi4iIiIg0KpVuITGzXgDOuWcqGvU3RRERiVSda0ccfLCvl5Gd7QuA/vznbGnRlv7b1sGdd8KgQTBqFDz4oG/bKiIiIiJST6pagfFvYASAmU11zp1fP1MSEZFoqNOWDzPfgvXII8loeiJjcr/gny1WwtSpsGCBHzffDMcfDxMm+PoZIiIiIiIxVFURTyt1uV+sJyIiIvFRptBnBQoCKXzY9wh4+mnYtAleecV3LUlLg1mz4MoroXNn/v7a3ZyxdA58912V5xMRERERqY2qVmC4Si6LiEgDd2C3kqq2m5SpVdG8OVxwgR87dsCrr/p6GbNmcepXczn1q7nQ6WE480y46CI4/XT/mFqIRY0M1d8QERERabiqWoFxuJntNLNdwGHhyzvNbJeZ7ayvCYqISPRV1K2kxtq2hZ/8BN59F9au5Yaz/pdfn3wdoba94OWXfcjRsaMv/Pn667B3bxRfgYiIiIgkm0pXYDjnUupzIiIiUn9q1a2kCqG85rw59BgApo46k8yWXxN87Rn45BPfovWFFyA9HcaP54Tv+/FhnyPq9HwiIiIiknyqWoEhIiKNVJ27lRygzIqOAkdWxikwf75vxTp5MowYAbt2wT//yaS5L3PtvH8RuuZm+M9/IC+vTs8tIiIiIsmhqhoYIiLSANW0vkNdwwuoYkVH375wyy1+rFhB6IXp/Ghnb/ICKTxekE/mVTcR3LvZFwW9+GLf1SRVfzWJiIiISHlagSEiInUW0YqOgQPJGncm+1KbUBhIIS+tCVkjT4Rt22DKFDj5ZOjaFa65Bt5+WyszRERERKQMBRgiIhJVVa3oKLNSIy2FjL/+ARYvht/+FgYPhm++gSeegFNPhU6d4LLLYNo0FQAVEREREQUYIiJSfypcqXHIIXDXXbB0KSxc6MOMQw7xbVqffRbGj/fdTC6+2Hc32b273Hn73Dq9uEVqJEq3kRURERGRhkEBhoiIxEW5lRpmcNhhPsxYvBiWLYO77/YFQHfv9uHFxRf7MOOcc+C553zIEaHSocWEJ7OiEmLUNDgRERERkdpTgCEikqRyJp9R44Kf9WrwYLj9dgiFfDeT+++Ho47y20lefx1+/GO/zeS007h44du02/Ntlacr0yklv7DMdRERERFJfAowREQk8fXtCzfdBB99BOvWwSOP+I4lBQXw1lvc+9bDfPLIRH/skUf8fQ5QaacUEREREWkQFGCIiEjD0q0bTJoEs2bBxo3wxBN80DdIgQXggw/gZz+DHj1g9Gi/BWXJEnAusk4pIiIiIpKwUuM9ARERaRzish2lY0e46iouz+5K6727+WJkHkyd6tuwfvKJH3fcAf37w/jxjNrYiVD3gxVeiIiIiDRACjBERKRexSro2NmsFaFxRxG89FLYswdmzvS1Mt54A1auhAce4BVga/PW8M35vrvJSSdBixYxmU8sFRUOTegaJiIiIiJRpi0kIiLSoFXYXaRFCx9QPPWU32by3//CzTeT07Yr7b/fCU8/7TuZdOjgfz79NGzZEpU5iIiIiEhsaAWGiIg0aBV1FymzRSQlBcaNg3HjOC7lOAZ+s5qZ/XfAv//tt5i8/rofgQCMHeuDj/HjYcCAKp/3wOBEdTVEREREYksrMEREpEGrUXcRM1Z07O3bs86fD2vXwmOPwSmn+KBjzhz45S9h4EAYNgxuvRU+/BDy88udSm1ZRUREROqXVmCIiEiDVpPuIuVqRnTvDtde68fOnfDWW35lxowZvnvJkiVw773Qrh2cdhqceaYPOw46qEbBiWpWiIiIiNSdAgwREWk06rSFo3VruOgiP/bv9ysv3nzTFwHNzobMTD9SUuDoowmeeSbdv23HutadorZ9REGHiIiISOUUYIiIiByoSRM44QQ/HngAvvrKhxlvvum3mcyeDbNn8xGQ07YrffZd5FdnHHOMf2w9CeVurzI4USAiIiIijYlqYIiIiFRn0CC46SaYNct3K3npJZg4kW3NW9NnxwZ46CHfkrVDB7jgAvjHP2Dz5phMpcKuKyIiIiJJQCswREREaqJt2+KtJiO7TGP4+q94tdc2vzpj0SKYOtUPMxg9Gs48k2Eb01nSuV9Unr7arisiIiIijZQCDBERkQNEuvWiMJDCpz0OhnvOgHvugdxcmD7dhxmzZsG8eTBvHm8CW1q0hQ1n+WKgJ58M7avollKFGnVdibJ4bknRdhgRERFRgCEiIg1ewnyp7d0brr/ej+++g/fegzff5O33F5HdoScZ731C8LnnIBDwqzNOPdUHGsGgLw4agZp0XRERERFpTBRgiIiIxELLlnD22YQOH8dPD/oIgGbmyFz6MsH/vARZWX7ceaevnXHyyRyTN5gFPYZWW5yziMILERERSSYq4ikiIlJLOZPPqHb1R9aqrb4ehhl5FiDrxt/Btm0wbRpcdx306QPffEPog0+Z33MYe1ObMOGR2YTu+BPMnQsFBXWepwp9ioiISGOgAENERCSGKqxZ0aoVnHUWPPoorFoFy5aRdd2t7E9JpTCQQh5G1nsL4KijoGNHuOQSeOYZWL8+4udVtxIRERFpbBRgiIiIVCIaX/qrrVlhBoMHkzHxbAoDKeAcaWmpZIweDAMGwPbtvm3r5ZdD9+5w6KHc8d4THLdyga+zUYmKupVEi8IQERERiYeYBRhm1tPM3jezpWa2xMx+ET7ezsxmmtmK8M+DwsfNzP5qZtlm9oWZjYjV3ERERCoTy5ULVdWsKL7NjMyfHkXwoT/AihV+PPwwnH66r6uxeDFXLXidf/zrTmjXDk44ASZPhlAICguLzxftbiWJsqJD4YmIiEjyiuUKjHzgZufcwUAGMMnMhgK3Au855wYC74WvA5wGDAyPa4DHYjg3ERGRCsVy5UKkygQdAwbADTf49qzbtsEHH8Dtt8OoUZCXB++/D7fdBiNHQqdOfrvJlCkEbVfxKaLRrSSefy6JEp6IiIhIfMWsC4lzbgOwIXx5l5ktBboD44Hjwnd7BvgAuCV8/FnnnAOyzKytmXUNn0dERKReRHvlQlQ1aQLHHuvH3XfD1q0waxbMnAnvvAO5uX67yUsvAfBeux7M6TOc4NhCaHccpKfX+qnj+edSUXiiDiwiIiLJp17aqJpZH+AIYB7QuSiUcM5tMLNO4bt1B9aUetja8LEyAYaZXYNfoUGvXr1iOm8REUk+1dasSCTt28OFF/rhHGRn+yBj5kyYNYv+29bSf9taOPtNSE2FMWPgpJPgBz/wKzbS0gDoc+t0gCo7qsTzzyWhQyURERGpNzEPMMysFTAVuNE5t9PMKr1rBcdcuQPOPQ48DjBy5Mhyt4uIiERLtL6kV9dqNSrMYOBAPyZNgrw8zr/iL4zL+YwbC3Ng3jyYM8eP3/7Wd0I59lg44QSGbkpjaac+ET9VfYc6DSpUEhERkZiJaYBhZmn48CLTOfdq+PCmoq0hZtYV2Bw+vhboWerhPYDI+8WJiIhIibQ0Qj2GEuoxlBsnnwE7dvh6GeHVGSxf7utqTJ/ODGBb89aw6mRfFPTEE30QUvk/OsSNwgsREZHkFbMAw/xSiynAUufcA6VumgZcBkwO/3y91PEbzOxF4EjgW9W/EBGRZBLTlRpt28K55/oBsHatDzRmzWLd1Ol037UF/vUvP8C3bD3hhJJAo2fPys9dB5FsXxERERGB2K7AGAtMBBaZ2efhY7fjg4uXzexKYDVwYfi2GcDpQDawB7gihnMTERFp9KoMBXr0gIkTYeJExnZ6k947NjA76OC99/wKjXXr4Lnn/AAYMIB7Wg7go96Hw+ZRvuOJiIiISD2KZReSD6m4rgXAiRXc3wGTYjUfERERqYQZuQd1I3TKUQSvuQYKC2HJkpIwY/ZsyM7mR2Tzo4VvwbR74dBD/eqM44+HceOgXbt4vwpAKzpEREQas3rpQiIiItKQJNOX31Du9uLLE57MKimSeeihftx4I+TnQyjkw4xZs+DDD2HRIj8eesjXyjjssJIWr8ccAx06xPFViYiISGOkAENERCSJZa3aWnw5L7+QrFVbyxfKTE2FI4/047bbYO9eyMryKzRmz/YdThYu9OOvf/WPGTasJNA49tg6bzlJplBJREREKqYAQ0REJIll9GtffDktNVDmeqWaNYPjjvMD4PvvfaAxe7Yfc+fC4sV+/O1v/j4HH1w20OjatcwpQ7nb1WFEREREqqQAQ0REJImVDg2Kt4/UVPPmvhbG8cf763v3wvz5JYHGxx/D0qV+/P3v/j6DBhE68RwCrcZSaIGy21fqgWpliIiINDyBeE9AREREEkPUwoNmzXwdjN/8Bt59F3bs8HUz7r4bTj4ZWraEr74i67Ov/f3NyNuXR9Z9j8PTT8OKFeBcdOYiIiIijYZWYIiIiEhsNWkCY8f6cfvtkJcHn35KxrvzSN1eSAGQVpBPxmv/gL8t84/p3BmOPrpkDB/ua3FESFtSREREGh+twBAREZGo6nPr9OItGhVKS4MjjyT465+zPzWNQozMk7sQ/OU1cN550LEjbNoEU6fC//wPjBoFbdvCSSfBXXf54qG7d5c77YEdVUpfFxERkYZPKzBEREQkrlwgQPCUo+CUo3xg4ZzfRvLhh37MmQPZ2X47yrvv+gelpMCIESUrNMaOJWvVruJzVtpRRURERBosBRgiIiKSWMxg0CA/fvITf2zjRvjoo5JQ47PP4JNP/HjwQQAyjvwBaeMmkRdIIS0lhYy+7ap9Km01ERERaTgUYIiIiCS5BtGJo0sXOP98PwB27YJ580oCjblzCc57lxfXrCWr16FkrF5EQxATzgAAIABJREFU8MltMGYMHHWU/zlyJDRvXm6rSX12P4k2dVMREZFkogBDREREGp70dPjBD/wAXxj088+ZcevjjFz7JcG8rbBlC7z+uh/gi4AOH07WMRMgdYDvfqKtJiIiIg2GiniKiIhITNRrEc20NBg1iimjzuG6c2/3RUC/+gqeeQauvRYOPxwKC2HBAjJe/D+a5e8npSCftLz9ZEy5H+6/Hz7+GPbujek0qy1wGmeJPj8REUluWoEhIiIiUVPT7Rkx2/pgBgMH+vHjH/tju3bB/PkE587lptefY1ezlhy3KkRw/TLIDD+uSRNfHHTMmJLRo0ds5igiIiI1ogBDREREoiZr1dbiywm3PSM9HU48EU48kXt2T8dcITc/cQfMnVsyliyBrCw/wsVB16V3pPtpx/sw48gj4YgjoFmzOL8YERGR5KMAQ0RERKImo1/74stpqYEy1+tDTVZ0OAvAkCF+XHGFP/jtt744aDjQ2PnBh3TftQVeftkP8LU0Dj/chxmjR/ufgwZBIH47c9VNRUREkoECDBEREYma0l+iG2R3jzZt4OST/QAOv+UNBnyzhpljmvhQY/58v0ojFPLj0UdLHjdqVEmgMXq075xSiWgEDo2pm4qIiEgkFGCIiIhITDSGL9POAqzo2BuuOgOuusof3LXLhxfz5/vVGvPnw9q18O67fhTp1csHGeFQI9R5QPFN1QUOkbRHjeV2Ha3oEBGRRKQAQ0RERJJOnYqHpqfDccf5ETZ60rMM3/AVjw8u8IHGJ5/A6tV+/OtfAGSNuYjA0RMoDKSQl1dA1sdLCHbP8FtSaiHa23W0okNERBKdAgwRERGRalS3ImFzenveSR8DfwwHIwUFsHx5yQqNefPIWLOYJgX55DlHWmEBGb+8Gq7M9UVBR470IxiEwYMjmlO0t+skdAFWERERFGCIiIiIVKhOKxJSUmDoUD/CBUKDe/bQ+4anSN+/h9tS1hBsug++/x4+/tiPIi1b8tJBfVjUZQBk7vDBxsCBVRYJjUbQEO8CrCIiItVRgCEiIiJRVaftGQkk6isSWrRgeae+AAQn/y/wMGzdCp9+CgsWlIzVqznyuyUcuXYJXPq6f2x6OowYUbJSY+RI6N+/Dq+uvAZfgFVERBo9BRgiIiIiFajNioQaF79s3x5OOsmPIlu2cM5Nz9Isfx+/cjkEP/qPLxI6e7YfRdq2JbNNLxZ1GQiv7PHbT/r2BbPIn78SkbyGSAqNioiIRJMCDBEREZEKRLoiIdrFL0N7Uvm8+xB/vrQjyLz3Xr/dJBTyKzRCIV8kdONGxu7YwdjcL+Ciqf7BBx3kV2occYQfI0b47ScpKbWej4iISKJQgCEiIiJSjaoCiWhvNanwfMcPgDPO8KPI+vUloUbR2LwZ3nvPjyItWsDhh5eEGkccAcOGQdOmtZ6jiIhIPCjAEBEREamDmmw1iWS7RcTn69bNj7PO8ted81tNPvus7Fi9GubO9aNIaqovMFo61Bg+vNq5iYiIxJMCDBEREZE6iHbxy1qfzwx69vTj7LNLjm/dWj7UWL4cvvjCj2eeKb5rTv/+Psy4Z2FJsNGlS51ej2pliIhItCjAEBEREYmSaHfuiMr52renz7v7gKHkPP8rf+y773x4UTrUWLQIVq7041//Knl8ly4lYcZhh/ntKAMGFN9c48KlIiIitaQAQ0RERKQSjXbVQMuWMGaMH0Xy8mDp0rKhxuefw8aN8J//+BEW6nMYTc6/k7xAKhMe+5DMI1sQPPYIaNeu1lPSSg0REamOAgwRERERgbQ0v8LisMPgssv8scJC+PprH2YsXOjHF1+Q1WUw+YEUXCBAXkE+WX/6P4LnvQLdu/sVGkXnOfxwUgoLKAioC4qIiNSdAgwRERGRBJMwqxACAejf348LLig+nLFkNe6ZhZgrJC0AGekFvtvJunV+zJhRfN/MHsOYMWQsoWtmETy4R0nA0aFDraaklRoiIslLAYaIiIhIHSXbl+ngIb1wgUUAZF53NMHe46GgAFat8rU1wis1Qmt3cvmxN7A/JZWXC/LJ/POvCa5f5k/StWtJmDFsGIds2k52ux5xfFUiIpLoFGCIiIiIJIlYFNwsPl9KCgwc6Mf55wOQ9X42e99aBmbkBQJkXXAlwflTfcHQDRv8eOstAH7fbQhzex1GaMzfCXZLh0MOgWHD/Bg40G9xERGRpKYAQ0RERKQRC+VuL7484cmsqLR6jVRGv/a+vSuQlpZCxnVXQu9f+toaRas1vviCUPZmftTlJPICKTQpyCfzxV8TfPXVkhOlpcGQITBsGNevS+OrDr1h5RDo29dvc6kDbUkREWk4FGCIiIiINGJZq7YWX87LLyRr1dZ6CzBKP0+Z4CQQ8K1YBwyA884j6/1s9r293M8xJYWsX/6B4OYQLFkCixf7QqKLFsGiRZzYbQgt939P6JgzCe5YDUOHlqzUKBrduhUHJ4kq2sGJghgRSQYKMEREREQasYx+7Ysvp6UGylyvi5p+Ua4qNCk3x/NOhN4lRUPZvRu+/JLQJ8u5OKc1hWY0KSwg84XbCS5YAAsWlD1h27Z80rIbX3XsxYQrToODD/ZBRx2CDQUEIiLxpwBDREREpBGrdBVEAql2jq1awejRZH3Xjvw1pVZqPPEywfQdfpVG6bFtG6N27GDUui/hF2+VnKd1ax9mFI2hQ+m8cyubW7WLWn0QBR0iIrGjAENEREQkSSRieHGgGq3UOKQH9D4Uxo0ruZNzsGkTE26cwuAtufx2gMGXX/qxbRvMm+cHEOo2hG8vuRvDMeGRD8jcPIvggE4lAcfAgdCkScxeq4iI1IwCDBERERGJmWiuRIhoNYkZdOnCR32G81Gf4fy26Pmdgy1bYOlSP778kqztrdifkkphIIW8AkfWV5sIPvtIyblSUnydjqFD+eXGNLLb94RPu8LgwdCyZVReU7Q7w8Si04yISKJQgCEiIiIiDU6kX9KLv9CbQadOfhx7LAAZudspfOxjcI60JqlkTDwLxvXxqzWWLvXFQ5cvh+XLGdNtCOYcobMmEFy/DHr39p1RhgzxgUbRcK7aOhvR7gwTz04zIiL1SQGGiIiISCPXEOoxRHOOkX6hLz5mRuY1R5W/z/ffw/LlhBZ8xcXLm1JgAZq6fDJf+g3B3CWQmwtvv13mIf/sfTjvDMwgdPnbBPu2Lwk2Bg0qXrUR7c4w8ew0IyJSnxRgiIiIiEijUpsv9BXe3rw5DB9O1vZW5GeHi4daE7Kemkqwh4Nly4pXaLB8OaHt+Vx12i3sT0nl5YJ8Mh//tV+tUaRHDxg8mIyDjySl2ZEUBAKkpQTI6Fu3sCFWnWZERBKNAgwRERERaTAiWakR7S/05c43sBP0PshvHykl6/1s9r61DMzICwTImjiJ4PL3fMCRnQ1r18LatQTfe4+Xuw0hq9ehZKxeRPCBHF8wtPRWlKLRpk2182sInWZERKJBAYaIiIiINCrR/kIf6fky+rUvrn+RlpZCxlUToPcN/sb8fMjJKV6tsfzZmQTXLSVYuAP27oVFi/w4UMeOvpDowIElo+h669ZVzrUyavUqIg2VAgwRERERabSq+0Jf0y/xVZ2vyqAjNdUHDwMGwBlncPvmwSXP/+23ZbaiFI8VK3znlC1bYO7c8k/YqVNxmDFpVSE5B3XzXVIGDoT09Bq9rooo6BCRRKMAQ0REREQkymoUnLRpA6NH+1FaYSGsX++DjOxs/7NorFwJmzf78fHH/KroMdPu9T87dy67WqPo8oABUXuNIiL1LWYBhpk9BZwJbHbODQsfawe8BPQBcoCLnHPbzcyAh4DTgT3A5c65T2M1NxERERFp3BrFqoFAgD6PLAQgZ/LVZW8rLIR168qHG9nZfmza5MeHH5Y77ZP9R/HegNGErnnfd0rp1w/69/c/27Wr8TQjXakR7RUdWiEiknxiuQLjH8AjwLOljt0KvOecm2xmt4av3wKcBgwMjyOBx8I/RURERETiLuG+JAcC0LOnHyecUPa2wkJfMPTAVRvZ2YT2pHDD+FvZn5LKawX5ZD5yQKeUtm2LA41Ld3Tgiy4DCU1tQnDEAP9cqVrALSLxE7PfQM65/5pZnwMOjweOC19+BvgAH2CMB551zjkgy8zamllX59yGWM1PRERERCTaEiLoCASgVy8/Dgg3smZ9xd63vyrplPKj6wiu/MBvSVm5EnbsgE8/JbRxD/+65G72p6Qy4eNdZP78dIKbs6F375LVGuGgY+imdeS27RKf1yoiSaW+I9TORaGEc26DmXUKH+8OrCl1v7XhY+UCDDO7BrgGoFevXrGdrYiIiIhIAgjlbo9Ke9SM/h3BVgDhTinXTYTeP/c3OucLhq5aRdbctezbmIazAHlA1rCxBN9ZVhJ0lPKHcEvY0CH3E0yn7JaUosvdukX9tcTqfCKSuBJlDZhVcMxVdEfn3OPA4wAjR46s8D4iIiIiIg1dKHd78eUJT2bFviWsme9s0qkTGV234x77GIC0pmlkPH4fdHrYt4JduRJWrYKVKwlt+I4f9TqDvEAKTQryyXzx1wTnvVD+tfQ5lKbn30VeIIUJj84hMz2H4IBO0KcP9O3rnzfcgjaS2hax+LNRLQ+RxFffAcamoq0hZtYV2Bw+vhboWep+PYD19Tw3EREREZGEkbVqa/HlvPxCslZtjepKg1q1hD34YD+K5vh+NvveXu7nmJJC1v97mGBgXZmQg1WryOoyhLxACoWBFPIK8sma/iHBrFdKnrBZMx9m9OnDHzYFWNumE7z8XUnA0aFDccABsf+zSVQKOyTZ1XeAMQ24DJgc/vl6qeM3mNmL+OKd36r+hYiIiIgks4x+7Ysvp6UGylyvT1UFA+XmeEIQev+g/P2WrcM9/Sk45+83dhh0K/ArOnJyYNs2WLYMli1jYtGDPvhHyQlatCgOOOjbl4zuBxMo7EmhBUhLTSGjb+XdU5LxS38yvmZJDrFso/oCvmBnBzNbC/wOH1y8bGZXAquBC8N3n4FvoZqNb6N6RazmJSIiIiLSEFS53SNBRDrH4JDuOPvc3+/asQR7n1X2Dt9+C7m5kJPD9VM+5vu0ptyQtpHgsvnw9df+9i+/9AMIAq+Ea29krF5E8KG1Jas1evcuKWLaqxeddm1lS6vI/+wiramhWh4i9S+WXUh+WMlNJ1ZwXwdMitVcREREREQasmh+sY3Vv8pHOscK79emDRx2GKE2PZlxcAoAc9MCZE6e7O+/Y0fJao2vv4acHLbNmMcJKz/h4L1bYfduWLzYjwPMB/YHUuGVnj7UOCDgoFcvQoE2xfevqqZGtGtv1PZ8CjskWSVKEU8REREREUkg8dh+UGlti7ZtYfhwP8Kubh7eJvHH02H79pJwY/VqP3JzYfVqZm2HpZ36+ZUas2dX/LwZFxIYd6mv0bE/n6y/PkOwW0HZoKNz56jX3qjJ+WoTdijokMZGAYaIiIiISIJqTDUMInkttar7YQbt2vkxYkSZm0K52/nJYx+DczRLDZB5RCrBb9eUhBzhkbFhKU0K8slzjrTCAjJe/D9Yv6zs8zRpQsbh42hyzCTyAqmkBRwZiz+CfdnQsyf06AGtW5cpNhrN1xtp2BGLDi3xoo4vciAFGCIiIiIiUmvR/DJYk7ofkTxv8Zd+M/IKHVnt+xG84OTyz1tYSLMbX6Jt/j7+Nro1we5Xl1vJwTffEPzkPV5Yt66k9saBIUerVj7I6NGjJNQ4cBx0UHHIUZPXG2nYkawdWiQ5KMAQEREREZGEE40v3RGvcAgE2NGitX/eyysJRvbsgTVrePC3r9Bt5xYmnTQE1qyBtWv9WLPG1+IId1OpVIsWZQKNX67Yy4b0DgQXOdgePt6+fbmVHJGGHYnSvUYkFhRgiIiIiIhIoxTVTi4tWsDgwXzY9wgA/vT7A4IO53y3lHCgcctDM+i66xtuHNy8JORYuxZ27oSvvvIDuKHo8e88WnKuZs0qXL1x0oqNbGzVnmDaEZCfDqnlv841hO41NaWOL1LEfAOQhmnkyJFuwYIF8Z6GiIiIiIgkqHjVPSh63qnXHVX2y/LOnWUDjaLVG6Wv79hR/RMEAtC5M3TvDt26lRmXvb2OTa3a8dafLqlwNceBc4zXn01VzxvK3c75j30MQLO0QFQ6vkTzfFI7ZhZyzo2s7eO1AkNERERERCSKqiyk2bo1DB3qR2V27y4fchQFHRs2wLp1sHmzv7xhQ7mHP1N04emfQZMm5QIOunWD7t0Zk7uaza3a+VCldetKpxOPoCOeHV8kcSnAEBERERGRRiuh2sFGqlUrGDLEj8rk5cHGjbB+vQ801q8vGaWv79jhW8zm5JQ7xQtFF568zj9nBSEH3boxas3XbGl5kA9WWrWK/HXUQbRrecSiNoi6mtQ/BRgiIiIiIiJRVC+FNNPSfKeTnj2rvt+ePWXDjVIBx3PLd5Lbtgun5YQIfr2wTG2O0l4puvDET6FlS+jSBbp2LfvzwMudOkFKSpVTq6oWRU1qeUQSJDTG2iDJSAGGiIiIiIhIFCXUl+UWLWDAAD9KCeVu5zfhmhD/HHsBmZccQtB2l1/BsW4dz6/cQ067bpySEyKY8wWsXOlHVQIB6NixXMARatebQGFXCi3AhMfnkjnxcIJDuld5qmj/+Sm8aLgUYIiIiIiIiMRIon5ZLrfNZfN+gseXr80Ryt3O7eGg49mjLyTzkmEEU/f42hsbN/pR0eUtW2DTJj9KP2/GhTDuUjAjb38eWVf8D8FFMypcyXHhF5vY2qItfNLJr+jo3Nl3aKlEvLqLqKtJ/VGAISIiIiIikmQi3eZSPujYR/D4wTB4cLn79rl1OhwR3sqRl+cLjR4QbGRs3ktKYaF/XldIxsbl8N13kJ3tRyn3FV2Y+vuSg+npPsjo3Bk6dSLUbTCpzTPIDwSY8PePyBzTkuDgbv72Nm3KdGCJZq2KKgu1SswowBAREREREYmyRC/sGOk2l1rX80hL80VAu5fdHhIE8sI1KzInjSN43+ewa1f5VRwbNvgApGgVR9HlXbv8CIcdWRkXUjhuDFiAvPx8sv74GMGscNWOJk1KVm5U97NDB0iN/OtxsnY1iXfhUgUYIiIiIiIiSayqL961qecR6ZaKMq1lW7eucFVHGc7Bt9+WCTUyVu8gsMHhCgtJM8hoA/Tv728v3Y62OmbQvn1JqNGpk6/hUdHo0IGMPm2LHxrNQq3xDggSnQIMERERERERqVZVoURNtlTU+su5GbRt60c47AgC+UUrOq47mmDv8SX337PHr9yoaCXHgce2boVvvvFjyZJqpxI04+k+I5jTZzhnNNtN8Msni8ONykIPmjat3euuQKRBR6wCkXjV/VCAISIiIiIikoQi/VIbyf0SYUtFuedr0QL69PGjOvn5PrzYtImJ90yj/Z5v+csJ3f2xLVvKj+3bOf7rEMd/HYp8gunplYcb4csDtuSyuVU7QsvW+VoepWp4xFsi1P1QgCEiIiIiIiJ1UutaGYkiNbW4C8qcvn7LycTzjqr8C3peHmzbVnG4sWVLueAjf/MWUovqd6xaVeEpQ92GsPaSu9mfksqEJ+aR+epdBPO2+q0t7dv7oCN8+cwv9/N1u+6EXphOsF/Jcdq08S1sDzx3FFZMJEJIpQBDRERERERE6qQ2tTKiJS7dRdLSSrqhRGDgLW/Qeu93LJw0otLQIyvQm/0pqRQGUsjDkdVlsC9Ium5d2Tl2G8K7RUHHgn1k/vJyguuX+RsDAWjXDtq3J9T3cNIO+SH5gRQmPDqHzPQcgp2blwlCikdaWrWvIRFCKgUYIiIiIiIiEjUNuRtHrFYZOAvwbfN0Qk07Ehw7qML7ZORup/CxjwFIa9qEjKcfhOb3+focReObb8jamMa+/Wm4QIA8M7LGnEowtNffvmtXcS2PrIMOoyAQ8PcryCdr+oclHVoOsKtJc7Y3b02vgT3LhxvhEWzXjv7frOablm156tKxBHu2qfOfS00pwBAREREREZE6awydM2KxyiDSVR2VrmLp3bvsHHO344qCjiapZNz/W+j9oL9x/36/tWXrVjKyN8NHuzBXSFpKgIyjD4WDW5UJQ9i6FbZtI33/96Tv/x4WbKrytbxXdOGv+PocBx3kV3xEOupIAYaIiIiIiIgIsdkKU5tVHbVubdukSXEtj+Ahh1A4t6hDy1EEe59V8QkLCxl681SaFOYz5ZSeBNlZEm4UjW3byo+dO0su1xMFGCIiIiIiIiIHiNZWmJqs6qjpKpZI51hlC9w137KnaQv2ABPm7iLzqjGRnTcvD3bsqDjcqGxkZ0f4yiqmAENEREREREQkLNpbYWKxqiMhWuCmpZW0gj1An1unQ0vIefiA569jW9jy/VVEREREREREJOoSscBpInQXiZRWYIiIiIiIiIjEUCIXOI1lC9xQ7vaonk8rMEREREREREQkKmHDgV1XSl+vK63AEBEREREREUli0VwhUuuaGhHQCgwRERERERERiYpY1tRQgCEiIiIiIiIiURHLmhraQiIiIiIiIiIiUROroqVagSEiIiIiIiIiCU8BhoiIiIiIiIgkPAUYIiIiIiIiIpLwFGCIiIiIiIiISMJTgCEiIiIiIiIiCU8BhoiIiIiIiIgkPAUYIiIiIiIiIpLwFGCIiIiIiIiISMJTgCEiIiIiIiIiCU8BhoiIiIiIiIgkPAUYIiIiIiIiIpLwFGCIiIiIiIiISMJTgCEiIiIiIiIiCU8BhoiIiIiIiIgkPAUYIiIiIiIiIpLwFGCIiIiIiIiISMIz51y851BrZrYLWB7veUi96wB8E+9JSL3T+56c9L4nL733yUnve3LS+56c9L4np8HOufTaPjg1mjOJg+XOuZHxnoTULzNboPc9+eh9T05635OX3vvkpPc9Oel9T05635OTmS2oy+O1hUREREREREREEp4CDBERERERERFJeA09wHg83hOQuND7npz0vicnve/JS+99ctL7npz0vicnve/JqU7ve4Mu4ikiIiIiIiIiyaGhr8AQERERERERkSSgAENEREREREREEl6DDTDM7FQzW25m2WZ2a7znI7FhZj3N7H0zW2pmS8zsF+Hjd5rZOjP7PDxOj/dcJbrMLMfMFoXf3wXhY+3MbKaZrQj/PCje85ToMbPBpT7Tn5vZTjO7UZ/3xsfMnjKzzWa2uNSxCj/f5v01/Pf9F2Y2In4zl7qo5H2/z8yWhd/b18ysbfh4HzP7vtTn/u/xm7nUVSXvfaW/283stvBnfrmZnRKfWUtdVfK+v1TqPc8xs8/Dx/WZbwSq+O4Wtb/jG2QNDDNLAb4CTgLWAp8AP3TOfRnXiUnUmVlXoKtz7lMzSwdCwDnARcBu59yf4zpBiRkzywFGOue+KXXsT8A259zkcHB5kHPulnjNUWIn/Ht+HXAkcAX6vDcqZnYMsBt41jk3LHysws93+EvNz4DT8f89POScOzJec5faq+R9PxmY5ZzLN7N7AcLvex/gzaL7ScNWyXt/JxX8bjezocALwGigG/AuMMg5V1Cvk5Y6q+h9P+D2+4FvnXO/12e+cajiu9vlROnv+Ia6AmM0kO2cW+Wc2w+8CIyP85wkBpxzG5xzn4Yv7wKWAt3jOyuJo/HAM+HLz+B/IUrjdCKw0jmXG++JSPQ55/4LbDvgcGWf7/H4//l1zrksoG34f5CkganofXfOveOcyw9fzQJ61PvEJOYq+cxXZjzwonNun3PuayAb///+0sBU9b6bmeH/QfKFep2UxFQV392i9nd8Qw0wugNrSl1fi77UNnrhZPYIYF740A3hpUZPaStBo+SAd8wsZGbXhI91ds5tAP8LEugUt9lJrF1C2f+p0ee98avs862/85PHT4D/lLre18w+M7PZZjYuXpOSmKrod7s+88lhHLDJObei1DF95huRA767Re3v+IYaYFgFxxreXhiJmJm1AqYCNzrndgKPAf2B4cAG4P44Tk9iY6xzbgRwGjApvAxRkoCZNQHOBl4JH9LnPbnp7/wkYGa/BvKBzPChDUAv59wRwE3A82bWOl7zk5io7He7PvPJ4YeU/YcKfeYbkQq+u1V61wqOVfl5b6gBxlqgZ6nrPYD1cZqLxJiZpeE/AJnOuVcBnHObnHMFzrlC4Am0tLDRcc6tD//cDLyGf483FS0rC//cHL8ZSgydBnzqnNsE+rwnkco+3/o7v5Ezs8uAM4EJLlycLbx9YGv4cghYCQyK3ywl2qr43a7PfCNnZqnAecBLRcf0mW88KvruRhT/jm+oAcYnwEAz6xv+l7pLgGlxnpPEQHh/3BRgqXPugVLHS++NOhdYfOBjpeEys5bhwj+YWUvgZPx7PA24LHy3y4DX4zNDibEy/yqjz3vSqOzzPQ34cbhSeQa+4NuGeExQos/MTgVuAc52zu0pdbxjuJgvZtYPGAisis8sJRaq+N0+DbjEzJqaWV/8ez+/vucnMfUDYJlzbm3RAX3mG4fKvrsRxb/jU6M853oRrlR9A/A2kAI85ZxbEudpSWyMBSYCi4raLAG3Az80s+H4JUY5wE/jMz2Jkc7Aa/53IKnA8865t8zsE+BlM7sSWA1cGMc5SgyYWQt8h6nSn+k/6fPeuJjZC8BxQAczWwv8DphMxZ/vGfjq5NnAHnxXGmmAKnnfbwOaAjPDv/OznHPXAscAvzezfKAAuNY5F2kRSEkwlbz3x1X0u905t8TMXga+xG8rmqQOJA1TRe+7c24K5etcgT7zjUVl392i9nd8g2yjKiIiIiIiIiLJpaFuIRERERERERGRJKIAQ0REREREREQSngIMEREREREREUl4CjBEREREREREJOEpwBARERERERGRhKcAQ0REREREREQSngIMERGRODGzAjP73MwWm9krZtbCzPqY2eIanudyM+sWpTmNNLO/RuNc1TzPDDNrW8vHnmNmQ6M9p9oys+PM7Fszm1Hq2EAze9PMVppZyMzeN7Njqjjz6wAfAAAF0klEQVRHSzPbamZtDjj+bzO7yMwuNrNsM3szlq9FREQkkSnAEBERiZ/vnXPDnXPDgP3AtbU8z+VAVAIM59wC59zPo3GuipgXcM6d7pzbUcvTnAMkTIARNsc5dzqAmTUDpgOPO+f6O+eCwM+AfpU92Dn3HfAO/rURPk8b4GjgTefcS8BVMZy/iIhIwlOAISIikhjmAAPCl1PM7AkzW2Jm75hZcwAzG25mWWb2hZm9ZmYHmdkFwEggM7yao7mZnWhmn5nZIjN7ysyahh+fY2Z3mdmn4duGHDiJ8GqCN8OX7ww//gMzW2VmPw8fv9fMri/1mDvN7GYza2Vm75U6//jw7X3MbKmZPQp8CvQMz6VD+PZ/h1cpLDGza0qdd7eZ3W1mC8Ovu7OZHQWcDdwXfr39D5j/P8zssfCKh1Vmdmz4NSw1s3+Uut9jZrYg/Jx3lTr+/9u7txArqzCM4/+HDCozO1BG0YkSo0kb81ASUVokkRTWqEFEFAhFJQVGgaEXRQTZ4aYDVDAFpaAFKolBBwsqQ0unMjLIVIhQjMkiTU2fLtbatpMZ3XbAUZ7f1Wbt9a31ft++mf3u913zuKSv6zOeXccm1SqZLkkftvBZ3gJ8YnthY8D2V7Y763r9a0zL6+d0Q502B7i5aZ2JwBLbW1vYMyIi4rCXBEZERMRBJqkfcC3wZR0aDDxruw34Gbipjr8KPGh7WJ07y/Z8YAVwi+12wEAnMMX2UKAfcFfTdpttXww8D0xvIbzzgfHAaGCWpCOBucCUpjmTgXnA78DEuv5Y4ElJqnOGAK/aHm57/V573FGrFEYC0ySdVMf7A8tsXwR8CEy1/TGwEHigVq9810PMJwDjgPuBRcDTQBswVFJ7nTPD9khgGHCFpGGSTqQkDdrqM360zp0JjK9xXN/CM2ujJGp6MwN4z/YoynN6QlJ/YAkwoun+b6YkNSIiIoIkMCIiIg6moyWtoiQgNgAv1/Hvba+qrz8Dzq7tBMfb/qCOvwL0dKbCkHr9t73Me7N53RZifMv2dtubgU3AINsrgVMknSbpIqDb9gZAwGOSvgDeAU4HBtV11tte1sse0yR1AcuAMygJHChtNY0zH1qNF2CRbVOSPBttf2l7N7C6aY3Jkj4HVlISDhcAv1CSMC9JuhFoVD58BHRKmgoc0WIMe9Rqma8kNZ79NcBD9bNfChwFnGl7ByU501GrU9opbSURERFB+VUmIiIiDo5ttWpij1qwsL1paBdw9AGsqf2831h7F639HbB3LI1r5gMdwKmUigworRMnAyNs75S0jvLlHOC3HoOVrgSuBsbY3ippadM1O2si4kDibY55917x7wb6STqHUn0yynZ3bS05yvYfkkYDV1GqH+4Bxtm+U9IlwHXAKknttn/ax/6raUoa2Z4oaSQwu3HbwE221/Rw7Rzg4Tpnge2dLd5zRETEYS8VGBEREYcA21uAbkmX16FbgUY1xq/AgPr6G0rFxnk9zPsvzaV8ye+gJDMABgKbavJiLHBWC+sMpFRwbK1nclzawjXN9/tPHEdJqGyRNIjSvoOkY4GBthcD91EqIJB0ru1Pbc8ENlOqRPbldeAySc3tJsc0vX4buLfRXiNpeNN771MqUO4m7SMRERF/kwqMiIiIQ8dtwAuSjgHWArfX8c46vg0YU8fn1bM1lgMv/NeB2F4taQDwg+0f6/BrwCJJK4BVlGTK/iwB7qxtJ2sobST7Mxd4sR4q2tHLORj7ir1L0kpKpcRaSosIlKTIApX/IiLKGRpQzqgYXMfeBbr2s/42SROApyQ9A2ykJF0aZ2o8AjwDfFGTGOuACfXa3ZLeACZRzv2IiIiISn9VZkZERETEgaptMNNtTzgc9omIiOir0kISERER8e/sAC6UtPj/2kDSFOA5oPv/2iMiIqKvSwVGRERERERERPR5qcCIiIiIiIiIiD4vCYyIiIiIiIiI6POSwIiIiIiIiIiIPi8JjIiIiIiIiIjo8/4ETepQs8BPY6AAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 1080x432 with 1 Axes>"
      ]
     },
     "execution_count": 20,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "fig.tight_layout()\n",
    "ax.legend()\n",
    "fig"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "metadata": {},
   "outputs": [],
   "source": [
    "if Saveplots:\n",
    "    fig.savefig(\"ExampleLikelihoodFit\", dpi=600)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Questions:\n",
    "---\n",
    " 0. Given the data (first without signal), see if you can write a ChiSquare fit to match it.\n",
    "    Remember to give the fit good starting values. Alternatively, consider the ChiSquare fit\n",
    "    that is provided.\n",
    "    Can you also write an unbinned likelihood fit? Remember, that contrary to the ChiSquare fit,\n",
    "    this should not include the binwidth in the normalisation of the PDF. And if you feel on top\n",
    "    of things, then perhaps a binned one?\n",
    "\n",
    " 1. Run the program without a signal peak and high statistics. Are the two/three fit\n",
    "    results similar, and was that what you would expect? And are the models good?\n",
    "\n",
    " 2. Now lower the number of background events from 25000 to 500 and ask yourself the above\n",
    "    questions. Also, think about how you _quantify_ if the models are good.\n",
    "\n",
    " 3. Switch on the signal of 50 events, both in the data and in the fit(s). How confident are\n",
    "    you, that you can see the signal in the data? Lower the number of signal events, until the\n",
    "    point where you are no longer certain, that you can see a peak by eye. Then consider the\n",
    "    significance (i.e. the number of sigma's) of Nsig. Does the size of this correspond well\n",
    "    with your observation?\n",
    "\n",
    " 4. Try to set the number of signal events to 20, but choose a random location of the signal\n",
    "    peak (mux) and set the initial value to say 125, i.e. pretending that you don't know, if or\n",
    "    where there is a peak.\n",
    "    See if the fit can find the position of the signal peak and hence converge correctly.\n",
    "\n",
    " 5. Also try to remove the signal from data, and see if the fit (including signal) still \"finds\"\n",
    "    any peaks. This is the pitfall of searching for peaks, that may very well not be there!!!\n",
    "\n",
    " 6. Given 500 background and 20 signal events, ask yourself what helps the most to discover\n",
    "    a signal peak:\n",
    "     - Removing half the background below the peak.\n",
    "     - Improving the resolution (sigma) by a factor two.\n",
    "     \n",
    "    Discuss first what you think (and why) with those around you, and then try it out with this fit.\n",
    "\n",
    "Advanced questions:\n",
    "---\n",
    " 0. In general, how far can the initial parameter values stray from the true ones,\n",
    "    and still find the correct signal peak and hence true minimum of the Chi2/LLH?\n",
    "    Note here, that also the width and the normalisation are generally unknown in real data!\n",
    "\n",
    " 1. Make a loop to repeat the fits many times, and see how the result of the\n",
    "    fit compares to the true values. Which fit is least biased? And does this depend\n",
    "    on the value of the input parameters?\n",
    "\n",
    " 2. The likelihood value does not in itself tell you much. But if repeated many times\n",
    "    (as above) it can be used to evaluate the fit at hand. Try to produce new data and\n",
    "    fit it 1000 times, and see the distribution of likelihood values obtained... and\n",
    "    then compare the \"original data LLH value\" with this distribution! From this, one\n",
    "    can actually obtain a probability for the LLH in the same way as for the ChiSquare.\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "executable": "/usr/bin/env python",
  "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"
  },
  "main_language": "python"
 },
 "nbformat": 4,
 "nbformat_minor": 2
}