{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "# Reproducing HERA Memo 4" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "In this demo, we are going to reproduce the results of \n", "[Pober2015](http://reionization.org/wp-content/uploads/2015/05/HERA4_sensecalc.pdf), a \n", "memo which estimated the sensitivity of the HERA array in its ultimate configuration. \n", "This memo was based on [Pober+2014](https://arxiv.org/pdf/1310.7031.pdf), which looked \n", "at the sensitivity of a \"concept HERA\" (which differs somewhat from the final \n", "instrument), as well as a few other well-known arrays.\n", "\n", "The purpose of the demo is firstly to show how to estimate the sensitivity of a \n", "realistic array in a realistic situation. Secondly, it is a regression test, to ensure\n", "that the results given by the new code are consistent with the original code used in\n", "that memo.\n", "\n", "The \"reproduction\" here is not expected to be exact -- there are several small tweaks \n", "made in the current version that are not present in the code used for the memo (not to \n", "mention the code used for Pober+2014, which was not released at the time), for example \n", "cosmological calculations are done with astropy instead of approximate fitting \n", "functions. Nevertheless, we expect reasonable consistency when the settings are \n", "tuned to match the original code. " ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "
\n", "\n", "A Note On The Source Data\n", "\n", "The data to which we compare in this notebook is precisely that used in Pober15. In that\n", "memo, the thermal noise estimates are not presented directly (either as a table or\n", "in plots), and only the ultimate SNR of the entire dataset is shown. We sourced the data\n", "for this notebook by installing and running the original code used for the memo with \n", "identical settings to those presented in the memo. This cannot be repeated easily\n", "since the original code requires Python <= 2. We include the resulting thermal \n", "sensitivities in this repo for posterity.\n", "\n", "
" ] }, { "cell_type": "code", "execution_count": 1, "metadata": {}, "outputs": [], "source": [ "import matplotlib.pyplot as plt\n", "import numpy as np\n", "from astropy import constants\n", "from astropy import units as un\n", "from astropy.cosmology import Planck15\n", "from astropy.cosmology.units import littleh\n", "\n", "import py21cmsense as p21c" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Original Data" ] }, { "cell_type": "code", "execution_count": 2, "metadata": {}, "outputs": [], "source": [ "memo_data = np.load(p21c.data.PATH / \"hera331drift_mod_0.135.npz\")" ] }, { "cell_type": "code", "execution_count": 3, "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAkMAAAHHCAYAAAC88FzIAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjkuMCwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy80BEi2AAAACXBIWXMAAA9hAAAPYQGoP6dpAABcfElEQVR4nO3dd1xWdf/H8dfFFgUUURBFceIGxW2u0lw5y5WaI63uuyyzrGyZLTNLu0u6G3dmmavMUZm5R5mluVdO3AIOhgwZ13V+f5D8IkQBgQNc7+fjwePhdc65zvlcHOC8Pec7LIZhGIiIiIjYKQezCxARERExk8KQiIiI2DWFIREREbFrCkMiIiJi1xSGRERExK4pDImIiIhdUxgSERERu6YwJCIiInZNYUhERETsmsKQSD565ZVXsFgseXrvnDlzsFgsnDx5Mn+L+puTJ09isViYM2dOgR0jt4piTXl19OhR7r77bry8vLBYLCxbtszskgpNx44dadiwodlliOSJwpAIcODAAYYNG0blypVxdXXF39+foUOHcuDAAbNLK3SPP/44FouFY8eOZbvNCy+8gMViYe/evYVYWdE3YsQI9u3bxxtvvMHcuXNp1qyZ2SVla+PGjVgslowvZ2dnatSowQMPPMCJEyfMLu+WFi1axLBhw6hduzYWi4WOHTvecLt/fs6/f/3222+FW7QUWU5mFyBitiVLljBkyBC8vb158MEHqV69OidPnuSzzz5j8eLFLFy4kH79+uVoXy+++CLPPfdcnuoYPnw4gwcPxtXVNU/vzy9Dhw7lgw8+YP78+bz88ss33GbBggU0atSIxo0b3/bxqlWrRlJSEs7Ozre9LzMlJSWxdetWXnjhBR577DGzy8mxxx9/nObNm5OamsrOnTv55JNPWLFiBfv27cPf39/s8rL13//+lx07dtC8eXMuX758y+2vf86/q1WrVkGVJ8WMwpDYtePHjzN8+HBq1KjB5s2bqVChQsa6J554gnbt2jF8+HD27t1LjRo1st1PQkICpUuXxsnJCSenvP1aOTo64ujomKf35qeWLVtSq1YtFixYcMMwtHXrVsLDw3nrrbdu6zhpaWnYbDZcXFxwc3O7rX0VBRcvXgSgbNmyt9z2+s9LUdCuXTvuu+8+AEaNGkWdOnV4/PHH+eKLL5g0aZLJ1WVms9lISUnBzc2NuXPnUrlyZRwcHHL0eO7vn1Pkn/SYTOza9OnTSUxM5JNPPskUhAB8fHz4+OOPSUhI4O23385Yfr1d0MGDB7n//vspV64cd9xxR6Z1f5eUlMTjjz+Oj48PHh4e9O7dm3PnzmGxWHjllVcytrtRm6HAwEDuuecefvnlF1q0aIGbmxs1atTgyy+/zHSMK1eu8PTTT9OoUSPKlCmDp6cn3bt3Z8+ePXn6vgwdOpQ///yTnTt3Zlk3f/58LBYLQ4YMISUlhZdffpnQ0FC8vLwoXbo07dq1Y8OGDZnec71d0DvvvMN7771HzZo1cXV15eDBgzdsM7R3715GjhxJjRo1cHNzw8/Pj9GjR2e5A3D9+33s2DFGjhxJ2bJl8fLyYtSoUSQmJmap/auvvqJFixa4u7tTrlw52rdvz+rVqzNts3LlStq1a0fp0qXx8PCgZ8+et3xc+sorr1CtWjUAJk6ciMViITAwMFONN/p5SUtL47XXXsv4fgQGBvL888+TnJycaf/Xfw42btxIs2bNKFWqFI0aNWLjxo1A+t3NRo0a4ebmRmhoKLt27bppvTdz5513AhAeHp6x7MMPP6RBgwYZj5AfffRRYmJibvj+HTt20KZNG0qVKkX16tX56KOPsmyTnJzM5MmTqVWrFq6urgQEBPDMM89k+dwWi4XHHnuMefPmZRz/p59+AiAgIAAHh9xdwq5evUpaWlqu3iP2QWFI7Nr3339PYGAg7dq1u+H69u3bExgYyIoVK7KsGzBgAImJibz55puMHTs222OMHDmSDz74gB49ejBt2jRKlSpFz549c1zjsWPHuO++++jSpQvvvvsu5cqVY+TIkZku0CdOnGDZsmXcc889zJgxg4kTJ7Jv3z46dOjA+fPnc3ys64YOHQqkB5+/s1qtfP3117Rr146qVasSFxfH//73Pzp27Mi0adN45ZVXuHjxIl27dmX37t1Z9vv555/zwQcf8NBDD/Huu+/i7e19w+OvWbOGEydOMGrUKD744AMGDx7MwoUL6dGjB4ZhZNl+4MCBXL16lalTpzJw4EDmzJnDlClTMm0zZcoUhg8fjrOzM6+++ipTpkwhICCA9evXZ2wzd+5cevbsSZkyZZg2bRovvfQSBw8e5I477rhpw/b+/fszc+ZMAIYMGcLcuXN57733Mm1zo5+XMWPG8PLLL9O0aVNmzpxJhw4dmDp1KoMHD85yjGPHjnH//ffTq1cvpk6dSnR0NL169WLevHk8+eSTDBs2jClTpnD8+HEGDhyIzWbLtt6bOX78OADly5cH0sPco48+ir+/P++++y733nsvH3/8MXfffTepqamZ3hsdHU2PHj0IDQ3l7bffpkqVKvzrX/9i9uzZGdvYbDZ69+7NO++8Q69evfjggw/o27cvM2fOZNCgQVnqWb9+PU8++SSDBg3iP//5T0bIzK1Ro0bh6emJm5sbnTp14o8//sjTfqSEMkTsVExMjAEYffr0uel2vXv3NgAjLi7OMAzDmDx5sgEYQ4YMybLt9XXX7dixwwCM8ePHZ9pu5MiRBmBMnjw5Y9nnn39uAEZ4eHjGsmrVqhmAsXnz5oxlUVFRhqurq/HUU09lLLt27ZphtVozHSM8PNxwdXU1Xn311UzLAOPzzz+/6Wc2DMNo3ry5UaVKlUz7/emnnwzA+Pjjjw3DMIy0tDQjOTk50/uio6MNX19fY/To0VmO6+npaURFRWWp8581JSYmZqlnwYIFWb4X17/ffz+WYRhGv379jPLly2e8Pnr0qOHg4GD069cvy/fJZrMZhmEYV69eNcqWLWuMHTs20/qIiAjDy8sry/J/uv45pk+fnml5dj8vu3fvNgBjzJgxmZY//fTTBmCsX78+Y9n1n4Nff/01Y9mqVasMwChVqpRx6tSpjOUff/yxARgbNmy4ab0bNmwwAGP27NnGxYsXjfPnzxsrVqwwAgMDDYvFYmzfvt2IiooyXFxcjLvvvjvT923WrFkZ772uQ4cOBmC8++67GcuSk5ONkJAQo2LFikZKSophGIYxd+5cw8HBwfj5558z1fPRRx8ZgLFly5aMZYDh4OBgHDhw4KafpUGDBkaHDh1uuG7Lli3Gvffea3z22WfG8uXLjalTpxrly5c33NzcjJ07d950v2I/dGdI7NbVq1cB8PDwuOl219fHxcVlWv7II4/c8hjXb+n/+9//zrR83LhxOa6zfv36me5cVahQgaCgoEw9flxdXTMeGVitVi5fvkyZMmUICgq64aOunBg2bBhnz55l8+bNGcvmz5+Pi4sLAwYMANLbObm4uADp/+O/cuUKaWlpNGvW7IbHvffee7M8jryRUqVKZfz72rVrXLp0iVatWgHccL//PBft2rXj8uXLGeds2bJl2Gw2Xn755SyPVq4/1lyzZg0xMTEMGTKES5cuZXw5OjrSsmXLLI/+cuufNf74448ATJgwIdPyp556CiDL3cj69evTunXrjNctW7YE0h9rVa1aNcvynPYIGz16NBUqVMDf35+ePXuSkJDAF198QbNmzVi7di0pKSmMHz8+0/dt7NixeHp6ZqnRycmJhx9+OOO1i4sLDz/8MFFRUezYsQOAb775hnr16lG3bt1M3+frj+f++X3u0KED9evXz9FnuZE2bdqwePFiRo8eTe/evXnuuef47bffsFgsRa5NlJhHDajFbl0POddDUXayC03Vq1e/5TFOnTqFg4NDlm1z04vl7xe668qVK0d0dHTGa5vNxn/+8x8+/PBDwsPDsVqtGeuuP+7IrcGDBzNhwgTmz59Px44duXbtGkuXLqV79+6UK1cuY7svvviCd999lz///DPTY5MbfX9y8j2D9DZQU6ZMYeHChURFRWVaFxsbm2X7f36PrtcXHR2Np6cnx48fx8HB4aYX1aNHjwL/32bmnzw9PXNUe3b++dmv/2z882fBz8+PsmXLcurUqUzL//kZvby8gPS2Mzda/vefj5t5+eWXadeuHY6Ojvj4+FCvXr2MTgDXawgKCsr0HhcXF2rUqJGlRn9//ywNw+vUqQOktxtr1aoVR48e5dChQ9mG4n+e75z+zORGrVq16NOnD0uWLMFqtRaJjgtiLoUhsVteXl5UqlTplmPl7N27l8qVK2e5GP797kVByu4PtfG3tjNvvvkmL730EqNHj+a1117D29sbBwcHxo8fn+e2IxUrVqRLly58++23hIWF8f3333P16tWM9kSQ3iB55MiR9O3bl4kTJ1KxYkUcHR2ZOnVqRtuTv8vp92zgwIH8+uuvTJw4kZCQEMqUKYPNZqNbt243/Dw5+R7dyvX9zp07Fz8/vyzr89pL8LrsPntOB+nM7jPe7mdv1KgRnTt3ztG2+cFms9GoUSNmzJhxw/X/DHcF9XsWEBBASkoKCQkJtx10pfhTGBK7ds899/Dpp5/yyy+/ZPTw+buff/6ZkydPZrr1nxvVqlXDZrMRHh5O7dq1M5bfbEDDvFi8eDGdOnXis88+y7Q8JiYGHx+fPO936NCh/PTTT6xcuZL58+fj6elJr169Mh23Ro0aLFmyJNNFffLkyXk+ZnR0NOvWrWPKlCmZuvZfv3OTFzVr1sRms3Hw4EFCQkKy3QbSQ2BhhIPrPxtHjx6lXr16GcsjIyOJiYnJ6J1mpus1HD58ONPQEikpKYSHh2f5Pp0/fz7LsAFHjhwByGj4XLNmTfbs2cNdd92V59Ha88OJEydwc3OjTJkyptUgRYfaDIldmzhxIqVKleLhhx/O0m37ypUrPPLII7i7uzNx4sQ87b9r165Aetfkv/vggw/yVnA2HB0ds9wJ+Oabbzh37txt7bdv3764u7vz4YcfsnLlSvr3759pTKDrdyX+fuzff/+drVu35vmYN9onkKV3Vm707dsXBwcHXn311Sx3lq4fp2vXrnh6evLmm29m6SUF/z+OUH7p0aMHkPVzXb9jkpsehwWlc+fOuLi48P7772c6H5999hmxsbFZakxLS+Pjjz/OeJ2SksLHH39MhQoVCA0NBdLv+p07d45PP/00y/GSkpJISEjI189wo/O2Z88evvvuO+6+++5cd8+Xkkl3hsSu1a5dmy+++IKhQ4fSqFGjLCNQX7p0iQULFmTcNcit0NBQ7r33Xt577z0uX75Mq1at2LRpU8b/lvPrf8b33HMPr776KqNGjaJNmzbs27ePefPm3XSgyJwoU6YMffv2zehi//dHZNePu2TJEvr160fPnj0JDw/no48+on79+sTHx+fpmJ6enrRv3563336b1NRUKleuzOrVqzONe5NbtWrV4oUXXuC1116jXbt29O/fH1dXV7Zv346/vz9Tp07F09OT//73vwwfPpymTZsyePBgKlSowOnTp1mxYgVt27Zl1qxZea7hn4KDgxkxYgSffPIJMTExdOjQgW3btvHFF1/Qt29fOnXqlG/HyqsKFSowadIkpkyZQrdu3ejduzeHDx/mww8/pHnz5gwbNizT9v7+/kybNo2TJ09Sp04dFi1axO7du/nkk08yRhgfPnw4X3/9NY888ggbNmygbdu2WK1W/vzzT77++mtWrVqVo2lMNm/enNG4/+LFiyQkJPD6668D6UNitG/fHoBBgwZRqlQp2rRpQ8WKFTl48CCffPIJ7u7utz1wqJQcCkNi9wYMGEDdunWZOnVqRgAqX748nTp14vnnn7/tySe//PJL/Pz8WLBgAUuXLqVz584sWrSIoKCgfBt5+fnnnychIYH58+ezaNEimjZtyooVK/I8NcjfDR06lPnz51OpUqUsjYtHjhxJREQEH3/8MatWraJ+/fp89dVXfPPNNxkDAubF/PnzGTduHGFhYRiGwd13383KlStva3qIV199lerVq/PBBx/wwgsv4O7uTuPGjRk+fHjGNvfffz/+/v689dZbTJ8+neTkZCpXrky7du0YNWpUno+dnf/973/UqFGDOXPmsHTpUvz8/Jg0adJtPWbMb6+88goVKlRg1qxZPPnkk3h7e/PQQw/x5ptvZplCpVy5cnzxxReMGzeOTz/9FF9fX2bNmpVpHC4HBweWLVvGzJkz+fLLL1m6dCnu7u7UqFGDJ554IqPB9a2sX78+y1hSL730EpD+mPZ6GOrbty/z5s1jxowZxMXFUaFCBfr3758x6KMIgMXITQtDEckXu3fvpkmTJnz11VdZ7raIiEjh0sNSkQKWlJSUZdl7772Hg4NDxv9eRUTEPHpMJlLA3n77bXbs2EGnTp1wcnJi5cqVrFy5koceeihLN2IRESl8ekwmUsDWrFnDlClTOHjwIPHx8VStWpXhw4fzwgsv3PbYNSIicvsUhkRERMSuqc2QiIiI2DWFIREREbFrarBwCzabjfPnz+Ph4WHq0PEiIiKSc4ZhcPXqVfz9/W850rjC0C2cP39ePX5ERESKqTNnzlClSpWbbqMwdAseHh5A+jdTMxuLiIgUD3FxcQQEBGRcx29GYegWrj8a8/T0VBgSEREpZnLSxEUNqEVERMSuKQyJiIiIXVMYEhEREbumMCQiIiJ2TWFIRERE7JrCkIiIiNg1hSERERGxawpDIiIiYtcUhkRERMSuKQyJiIiIXVMYEhEREbumMCQiIiJ2TWFIREREzGEYHN78DYbNamoZJT4MxcTE0KxZM0JCQmjYsCGffvqp2SWJiIgIsGfV5wStH8Ohtztjs5oXiJxMO3Ih8fDwYPPmzbi7u5OQkEDDhg3p378/5cuXN7s0ERERu5UQe5nKv00B4IpPKA6OjqbVUuLvDDk6OuLu7g5AcnIyhmFgGIbJVYmIiNi3Q189jQ8xnLb4E3r/q6bWUuTD0ObNm+nVqxf+/v5YLBaWLVuWZZuwsDACAwNxc3OjZcuWbNu2LdP6mJgYgoODqVKlChMnTsTHx6eQqhcREZF/OrF7E02jlgJwueNblPrrpoVZinwYSkhIIDg4mLCwsBuuX7RoERMmTGDy5Mns3LmT4OBgunbtSlRUVMY2ZcuWZc+ePYSHhzN//nwiIyMLq3wRERH5G2taKvwwHgeLwTaPLjTp0Mfskop+GOrevTuvv/46/fr1u+H6GTNmMHbsWEaNGkX9+vX56KOPcHd3Z/bs2Vm29fX1JTg4mJ9//jnb4yUnJxMXF5fpS0RERPLHzsXTqJF2glijNIFDZ5pdDlAMwtDNpKSksGPHDjp37pyxzMHBgc6dO7N161YAIiMjuXr1KgCxsbFs3ryZoKCgbPc5depUvLy8Mr4CAgIK9kOIiIjYiUvnT1D/zw8AOFD/SSr6FY1rbLEOQ5cuXcJqteLr65tpua+vLxEREQCcOnWKdu3aERwcTLt27Rg3bhyNGjXKdp+TJk0iNjY24+vMmTMF+hlERETsxbn5j1Oaaxx0qkfL+yaYXU6GEt+1vkWLFuzevTvH27u6uuLq6lpwBYmIiNihAxsWEhz/M2mGA069/4OjiV3p/6lY3xny8fHB0dExS4PoyMhI/Pz8TKpKRERE/u5aQhzlN78IwO++Q6jTuKXJFWVWrMOQi4sLoaGhrFu3LmOZzWZj3bp1tG7d2sTKRERE5Lp9857Hz7jIBXxoPOxNs8vJosg/JouPj+fYsWMZr8PDw9m9ezfe3t5UrVqVCRMmMGLECJo1a0aLFi147733SEhIYNSoUbd13LCwMMLCwrCaODy4iIhIcXf60HZCzs0HC5xr/SrNPMuaXVIWFqOID8e8ceNGOnXqlGX5iBEjmDNnDgCzZs1i+vTpREREEBISwvvvv0/LlvlzCy4uLg4vLy9iY2Px9PTMl32KiIjYA8Nm5ejUttRJPcSOUm1p+swKLBZLoRw7N9fvIh+GzKYwJCIikjc7l86k6Z5XSDDciB29Bf9qtQrt2Lm5fhfrNkMiIiJSNMVePEfNPe8AsKvWvws1COWWwpCIiIjku+PzxuNFPMccqtNi0CSzy7kphaFshIWFUb9+fZo3b252KSIiIsXK4a0/0DRmNTbDQkr3Gbi4uJhd0k2pzdAtqM2QiIhIzqUmJxExLZQA2zl+9e5Hm8fnmFKH2gyJiIiIKXYveIUA2zkuUpZ6w94xu5wcKfLjDImIiEjRlhAdyfGty7EdWUVw9AawwPGmL9DK28fs0nJEYUhERERyxzA49+fvRGz/Ds+zG6iZfIjGlr9a3Vhgu3s7Wt4zxtwac0FhSERERG4pOSGG47+v4NqBlVS98guVjWgqX19pgWOWQM5XbIdHox6EtOyCxaH4tMRRGBIREZEbunTqIKd/W4LbyXXUTtxDfcv/T1GVaLhysFQTkgI7U6VFH2pWr02tQhpdOr8pDGVDc5OJiIi9MayphO9cz5Vdy/GL2EQV21kyWv1Y4DR+nCp/B271ulG3VXealSljZrn5Rl3rb0Fd60VEpCRLir3MsV+XYj28khoxW/EkIWNdiuHIQZdGxAXciV+zPtSqG4yDQ/G4+5Ob67fuDImIiNgTwyDq5D7O/LaUMifXUvPafhpZbBmrrxgeHPZohVGnO3Xa9CHEp3j0CLsdCkMiIiIlnGFNJXzXBqJ3LcPvwnoq2y5Q8fpKCxy3VOVcxQ54NO5F/eadaF3ER4zObwpDIiIiJVBKYhzHti4n+cAPVL+yhRpczViXbDhx0DWYq1XvonKLftSoXY+axbTxc35QGBIRESkh4qJOc2LLYpyP/kTthB3Ut6RlrIsxynDQozUE9SCobR+aeJc3sdKiRWFIRESkuDIMIo7t4vxvi/E6vYaaqUcIub7OAmfw46RPB9wb9aZBqy60cXU1sdiiS2EoG+paLyIiRZFhs3J672Yubf+WShfW4m+7gN9f62yGhT+d6nDR/04qNOtP3YbNCHAsPoMfmkVd629BXetFRMRstpRrHN/+I/G7l1P10kbKGzEZ65INZ/a7NSWx+t0Etr6XgGrVzSu0CFHXehERkWIuOSGaY1uWkXbge2rG/kptkjLWxRnuHCjTGiOoJ0F39CVU7X9ui8KQiIhIEZEQHcHxzV/jeOQHasfvoMHfGkBHUY4jZdvj0rA3Ddr0oLW7u4mVliwKQyIiIiaKizrNiZ8X4nZ0BbWT9mSa/f0U/pys0AnPJv1o0LwTdzjrsl0Q9F0VEREpZDHnj3Hy5wWUPvEjtZMPZuoBdtihBhH+XajQYgB1GzajWjGZ/qI4UxgSEREpBJdO7uf0loWUPbmSGqnH/j8AAQcc63KpSlf8Ww+kTlADgux4AEQzKAyJiIgUkIsndnPmlwX4nFpJVeupjBngrYaFA84Nia7WjaptB9GgRm1T67R3CkPZ0DhDIiKSF9cDUIVTPxJgPU2Fv5anGo7scw0hvnp3atwxiMYBVU2tU/6fxhm6BY0zJCIit/LPAHRdiuHIPrdQEmvdQ+32A/HzrWRilfZF4wyJiIgUsOzuAF0PQAm1ehHUfhChvr6m1im3pjAkIiKSQ5dPHeDU5q+ocPKHLAFor2soibUVgIojhSEREZGbiD1/nPBNc/E68T3VU49xfaznjABUqxd1Ogykma/fTfcjRZfCkIiIyD/EXzrD8Y3zKH10ObX+Ng5QmuHAHtemxNfqTVCHQQpAJYTCkIiICJAUE8XRTfNxObSUOkl7CP5rJGibYWGfcyOia/SidochhFYOMLlSyW8KQyIiYrdSE2M5snkRln2LqR3/B40tfw2nYoEDDkFcDLyHau2GEly9prmFSoFSGBIREbtiS03m2G/LubZjIbVjfqYBKekrLHDYUoMLAT2o3HYo9evUw6KRoO2CwpCIiJR8Nhun96znytavqB61hjrEZ6w6RSXCK/WgYpuh1GvYVFNh2CGFoWxoBGoRkeIv8ugOzm3+kspnV1DVuMj1MZ8vGmU5WL4Lni3up1GzDlRzcjS1TjGXRqC+BY1ALSJSvMRFhHNi/ed4n1hO1bSTGcuvGqXY49EBp5CBhLTrhZuri3lFSoHTCNQiImJXUhJiObJxPs77F1I7cQ8hf/UESzac2OPWgpT699Kw4wDu8PIyuVIpihSGRESkWDKsqYRvW0H8tq+oHb2Jhn9rCL3XsSFXavUjqNNQWvhpPjC5OYUhEREpViKPbOfC5s8JOPcjNYzojOUn8SfcvxdVOoygcVADEyuU4kZhSEREirz4i6c5vn4OZY9+S7W0k1yf+Sva8GBvuc6UaTGM4BadCFRDaMkDhSERESmSbClJHPn5a2w75xEUvy1jROhkw4ldpVqT1nAgwZ3upUPp0iZXKsWdwpCIiBQdhkHEoa1EbP6MGhE/Uff6eEAW2OdYj0s1+hN013BaqR2Q5COFIRERMV3i5XMcXTsb76PfEJB2iuvTn0ZQnkMVe+LbbhQNGzbRiNBSIBSGRETEFEZaMkd/Xkzqzq8IivuNYIsNgGuGMztL34GlyVCatO9DJ40HJAVMYUhERArVxeM7ObfuEwLP/0AdrqYv/Gti1Kia91H3rgdo4+d3852I5COFIRERKXCpibEcXvsF7vvnUSPlTyr8tTzKKMeBCj2ocMcoGgQ3o4Eeg4kJFIZERKRgGAbn9m3k0s+fUfviGhpyDYBUw5GdpVqR1ngYTe7sTyc3N3PrFLunMJQNTdQqIpI3SdGRHF37Kd6HF1El7TSV/1oeTmXCq/anxl1jaFkt0MwSRTLRRK23oIlaRURywGbl5LYVxG+dTVDsZpxJ/49kouHKTo+OODcfQdO23XDWoIhSSDRRq4iIFIrEy2c5uuoj/I4tItAWlbH8oEMtImsOpH6X0dxRscJN9iBiPoUhERHJHZuNk3+sIH7LpwTF/EKwJf0uUKxRml3luuLVdjTBoXdQ30GNoaV4UBgSEZEcSYqO4Oiqj6h4ZAGBtoj0hRbY71CXS3WHEnz3CDqW9TK3SJE8UBgSEZHsGQZndq0i9udPCIreSOO/2gLFGe7sLNeVsneMJTi0jUaGlmJNYUhERLK4FnuRo6s/xvvP+QRYzxHw1/KDDrWJrD2ERl1H0dHb29QaRfKLwpCIiGSIOPQrUetmUefSahqRCkC84cYOry6UaTuWJs3bqy2QlDgKQyIids6WksSR9XNx2fUZNZL/zJgk9YilOmdrDqFhtwfp4ONjao0iBUlhSETETsVdOE74T+9T7dQS6hIHQLLhxB+l2+PU6iFC29xNHY0LJHZAYUhExJ7YbJza/j0Jv3xM3bhfCbakj7t7wSjPwcr3UbPrv2mr0aHFzigMiYjYgZT4aA7/9F98Ds2lmvV8+kIL7HQK4WrjkTTrMoS7SmmOMLFPCkMiIiVY9OkDnFk5k1oXvqfRXxOlXjVK8Ue57pTv+G+aBDdTt3ixewpDIiIlzV+PwpJ+DqNu/O+U+2vxMapyquZQGncfSyef8qaWKFKUKAyJiJQQ1mtXOfzTJ3jtm00161kAbIaF7a4tSGv+MM079qWWsxpEi/yTwpCISDEXd+EYJ3+cSfUzS6hPIpD+KGx7uZ5UuGscLRoG61GYyE0oDImIFEeGwbk9a4nd8D5BMT/T+K9eYaeoxJHAoTTq+Qh3VtBs8SI5oTAkIlKMGGkpHNv4Fc7bPiQw5SiVASzwh1MTEkLG0PLuQVRzcTa7TJFiRWFIRKQYSE2I5vCPs/A9NIfatksAJBku/O55N14dxxHatKUehYnkkcJQNsLCwggLC8NqtZpdiojYsfjIcMJXvEON09/SkCQALhpe7K40kNo9H6djQFWTKxQp/iyGYRhmF1GUxcXF4eXlRWxsLJ6enmaXIyJ24uLhrUStepegy+twstgAOE4VTtQaRdN7xlK+rJfJFYoUbbm5fuvOkIhIUWGzceq3JaT+/D61kvZQAcACOxwbc7Xpv2jVZSA1XfRnWyS/6bdKRMRkRuo1jq6bTZk/wqiWlj4+UKrhyFb3jri2G0fzVh1xcFB7IJGCojAkImISa1Ich3/8AN/9/6OOcQWAOMOdbeX7ULnreNoH1TW5QhH7oDAkIlLIkmMjOfr9O1Q7No/6JAAQaZRjd8AwGvZ6nM6+FU2uUMS+KAyJiBSS+IjjnPxhGrXOLqUhKQCE48/RWg/SrNfDdPXyMLlCEfukMCQiUsCiw3dxfsU0gi6toiHpPcMOWGpxoeEjtO7xANVLuZpcoYh9UxgSESkgkfs3ErN6GkFxv2bMHL/dMYSE5uNoc1c/GmjSVJEiQWFIRCQ/GQYRu1eRsGYqNRN340v6zPG/ut6BQ7vxtGx7F47qGSZSpCgMiYjkB8Pg/B/fkbx+GtWTDgCQYjjyS5m7KXvXU7Rt0kzTZYgUUQpDIiK3wzA4+9tirBunUy35MADJhjM/e/bAr/uz3Fm/gckFisitKAyJiOSFzcbpLQuw/PwOASknAEg0XPmlbC8Cej5L5zp1TC5QRHJKYUhEJDesaZzcPBeXX2dQNfU0APGGG1u8+1O91zPcXaO6yQWKSG4pDImI5IQ1jfD1syn1+0wC084D6aNF/+ozgDp9nqZrVc0eL1JcKQyJiNyMzcqpTV/ismU61dPOARBtlOE338HU6/M03SpXMrlAEbldCkMiIjdis3H6l/k4/TyNan89DrtilOH3SsNo1PdpuvtVMLlAEckvCkMiIn9nGJz77RuMDVOp+lfD6BijNL/63k+jfs/QvZLmDRMpaRSGREQgfZyg7cuwrnuDgOSjAMQZpdhSYRB1+z5Ljyr+JhcoIgVFYUhE7JthELHzR5LXvEa1a4eA9N5hP3vfS+1+k+heNcDkAkWkoCkMiYjdurhvLfErX6F64j4gfZygTWX7Ur3PJLqri7yI3VAYEhG7E3viDy4ue55acb9TAbhmOLPRszcBvSbRvU5ts8sTkUKmMCQidiPxwhHOLHmBoIur8QJSDUc2lO6BX68X6FavntnliYhJFIZEpMRLiT5P+LcvU/PsEoKwArDBpQOlu02mS5OmmkBVxM4pDIlIiWVLjObYsjeoeuQLgkgB4DfHplxr/yId2t2Jg4NCkIjYQRg6c+YMw4cPJyoqCicnJ1566SUGDBhgdlkiUoCMlERO/PgeFfeEUceIB2APQUQ0f4ZOXfvj4uRgcoUiUpRYDMMwzC6iIF24cIHIyEhCQkKIiIggNDSUI0eOULp06Ry9Py4uDi8vL2JjY/H09CzgakXktljTOLXhf5T5dTrlbZcAOGZU4XCD8XTsPZLSbs4mFygihSU31+8Sf2eoUqVKVKqUPneQn58fPj4+XLlyJcdhSESKAcMgcuf3WFe9SLWUUwCcN8rzR/V/cce9j9HTo5TJBYpIUVbk7xVv3ryZXr164e/vj8ViYdmyZVm2CQsLIzAwEDc3N1q2bMm2bdtuuK8dO3ZgtVoJCNAgaiIlxdVTuzgxswu+3w/HP+UUV4wyLPd9DGPcDnqPnIi3gpCI3EKRD0MJCQkEBwcTFhZ2w/WLFi1iwoQJTJ48mZ07dxIcHEzXrl2JiorKtN2VK1d44IEH+OSTTwqjbBEpYKkx5zj66QhKf96JGnHbSTacWOFxH5dG/0aff71BZZ9yZpcoIsVEsWozZLFYWLp0KX379s1Y1rJlS5o3b86sWbMAsNlsBAQEMG7cOJ577jkAkpOT6dKlC2PHjmX48OE3PUZycjLJyckZr+Pi4ggICFCbIZEiwkiO58R3b+F/4BNKkf67usGpHW7dp9CqqbrJi0g6u2kzlJKSwo4dO5g0aVLGMgcHBzp37szWrVsBMAyDkSNHcuedd94yCAFMnTqVKVOmFFjNIpJHNitnN87G/Zep1LRdBmAPdbjQ6iU6d7kHJ8cif6NbRIqoYv3X49KlS1itVnx9fTMt9/X1JSIiAoAtW7awaNEili1bRkhICCEhIezbty/bfU6aNInY2NiMrzNnzhToZxCRW4vev5rzb7egyuan8bZd5oxRkWW136T6s1vo1q23gpCI3JZifWcoJ+644w5sNluOt3d1dcXV1bUAKxKRnLp24U/Of/0UNaJ/oRwQZ7izvuIDhA54jr4V1SZIRPJHsQ5DPj4+ODo6EhkZmWl5ZGQkfn5+JlUlIrfLuBbLiW9fodrRL6iBlVTDkdXuPanS7xX61qlpdnkiUsIU63vLLi4uhIaGsm7duoxlNpuNdevW0bp1axMrE5E8sdm4sGk2sW8HU/PobJywssUhlC13/0CPZ74kWEFIRApAkb8zFB8fz7FjxzJeh4eHs3v3bry9valatSoTJkxgxIgRNGvWjBYtWvDee++RkJDAqFGjbuu4YWFhhIWFYbVab/cjiEgOxJ/YRvS34wlIOABAuOHHngbP0a3fA7g5O5pcnYiUZEW+a/3GjRvp1KlTluUjRoxgzpw5AMyaNYvp06cTERFBSEgI77//Pi1btsyX42s6DpGCZbsaxcmvnyXwzFIcMIg33Fhd/gGaD3mBgAplzS5PRIqp3Fy/i3wYMpvCkEgBsaZydvX7lN32LmWMBABWO3XEq9cbtAxuaHJxIlLc2c04QyJSPMXsX8O17ydSJTkcgANGdU40e5luPfrirG7yIlLIFIZEpNCkXjnF2YUTqB61FoDLhgfr/B+h4+AnaeClyZNFxBwKQ9lQA2qRfGRN5exPM/DZ/i7VSSbNcGBlqZ4E9H+dgXUCza5OROyc2gzdgtoMidye+GO/cnXxY1S6dhyAndTlwh1v0O3Ou3B00DxiIlIw1GZIRExnJEZz6uuJBJ78hjJAtFGG1f6P0vn+J2nqUcrs8kREMigMiUj+Mgwu/ToXl3UvEmiLBWCl01349JvGoAa1TS5ORCQrhSERyTcpkYeJmv9vqsT+AcAxozK7gyfTq/d9uDpp4EQRKZoUhkTk9qVe4+z3r1Nx73+pQhrXDGeWew2jxdDJ3OerCVVFpGhTGMqGepOJ5EzcgdUkLx9PlZRzAGyxNCGh81sMbNsSi0UNpEWk6FNvsltQbzKRGzMSr3BmwXiqnlkOQIRRjnXVJnDPoEfwKu1icnUiYu/Um0xEClT0H4txWPk0Va3R2AwLy13vIXDAGwytXc3s0kREck1hSERyzLgawdl5jxEQsQZIbyC9M+Q1+vXup2k0RKTYUhgSkVszDK5s/QKXtS8SYLtKquHIEvcBNBn2BgMr+5hdnYjIbVEYEpGbskWf5sK8R6h8aQsAB4xADjWfyr3du+Gku0EiUgIoDInIjdlsXN70Ee6bX6WykUSy4cw3HsNoPWwy9/mpu7yIlBwKQ9lQ13qxZ9aLx4iaN5ZKMTsB2GnU4VTbadzfuRMOmk9MREoYda2/BXWtF7tis3Jp7Uw8fp2GKykkGq58U3Y0nYa9QNUKHmZXJyKSY+paLyK5ZrscTtSXo/CL3QXAVqMRFzu9zQMdWmvwRBEp0fIUhr777rtcv6dLly6UKqWZqkWKHMMgZstnuK57ET8jiXjDja/L/4tuw5+hdTl3s6sTESlweQpDffv2zdX2FouFo0ePUqNGjbwcTkQKSnwUF756iEoRGwD4w6jLmQ4zGNWpje4GiYjdyPNjsoiICCpWrJijbT081NZApKhJ2L0M2/ePU8kaS7LhxPwyI+jwwGSa+XqZXZqISKHKUxgaMWJErh55DRs2TI2PRYqKa7FEfj0e3xNLADhoq8aO0LcYfo/GDRIR+6TeZLeg3mRSkqQc20Ti12MpmxKJ1bCw0KU/DYdOJTjQ1+zSRETylXqTiUhmqde4tPwFfPb/DxfglK0iq+q8wrCBg3B30Z8BEbFvub4nnpSUxLlz57IsP3DgQL4UVFSEhYVRv359mjdvbnYpIrfFem43V2a2xmf//wBYaunMqYGreWjYUAUhERFy+Zhs8eLFjB8/Hh8fH2w2G59++iktW7YEoGnTpuzcubPACjWLHpNJsWWzEbvxP5Te/DpOpHHR8GJhpWcYOvwhvEu7mF2diEiBKrDHZK+//jo7duzA19eXHTt2MGLECJ5//nnuv/9+1PRIpAhJuMTFr0ZT4cImANYaLUjsOoPHWjdUl3kRkX/IVRhKTU3F1ze9oWVoaCibN2+mX79+HDt2TH9gRYqI1GObSVo0igqpl0g2nJldZiw9Rr1ANZ8yZpcmIlIk5arNUMWKFdm7d2/Ga29vb9asWcOhQ4cyLRcRE1jTiFkxGceveuOZeoljNn/mNJjNg0++oSAkInITuWozdPbsWZycnPDz88uybsuWLbRt2zbTsvj4eMqUKd5/hNVmSIqF2LNc/vIByl/eAcBSOuHZfyZ3Na5ucmEiIubIzfU7V3eGqlSpkhGEZs6cmWndP4PQ1atX6dq1a252LyJ5kLL/exLfb0X5yzu4apRipueztBi/QEFIRCSH8tyv9vnnn6d8+fI88MADWdYlJCTQrVs3Ll++fFvFichNpF4j5rtJlN03Gxdgj60G20KnM+6eOzWStIhILuQ5DM2dO5fhw4dTtmxZevfunbE8ISGBrl27cvHiRTZt2pQvRYrIP1w6SvSXwygX9ycAcy29CBwyjbH1KptcmIhI8ZPnMHTfffcRExPDkCFDWLFiBR07dsy4IxQZGcmmTZuoVKlSftYqIkDyH/NgxVOUM5K4bHjwmc8zjBzxEBU93cwuTUSkWLqt4WfHjBnDlStX6NOnD8uXL+fll1/m/PnzbNq0CX9///yq0RRhYWGEhYVhtVrNLkUkXVoysUufxuvAlwD8aqvPn63e5alurXF00NAWIiJ5lS8TtT733HNMnz6dwMBANm7cSEBAQH7UViSoN5kUCbHniPliCGWv7MFmWPjMcQCN73+DlrUqml2ZiEiRVCgTtfbv3z/Ta2dnZ3x8fHjiiScyLV+yZEleDyEigPX4JpIXjKBsWjSxhjsfej/HmNGPUMHD1ezSRERKhDyHIS8vr0yvhwwZctvFiMjfGAZJG2fiuuk13LFxwFaNtY3e4en+nXFWbzERkXyT5zD0+eef52cdIvJ31+KIW/QQnuErAVhma49Drxk80by2yYWJiJQ8t9WA+rpr166xd+9eoqKisNlsGcstFgu9evXKj0OI2I+Lh7n6xWA840+QYjjyvssYuo98ngaVy5pdmYhIiXTbYeinn35i+PDhNxxg0WKxqDeWSC5Y9y0hbemjeNgSuWB487HvZJ4YcT/lSruYXZqISIl12w0Pxo0bx8CBA7lw4QI2my3Tl4KQSA5Z00j8YRKO347C1ZbIVmt9vg39ipceGaEgJCJSwG77zlBkZCQTJkzA19c3P+oRsT/xF4mfN4wyF34DYLbRi8oDpvJY45IzRIWISFF223eG7rvvPjZu3JgPpYjYoYh9JM66gzIXfiPecOMVt2do9+//0lVBSESk0Nz2oIuJiYkMGDCAChUq0KhRI5ydnTOtf/zxx2+rQLNp0EUpKLZDK0j75kFcbEmcsPnxecCbTBzeG08351u/WUREbqpQBl28bsGCBaxevRo3Nzc2btyIxfL/0wJYLJZiH4ZE8p1hkLr5PRw3TMEFg1+sDdjd+j9M6dYMB02rISJS6G77zpCfnx+PP/44zz33HA4OJW8gON0ZknyVlkLS0nGUOrAQgPnWzpTq8w79mlU3uTARkZKlUO8MpaSkMGjQoBIXhDRRq+S7hMskzh2Me8Q2rIaFdxxG0WnUi7SoUd7sykRE7Npt3xl68sknqVChAs8//3x+1VSk6M6Q5IuoP0n64l5KJZwlzijFG6We4V9jHibQp7TZlYmIlEiFemfIarXy9ttvs2rVKho3bpylAfWMGTNu9xAixdvRtaQsGkGptHhO2Sryge/rvDSqP17uaigtIlIU3HYY2rdvH02aNAFg//79mdb9vTG1iN0xDKy/fYRl1fO4YON3W11W1n+bN+9rh4tTyXqsLCJSnN12GNqwYUN+1CFSslhTSfnhaVx2zQHg67QORN85jcmd6uo/CSIiRUy+TNQqIn+TFM21+cNxO/MzNsPCO8b9NBr4IgMb+5tdmYiI3ECe7tXv3bs30+z0t3LgwAHS0tLyciiR4iX2LNc+7ozbmZ9JMFyZ6PQsXce+SXcFIRGRIitPYahJkyY3nKU+O61bt+b06dN5OZRI8RF1iOSP7sQt5hgXDG+e9pzOhHHjCQ4oa3ZlIiJyE3l6TGYYBi+99BLu7u452j4lJSUvhxEpPk5tJfWrgbimxnHUVpkP/N9i+qgelHHVk2gRkaIuT3+p27dvz+HDh3O8fevWrSlVqlReDiVS9B36gbRvRuNsS2aHrTYLak1n+v3tcXVyNLsyERHJgTyFIc1SL5LO+ONzjB8m4ISNNdambAl5m2n9muGoOcZERIoN3cMXyQvDwLbxLRw2vYUFWJjWkQvt3mTy3fXVdV5EpJhRGBLJLZsV6w8TcNw5B4D30/ri0W0yT95Rw9y6REQkTxSGRHIjNYm0b0bjdORHbIaFKdaRNLl3In2bVDa7MhERySOFIZGcSoombd5gnM7+RrLhxETbOPoN/zedgiqaXZmIiNwGhSGRnIg9R+qX/XG+/CdxhjtPWJ7hsTEjCK3mbXZlIiJymwp1tsjff/+9MA8nkj8uHiHt0844X/6TCKMcjzi/znOPjFEQEhEpIQo1DA0YMKAwDydy+yIPkja7O07x5zluq8QT7tOY9u8hBPl5mF2ZiIjkk3x/TDZw4MAbLjcMgytXruT34QpMWFgYYWFhWK1Ws0sRs1zYQ9qcPjglR7PfFsgb3m8QNqYLPmVcza5MRETykcUwDCM/d+jt7c3cuXMpU6ZMpuWGYTBo0CAiIyPz83AFLi4uDi8vL2JjY/H09DS7HCksZ3dg/bIvjilx7LbV5J2Kb/LhmLvwdHM2uzIREcmB3Fy/8/3OUMeOHfHw8KB9+/ZZ1jVu3Di/DyeS/07/jnVufxxT4/nDVocZFd/k4zEd8VAQEhEpkfL9zlBJoztDdubkL1i/GoBjWiK/2erxn4qv88mYDgpCIiLFTG6u37luQJ2UlMS5c+eyLD9w4ECuthEpco5vwDb3XhzTEtlsbcR7vm8oCImI2IFchaHFixdTu3ZtevbsSePGjTN1lR8+fHiOtxEpco6sxjZ/EA7Wa6y3hvBhpdf59MH2CkIiInYgV2Ho9ddfZ8eOHezevZvPP/+cBx98kPnz5+d4Gz2RkyLpzxXYFt6PgzWZ1dZQPq40hU9Ht1UQEhGxE7lqQJ2amoqvry8AoaGhbN68mX79+nHs2LEcbaPZvKXIObAUY/EYHIw0frC25MtKLzL7wTaUcdXg7CIi9iJXd4YqVqzI3r17M157e3uzZs0aDh06lLE8J9uIFAl7v8ZYPBqLkcYS6x18WeklBSERETuUq95kZ8+excnJCT8/vyzrtmzZQtu2bTl79izOzs74+vqyf/9+GjZsmGWb4kS9yUqoXV9hLH8MCwZfp3Xg28oT+Wx0awUhEZESosB6k1WpUuWGQQjICDleXl4sX76cli1bEhIScsNtREy1ax4sfxQLBl+l3cXiys8oCImI2LF8m5ts8+bNjBgxgkqVKvHiiy9SpUoVNZiWoufAUozvHgPg87SuLK/8FLNHt1IQEhGxY7cVhiIiInjrrbeoXbs2PXr0IC0tja+//prz588zZcqU/KpRJH8cWY3x7Rgsho35aZ34sfITzBndUkFIRMTO5fkq0KtXL9atW0enTp145ZVX6Nu3L6VLl85Yr55jUqSE/4yxaDgWWxrLrW1YWHE880a1oLSCkIiI3cvzlWDFihXcf//9jB8/nmbNmuVnTSL56+wfGPMHYbFeY401lPc9JrBodGuNIyQiIsBtPCb79ddfKVWqFHfeeSdBQUG8+uqrHD9+PD9rE7l9EfsxvroXS2oCv1gbMMX1aeaMaYtPGVezKxMRkSIiz2GoVatWfPrpp1y4cIFnn32W1atXU6dOHVq1asUHH3xAZGRkftYpknuXjmHM7YvlWgw7bLV5yvFZ/jfmDgK83c2uTEREipB8nbX+8OHDfPbZZ8ydO5fIyEgsFgtWqzW/dm8KjTNUTMWcxpjdHUvcWQ7YqjHSeJn/PngnzQK9za5MREQKQYHOWn8zQUFBvP3225w9e5YlS5bQs2fP/Ny9SM5cjcD4sg+WuLMcs/kzMu15pg1tpyAkIiI3lK93hkoi3RkqZhKvwJyeEHWQM7YKDEh5mWcG3kn/plXMrkxERAqRaXeGREx1LQ6+uheiDhJplOX+1OcZ0/MOBSEREbkphSEpGVISYcFgOL+Ty4YHQ1Oe554ObRjTrobZlYmISBGnMCTFX1oKfD0cTm3hqlGKB1KeIzS0Nc90DTK7MhERKQYUhqR4Mwz4/gk4tpZEw5WRKc9QuV4r3ujXUKOgi4hIjigMSfG26W3YM580HPh36hM4Bbbm/SFNcHLUj7aIiOSMJmaS4mv3Atj4JgAvpY4iyrc9C0c0w83Z0eTCRESkOLGL/z7369ePcuXKcd9995ldiuSX8M0Y340D4MO03qxz78Hskc3x1HxjIiKSS3YRhp544gm+/PJLs8uQ/BL1JywchsWWyvfWVvyHwXzyQDP8vNzMrkxERIohuwhDHTt2xMPDw+wyJD9cjYR5AyA5lu22Ojyd+ghv3xdCSEBZsysTEZFiqsiHoc2bN9OrVy/8/f2xWCwsW7YsyzZhYWEEBgbi5uZGy5Yt2bZtW+EXKgUvJQEWDILY04QbfjyUMoGxnerTJ6Sy2ZWJiEgxVuTDUEJCAsHBwYSFhd1w/aJFi5gwYQKTJ09m586dBAcH07VrV6Kiogq5UilQNit8OwbO7yIaT0amPEPz+rWZ0KWO2ZWJiEgxV+R7k3Xv3p3u3btnu37GjBmMHTuWUaNGAfDRRx+xYsUKZs+ezXPPPZfr4yUnJ5OcnJzxOi4uLvdFS/5b9Twc/pEUnHkweQLufnWYOSgEBweNJSQiIrenyN8ZupmUlBR27NhB586dM5Y5ODjQuXNntm7dmqd9Tp06FS8vr4yvgICA/CpX8uq3/8LvHwHwRMq/OV26IZ8+EEpp1yKf5UVEpBgo1mHo0qVLWK1WfH19My339fUlIiIi43Xnzp0ZMGAAP/74I1WqVLlpUJo0aRKxsbEZX2fOnCmw+iUHDv0AP00C4I3U+1lnac3Hw0OpUs7d5MJERKSksIv/Wq9duzbH27q6uuLq6lqA1UiOnd2R3k4Ig7lpnfnU2pN3BjQitJq32ZWJiEgJUqzvDPn4+ODo6EhkZGSm5ZGRkfj5+ZlUleSLK+EwfyCkJbHB1oRX0kbwcPua3BdaxezKRESkhCnWYcjFxYXQ0FDWrVuXscxms7Fu3Tpat25tYmVyW67FpgehxEv8SXUeTRlHh7qVeKZbXbMrExGREqjIPyaLj4/n2LFjGa/Dw8PZvXs33t7eVK1alQkTJjBixAiaNWtGixYteO+990hISMjoXZZXYWFhhIWFYbVab/cjSG7YbLDkYbh0hIsOPjyQ+DRVfH34z+AQHNVzTERECoDFMAzD7CJuZuPGjXTq1CnL8hEjRjBnzhwAZs2axfTp04mIiCAkJIT333+fli1b5svx4+Li8PLyIjY2Fk9Pz3zZp9zExmmw8U1SLc70uzaZc6WCWP7oHVQtrwbTIiKSc7m5fhf5MGQ2haFCdGQVzB8EGDyd+jDLjI58NaYlrWqUN7syEREpZnJz/S7WbYakBLl8HL4dCxh8Ze3CYmsHXuhZT0FIREQKnMKQmC85HhYOheRY9jrUY0rqcLo28GVkm0CzKxMRETugMJSNsLAw6tevT/Pmzc0upWQzDPjuMbh4iBjH8jyY+BgVy3rw9r3BWCxqMC0iIgVPbYZuQW2GCtiW92HNS1gtTgy49iJ7LUF880hrmlQtZ3ZlIiJSjKnNkBQPJzbC2skAvJr2ADuNOjzbra6CkIiIFCqFITFHzGn4ZhQYNlY63cUXqXdxZ92KPHhHdbMrExERO6MwJIUvNQkWDYOkK5x2C2J8/HAqeZXi3QHBOGhgRRERKWQKQ1K4DAN+mAAX9nDNpRyDYx4lzcGV94c0oVxpF7OrExERO6QwlA31Jisg2/8He+ZjWBx4OOlRzuPDhC51aB6omehFRMQc6k12C+pNlo9ObYUv7gFbGh+7jmJqbBfa1fbhi1Et9HhMRETylXqTSdETdwG+GQG2NHZ73cXU2M5U9HBl5qAQBSERETGVwpAUPGtqehCKjyTWsw5DIofiYLHwn8FN8CnjanZ1IiJi5xSGpOCtfx3O/I7VxZOBMY+ShBuP31Wb1jU175iIiJhPYUgK1rG1sOU9AN50epTDKRVoU7M84+6sbW5dIiIif1EYkoJzNQKWPAzAb+X78dmVRviUceG9QSE4qp2QiIgUEQpD2VDX+ttks8KSsZB4iatl6zLiXB8sFpg5KISKnm5mVyciIpJBYSgbjz76KAcPHmT79u1ml1I8/TwDwjdjOJdmdMK/ScaF0W2r0652BbMrExERyURhSPLfqV9h45sALKjwBNuv+lDDpzRP3x1kcmEiIiJZKQxJ/kq8At+OAcPGhWp9eP5EQxwsMH1AMKVcHM2uTkREJAuFIck/hgHLH4W4c1jL1WTQuQEAjGlXg9Bq5UwuTkRE5MYUhiT//P4RHP4RHF2ZUXYSp+MdqFmhNBO61DG7MhERkWwpDEn+OL8LVr8EwMHGzxJ2yB0HC7wzIBg3Zz0eExGRokthSG7ftTj4ZhTYUkmp3ZMH9jYC4KH2NWlSVY/HRESkaFMYkttjGPDDkxAdDl4BvGg8wqWEVGpXLMP4zhplWkREij6FoWxo0MUc2vUV7F8MFke2NpnG1/uv4uhg0eMxEREpNiyGYRhmF1GUxcXF4eXlRWxsLJ6enmaXU7RE/QmfdIS0JBLavUj7X0O4nJDCY51q8XRXjSkkIiLmyc31W3eGJG9Sk2DxKEhLgpp38uyFTlxOSCHI14Nxd9UyuzoREZEcUxiSvPlpEkQdhNIVWR00hR/2R+LoYOHdgcG4OunxmIiIFB8KQ5J7B5bBjs8BCzHdw3j2pwgAHu1Ui4aVvUwtTUREJLcUhiR3Yk7D948DYNzxJM/u8iY6MZV6lTx5rJMej4mISPGjMCQ5Z02Db8fCtVio3IzvvUew6kAkTg4W3hnQGBcn/TiJiEjxo6uX5NymaXDmN3D15FK3D3np+yMAjLuzNg389XhMRESKJ4UhyZnwn2Hz9PR/93qP5zdcJTYplQb+nvy7U01zaxMREbkNCkNyawmXYclYwIAmw1htacvqg9cfjwXj7KgfIxERKb50FcuGRqD+i2HA8kfh6gUoX5uEO9/kle8OADC2fQ3qVdJAlCIiUrxpBOpbsPsRqH//BFZOBEcXGLOOqbuc+XjzCSqXLcXaCR0o5aIxhUREpOjRCNSSPyL2weoX0//d5TX+tATy2S/hALzap4GCkIiIlAgKQ3JjKQmweDRYk6FON2zNH+LFpftJsxncXd+Xu+r5ml2hiIhIvlAYkhv76Tm4dATK+EGfD1m88xx/nIrG3cWRyb0bmF2diIhIvlEYkqz2fws7vwQs0P8TruDBmysPAfBk5zpULlvK3PpERETykcKQZBZ9Er4fn/7vdhOgRgem/niImMRU6vp5MLJtoInFiYiI5D+FIfl/1lT4dgwkx0GVFtBxEtvCr/DNjrMAvNGvocYUEhGREkdXNvl/G6fC2e3g6gn3/o9UHHlx2T4ABjcPILSat8kFioiI5D+FIUl3YhP8PCP9373+A+Wq8dkv4RyJjMe7tAvPdqtrbn0iIiIFRGFI4NIx+GYEYEDTB6Bhf85cSeS9tekTsU7qXpdypV3MrVFERKSAKAzZu4RLMO9eSIoG/6bQbRoAU74/wLVUGy2qe3NfaBWTixQRESk4CkP2LDUJFgxO70FWthrcvwhc3Fl9IIK1h6JwcrDwRt+GWCwWsysVEREpMApD2SjxE7XarOk9x85uB7eyMHQxlKlIQnJapolYa/t6mFuniIhIAVMYysajjz7KwYMH2b59u9mlFIzVL8KfP6RPwDpkAVSoA8D7645yPvYaVcqV4vE7a5tcpIiISMFTGLJHv30Ev32Y/u++/4VqbQD4MyKO//01EeuU3pqIVURE7IPCkL059EP6vGMAnV+BRvcBYLMZvLh0P1abQdcGmohVRETsh8KQPTm7I72dEAaEjoK24zNWfbPjzP9PxNpLE7GKiIj9UBiyF1fCYf5ASEuC2ndDj3fgr15i0QkpvLXyTyB9IlZ/TcQqIiJ2RGHIHiRegXn3QeIlqBQM930Ojk4Zq99e9SfRiakE+WoiVhERsT8KQyVd6jVYeD9cPgZeAXD/1+BaJmP1rtPRLNx+BoDX+moiVhERsT+68pVkNhss+xec3gquXjD0G/Dwy1httRm8tHw/hgH9m1amRXVNxCoiIvZHYagkWzcFDiwBB2cY/BVUrJdp9bzfT7H/XBwebk5M6l4vm52IiIiUbE633kSKHcOALf+BLe+lv+4zC6q3z7TJxavJTF91GICJXYOo4OFayEWKiIgUDQpDJY01DX56Frb/L/11pxcheHCWzab+eIir19JoWNmToS2rFXKRIiIiRYfCUEmSHA+LR8PRVYAFur4Jrf6VZbPfT1xmya5zWCzwWp+GODpoIlYREbFfCkMlRdyF9HGEIvaCkxv0/xTq986yWarVxkvL9wMwuHkATaqWK+xKRUREihSFoZIg8gDMGwhxZ8HdB+5fBFWa3XDTOVtOciQynnLuzjzTtW4hFyoiIlL0qDeZWZLj4ZuRcHRtehf4vDq+AWZ3Sw9C5WvDmLXZBqGI2Gu8t/YIAM91r0u50i55P66IiEgJoTtDZtn3DRxYmv5VLjB9rrAmw6F0+ZzvY+dc+GE82NKgWlsY9BW4Zz9W0GsrDpKQYqVp1bIMCA247Y8gIiJSEujOkFmqt4eW/0ofDDH6JKydDDPqwZKH4PTv6d3js2MYsP51+O6x9CDUaAAMX3rTIPTL0Uus2HsBB0v6SNMOajQtIiICKAxlKywsjPr169O8efOCOUD5mtD9LXjqEPT+ACqFgDUZ9i6C2XfDR3fA9s8g+Wrm96UlpwemzdPTX7efmN5Y2in7cYKS06y8/Fej6QdaB9LA36tgPpOIiEgxZDGMm92CkLi4OLy8vIiNjcXT07NgD3ZuB2yfDfsXQ9q19GUuHhA8CJo9CJ6VYOEwOPULWByh13vQ9IFb7jZswzGmrzqMTxlX1j/dAU8354L9HCIiIibLzfVbYegWCjUMXZd4BfYsgD9mp0+wep2rFyTHpgekgV9ArbtuuaszVxLpMnMT11JtzBwUTL8mVQqwcBERkaIhN9dvPSYrity9ofWj8Ngf8MByqNc7/U5Qcix4VoYHV+UoCAG8+sNBrqXaaFndm74hlQu4cBERkeJHvcmKMosFanRM/4o7D8fWQZ2uUKZijt6+/s9I1hyMxMnBwmt9G2KxqNG0iIjIPykMFRee/tB0eI43v5ZqZfJ3BwAYfUd16vh6FFRlIiIixZoek5VQH286wZkrSfh5uvHEXbXNLkdERKTIUhgqgS7HJ/PJ5uMAPN+zHqVddQNQREQkOwpDJdCHG4+TkGKlYWVP7mlUyexyREREijSFoRLmXEwSc387BcDErnU10rSIiMgtKAyVMP9Ze4SUtPSu9O1r+5hdjoiISJGnMFSCHIuKZ/GOswA8062uutKLiIjkgMJQCTJjzWFsBnSuV5HQauXMLkdERKRYUBgqIfadjeXHfRFYLPB01yCzyxERESk2FIZKiLdX/QlA35DK1PUrpDnURERESgCFoRJg6/HL/Hz0Ek4OFp7sXMfsckRERIoVhaFizjCMjLtCQ1pUpWp5d5MrEhERKV4Uhoq5NQcj2XU6BjdnB8bdWcvsckRERIodhaFizGozeGf1YQBGta1ORU83kysSEREpfhSGirHlu89xJDIeTzcnHmlf0+xyREREiiWFoWIqJc3GzLVHAHikY0283J1NrkhERKR4UhgqphZuP82ZK0lU8HBlVJvqZpcjIiJSbCkMFUOJKWm8v+4YAI/fWYtSLo4mVyQiIlJ8KQwVQ59vOcml+GQCvEsxqHlVs8sREREp1hSGipmYxBQ+2nQcgKe6BOHipFMoIiJyO3QlLWY+2nSCq9fSqOvnQe9gf7PLERERKfbsIgz98MMPBAUFUbt2bf73v/+ZXU6eRcVdY86v4QA8fXcQDg4WkysSEREp/pzMLqCgpaWlMWHCBDZs2ICXlxehoaH069eP8uXLm11arr2//ijXUm00rVqWu+pVNLscERGREqHE3xnatm0bDRo0oHLlypQpU4bu3buzevVqs8vKtVOXE1i47QwAz3Sri8Wiu0IiIiL5ociHoc2bN9OrVy/8/f2xWCwsW7YsyzZhYWEEBgbi5uZGy5Yt2bZtW8a68+fPU7ly5YzXlStX5ty5c4VRer4xDIPJ3x0gzWbQvk4FWtUofne1REREiqoiH4YSEhIIDg4mLCzshusXLVrEhAkTmDx5Mjt37iQ4OJiuXbsSFRVVyJUWnO/2nGfj4Yu4ODrw8j31zS5HRESkRCnyYah79+68/vrr9OvX74brZ8yYwdixYxk1ahT169fno48+wt3dndmzZwPg7++f6U7QuXPn8PfPvhdWcnIycXFxmb7MdCUhhSnfHwRg3J21qFWxjKn1iIiIlDRFPgzdTEpKCjt27KBz584ZyxwcHOjcuTNbt24FoEWLFuzfv59z584RHx/PypUr6dq1a7b7nDp1Kl5eXhlfAQEBBf45bub1Hw5yJSGFIF8PHu6gyVhFRETyW7EOQ5cuXcJqteLr65tpua+vLxEREQA4OTnx7rvv0qlTJ0JCQnjqqadu2pNs0qRJxMbGZnydOXOmQD/DzWw6cpElu85hscBb9zbSAIsiIiIFoMR3rQfo3bs3vXv3ztG2rq6uuLq6FnBFt5aQnMbzS/YBMLJNIE2qljO5IhERkZKpWN9q8PHxwdHRkcjIyEzLIyMj8fPzM6mq/DFjzRHOxSRRuWwpnr47yOxyRERESqxiHYZcXFwIDQ1l3bp1GctsNhvr1q2jdevWJlZ2e3afieHzLekjTb/RryGlXe3iBp6IiIgpivxVNj4+nmPHjmW8Dg8PZ/fu3Xh7e1O1alUmTJjAiBEjaNasGS1atOC9994jISGBUaNG3dZxw8LCCAsLw2q13u5HyJVUq43nvt2LzYB+TSrTMUgjTYuIiBQki2EYhtlF3MzGjRvp1KlTluUjRoxgzpw5AMyaNYvp06cTERFBSEgI77//Pi1btsyX48fFxeHl5UVsbCyenp75ss+bCdtwjOmrDuNd2oW1EzrgXdqlwI8pIiJS0uTm+l3kw5DZCjMMHb8YT/f//ExKmo33BoXQt0nlW79JREREssjN9btYtxkqSWw2g0nf7iMlzUaHOhXoE5L9wJAiIiKSfxSGiogF20+z7eQV3F0ceaNfQ03EKiIiUkgUhrIRFhZG/fr1ad68eYEfKyL2Gm/9+CcAT98dRJVy7gV+TBEREUmnNkO3UNBthgzD4KG5O1hzMJLggLIs+VcbHB10V0hEROR2qM1QMfLT/gjWHIzEycHCtHsbKQiJiIgUMoUhE8UmpvLydwcA+HfHmtT1K/iu+yIiIpKZwpCJpq48xMWrydSsUJpH76xldjkiIiJ2SWHIJL8ev8TC7WcAeOvexrg6OZpckYiIiH1SGMpGQfcmi0tKxauUM8NaVaV5oHeBHENERERuTb3JbqEge5NFXb1GKWdHPNyc83W/IiIi9i431+8iP1FrSVbRw83sEkREROyeHpOJiIiIXVMYEhEREbumMCQiIiJ2TWFIRERE7JrCUDYKc6JWERERMY+61t9CQU/UKiIiIvlPE7WKiIiI5JDCkIiIiNg1hSERERGxawpDIiIiYtcUhkRERMSuKQyJiIiIXVMYyobGGRIREbEPGmfoFmJjYylbtixnzpzROEMiIiLFRFxcHAEBAcTExODl5XXTbZ0KqaZi6+rVqwAEBASYXImIiIjk1tWrV28ZhnRn6BZsNhvnz5/Hw8MDi8VidjklxvXErjtu5tO5KDp0LooOnYuiI6/nwjAMrl69ir+/Pw4ON28VpDtDt+Dg4ECVKlXMLqPE8vT01B+aIkLnoujQuSg6dC6Kjryci1vdEbpODahFRETErikMiYiIiF1TGBJTuLq6MnnyZFxdXc0uxe7pXBQdOhdFh85F0VEY50INqEVERMSu6c6QiIiI2DWFIREREbFrCkMiIiJi1xSGRERExK4pDEmBCQsLIzAwEDc3N1q2bMm2bduy3XbOnDlYLJZMX25uboVYbcm1efNmevXqhb+/PxaLhWXLlt3yPRs3bqRp06a4urpSq1Yt5syZU+B1lnS5PQ8bN27M8jthsViIiIgonIJLsKlTp9K8eXM8PDyoWLEiffv25fDhw7d83zfffEPdunVxc3OjUaNG/Pjjj4VQbcmWl3NRENcLhSEpEIsWLWLChAlMnjyZnTt3EhwcTNeuXYmKisr2PZ6enly4cCHj69SpU4VYccmVkJBAcHAwYWFhOdo+PDycnj170qlTJ3bv3s348eMZM2YMq1atKuBKS7bcnofrDh8+nOn3omLFigVUof3YtGkTjz76KL/99htr1qwhNTWVu+++m4SEhGzf8+uvvzJkyBAefPBBdu3aRd++fenbty/79+8vxMpLnrycCyiA64UhUgBatGhhPProoxmvrVar4e/vb0ydOvWG23/++eeGl5dXIVVnvwBj6dKlN93mmWeeMRo0aJBp2aBBg4yuXbsWYGX2JSfnYcOGDQZgREdHF0pN9iwqKsoAjE2bNmW7zcCBA42ePXtmWtayZUvj4YcfLujy7EpOzkVBXC90Z0jyXUpKCjt27KBz584ZyxwcHOjcuTNbt27N9n3x8fFUq1aNgIAA+vTpw4EDBwqjXPmHrVu3Zjp3AF27dr3puZOCExISQqVKlejSpQtbtmwxu5wSKTY2FgBvb+9st9HvReHIybmA/L9eKAxJvrt06RJWqxVfX99My319fbNt7xAUFMTs2bNZvnw5X331FTabjTZt2nD27NnCKFn+JiIi4obnLi4ujqSkJJOqsj+VKlXio48+4ttvv+Xbb78lICCAjh07snPnTrNLK1FsNhvjx4+nbdu2NGzYMNvtsvu9UBuu/JPTc1EQ1wvNWi9FQuvWrWndunXG6zZt2lCvXj0+/vhjXnvtNRMrEzFHUFAQQUFBGa/btGnD8ePHmTlzJnPnzjWxspLl0UcfZf/+/fzyyy9ml2L3cnouCuJ6oTtDku98fHxwdHQkMjIy0/LIyEj8/PxytA9nZ2eaNGnCsWPHCqJEuQk/P78bnjtPT09KlSplUlUC0KJFC/1O5KPHHnuMH374gQ0bNlClSpWbbpvd70VO/6bJzeXmXPxTflwvFIYk37m4uBAaGsq6desyltlsNtatW5cpzd+M1Wpl3759VKpUqaDKlGy0bt0607kDWLNmTY7PnRSc3bt363ciHxiGwWOPPcbSpUtZv3491atXv+V79HtRMPJyLv4pX64X+docW+QvCxcuNFxdXY05c+YYBw8eNB566CGjbNmyRkREhGEYhjF8+HDjueeey9h+ypQpxqpVq4zjx48bO3bsMAYPHmy4ubkZBw4cMOsjlBhXr141du3aZezatcsAjBkzZhi7du0yTp06ZRiGYTz33HPG8OHDM7Y/ceKE4e7ubkycONE4dOiQERYWZjg6Oho//fSTWR+hRMjteZg5c6axbNky4+jRo8a+ffuMJ554wnBwcDDWrl1r1kcoMf71r38ZXl5exsaNG40LFy5kfCUmJmZs88+/UVu2bDGcnJyMd955xzh06JAxefJkw9nZ2di3b58ZH6HEyMu5KIjrhcKQFJgPPvjAqFq1quHi4mK0aNHC+O233zLWdejQwRgxYkTG6/Hjx2ds6+vra/To0cPYuXOnCVWXPNe7aP/z6/r3f8SIEUaHDh2yvCckJMRwcXExatSoYXz++eeFXndJk9vzMG3aNKNmzZqGm5ub4e3tbXTs2NFYv369OcWXMDc6D0Cmn/N//o0yDMP4+uuvjTp16hguLi5GgwYNjBUrVhRu4SVQXs5FQVwvLH8VIyIiImKX1GZIRERE7JrCkIiIiNg1hSERERGxawpDIiIiYtcUhkRERMSuKQyJiIiIXVMYEhEREbumMCQipurYsSPjx4/P8fYbN27EYrFgsVjo27dvnvdjlo4dO2bUv3v3brPLEREUhkSkmDp8+DBz5szJ03urV6/O2rVrM4JVuXLluHbtWqZttm/fnhFa8tOSJUvYtm1bvu5TRG6PwpCIFEsVK1akbNmyuX7f3r17iY6OpkOHDhnLPDw8WLp0aabtPvvsM6pWrXq7ZWbh7e1NhQoV8n2/IpJ3CkMiUqSsWLECLy8v5s2bl+v32mw2nnnmGby9vfHz8+OVV17Jss3y5cvp1q0bzs7OGctGjBjB7NmzM14nJSWxcOFCRowYkem9c+bMoWzZsixbtozatWvj5uZG165dOXPmTKbtvv/+e5o3b46bmxs+Pj7069cv159FRAqPwpCIFBnz589nyJAhzJs3j6FDh+b6/V988QWlS5fm999/5+233+bVV19lzZo1mbb57rvv6NOnT6Zlw4cP5+eff+b06dMAfPvttwQGBtK0adMsx0hMTOSNN97gyy+/ZMuWLcTExDB48OCM9StWrKBfv3706NGDXbt2sW7dOlq0aJHrzyIihcfJ7AJERADCwsJ44YUX+P777zM9wsqNxo0bM3nyZABq167NrFmzWLduHV26dAHg3Llz7N27l+7du2d6X8WKFenevTtz5szh5ZdfZvbs2YwePfqGx0hNTWXWrFm0bNkSSA9g9erVY9u2bbRo0YI33niDwYMHM2XKlIz3BAcH5+nziEjh0J0hETHd4sWLefLJJ1mzZk2egxCkh6G/q1SpElFRURmvv/vuO+64444btjUaPXo0c+bM4cSJE2zdujXbO1NOTk40b94843XdunUpW7Yshw4dAmD37t3cddddef4MIlL4FIZExHRNmjShQoUKzJ49G8Mw8ryfv7cDArBYLNhstozX3333Hb17977he7t3705SUhIPPvggvXr1onz58nmqoVSpUnl6n4iYR2FIRExXs2ZNNmzYwPLlyxk3blyBHCM+Pp4NGzZkaS90nZOTEw888AAbN27M9hEZQFpaGn/88UfG68OHDxMTE0O9evWA9LtT69aty9/iRaRAKQyJSJFQp04dNmzYwLffflsggyf+9NNP1KlTh8DAwGy3ee2117h48SJdu3bNdhtnZ2fGjRvH77//zo4dOxg5ciStWrXKaCQ9efJkFixYwOTJkzl06BD79u1j2rRp+f1xRCQfKQyJSJERFBTE+vXrWbBgAU899VS+7nv58uXZPiK7zsXFBR8fn5sOtOju7s6zzz7L/fffT9u2bSlTpgyLFi3KWN+xY0e++eYbvvvuO0JCQrjzzjs1yKJIEWcxbucBvYhIIdu4cSOdOnUiOjo6x4MupqWl4evry8qVK2+rm/ucOXMYP348MTExed4HwMmTJ6levTq7du0iJCTktvYlIrdPd4ZEpFiqUqUKQ4YMydG2V65c4cknn8zUC8ws3bt3p0GDBmaXISJ/oztDIlKsJCUlce7cOQDKlCmDn59foR07P+4MnTt3jqSkJACqVq2Ki4tLPlUnInmlMCQiIiJ2TY/JRERExK4pDImIiIhdUxgSERERu6YwJCIiInZNYUhERETsmsKQiIiI2DWFIREREbFrCkMiIiJi1xSGRERExK79H9warZwO0uFlAAAAAElFTkSuQmCC", "text/plain": [ "
" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "plt.plot(memo_data[\"ks\"], memo_data[\"T_errs\"], label=\"Thermal Variance\")\n", "plt.plot(memo_data[\"ks\"], memo_data[\"errs\"], label=\"Sample+Thermal Variance\")\n", "plt.yscale(\"log\")\n", "plt.title(\"Original Variance from Pober15\")\n", "plt.ylabel(r\"$\\Delta^2_{21}$ [mK$^2$]\")\n", "plt.xlabel(\"k [h/Mpc]\");" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## HERA" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "P15 used a concept version of HERA that was a perfect hexagon (with no split-core), \n", "with 331 elements and no outriggers. It also used an approximation for the separation\n", "between hexagon rows of 12.12 m, which is close (but not exactly) the same as\n", "$14.0 \\sqrt{3}/2$. Let's create such an array:" ] }, { "cell_type": "code", "execution_count": 4, "metadata": {}, "outputs": [], "source": [ "hera_ants = p21c.antpos.hera(hex_num=11, row_separation=12.12 * un.m)" ] }, { "cell_type": "code", "execution_count": 5, "metadata": {}, "outputs": [], "source": [ "assert hera_ants.shape == (331, 3)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "We will reproduce the results at $\\nu=135$ MHz (the default case in P15). \n", "Like P15, we use a dish size of *7 wavelengths at 150 MHz*, which is close to, but not\n", "exactly, 14 m:" ] }, { "cell_type": "code", "execution_count": 6, "metadata": {}, "outputs": [], "source": [ "beam = p21c.GaussianBeam(\n", " frequency=135 * un.MHz, dish_size=7 * (constants.c / (150 * un.MHz)).to(\"m\")\n", ")" ] }, { "cell_type": "code", "execution_count": 7, "metadata": {}, "outputs": [ { "data": { "text/latex": [ "$13.990315 \\; \\mathrm{m}$" ], "text/plain": [ "" ] }, "execution_count": 7, "metadata": {}, "output_type": "execute_result" } ], "source": [ "beam.dish_size" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Now, let's create the HERA Observatory model. Here, we set the receiver temperature to 100K,\n", "as per Table 1 of P14. We also set the latitude of the instrument to that of Green Bank (38:25:59.24),\n", "which was the assumed location of HERA in P14 & P15. We also set the keyword `beam_crossing_time_incl_latitude` to `False`, which means that the beam crossing time is solely determined by the beam's size, and not the latitude of the telescope (this setting is purely available for backwards compatibility, and its default is `True`):" ] }, { "cell_type": "code", "execution_count": 8, "metadata": {}, "outputs": [], "source": [ "hera = p21c.Observatory(\n", " antpos=hera_ants,\n", " beam=beam,\n", " latitude=0.6707845 * un.rad,\n", " Trcv=100 * un.K,\n", " beam_crossing_time_incl_latitude=False,\n", ")" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Now, set up the observation itself. Here we assume a sky temperature model of \n", "$60{\\rm K} (\\nu/300 {\\rm MHz})^{-2.6}$, which is what P15 used (though there is \n", "no record of it in that memo, and this is slightly different from the sky model\n", "reference in Pober+14). \n", "\n", "Furthermore, we set the number of days of observation to 180, with 6 hours per night, \n", "as described in that memo. \n", "\n", "The coherent observation time in P15 is the \"beam crossing time\", i.e. the size\n", "of the FWHM of the beam, *at 150 MHz*. In the modern version of the code, the observation\n", "duration is taken to be the FWHM at the *observation* frequency. For the sake of this\n", "comparison, we scale the observation duration back to the value at 150 MHz. \n", "\n", "The integration time in 21cmSense plays a negligible role to first order: it controls\n", "how many times the UVWs are phased within the observation duration. The baselines are\n", "gridded in the UVW plane using a delta-function kernel (i.e. nearest-neighbours), and\n", "so changing the integration time merely changes how resolved their evolution across \n", "the UV-plane is. In the original 21cmSense, this value was hard-coded to 1 minute, \n", "so we set this value here.\n", "\n", "**Note:** the phasing of UVWs has received a significant accuracy upgrade in the new\n", "21cmSense, as we have moved to a pyuvdata-backed method. This will cause some small \n", "differences in the final sensitivity (sub-percent level).\n", "\n", "Finally, as per S3.1.2 of P14, we use 8MHz of bandwidth, with 82 channels (providing\n", "~0.1 MHz channel width)." ] }, { "cell_type": "code", "execution_count": 9, "metadata": {}, "outputs": [], "source": [ "obs = p21c.Observation(\n", " observatory=hera,\n", " tsky_amplitude=60 * un.K,\n", " tsky_ref_freq=300 * un.MHz,\n", " spectral_index=2.6,\n", " n_days=180,\n", " time_per_day=6 * un.hour,\n", " bandwidth=8 * un.MHz,\n", " n_channels=82,\n", " integration_time=60 * un.s,\n", " lst_bin_size=beam.at(150 * un.MHz).fwhm.value * 12 / np.pi * 3600 * un.s,\n", " use_approximate_cosmo=True,\n", " cosmo=Planck15.clone(H0=70.0, Om0=0.266),\n", ")" ] }, { "cell_type": "code", "execution_count": 10, "metadata": {}, "outputs": [ { "data": { "text/latex": [ "$2077.3813 \\; \\mathrm{s}$" ], "text/plain": [ "" ] }, "execution_count": 10, "metadata": {}, "output_type": "execute_result" } ], "source": [ "obs.lst_bin_size.to(un.s)" ] }, { "cell_type": "code", "execution_count": 11, "metadata": {}, "outputs": [ { "name": "stderr", "output_type": "stream", "text": [ "finding redundancies: 100%|██████████| 330/330 [00:00<00:00, 1139.98ants/s]\n", "gridding baselines: 100%|██████████| 1260/1260 [00:00<00:00, 7322.08baselines/s]\n" ] }, { "data": { "text/plain": [ "(41, 41)" ] }, "execution_count": 11, "metadata": {}, "output_type": "execute_result" } ], "source": [ "obs.uv_coverage.shape" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Now construct the sensitivity calculation. Here we use the \"moderate\" foreground model,\n", "with a horizon buffer of 0.1 h/Mpc, as mentioned in P15. The `Legacy21cmFAST` theory \n", "model is the same model used in the original code." ] }, { "cell_type": "code", "execution_count": 12, "metadata": {}, "outputs": [], "source": [ "sense_moderate = p21c.PowerSpectrum(\n", " observation=obs,\n", " foreground_model=\"moderate\",\n", " horizon_buffer=0.1 * littleh / un.Mpc,\n", " theory_model=p21c.theory.Legacy21cmFAST(),\n", ")" ] }, { "cell_type": "code", "execution_count": 13, "metadata": {}, "outputs": [ { "name": "stderr", "output_type": "stream", "text": [ "calculating 2D sensitivity: 100%|██████████| 569/569 [00:01<00:00, 405.88uv-bins/s]\n", "averaging to 1D: 100%|██████████| 159/159 [00:01<00:00, 118.37kperp-bins/s]\n", "averaging to 1D: 100%|██████████| 159/159 [00:01<00:00, 124.37kperp-bins/s]\n" ] } ], "source": [ "sense1d = sense_moderate.calculate_sensitivity_1d(thermal=True, sample=True)\n", "sense1d_t = sense_moderate.calculate_sensitivity_1d(thermal=True, sample=False)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Plot the resulting comparison" ] }, { "cell_type": "code", "execution_count": 14, "metadata": {}, "outputs": [ { "name": "stderr", "output_type": "stream", "text": [ "/tmp/ipykernel_8300/2451292649.py:12: RuntimeWarning: invalid value encountered in subtract\n", " ax[1].plot(memo_data['ks'], 100*(sense1d[:len(memo_data['ks'])].value- memo_data['errs'])/memo_data['errs'])\n", "/tmp/ipykernel_8300/2451292649.py:13: RuntimeWarning: invalid value encountered in subtract\n", " ax[1].plot(memo_data['ks'], 100*(sense1d_t[:len(memo_data['ks'])].value- memo_data['T_errs'])/memo_data['T_errs'])\n" ] }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAABLsAAAJjCAYAAADkuxODAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjkuMCwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy80BEi2AAAACXBIWXMAAA9hAAAPYQGoP6dpAADLHElEQVR4nOzdd3xUVfrH8c+dTHomvZIGUqRIb4KFIoqgKCquZUWk2AgosPoTBQW7uytiIZZVaaKIfV0BUZEiiIggIITeCSmQkElvM/f3BzIyBjCUZCjf9+uVl8w5957z3BmM5sk5zzFM0zQRERERERERERE5B1g8HYCIiIiIiIiIiMjpomSXiIiIiIiIiIicM5TsEhERERERERGRc4aSXSIiIiIiIiIics5QsktERERERERERM4ZSnaJiIiIiIiIiMg5Q8kuERERERERERE5Z1g9HcCZyOl0sm/fPmw2G4ZheDocEREREREREZHznmmaFBQUUKdOHSyWY6/fUrLrKPbt20diYqKnwxARERERERERkT/Zs2cPCQkJx+xXsusobDYbcOjNCw4O9nA0IiIiIiIiIiKSn59PYmKiK29zLEp2HcXhrYvBwcFKdomIiIiIiIiInEH+quSUCtSLiIiIiIiIiMg5Q8kuERERERERERE5ZyjZJSIiIiIiIiIi5wzV7DoFDoeDiooKT4chIoK3tzdeXl6eDkNERERERMTjlOw6CaZpkpmZSV5enqdDERFxCQ0NJTY29i+LNYqIiIiIiJzLlOw6CYcTXdHR0QQEBOgHSxHxKNM0KS4uJjs7G4C4uDgPRyQiIiIiIuI5SnadIIfD4Up0RUREeDocEREA/P39AcjOziY6OlpbGkVERERE5LylAvUn6HCNroCAAA9HIiLi7vD3JdUSFBERERGR85lWdp0kbV0UkTONvi+JiIiIiMhhlZWVbE1bRXlJAbF1mxAdU8fTIdUaJbtERERERERERM4Rpmmy9LNUEtel0tjcB0CZaeUnWxeSb5tIXHyyhyOsedrGKCIiIiIiIiJyjlgybRyX/jYGe1B9Nl81jX23fMO6Jg/SoHAVjrevJDtzj6dDrHFKdslptXDhQgzDIC8vz9OhiIiIiIiIiJxXstJ3cvGOSfxS5w5a/OMrGnXuS50mHWl76xOYg78hiCK2fTLe02HWOCW7PCC/tIIMe8lR+zLsJeSX1kxx6bvuugvDMLjvvvuq9KWkpGAYBnfddVeNzF0b1qxZw2233UZiYiL+/v40adKEV155xe2ajIwMbr/9dho1aoTFYmHEiBG1EpvD4eCFF16gcePG+Pv7Ex4eTseOHXnnnXdqZX4RERERERE592379m0qsNL4lqer9EUlNmJTfD8u2j+b8rIyD0RXe87ZZFdeXh7t2rWjVatWXHTRRbz99tueDgk4lOgaMPlnbnnrJ/bluSe89uWVcMtbPzFg8s81lvBKTEzkww8/pKTkj7lLS0v54IMPSEpKqpE5T1R5eflJ3bdy5Uqio6OZMWMG69evZ8yYMTz66KNMmjTJdU1ZWRlRUVGMHTuWli1bnq6Q/9KTTz7JxIkTefrpp0lLS2PBggXcc889WgEnIiIiIiIip401fxfp3okEhYQftd+3bkdsRgn5B7NqObLadc4mu2w2G4sXL2b16tUsX76c5557jpycHE+HRVFZJTmF5ezOLebW//yR8NqXV8Kt//mJ3bnF5BSWU1RWWSPzt2nThsTERD777DNX22effUZSUhKtW7d2u7asrIwHHniA6Oho/Pz8uPTSS1mxYoXbNXPmzKFRo0b4+/vTrVs3du7cWWXOJUuWcNlll+Hv709iYiIPPPAARUVFrv66devy9NNPc+eddxIcHMw999zD1KlTCQ0NZd68eTRp0oSgoCCuvvpqMjIyjvlsgwYN4pVXXqFLly5ccMEF3HHHHQwcONDtWevWrcsrr7zCnXfeSUhIyDHHmjx5Ms2aNcPX15e4uDiGDRvm6jMMg7feeotrr72WgIAAmjRpwrJly9i6dStdu3YlMDCQzp07s23bNtc9X375JUOHDuXmm2+mXr16tGzZksGDB/PQQw+5rnE6nTz//PPUq1cPf39/WrZsySeffOLqP7xFdP78+bRr146AgAA6d+7Mpk2bXNesWbOGbt26YbPZCA4Opm3btvzyyy/V/ixERERERETk7GX6hhJWeeCYi0hK9m+n0rQQFBxWy5HVrnM22eXl5UVAQABwKGljmiamaXo4KogL8efDey4mKTzAlfBauSvXlehKCg/gw3suJi7Ev8ZiGDRoEFOmTHG9njx5MgMHDqxy3f/93//x6aefMm3aNFatWkWDBg3o2bMnubm5AOzZs4cbb7yRPn36sHr1aoYMGcLo0aPdxti2bRtXX301N910E2vXrmXWrFksWbLELXkE8OKLL9KyZUt+/fVXHn/8cQCKi4t58cUXee+991i8eDG7d+92Sw5Vh91uJzz86BntY3njjTdISUnhnnvu4bfffuPLL7+kQYMGbtccTs6tXr2axo0bc/vtt3Pvvffy6KOP8ssvv2CaptszxsbG8v3337N///5jzvv8888zffp03nzzTdavX8/IkSO54447WLRokdt1Y8aMYcKECfzyyy9YrVYGDRrk6vv73/9OQkICK1asYOXKlYwePRpvb2+g+p+FiIiIiIiInJ1iLrmDSPJYNbvq7rbi4kLqbJnJOtsl+AXYPBBdLTLPUIsWLTKvvfZaMy4uzgTMzz//vMo1kyZNMpOTk01fX1+zQ4cO5vLly936Dx48aLZo0cL09/c3J02aVO257Xa7CZh2u71KX0lJiZmWlmaWlJSc8DMdKf1gsXnZP783kx/5yvV12T+/N9MPFp/SuMczYMAA8/rrrzezs7NNX19fc+fOnebOnTtNPz8/c//+/eb1119vDhgwwDRN0ywsLDS9vb3N999/33V/eXm5WadOHfNf//qXaZqm+eijj5pNmzZ1m+ORRx4xAfPgwYOmaZrm4MGDzXvuucftmh9++MG0WCyu9zA5Odns27ev2zVTpkwxAXPr1q2uttTUVDMmJqbaz7t06VLTarWa8+bNO2p/ly5dzAcffLBKe506dcwxY8Ycc1zAHDt2rOv1smXLTMB89913XW0zZ840/fz8XK/Xr19vNmnSxLRYLGbz5s3Ne++915wzZ46rv7S01AwICDB//PFHt7kGDx5s3nbbbaZpmuaCBQtMwPzuu+9c/bNnzzYB13tps9nMqVOnHjXu6nwWcnY7Xd+fRERERETk7PXri9eZxU9Emovff87MO3jQdDqd5rpffjBXP9vFLHkiwty+domnQzxpx8vXHOmMXdlVVFREy5YtSU1NPWr/rFmzGDVqFOPGjWPVqlW0bNmSnj17kp2d7bomNDSUNWvWsGPHDj744AOyss6cPal1Qv2ZeIt7zaiJt7SkTmjNreg6LCoqimuuuYapU6cyZcoUrrnmGiIjI92u2bZtGxUVFVxyySWuNm9vbzp06MCGDRsA2LBhAx07dnS7r1OnTm6v16xZw9SpUwkKCnJ99ezZE6fTyY4dO1zXtWvXrkqcAQEB1K9f3/U6Li7O7fM9nnXr1nH99dczbtw4rrrqqmrdA5Cdnc2+ffu44oorjntdixYtXH+OiYkBoHnz5m5tpaWl5OfnA9C0aVPWrVvHTz/9xKBBg8jOzqZPnz4MGTIEgK1bt1JcXMyVV17p9l5Nnz7dbTvkn+eOi4tzxQ0watQohgwZQo8ePXjhhRfc7q3uZyEiIiIiIiJnr6ZDP2BjeHcu2fRP/CY2IG98Is3+dw11Knaz++qp1Gt+yV8PcpazejqAY+nVqxe9evU6Zv9LL73E3Xff7dp+9+abbzJ79mwmT55cZStdTEwMLVu25IcffqBfv35VxiorK6PsiJMIDicoatK+vBJGzlrj1jZy1ho+vOfiWkl4DRo0yLV97VgJxdOhsLCQe++9lwceeKBK35EF8QMDA6v0H95+d5hhGNXaipqWlsYVV1zBPffcw9ixY08oXn//6r33R8ZmGMYx25xOp6vNYrHQvn172rdvz4gRI5gxYwb9+/dnzJgxFBYWAjB79mzi4+Pd5vL19f3LuQ/PM378eG6//XZmz57N3LlzGTduHB9++CE33HBDtT8LEREREREROXv5+AfS+sFZHNi7md0/foKzrIjAhIu48LJ+RFm9/3qAc8AZm+w6nvLyclauXMmjjz7qarNYLPTo0YNly5YBkJWVRUBAADabDbvdzuLFi7n//vuPOt7zzz/Pk08+WSuxg3sx+qTwACbe0pKRs9a4anjVRsLr6quvpry8HMMw6NmzZ5X++vXr4+Pjw9KlS0lOTgagoqKCFStWMGLECACaNGnCl19+6XbfTz/95Pa6TZs2pKWlVal5VVPWr19P9+7dGTBgAM8+++wJ32+z2ahbty7z58+nW7duNRDhH5o2bQocWsXYtGlTfH192b17N126dDmlcRs1akSjRo0YOXIkt912G1OmTOGGG26o9c9CREREREREPCcyoRGRf3vM02F4xBm7jfF4Dhw4gMPhcG0fOywmJobMzEwAdu3axWWXXUbLli257LLLGD58uNs2syM9+uij2O1219eePXtqLPYMe0mVYvRtk8OrFK3PsJfUWAxwqID/hg0bSEtLw8vLq0p/YGAg999/Pw8//DBff/01aWlp3H333RQXFzN48GAA7rvvPrZs2cLDDz/Mpk2b+OCDD5g6darbOI888gg//vgjw4YNY/Xq1WzZsoX//ve/NVIUfd26dXTr1o2rrrqKUaNGkZmZSWZmZpWi8KtXr2b16tUUFhayf/9+Vq9eTVpamqt//PjxTJgwgVdffZUtW7awatUqXnvttVOKrV+/fkycOJHly5eza9cuFi5cSEpKCo0aNaJx48bYbDYeeughRo4cybRp09i2bZtr3mnTplVrjpKSEoYNG8bChQvZtWsXS5cuZcWKFTRp0gSo3c9CRERERERExFPOypVd1dGhQwdWr15drWt9fX2rbBWrKYG+ViKCfADcVnDVCT10SuOt//mJiCAfAn1r/qMJDg4+bv8LL7yA0+mkf//+FBQU0K5dO+bNm0dY2KEjSpOSkvj0008ZOXIkr732Gh06dOC5555zOx2wRYsWLFq0iDFjxnDZZZdhmib169fnlltuOe3P88knn7B//35mzJjBjBkzXO3Jycns3LnT9bp169auP69cuZIPPvjA7ZoBAwZQWlrKxIkTeeihh4iMjDzq9tcT0bNnT2bOnMnzzz+P3W4nNjaW7t27M378eKzWQ5/1008/TVRUFM8//zzbt28nNDSUNm3a8Nhj1cvEe3l5kZOTw5133klWVhaRkZHceOONrlWLtflZiIiIiIiIiHiKYVanCJKHGYbB559/Tt++fYFD2xgDAgL45JNPXG1wKEmRl5fHf//731OaLz8/n5CQEOx2e5WEUGlpKTt27KBevXr4+fmd3PilFRSVVRIXUnWrYoa9hEBfK8F+58c+WhE5fU7H9ycREREREZEz1fHyNUc6K7cx+vj40LZtW+bPn+9qczqdzJ8/v8ppgCciNTWVpk2b0r59+9MR5jEF+3kfNdEFEBfir0SXiIiIiIiIiMhJOmO3MRYWFrJ161bX6x07drB69WrCw8NJSkpi1KhRDBgwgHbt2tGhQwdefvllioqKXKcznoyUlBRSUlJcmUIRERERERERETm7nLHJrl9++cXtNLxRo0YBh7YqTp06lVtuuYX9+/fzxBNPkJmZSatWrfj666+rFK0XEREREREREZHzxxmb7OratSt/VU5s2LBhOklORERERERERERczthkl4iIiIiIiIiI/MF+cD8bl8+jcuM8/EqyKfOLwK/t7bS6pDf7M3bgHxRKcGiEp8P0OCW7RERERERERETOcLn7MyC1PR0pINuIICugEbEFq4j/fjarl7QhqnwvGdZwGDbnvE94nZWnMdaU2jqNUURERERERETkRGz88DHCzAIAKkwvom6dRJ2xaaxq8QQtylYRTzZBjoOUFOZ5NtAzgJJdR0hJSSEtLY0VK1Z4OhQREREREREROY9UVjpYtWQOSz97nRXffUxJSamrLz//IM0PzGVF1I2kGzHEk41jcm82/TKfmN/ewmKAaULR9ZOJSajvwac4M2gbo4iIiIiIiIiIB/00eyoJK56jDVmutuwlYaxs+iANOvUhff0S2hol1Ok5Ai+/MaRP7k28mQVz+gGQQRRxxn7se9Kg1aWeeowzhlZ2iYiIiIiIiIh4yIp5M+j484NEmAdZ1/lVzMf2se/Wb9kX2p5L08bj825Xopb/CwCfgFBiExtQ0CvVbYycK1859AdHWW2Hf0bSyi4PSc8rYeby3azfZ8fby8IVTaK5rmU8/j5eng5NRERERERERGqB6XQSufyflOGNv1FO2LJnyWpyCXUad8AS8ByF7y4m3CigzPTBaRps+/FTnBdfh21uits4sd/eD0B0ow6eeIwzjlZ2HaG2CtR/9MseuvxrAVN/3ImXxcBeUsHoz37jigkL2ZJVUGPzdu3aleHDhzNixAjCwsKIiYnh7bffpqioiIEDB2Kz2WjQoAFz58513bNu3Tp69epFUFAQMTEx9O/fnwMHDpzSmACLFi2iQ4cO+Pr6EhcXx+jRo6msrKyxZxcRERERERGpTfl5OWxYtZjFb41g9/gLOTguno1PteOHmf9m97Z15OflsD3tF+o5d7Ou9bhDtbjMLByTe7Px529xTLmGIONQ3a6dTe7nt6BONFz/Cua7PYk3s0g3YtjY+xMyiCQSO6WmN7628/sUxsOU7DpCbRSo/3lHLo98upab2yWw/LEreGdAe2bd24kF/+iKzc+bu6asoLTCUWPzT5s2jcjISH7++WeGDx/O/fffz80330znzp1ZtWoVV111Ff3796e4uJi8vDy6d+9O69at+eWXX/j666/Jysrib3/720mPCZCenk7v3r1p3749a9as4Y033uDdd9/lmWeeqbHnFhEREREREakt+Xk5ZL/anYb/vZ62GTM5ENWR7Q0HUREQReeNzxI9vSsZr11N/r5NANRpdRVeg+a4El6N5/RzJbTKTCuGs4Kgqx4jxCwgjgMUmP5sS+zHgZ8/ItAsptK04GdU4Jjcm6y92zz89J6nZFcte/uH7VwYY+PZvs0J9P1jF2ndyEDeuKMN6XklfLU2o8bmb9myJWPHjqVhw4Y8+uij+Pn5ERkZyd13303Dhg154oknyMnJYe3atUyaNInWrVvz3HPP0bhxY1q3bs3kyZNZsGABmzdvPqkxAV5//XUSExOZNGkSjRs3pm/fvjz55JNMmDABp9NZY88uIiIiIiIiUhuK8nNJcuzBajjJJ4g61z1B2zueJuqW18gxQvEzKqhTuRefwGAAMjavPGotrvS2D+NrVOIXHk9UclN2WuthJ5Byw5fLd6fS9MA8NifcyK4bvyTdiKHQKwz/oFAPPPGZRTW7apFpmizatJ+HejbCYjGq9F8QFUTb5DAWbMqmX9uEGomhRYsWrj97eXkRERFB8+bNXW0xMTEAZGdns2bNGhYsWEBQUFCVcbZt20ajRo1OeEyADRs20KlTJwzjj/fgkksuobCwkL1795KUlHQ6HlVERERERESkRpimyfpVS8nd/BOG1Zukdr1IrtfI1b9/1wbiDAf7CSOOA6RP7s3GXqnY5qYQy0HyCcCCSd2W3di6oAH+P7/GnvqtqtTiarpiLHYjkKZdbsbHLwCGf0NJYR4xCfXBNAk3DMJ/vzYrYi5xQaEEh2oro5Jdtcg0ocLpJMDn2G97gI8XFZU1t7rJ29vb7bVhGG5thxNQTqeTwsJC+vTpwz//+c8q48TFxZ3UmCIiIiIiIiJns7Rfl+D1vwe5yLkVp2lgMUwq141nme0KYq8dS0RsIkXbfyKfQBwDvyF96rXEm1kwpx8A6UYMBy59kpY/3MeWrWtw9niaRnP/TuX0y/EzKthHFDuapdBq3fMEGaXkYiN3/z5iExsQHBrxRzLLcF9EE5NQv7bfijOWtjHWIovFoHl8CPM3ZB21315cwYqdubRICKnlyI6uTZs2rF+/nrp169KgQQO3r8DAwJMet0mTJixbtgzTNF1tS5cuxWazkZBQMyvaRERERERERE7Vjo2rafDF9Vzg2MGvrZ7EeDyb0n/s5Lfmo2lZuJjYmT3Y91ovDEclXqaD8Li6VbYmFvRKxRoQBoCX1ZvQhAs5aATjZ1QAUIf9XLJ+PKWGP7kEE06BanGdICW7aln/i5NZsGk//12d7tZe6XAy7st1OE24pf2ZsY0vJSWF3NxcbrvtNlasWMG2bduYN28eAwcOxOE4+SL6Q4cOZc+ePQwfPpyNGzfy3//+l3HjxjFq1CgsFv2VFBERERERkTNTxrevYMWBt+Eges3rZGXsxs8WRlzHmyjGH3+jgghHFkH1OxBolPLLF69W2Zpom5tC4Y/vsJ8wki9sg39QKLnWWNKJZv1lr7Ou00ts7fU+YWO3Uj54gWpxnQRtYzxCamoqqampp5TI+Ss3tUlg2fYcHvxwNZ+s3Ev3xtEUlFby+a/p7MktZuItrYiy+dbY/CeiTp06LF26lEceeYSrrrqKsrIykpOTufrqq08pKRUfH8+cOXN4+OGHadmyJeHh4QwePJixY8eexuhFREREREREqic/L4f07Wlkr/of4ZlLsJgOckMuIqb7UEJCQ/H/vRZW/ZxFrAq7mjj7qkOnJbrV4sqj3PRid+BFtO3Sj01LJ9Bh/bNYDSfpRgwFvVKxzRlKvJlFnfxvWZFwJ1HePoe2JQ6bQ0lhHs3+tBUxNrEBWYNVi+tEGeaRe8kEgPz8fEJCQrDb7QQHB7v1lZaWsmPHDurVq4efn99Jje90mnyxOp33ftrF+vR8fKwWujeOZvCl9WiZGHoankBEzken4/uTiIiIiMj5Jj8vh9xXupDk3Es5VjaHXIrT6kty7o8EmwXYjSD2W+OJGzYHv4mNWNX4H9S99G84Jvc+VIvrd+lGDAd8kzBwEnPH2zjfuYo4DgCwwetC8gLqkliwmgQO3ZNONNYhX6vW1gk4Xr7mSFrZ5QEWi8GNbRK4sY3qU4mIiIiIiIh4UtaerdR37sVimBwklOibXyQ2sQH7tq/He/oVhFNAZWUWJYV5FFgiIfM3YhMfY2OvVFfReQB7z1cJnzeCzNC2+AeFkmGNxFlpIbPpYHx3LSCibA8HwlqR1/gawpY9Q6E1nDhtTawRSnaJiIiIiIiIyDnLNE2yD+Ri4iQ6IhKLxf0Uw+xVX5GMFzmEEcd+t62JQZRQYXqR7t+I1gn1+bnuTbTc/jZrFv+PyAX/cBsn5uu7iSCfkksGum1NbJtQHxjtdm1Ws0u1NbEGKdklIiIiIiIiIucc+8H9rJn9JjHbP+dC56GTDLcaddnX5C4aXnwtgcHhBIdGYMtYyobAdsTclkr64a2Jv6/YSjdi2BPTjcZZXwHQrO9D5LzyORfN74+XYZJBJDtajOKiNc8QYeRTZPoRHJ0MQHBoxDGTWdq6WLN09J2IiIiIiIiInFPseQcoe6U9l299EdMniN/a/5PfOv6b0sA6XJ42nqB3L2Pfa73Iz8vBYjqpxJvYxAYU9Ep1G6egVyqmfwQWnAAUFhzEy1mBl2FimhDHATqvfQyLYZJPIIFGKY4p15C1d5snHlt+p5VdIiIiIiIiInJO2bTsKzpwEABbWSZeLXoQm9iAzIu6c/DdroQZBUQ4DtXhKo5uS7PdM9i2YSW2uSlu49jmDMXH4scu/4toDr/X4orA6bBQceNkyvOz8fL2IbF5F3JzsiiY3JtCrzDV4vIwJbtERERERERE5KyybesmDqZvwS84kiYtOuDl5b5xzUj7kn1GDKZpEm9mudXhCqOActNKht8FtEioj1evFLzfmkrsh70INMpIN2Io6JWKbc5Q4skGJ6xqfB+AWy2u+D9tRYwNsJE1eK5qcZ0BlOw6QmpqKqmpqTgcDk+HIiIiIiIiIiJ/svanb+G78bSoXOdq2/plMvvbP8QFLS7B//dEU2TRNvZEXELy9WOOWodrb0QnEnKXAeBwmhQYgURgp9K0sDOkA5W/zsbX8AHz0Bwxq18hq+3VxCTUVy2us4Bqdh0hJSWFtLQ0VqxY4elQREREREREROQIvy3/lqZz/0aTig2savwwOXctYcuVUyj3j6Lj8gewvtOVjEm9yc/LocLii1fpwWPX4XJUUGH4Aoe2Jh6wxpFFBOtCu9MwfxnNMr+gMKguv138EulGDIVeYfhra+JZQyu7REREREREROTMN/8pTAy8DQcxm6ZTccmtNLzkRjITmlEy+XIiyKe00o+Swjzykq6kxba32LDqB4L/VIcreM79JJoFrE+6nbq4b01sdcTKrMjf/5l1UVdtTTzLaGWXnFYLFy7EMAzy8vI8HYqcRpdffjkffPBBta8fP348rVq1qrmAasBdd91F3759a22+rl27MmLEiGpfP3r0aIYPH15zAYmIiIiIeNDOndtY9vUH/PTtx+Tk5lTpz9i7nWZlv7G60QOkGzHEm1k4Jvdm48/f4pjah0CjFIBdFw4iJqE+jXsPo9Tw44L/9j1Us8uIYWPvT8gggjrsx48yItr0dY0fHBpxzC2Ih7cuytlDya7zyF133YVhGNx3331V+lJSUjAMg7vuuqv2AztN1qxZw2233UZiYiL+/v40adKEV155xe2ajIwMbr/9dho1aoTFYjmhZMOpcDgcvPDCCzRu3Bh/f3/Cw8Pp2LEj77zzTq3Mfyq+/PJLsrKyuPXWW13JzON9LVy40NMh1xiHw8HEiRNp3rw5fn5+hIWF0atXL5YuXVrjcz/00ENMmzaN7du31/hcIiIiIiK1Zcemtfzy/FUkTGlHp5/u5+KlQ/B5pRmLUoeyZ0ca+XmHEl/2jO1YDJPI1tfgNWiOK+HVeE4/VzIrn0CM8gIASkqLKcMbX6MSgDJ8ML7+P+LIwWFa8DJM/L68m6y92zz27FJzlOzyhFI72NOP3mdPP9RfQxITE/nwww8pKSn5I5zSUj744AOSkpJqbN4TUV5eflL3rVy5kujoaGbMmMH69esZM2YMjz76KJMmTXJdU1ZWRlRUFGPHjqVly5anK+S/9OSTTzJx4kSefvpp0tLSWLBgAffcc89ZsQLu1VdfZeDAgVgsFjp37kxGRobr629/+xtXX321W1vnzp1rJA7TNKmsrKyRsas7/6233spTTz3Fgw8+yIYNG1i4cCGJiYl07dqVL774okbnj4yMpGfPnrzxxhs1Oo+IiIiISG3Zu3MzsR9cQevSn1nb4D6KUtYeqsOVfCuXZM8kfGpX9r3Wi/y8HGxh0QDk7tl81Dpc2Zc+RYBZgnfQoc2H/kGh5FmjSSeaX5uNxh7RiuLwpvzW7nky71ykOlznOCW7alupHWbcBFN7g32ve59976H2GTfVWMKrTZs2JCYm8tlnn7naPvvsM5KSkmjdurXbtWVlZTzwwANER0fj5+fHpZdeWqV4/5w5c2jUqBH+/v5069aNnTt3VplzyZIlXHbZZfj7+5OYmMgDDzxAUVGRq79u3bo8/fTT3HnnnQQHB3PPPfcwdepUQkNDmTdvHk2aNCEoKMiVVDmWQYMG8corr9ClSxcuuOAC7rjjDgYOHOj2rHXr1uWVV17hzjvvJCQk5JhjTZ48mWbNmuHr60tcXBzDhg1z9RmGwVtvvcW1115LQEAATZo0YdmyZWzdupWuXbsSGBhI586d2bbtj98QfPnllwwdOpSbb76ZevXq0bJlSwYPHsxDDz3kusbpdPL8889Tr149/P39admyJZ988omr//Cqqvnz59OuXTsCAgLo3LkzmzZtcl2zZs0aunXrhs1mIzg4mLZt2/LLL79U+7P4s/379/P999/Tp08fAHx8fIiNjXV9+fv74+vr69bm4+Pjuv+9996jbt26hISEcOutt1JQUHDCzzt37lzatm2Lr68vS5YsoWvXrgwfPpwRI0YQFhZGTEwMb7/9NkVFRQwcOBCbzUaDBg2YO3euayyHw8HgwYNdc1144YVVVv39lY8++ohPPvmE6dOnM2TIENfn+J///IfrrruOIUOGuN7Lw9s4j/f8R3rqqae46KKLqrS3atWKxx9/3PW6T58+fPjhhycUt4iIiIjImWrX16/gSzlehknM9k8pKK0gom5z6lw5DLsRRKBRRrgjm5LCPOLrN2ebtQE+v7xJ+o6N2P5Uh+uCH0bhxEKjbn8HDm1LjBs2B+uQr2l986O0HjaD1sPep/m1Q4mv3wLr4LnEDZuj7YnnKCW7altZIRTth4M7Yeo1fyS87HsPvT6481B/WWGNhTBo0CCmTJniej158mQGDhxY5br/+7//49NPP2XatGmsWrWKBg0a0LNnT3JzcwHYs2cPN954I3369GH16tUMGTKE0aNHu42xbds2rr76am666SbWrl3LrFmzWLJkiVvyCODFF1+kZcuW/Prrr64f7ouLi3nxxRd57733WLx4Mbt373ZLDlWH3W4nPDz8hO554403SElJ4Z577uG3337jyy+/pEGDBm7XHE7OrV69msaNG3P77bdz77338uijj/LLL79gmqbbM8bGxvL999+zf//+Y877/PPPM336dN58803Wr1/PyJEjueOOO1i0aJHbdWPGjGHChAn88ssvWK1WBg0a5Or7+9//TkJCAitWrGDlypWMHj0ab29voPqfxZGWLFniSuidqG3btvHFF1/w1Vdf8dVXX7Fo0SJeeOGFE37e0aNH88ILL7BhwwZatGgBwLRp04iMjOTnn39m+PDh3H///dx888107tyZVatWcdVVV9G/f3+Ki4uBQ4m1hIQEPv74Y9LS0njiiSd47LHH+Oijj6r9PB988AGNGjVyJf6O9I9//IOcnBy+/fbbaj//kQYNGsSGDRvcksm//vora9eudft3s0OHDuzdu/eoSWURERERkTPJ3r27WbFoNquXL6S07Oi7d5Ky5rMmuGvVGlyTexNBPhWmF7sDWxyqpWUYVHR5jKblvxE29XLX1sUVnd/AbgYSQhGlhi/FRX/8gll1uM5jplRht9tNwLTb7VX6SkpKzLS0NLOkpOTkJ8jbY5ovtzDNccGH/rnrJ/fXeXtOIfpjGzBggHn99deb2dnZpq+vr7lz505z586dpp+fn7l//37z+uuvNwcMGGCapmkWFhaa3t7e5vvvv++6v7y83KxTp475r3/9yzRN03z00UfNpk2bus3xyCOPmIB58OBB0zRNc/DgweY999zjds0PP/xgWiwW13uYnJxs9u3b1+2aKVOmmIC5detWV1tqaqoZExNT7eddunSpabVazXnz5h21v0uXLuaDDz5Ypb1OnTrmmDFjjjkuYI4dO9b1etmyZSZgvvvuu662mTNnmn5+fq7X69evN5s0aWJaLBazefPm5r333mvOmTPH1V9aWmoGBASYP/74o9tcgwcPNm+77TbTNE1zwYIFJmB+9913rv7Zs2ebgOu9tNls5tSpU48ad3U+iz+bOHGiecEFFxzzvTj8d+rPxo0bZwYEBJj5+fmutocfftjs2LHjCT/vF1984XZNly5dzEsvvdT1urKy0gwMDDT79+/vasvIyDABc9myZceMPSUlxbzpppv+8lkOa9y48TH7c3NzTcD85z//Wa3nP/wcR/7969Wrl3n//fe7Xg8fPtzs2rWr2zyHvzctXLjwqHGclu9PIiIiIiKnYPP6Veby53uaZU+EHfoZd1ywmT6uvrlgxvPmvl2bTfvBA65rS56IMH/64GkzY/cWc+/4hq7rzXHB5t7xDc1fX7jSXPtcF9f1mXu2mgfGJbquKXgi2qx8IsQsfSLCzBsX57ovc8/Wo0Qm54Lj5WuOZPVMiu08F5IAd83+YyXX5KsOtYfVPdQeklCj00dFRXHNNdcwdepUTNPkmmuuITIy0u2abdu2UVFRwSWXXOJq8/b2pkOHDmzYsAGADRs20LFjR7f7OnXq5PZ6zZo1rF27lvfff9/VZpomTqeTHTt2uFYMtWvXrkqcAQEB1K//RxY+Li6O7Ozsaj3junXruP766xk3bhxXXXVVte4ByM7OZt++fVxxxRXHve7wCiOAmJgYAJo3b+7WVlpaSn5+PsHBwTRt2pR169axcuVKli5dyuLFi+nTpw933XUX77zzDlu3bqW4uJgrr7zSbZ7y8vIq20uPnDsuLs4Vd1JSEqNGjWLIkCG899579OjRg5tvvtn1Hlb3szhSSUkJfn5+x30vjqVu3brYbDa3WA9/fifyvEf7u3Hke+Dl5UVERESV9x9w+/uSmprK5MmT2b17NyUlJZSXl5/wiZGmaVb72uM9/9HcfffdDBo0iJdeegmLxcIHH3zAxIkT3a7x9/cHcK1YExERERE5k+zY/BtJs66kHg7WXjCE5MvvoPBgFnlL3qHrluexb36VDO9E+H37YJ4lFOf+LcQmNmBjr1SY0881Vv7VrxE872EOBjV0tfkHhZJhjaO00pf9Le/HLCvAGhhB/ctvoSTfTuHk3hR6hRGnOlznPSW7jpCamkpqaioOh6PmJwtJgBv+80eiCw69ruFE12GDBg1ybV9LTU39i6tPXmFhIffeey8PPPBAlb4jC+IHBgZW6T+8/e4wwzCqlWxIS0vjiiuu4J577mHs2LEnFO/hZMJfOTI2wzCO2eZ0Ol1tFouF9u3b0759e0aMGMGMGTPo378/Y8aMobDw0LbV2bNnEx8f7zaXr6/vX859eJ7x48dz++23M3v2bObOncu4ceP48MMPueGGG6r9WRwpMjKSgwcP/sW7cXRH+/wOx3kiz1vdvxvHe18+/PBDHnroISZMmECnTp2w2Wz8+9//Zvny5dV+nkaNGrkSvX92uL1Ro0bHjfHIvw9/1qdPH3x9ffn888/x8fGhoqKCfv36uV1zeAtxVFRUteMWEREREaktmXNeIInKQzW4dv6Xiu73ktzmKnyjLsD+7hJCjCKKKw9QUphHcGgEuxP60Hz3TNav+J7QP9Xgipp7D5Hkkdb+SVdbcGgEDJtDSWEerf60PTEgJIqswXOJCwrV9kRRsutIKSkppKSkkJ+ff9zi5aeFfS98fo972+f31MrKLoCrr76a8vJyDMOgZ8+eVfrr16+Pj48PS5cuJTk5GYCKigpWrFjBiBEjAGjSpAlffvml230//fST2+s2bdqQlpZWpeZVTVm/fj3du3dnwIABPPvssyd8v81mo27dusyfP59u3brVQIR/aNq0KQBFRUU0bdoUX19fdu/eTZcuXU5p3EaNGtGoUSNGjhzJbbfdxpQpU7jhhhtO6rNo3bo1mZmZHDx4kLCwsFOK60in83mrY+nSpXTu3JmhQ4e62o48QKA6br31Vm6//Xb+97//VanbNWHCBCIiIqqsVDsRVquVAQMGMGXKFHx8fLj11lurJF/XrVuHt7c3zZo1O+l5REREREROhr2whOKSIsJDQ/H1rppKqKyooOnBBayKuo46OT8dqqk1uTcbe6Vim5tCCEU4TYMd0d3p/HuiqvF1D1OS+gWNvuqHt+EgnWjSL32OhktGEkkeJaYP4UlN3eYJDo04ZjLrWPW55PyjZJcnHFmMPqzuoRVdn9/zR9H6Wkh4eXl5uVajeHl5VekPDAzk/vvv5+GHHyY8PJykpCT+9a9/UVxczODBgwG47777mDBhAg8//DBDhgxh5cqVTJ061W2cRx55hIsvvphhw4YxZMgQAgMDSUtL49tvv2XSpEmn9ZnWrVtH9+7d6dmzJ6NGjSIzM9P1fEeuhFm9ejVwaHXR/v37Wb16NT4+Pq7k0/jx47nvvvuIjo6mV69eFBQUsHTpUoYPH37SsfXr149LLrmEzp07Exsby44dO3j00Udp1KgRjRs3xmq18tBDDzFy5EicTieXXnopdrudpUuXEhwczIABA/5yjpKSEh5++GH69etHvXr12Lt3LytWrOCmm24CTu6zaN26NZGRkSxdupRrr732pJ//z2w22yk/74lo2LAh06dPZ968edSrV4/33nuPFStWUK9evWqPceutt/Lxxx8zYMAA/v3vf3PFFVeQn59PamoqX375JR9//PFRV6GdiCFDhri2ky5durRK/w8//OA6TVNEREREpKbl5+Wwcfm3OFe/T9vipYQYDrIIJy3uRpIuv4OouCRX4qmoII8Qowiv+t3wuv4J0if3Jt7Mcm1NTDdiKPay4V1ud41fUlaCiQVv49DuqniyiV8yhErTQhF+BBql5E7tQ9bguUpkyQlRsqu22dPdE12HE1tH1vCaeg3cNQdC4v9isFMTHBx83P4XXngBp9NJ//79KSgooF27dsybN8+1wicpKYlPP/2UkSNH8tprr9GhQweee+45t9MBW7RowaJFixgzZgyXXXYZpmlSv359brnlltP+PJ988gn79+9nxowZzJgxw9WenJzsdnrdkTWhVq5cyQcffOB2zYABAygtLWXixIk89NBDREZGVtlOdqJ69uzJzJkzef7557Hb7cTGxtK9e3fGjx+P1XroX8Onn36aqKgonn/+ebZv305oaCht2rThscceq9YcXl5e5OTkcOedd5KVlUVkZCQ33ngjTz55aNnvyXwWXl5eDBw4kPfff/+0Jrvg1J/3RNx77738+uuv3HLLLRiGwW233cbQoUOZO3dutccwDIOPPvqIl19+mYkTJzJ06FD8/Pzo1KkTCxcudKtvd7IaNmxI586dyc3NrVIPDw5txxw/fvwpzyMiIiIi8lfy83LIfeVy2jv3kukVQ1qTEXiH1aF82w902jcN48Op7LLWg+HzCA6NICDIRrlppSRrG7G9BlapwZV31StEzhtKrt8fNXkP1eCKxOHwoqjHvyjL3YPFy5vEdr0pLisnTzW45CQZ5olUXD5PHN7GaLfbqySESktL2bFjB/Xq1Tu5wt2ldphxExTtr7qC6/CKr8AouONT8KvhrZQi1ZCZmUmzZs1YtWqVa0ur1AzTNGnYsCFDhw5l1KhRbn1z587lH//4B2vXrnUlSP/slL8/iYiIiIj8bteWNSTM6IKXYZJONF6D5xKb2IDMPVsx3r2SGHLJNwMouXuJa9XVqon9iLGvpfTWWfjNuuXQyq7f5RBMBPnsuPEr6rW4zNWen5dDSWHeUVduZe3dhr9qcMkRjpevOZKlFmMSOJTAuuPT31du/WmrYkjCoXYluuQMEhsby7vvvsvu3bs9Hco5bf/+/UyaNInMzEwGDhxYpb+oqIgpU6YcM9ElIiIiInI6ZaxbgpdhkkkE8WTjmNybjT9/i2Ny70OJLgIwMQiJiHPdE3vNY4SRR/zMKw/V7CKaVV2mcBAbEeRTbPriHxbnNk9waMQxtyjGJNRXoktOin5q8gS/kGMns2p466LIyejbt6+nQzjnRUdHExkZyX/+85+jHgZwqltpRUREROT8lp+Xg/1gNrs3/Iqx7mP8KvIo8o/H/+K7iE+sS6AtzD2xlPkbe4w6eA/631Hrb2W0eYh2Kx9m774dJNQ/dICSV0AIhQQQbRw60T3UzCdq4RC8cLpqcB2c3Fs1uKTGKdklInIG0I5yEREREakp+Xk5ZL52FUmVu0g0KthlrUdBQBLxhT8R883/KDT92WtNguFz/0h4efsTaBZii01m25/qbxX0SsXcf2jnh59/gKv9UA2uGCoqvcnt8DCOgiwMb3+SO91IaaVTNbik1ijZJSIiIiIiInIOKynMI6FyD35GBfsJxffOj0lOakjmzo0cnNqDMKOAWEc6JYV5rmRXTPsbCN87mSWfvUq9jW+5jWebM5RKazTbrA2pH1fX1R4cGgHD5lBSmEfzo9XgGjyXONXgklqgml0iIiIiIiIi57CS4mICjDJyCSaKPBxTrjlUf2vadYRRQAH+h+pvRdZx3VOv5eVs8mlGp7SnXfW3Nvb+hH1EEU82F1Wuw97kNjAMt7lUg0vOBEp2iYiIiIiIiJyF8vLz+fmHeSxfNJvs/fuPeV3mqjmUm14U3zGbdCOGeDOLxnP6HUpiGTHsuHQiYRSw87cfXfdk7d2GrTwbL+NQuQ0TKJ83nlAzn8MVOGLWvUXW3m01+YgiJ0XbGEVERERERETOIvuz95E283Fa5M6jg1EAQNH3viwKvZqEniOJqpPkvoLKWY4DL+LqXcSWo9TfsgWFA2BWlrvaD9XfiiC90iDv4kep2L4YS0Uhm8IvI7RFL3y+uFv1t+SMpWSXiIiIiIiIyFniYO5+rK+353KzkDVRfSjvMRSrty97f/yY9tsmY5k1h13WejD8a1fCK7LhxfhvfpFlX/2HpNUvuY1nm5vCjpiriDO9SWjc1tV+ZP2tZgn1gUFu92VFqP6WnLm0jfEIqampNG3alPbt23s6FBEREREREZEqNv7wKWEUYhgQlfMzZmA0EfXbEtN1CIUE4GdUEOXIoKQwz3VPg3Y92G1JpO2qx1xbFzf2/sS1pbFTxgw2hlyKLSzabS7V35KzlZJdR0hJSSEtLY0VK1Z4OhQ5jv79+/Pcc88d95qpU6cSGhpaOwGdJuPHj6dVq1a1Nt9dd91F3759q339m2++SZ8+fWouIBERERER+UvWLV+zxxJPOtHEm1k4Jvc+VGx+cm+iOUix6ctBS6RbkiorfTu+zmJ8DAcO02BHUFv2bVzOHt+GmCZ4GSbR+etUf0vOGUp2yVllzZo1zJkzhwceeMDVVrduXV5++WXPBVXLpk2bRvv27QkICMBms9GlSxe++uqrGp930KBBrFq1ih9++KHG5xIREREROR/kl1awJ6eQnxb8j4WvDGbxi7cxf/IT7Nm7mwx7CfmlFVXuCS7NIDO4JV6D5x612PyauH6EO3Pc7vEPCiXPGkUGkayNvp4mhT9x6faJJFVs59fkgewjigJrBP6qvyXnCNXs8pS9K2HFO5CxBqw+0OhqaHsX2GI9HdkZ7bXXXuPmm28mKCjII/NXVFTg7e3tkbkBHnroISZNmsQzzzxD3759qaioYMaMGVx//fW88sorDBs2rMbm9vHx4fbbb+fVV1/lsssuq7F5RERERETOB/mlFTz4xn8ZfvAFLrZsIcsSTZE1jPjd38Lbr/OacRu/RvbhjSHdCPb742eQMp8Q/Iv2EJvYgI1HKTZvLJ1MocVG2BFzHVl/q/URK77q/P6VtXe46m/JOUUruzxh8b/hne6waykkd4LIRrD0VZjUAXb/VGPTdu3aleHDhzNixAjCwsKIiYnh7bffpqioiIEDB2Kz2WjQoAFz5851u2/dunX06tWLoKAgYmJi6N+/PwcOHDjlcRctWkSHDh3w9fUlLi6O0aNHU1lZecz4HQ4Hn3zyidtWuq5du7Jr1y5GjhyJYRgYhuF2z7x582jSpAlBQUFcffXVZGRkuPW/8847NGnSBD8/Pxo3bszrr7/u6tu5cyeGYTBr1iy6dOmCn58f77//vmv733PPPUdMTAyhoaE89dRTVFZW8vDDDxMeHk5CQgJTpkxxm+uRRx6hUaNGBAQEcMEFF/D4449TUVH1NzXH8tNPPzFhwgT+/e9/89BDD9GgQQOaNGnCs88+y4gRIxg1ahR79uwB/tjG+VfPf9j06dOJiIigrKzMrb1v377079/f9bpPnz58+eWXlJSUVDtuERERERGpqijvAK/mDaOlsZXxXsNwDF/DBY/9TOagVcw3OvAQ03n6wAiK83Pd7jMv6sdFFb+x/OuZ2OamuPUFz7mPlnnfkZl8fZX5VH9LzidKdtW2zfPg+2eg66PwwK9wzQS48T8waj3EXgQzb4XS/Bqbftq0aURGRvLzzz8zfPhw7r//fm6++WY6d+7MqlWruOqqq+jfvz/FxcUA5OXl0b17d1q3bs0vv/zC119/TVZWFn/7299Oadz09HR69+5N+/btWbNmDW+88QbvvvsuzzzzzDFjX7t2LXa7nXbt2rnaPvvsMxISEnjqqafIyMhwS+YUFxfz4osv8t5777F48WJ2797NQw895Op///33eeKJJ3j22WfZsGEDzz33HI8//jjTpk1zm3f06NE8+OCDbNiwgZ49ewLw/fffs2/fPhYvXsxLL73EuHHjuPbaawkLC2P58uXcd9993Hvvvezdu9c1js1mY+rUqaSlpfHKK6/w9ttvM3HixGp/djNnziQoKIh77723St8//vEPKioq+PTTT6v9/Ee6+eabcTgcfPnll6627OxsZs+ezaBBf5y60q5dOyorK1m+fHm14xYRERERkaryNi7EZpTgZZgMrPiIB976ipW7cvnHBz/RzLkZgHiyifV1XxDQrMed7PJKpt2y+w9tXSSaNVd9yH5CqcMBvKkkqm1fDzyRyJlDya7atiwVEtpDl0fA4vVHu38Y3PTOoUTX2lk1Nn3Lli0ZO3YsDRs25NFHH8XPz4/IyEjuvvtuGjZsyBNPPEFOTg5r164FYNKkSbRu3ZrnnnuOxo0b07p1ayZPnsyCBQvYvHnzSY/7+uuvk5iYyKRJk2jcuDF9+/blySefZMKECTidzqPGvmvXLry8vIiO/uOEkPDwcLy8vLDZbMTGxhIb+8c20IqKCt58803atWtHmzZtGDZsGPPnz3f1jxs3jgkTJnDjjTdSr149brzxRkaOHMlbb73lNu+IESNc18TFxbnmffXVV7nwwgsZNGgQF154IcXFxTz22GOu98DHx4clS5a4xhk7diydO3embt269OnTh4ceeoiPPvqo2p/d5s2bqV+/Pj4+PlX66tSpQ3BwsNtn8lfPfyR/f39uv/12t9VoM2bMICkpia5du7raAgICCAkJYdeuXdWOW0REREREqjq4cQm5hFAZnESyJZsJxWN49s1pTCgeQ7JlPyXeYfhSQW65+4/tOdnp+DiK8TJMAKLMXBrNu5Mo8igzrVgNJ96f9lexeTmvqWZXbTLNQ1sXr3oW/rTdDoDgOpDcGXYshg5310gILVq0cP3Zy8uLiIgImjdv7mqLiYkBDq3qgUMF4RcsWHDUGlnbtm2jUaNGJzXuhg0b6NSpk9u2w0suuYTCwkL27t1LUlJSlflKSkrw9fWtslXxWAICAqhf/49lunFxca75i4qK2LZtG4MHD+buu/94rysrKwkJCXEb58iVZIc1a9YMi+WP/+jExMRw0UUXVXkPDs8HMGvWLF599VW2bdtGYWEhlZWVBAcHV+tZDjNNs9rXHu/5j+buu++mffv2pKenEx8fz9SpU7nrrruqvN/+/v6uFXoiIiIiIgL5eTls+PVHitd8Tqx9NQAZ4R2o02MYYTZ//I9SD8vLUYrdEkr44DmUvdOL5ILdfOY7HoAyWxLbm6TQ7OdHqChz/39v/6BQMqwROB0WSnpOpCT9N0zTQWTjy7EEx+KYcg2FXmHEqdi8nMeU7KptpgnGcRbUGQZQ/YTGifpzcXXDMNzaDic2Dq+uKiwspE+fPvzzn/+sMtbhVU4nM+7JiIyMpLi4mPLy8qOubvqzo8V0OFlUWFgIwNtvv03Hjh3drvPy8nJ7HRgYWK2xj9Z2+HmXLVvG3//+d5588kl69uxJSEgIH374IRMmTPjL5zisUaNGLFmy5KjPv2/fPvLz813Jx2PFeLxkWevWrWnZsiXTp0/nqquuYv369cyePbvKdbm5uURFRVU7bhERERGRc1l+Xg75L3emg5lJvhHErqhuYDpof+BLfD78mHwjkAxrHAyb45bwssY2JTH7M1bvzec/5UN5ndGuvpHlQ/nb1l/Iw0ZEdILbfEcWm2+QUB+4xq0/a/BcFZuX8562MdYmw4CkTrD+86P3F2TBrh8hqXPtxnUcbdq0Yf369dStW5cGDRq4fR0tCVRdTZo0YdmyZW7Jl6VLl2Kz2UhISDjqPa1atQIgLS3Nrd3HxweHw3FC88fExFCnTh22b99e5bnq1at3Yg9TDT/++CPJycmMGTOGdu3a0bBhwxPeCnjrrbdSWFhYZZslwIsvvoi3tzc33XTTKcU5ZMgQpk6dypQpU+jRoweJiYlu/du2baO0tJTWrVuf0jwiIiIiIueKLasWEm9mYhhQSADR142nxbCZFN7xNQ68iMJOcOUBSgrz3O5retUgSg1fcmel8EiJ+y/Bx5b8k/Y5X5IWez1WH98qc6rYvMjxKdlV2y6+D3b/eOj0xSNX2ZQVwhf3gXcAtLrdc/H9SUpKCrm5udx2222sWLGCbdu2MW/ePAYOHHjCCaYjDR06lD179jB8+HA2btzIf//7X8aNG8eoUaPctgceKSoqijZt2rjVwQKoW7cuixcvJj093e2UyL/y5JNP8vzzz/Pqq6+yefNmfvvtN6ZMmcJLL7100s91LA0bNmT37t18+OGHbNu2jVdffZXPPz9G0vMYOnXqxIMPPsjDDz/MhAkT2LZtGxs3bmTs2LG88sorTJgwoUpy6kTdfvvt7N27l7ffftutMP1hP/zwAxdccIHb9kgRERERkfNZ2ZrPyDHCSCeaeLJxTO7Nxp+/xfl+PwKMMhymwY6wS6okp/Kcfrxm6U83YxXJlv0UWUP4tdEDFBsB1LEcxIqD1NyOZNh1ErrIiVKyq7Y16QOXjoJvH4fXL4ZvHof/PQgTm8Hu5XDLDPAP9XSULnXq1GHp0qU4HA6uuuoqmjdvzogRIwgNDT1mUqo64uPjmTNnDj///DMtW7bkvvvuY/DgwYwdO/a49w0ZMoT333/fre2pp55i586d1K9f/4S21w0ZMoR33nmHKVOm0Lx5c7p06cLUqVNrZGXXddddx8iRIxk2bBitWrXixx9/5PHHHz/hcV5++WVef/11Zs6cyUUXXUS7du1YvHgxX3zxBcOHDz/lOENCQrjpppsICgqib9++VfpnzpzpVuNMREREROR8l5y/gu0xPfEaPJd0I4Z4M4vGc/odOinRiGGN7XIi89Oq3BdUls1AvsQwoMzwxb8in9abX6UcK6VeNnyNSiZUPktQ2bHr7orI0RnmiVS8Pk/k5+cTEhKC3W6vUkC8tLSUHTt2UK9ePfz8/E5+kh0/wIp3IGM1ePnChVdDu8EQlnxqwZ/jSkpKuPDCC5k1axadOnXydDjnpCuuuIJmzZrx6quvurWvX7+e7t27s3nz5ipF/OXMcNq+P4mIiIich/LzctiatpLcVV8Sf2ApVirYF9gM/873kJSUTIAt7KjbA/c92Yid0VfQ+f432PjztzSe08/Vt7H3J+T9OIXIwk00GLvS/cZSO5XTb4Si/VgHzQH/MMzyYoyAcCjIoHJybwiMwnrnZ+Cn//8WgePna46kAvWeUu+yQ19yQvz9/Zk+ffoJbVeU6jl48CALFy5k4cKFvP7661X6MzIymD59uhJdIiIiInLOyc/LIeeVLrRy7qXE8GV75BVUeAfQJGsR0d/eSJ4ZxD7vhCpF5gEyQ1uRlDWf9B0bsc1NceuzzRlKjFnMtjj3IvIA+IUcSmSVFUJIPACGz+91kUMSsA6aC75BSnSJnAQlu+Ss07VrV0+HcE5q3bo1Bw8e5J///CcXXnhhlf4ePXp4ICoRERERkZqXuXsLDZx7sRgmdmxEXf8UsYkNyNy1ifwp3Qg1CqmozKSkMK9Ksiui23ASPruW/KndCDaKSTdiKOiVim3OUOLJxgRsra4/+sR+IcdOZv2eABORE6dkl4gAsHPnTk+HICIiIiLiEVmr5lAXL3IIow77SZ/cm429UrHNTSGYIspNKxk+ybQ4ygmIfuF1OIiNMKOACtPC1uD2OJfNItwsxwQMA4LmjSQrofExT1AUkdNLBepFRERERETkvGbL+olNAW0wBn991CLzvyb8nfrlG496r39QKNnWeLIJZ1PIZVxY+AtN7D+QGdOFbT2nkW7EUOgVhn9QaO0+lMh5TCu7TpLq+ovImUbfl0REROR8l5+Xw8bffsG+dg7+BTuotAZCsxu4+IobsGfvxj8o9KhF5i2YOEyITWzAxl6pcESR+YJeqTg2/oBxjDmDQyNg2BxKCvO46IiVW7G//zMrqTlxx5hXRGqGkl1HSE1NJTU1FYfDccxrvL29ASguLsbf37+2QhMR+UvFxcXAH9+nRERERM4n9rwDFL3ckfbmAUoMP9IDmhBYuo06y+ew4ed/EmLaybBGHLXIfEX8xTTb8gabVi87apF5qyWA7QEtuegYcweHRhwzmaWtiyK1T8muI6SkpJCSkuI6yvJovLy8CA0NJTs7G4CAgAAM41g5fhGRmmeaJsXFxWRnZxMaGoqXl5enQxIRERGpdb/OfpuuHAAD8rBhu+0dYhPqs+m7KTRaMgovwyTdYTlqkflGve7H3PoWSZ9fh79RXqXIPE5Ye9H9HnoyETlRSnadhNjYQwtSDye8RETOBKGhoa7vTyIiIiLnE6fDSfy2WWz0bYmtLIN4sv8oMv/jC3gZh8o97O04no5HWWlVXFKCBX+ijDwcpsEO/2aULfuIRkdcE/HLS2S1vFIrtUTOAkp2nQTDMIiLiyM6OpqKigpPhyMigre3t1Z0iYiIyHkrc+92Gjp3sLbta3g16Uz65N7Em1mu2lvpRBNEMY49K4Dbq9zvHxRKhjUWR6UXWbGXUTfnV7zKN5IT3pLSFjcQsOhJCr3CiFOReZGzgpJdp8DLy0s/XIqIiIiIiJxGTqfJr2t/JT99Mz62cFp37EqAr89x7zErSw79wT/06EXme7+OOXcYFmfZUe8/ssh8yyNWbsX9/s+sBheryLzIWUTJLhERERERETkjLF80F/9F42nr3Ohq2zU/lk1Nh9Oi89UE2MKOmnCKTqiPnSDsv80ls27TKkXmQ+fcSyw55CS2OubcKjIvcu6weDoAERERERERkVVLv6Ht97fT2LGFtJaP4Ry+muwbP6UkrBFXbRiDz7tdyZjUm/y8nCr3evsGsDW+L22zPsV490rizSzSjRg29v6EdKKJJQeHaRB5YScPPJmI1DYlu0RERERERMTjvBY+i4mBj+EgZO27ZBc7iG7Rg9AbX6IQf8IpwFaZQ0lh3lHvj+t2N144iCGXQtOPtMTb2b36e0zAaXKoSP2Mm8jau602H0tEPEDJLhEREREREfGojPRdNC9fw5pGw0g3Yog3s3BM7s3Gn7/FMeUagijBNGH7BXccc0thUGQiO60XYCeQMksAV+6eSNd9/8Ee3optV00l3Yih0CsMfxWZFznnqWaXiIiIiIiIeFRh9i7iDJPQi3ridfmdVU9TNGLwpwRrZdExxwgOjYDhX1NSmEdM/AVQWYqPlw/NLIcOFctKbqEi8yLnCSW7RERERERE5JTl5+WwfdsmMlZ/g3/2akzDi7LkrrS5+k7MogP4HyfRFBoZC0D2rjQatLykymmKmZc8RfMfhuJtiz5uDG5F5r393fpUZF7k/GGYpml6OogzTX5+PiEhIdjtdoKDgz0djoiIiIiIyBktPy+H3FcuI8m5D6dhsMe/GV7OMpLKNpNphoNhYLdGETdszjETXpuf70xlZSUBt07B+4MbDq3s+p3dDMSfMsofXEdQeFxtPZaInGGqm69RzS4RERERERE5JZt+/YEkZzoWw2Q/4fjfPpWkR1ew69qZhGMnlpzjFpcH8LlyHI0qtxAzo6vrNMVFbV/FbgYQYhRRavhSWHTsbYwiIocp2SUiIiIiIiKnpGz1Jxw0QkknmjgOuIrLW2ePwMdwYJqwOeGm424l9I9tQJ4RjL9RDkCkM4cuKx/A16jETiDBFOGY3FunKYrIX1LNLhERERERETklF9qXsD3+OhKvfvCoxeXzfaMJ27/iuGP4B4WSYY2lrNKHok4P4yzNxycojORON5GXc4DCyb0p9AojTqcpishfULJLRERERERETokv5ZR7BxOb2KBKcfmCXqkULn6dwPIDxx0jODQChs2hpDCPRn9aARYbEErW4Lk6TVFEqkXJLhERERERkfNYfl4OGRn72PrrQsK3fk6QI4887xiKmt1K45YdCA+P/ssE0z7f+tj2LiJz9xZsc1Pc+mxzhhJtlrAjrvdfxuJ2muKf6DRFEakuJbtERERERETOU/l5OWS8ehV1HTu50Khke2AryoObUzd3PYlrHqB4tQ+7rXVh+NfHTXiZ7QbTYumD5EzuTiz5pBsxFPRKxTZnKPFkgwE5ra6vvQcTkfOaCtSLiIiIiIicp0oK80h07MbXqCSbcAJufZfG907D++8zOUAIAUY5dSr3HvcURYDwxpdShB8R5FNs+vBbnb+xY/1y7Nhc1wTMG6ni8iJSK5TsEhEREREROU/Z8wsJMMrJJZhocl2nKDqmXEMkdvIJwIIT78Cw447jbwtjr1cyB7GR5x3L1emv0WvXv/EJCmXTxf8k3Yih0CsMfxWXF5FaoG2MIiIiIiIi56nc9d/hMA1K7pxN+ox+VU5RLOjxbxp/eye/rf+R8EuvPeY4waERMHwuJYV51EmoDxWlYFhoYPUBIOuiK1VcXkRqjVZ2iYiIiIiInKe8cGJiEF6nEQW9Ut36Cnql4h1R99B1hvMvxwoOjfijiLy3H/ye6IJDxeWV6BKR2qKVXSIiIiIiImeh/Lwc9mbsY9uvi/Hf9T0WZwUFEc1pcOU9RPpW4l+NlVSJLbthXf8MSz5/lYZb3nHrs81NYWtkD+Lxpt5FnWvyUURETivDNE3T00GcafLz8wkJCcFutxMcHOzpcERERERERNzk5+WQ/eoVJDt242042OvbgEprAPFFaVSYFooMf3KtccQNm/OXCa9tz3UkuWwzVsP5xymKc1OIN7Nwmgbrwq6gxYhPa+nJRESOrbr5Gm1jFBEREREROcvsz9xLsmMX3oaDLMKx3jGLug//QM6ts8EwiMJOaGX2X56imLV3G4HlB7AaTpwm7LBewLaf55BuxGGaYDFMIvPW6BRFETmrnLPJrj179tC1a1eaNm1KixYt+Pjjjz0dkoiIiIiIyGmxd+0CvA0nWUQQc+QpirPuJIAyKk0LO/2a/FFD6xj8g0KxWyPZRxQb6w+mnpHJJbmfE+tTysZmI0gnmgJrhE5RFJGzyjm7jTEjI4OsrCxatWpFZmYmbdu2ZfPmzQQGBv7lvdrGKCIiIiIiZ7KV/7qGoMo8QvpPwzG596FTFH+XbsSwL7Y7jfd9ge3JzL8cKz8vh5LCvKMmxrL2bqtW7S8Rkdpw3m9jjIuLo1WrVgDExsYSGRlJbm6uZ4MSERERERE5DbzNcoq8bMQmNjjqKYoVQYn4UlGtsdxOUfwTnaIoImejMzbZtXjxYvr06UOdOnUwDIMvvviiyjWpqanUrVsXPz8/OnbsyM8//3zUsVauXInD4SAxMbGGoxYREREREflrm/ZkMfebOXz97Tx2ZR884fsro5pRv3gtW9J+xTY3xa3PNjeFgB1fs9fngtMVrojIWeWMTXYVFRXRsmVLUlNTj9o/a9YsRo0axbhx41i1ahUtW7akZ8+eZGdnu12Xm5vLnXfeyX/+859jzlVWVkZ+fr7bl4iIiIiIyOm2bfs25r1wK3HvtKLXj7dx9dK/EZDagi9eGcG2zWnk5+VUa5wLeg4liGJiZ/Ui3swi3YhhY+9PSDdiiDezaFmxloIGfWr4aUREzkxnRc0uwzD4/PPP6du3r6utY8eOtG/fnkmTJgHgdDpJTExk+PDhjB49GjiUxLryyiu5++676d+//zHHHz9+PE8++WSVdtXsEhERERGR0yU9fTfB/+lAICVsSbqF+lfdg6OynD0Lp1Nv1yxKTW/2WOtSZ/jcv9w6mLV3Gz7vdCGMAspMK99FDcARHE/Uztlc7FiJYUA60ViHfP2XRepFRM4W53TNrvLyclauXEmPHj1cbRaLhR49erBs2TIATNPkrrvuonv37sdNdAE8+uij2O1219eePXtqNH4RERERETn/bPrmbYIowWJA0N6FHCAU33qdCerxEHkEEWiUEeHIoqQw7y/H8g8KJdsazwFCSQ9qRq8Dk7lu+1PUtR5gc7MHSSeaQmu4TlEUkfOS1dMBnIwDBw7gcDiIiYlxa4+JiWHjxo0ALF26lFmzZtGiRQtXva/33nuP5s2bVxnP19cXX1/fGo9bRERERETOXzF7v2Wzf0uCyjIPbT2c3JuNvVKxzU0hlnzKTS92+11I22qsxAoOjYBhcygpzOOChPrgqABnJXHe/sQBWXsHEKdTFEXkPHVWJruq49JLL8XpdHo6DBEREREREQBiHBlsif8bIVfeS/rk3sSbWTCnHwDpRgwZtmaElFR/l0lwaMQfySwv70Nfh+fS1kUROY+dldsYIyMj8fLyIisry609KyuL2NhYD0UlIiIiIiJybCVeNpy5u4hNbEBBL/eDuPJ7TcKrNJdy7xAPRScicu44K1d2+fj40LZtW+bPn+8qWu90Opk/fz7Dhg076XFTU1NJTU3F4XCcpkhFRERERORsll9awfbsfLb9NBtnVhp4+xPZ+joua9eS/YVlBPpaCfbz/uuBgNx619F6yzssX/odCd+luPWFz7mXGHJJa/1CTTyGiMh55Yw9jbGwsJCtW7cC0Lp1a1566SW6detGeHg4SUlJzJo1iwEDBvDWW2/RoUMHXn75ZT766CM2btxYpZbXiapudX8RERERETl35ZdWMO6l1xhVmkqi5QClhh9WswJMk298ujPZuAHf4GjeGNKtWgmv0rwsil7pSLAzH2/DwV6i2X7pv2m65AEisVNuWskZsJC4C6rWGRYRkerna87YlV2//PIL3bp1c70eNWoUAAMGDGDq1Knccsst7N+/nyeeeILMzExatWrF119/fcqJLhEREREREYCVS+bxYtlTVBhWRnuP5oGhD1InwMGmr9/kilX/pAs/sCMnmeL8rwj2++ufQ+yFhZimBW/j0E6SELOANj/cS5BRSine+BkVON+7iazBc1VzS0TkFJyxya6uXbvyV4vOhg0bdkrbFkVERERERI4l6JdJmIYFPyq4v2wyD7xVn0dv68HzG5J50/Qn0pLPBb52AnwrqzWef1AoGdYoHJVeOK98isr9WzEMLyytepJvCcExuTeFXmHEBYXW7IOJiJzjzthkl4iIiIiIiKcUFxfRsuRnNjUYxIUHviHZvosJxWMY+WY6E71fJ9KSjwML6wM60j4kvlpjBodGwLA5lBTmVVm5FQBkDZ5LXFDoHycsiojISTkrT2OsKampqTRt2pT27dt7OhQREREREfGgypJCfAwHJZEXYR00hzJbEsmWbD7zHU+yJZsyWxI7rBeAs+KExg0OjTjmFsWYhPpKdImInAZKdh0hJSWFtLQ0VqxY4elQRERERETEg2yh4RQQQP7W5ewzIxhZPtSt/6HSgcRUpOMdUdczAYqIyDFpG6OIiIiIiJwz9uQUsWzx11TaM7AGx3Lx5VeTFBl0wuMYXt7sSLie1nu+ZNRrFzG+4mW3pQJPl/+bAEoJvnjA6QteREROC8P8qyrw56HqHmUpIiIiIiJnBqfT5OPpr3Lx9kkkW7Jd7TucsfzUcCRdu/Yg0BZ2QtsE9+7dheXtbkRzEKvhpMAnhlXJA2m/5WUCKCXPDGCI/yu8dn8f4kL8a+KxRETkCNXN12gbo4iIiIiInPVmf/g6f9vxBDGWPDKvfAP+bwdld3yFV1R9btk6Gu93u5MxqTf5eTnVHjPEB7wsFqyGExOwlWfRZcsLlFltlPuEEWoUM6n8cYLKsv9yLBERqT3axigiIiIiIme1krIKWm9+mVJ88Kccx3fjyUzqQGyDy/C5PoayyZcQiZ2ySl9KCvOqvbrLFhyKf1wylUUBWG+ZTmVZIV4+gYTFtYCCDCon9yYyMAprcGjNPqCIiJwQJbuOkJqaSmpqKg6Hw9OhiIiIiIhINW34+TvakM3GLq9h++EZ4s0s0if3ZmOvVGxzU/CnHICdrR7mkmOchHhUfiFY7/wMygohJN79h6eQBKyD5oJvEPiFnNbnERGRU6NtjEfQaYwiIiIiImehggwA4tpcg9egOaQbMcSbWTSe0494M4t9RAPgazmJX2r7hUBI/NH7QuKV6BIROQMp2SUiIiIiIme1yNgEANJW/0xsYgMKeqW69W+68D4AImITaz02ERGpfUp2iYiIiIhIrTNNk5zCMg4UlnGqB8QntbqCLEsM5pKX2LplA7a5KW79HTf+kxwjnHptrz6leURE5Oygml0iIiIiIlJr7HkHmP/DUkrWfUWjktUYwHe+zfDuOITOF8YSaAurdgF5F4sXjh5P0vmb+yia0YVAo4y9RLM8cQhX755AoFFGielL5r6dxCY2qJHnEhGRM4dhnuqvUc5B+fn5hISEYLfbCQ4O9nQ4IiIiIiLnBHveAXJfvoxkcx/lFj9y4ruD4UVY+gK8HSXYCeKAdxx1hs054YRX1t5tWN/pRgR29zm9wql0VBJBPulGDNbBc4k5kSL1IiJyxqhuvkYru0REREREpFasX/0zF5v7sBiQa9rw6vk0sYkNyNzxG7ZpPYjETmWllZLCvBNOdvkHhZJhjaO00hdrj7H4Wpz4h8cT0rA7mft2kj65N4VeYcQFhdbMw4mIyBlDya4jpKamkpqaisNxEqe0iIiIiIjIcRWt/R8lhj92gqjDftIn92Zjr1Rsc1MIpJRK08KWgJZcdhIrr4JDI2DYHEoK86qs3IpNbEDW4LnEBYWe+BZJERE562gb41FoG6OIiIiIyOm3+el2FAXXJ+7GZ3FM7k28meXqSzdiyAxrR/TBVSSO2+jBKEVE5ExV3XyNTmMUEREREZFaYTWcFDksxCY2oKBXqltfQa9UCi02rDg9FJ2IiJwrtI1RRERERERqRUlMOy7c+zVp69YQMjfFrc82Zyg200leVDviPBSfiIicG7SyS0REREREqjBNk105RWzOKqC4vPK0jJnU8wHCjXwSP+5JvJlFuhHDxt6fsJdo4skm3jiAX8ubTstcIiJy/jqplV1ffvnlCd9z5ZVX4u/vfzLTiYiIiIhILcnPy+G7ZSvZv3oOrUqW4UsFs4265DT+Oz07XERERORJF3kvtvhTgY1wI58K08JyrzaUfPcRlzpx/Rre5/uxZNVrW6XIvIiISHWdVLKrb9++J3S9YRhs2bKFCy644GSmExERERGRWpCfl0PWy1253tyDaXiRm9gNIyCMC3YtInjTIPI2BrLPOxGGzTmphJd/UCgZ1jpUVlopi23DFbm/YJgmBYntyG5xExVfj6XQK4y4oNDT/3AiInLeOOmaXZmZmURHR1frWpvNdrLTiIiIiIhILdm1exdNzL14GSYZhGJc/QKxiQ3I3LURc0p3Qo0iKiozKCnMO6lkV3BoBAybQ0lhHolHrNw6fJ5WVnxb4oJCT3rlmIiICJxkza4BAwac0JbEO+6447hHQp4pUlNTadq0Ke3bt/d0KCIiIiIitW7v2gVYDSeZRBDHARyTe7Px529xTL2OEIooNb3JJJro+JPfsREcGnHMLYoxCfWV6BIRkVNmmKZpejqIM01+fj4hISHY7fazIkknIiIiInI6LJ/Qj4jS3QTd+QGOyb2JN7NcfelGDBkX/I12216jdHQWfn5+HoxURETOR9XN1+g0RhERERERAcDXC0ocFmIS6lPQK9Wtr6BXKvleoQD46KcIERE5g53wf6ZKSkpIT0+v0r5+/frTEpCIiIiIiHhGcMNLaeLYxMLv52Kbm+LWZ5ubQsiWz9njdyEWH63qEhGRM9cJJbs++eQTGjZsyDXXXEOLFi1Yvny5q69///6nPTgREREREak99a4YSJnFj4sX30m8mUW6EcPG3p+wlxjizSzamusob3qTp8MUERE5rhM6jfGZZ55h5cqVxMTEsHLlSgYMGMBjjz3G7bffjkp/iYiIiIjUrPzSCtan21m2+jfK0teB1Y+kFpdzffv6FJRWEOhrJdjP+6THzz5wAIvpS5RRjMM0WGE2JmfeR3Ry+rp+Te636h2y2vQ7ZpF5ERERTzuhZFdFRQUxMTEAtG3blsWLF3PDDTewdetWDMOokQBFRERERORQouu+l97n3pJ3eNBrHV4c+mVz7rdBTPv+OuZZLscWGsEbQ7qddMLLPyiUDGsMjkoLZv3udNv3IxZnGiWRTci5aASl81+g0BpOXFDoaXwyERGR0+uEkl3R0dGsXbuWFi1aABAeHs63337LgAEDWLt2bY0EKCIiIiIi8OPPPzOl7B94WZy8bb2DG+58kBifUioXvc39adPp7/ic3TmJFOd/RbBfzEnNERwaAcPmUFKYR+wRK7dsv/8zq15X4oJCD10nIiJyhjqhml3vvfce0dHRbm0+Pj7MnDmTRYsWVbm+sLDw1KITEREREREArD9NwttwYDWc9CqfR8oHv7KyLJ6hOy4lzwwk2Cihro+dWN/KU5onODTimFsUYxLqK9ElIiJnvBNKdiUkJBAbGwvAxIkT3fouueQSt9cFBQX07NnzFMOrXampqTRt2pT27dt7OhQREREREZeKSgdtixazNeFGKkOSSbZkM6F4DM++OY0JxWMItxTiwMJa/w4QEu/pcEVERDzqhJJdR3rssceYPn36UfuKioq4+uqrycnJOenAPCElJYW0tDRWrFjh6VBERERERFycFcWEGYXkRrbDOmgOZbYkki3ZfOY7nmRLNmW2JLZ4N8bPoZ0VIiIiJ53seu+997j33nv58ssv3dqLioro2bMn+/fvZ8GCBaccoIiIiIjI+c7XN5Bi/Nm/cz37zAhGlg9163+w9F5s5dn4hcZ6KEIREZEzxwkVqD9Sv379yMvL47bbbmP27Nl07drVtaIrKyuLRYsWERcXdzpjFRERERE5P1ks7Eu6lo67vmLYK635t+Mlt19bP132T6Is+VS2v91zMYqIiJwhTnplF8CQIUMYN24c119/PQsXLqRXr17s27ePBQsWUKdOndMVo4iIiIjIeS/oiofwpYIZzkdItmRz0KcOH144kTyCiLLkU2T68o85mWTYSzwdqoiIiEed9Mquw/7v//6P3NxcrrjiCurWrcvChQtJSEg4HbGJiIiIiJzx7CXlLNuew5I1myndvwPDL4T2bdpyXat4DhaXE+hrJdjP+5TnCfT3p8ziR4hZBEBIeQa3bhqJAyuV1iACKwuZVP44gWUdgeRTnk9ERORsddLJrhtvvNHttbe3N5GRkTz44INu7Z999tnJTiEiIiIickazl5Rz97+nc3fF+4z3WoMVBwDr/leXJ7+7nTTqEmQL440h3U454WULDsU/rh6VhYFYe/8TS2EWWP3wanAFOMqpnNybyMAorMGhp+HJREREzl4nnewKCQlxe33bbbedcjAiIiIiImeT/307nxmO/8NiMZnm9Teuu2UQUc4c4hdO4vnM5zjgtJFZHkdx/lcE+8Wc2mR+IVjv/AzKCiEkvkq3ddBc8A0Cv5Cj3CwiInL+OOlk15QpU05nHCIiIiIiZxWH0yRx9US8DBMvnPSo+J77P7uCR2/rwfMHC5hqriTSUkCAfzABvpWnZ1K/kGMns46SABMRETkfnXLNLoDS0lLWrl1LdnY2TqfT1W4YBn369DkdU4iIiIiInFH27NvHxY6V7L3oXuL3ziHZvosJxWMY+WY6E71fJ8hSign8En0TlysRJSIiUmtOOdn19ddf079/f3Jycqr0GYaBw+E41SlERERERM44loJ9+BqV5Cd0J/mqYZS904vkgt185jsegFJbEvaCInwq8j0bqIiIyHnGcqoDDB8+nL/97W9kZGTgdDrdvpToEhEREZFzVWxMHACbNq1nnxnByPKhbv2jigcSZBYSGRntifBERETOW6e8sisrK4tRo0YRE3OKBTdFRERERM4iPuEJ7LW1pN72D7jvtRheq3jJ7VfJz1RMwN8oJ6B1P88FKSIich465ZVd/fr1Y+HChachFM9LTU2ladOmtG/f3tOhiIiIiMhZwPuKx2hlbGFm5QiSLdkc9I1n8gUTOWgGEW4ppAhfHvhwDRn2Ek+HKiIict4wTNM0T2WA4uJibr75ZqKiomjevDne3t5u/Q888MApBegJ+fn5hISEYLfbCQ4O9nQ4IiIiInKC8orLOVBYRnigL+GBPjU2T0H2LirfuJwwMw+AXIIJpARvHJjeAXhVFJJpiSXwvm+wRSfXWBwiIiLng+rma055G+PMmTP55ptv8PPzY+HChRiG4eozDOOsTHaJiIiIyNlp9Z6DpH6zDr8d3xBt5nKAEMrrX01Kz5ZEBPkQ6Gsl2M/7rweqJltwKJVxDagszMZ66XDCiw+CXzA0uQ4wqZzcm8jAKKzBoadtThERETm+U17ZFRsbywMPPMDo0aOxWE55V+QZQSu7RERERM4+K3bm8tXbT/IP60cEG8U4rAF4VRZTRAAvO/vxk/clBIeF88aQbqc14UWpHcoKISS+ap89HXyDwC/k9M0nIiJynqq1lV3l5eXccsst50yiS0RERETOTj99/BJPek+l0PRjmM+zPDZ0MHU4gGP+i4z5bSq5FZ+xL6cOxflfEex3Gg9X8gs5djLraAkwERERqVGnnKEaMGAAs2bNOh2xiIiIiIiclC0ZefQrfJ8Kiz9BRikPl77CA2/+j5X2QAZtuYQC049wSyEX+BcT61vp6XBFRESkBp3yyi6Hw8G//vUv5s2bR4sWLaoUqH/ppZdOdQoRERERkePK3bSEhkYu+X2mYSx8gmT7LiYUj2Hkm+lM9H4dm6UUgHUN7qaDVluJiIic00452fXbb7/RunVrANatW+fWd2SxehERERGRmhLitAOwL6QljQfNoeydXiQX7OYz3/EAFAckEFC8F3/vU/7fXxERETnDnfJ/7RcsWHA64hAREREROWn1GzWDxfDjom8JvqE/z5QP5XVGu/r/VdyL8bxNg0ZNPRiliIiI1Ab9aktEREREznre8S3JsTWm5Y7/MORlH95wvuRWnfYRcxr7jEiI7oi/58IUERGRWnBSBerXrl2L0+ms9vXr16+nslKFQEVERESkhhgG5tUv0Nyyg8/Nf5BsySaTCF503Eap6Y2/UY6Xs4IH3ppNhr3E09GKiIhIDTqpZFfr1q3Jycmp9vWdOnVi9+7dJzOViIiIiEi1+EZdQIERjK9x6JesseTwkNdMvKIa4AiMJsZiZ1L54wSVZXs4UhEREalJJ7WN0TRNHn/8cQICAqp1fXl5+clMIyIiIiJnmfzSCrZmFTAvLYsft+bgcJq0TQ7jzk7JBPlZCfS1Euzn/dcDnQRbcCiVcRdQWWTDesPrYDohMBrv6MZg30vl5N5EBkZhDQ6tkflFRETkzGCYpmme6E1du3Y94ZMWP/jgA+Li4k50Ko/Iz88nJCQEu91OcHCwp8MREREROSvkl1bw90nfEJ6zmht8fqaLTxoW08kKR0PeLuuB3S+BiPBw3hjSrcYSXpTaoawQQuKr9tnTwTcI/EJqZm4RERGpUdXN15xUsutcp2SXiIiIyIlbs3UXIdN7kGxkk21EENDuNmwB/pT/9gU+B7dw0BnIXq94ou7/itjoGE+HKyIiImeZ6uZrTqpml4iIiIjIn61f9SPJRjaGAaVOCwPXtWRl/aHcXvwweWYAYZYiEi05xPrq4CIRERGpOSdVs0tERERE5M+CdnxNsZcN36AwkvN3M6F4DCPfTGei9+uEWopxYGWdWY9Lj7bFUEREROQ00couERERETktLqzYQJqtM9bBcymzJZFsyeYz3/EkW7IpsyWxOuYGGjm3ezpMEREROccp2XWE1NRUmjZtSvv27T0dioiIiMhZJ9DPhwP2QraXhzKyfKhb34iy+9mUU47VWkOF6UVERER+pwL1R6EC9SIiIiInzj77Sbx+fp3bzGeYxL9ItmS7+nY5o/A1Kqms152EuyZ7MEoRERE5W9VYgfqSkhLS09OrtK9fv/5EhxIRERGRc0hpq/54GSazjMdItmRzwLsOUxu/RTrRJFv2E8NBXt7XlAx7iadDFRERkXPYCSW7PvnkExo2bMg111xDixYtWL58uauvf//+pz04ERERETl7BHh7UWIEEGCUY2KwyRlP7JaZhFmKMQHDgIcq3yaoLPsvxxIRERE5WSeU7HrmmWdYuXIlq1evZsqUKQwePJgPPvgAAO2GFBERETm/2YJDCY6rT6UtHqPjfVxSx4ur44oI6HgXxsC5VIYkExmbgC041NOhioiIyDnMeiIXV1RUEBMTA0Dbtm1ZvHgxN9xwA1u3bsUwjBoJUERERET+Wn5pBXsPFrNo0wHW77PjY7XQvXE0PZvFcqCwjEBfK8F+NVwc3i8E652fQVkhhMRX6bYOmgu+QeAXUrNxiIiIyHnthJJd0dHRrF27lhYtWgAQHh7Ot99+y4ABA1i7dm2NBCgiIiIix5dfWsF1ry1hd04hl3lv4MqwTAorrfz71wv5V9gFVDqcxIT4MW1Qh1pJeB0zmXWUBJiIiIjI6XZCya733nsPq9X9Fh8fH2bOnMmwYcOqXL9u3TouuuiiU4tQRERERI5r+fYconNXMtXnLepasnGWBmFxVHCfbxlzCzvy7/J+lBhxFJVV1nyyS0RERMTDTijZlZCQcMy+Sy65BICCggJmzpzJu+++y8qVK6msrDy1CEVERETkuJb/uIiZvs/iwMK9ZSPY4N+Fl29pzDezUhlW9g5dfVdR4NOYaN/LAX9PhysiIiJSo06oQP3xLF68mAEDBhAXF8fYsWNJSEhQ0XoRERGRGmaaJhftmo4B+FDJE34zqTi4lxvf/pX/5l9IqeGHv1GBb9G+Q7W0RERERM5xp5TsyszM5IUXXqBhw4b07t2byspKPvroI/bt28eTTz55umIUERERkWNwOk26W1ayJflWCKtLvJnFhz5P08bYzIc+TxNJHpV4sdm/hWpmiYiIyHnhhLYxHqlPnz7Mnz+fbt26MX78ePr27UtgYKCrX6czioiIiNQ8L8Mk2Cjhi+I4Qm78FMfk3iRbsvjMdzwA6UYMe50RRPpqxb2IiIicH056Zdfs2bO58cYbefLJJ/n73//ulugSERERkVpi8aLILxZr5q9cO20Hw0vvd+seWXo3F5BOcGx9DwUoIiIiUrtOOtn1448/4u/vT/fu3bnwwgt56qmn2LZt2+mMTURERESqobLVHfS1LKVB0a+87POGW9+b3i8TZdj5x7aWZNhLPBShiIiISO056WTXxRdfzNtvv01GRgaPPPII33zzDY0aNeLiiy/mtddeIysr63TGKSIiIiLHYFw8lFxrFO/7PkeSkcU+ZzhjKwZSYvgTbimkgAB8AkMI9D3pChYiIiIiZw3DPI1HJm7atIl3332X9957j6ysLAzDwOFwnK7ha01+fj4hISHY7XaCg4M9HY6IiIjI8dnTcb7bE0v+HsAAfv/fu6BYcJRDSS7O0LpYBs5RkXoRERE5a1U3X3Naf7134YUX8q9//Yvnn3+e//3vf0yePPl0Di8iIiIiR+MbhCU4Fry84NYPoKwQrD4Q0xwKM2HqNVgCo8A3yNORioiIiNS407qy61yhlV0iIiLyV7bvL2T22gzsJRUkRwRwXat4Qvy9PRdQqf1QkutoK7fs6YcSXX4htR+XiIiIyGnikZVdIiIiIue6A4VlPPHfdcz5LRObn5WoIF925xbz3JyNPNijAbd3TCbYzwNJL7+QYyeztHVRREREziNKdomIiIhUU35pBTe9NBdb8R4+a7iVVvnfYym1U5GQzGfGFUyam85XP8fwwfCenkl4iYiIiIiSXSIiIiLVtX1POm9UjuVC372U7vGjpMUtBEbXo3L7Mvptf5nrfa1sLEiiyN6OYL9YT4crIiIicl6yeDoAERERkbPFhp37aGjsw8swyXUGMmDzJaxMvJP++25ivzMEP6OCesY+yovzPR2qiIiIyHlLK7tEREREqsnPvh1vw0FlYAyJRVlMKB7DyDfTmej9OrGWPCq8bVjKKykwdMCNiIiIiKdoZZeIiIhINTVxbibXDGLvDV9SZksi2ZLNZ77jSbZkU2ZLYkGTp7AZJSQ593o6VBEREZHzlpJdIiIiItXUICYEH8PJEwtyGVF2v1vf8LL7+HSdHYDgQH9PhCciIiIinOPJrhtuuIGwsDD69evn6VBERETkHGBteAVBFFNnx6eMLp3o1jem9GWurZhHeWAcRF7ooQhFRERE5JxOdj344INMnz7d02GIiIjIOSIj8EJ+s1zIM96TSbZks8sZzY1l49ltRpFsyeZayzJmll1KRmGFp0MVEREROW+d08murl27YrPZPB2GiIiInCOCyrKJMXOxGk4A6tSJZ3rztSSE+AFgGNDTsZigsmxPhikiIiJyXjtjk12LFy+mT58+1KlTB8Mw+OKLL6pck5qaSt26dfHz86Njx478/PPPtR+oiIiInDdswaGExyVTGZwEvf6Fd0gcQWVZWC7oArfMoDI4ici4RGzBoZ4OVUREROS8ZfV0AMdSVFREy5YtGTRoEDfeeGOV/lmzZjFq1CjefPNNOnbsyMsvv0zPnj3ZtGkT0dHRHohYREREznl+IVjv/AzKCiEkHjre69ZtrdMGfIPAL8RDAYqIiIjIGZvs6tWrF7169Tpm/0svvcTdd9/NwIEDAXjzzTeZPXs2kydPZvTo0Sc0V1lZGWVlZa7X+fn5Jxe0iIiInBb78kpYu9eOl8Wgfd0wQgN8PB3SH/xCjp3MComv3VhEREREpIozNtl1POXl5axcuZJHH33U1WaxWOjRowfLli074fGef/55nnzyydMZooiIiJyEnTlFPPnlehZt3o/TPNTma7Vwc7sEhlxWj/BAX4L9vD0bpIiIiIic0c7KZNeBAwdwOBzExMS4tcfExLBx40bX6x49erBmzRqKiopISEjg448/plOnTlXGe/TRRxk1apTrdX5+PomJiTX3ACIiIlLFvrwSrn55MeWVTh66siG3NPOj3GHw2cZSXv1+K5+s3EvjWBvTB3dUwktEREREjumsTHZV13fffVet63x9ffH19a3haEREROR4Zv68G78KO0O95tL3h2VELM4C4J6oi8jy6c43xQ3Iyy2nqKxSyS4REREROaazMtkVGRmJl5cXWVlZbu1ZWVnExsZ6KCoRERE5FT/+tpUlAQ8R6CxgbmUHJvkPZEinOmT+MJ2nnK/yD98A9luTifPtCvh7OlwREREROUNZPB3AyfDx8aFt27bMnz/f1eZ0Opk/f/5RtylWV2pqKk2bNqV9+/anI0wRERE5Ae0KviXIWYABtLTuYkF+PN3nRfKP4gHYCSLEKCa0IvvQSYgiIiIiIsdwxia7CgsLWb16NatXrwZgx44drF69mt27dwMwatQo3n77baZNm8aGDRu4//77KSoqcp3OeDJSUlJIS0tjxYoVp+MRRERE5ARc47Wc7f7NIawu8WYWH/o8TRtjMx/6PE0IhTiwsD2glU48FBEREZHjOmO3Mf7yyy9069bN9fpwAfkBAwYwdepUbrnlFvbv388TTzxBZmYmrVq14uuvv65StF5ERETODo2MdFILe3HFNcOI+vQmki1ZfOY7HoC9xPCbI4lOfgc9G6SIiIiInPEM0zRNTwdxpsnPzyckJAS73U5wcLCnwxERETkvOCc05dOS1owu+jst2eRKdAHcWDaekd6f0a5uOP6D/uu5IEVERETEY6qbrzljtzGKiIjI+aW4wTX0qFxMXXMvE71fd+ub5PMalxi/8VpmUzLsJR6KUERERETOBkp2HUEF6kVERDzH7HgvfkY5c/weI9mSTVFAIksuex9HYCx1jBychoWdQa0I9D1jqzCIiIiIyBlA2xiPQtsYRUREPMCejvPt7lgKMw+9jmgATgcc3AEWb3BW4Ayti2XgHBWpFxERETkPVTdfo1+NioiIyJnBNwhLaCJYfeDioZCzDSxWqN8NohrD9OuwBEaBb5CnIxURERGRM5iSXSIiInJm8AuBOz6FssKjr9y6a86hRJdfSO3HJiIiIiJnDSW7RERE5MzhF3LsZJa2LoqIiIhINSjZJSIicg7KLSrjf2v28W1aNnsPFhMa4MP1repwc7tECkorCPS1Euzn7ekwRUREREROOyW7jpCamkpqaioOh8PToYiIiJy07IJSrpiwiILSStokhdLzolh2HSjm2dkbmLJ0J5VOJzHBfkwb1EEJLxERERE55yjZdYSUlBRSUlJc1f1FRETORi99s5nC0koADhSWM6BTXeqE+rN8ew63v7Mch9PEarFQVFapZJeIiIiInHOU7BIRETmHFJVV8tXaDO7uGMGBDT/SPP8ntr/2HD4XJPPh7otwOBsQSy7P9mpMXIi/p8MVERERETntlOwSERE5h2zIyMdSlsfDm1PwLt9PrjWEX8svIG/TUiZaPuYu/yaEO/fj+10daPiVTjYUERERkXOOkl0iIiLnmLu9ZmMt2Q9AYFAwYw8MIoNwbrT8wASftzAMk/xyPygrVLJLRERERM45Fk8HICIiIqdP0xh/bvVeRFrE1VSGJONbuIcPfZ6mjbGFB62fYWBimpBz1WsQEu/pcEVERERETjslu0RERM4hATnriSKPcVmX0LfwMXY5o0m2ZPOZ73iSLdnsdkbhwEJY/kZPhyoiIiIiUiOU7DpCamoqTZs2pX379p4ORURE5ORUlgFQbASwrsjGI84Ut+5/OFIox5v3lmwhw17iiQhFRERERGqUkl1HSElJIS0tjRUrVng6FBERkZMT3QTTy4ebg9fTLDCfV3zfdOueGfI6AUYZmYEXEuir0p0iIiIicu5RsktERORcEhCOcdFNDHB8zn/9nyLGkQFhdWHQNxCSiLU4G9PizejbehLs5+3paEVERERETjslu0RERM41nYZhKS/EWrgPfIKg7V2waQ5UFINhwXBWYPuwL9jTPR2piIiIiMhpp/0LIiIi55rQRIi9CHJ3gLc/fDce/MOg+d/gopvg83sgMAp8gzwdqYiIiIjIaadkl4iIyLnGLwTu/C+UFUJIPDidYDliMfddcw4luvxCPBejiIiIiEgNUbJLRETkXOQX8kcyy/KnqgUh8bUfj4iIiIhILVGyS0RE5HcFpRWs2p2H02nSLD6YaJufp0MSEREREZETpGTXEVJTU0lNTcXhcHg6FBERqUU5hWW8+M0m/rt6H8Xlh/4bYLUY9G4ex9Bu9akT6q+TC0VEREREzhKGaZqmp4M40+Tn5xMSEoLdbic4ONjT4YiISA3KKy6ny78XYi+pYGDnuvTvlIyvtxffrM/kle+2UFReSZO4YGYM6aiEl4iIiIiIB1U3X2M5Zo+IiMh5YOGmbOwlFQDM35iNn7cX8aH+9GwWi5+3FxUOk925xRSVVXo4UhERERERqQ4lu0RE5Lw2f+N+6kcFkhjmz+7cYm79z0+s3JXLrf/5icz8UgJ8vAgP9CEuxN/ToYqIiIiISDUo2SUiIue1THsJLRJCmXV3B/4WksbA/DdY9/Y9tMn7hgZhVu7qXJcDBWWeDlNERERERKpJBepFROS8FhHoi3PfWurMGMS/yrax0xJDGd4MsH5LBR8xI30EyYEXeDpMERERERGpJiW7RETkvHZrE286/W8kTguM8RnNzPzmgEE9I4PXjDcZsPtxrg1qBKWXgl+Ip8MVEREREZG/oG2MIiJyXru89Ht8jEosZiX3lU6mTUgxH97TEX//QIIdB7EYJuEV+6Cs0NOhioiIiIhINWhll4iInNfKty5kk9GEUEcOyZZsJpaMYeTbQ3nD+3WSLAfIN/0pK/emknDiPB2siIiIiIj8Ja3sEhGR85q1soQ8n1j+EfAsZbYkki3ZfOY7nv9v787jo6oP9Y8/Z5bMZB0SIIGQhX1TNtlEKkhFEa243Yq1KlqvrYqtiktrl6v86i1eaYGqVFzqRtWCdV9QFEEU2RFcWJQtgbCTfZ/l/P6YJCQkQIAkZ2byeb9e05k558yZZwaK8OT7/Z5M2wH5PJl6xzlebqNSsS5+PgQAAACEA8quWmbPnq2+fftq6NChVkcBALQQR0ofnRu1RY//arxcP3227r6rntFPUw8putMZSnA7LUoIAAAA4GQYpmmaVocINYWFhfJ4PCooKFBCQoLVcQAAzWnvBumpUdLIu6SNb0l5O4/si02WSg5IVzwlDbjGooAAAAAApMb3NYzsAgC0bh0HSMN+KS2bFSy6YttL5z8kRcUGiy5njJQxwuKQAAAAABqLsgsA0LoV5Eg/LAw+tkdJJQelRQ9J7jZSdKLkLZVemhA8DgAAAEDIY7VdAEDr5ooLjuaSpEnvBUdyBXzBbUV7pBcuCT52xVmbEwAAAECjUHYBAOooLPeqpMInQ4be3bBHh4orlJLg1oSBqfL6A4p1OSJrsXa3R7rudamiWPJ0qrvPkybd+EGw6HJ7rMkHAAAA4KRQdgEAahSWe3XDP1dq+8ESlVT4FeWwqYPHrZz8Mk1bsEkxUQ51aRejl24eHnmF17HKrKMLMAAAAAAhjbILAFCjpMKnHYdKVVjukyfaqf/cOkI9UuK1eW+hrn5quQrKvMo6XKqSCl9klV0AAAAAIgYL1AMAasS7nfL6A4p3O1RQ5tXNL67R2qxc/XLuWhWW+xQbZZfNMNQuzmV1VAAAAABoEGUXAKDG0u8PqrTSrxduGqqMpBhl55bqqieXKzu3VBlJMXri2rN0uKRSq3fmWh0VAAAAABpE2QUAqFFc4ZMkDUhro5kTB9TZN3PiAPVLC65rVVzua/FsAAAAANAYlF0AgBpd28VKkj74Zq/unrehzr67523Qh9/uCx7XPrbFswEAAABAY1B2AQBqDM5MVJd2sbr3P1/XTF18/bYRNVMap777nQakedQ9Od7qqAAAAADQIK7GWMvs2bM1e/Zs+f1+q6MAgCX2FZbLKC9QF1+OLnGu0aTACsXNzdPbznZ61nGOPvSdpeKCjtpbUKaOnmir4wIAAABAPYZpmqbVIUJNYWGhPB6PCgoKlJCQYHUcAGgxhfmHdfDxC5Tp3yG/4dQbvnO0NdBJfW3ZmmBfLimgLHsXJf96oRLatLU6LgAAAIBWpLF9DSO7AAA1EoxyxRs5MhSQI6G9xl8zU4fs7ZUcOCTnKxdLhTnqZtsjwyi3OioAAAAANIg1uwAAR1QUyfCVS7HtpYLdajP/cnUv/04J8y6TCnOk6CQZvgrJMKxOCgAAAAANouwCAByxZ13w/hcfSYmdpbyd0nMXBu8TO0vX/Ucy/dLeDcc+BwAAAABYiLILQKtUWO7VnvxSrd+Vr/mrd+mdDXtUUOqVJO0tKFNhudfihBaxVc1uj06Urni67r4rnpbcbeoeBwAAAAAhhn+tAGh1Csu9+umcL7X9YIm8/iPX6HA7bbpyUJo+/+Gg2sW79OIvhinB7bQwqQU6nxssslbOkb6eV3ffm7+Uel4kRcVJ6cOtyQcAAAAAJ0DZBaDV2bSnUN/vL5ZpSu3jXHr99hFyO+x6aul2/fOLHZIkwzBUUuFrfWVXQkep90+kzx6VZEptMqUrn5HeuCU4lXHlHGnwjZKbK9UCAAAACE1MYwTQ6sxbvUvJ8S6ltYnWweIKXffsKu3KK9XHG/fXHPPIlf3U0RNtYUqLFORUrdtVNeLN6Q4WXDb7kWO2LQ4eBwAAAAAhiLILQKtS4fPrvW/2atI5nTX/1hHKSIpRdm6prnpyubJzS5WeGK020U4t23bI6qjWcMVJcSnBEV2Xz5FSzpRKDkppw6T/ei64PS45eBwAAAAAhCCmMUaownKvSip8ahvr0t6CMjntNnX0uGUYhvYWlCnW5Wh907MASUXlPlX6AurePk6pbaI1c+IAXfXk8pr9s64ZqP/33iYdKqq0MKWF3B7putelimLJ00ka+LO6+9PPDhZdbo81+QAAAADgBCi7IlBhuVfX/3Oldhwskc1mKL/qCnM9U+L008Hpemn5zta7+DZavQS3UzFRdm3cW6gzO3l097wNdfbf+e/1Kij16rye7S1KGALcnmOXWZ5OLZsFAAAAAE4SZVcEyiup1Oa9RarwBRTncujv1wxUtNOul1dm638/2CSpFS++jVYvymHTZQM7ae7yLL22Zrdy8suUkRSjmRMH6O55G5SdWypJ+lGPdhYnBQAAAACcCtbsikBLthyUz28qOd6l4gqf/rbwe7WNi9KOQyU1x/zlijNb5+LbgKSrh6Qpr7RSOfllahsbpdnXDlJMlENDOyfWHHPP/A3aW1BmYUoAAAAAwKmg7IpAr67K1rgzU/TW5JENLr6dFBulTzYdsDomYJluyXHq3SFebodNh0sqdekTyzT+759r8ZaDuu28bkpPjFbbuCjFuhj8CgAAAADhhn/JRaDtB0v0s2EZx1x8+5mlO7S91igvhKbqiwwUlfu0ckeuZJoa1qWtenWI5yIDpynB7dS/fzVCJRXBxep/2F+s6Ci7Bmcmyu2064YRmXy/AAAAABCmKLsiUJzboX2F5dqTX1Zv8e27521QvNuhzm1jLUqHxigs9+raZ1Zo6/5ilfsCctgMSZIvYGpwZqL2FZQrOYGLDJyOBLez5rvLPOr/D0zxBQAAAIDwRdlVy+zZszV79mz5/X6ro5yWi/t10GtrdundDXu0O6/hxbd/NizD4pQ4nkNFFdqyr0hev6m2sVF64/Zz1NETrX+vztbUdzbKb5oyDHGRAQAAAAAAjsKaXbVMnjxZGzdu1OrVq62OclomDEjV4eJK7c4rU4cEt/79y7M1KD1Rd43tIbsRHCH01GfbWHw7hC3fflj+gKkOCW4dLqnU9f9cpW9y8vXs5zvkN01J0sQh6YxAaioBv1T1vQIAAAAAwhsjuyJQ744J6p4Sp+0HSrSvsFzXPL1CpZV+HSqu0JmdEpRX4lW7eBeLb4ewt77K0Zheyfrz5WfqmqdX1FxkQJIykmKU2TZGn245oF+f38PipGGs5KC06p/Sd29Ih7ZIDrfU+xLpnN9Ise0lV5zk9lidEgAAAABwkmg7IlCC26nXbztH+aWV+mZ3oTbszpfTbujHvZN1Vkai9hWWs/h2iMstqVT/tDYNXmRg5sQBWvjdfn303T4LE4a54gPS44OlikKp53hpxGSp9JD01cvSP8dK7kQpMVO67nUKLwAAAAAIM5RdEap68e2MpFhd0r9jnX1MfQt9qW2i9c3ugmNeZKBdXJQ6JfLreMqWzw4WXZJ0cJN0yV8lT5p0xhXSkyOlkgOS0y1VFFN2AQAAAECYYc0uIARNHJquVTtzdfnsZcrOLVVGUoxev22EMpJilJ1bqnXZ+bqgb4rVMcNTICB9+7p0xpVSYmcpb6f0wiVS9kpp7hWSN3gRB501SfJ0sjIpAAAAAOAUUHYBIWhAukduh00HiiqU4Hbof684U/Fup8b2SZZRdcw/P9/BRQZORVmeVLBLOuNy6cb3jxRez10YvE/sLCWfIeVnWxoTAAAAAHBqKLuAEOSJjlKvjvGKdznkC5i6/p+rdOHMpZq3epeuGpymtMRoLjJwqhxRwfuKouDUxSuerrv/8qckX5nkcLV8NgAAAADAaeNfykAISnA7Nffm4Sqp8Cne7dS3OQUyTenMTgmKdzu1t6CMiwycKle8lDlS+upfUpfR0pu/rLv/tUlS8T6p5zhr8gEAAAAATgsju4AQleB2qqMnWnEuh87u2lYjurVVfFW51dETTdF1OkbeKWUvl54ccWTq4i8WSvEdg0WXPUpq19PqlAAAAACAU0DZBaD1STlTik4KTmWUISV0kj64VyraK9mckr9SevFSqSDH6qQAAAAAgJPENEYArY8rTmrbTSqKkXr/RCrMCa7fNeo+qWN/6aXLpNj2weMAAAAAAGGFsgtA6+P2SNe9LlUUS55O9fff+EGw6HJ7Wj4bAAAAAOC0UHYBaJ3cnmOXWQ0VYAAAAACAsMCaXQAAAAAAAIgYlF0AAAAAAACIGJRdAAAAAAAAiBiUXQAAAAAAAIgYlF0AAAAAAACIGJRdAAAAAAAAiBiUXQAAAAAAAIgYlF0AAAAAAACIGJRdAAAAAAAAiBgOqwMAOIGyfGnvesk0pdRBUnQbiwMBAAAAABC6KLuAUFW0T/pkqrTxLclbGtzmiJYG/kw6+3YpLllyeyyNCAAAAABAqInYaYzvvfeeevXqpR49eujZZ5+1Og5wckoPS08MlTa8Kg25SbpjjXTHWunce4Lb5vxImnuFVF5gdVIAAAAAAEJKRJZdPp9PU6ZM0aeffqqvvvpK06dP1+HDh62OBTTet29IFYWSTGnz+5IzWmrXPTiqy+WRfOVSXpZUUWx1UgAAAAAAQkpEll2rVq3SGWecoU6dOikuLk7jx4/XwoULrY4FNN7m96W0oVJiZylvp/TCJVL2yuB98T7J4ZbadpM8naxOCgAAAABASAnJsmvp0qW69NJLlZqaKsMw9NZbb9U7Zvbs2ercubPcbreGDx+uVatW1ezbs2ePOnU6UgJ06tRJOTk5LREdaBqFOcGy68b3jxRez10YvE/sLA34mVRy0NqMAAAAAACEoJAsu0pKSjRgwADNnj27wf3z5s3TlClT9OCDD2rdunUaMGCAxo0bpwMHDrRwUqCZxLSTDm+TPGnSFU/X3XfF01Lx/uAxAAAAAACgjpAsu8aPH6+HH35YV1xxRYP7Z8yYoVtuuUU33XST+vbtqzlz5igmJkbPPfecJCk1NbXOSK6cnBylpqYe8/0qKipUWFhY5wZYqv/V0taPpW2fSm/+su6+/9woff+hNGCiJdEAAAAAAAhlIVl2HU9lZaXWrl2rsWPH1myz2WwaO3asli9fLkkaNmyYvv32W+Xk5Ki4uFgLFizQuHHjjnnOadOmyePx1NzS09Ob/XMAx9V/otS2u/Svq4JTF9tkSJPeC47mKtwjGXapyyirUwIAAAAAEHLCruw6dOiQ/H6/UlJS6mxPSUnRvn37JEkOh0N/+9vfNGbMGA0cOFD33HOP2rZte8xzPvDAAyooKKi57dq1q1k/A3BCZXmSt0wyA8Hn+bukFy+VSg8FF6cPeKWXfyoVsBYdAAAAAAC1OawO0FwmTJigCRMmNOpYl8sll8vVzImAk+CKk+I7SDa7NOEJKXebZJpS5kgpKiZ4VcbY9sHjAAAAAABAjbAru9q1aye73a79+/fX2b5//3516NDBolRAE3N7pOtelyqKJU8nqcu5dfff+EGw6HJ7rMkHAAAAAECICrtpjFFRURo8eLAWLVpUsy0QCGjRokUaMWKEhcmAJub2BIuuhng6UXQBAAAAANCAkBzZVVxcrK1bt9Y837Fjh9avX6+kpCRlZGRoypQpmjRpkoYMGaJhw4Zp1qxZKikp0U033XRa7zt79mzNnj1bfr//dD8CAAAAAAAALGCYpmlaHeJoS5Ys0ZgxY+ptnzRpkl544QVJ0hNPPKHp06dr3759GjhwoB577DENHz68Sd6/sLBQHo9HBQUFSkhIaJJzAgAAAAAA4NQ1tq8JybLLapRdAAAAAAAAoaWxfU3YrdkFAAAAAAAAHAtlFwAAAAAAACJGSC5Qb5WIXKB+/3fSnq8ke5TUZbQUn2J1IgAAAAAAgGbDml0NiIg1u/Z9I717p5Sz9sg2m0MacI008m4prr3k9liXDwAAAAAA4CQ0tq9hZFckOrxVeubHkhmQLpkhDbpOqiyRNrwqLfqz9M3rUsoZ0vVvUHgBAAAAAICIwppdkWjlU5LfKwV80pePSSUHpZgkqe9lkite8pVJBbukimKrkwIAAAAAADQpyq5ItPn94HTFxM5S3k7phUuk7JXB+5IDks0ppQ+XPJ2sTgoAAAAAANCkKLsiUfEBqdNg6cb3jxRez10YvE/sLHUdLVUyqgsAAAAAAEQeyq5aZs+erb59+2ro0KFWRzk9CanSvq8lT5p0xdN1910+R8rdHjwGAAAAAAAgwnA1xgaE/dUYP3tU+mKmdN1/pLduD47oqhbbTio5JP3iIynjbMsiAgAAAAAAnIzG9jWM7IpEw26R4lKkFy4NFl2eNOmqfwavvFhySHLGSgms1wUAAAAAACIPZVckqiwNXonR9AefF+yWXr85uM2VIHlLpBd/IhXkWJsTAAAAAACgiTmsDoBm4IqT4jtINntwRFfJQcnulNLPlsrzg1dljG0fPA4AAAAAACCCUHZFIrdHuu51qaJY8hw1XdEVJ934QfDe7bEmHwAAAAAAQDNhGmMtEXM1RilYZB1ddFXzdKLoAgAAAAAAEYmrMTYg7K/GCAAAAAAAEGG4GiMAAAAAAABaHcouAAAAAAAARAzKLgAAAAAAAEQMyi4AAAAAAABEDMouAAAAAAAARAzKLgAAAAAAAEQMyq5aZs+erb59+2ro0KFWRwEAAAAAAMApMEzTNK0OEWoKCgrUpk0b7dq1SwkJCVbHAQAAAAAAaPUKCwuVnp6u/Px8eTyeYx7naMFMYaOoqEiSlJ6ebnESAAAAAAAA1FZUVHTcsouRXQ0IBALas2eP4uPjZRiG1XFajaFDh2r16tVWx0Arxe+/8MevYfPi+z221vjdROJnDtfPFOq5q38Cz4wJAKcq1P+cQ8syTVNFRUVKTU2VzXbslbkY2dUAm82mtLQ0q2O0Ona7nb8EwTL8/gt//Bo2L77fY2uN300kfuZw/UzhkjshISEscgIIPeHy5xxazvFGdFVjgXqEjMmTJ1sdAa0Yv//CH7+GzYvv99ha43cTiZ85XD9TuOYGgMbizzmcCqYxAgAAAGgWhYWF8ng8KigoYGQGAKDFMLILAAAAQLNwuVx68MEH5XK5rI4CAGhFGNkFAAAAAACAiMHILgAAAAAAAEQMyi4AAAAAAABEDMouAAAAAAAARAzKLgAAAAAAAEQMyi4AAAAAAABEDIfVAUJRIBDQnj17FB8fL8MwrI4DAAAAAADQ6pmmqaKiIqWmpspmO/b4LcquBuzZs0fp6elWxwAAAAAAAMBRdu3apbS0tGPup+xqQHx8vKTgl5eQkGBxGgAAAAAAABQWFio9Pb2mtzkWyq4GVE9dTEhIoOwCAAAAAAAIISdacooF6gEAAAAAABAxKLsAAAAAAAAQMSi7AAAAAAAAEDEouwAAAAAAABAxKLtagcWbD6io3Gt1DAAAAAAAgGbH1RgjXNbhEv1y7hrFu5369Y+769rhGXI57FbHAgAAAAAAaBaM7Ipwh0sqlZ4Yo9ySSk19d6PGzvhMb6/PUSBgWh0NAAAAAACgyRmmadJ6HKWwsFAej0cFBQVKSEiwOs5p8/kDmr9mt2Z98r0OFFVIkvp2TNBvx/fWqB7tZBiGxQkBAAAAAACOr7F9DWVXAyKt7KpWVunXc8t2aM6SbSqq8EmSzunWVr+9qLcGpLexNhwAAAAAAMBxUHadhkgtu6rllVRq9uKteml5lir9AUnSJf066t5xvdSlXazF6QAAAAAAAOoLu7IrEAjos88+0+eff66srCyVlpaqffv2GjRokMaOHav09PRGn2vatGl64403tHnzZkVHR+ucc87R//3f/6lXr16Nen2kl13VdueVaubHP+iNr3bLNCWHzdA1w9L1m/N7KDnebXU8AAAAAACAGo3tayxfoL6srEwPP/yw0tPTdfHFF2vBggXKz8+X3W7X1q1b9eCDD6pLly66+OKLtWLFikad87PPPtPkyZO1YsUKffzxx/J6vbrwwgtVUlLSzJ8mvKQlxuhvVw/QgjvP1Y97J8sXMPWvFdka/egSzVi4RUXlXqsjAgAAAAAAnBTLR3alp6drxIgRuvHGG3XBBRfI6XTWOyYrK0uvvPKKnnrqKf3hD3/QLbfcclLvcfDgQSUnJ+uzzz7TqFGjTnh8axnZdbQV2w/rkQWbtX5XviQpKTZKd4zprp+fnSGXw25tOAAAAAAA0KqFzTTGTZs2qU+fPo061uv1Kjs7W926dTup99i6dat69Oihb775RmeeeWa9/RUVFaqoqKh5XlhYqPT09FZXdkmSaZr66Lv9evSjzdp+MDgSLi0xWvde2EsTBqTKZuPKjQAAAAAAoOWFTdnV3AKBgCZMmKD8/Hx98cUXDR7z0EMPaerUqfW2t8ayq5rPH9Bra3dr5sff60BRsAjs0zFBv72ol0b3bC/DoPQCAAAAAAAtJ6zLLp/Pp6eeekpLliyR3+/XyJEjNXnyZLndJ79o+m233aYFCxboiy++UFpaWoPHMLLr2Moq/Xpu2Q7N+Wybisp9kqQRXdvqd+N7a0B6G2vDAQAAAACAViOsy67bb79d33//va688kp5vV699NJL6tmzp1599dWTOs8dd9yht99+W0uXLlWXLl0a/brWumbX8eSVVOofS7bqxS+zVOkPSJIu6ddR947rpS7tYi1OBwAAAAAAIl1YlV1vvvmmrrjiiprn3bt315YtW2S3BxdF37x5s84++2zl5+c36nymaerXv/613nzzTS1ZskQ9evQ4qTyUXce2O69UMz/+QW98tVumKdlthq4Zmq47z++h5ISTH3kHAAAAAADQGGFVdl166aWy2+36xz/+odTUVF199dXyeDy66qqr5PV69cwzz6isrEwff/xxo853++2365VXXtHbb7+tXr161Wz3eDyKjo4+4espu05s875CTf9wixZtPiBJinba9d/ndtEvR3VVvLv+FTUBAAAAAABOR1iVXZI0b948/elPf9Kvf/1rXX/99frzn/9cZ82uhx56SO3bt2/UuY61ePrzzz+vG2+88YSvp+xqvJXbD+uRDzfrq+x8SVJijFN3/LiHrjs7Qy6H3dpwAAAAAAAgYoRd2SVJ+fn5uv/++7VhwwbNmTNHgwYNsiQHZdfJMU1TCzfu16Mfbta2gyWSpLTEaN1zYU9dNqCTbDau3AgAAAAAAE5PWJZd1ZYuXarJkyfroosu0p///OdTugrj6aDsOjU+f0D/WbtbMz/5XvsLg1e37NMxQfdf1Evn9Wx/zBF3AAAAAAAAJ9LYvsbWgpmOKTs7W1dffbX69eunn//85+rRo4fWrl2rmJgYDRgwQAsWLLA6IhrBYbfpmmEZWnLvGN1/US/Fux3atLdQNz2/Wj97ZoXW78q3OiIAAAAAAIhwITGy67zzzlOHDh1044036qOPPtK2bdv0zjvvSJI2bdqkX/3qV+rQoYPmz5/fInkY2dU08ksr9Y8l2/TClztV6QtIki7u10H3XthLXdvHWZwOAAAAAACEk7CaxhgXF6cNGzaoW7duMk1TXbp00c6dO+sc8/TTT+uXv/xli+Sh7GpaOfllmvnx93p93W6ZpmS3GZo4NF13nd9DyQktO0UVAAAAAACEp7Aqu0aPHq20tDRNmjRJn3zyiTZt2qR3333XsjyUXc1jy74iTf9osz7ZdECSFO206+YfddEvR3dVgttpcToAAAAAABDKwqrsysrK0j333KNNmzZp4MCBmj59ulJTUy3LQ9nVvFbtyNUjCzZpXXa+JCkxxqnJY7rr+hGZcjns1oYDAAAAAAAhKazKrlBD2dX8TNPUwo379eiHm7XtYIkkqVObaN1zYU9dNrCT7Dau3AgAAAAAAI4Im7KrpKREsbGxzXb8qaDsajk+f0Cvr9utGR9/r/2FFZKk3h3i9dvxvXVez/YyDEovAAAAAADQ+L7G1oKZGtS9e3c98sgj2rt37zGPMU1TH3/8scaPH6/HHnusBdOhuTnsNk0cmqEl947Rby/qrXi3Q5v3Femm51frmqdX6KvsPKsjAgAAAACAMGL5yK4tW7bo97//vd5//30NGDBAQ4YMUWpqqtxut/Ly8rRx40YtX75cDodDDzzwgH71q1/Jbm/edZ0Y2WWd/NJKPblkm57/cqcqfQFJ0vgzO+jecb3UrX2cxekAAAAAAIBVwmYaY7Xs7Gy99tpr+vzzz5WVlaWysjK1a9dOgwYN0rhx4zR+/PhmL7mqUXZZb09+mWZ+/L1eX7dbAVOy2wxNHJquu87voeQEt9XxAAAAAABACwu7siuUUHaFji37ijT9o836ZNMBSZLbadPNP+qiX43upgS30+J0AAAAAACgpVB2nQbKrtCzemeuHlmwWWuzgmt4tYlx6o4x3XXd2ZlyO1tmxB8AAAAAALAOZddpoOwKTaZp6uON+/XoR1u09UCxJKlTm2hNuaCnLh/USXYbV24EAAAAACBSUXadBsqu0ObzB/T6ut2a+fEP2ldYLknq3SFev72ot87r1V6GQekFAAAAAECkoew6DZRd4aHc69cLX+7UPxZvVWG5T5I0rEuSfje+t87KSLQ4HQAAAAAAaEqUXaeBsiu8FJR69Y/Ptur5ZTtV6QtIki46o4Puu6iXurWPszgdAAAAAABoCo3ta2wtmKlRPv/8c1133XUaMWKEcnJyJElz587VF198YXEyhCpPjFMPjO+jJfeep6uHpMlmSB9+t08XzlyqB974RvurpjoCAAAAAIDIF1Jl1+uvv65x48YpOjpaX331lSoqKiRJBQUF+stf/mJxOoS61DbRevS/BujDu0ZpbJ8U+QOmXl2VrdHTF+vRDzeroMxrdUQAAAAAANDMQqrsevjhhzVnzhw988wzcjqdNdtHjhypdevWWZgM4aRnSryenTRE/7l1hAZnJqrcG9A/lmzT6OmL9ezn21Xu9VsdEQAAAAAANJOQKru2bNmiUaNG1dvu8XiUn5/f8oEQ1oZ0TtJ/bh2hZ24Yoh7Jccov9erh9zfp/L99pv+s3S1/gOXqAAAAAACINCFVdnXo0EFbt26tt/2LL75Q165dLUiEcGcYhi7om6IFd56rR6/qr44et3Lyy3Tvaxt08d8/16eb94trNAAAAAAAEDlCquy65ZZbdOedd2rlypUyDEN79uzRyy+/rHvvvVe33Xab1fEQxhx2m64emq7F956nB8b3VoLboS37i/SLF9Zo4tMrtC47z+qIAAAAAACgCRhmCA1rMU1Tf/nLXzRt2jSVlpZKklwul+699179+c9/brEcjb2UJcJXQalX//hsq15YtlMVvoAkadwZKbpvXG91T46zOB0AAAAAADhaY/uakCq7qlVWVmrr1q0qLi5W3759FRfXsuUDZVfrsbegTLM+/kGvrd2lgCnZDGni0HTdeX5PdfC4rY4HAAAAAACqhGXZVVBQIL/fr6SkpDrbc3Nz5XA4Wqx4ouxqfX7YX6RHP9qijzfulyS5HDb94kdddOvobvJEO0/wagAAAAAA0Nwa29eE1Jpd11xzjf7973/X2z5//nxdc801FiRCa9EjJV7P3DBE/7l1hIZkJqrCF9CTS7Zp1KOL9czS7Sr3+q2OCAAAAAAAGiGkRnYlJSVp2bJl6tOnT53tmzdv1siRI3X48OEWycHIrtbNNE0t2nRA//fhZv1woFiSlOpx6+4LeurKs9JktxkWJwQAAAAAoPUJy5FdFRUV8vl89bZ7vV6VlZVZkAitkWEYGts3RR/eNUqP/ld/dfS4taegXPf952uN//tSLdq0XyHUEQMAAAAAgFpCquwaNmyYnn766Xrb58yZo8GDB1uQCK2Z3Wbo6iHpWnzvefr9xb3liXbq+/3FuvnFNZr41AqtzcqzOiIAAAAAADhKSE1jXLZsmcaOHauhQ4fq/PPPlyQtWrRIq1ev1sKFC3Xuuee2SA6mMaIhBaVePfnZNj2/bIcqfAFJ0oV9U3T/Rb3UPTne4nQAAAAAAES2sLwaoyStX79e06dP1/r16xUdHa3+/fvrgQceUI8ePVosA2UXjmdvQZlmffyDXlu7SwFTshnSTwen664LeqijJ9rqeAAAAAAARKSwLbtCAWUXGmPrgSI9+uEWLdy4X5Lkcth008guum10N3linBanAwAAAAAgsoRt2RUIBLR161YdOHBAgUCgzr5Ro0a1SAbKLpyMtVm5emTBZq3eGVzDyxPt1OQx3XTDiM5yO+0WpwMAAAAAIDKEZdm1YsUKXXvttcrKyqp3tTvDMOT3+1skB2UXTpZpmlq06YAe/Wizvt9fLEnq6HHr7gt66qqz0mS3GRYnBAAAAAAgvIVl2TVw4ED17NlTU6dOVceOHWUYdQsCj8fTIjkou3Cq/AFTb6zbrZkff689BeWSpJ4pcbp/XG+d3ye53u9pAAAAAADQOGFZdsXGxmrDhg3q3r27pTkou3C6yr1+zV2epScWb1VBmVeSNLRzon43vrcGZyZZnA4AAAAAgPDT2L7G1oKZTmj48OHaunWr1TGA0+Z22nXLqK5aev8Y3XZeN7kcNq3emaernlyuW15aox/2F1kdEQAAAACAiBRSI7vefPNN/fGPf9R9992nfv36yemse0W7/v37t0gORnahqe0rKNesT77X/DW7FDAlmyH9dHC67rqghzp6oq2OBwAAAABAyAvLaYw2W/2BZoZhyDRNFqhHRNh6oEjTP9qij77bL0lyOWy6cWRn3T66uzwxzhO8GgAAAACA1issy66srKzj7s/MzGyRHJRdaG5rs/L0fws2a9XOXElSgtuhyWO6a9I5neV22i1OBwAAAABA6AnLsitUUHahJZimqcVbDuj/FmzRlqo1vDp63Lp7bE9dNThNdhtXbgQAAAAAoFpYLlAvSXPnztXIkSOVmppaM9Jr1qxZevvtty1OBjQtwzD0494p+uDOc/XXnw5QqsetvQXluv/1r3XRrKX6eON+0UUDAAAAAHByQqrsevLJJzVlyhRdfPHFys/Pr1mjq02bNpo1a5a14YBmYrcZ+q/Bafr03vP0h4v7qE2MUz8cKNYtL63RT+cs15qqqY4AAAAAAODEQmoaY9++ffWXv/xFl19+ueLj47VhwwZ17dpV3377rc477zwdOnSoRXIwjRFWKijz6qnPtum5ZTtU7g1Iksb2SdFvL+qlHinxFqcDAAAAAMAaYTmNcceOHRo0aFC97S6XSyUlJRYkAlqeJ9qp+y/qrSX3jtHPhqXLZkifbNqvcbOW6g9vfqPDxRVWRwQAAAAAIGSFVNnVpUsXrV+/vt72Dz/8UH369Gn5QICFOnjcmnZlfy28e7Qu7JuigCm9vDJb501foqc+26YKn9/qiAAAAAAAhByH1QFqmzJliiZPnqzy8nKZpqlVq1bp1Vdf1bRp0/Tss89aHQ+wRPfkOD19wxAt33ZYD7+/Ud/tKdS0BZv18spsPTC+ty46s4MMgys3AgAAAAAghdiaXZL08ssv66GHHtK2bdskSampqZo6dapuvvnmFsvAml0IVf6AqdfX7db0j7boYFFwOuOwLkn60yV91S/NY3E6AAAAAACaT2P7mpApu3w+n1555RWNGzdOKSkpKi0tVXFxsZKTk1s8C2UXQl1JhU9PfbZNTy3drgpfQIYhXTkoTfdf1EspCW6r4wEAAAAA0OTCruySpJiYGG3atEmZmZmW5qDsQrjYk1+mRz/crLfW75EkRTvtunV0N/1yVFdFR9ktTgcAAAAAQNMJy6sxDhs2TF999ZXVMYCwkdomWrOuGaQ3bz9HZ2W0UZnXr5mffK8f/22J3vxqtwKBkOmyAQAAAABoESFVdt1+++2655579MQTT2j58uX6+uuv69waa+nSpbr00kuVmpoqwzD01ltvNV9oIAQMykjU67edo8d/Nkid2kRrb0G57p63QVc8+aXWZuVaHQ8AAAAAgBYTUtMYbbb63ZthGDJNU4ZhyO/3N+o8CxYs0LJlyzR48GBdeeWVevPNN3X55Zc3OgfTGBHOyr1+/fOLHfrH4q0qqQz+f+aS/h31u4t6Kz0pxuJ0AAAAAACcmrBcsysrK+u4+09lLS/DMCi70CodKCrXjIXfa96aXTJNKcph080/6qLbz+umeLfT6ngAAAAAAJyUxvY1jhbMdEJWLUxfUVGhioqKmueFhYWW5ACaUnK8W49c1V83jOish9/fqC+3HdaTS7bptTW7dM+FvXT1kHTZbYbVMQEAAAAAaFIhtWaXJM2dO1cjR45UampqzUivWbNm6e23326295w2bZo8Hk/NLT09vdneC2hpfVMT9PJ/D9czNwxRl3axOlRcqQfe+EaXPPa5lm09ZHU8AAAAAACaVEiVXU8++aSmTJmiiy++WPn5+TVrdLVp00azZs1qtvd94IEHVFBQUHPbtWtXs70XYAXDMHRB3xR9dNco/eknfZXgdmjzviL9/NmV+u8XV2v7wWKrIwIAAAAA0CRCqux6/PHH9cwzz+gPf/iD7HZ7zfYhQ4bom2++abb3dblcSkhIqHMDIlH1ul2f3TdGN57TWXaboU82HdCFM5dq6rvfKb+00uqIAAAAAACclpAqu3bs2KFBgwbV2+5yuVRSUmJBIiAyJcZG6aEJZ+iju0bpx72T5QuYen7ZTo2evkTPL9shrz9gdUQAAAAAAE5JSJVdXbp00fr16+tt//DDD9WnT59Gn6e4uFjr16+vOdeOHTu0fv16ZWdnN1FSIDJ0T47TczcO1dybh6lXSrwKyrya+u5GjZu1VIs27VcIXawVAAAAAIBGCamrMU6ZMkWTJ09WeXm5TNPUqlWr9Oqrr2ratGl69tlnG32eNWvWaMyYMXXOK0mTJk3SCy+80NSxgbB3bo/2ev83bTVvzS7NWPi9th8s0c0vrtGPurfTH3/SR707MLUXAAAAABAeDDPEhm68/PLLeuihh7Rt2zZJUmpqqqZOnaqbb765xTIUFhbK4/GooKCA9bvQ6hSWe/WPxdv03Bc7VOkPyGZIE4dmaMoFPdU+3mV1PAAAAABAK9XYvsbysuudd97R+PHj5XQ662wvLS1VcXGxkpOTWzwTZRcg7cot1SMLNuv9b/ZKkuJcDk0e0103jewst9N+glcDAAAAANC0wqbsstvt2rdvn9q3by+73a69e/daUnDVRtkFHLF6Z67+/N5Gfb27QJKUlhitB8b30cX9OsgwDIvTAQAAAABai8b2NZYvUN++fXutWLFCkmSaJv94BkLM0M5Jeuv2kZpx9QB1SHBrd16ZJr+yTj+ds1wbduVbHQ8AAAAAgDosL7tuvfVWXXbZZbLb7TIMQx06dJDdbm/wBsAaNpuhK89K06f3jtZdY3so2mnXmqw8XTZ7me6et157C8qsjggAAAAAgKQQmMYoSZs3b9bWrVs1YcIEPf/882rTpk2Dx1122WUtkodpjMDx7Sso16MfbdYb63IkSW6nTb8c1U23ju6qmKiQusgrAAAAACBChM2aXbUXqJ86daruu+8+xcTEWBmJsgtopK935+vP723U6p15kqSUBJfuG9dbVw7qJJuNKckAAAAAgKYTNmUXC9QD4c00TS34dp+mLdikXbnB6Yz9Onn0p5/01bAuSRanAwAAAABEChaoB9AiDMPQxf066uO7R+t343srzuXQNzkFuvqp5brtX2uVfbjU6ogAAAAAgFbE8rKLBeqByOB22nXr6G5act95unZ4hmyGtODbfRo74zNN+2CTCsu9VkcEAAAAALQClk9jlFigHohEm/cV6n/f36TPfzgkSWobG6W7L+ipa4amy2G3vGcHAAAAAISZsFmzqzYWqAcii2maWrzlgB5+f5O2HyyRJPVMidMfL+mrUT3bW5wOAAAAABBOwrLsChWUXUDT8voDenlFlmYt+kH5pcHpjGN6tdcfLumj7snxFqcDAAAAAISDsCm7zjrrLC1atEiJiYkaNGjQcReoX7duXYtkouwCmkdBqVd/X/SDXlq+U76AKbvN0HXDM3TX2J5KjI2yOh4AAAAAIIQ1tq9xtGCmBl122WVyuVySpMsvv9zaMACalSfGqf+5tK+uOztDf/lgsz7ZtF8vLs/Sm1/l6Dfn99ANIzorysF6XgAAAACAU2f5yK5QxMguoGV8ufWQ/t97G7V5X5EkqUu7WD0wvrcu6Jty3FGeAAAAAIDWJ2ymMdZmmqbWrl2rnTt3yjAMdenS5YRTG5sDZRfQcvwBU6+t2aW/Lvxeh4orJEkjurbVH3/SR2ekeixOBwAAAAAIFWFXdi1evFg333yzsrKyVB2puvB67rnnNGrUqBbLQtkFtLziCp/+sXirnv1ihyp9ARmG9NPBabr3wl5KTnBbHQ8AAAAAYLHG9jUhsTjO1q1b9ZOf/ESdO3fWG2+8oU2bNmnjxo167bXXlJaWposvvljbt2+3OiaAZhTncuj+i3pr0ZTR+kn/jjJNaf6a3Trvr0v0xKc/qNzrtzoiAAAAACAMhMTIrjvuuEObNm3SokWL6u0zTVNjx45V37599fjjj7dIHkZ2AdZbm5Wr//feJm3YlS9J6tQmWvdf1EsTBqSynhcAAAAAtEJhNbJryZIluuuuuxrcZxiG7rrrLi1evLhlQwGw1ODMJL152zmaNXGgOnrcyskv053/Xq8rn/xS67LzrI4HAAAAAAhRIVF2ZWdnq1+/fsfcf+aZZyorK6sFEwEIBTabocsHddKn95ynKRf0VLTTrq+y83XlP77Ub179Sjn5ZVZHBAAAAACEmJAou4qLixUTE3PM/TExMSotLW3BRABCSXSUXb85v4eW3Heefjo4TYYhvbNhj3781yX660dbVFLhszoiAAAAACBEhMSaXTabTZ9++qmSkpIa3H/o0CFdcMEF8vtbZoFq1uwCQtu3OQX683sbtXJHriSpfbxL913YS1cNTpPdxnpeAAAAABCJGtvXhEzZZRiGGopSvd0wDMouADVM09RH3+3XtAWblHU4OPKzb8cE/eknfTWiW1uL0wEAAAAAmlpYlV2NXY8rMzOzmZMEUXYB4aPC59dLX2bpsU9/UFF5cDrjhX1T9PuL+6hzu1iL0wEAAAAAmkpYlV2hhrILCD+Hiys065Mf9MqqbPkDppx2Q5NGdNavz+8hT7TT6ngAAAAAgNNE2XUaKLuA8PXD/iI9/P4mffb9QUlSYoxTd1/QU9cOy5DDHhLX5AAAAAAAnALKrtNA2QWEvyVbDuh/39+kHw4US5K6J8fpD5f00ZheyRYnAwAAAACcCsqu00DZBUQGnz+gV1dla8bH3yuv1CtJGtWzvf54SR/1TIm3OB0AAAAA4GRQdp0Gyi4gshSUefXEpz/ohS93yus3ZTOknw3L0JQLeqptnMvqeAAAAACARqDsOg2UXUBk2nmoRNMWbNJH3+2XJMW7HLrjx91148jOcjnsFqcDAAAAABxP2JRdgwYNkmEYjTp23bp1zZwmiLILiGzLtx3Ww+9v1Hd7CiVJGUkxuvP8Hrqkf0e5nZReAAAAABCKGtvXOFowU4Muv/xyqyMAaGVGdGurd+74kV5ft1vTP9qi7NxS3fPaBv3vB5s0cWi6fj48Q2mJMVbHBAAAAACcAstHdoUiRnYBrUdJhU8vfLlT/1qRpb0F5ZIkmyH9uHeKbhiRqR91byebrXGjTwEAAAAAzSdspjGGIsouoPXx+QP6ZNMBzV2xU8u2Hq7Z3qVdrK47O1P/NThNnminhQkBAAAAoHULy7LL7/dr5syZmj9/vrKzs1VZWVlnf25ubovkoOwCWretB4r1rxVZen3tbhVV+CRJ0U67Lh+UquvP7qy+qfy5AAAAAAAtrbF9ja0FM53Q1KlTNWPGDE2cOFEFBQWaMmWKrrzyStlsNj300ENWxwPQSnRPjtNDE87Qit+fr4cvP1O9UuJV5vXr1VW7dPFjn+u/nvxSb6/PUaUvYHVUAAAAAMBRQmpkV7du3fTYY4/pkksuUXx8vNavX1+zbcWKFXrllVdaJAcjuwDUZpqmVu3I1UsrsvTRt/vkCwT/2GwX59LPhqXr2uEZ6uiJtjglAAAAAES2sJzGGBsbq02bNikjI0MdO3bU+++/r7POOkvbt2/XoEGDVFBQ0CI5KLsAHMv+wnK9uipbr6zM1oGiCkmS3Wbowr4pun5EpkZ0bSvDYEF7AAAAAGhqYTmNMS0tTXv37pUUHOW1cOFCSdLq1avlcrmsjAYAkqSUBLfuGttTy373Y82+9iwN75Ikf8DUgm/36dpnVurCmUv10vKdKir3Wh0VAAAAAFqlkBrZ9bvf/U4JCQn6/e9/r3nz5um6665T586dlZ2drbvvvluPPPJIi+RgZBeAk7FlX5FeWr5Tb36Vo9JKvyQpNsquK89K0/UjMtUzJd7ihAAAAAAQ/sJyGuPRli9fruXLl6tHjx669NJLW+x9KbsAnIrCcq/eWLtbc1dkadvBkprtZ3dN0g0jOuuCvily2kNqQC0AAAAAhI2IKLusQtkF4HSYpqkvtx3WS8t36uON+1W1nr1SEly6dlimfjYsXckJbmtDAgAAAECYCduy64cfftDixYt14MABBQKBOvv+53/+p0UyUHYBaCp78sv0ysps/Xt1tg4VV0qSHDZDF53ZQTeM6KyhnRNZ0B4AAAAAGiEsy65nnnlGt912m9q1a6cOHTrU+QegYRhat25di+Sg7ALQ1Cp8fn347T69tDxLa7Pyarb37hCv60dk6vKBnRTrcliYEAAAAABCW1iWXZmZmbr99tv129/+1tIclF0AmtN3ewo0d3mW3lqfo3JvcARrvMuhqwYHF7Tv1j7O4oQAAAAAEHrCsuxKSEjQ+vXr1bVrV0tzUHYBaAkFpV69tnaX/rUiSzsPl9Zs/1H3drp+RKbO750sBwvaAwAAAICkMC27br75Zg0dOlS33nqrpTkouwC0pEDA1OdbD2nu8p1atPmAqv9U7tQmWtcOz9A1Q9PVNs5lbUgAAAAAsFhYll3Tpk3TjBkzdMkll6hfv35yOp119v/mN79pkRyUXQCssiu3VC+vzNa81dnKK/VKkqLsNl3Sv6OuH5GpQeltWNAeAAAAQKsUlmVXly5djrnPMAxt3769RXJQdgGwWrnXr/e/3quXlu/Uht0FNdvPSE3QpBGddemAVEVH2S1MCAAAAAAtKyzLrlBB2QUglGzYla+Xlmfp3a/3qNIXXNDeE+3U1UPSdN3ZmcpsG2txQgAAAABofmFfdlXHsmK6DmUXgFCUW1Kp+WuCC9rvziuTJBmGNLpne90wIlOjeybLbmOKIwAAAIDI1Ni+JuQu8/XSSy+pX79+io6OVnR0tPr376+5c+daHQsALJcUG6VbR3fTZ/eN0T8nDdHonu1lmtKSLQf1ixfW6Ly/LtZTn21TXkml1VEBAAAAwDIhVXbNmDFDt912my6++GLNnz9f8+fP10UXXaRbb71VM2fOPOnzzZ49W507d5bb7dbw4cO1atWqZkgNAC3LbjN0fp8UvfiLYVpy73n67x91UYLboV25ZZq2YLPOnrZI9762QV/vzrc6KgAAAAC0uJCaxtilSxdNnTpVN9xwQ53tL774oh566CHt2LGj0eeaN2+ebrjhBs2ZM0fDhw/XrFmz9Nprr2nLli1KTk4+7muZxggg3JRV+vXOhhy9tDxL3+0prNk+IL2Nbjg7U5f07yi3kwXtAQAAAISvsFyzy+1269tvv1X37t3rbP/hhx/Ur18/lZeXN/pcw4cP19ChQ/XEE09IkgKBgNLT0/XrX/9av/vd7477WsouAOHKNE2ty87X3OU79cE3+1TpDy5onxQbpauHpOvnwzOUnhRjcUoAAAAAOHlhuWZX9+7dNX/+/Hrb582bpx49ejT6PJWVlVq7dq3Gjh1bs81ms2ns2LFavnx5veMrKipUWFhY5wYA4cgwDA3OTNSsawbpywd+rPvG9VKqx63ckkrN+WybRk1frP9+cbU++/6gAoGQ+VkHAAAAADQZh9UBaps6daomTpyopUuXauTIkZKkZcuWadGiRQ2WYMdy6NAh+f1+paSk1NmekpKizZs31zt+2rRpmjp16umFB4AQ0y7OpcljuutXo7pq0eYDmrs8S19sPaRPNh3QJ5sOqEu7WP18eIZ+Ojhdnhin1XEBAAAAoEmE1Miuq666SitXrlS7du301ltv6a233lK7du20atUqXXHFFc32vg888IAKCgpqbrt27Wq29wKAluaw2zTujA76138P1ydTRuvGczor3uXQjkMlevj9TRo+7RM98MbX2riHUa0AAAAAwl9IrdnVVCorKxUTE6P//Oc/uvzyy2u2T5o0Sfn5+Xr77beP+3rW7AIQ6UoqfHrzqxzNXZ6lLfuLarYPyUzU9SMyNf7MjopyhNTPQwAAAAC0cmGzQH1hYWFNwBOtlXUyxdPw4cM1bNgwPf7445KCC9RnZGTojjvuYIF6AKhimqZW7cjVSyuy9NG3++SrWserXVyUfjYsQ9cOz1BHT7TFKQEAAAAgjMouu92uvXv3Kjk5WTabTYZh1DvGNE0ZhiG/39/o886bN0+TJk3SU089pWHDhmnWrFmaP3++Nm/eXG8tr6NRdgFojQ4UluuVVdl6ZWW2DhRVSJLsNkNj+yTrvF7JGpKZqG7t42Sz1f9zGgAAAACaW2P7GssXqP/000+VlJQkSVq8eHGTnXfixIk6ePCg/ud//kf79u3TwIED9eGHH56w6AKA1io5wa27xvbU5DHdtfC7/Xpp+U6t3JGrj77br4++2y9JSnA7dFZmooZkJuqszEQNTG+jmCjL/1MCAAAAADUsH9lVW3Z2ttLT0+uN7jJNU7t27VJGRkaL5GBkFwAEbdlXpHc37NGarFxt2FWgMm/dEbZ2m6G+HRM0ODOx5pbahmmPAAAAAJpe2ExjrK32lMbaDh8+rOTk5JOaxng6KLsAoD6vP6BNewu1Niuv5ra3oLzecaked83or8GZSerTMV4OO4vdAwAAADg9YTONsbbqtbmOVlxcLLfbbUEiAEA1p92m/mlt1D+tjW4a2UWStCe/TGuy8rQuK09rsnK1aW+R9hSUa8/Xe/Xe13slSdFOuwamtwmO/OqcqLPSE+WJcVr5UQAAAABEsJAou6ZMmSJJMgxDf/rTnxQTE1Ozz+/3a+XKlRo4cKBF6QAAx5LaJloT2kRrwoBUSVJJhU8bdudr7c48rc0OlmCF5T4t335Yy7cfrnldz5Q4Dc5M1FkZiRrSOUmd28Y0+MMOAAAAADhZIVF2ffXVV5KCI7u++eYbRUVF1eyLiorSgAEDdO+991oVDwDQSLEuh87p1k7ndGsnSQoETG09WKy1WXlaszNP67LztONQib7fX6zv9xfr1VW7JEltY6N0VtWaX0MyE3VmJ4/cTruVHwUAAABAmAqpNbtuuukm/f3vf7d8nSzW7AKA5nOouELraq379XVOgSp9gTrHOO2GzuzkqVr3K3jlx+R4prMDAAAArVlYLlBfUFAgv9+vpKSkOttzc3PlcDharHii7AKAllPh8+vbnMKadb/WZuXpUHFlveMykmI0pKr4GtI5UT2S42W3MfURAAAAaC3CsuwaP368Lr30Ut1+++11ts+ZM0fvvPOOPvjggxbJQdkFANYxTVPZuaV1rvq4ZX+Rjv6vVbzLoYEZbTQkM0mDMxM1MKON4lwhMTsfAAAAQDMIy7IrKSlJy5YtU58+feps37x5s0aOHKnDhw8f45VNi7ILAEJLYblXX2Xna23VlR+/ys5TSaW/zjE2Q+rdISG47lfn4OL3aYnRLHwPAAAARIjG9jUh9SPwiooK+Xy+etu9Xq/KysosSAQACAUJbqdG92yv0T3bS5J8/oC27C+qGfm1ZmeecvLLtHFvoTbuLdTcFVmSpJQEV52rPvbtmKAoh83KjwIAAACgmYXUyK4xY8bozDPP1OOPP15n++TJk/X111/r888/b5EcjOwCgPCzr6Bc67KDxdfa7Dx9l1MgX6Duf+JcDpsGpLepuerjWRmJSoyNOsYZAQAAAISSsJzGuGzZMo0dO1ZDhw7V+eefL0latGiRVq9erYULF+rcc89tkRyUXQAQ/soq/fp6d77WVE19XJudp/xSb73jurWP1eCqqz4OzkxSt/axTH0EAAAAQlBYll2StH79ek2fPl3r169XdHS0+vfvrwceeEA9evRosQyUXQAQeUzT1LaDJXWu+rjtYEm949rEONWvk0ed28Yqs22MMtvGqnPbGKUnxcjttFuQHAAAAIAUxmVXKKDsAoDWIa+kUuuyq9b9ysrThl35qvAFjnl8R49bGUkx6tw2VhltY2oVYjGKdztbMDkAAADQ+oR92VVeXq7Kyso621qqeKLsAoDWqdIX0Ma9hfp+X5F2Hi5RVm6psg6XKOtQqYoq6l9Apbak2ChlVhVgGUkx6twuRhlJwVFhSbFRTI0EAAAATlNYll2lpaW6//77NX/+fB0+fLjefr/f38Crmh5lFwCgNtM0lVfqDRZfh0u183CJsqvvc0t1qLjyuK+PczlqRoBlto1VZlLV9Mh2MUqJd8tmowgDAAAATqSxfY2jBTOd0H333afFixfrySef1PXXX6/Zs2crJydHTz31lB555BGr4wEAWinDMJQUG6Wk2CgNykist7+o3Kvs3FJlHa6+ldTc7y0sV3GFT9/tKdR3ewrrvTbKYasqv6qKsFqFWKfEaDnttpb4iAAAAEDECKmRXRkZGXrppZd03nnnKSEhQevWrVP37t01d+5cvfrqq/rggw9aJAcjuwAATaXc69fuvFLtPFR6ZFpkVRG2O69MvsCx/zNstxnq1Ca6ZlRYcI2wYCGWwYL5AAAAaGXCcmRXbm6uunbtKim4Pldubq4k6Uc/+pFuu+02K6MBAHBK3E67uifHq3tyfL19Pn9Ae/LLlZVbop2HS5V16Mg6Ydm5pSr3BpSdW6rs3FJ9/kP9c3dIcNedHtn2yOL5CSyYDwAAgFYqpMqurl27aseOHcrIyFDv3r01f/58DRs2TO+++67atGljdTwAAJqUw25TRtsYZbSN0bk96u4LBEwdKKo4MhKsqhCrXiusqNynfYXl2ldYrpU7cuudOyk2qurKkTHKaBurDglutYuLUvt4l9rFudQ+3sXIMAAAAESkkJrGOHPmTNntdv3mN7/RJ598oksvvVSmacrr9WrGjBm68847WyQH0xgBAKHMNE3ll3qDV4ysvU5Y1aiwEy2YXy3e5ahTflWXYXW3BW9RDtYOAwAAgLXC8mqMR8vKytLatWvVvXt39e/fv8Xel7ILABDOiit8wamQh0uDo8FyS3WwqFwHiyt1qKhCB4srVOkLnNQ5PdHOWoVY/VFi7avuk2KjWFQfAAAAzSLsyi6v16uLLrpIc+bMUY8ePU78gmZE2QUAiGSmaaqw3KdDxRU6WBS8VT8+cl9Z8/x4i+gfzTCkxJgotY9zqV181X1c/dFi1cWY3WY04ycFAABAJAm7BeqdTqe+/vprq2NEnpLD0sI/SjFJUnRi1f1Rj2OSJGe01UkBAC3EMAx5op3yRDvVrX3ccY8NBEwVlHmPFGMNlGHV2w8XVyhgSrkllcotqdSW/cfPYTOkpNha0yePKsPax7uUGBOlhGiHEqKdiotyyEY5BgAAgBMImZFdknT33XfL5XLpkUcesTRHRI3s2vetNGfkiY9zRB8pxI4uwqITaz2udR/dRrKxuDEAIMgfMJVXWnncUWLV94dLKnWyfwOxGVK826mEaIc80U4luIM3T3RwW4LbKU9M1fbax1QVey6HTYZBWQYAABCuwm5klyT5fD4999xz+uSTTzR48GDFxsbW2T9jxgyLkoWx2HbS+Q9KZblSaZ5Ullf1ODd4X5YnBXySr0wqzAneTobb00ARllj38dGjyaLigvNcAAARxW4zaha0793h+Mf6/AHlllQec6RY9X1eqVeFZV5V+gMKmFJBmVcFZV7tUtlJ54uy25RQuxiLDhZhCW5HrcdHlWdV2+PdDtYiAwAACBMhNbJrzJgxx9xnGIY+/fTTFskRUSO7TsQ0pYrCYOlVXYCV5h0pwmq21b7PlyoKTv09bc4GyrGjR5M1UJQ5oprsYwMAwku516/CMq8Ky70qKPPVPC6sKr8Ky321HntVWOar9dirk1h27Jhio+wNFmIJ0U7FRNkV63IE76Mcio6yK9ZlV0yUQ7FRDsW47IqJqn5ul4PiDAAA4KSF1QL127dvV5cuXUJmakGrKrtOld8bLL2OHiV2dDlWll93m7/i1N8zKr6qAGtoWuXR26qeuxIkG/+gAIDWLBAwVVLpU2G5TwWlJy7Jji7SSir9TZ4pymFTbFX5FRNlV4zLUef5kaIsuK92UVZ9bHRVsRbjqirYnHbWNAMAABEtrMouu92uvXv3Kjk5WZI0ceJEPfbYY0pJSbEkD2VXMzFNyVtaVX4dNZ3yeKPJyvIlneJvU8MeXFuswXLsOGuTsWA/AKCKzx+oKcWCI8uqSrGax16VVvpVWumruverpCL4uKTSp7Kq5yWVfvmbYojZcUQ7jxRl0U67XE6b3I6qe6ddLkfw3l17u8Nes81V5xi73I7gNvfRxzmCx1GuAQCAlhRWZZfNZtO+fftqyq74+Hht2LBBXbt2tSQPZVeICfil8oIGirCjC7OjRpN5S079PWsW7K9aiP/oKZYNjSZjwX4AwHGYpqlKfyBYflX6VVpVgJXWKsZqF2XVBVpJRdV9pV9l9Z4HX2fV3+aiHLYGC7HapZrLYZfTbijKYZPTblOUw6aoqntn7fujjmnwWLtNUQ6j4WPslG8AAES6sFygHmiQzR4slGKSpLbdGv86b/lRhVhDo8macsF+o2rB/gbWHzve2mRRsSzYDwCtgGEYVSOi7GoT03TnNU1T5d5AsCyr8KvU61NJhU/l3oDKvX6VewOq8PmPPPf5VeENHLn3+lXu9avCd+T4ho6p3u+rNTqt0hdQpS8glfua7gOdBrvNUJTdVlOuRdltcthtctgNOW3Be4fdJqfNCG6z2+SwGbLbgq+pve/I47qvr36No+p9HDXbg4+r7+12o+rcVc9ttZ7X7LPJbhhHHVv/NZR4AACcnJAouwzDqLdeV6is34Uw5nRLzo5SQsfGv+Z4C/YfczRZftWC/aZUnh+85e1o/Hvao+pfwbKhUswZXasUq7o/5vPGHHMar6nzsKnPezp5T3SO45zXZg9Oe7U5guu82Ry1njtY+w1AyDIMQ9FVa3gprvnfz+cPqNwXUIXXr3Jf/bKs4qhSrcIXkNd/5L6y6t7rN+ttq/QFVFn7uT8gr8+s+3p/QN6q47z+ukPa/AFTZQG/yrzN/z20JMPQMYuzOjcjWIzZjSPbgs8lh80mmy1YCNqMI6+3GUeOcxx9DvuRc9kMQzZDNcfaDMle9Xf46vcyqrYFHwfft+Zx1TkN40iG6vexGao6Z9Vj48h72KpeU5PTCP6etxlHzhc8vuq+9jmq91e/TrVfX/e+5jxV91Ld54bBv08AIJyERNllmqZuvPFGuVwuSVJ5ebluvfVWxcbG1jnujTfesCIeWhOjanSW2yMldm786xq7YP/Ro8n8FZK/UireH7whtNUpwOxHFWQNPXdIhq1WYXb0tmM9tx85X0PvWWfb0cXciXLU2l7vmIbet+qm1vYX/AbmhB1zntgxtp/svLKTOf/JZjklJ/FrftL/AGzOc0eAJvvMzfPdORTs1Gp6NXvV7YRZGnXgSTElVfoD8vlN+fzVZViwHPMFqoqygOT3B+Qzg8f4/Kb8pimv35Q/EJAvYAa3BQJHjg2Y8gUC8vtN+QKqdVz1vqrXmEfe+8g5g/sD1fdmoOocpvxV7+0PBN/Pb1Ztr74d78+NQNWtwe9W8lfdUFdTT/E1apdkqi7Baj+uuldVSaYj2yRVlWtH9qnmmCPnOLKv+j3rvl9wVwPvXfUaW9Xj6nKu+pjqUzaUt/p9juw/8rzmHY0jr69zrurzH+M8OupctT9H7XzV5zr6/VQ7/1GfrSZI9XlV+z2P+gxH56zzWYw6eVUn75H3b/C1DfzMtt5nbCDr0a+t/p7qHH/U43r76g0YUX21NhoNHFfnV7mhzCd4TYMHH/3aYz1u4Ltr6Pwn+s9ive/heAEaOP/x3uNYeY/e25jXH+/9G3xtvd8fDZ8juV2SumdmNOqcrUlIlF2TJk2q8/y6666zKAlwiuxOKa598NZYdRbsP7oca2DBfm959QuPvP5knp/Ka+r8DdGsc9fy73v0MSfYf7LvawaCU1iPJ+CT5Du9q4oCAJqMIclVdQtbtqobIk8DfzUBgKa2st2V6n7H81bHCDkhUXY9/zy/MGiFDCO4XldUrNQm3eo0qBaoKr1Mf/A+4Gtgmz94q/Pcd6Qwq3nuP/K43jb/Uec80Xsc631PJketczWY4zi5WuPomoZ+6nbM7+FYP8472eObMctJOYl/mZ30sImTOfdJntqqk1p/rZ8GhFCmkPx+mlMEft4I/DVsik9knsqZGnH46WU7+oeCp3WWkGTW/E+jjw5JoZssxITYF9VQHMPikNHusP6RT7MJibILAEKGzSbZoqxOAQAAmlFT/FigNf4YCEDo6W91gBDFoGkAAAAAAABEDMouAAAAAAAARAzKLgAAAAAAAEQMyi4AAAAAAABEDMouAAAAAAAARAzKLgAAAAAAAEQMh9UBQpFpmpKkwsJCi5MAAAAAAABAOtLTVPc2x0LZ1YCioiJJUnp6usVJAAAAAAAAUFtRUZE8Hs8x9xvmieqwVigQCGjPnj2Kj4+XYRhWx2k1hg4dqtWrV1sdA60Uv//CH7+GzYvv99ha43cTiZ85XD9TqOcuLCxUenq6du3apYSEBKvjAAhDof7nHFqWaZoqKipSamqqbLZjr8zFyK4G2Gw2paWlWR2j1bHb7fwlCJbh91/449ewefH9Hltr/G4i8TOH62cKl9wJCQlhkRNA6AmXP+fQco43oqsaC9QjZEyePNnqCGjF+P0X/vg1bF58v8fWGr+bSPzM4fqZwjU3ADQWf87hVDCNEQAAAECzKCwslMfjUUFBASMzAAAthpFdAAAAAJqFy+XSgw8+KJfLZXUUAEArwsguAAAAAAAARAxGdgEAAAAAACBiUHYBAAAAAAAgYlB2AQAAAAAAIGJQdgEAAAAAACBiUHYBAAAAAAAgYlB2AQAAALDEe++9p169eqlHjx569tlnrY4DAIgQhmmaptUhAAAAALQuPp9Pffv21eLFi+XxeDR48GB9+eWXatu2rdXRAABhjpFdAAAAAFrcqlWrdMYZZ6hTp06Ki4vT+PHjtXDhQqtjAQAiAGUXAAAAgJO2dOlSXXrppUpNTZVhGHrrrbfqHTN79mx17txZbrdbw4cP16pVq2r27dmzR506dap53qlTJ+Xk5LREdABAhKPsAgAAAHDSSkpKNGDAAM2ePbvB/fPmzdOUKVP04IMPat26dRowYIDGjRunAwcOtHBSAEBrQ9kFAAAA4KSNHz9eDz/8sK644ooG98+YMUO33HKLbrrpJvXt21dz5sxRTEyMnnvuOUlSampqnZFcOTk5Sk1NbZHsAIDIRtkFAAAAoElVVlZq7dq1Gjt2bM02m82msWPHavny5ZKkYcOG6dtvv1VOTo6Ki4u1YMECjRs3zqrIAIAI4rA6AAAAAIDIcujQIfn9fqWkpNTZnpKSos2bN0uSHA6H/va3v2nMmDEKBAK6//77uRIjAKBJUHYBAAAAsMSECRM0YcIEq2MAACIM0xgBAAAANKl27drJbrdr//79dbbv379fHTp0sCgVAKC1oOwCAAAA0KSioqI0ePBgLVq0qGZbIBDQokWLNGLECAuTAQBaA6YxAgAAADhpxcXF2rp1a83zHTt2aP369UpKSlJGRoamTJmiSZMmaciQIRo2bJhmzZqlkpIS3XTTTRamBgC0BoZpmqbVIQAAAACElyVLlmjMmDH1tk+aNEkvvPCCJOmJJ57Q9OnTtW/fPg0cOFCPPfaYhg8f3sJJAQCtDWUXAAAAAAAAIgZrdgEAAAAAACBiUHYBAAAAAAAgYlB2AQAAAAAAIGJQdgEAAAAAACBiUHYBAAAAAAAgYlB2AQAAAAAAIGJQdgEAAAAAACBiUHYBAAAAAAAgYlB2AQAAAAAAIGJQdgEAAAAAACBiUHYBAAAAAAAgYlB2AQAAAAAAIGL8fwocmUCZ9eN5AAAAAElFTkSuQmCC", "text/plain": [ "
" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "fig, ax = plt.subplots(\n", " 2,\n", " 1,\n", " sharex=True,\n", " figsize=(12, 6),\n", " constrained_layout=True,\n", " gridspec_kw={\"height_ratios\": [0.65, 0.35]},\n", ")\n", "ax[0].scatter(sense_moderate.k1d, sense1d, label=\"Modern 21cmSense\", marker=\"x\", color=\"C0\")\n", "ax[0].scatter(memo_data[\"ks\"], memo_data[\"errs\"], label=\"memo\", color=\"C0\", lw=1, facecolor=\"none\")\n", "\n", "ax[0].scatter(\n", " sense_moderate.k1d, sense1d_t, label=\"Modern 21cmSense (Thermal Only)\", marker=\"x\", color=\"C1\"\n", ")\n", "ax[0].scatter(\n", " memo_data[\"ks\"], memo_data[\"T_errs\"], label=\"memo (thermal Only)\", color=\"C1\", facecolor=\"none\"\n", ")\n", "\n", "ax[0].set_yscale(\"log\")\n", "ax[0].set_xscale(\"log\")\n", "ax[0].set_ylabel(r\"$\\Delta^2$ [mK$^2$]\")\n", "ax[1].set_ylabel(\"Fractional Difference (%)\")\n", "ax[1].plot(\n", " memo_data[\"ks\"],\n", " 100 * (sense1d[: len(memo_data[\"ks\"])].value - memo_data[\"errs\"]) / memo_data[\"errs\"],\n", ")\n", "ax[1].plot(\n", " memo_data[\"ks\"],\n", " 100 * (sense1d_t[: len(memo_data[\"ks\"])].value - memo_data[\"T_errs\"]) / memo_data[\"T_errs\"],\n", ")\n", "ax[0].legend();" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Tests" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The following are just assertions that the code reproduces the memo data within some\n", "tolerances, as can be seen already by the above plots. This is just for CI testing \n", "purposes." ] }, { "cell_type": "code", "execution_count": 15, "metadata": {}, "outputs": [], "source": [ "mask = ~np.isinf(memo_data[\"errs\"])\n", "\n", "assert np.allclose(\n", " sense1d[: len(memo_data[\"ks\"])][mask][:-1].value, memo_data[\"errs\"][mask][:-1], rtol=1e-1\n", ")\n", "assert np.allclose(\n", " sense1d_t[: len(memo_data[\"ks\"])][mask][:-1].value, memo_data[\"T_errs\"][mask][:-1], rtol=1e-2\n", ")" ] } ], "metadata": { "kernelspec": { "display_name": "Python 3.10.6 ('sense')", "language": "python", "name": "python3" }, "language_info": { "codemirror_mode": { "name": "ipython", "version": 3 }, "file_extension": ".py", "mimetype": "text/x-python", "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", "version": "3.12.4" }, "orig_nbformat": 4, "vscode": { "interpreter": { "hash": "6fa29f5e4f57cd8c712524a97e07d2162cf815d6c48daa40d8ca34da05a68238" } } }, "nbformat": 4, "nbformat_minor": 2 }