From fd6abdbdfb15ef0636a5797dc91fc072bc5f2c96 Mon Sep 17 00:00:00 2001 From: thea-mckenna Date: Thu, 20 Jul 2023 15:45:09 -0400 Subject: [PATCH] nested sampler tutorial! --- docs/tutorials/dynesty_tutorial.ipynb | 185 +++++++++++++++----------- 1 file changed, 110 insertions(+), 75 deletions(-) diff --git a/docs/tutorials/dynesty_tutorial.ipynb b/docs/tutorials/dynesty_tutorial.ipynb index 3a415af6..74f2ebf1 100644 --- a/docs/tutorials/dynesty_tutorial.ipynb +++ b/docs/tutorials/dynesty_tutorial.ipynb @@ -1,132 +1,165 @@ { "cells": [ { + "attachments": {}, "cell_type": "markdown", "metadata": {}, "source": [ - "# Nested Sampler Introduction" - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "***Mention orbit fraction and amount of constraint as motivation for having nested sampling as an option once tests are done and we have more information***\n", + "# Nested Sampler Introduction\n", + "Here, we will explain how to sample an orbit posterior using nested sampling from the `dynesty` package. An advantage of nested sampling is that it computes the evidence and posterior at the same time. To compute the evidence, Dynesty first uses live points to make iso-likelihood contours with probability greater than that of the previous contour. It then takes the integral of these contour shells and calculates prior “volumes” of these shells. As each contour is made, the live point with the least probability is removed (it “dies”), and it is replaced with a new live point sampled from the prior, which must have a higher probability. `dynesty` then uses these “dead” points to approximate the evidence integral by weighting each point and summing them. To estimate the posterior, dynesty then uses the calculated “dead” point weighting and evidence to get its importance weight, or the probability of the parameter set. For more info on how nested sampling works in `dynesty`, go to their [website](https://dynesty.readthedocs.io/en/stable/index.html).\n", "\n", - "Here, we will explain how to sample an orbit posterior using nested sampling from the dynesty package. An advantage of nested sampling is that it computes the evidence and posterior at the same time. To compute the evidence, Dynesty first uses live points to make iso-likelihood contours with probability greater than that of the previous contour. It then takes the integral of these contour shells and calculates prior “volumes” of these shells. As each contour is made, the live point with the least probability is removed (it “dies”), and it is replaced with a new live point sampled from the prior, which must have a higher probability. Dynesty then uses these “dead” points to approximate the evidence integral by weighting each point and summing them. To estimate the posterior, dynesty then uses the calculated “dead” point weighting and evidence to get its importance weight, or the probability of the parameter set.\n" + "In this tutorial, we will demonstrate how to generate a synthetic data set with a user-controlled orbit fraction (fraction of orbit covered by the synthetic data points) to feed into the nested sampler. We will also show plots using the sampler results from `dynesty` and save these results.\n" ] }, { + "attachments": {}, "cell_type": "markdown", "metadata": {}, "source": [ - "# Read in Data and Set up Sampler" - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "desscription here" + "# Generate Synthetic Data\n", + "We generate synthetic data using the `generate_synthetic_data` function in `orbitize.system` to feed into the nested sampler. This function sets ground truth values for all of the parameters, the desired fractional orbit coverage, the number of observations, and the desired observational uncertainty. The observations are uniformly spaced synthetic data points covering the years 2001-2003 such that the orbital period and semi-major axis are determined by selecting the orbital fraction. A three year period is used for all orbits, so that for smaller orbit fraction, the orbit has a larger semi-major axis. The function then adds Gaussian noise to the observations to make them more realistic. This approach simulates similar observing strategies for orbits of varying lengths." ] }, { "cell_type": "code", - "execution_count": 2, + "execution_count": 3, "metadata": {}, - "outputs": [], + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + " epoch object quant1 ... quant_type instrument\n", + "------------------ ------ ------------------- ... ---------- ----------\n", + " 51550.0 1 11.690679585649285 ... radec defrd\n", + " 51587.93103448276 1 36.7951190350705 ... radec defrd\n", + "51625.862068965514 1 53.44625230418348 ... radec defrd\n", + "51663.793103448275 1 57.74172150011931 ... radec defrd\n", + "51701.724137931036 1 43.476583812616035 ... radec defrd\n", + " 51739.65517241379 1 17.60198047834588 ... radec defrd\n", + " 51777.58620689655 1 -10.554719202310748 ... radec defrd\n", + " 51815.51724137931 1 -42.902152497869935 ... radec defrd\n", + " 51853.44827586207 1 -65.65070346148586 ... radec defrd\n", + "51891.379310344826 1 -92.15759521127484 ... radec defrd\n", + " ... ... ... ... ... ...\n", + " 52270.68965517241 1 -173.0917550491587 ... radec defrd\n", + "52308.620689655174 1 -170.18957371562144 ... radec defrd\n", + " 52346.55172413793 1 -161.916220629315 ... radec defrd\n", + " 52384.48275862069 1 -153.65792545220265 ... radec defrd\n", + " 52422.41379310345 1 -145.00885130393408 ... radec defrd\n", + " 52460.34482758621 1 -129.25073055046082 ... radec defrd\n", + "52498.275862068964 1 -110.71976056947764 ... radec defrd\n", + "52536.206896551725 1 -92.64517427318272 ... radec defrd\n", + "52574.137931034486 1 -74.64835433878045 ... radec defrd\n", + " 52612.06896551724 1 -53.65409077389522 ... radec defrd\n", + " 52650.0 1 -27.50155509020929 ... radec defrd\n", + "Length = 30 rows\n" + ] + } + ], "source": [ "import orbitize\n", "from orbitize import read_input, system, priors, sampler\n", "import matplotlib.pyplot as plt\n", + "from orbitize.system import generate_synthetic_data\n", + "import numpy as np\n", "\n", + "# generate data\n", + "plx = 60.0 # [mas]\n", + "mtot = 1.2 # [M_sol]\n", + "orbit_frac = 95.\n", + "data_table, sma = generate_synthetic_data(orbit_frac, mtot, plx, num_obs=30)\n", "\n", - "data_table = read_input.read_file('{}/GJ504.csv'.format(orbitize.DATADIR))\n", - "\n", - "# number of secondary bodies in system\n", - "num_planets = 1\n", - "# total mass & error [msol]\n", - "total_mass = 1.22\n", - "mass_err = 0.08\n", - "# parallax & error[mas]\n", - "plx = 56.95\n", - "plx_err = 0\n", - "sys = system.System(\n", - " num_planets, data_table, total_mass,\n", - " plx, mass_err=mass_err, plx_err=plx_err\n", - ")\n", - "# alias for convenience\n", + "# initialize orbitize System object\n", + "num_secondary_bodies = 1\n", + "sys = system.System(num_secondary_bodies, data_table, mtot, plx)\n", + "print(data_table)\n", "lab = sys.param_idx\n", "\n", - "mu = 0.2\n", - "sigma = 0.05\n", - "\n", - "sys.sys_priors[lab['ecc1']] = priors.GaussianPrior(mu, sigma)\n", - "sys.sys_priors[lab['inc1']] = 2.5" - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "# Running the Dynesty Sampler" + "# other ground truth values\n", + "ecc = 0.5 \n", + "inc = np.pi/4 \n", + "aop = np.pi/4\n", + "pan = np.pi/4\n", + "tau = 0.8" ] }, { + "attachments": {}, "cell_type": "markdown", "metadata": {}, "source": [ - "description" + "# Running the Dynesty Sampler\n", + "`dynesty's` nested sampler has multiple prior bound shape options, as well as a static and dynamic nested sampler, more details of which can be found [here](https://dynesty.readthedocs.io/). These can also be configured in `run_sampler` keywords. Here we stick with the default options of a static sampler and 'multi' prior bound shape." ] }, { "cell_type": "code", - "execution_count": 1, + "execution_count": 5, "metadata": {}, "outputs": [ { - "ename": "NameError", - "evalue": "name 'sampler' is not defined", - "output_type": "error", - "traceback": [ - "\u001b[0;31m---------------------------------------------------------------------------\u001b[0m", - "\u001b[0;31mNameError\u001b[0m Traceback (most recent call last)", - "Cell \u001b[0;32mIn[1], line 1\u001b[0m\n\u001b[0;32m----> 1\u001b[0m nested_sampler \u001b[38;5;241m=\u001b[39m \u001b[43msampler\u001b[49m\u001b[38;5;241m.\u001b[39mNestedSampler(sys)\n\u001b[1;32m 3\u001b[0m \u001b[38;5;66;03m# number of orbits to accept\u001b[39;00m\n\u001b[1;32m 4\u001b[0m n_orbs \u001b[38;5;241m=\u001b[39m \u001b[38;5;241m500\u001b[39m\n", - "\u001b[0;31mNameError\u001b[0m: name 'sampler' is not defined" + "name": "stderr", + "output_type": "stream", + "text": [ + "iter: 16803 | batch: 8 | bound: 6 | nc: 1 | ncall: 66312 | eff(%): 25.066 | loglstar: -130.534 < -126.270 < -126.777 | logz: -137.282 +/- 0.082 | stop: 0.902 " ] } ], "source": [ - "nested_sampler = sampler.NestedSampler(sys)\n", + "# fix some parameters for faster convergence\n", + "sys.sys_priors[lab['sma1']] = sma\n", + "sys.sys_priors[lab['aop1']] = np.pi/4 \n", + "sys.sys_priors[lab['pan1']] = np.pi/4\n", + "sys.sys_priors[lab['tau1']] = 0.8 \n", + "sys.sys_priors[lab['plx']] = plx\n", + "sys.sys_priors[lab['mtot']] = mtot\n", "\n", - "# number of orbits to accept\n", - "n_orbs = 500\n", - "\n", - "_ = nested_sampler.run_sampler(n_orbs, static = False, bound = 'single')\n", - "nested_sampler.results.save_results('test2.hdf5')" - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "# Plotting" + "# run the nested sampler using 2 threads for parallelization\n", + "threads = 2\n", + "nested_sampler = sampler.NestedSampler(sys)\n", + "samples, num_iter = nested_sampler.run_sampler(num_threads = threads)" ] }, { + "attachments": {}, "cell_type": "markdown", "metadata": {}, "source": [ - "description" + "# Plotting\n", + "Full plotting capabilities of orbitize (including advanced corner and orbit plots) can be found in the [orbitize plotting tutorial](https://orbitize.readthedocs.io/en/latest/tutorials/Plotting_tutorial.html).\n", + "\n", + "Here, we will plot the posterior samples from dynesty against the ground truth variables used to generate the synthetic data set." ] }, { "cell_type": "code", - "execution_count": null, + "execution_count": 6, "metadata": {}, - "outputs": [], + "outputs": [ + { + "data": { + "image/png": "iVBORw0KGgoAAAANSUhEUgAAC6UAAAHvCAYAAAAWvpggAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8qNh9FAAAACXBIWXMAAA9hAAAPYQGoP6dpAACwlElEQVR4nOzdeZxd8/0/8NdNIpOELI0lkciCxhIigtiixLe2UNRSSxFLF62lSFualjZaBN8qRfHVqqXWtgj9am0liSWUEvtWQlJEfIlEEhKS+/vDL1MjM5GZzMydnDyfj8d9PHLPOZ9z3/fMOa/PNd73TKlcLpcDAAAAAAAAAAAAAAAN0KrSBQAAAAAAAAAAAAAAsOzSlA4AAAAAAAAAAAAAQINpSgcAAAAAAAAAAAAAoME0pQMAAAAAAAAAAAAA0GCa0gEAAAAAAAAAAAAAaDBN6QAAAAAAAAAAAAAANJimdAAAAAAAAAAAAAAAGkxTOgAAAAAAAAAAAAAADaYpHQAAAAAAAAAAAACABtOUDgAAAAAtxKuvvppSqZS+fftWupQW4YknnshXvvKVdO3aNa1atUqpVMrYsWMrXdYyadSoUSmVShk1alSlSyHJ0KFDnc9LwHnb9A477LCUSqVcccUVlS6lSZRKpZRKpUqX0WBF//kAAAAAUCya0gEAAABosAceeCDf/va3s95666Vz586pqqpKz54985WvfCW/+93vMnv27EqXyDJq2rRp2X777XPbbbelQ4cO2WqrrTJkyJB07tz5c8f27du3uhFx4aN9+/ZZe+21c8QRR+SZZ55phndQ03vvvZdRo0blvPPOa/bXbsnuuOOO7L333unZs2fatm2bLl26ZN11183uu++ec845Jy+//HKlS2QZ9tFHH+WKK67IXnvtlT59+qRDhw7p0KFD+vTpkz322CO/+c1v8vbbb1e6zGXGwi8J1PfRVF++GDNmTEaNGpWJEyc2yf4/68c//nFKpVKGDBmyRNvfdNNNKZVK6dq1a+bNm9fE1QEAAABA5bWpdAEAAAAALHvmzJmTww8/PH/84x+TJO3atcvaa6+d9u3b5/XXX89tt92W2267LT/96U9zxx13ZMCAARWueNmwwgorZN11103Pnj0rXUrFXX/99Zk+fXr23HPP3HTTTWnVqv731+jXr19WW221JJ80hb/00ku5/PLLc+211+ZPf/pTdt9998Yuu07vvfdeTj311PTp0yfHH398s73uQqusskrWXXfdrLLKKs3+2nU55phj8pvf/CZJsuKKK6Zfv37p0KFDXnvttfzv//5v/vd//zdvvvlmfvnLX1a4UpZFjz32WL72ta/llVdeSZJ07do166yzTlq3bp3XX389f/nLX/KXv/wlJ554Yi688MIcfvjhFa645evdu3etDdlPPfVUZs6cmV69eqV3796LrF+SLxM1xJgxY3LllVemb9++2XjjjZvkNT5t+PDhGT16dB588MG88sorWWuttRa7/dVXX50k2W+//dK2bdsmrw8AAAAAKk1TOgAAAAD18tFHH2WnnXbKAw88kO7du+ess87K1772tbRv3756m2effTbnn39+Lrvssrz88sua0pdQz5498/zzz1e6jBZh4XHYeeedG9SQnnxyV9vDDjus+vlbb72Vgw8+OHfffXcOP/zwvPrqq1lppZUao9wW75hjjskxxxxT6TKqXXfddfnNb36TVq1a5dxzz82RRx6Zqqqq6vXPPvtsrr322nTq1KmCVbKs+uc//5ltt902c+bMyY477phf/OIX2XzzzVMqlaq3ef755/P73/8+F198cR5++GFN6UvgiCOOyBFHHLHI8qFDh2bcuHE54ogjMmrUqOYvrJmst9562WyzzfLoo4/mmmuuySmnnFLntu+9917++te/JkkOOeSQ5ioRAAAAACqqYf83BwAAAIDl1qmnnpoHHngg3bp1y4QJEzJ8+PAaDelJ0r9//1xyySW59957q+9UDfXxwQcfJMki59bS6NatW/7whz+kqqoq77zzTu66665G2zf1c+WVVyb5pMn1e9/7Xo2G9OSTDDnttNNy4oknVqI8lmFz587N1772tcyZMyfDhw/P7bffni222KJGQ3rySYPx2WefnaeffjpbbbVVhaplWbOwwfyaa65Z7HZ//OMfM3fu3Ky11lq13l0eAAAAAIpIUzoAAAAAS2zGjBk5//zzkyTnnXde+vbtu9jtt9lmm2y99daLLL/tttuyyy67ZJVVVklVVVXWXHPNHHXUUZkyZUqt++nbt29KpVJeffXVjBs3LjvssEO6dOmSrl27Zq+99spLL71Uve2tt96aL33pS+nUqVO+8IUv5MADD8wbb7yxyD7Hjh2bUqmUoUOH5qOPPsqpp56addZZJ+3atUvPnj1z9NFH59133621noceeignnnhiNttss6y22mqpqqpKr169csghh+SZZ56pdcyoUaNSKpUyatSovP322znmmGPSt2/frLDCCtV383711VdTKpVqPa6vvfZajjzyyKy11lqpqqpKx44ds9Zaa2WvvfbK9ddfX+trPvjgg9l7773TrVu3tG3bNmussUaGDx+e5557rtbthw4dmlKplLFjx+b555/P1772tayyyipp3759Nt100/zxj3+sddznKZfLufrqq7PddtulS5cuad++fdZbb72cdNJJixzjhcfpiiuuSJIcfvjhKZVK1T+rpdW9e/f069cvSWqcN0n9z8t33nknP/jBD7LeeuulXbt2WXHFFdO3b9/ssssuueiii6q3O+yww7Lmmmsm+eTnuPD9LHx81vPPP58jjjgiffv2TVVVVVZeeeXstttuueeee2qt49PXx7333pthw4ZllVVWqf5ZJjXPv9oszbkyceLE7LvvvunWrVtatWpV/bNbnFdeeSVJsvHGG3/utp9111135ZhjjsnAgQPTtWvXtGvXLmuvvXa++93vZvLkybWOOeyww6rPq9deey0HH3xwunXrlpVWWilbbbVVjS8oPPXUU9lnn32y2mqrpUOHDtl2223z0EMP1brfT/8Mr7322my++eZZaaWV0rVr13z1q1/N008/Xe/3lyT/+Mc/csABB6Rnz55p27ZtunXrlq997Wt5/PHHa92+IflQl6U9vm+88UaOOOKIrL766mnXrl022GCD/OY3v6nz9T7++OOcffbZ1ddRz549861vfStvvfVWvepe6A9/+EMmTZqUbt265aKLLvrcv7TQp0+fHHrooTWWfTaLf/vb32bw4MHp2LHjItdsfa+dT1+vtfn0tVXX8vrm8+zZszNy5MisueaaadeuXfr27Zvvf//7mTVr1mKPTWP49DVy4403Ztttt02XLl1qHIP6HpOFP5+FX2759DyxuJz729/+lm233TYdO3ZM586dM2zYsDqvqboceOCBadOmTV544YU8+uijdW539dVXJ0kOPvjgJMn8+fNzyy235IgjjsgGG2yQzp07p0OHDll//fVz4okn5v/+7//qVcenr7nafF7m13eeWdL5DgAAAIDlm6Z0AAAAAJbYbbfdlvfffz+rrrpq9t133wbtY+TIkfnKV76SO+64I+3bt8+AAQMybdq0XHzxxRk4cOBim7xuvvnmfPnLX85TTz2VtddeO/PmzcuYMWOy3XbbZerUqTn33HOz55575tVXX81aa62VDz74INdff33+67/+Kx9++GGt+yyXy9lrr72qG7fWX3/9TJs2LRdddFG22GKLTJs2bZExBx98cP77v/87r776arp165b1118/77//fq6++uoMHjx4kWbCT3v77bez2Wab5ZJLLknnzp3Tv3//tG7derHH7NVXX81mm22WSy+9NG+99VbWXXfdfPGLX8yMGTMyZsyYnHnmmYuMufjii7PNNtvk5ptvTpIMHDgws2fPzh/+8Idssskmue222+p8vX/+858ZPHhw7rjjjvTt2zcdO3bMY489lv3337+60W5JlcvlHHzwwTnkkEMyfvz4rLzyyunfv38mTZqUs88+O5tsskl1g3KS9O7dO0OGDKm+w36/fv0yZMiQDBkyJAMGDKjXay+ups+q73k5Y8aMbLHFFjnnnHMyadKkrL322llvvfXywQcf5M4778yPf/zj6m3XWWedbLbZZkmSqqqq6vez8PFpf/zjHzNw4MBcfvnleffdd9O/f/+0bds2f/3rX7PDDjvkggsuqPN9XXfdddlhhx3y8MMPZ6211soaa6yxRMdjac6V8ePHZ8stt8wdd9yRXr16VTfff55OnTol+aT5ur6GDRuWiy66KFOnTk2fPn3Sr1+/vPXWW7nkkkuyySab5Nlnn61z7KRJk7LZZptlzJgx6dWrV9q3b5+HHnoou+66a+65557cf//92WqrrXLPPfekd+/eadu2be677758+ctfrvMLJ0ly9tln56CDDsqUKVOy/vrr5+OPP84tt9ySzTffPPfff3+93t+5556bLbfcMjfccEM+/PDDbLjhhpk/f37+/Oc/Z4sttshNN91UY/uG5MPiLM3xfe2117LpppvmuuuuS48ePbLyyivn2WefzTHHHJPTTz99ke3nz5+fvffeOyeddFJeeOGF9OnTJ6uuumouv/zybL755nnnnXfqVXuS6ubs4cOHZ8UVV6z3+M/67ne/m29/+9t56623st5666VLly7V65bm2mmo+ubz7Nmz81//9V8588wz89prr6Vfv35ZccUVc+6552a77bbL3LlzG73G2px11lnZd9998+KLL2adddbJqquu2uB9tWvXrs55YsiQIendu/ciYy655JLstttu+de//pV11lkn8+fPz+23355tt902zz///BK/9qqrrpqdd945SeqcD1977bXq635hU/qbb76Zr371q7nyyiszffr0fPGLX0yfPn3y6quv5r//+78zePDgBn8Ro77qO8/UZ74DAAAAYDlXBgAAAIAldPTRR5eTlL/61a82aPxf/vKXcpJymzZtyldffXX18hkzZpT32muvcpJy3759y3PmzKkxrk+fPuUk5RVWWKF8zjnnlOfPn18ul8vl6dOnl7fccstykvJuu+1W7tChQ/maa66pHjd58uTyWmutVU5Svuiii2rs8957762upVOnTuV77rmnet1rr71WHjhwYDlJed99913kfVx55ZXll19+ucayjz76qPy73/2u3KZNm/Jaa61VXeNCP/vZz8pJyq1bty5vtdVW5SlTplSv++CDD8rlcrk8adKkcpJynz59aow95phjyknKhx56aPn999+vse65554r/8///E+NZY8//ni5TZs25STls88+u7qWDz/8sHzUUUeVk5Q7d+5cfuONN2qM22677aqP8zHHHFNd14IFC8onnXRSOUm5R48e5Y8//niRY1KXCy64oJyk3LFjx/Kdd95ZvfzNN98sDxkypJykvMUWWywy7tBDDy0nKV9++eVL/FoLLTxfahv75ptvlquqqspJyjfeeGO5XG7YefnLX/6ynKS80047ld95550ar/Haa6+Vzz333BrL6vrZftoTTzxRrqqqKrdr16586aWX1jiHbr311nKnTp3KrVu3Lk+cOLHW99u6devyqaeeWv7oo4/K5fInP7cPP/ywXC7/5/z72c9+VmPs0p4rrVu3Ln/7298uz549u3rdZ6/f2vzkJz8pJymXSqXyd77znfI//vGPJT6v/ud//qf8+uuv11g2Z86c8umnn15OUh46dOgiYxaeTyussEL5gAMOKM+cObNcLpfL8+fPr36fAwcOLPft27c8YsSI8ty5c6uPw+67715OUt5vv/0W2W+SWrNp9uzZ5YMOOqj6Z/7ZY7Lw+N177701lv/tb38rl0ql8iqrrFJ9fi60MF86duxY4+dR33z4PEt7fPfdd9/y9OnTq9dddNFF5STldu3a1VheLpfLv/71r8tJyl/4whfK9913X/XySZMmlTfccMPyCiusUOt5uzgrrbRSOUn5lltuWeIxn7Xwem3dunV5xRVXrLGvhT/Lhl47C6/XSZMm1fradZ0bDc3nE044ofo8fPrpp6uXT5w4sdyzZ8/qY9yQrP1sbbX9nBZeI23bti1feuml5QULFpTL5U/mzIVZ1dBjsiTzxMLX79ChQ43tZs6cWf7yl79cTlLef//96/N2yzfccEM5Sblbt2615tbCa2XLLbesXvbee++Vr7jiikXmi+nTp1dfw4cddtgi+6rrPX7ee68r8xsyz9R3vgMAAABg+aUpHQAAAIAl9tWvfrWcpHzCCSc0aPzCJuTjjjtukXWzZ88ur7LKKuUk5csuu6zGuoUNa3vuueci4+64447qprPa9nvJJZeUk5T32GOPGssXNqUnKf/qV79aZNwTTzxR3TT72Qb0xTn44IPLScoPPPBAjeULG8SqqqoWafhcqK7G5Z133rmcpPzEE08sUQ0Lm2FrO14LFiwob7DBBuUk5VNOOaXGuoWNfwMHDlykqX7evHnl7t27l5OUH3vssSWqY8GCBeVevXqVk9TatPbvf/+73LZt23KS8t///vca65qiKf2tt94q77DDDtVNsAsbkxtyXh555JH1anxdkqb0vffeu5yk/Otf/7rW9Qsb/I844oha3+/uu+9e577ralBsinNlScycObO86aabVl+DC5tGhwwZUj7ppJPKjzzySL33WS6Xy9tss005Sfnf//53jeULz6fVV1+9RgN9ufxJs2i7du3KScqDBg2qbppd6Pnnny8nKXfq1GmR11tY+2fzpVwul+fOnVt9zfz+97+vsa6uJttNNtlksefV97///XKS8s9//vPqZfXNh6Xxece3e/fu5VmzZi0ybuH7uummm6qXLViwoNy7d+9ykvJvfvObRcb885//rD6+S9qU/t5771WPWZrjsfB6TVI+55xzat2modfO0jal1yefZ86cWe7QoUM5Sfm2225b5LVuuumm6vfZ1E3pxx57bJ3jm6MpvbbXf/LJJ6u/PFAfH3zwQblz587lJOXbb799kfX9+/ev87yuS69evcodOnSobtRfqLGb0hsyz9R3vgMAAABg+dUqAAAAALCE3n///STJiiuuWO+xs2bNyoQJE5Ikxx577CLrO3TokG9961tJkjvvvLPWfXzjG99YZNnGG2+82PWDBg1Kkrzyyiu17rNt27b55je/ucjyjTbaKNtss03K5XKt9Tz//PP52c9+lr333jtDhw7NNttsk2222Sbjxo1LkjzxxBO1vt4OO+yQHj161LquLr169UqS/PnPf065XP7c7RfWW9txLpVK+d73vldju8864ogj0qpVzV8drrDCChk4cGCSuo/lZz333HOZMmVK2rVrV/2z/bSePXtmn332WWwtS+OMM86o/rlsuOGG6dWrV+6+++6ssMIK+e1vf5uOHTs2+Lxc+DO5+eab8/HHHy91rfPmzctf//rXtG7dOocddlit2+yxxx5JUn2Ofdbw4cPr/bpLe64cfPDBi5wrS6Jjx465//77c84552T99ddPksyZMycPPPBAzjrrrAwePDhf/epX895779U6/tFHH82PfvSj7LHHHtluu+2qf84vvvhikuTJJ5+sddyBBx6YDh061FjWuXPnrLnmmkmSww8/PKVSqcb6ddddN+3bt8/MmTPzzjvv1Lrfo48+epFln86WO+64o44j8R+vvfZaHnvssay22mrVP+vPqu0cqG8+LImlOb61zQ+DBw9OUjM7nnvuuUyePDnt2rWr9ZzfZJNNsuWWW9ar7oVzVFL3PLXLLrukVCot8qhLXdfV0l47DVWffL7vvvsyZ86c9OnTJ8OGDVtkX3vuuWd69uzZqPXVpSH51Jhqm+cHDBiQdu3aZcaMGXVe27Vp165d9t133yTJ1VdfXWPdY489lmeffTYrrLBC9t9//0XG3nPPPTnhhBOy2267Zdttt62+tmbMmJE5c+bkpZdequc7W3INnWcae74DAAAAoLjaVLoAAAAAAJYdHTt2TJLMnj273mP/9a9/ZcGCBamqqspaa61V6zYbbLBBklQ3Pn7W2muvvciyVVdddYnWz5o1q9Z9rrHGGtXv67PWX3/93H///YvUM3r06Jx88slZsGBBreOS5N13361zn/V19NFH58orr8wvfvGLXHXVVdlll13ypS99Kdtvv/0iDe7vvfde3n777SRJ//79a91fQ45zkqy22mpJ6j6Wn7Vw/717966zQfTzalkaL730UnWDX9u2bdO9e/dsu+22+f73v1/9ZYaGnpeHH354/vu//ztXXHFF/va3v9X4mdS1n8V58cUX8+GHH6Zt27bZdddda91mYcPx66+/Xuv6+p5bjXGuNOR8Xqhdu3YZMWJERowYkTfeeCMPP/xw7rvvvowZMyaTJk3KLbfckr333jv33HNP9ZhyuZxjjjkmF1100WL3Xdf1V9e5veqqq+a5555b7PrJkydn1qxZWXnllRdZX9dxWLh8Sc7vp556Kkny4YcfZptttql1mw8//DBJzXOgPvnweZrq+NaWHQuPSZ8+fRb5osBC66+/fh566KHPrXuhT2d5XfPUhhtuWF3HvHnz8sgjj9S5v1VWWSWrrLLKIssb49ppqIYc4/XWW6/WxvtWrVplnXXWqTNTGtPSZEVjWNy1PWXKlDqv7boccsghueyyyzJmzJjMmTOn+hxe2KS+66671tjfvHnzsv/++2fMmDGL3W9d11ZjaOg809jzHQAAAADF5U7pAAAAACyxhXdUnTRpUr3HLmyUW3XVVeu8K223bt2S1Lzb7afV1rj46X0tbn1ddxBe2Mi3pPWMHz8+P/7xj1MqlTJ69Og888wzmTVrVhYsWJByuZyf/OQnSZKPPvqo1n025C7zG2+8ccaPH5+ddtopr7/+ev7nf/4nBx98cNZYY43svPPOee6556q3/XRDYl3v7fOOc101Lrw775LejXlhLfU9xo3l8ssvT7lcTrlczty5c/Paa6/lD3/4Q4276zf0vOzRo0cmTJiQffbZJzNmzMiVV16Zb37zm1l77bWz1VZbVd99fUnNmDEjySeNiw888ECtjwcffDDJfxqTP6u+51ZTniv11aNHj+y111751a9+lRdffDHf//73kyT33ntvHnjggert/vCHP+Siiy7KiiuumIsuuigvvfRS5syZU/1zPuigg5LUff3V1fy88Gf/eevrmyP1Ob8XngMzZ86s8xz45z//mST54IMPqsfVJx8+z9Ie3/pkx6evvbosPH5LqnPnzllppZWSJK+++mqt2/zyl7/M/fffn/vvvz833XTTYvdX1/tpjGunoSp9jBuqsbKisV+/vvPaQttuu2369OmTWbNmVTeaz58/P9dff32ST5rWP+3MM8/MmDFj0r1791x11VV59dVX8+GHH1ZfW0OGDElS97XVGBo6zzT2fAcAAABAcWlKBwAAAGCJbb311kmSBx98MB9//HG9xi5sFHz77bfrbP566623kqTOO5c3hYV3u63NtGnTktSs55prrkmS/PCHP8yPfvSj9O/fPyuuuGJ10+qUKVOapM4tt9wyd9xxR6ZPn57bb789J510UtZYY43ceeed2XHHHfPee+8l+c9x/nT9n9Vcx3lhLXXV0Zy11GVpzsv1118/f/7zn/Pee+/l3nvvzahRo7LeeuvloYceyk477VRnU+zi6ujZs2d1k+LiHo2hJZ0rn9amTZucffbZ6d69e5LkH//4R/W6hdffOeeck+9+97v54he/mPbt21evb6rr7/PUlSO1ZUhdFv48hgwZ8rk//8+eW0uaD5+nOY/vp6+9uiwuO+qy5ZZbJknuu+++hhW2BJbm2vm8Lzg05C+R1KWpjnFja85j0lhKpVL1lzQW3h3973//e95888106dIlX/nKV2psv/DauuKKK3LIIYekT58+qaqqql5f32urIcdsaeaZxpzvAAAAACguTekAAAAALLFdd901K620UqZNm5Y///nP9Rr7xS9+Ma1atcrcuXPzyiuv1LrNM888kyRZZ511lrrWJTVlypQad739tIV3GP50PQsbrxY26H/WE0880bgFfsZKK62UnXfeOWeeeWaef/75rL322nn99dfzt7/9LUnSpUuX6rviPvvss7Xuo7mO88L9T548uc5jXImf+ac1xnlZVVWVoUOH5mc/+1mefvrpDBkyJLNmzcp1111XvU1dd2FfqF+/fllhhRXy5ptv5t13323gu6mflnSufFarVq3Sp0+fJJ/c1XehxV1/H330Ub3uCt6Y6nrd2jKkLv37968es2DBggbV8Xn58Hma8/h+Oh/mzJlT6zYNeb399tsvySd3fW+qZualuXYW3rG7rkbxl19+ubHKrH7tF154odbm5QULFuSFF15otNdrqIYek8/L1aY2fPjwJMldd92Vt99+u7o5fb/99qvRcJ4s/tp655138vrrr9frtT/vmP3rX/9aZFljzDNLMt8BAAAAsPzSlA4AAADAEuvSpUuOPfbYJMnxxx//uXfGfOCBB/Lggw8m+aRZcmEz1gUXXLDIth988EF+97vfJUl23nnnRqx68ebNm5fLLrtskeVPP/107rvvvpRKpey4447VyxfeNXjhXXA/7c4772zypvRP69ChQwYMGJAkeeONN6qXLzx+tR3ncrlcvbypj/P666+f3r1758MPP6z+2X7aG2+8kRtvvLFZaqlLY5+XrVu3zuDBg5PU/JksPG8++OCDWsd16NAhO++8cxYsWJDzzz+/fm9iKVTqXPm8uzO/99571c2+/fr1q16+uOvv8ssvX+wdoZvSRRddtMiyT2fLTjvt9Ln76NevXzbccMO8++67ueqqq5a6prryYXGa8/iut9566dWrVz744INa3+/EiRMzYcKEeu93+PDh6du3b956660cddRRDW7w/zwNvXbWWmutJMkjjzyyyLgbb7wx06dPb7Qat9lmm3To0CGvvvpq7rjjjkXW33rrrfVuhm4KDT0mn5erTW3dddfN4MGD8/HHH+f3v/99br755iTJIYccssi2i7u2zjnnnMyfP79er724Y/bvf/+71p93Y88zdc13AAAAACy/NKUDAAAAUC+jRo3KVlttlbfeeitbbbVV/vCHP+TDDz+ssc2LL76Yo48+OkOHDq3RfHrSSScl+aSB89prr61e/v7772f48OF5++2307dv3xxwwAHN82aStGnTJj/72c8ybty46mX//ve/q++Auvfee2fttdeuXrfNNtskSc4888xMmjSpevkjjzySI444Iu3atWv0Gr/73e/mhhtuWORuwuPHj8/f//73JMkmm2xSvfz73/9+2rRpk1tuuSXnnHNOdVPmvHnzctxxx+Xpp59O586d893vfrfRa/20UqmUH/7wh0mSn/3sZ9W1Jp805h1wwAGZN29ettxyy2y//fZNWsviNOS8/MlPfpLLLrss7733Xo19Pf300/njH/+YpObPZNVVV03Hjh0zbdq0Ou/+/Itf/CJVVVU57bTTcuaZZy7SaPnmm2/m17/+dS655JKler+fVqlzZdddd81BBx2Ue+65Jx999FGNdRMnTsyee+6Z999/P6uvvnqNpt6F19/JJ59co0H69ttvzw9/+MMmuf6WxG233ZZf//rX1Xej/uCDD/Ktb30rb7zxRnr16rXEmXbWWWelVCrl6KOPzu9+97t8/PHHNda/8sorOf3003PTTTdVL6tvPixOcx7fVq1aZcSIEUk+uZ4WfoEpSV577bUceuihWWGFFeq936qqqtxwww1p3759rrrqquy888556KGHFrlT+NSpU5fqWmrotTNs2LAkydlnn52XXnqpevkjjzyS733vew16z3Xp1KlTvvWtbyVJjjrqqBrZ8+STTzb66zVUQ4/Jwsbs8ePH13on+OawsAF91KhRmTVrVtZcc80MGTJkke0WXlvf//73q/9ySLlczlVXXZVf/vKX9b62Fh6zMWPG5K9//Wv18jfffDMHHXTQItmxUEPmmfrOdwAAAAAsvzSlAwAAAFAvbdu2zZ133pl99tknU6dOzfDhw9O1a9cMGDAgm2++edZYY42su+66ueiii9K9e/d88YtfrB77la98JT/60Y/y0Ucf5aCDDkrv3r0zePDgrL766vnzn/+cL3zhC/njH/9YfUfR5rD11ltnm222ydChQ7Puuutmk002yZprrpnHH388a621Vi688MIa23/729/OWmutlZdffjnrrbdeNtpoo6y33nrZfPPN07lz5xx11FGNXuOECRNywAEHpHPnzunfv3+22GKL9O3bN9ttt13ef//9HHzwwTWaujfeeOOcf/75KZVK+cEPfpAePXpk8803T7du3XLBBRekqqoq11xzTbp3797otX7W0Ucfna9//euZOXNmdthhh/Tr1y+bbrppevfunfvuuy+9e/fONddc0+R1LE5Dzstnnnkm3/zmN7PyyiunX79+2WKLLdKvX78MGDAgb7zxRrbffvsad8stlUr52te+luST5r3Bgwdn6NChGTp0aPU2G2+8ca677rpUVVVl5MiR6dq1awYNGpQtttgivXv3To8ePZboLxTUR6XOlfnz5+faa6/Nl7/85XTq1CkDBgzI4MGD07NnzwwaNCjjx49Ply5d8sc//jErrrhi9bgTTzwxXbt2zcMPP5w+ffpk0KBBWXPNNTNs2LBsuumm2WeffRq1ziV12mmn5fjjj68+ft27d89VV12Vdu3a5eqrr06HDh2WaD+77rprLrjggsydOzff+ta30rVr12y22WYZPHhwunfvnrXXXjsnn3xyjS/71DcfFqe5j++xxx6bXXfdNe+++26GDBmS9ddfP4MGDcoXv/jFvPfeeznyyCMbtN/NN98848aNS9++fXP33Xdnq622ysorr5xBgwZl0003Tc+ePdOzZ8+cfvrpad++fU455ZR6v0ZDr53DDz88G2ywQSZPnpz+/ftnwIABWXfddbP55ptn2223rf7LDY3ltNNOy6abbppJkyZlgw02yEYbbZQBAwZk4403zqqrrlqxa+bTGnpM9tprr7Rt2zbXX3991lxzzWy77bYZOnRorrjiimar/YADDkibNm2qv5x38MEHp1QqLbLdqaeemqqqqtx6663p2bNnNttss6yxxho59NBDc8ABB2SLLbao1+uuv/76+cY3vpGPP/44u+22W9Zaa60MGjQovXv3zrRp0+r8LNKQeaa+8x0AAAAAyy9N6QAAAADU20orrZQ///nPGT9+fL7xjW+kV69eefXVV/PEE0+kXC5nt912y2WXXZYXX3wxG264YY2xo0ePzl/+8pfsuOOOmTVrVp588smsssoq+c53vpMnnngigwcPbtb3UiqVcvPNN2fUqFFZsGBBnn322ay66qr57ne/m4cffniRhsJOnTrl/vvvz/Dhw9OpU6e88MILmTdvXkaMGJEJEyakY8eOjV7jueeem+OOOy4bbbRR/u///i8TJ05Mkuy888659dZbc9VVVy0y5rvf/W7uu+++fPWrX82CBQsyceLEdOjQIQcffHAee+yx7Lbbbo1eZ21KpVKuvvrqXHXVVfnSl76UadOm5ZlnnkmfPn3ywx/+MI899lj13W4rqb7n5cknn5wf/ehHGTx4cGbNmpWJEyfmgw8+yHbbbZerrroqd955Z9q0aVNjzK9//escd9xx6d69e5544omMGzeuxh36k0+aLJ999tkcd9xx6du3b1544YU8++yz6dChQ/baa69ceeWV+dGPftSo770S58qdd96ZP/zhD/n617+eddZZJ1OnTs3EiRMzZ86cbLHFFvnpT3+aF154ofruwgv17t07EyZMyN577522bdvm+eefT7t27XLqqafm9ttvX+SYN5cTTzwx11xzTXr16pVnnnkmpVIpe+yxRx5++OFsu+229drX0UcfnYkTJ+ab3/xmVl111TzzzDN56aWXssoqq+TAAw/Mn/70p+q/5JA0LB/q0tzHt3Xr1hkzZkxGjx6dddZZJ6+88kreeuutHHroofnHP/6RlVdeucH7Hjx4cF544YVcdtll2WOPPbLiiivm+eefz7PPPpvWrVtn1113za9//etMnjw5P//5zxv0Gg25dtq1a5d77rkn3/jGN9K1a9e89NJLadWqVX75y182yRd0VlpppYwdOzYnnXRSevfunRdeeCHvv/9+TjjhhIwbNy5VVVWN/pr11dBjsvbaa+cvf/lLtttuu0yfPj33339/xo0b16hf3Pk8q666anbZZZfq5wcffHCt22266aYZP358dtxxxyxYsCDPP/98VltttZx//vm58sorG/Tal1xySX7+859n7bXXzuuvv5633347Rx55ZCZMmJAuXbrUOa6+80xD5jsAAAAAlk+lcqX+piEAAAAAVNDYsWOz/fbbZ7vttsvYsWMrXQ6wDFp4R2S/ZgcAAAAAYHnnTukAAAAAAAAAAAAAADSYpnQAAAAAAAAAAAAAABpMUzoAAAAAAAAAAAAAAA2mKR0AAAAAAAAAAAAAgAYrlcvlcqWLAAAAAAAAAAAAAABg2eRO6QAAAAAAAAAAAAAANJimdAAAAAAAAAAAAAAAGkxTOgAAAAAAAAAAAAAADaYpHQAAAAAAAAAAAACABtOUDgAAAAAAAAAAAABAg2lKBwAAAAAAAAAAAACgwTSlAwAAAAAAAAAAAADQYJrSAQAAAAAAAAAAAABoME3pAAAAAAAAAAAAAAA0mKZ0AAAAAAAAAAAAAAAaTFM6AAAAAAAAAAAAAAANpikdAAAAAAAAAAAAAIAG05QOAAAAAAAAAAAAAECDaUoHAAAAAAAAAAAAAKDBNKUDAAAAAAAAAAAAANBgmtIBAAAAAAAAAAAAAGgwTekAAAAAAAAAAAAAADSYpnQAAAAAAAAAAAAAABqsTaULaGoLFizIG2+8kY4dO6ZUKlW6HKDgyuVy3n///fTo0SOtWvneT1OS70Bzku/NQ7YDzUm2Nx/5DjQn+d585DvQXGR785HtQHOS781HvgPNSb43D9kONKf6ZHvhm9LfeOON9OrVq9JlAMuZKVOmZI011qh0GYUm34FKkO9NS7YDlSDbm558BypBvjc9+Q40N9ne9GQ7UAnyvenJd6AS5HvTku1AJSxJthe+Kb1jx45JPjkYnTp1qnA1LBNmz0569Pjk32+8kay4YmXrYZkyc+bM9OrVqzp7aDryvQHkGzSYfG8esp0GMb/RQLK9+ch36k22sxTke/OR79SbfKeBZHvzke00iHyngeR785Hv1JtsZynI9+Yh22kQ+U4D1SfbC9+UvvDPU3Tq1EkAs2Rat/7Pvzt1Er40iD+N0/TkewPIN1hq8r1pyXYaxPzGUpLtTU++U2+ynUYg35uefKfe5DtLSbY3PdlOg8h3lpJ8b3rynXqT7TQC+d60ZDsNIt9ZSkuS7a2aoQ4AAAAAAAAAAAAAAApKUzoAAAAAAAAAAAAAAA2mKR0AAAAAAAAAAAAAgAbTlA4AAAAAAAAAAAAAQINpSgcAAAAAAAAAAAAAoME0pQMAAAAAAAAAAAAA0GCa0gEAAAAAAAAAAAAAaDBN6QAAAAAAAAAAAAAANJimdAAAAAAAAAAAAAAAGkxTOgAAAAAAAAAAAAAADaYpHQAAAAAAAAAAAACABtOUDgAAy5nRo0dn8ODB6dixY1ZbbbV89atfzQsvvFBjm3K5nFGjRqVHjx5p3759hg4dmmeeeeZz933jjTemf//+qaqqSv/+/XPzzTc31dsA4FNkO0AxyXeA4pHtAMUk3wGKSb4D1I+mdAAAWM6MGzcuRx99dB566KHcdddd+fjjj7PTTjtl9uzZ1ducffbZ+dWvfpULL7wwjzzySLp3754dd9wx77//fp37nTBhQvbff/8ccsgheeKJJ3LIIYdkv/32y8MPP9wcbwtguSbbAYpJvgMUj2wHKCb5DlBM8h2gfkrlcrlc6SKa0syZM9O5c+fMmDEjnTp1qnQ5LAtmz05WWumTf8+alay4YmXrYZkic5qPY90A8g0arOiZ8/bbb2e11VbLuHHjsu2226ZcLqdHjx45/vjjc9JJJyVJ5s6dm27duuWss87KkUceWet+9t9//8ycOTN/+9vfqpftsssu+cIXvpDrrrvuc+so+nGmiZjfaKCiZ05Lyfak+MeaJiDbWQpFzxz5zjJNvtNARc8b2c4yT77TQEXPHPnOMk22sxSKnjktJd+LfpxpIvKdBqpP5rhTOgAALOdmzJiRJOnatWuSZNKkSZk6dWp22mmn6m2qqqqy3Xbb5cEHH6xzPxMmTKgxJkl23nnnxY4BoGnIdoBiku8AxSPbAYpJvgMUk3wHWLw2lS4AWrJn3piRcoePl3j7L6zYNj27tG/CigAaR33yTbZBsZXL5YwYMSLbbLNNNtxwwyTJ1KlTkyTdunWrsW23bt3y2muv1bmvqVOn1jpm4f4+a+7cuZk7d27185kzZzboPcBCPr/DJyqZ7Yl8p3HJdvgP+U6RyHf4hGynaOQ7fEK+UySyHf7D/1elSOQ7TUVTOizGvhdPyAdt2y3x9u1XaJ27v7+dAAZapDfe+yA9/v+/65Nvsg2K7ZhjjsmTTz6Z+++/f5F1pVKpxvNyubzIsqUZM3r06Jx66qn1rBjq5vM7fKKS2Z7IdxqXbIf/kO8UiXyHT8h2ika+wyfkO0Ui2+E//H9VikS+01Q0pcNinL3vRlmzz2pLtO2/ps3K8TdMzPTZ84Qv0CJNnzOvuin9z9/dKuUOK37uGNkGxXbsscfm1ltvzfjx47PGGmtUL+/evXuST76hv/rqq1cvnzZt2iLf2P+07t27L/Lt/cWNGTlyZEaMGFH9fObMmenVq1eD3gskPr9DUvlsT+Q7jUu2wyfkO0Uj30G2U0zyHeQ7xSPb4ROVznfZTmOT7zQVTemwGGutumI26Nm50mUANLoNenROVvz8pnSgmMrlco499tjcfPPNGTt2bNZcc80a69dcc8107949d911VwYNGpQkmTdvXsaNG5ezzjqrzv1utdVWueuuu3LCCSdUL7vzzjuz9dZb17p9VVVVqqqqGuEdwSd8fmd51lKyPZHvNC7ZzvJOvlNU8p3lmWynyOQ7yzP5TlHJdpZ3LSXfZTuNTb7TVDSlAwDAcuboo4/Otddem1tuuSUdO3as/hZ+586d0759+5RKpRx//PE544wz0q9fv/Tr1y9nnHFGOnTokK9//evV+xk+fHh69uyZ0aNHJ0mOO+64bLvttjnrrLOy55575pZbbsndd99d65+wA6BxyXaAYpLvAMUj2wGKSb4DFJN8B6gfTekAALCcufjii5MkQ4cOrbH88ssvz2GHHZYkOfHEE/PBBx/kqKOOyvTp07PFFlvkzjvvTMeOHau3nzx5clq1alX9fOutt87111+fk08+OaecckrWXnvt3HDDDdliiy2a/D0BLO9kO0AxyXeA4pHtAMUk3wGKSb4D1I+mdAAAWM6Uy+XP3aZUKmXUqFEZNWpUnduMHTt2kWX77rtv9t1336WoDoCGkO0AxSTfAYpHtgMUk3wHKCb5DlA/rT5/EwAAAAAAAAAAAAAAqJ2mdAAAAAAAAAAAAAAAGkxTOgAAAAAAAAAAAAAADaYpHQAAAAAAAAAAAACABtOUDgAAAAAAAAAAAABAg2lKBwAAAAAAAAAAAACgwTSlAwAAAAAAAAAAAADQYJrSAQAAAAAAAAAAAABoME3pAAAAAAAAAAAAAAA0mKZ0AAAAAAAAAAAAAAAaTFM6AAAAAAAAAAAAAAANpikdAAAAAAAAAAAAAIAG05QOAAAAAAAAAAAAAECDVbQpffz48dl9993To0ePlEqljBkzpnrdRx99lJNOOikDBgzIiiuumB49emT48OF54403KlcwAAAAAAAAAAAAAAA1VLQpffbs2Rk4cGAuvPDCRdbNmTMnjz32WE455ZQ89thjuemmm/Liiy9mjz32qEClAAAAAAAAAAAAAADUpk0lX3zYsGEZNmxYres6d+6cu+66q8ayCy64IJtvvnkmT56c3r17N0eJAAAAAAAAAAAAAAAsRkWb0utrxowZKZVK6dKlS53bzJ07N3Pnzq1+PnPmzGaoDAAAAAAAAAAAAABg+dSq0gUsqQ8//DA/+tGP8vWvfz2dOnWqc7vRo0enc+fO1Y9evXo1Y5UAAAAAAAAAAAAAAMuXZaIp/aOPPsoBBxyQBQsW5KKLLlrstiNHjsyMGTOqH1OmTGmmKgEAAAAAAAAAAAAAlj9tKl3A5/noo4+y3377ZdKkSbnnnnsWe5f0JKmqqkpVVVUzVQcAAAAAAAAAAAAAsHxr0U3pCxvSX3rppdx7771ZeeWVK10SAAAAAAAAAAAAAACf0qqSLz5r1qxMnDgxEydOTJJMmjQpEydOzOTJk/Pxxx9n3333zaOPPpprrrkm8+fPz9SpUzN16tTMmzevkmUDsJRGjRqVUqlU49G9e/dKlwUAAAAAAAAAAAA0QEXvlP7oo49m++23r34+YsSIJMmhhx6aUaNG5dZbb02SbLzxxjXG3XvvvRk6dGhzlQlAE9hggw1y9913Vz9v3bp1BasBAAAAAAAAAAAAGqqiTelDhw5NuVyuc/3i1gGwbGvTpo27owMAAAAAAAAAAEABtKp0AQAsn1566aX06NEja665Zg444IC88sorlS4JAAAAAAAAAAAAaICK3ikdgOXTFltskauuuirrrLNO3nrrrZx22mnZeuut88wzz2TllVeudczcuXMzd+7c6uczZ85srnIBAAAAAAAAAACAxXCndACa3bBhw7LPPvtkwIAB2WGHHXLbbbclSa688so6x4wePTqdO3eufvTq1au5ygUAAAAAAAAAAAAWQ1M6ABW34oorZsCAAXnppZfq3GbkyJGZMWNG9WPKlCnNWCEAAAAAAAAAAABQlzaVLgAA5s6dm+eeey5f+tKX6tymqqoqVVVVzVgVAAAAAAAAAAAAsCTcKR2AZveDH/wg48aNy6RJk/Lwww9n3333zcyZM3PooYdWujQAAAAAAAAAAACgntwpHYBm9+9//zsHHnhg/u///i+rrrpqttxyyzz00EPp06dPpUsDAAAAAAAAAAAA6klTOgDN7vrrr690CQAAAAAAAAAAAEAjaVXpAgAAAAAAAAAAAAAAWHZpSgcAAAAAAAAAAAAAoME0pQMAAAAAAAAAAAAA0GCa0gEAYDkzfvz47L777unRo0dKpVLGjBlTY32pVKr18d///d917vOKK66odcyHH37YxO8GgIXkO0DxyHaAYpLvAMUk3wGKR7YD1I+mdAAAWM7Mnj07AwcOzIUXXljr+jfffLPG4/e//31KpVL22Wefxe63U6dOi4xt165dU7wFAGoh3wGKR7YDFJN8Bygm+Q5QPLIdoH7aVLoAAACgeQ0bNizDhg2rc3337t1rPL/llluy/fbbZ6211lrsfkul0iJjAWg+8h2geGQ7QDHJd4Biku8AxSPbAerHndIBAIA6vfXWW7ntttvyjW9843O3nTVrVvr06ZM11lgjX/nKV/L44483Q4UANIR8Byge2Q5QTPIdoJjkO0DxyHYATekAAMBiXHnllenYsWP23nvvxW633nrr5Yorrsitt96a6667Lu3atcuQIUPy0ksv1Tlm7ty5mTlzZo0HAM1DvgMUT1NmeyLfASrFZ3eAYpLvAMUj2wE0pQMAAIvx+9//PgcddFDatWu32O223HLLHHzwwRk4cGC+9KUv5Y9//GPWWWedXHDBBXWOGT16dDp37lz96NWrV2OXD0Ad5DtA8TRltifyHaBSfHYHKCb5DlA8sh1AUzoAAFCH++67Ly+88EK++c1v1ntsq1atMnjw4MV+o3/kyJGZMWNG9WPKlClLUy4AS0i+AxRPU2d7It8BKsFnd4Biku8AxSPbAT7RptIFAAAALdNll12WTTfdNAMHDqz32HK5nIkTJ2bAgAF1blNVVZWqqqqlKRGABpDvAMXT1NmeyHeASvDZHaCY5DtA8ch2gE9oSgcAgOXMrFmz8q9//av6+aRJkzJx4sR07do1vXv3TpLMnDkzf/rTn3LOOefUuo/hw4enZ8+eGT16dJLk1FNPzZZbbpl+/fpl5syZOf/88zNx4sT85je/afo3BEAS+Q5QRLIdoJjkO0AxyXeA4pHtAPWjKR0AAJYzjz76aLbffvvq5yNGjEiSHHroobniiiuSJNdff33K5XIOPPDAWvcxefLktGrVqvr5e++9l29/+9uZOnVqOnfunEGDBmX8+PHZfPPNm+6NAFCDfAcoHtkOUEzyHaCY5DtA8ch2gPrRlA4AAMuZoUOHplwuL3abb3/72/n2t79d5/qxY8fWeH7uuefm3HPPbYzyAGgg+Q5QPLIdoJjkO0AxyXeA4pHtAPXT6vM3AQAAAAAAAAAAAACA2mlKBwAAAAAAAAAAAACgwTSlAwAAAAAAAAAAAADQYJrSAQAAAAAAAAAAAABoME3pAAAAAAAAAAAAAAA0mKZ0AAAAAAAAAAAAAAAaTFM6AAAAAAAAAAAAAAANpikdAAAAAAAAAAAAAIAG05QOAAAAAAAAAAAAAECDaUoHAAAAAAAAAAAAAKDBNKUDAAAAAAAAAAAAANBgmtIBAAAAAAAAAAAAAGgwTekAAAAAAAAAAAAAADSYpnQAAAAAAAAAAAAAABpMUzoAAAAAAAAAAAAAAA2mKR0AAAAAAAAAAAAAgAbTlA4AAAAAAAAAAAAAQINpSgcAAAAAAAAAAAAAoME0pQMAAAAAAAAAAAAA0GCa0gEAAAAAAAAAAAAAaDBN6QAAAAAAAAAAAAAANJimdAAAAAAAAAAAAAAAGkxTOgAAAAAAAAAAAAAADaYpHQAAAAAAAAAAAACABtOUDgAAAAAAAAAAAABAg1W0KX38+PHZfffd06NHj5RKpYwZM6bG+nK5nFGjRqVHjx5p3759hg4dmmeeeaYyxQIAAAAAAAAAAAAAsIiKNqXPnj07AwcOzIUXXljr+rPPPju/+tWvcuGFF+aRRx5J9+7ds+OOO+b9999v5koBAAAAAAAAAAAAAKhNm0q++LBhwzJs2LBa15XL5Zx33nn5yU9+kr333jtJcuWVV6Zbt2659tprc+SRRzZnqQAAAAAAAAAAAAAA1KKid0pfnEmTJmXq1KnZaaedqpdVVVVlu+22y4MPPljnuLlz52bmzJk1HgAAAAAAAAAAAAAANI0W25Q+derUJEm3bt1qLO/WrVv1utqMHj06nTt3rn706tWrSesEAAAAAAAAAAAAAFietdim9IVKpVKN5+VyeZFlnzZy5MjMmDGj+jFlypSmLhEAAAAAAAAAAAAAYLnVptIF1KV79+5JPrlj+uqrr169fNq0aYvcPf3TqqqqUlVV1eT1AQAAAAAAAAAAAADQgu+Uvuaaa6Z79+656667qpfNmzcv48aNy9Zbb13BygAAAAAAAAAAAAAAWKiid0qfNWtW/vWvf1U/nzRpUiZOnJiuXbumd+/eOf7443PGGWekX79+6devX84444x06NAhX//61ytYNQAAAAAAAAAAAAAAC1W0Kf3RRx/N9ttvX/18xIgRSZJDDz00V1xxRU488cR88MEHOeqoozJ9+vRsscUWufPOO9OxY8dKlQwAAAAAAAAAAAAAwKdUtCl96NChKZfLda4vlUoZNWpURo0a1XxFAQAAAAAAAAAAAACwxFpVugAAAKB5jR8/Prvvvnt69OiRUqmUMWPG1Fh/2GGHpVQq1XhsueWWn7vfG2+8Mf37909VVVX69++fm2++uYneAQC1ke8AxSPbAYpJvgMUk3wHKB7ZDlA/mtIBAGA5M3v27AwcODAXXnhhndvssssuefPNN6sff/3rXxe7zwkTJmT//ffPIYcckieeeCKHHHJI9ttvvzz88MONXT4AdZDvAMUj2wGKSb4DFJN8Byge2Q5QP20qXQAAANC8hg0blmHDhi12m6qqqnTv3n2J93neeedlxx13zMiRI5MkI0eOzLhx43LeeefluuuuW6p6AVgy8h2geGQ7QDHJd4Biku8AxSPbAerHndIBAIBFjB07NquttlrWWWedfOtb38q0adMWu/2ECROy00471Vi2884758EHH6xzzNy5czNz5swaDwCalnwHKJ7myPZEvgM0N5/dAYpJvgMUj2wH+A9N6QAAQA3Dhg3LNddck3vuuSfnnHNOHnnkkfzXf/1X5s6dW+eYqVOnplu3bjWWdevWLVOnTq1zzOjRo9O5c+fqR69evRrtPQCwKPkOUDzNle2JfAdoTj67AxSTfAcoHtkOUFObShcAAAC0LPvvv3/1vzfccMNsttlm6dOnT2677bbsvffedY4rlUo1npfL5UWWfdrIkSMzYsSI6uczZ870CxSAJiTfAYqnubI9ke8Azclnd4Biku8AxSPbAWrSlA4AACzW6quvnj59+uSll16qc5vu3bsv8u39adOmLfIt/0+rqqpKVVVVo9UJQP3Id4DiaapsT+Q7QCX57A5QTPIdoHhkO7C8a1XpAgBg9OjRKZVKOf744ytdCgC1eOeddzJlypSsvvrqdW6z1VZb5a677qqx7M4778zWW2/d1OUB0EDyHaB4ZDtAMcl3gGKS7wDFI9uB5Z07pQNQUY888kguvfTSbLTRRpUuBWC5MWvWrPzrX/+qfj5p0qRMnDgxXbt2TdeuXTNq1Kjss88+WX311fPqq6/mxz/+cVZZZZXstdde1WOGDx+enj17ZvTo0UmS4447Lttuu23OOuus7Lnnnrnlllty99135/7772/29wewvJLvAMUj2wGKSb4DFJN8Byge2Q5QP+6UDkDFzJo1KwcddFB++9vf5gtf+EKlywFYbjz66KMZNGhQBg0alCQZMWJEBg0alJ/+9Kdp3bp1nnrqqey5555ZZ511cuihh2adddbJhAkT0rFjx+p9TJ48OW+++Wb186233jrXX399Lr/88my00Ua54oorcsMNN2SLLbZo9vcHsLyS7wDFI9sBikm+AxSTfAcoHtkOUD/ulA5AxRx99NHZbbfdssMOO+S0005b7LZz587N3Llzq5/PnDmzqcsDKKyhQ4emXC7Xuf6OO+743H2MHTt2kWX77rtv9t1336UpDYClIN8Bike2AxSTfAcoJvkOUDyyHaB+NKUDUBHXX399HnvssTzyyCNLtP3o0aNz6qmnNnFVAAAAAAAAAAAAQH21qnQBACx/pkyZkuOOOy5XX3112rVrt0RjRo4cmRkzZlQ/pkyZ0sRVAgAAAAAAAAAAAEvCndIBaHb//Oc/M23atGy66abVy+bPn5/x48fnwgsvzNy5c9O6desaY6qqqlJVVdXcpQIAAAAAAAAAAACfQ1M6AM3uy1/+cp566qkayw4//PCst956OemkkxZpSAcAAAAAAAAAAABaLk3pADS7jh07ZsMNN6yxbMUVV8zKK6+8yHIAAAAAAAAAAACgZWtV6QIAAAAAAAAAAAAAAFh2uVM6AC3C2LFjK10CAAAAAAAAAAAA0ADulA4AAAAAAAAAAAAAQINpSgcAAAAAAAAAAAAAoME0pQMAAAAAAAAAAAAA0GBL3ZQ+c+bMjBkzJs8991xj1APAMsQcAFA5MhigmOQ7QPHIdoBiku8AxSTfAYpHtgM0n3o3pe+333658MILkyQffPBBNttss+y3337ZaKONcuONNzZ6gQC0HOYAgMqRwQDFJN8Bike2AxSTfAcoJvkOUDyyHaBy6t2UPn78+HzpS19Kktx8880pl8t57733cv755+e0005r9AIBaDnMAQCVI4MBikm+AxSPbAcoJvkOUEzyHaB4ZDtA5dS7KX3GjBnp2rVrkuT222/PPvvskw4dOmS33XbLSy+91OgFAtBymAMAKkcGAxSTfAcoHtkOUEzyHaCY5DtA8ch2gMqpd1N6r169MmHChMyePTu33357dtpppyTJ9OnT065du0YvEICWwxwAUDkyGKCY5DtA8ch2gGKS7wDFJN8Bike2A1ROm/oOOP7443PQQQdlpZVWSu/evTN06NAkn/zZiwEDBjR2fQC0IOYAgMqRwQDFJN8Bike2AxSTfAcoJvkOUDyyHaBy6t2UftRRR2XzzTfPlClTsuOOO6ZVq09utr7WWmvltNNOa/QCAWg5zAEAlSODAYpJvgMUj2wHKCb5DlBM8h2geGQ7QOXUuyk9STbbbLNstNFGmTRpUtZee+20adMmu+22W2PXBkALZA4AqBwZDFBM8h2geGQ7QDHJd4Biku8AxSPbASqjVX0HzJkzJ9/4xjfSoUOHbLDBBpk8eXKS5Hvf+17OPPPMRi8QgJbDHABQOTIYoJjkO0DxyHaAYpLvAMUk3wGKR7YDVE69m9JHjhyZJ554ImPHjk27du2ql++www654YYbGrU4AFoWcwBA5chggGKS7wDFI9sBikm+AxSTfAcoHtkOUDlt6jtgzJgxueGGG7LlllumVCpVL+/fv39efvnlRi0OgJbFHABQOTIYoJjkO0DxyHaAYpLvAMUk3wGKR7YDVE6975T+9ttvZ7XVVltk+ezZs2uEOADFYw4AqBwZDFBM8h2geGQ7QDHJd4Biku8AxSPbASqn3k3pgwcPzm233Vb9fGFQ//a3v81WW23VeJUB0OKYAwAqRwYDFJN8Byge2Q5QTPIdoJjkO0DxyHaAymlT3wGjR4/OLrvskmeffTYff/xxfv3rX+eZZ57JhAkTMm7cuKaoEYAWwhwAUDkyGKCY5DtA8ch2gGKS7wDFJN8Bike2A1ROve+UvvXWW+eBBx7InDlzsvbaa+fOO+9Mt27dMmHChGy66aZNUSMALYQ5AKByZDBAMcl3gOKR7QDFJN8Bikm+AxSPbAeonHrfKT1JBgwYkCuvvLKxawFgGWAOAKgcGQxQTPIdoHhkO0AxyXeAYpLvAMUj2wEqY4ma0mfOnLnEO+zUqVODiwGg5TEHAFSODAYoJvkOUDyyHaCY5DtAMcl3gOKR7QAtwxI1pXfp0iWlUmmx25TL5ZRKpcyfP79RCgOgZTAHAFSODAYoJvkOUDyyHaCY5DtAMcl3gOKR7QAtwxI1pd97771NXQcALZQ5AKByZDBAMcl3gOKR7QDFJN8Bikm+AxSPbAdoGZaoKX277bZr6joAaKHMAQCVI4MBikm+AxSPbAcoJvkOUEzyHaB4ZDtAy7BETemfNX369Fx22WV57rnnUiqVsv766+fwww9P165dG7s+AFoYcwBA5chggGKS7wDFI9sBikm+AxSTfAcoHtkOUBmt6jtg3Lhx6du3b84///xMnz497777bs4///ysueaaGTduXFPUCEALYQ4AqBwZDFBM8h2geGQ7QDHJd4Biku8AxSPbASqn3ndKP/roo7P//vvn4osvTuvWrZMk8+fPz1FHHZWjjz46Tz/9dKMXCUDLYA4AqBwZDFBM8h2geGQ7QDHJd4Biku8AxSPbASqn3ndKf/nll/P973+/OrCTpHXr1hkxYkRefvnlRi0OgJbFHABQOTIYoJjkO0DxyHaAYpLvAMUk3wGKR7YDVE69m9I32WSTPPfcc4ssf+6557Lxxhs3Rk0AtFDmAIDKkcEAxSTfAYpHtgMUk3wHKCb5DlA8sh2gcurdlP69730vxx13XH75y1/m/vvvz/33359f/vKXOeGEE3L88cfnySefrH4srY8//jgnn3xy1lxzzbRv3z5rrbVWfv7zn2fBggVLvW8A6q855wAAamrMDB4/fnx233339OjRI6VSKWPGjKle99FHH+Wkk07KgAEDsuKKK6ZHjx4ZPnx43njjjcXu84orrkipVFrk8eGHHy7tWwcoNPkOUDyyHaCY5DtAMcl3gOKR7QCV06a+Aw488MAkyYknnljrulKplHK5nFKplPnz5y9VcWeddVYuueSSXHnlldlggw3y6KOP5vDDD0/nzp1z3HHHLdW+Aai/5pwDAKipMTN49uzZGThwYA4//PDss88+NdbNmTMnjz32WE455ZQMHDgw06dPz/HHH5899tgjjz766GL326lTp7zwwgs1lrVr125J3h7Acku+AxSPbAcoJvkOUEzyHaB4ZDtA5dS7KX3SpElNUUetJkyYkD333DO77bZbkqRv37657rrrPje0AWgazTkHAFBTY2bwsGHDMmzYsFrXde7cOXfddVeNZRdccEE233zzTJ48Ob17965zv6VSKd27d2+0OgGWB/IdoHhkO0AxyXeAYpLvAMUj2wEqp95N6X369GmKOmq1zTbb5JJLLsmLL76YddZZJ0888UTuv//+nHfeec1WAwD/0ZxzAAA1VTKDZ8yYkVKplC5duix2u1mzZqVPnz6ZP39+Nt544/ziF7/IoEGD6tx+7ty5mTt3bvXzmTNnNlbJAMsM+Q5QPEXM9kS+AxQx32U7gHwHKCLZDlA59W5KT5LXX389DzzwQKZNm5YFCxbUWPe9732vUQpLkpNOOikzZszIeuutl9atW2f+/Pk5/fTTq//ERm0EMEDTaq45AIBFVSKDP/zww/zoRz/K17/+9XTq1KnO7dZbb71cccUVGTBgQGbOnJlf//rXGTJkSJ544on069ev1jGjR4/Oqaee2iR1AyxL5DtA8RQt2xP5DpAUL99lO8An5DtA8ch2gMqod1P65Zdfnu985ztp27ZtVl555ZRKpep1pVKpUUP7hhtuyNVXX51rr702G2ywQSZOnJjjjz8+PXr0yKGHHlrrGAEM0HSacw4AoKZKZPBHH32UAw44IAsWLMhFF1202G233HLLbLnlltXPhwwZkk022SQXXHBBzj///FrHjBw5MiNGjKh+PnPmzPTq1atxigdYRsh3gOIpYrYn8h2giPku2wHkO0ARyXaAyql3U/pPf/rT/PSnP83IkSPTqlWrpqip2g9/+MP86Ec/ygEHHJAkGTBgQF577bWMHj26zqZ0AQzQdJpzDgCgpubO4I8++ij77bdfJk2alHvuuWex3+avTatWrTJ48OC89NJLdW5TVVWVqqqqpS0VYJkm3wGKp4jZnsh3gCLmu2wHkO8ARSTbASqn3qk7Z86cHHDAAc0S2HPmzFnkdVq3br3In9T4tKqqqnTq1KnGA4DG0ZxzAAA1NWcGL/zFyUsvvZS77747K6+8cr33US6XM3HixKy++upNUCFAcch3gOKR7QDFJN8Bikm+AxSPbAeonHon7ze+8Y386U9/aopaFrH77rvn9NNPz2233ZZXX301N998c371q19lr732apbXB6Cm5pwDAKipMTN41qxZmThxYiZOnJgkmTRpUiZOnJjJkyfn448/zr777ptHH30011xzTebPn5+pU6dm6tSpmTdvXvU+hg8fnpEjR1Y/P/XUU3PHHXfklVdeycSJE/ONb3wjEydOzHe+851GqRmgqOQ7QPHIdoBiku8AxSTfAYpHtgNUTpv6Dhg9enS+8pWv5Pbbb8+AAQOywgor1Fj/q1/9qtGKu+CCC3LKKafkqKOOyrRp09KjR48ceeSR+elPf9porwHAkmvOOQCAmhozgx999NFsv/321c9HjBiRJDn00EMzatSo3HrrrUmSjTfeuMa4e++9N0OHDk2STJ48ucbdBd577718+9vfztSpU9O5c+cMGjQo48ePz+abb16ftwmw3JHvAMUj2wGKSb4DFJN8Byge2Q5QOfVuSj/jjDNyxx13ZN11102SlEql6nWf/ndj6NixY84777ycd955jbpfABqmOecAAGpqzAweOnRoyuVynesXt26hsWPH1nh+7rnn5txzz61XHQDId4Aiku0AxSTfAYpJvgMUj2wHqJx6N6X/6le/yu9///scdthhTVAOAC2ZOQCgcmQwQDHJd4Dike0AxSTfAYpJvgMUj2wHqJxWn79JTVVVVRkyZEhT1AJAC2cOAKgcGQxQTPIdoHhkO0AxyXeAYpLvAMUj2wEqp95N6ccdd1wuuOCCpqgFgBbOHABQOTIYoJjkO0DxyHaAYpLvAMUk3wGKR7YDVE6b+g74xz/+kXvuuSf/+7//mw022CArrLBCjfU33XRToxUHQMtiDgCoHBkMUEzyHaB4ZDtAMcl3gGKS7wDFI9sBKqfeTeldunTJ3nvv3RS1ANDCmQMAKkcGAxSTfAcoHtkOUEzyHaCY5DtA8ch2gMqpd1P65Zdf3hR1ALAMMAcAVI4MBigm+Q5QPLIdoJjkO0AxyXeA4pHtAJXTqtIFAAAAAAAAAAAAAACw7Kr3ndKT5M9//nP++Mc/ZvLkyZk3b16NdY899lijFAZAy2QOAKgcGQxQTPIdoHhkO0AxyXeAYpLvAMUj2wEqo953Sj///PNz+OGHZ7XVVsvjjz+ezTffPCuvvHJeeeWVDBs2rClqBKCFMAcAVI4MBigm+Q5QPLIdoJjkO0AxyXeA4pHtAJVT76b0iy66KJdeemkuvPDCtG3bNieeeGLuuuuufO9738uMGTOaokYAWghzAEDlyGCAYpLvAMUj2wGKSb4DFJN8Byge2Q5QOfVuSp88eXK23nrrJEn79u3z/vvvJ0kOOeSQXHfddY1bHQAtijkAoHJkMEAxyXeA4pHtAMUk3wGKSb4DFI9sB6icejeld+/ePe+8806SpE+fPnnooYeSJJMmTUq5XG7c6gBoUcwBAJUjgwGKSb4DFI9sBygm+Q5QTPIdoHhkO0Dl1Lsp/b/+67/yl7/8JUnyjW98IyeccEJ23HHH7L///tlrr70avUAAWg5zAEDlyGCAYpLvAMUj2wGKSb4DFJN8Byge2Q5QOW3qO+DSSy/NggULkiTf+c530rVr19x///3Zfffd853vfKfRCwSg5TAHAFSODAYoJvkOUDyyHaCY5DtAMcl3gOKR7QCVU++m9FatWqVVq//cYH2//fbLfvvt16hFAdAymQMAKkcGAxSTfAcoHtkOUEzyHaCY5DtA8ch2gMpp9fmb1HT77bfn/vvvr37+m9/8JhtvvHG+/vWvZ/r06Y1aHAAtizkAoHJkMEAxyXeA4pHtAMUk3wGKSb4DFI9sB6icejel//CHP8zMmTOTJE899VRGjBiRXXfdNa+88kpGjBjR6AUC0HKYAwAqRwYDFJN8Byge2Q5QTPIdoJjkO0DxyHaAymlT3wGTJk1K//79kyQ33nhjdt9995xxxhl57LHHsuuuuzZ6gQC0HOYAgMqRwQDFJN8Bike2AxSTfAcoJvkOUDyyHaBy6n2n9LZt22bOnDlJkrvvvjs77bRTkqRr167V3zACoJgaaw64+OKLs9FGG6VTp07p1KlTttpqq/ztb39rkpoBisLncIBiku8AxSPbAYpJvgMUk3wHKB7ZDlA59b5T+jbbbJMRI0ZkyJAh+cc//pEbbrghSfLiiy9mjTXWaPQCAWg5GmsOWGONNXLmmWfmi1/8YpLkyiuvzJ577pnHH388G2ywQZPUDrCs8zkcoJjkO0DxyHaAYpLvAMUk3wGKR7YDVE6975R+4YUXpk2bNvnzn/+ciy++OD179kyS/O1vf8suu+zS6AUC0HI01hyw++67Z9ddd80666yTddZZJ6effnpWWmmlPPTQQ01VOsAyz+dwgGKS7wDFI9sBikm+AxSTfAcoHtkOUDn1vlN6796987//+7+LLD/33HMbpSAAWq6mmAPmz5+fP/3pT5k9e3a22mqrOrebO3du5s6dW/3cn1QCljc+hwMUk3wHKB7ZDlBM8h2gmOQ7QPHIdoDKqfed0gGgMTz11FNZaaWVUlVVle985zu5+eab079//zq3Hz16dDp37lz96NWrVzNWCwAAAAAAAAAAANRFUzoAFbHuuutm4sSJeeihh/Ld7343hx56aJ599tk6tx85cmRmzJhR/ZgyZUozVgsAAAAAAAAAAADUpU2lCwBg+dS2bdt88YtfTJJsttlmeeSRR/LrX/86//M//1Pr9lVVVamqqmrOEgEAAAAAAAAAAIAlsER3Sn/yySezYMGCpq4FgBaoueaAcrmcuXPnNvnrACxLfA4HKCb5DlA8sh2gmOQ7QDHJd4Dike0ALcMSNaUPGjQo//d//5ckWWuttfLOO+80aVEAtBxNMQf8+Mc/zn333ZdXX301Tz31VH7yk59k7NixOeigg5Z63wBF4nM4QDHJd4Dike0AxSTfAYpJvgMUj2wHaBmWqCm9S5cumTRpUpLk1Vdf9a0igOVIU8wBb731Vg455JCsu+66+fKXv5yHH344t99+e3bcccel3jdAkfgcDlBM8h2geGQ7QDHJd4Biku8AxSPbAVqGNkuy0T777JPtttsuq6++ekqlUjbbbLO0bt261m1feeWVRi0QgMpqijngsssua8wSAQrL53CAYpLvAMUj2wGKSb4DFJN8Byge2Q7QMixRU/qll16avffeO//617/yve99L9/61rfSsWPHpq4NgBbAHABQOTIYoJjkO0DxyHaAYpLvAMUk3wGKR7YDtAxL1JSeJLvsskuS5J///GeOO+44oQ2wHDEHAFSODAYoJvkOUDyyHaCY5DtAMcl3gOKR7QCVt8RN6Qtdfvnl1f/+97//nVKplJ49ezZqUQC0TOYAgMqRwQDFJN8Bike2AxSTfAcoJvkOUDyyHaByWtV3wIIFC/Lzn/88nTt3Tp8+fdK7d+906dIlv/jFL7JgwYKmqBGAFsIcAFA5MhigmOQ7QPHIdoBiku8AxSTfAYpHtgNUTr3vlP6Tn/wkl112Wc4888wMGTIk5XI5DzzwQEaNGpUPP/wwp59+elPUCUALYA4AqBwZDFBM8h2geGQ7QDHJd4Biku8AxSPbASqn3k3pV155ZX73u99ljz32qF42cODA9OzZM0cddZTQBigwcwBA5chggGKS7wDFI9sBikm+AxSTfAcoHtkOUDmt6jvg3XffzXrrrbfI8vXWWy/vvvtuoxQFQMtkDgConMbM4PHjx2f33XdPjx49UiqVMmbMmBrry+VyRo0alR49eqR9+/YZOnRonnnmmc/d74033pj+/funqqoq/fv3z80331yvugCWR/IdoHhkO0AxyXeAYpLvAMUj2wEqp95N6QMHDsyFF164yPILL7wwAwcObJSiAGiZzAEAldOYGTx79uw695ckZ599dn71q1/lwgsvzCOPPJLu3btnxx13zPvvv1/nPidMmJD9998/hxxySJ544okccsgh2W+//fLwww/XqzaA5Y18Byge2Q5QTPIdoJjkO0DxyHaAymlT3wFnn312dtttt9x9993ZaqutUiqV8uCDD2bKlCn561//2hQ1AtBCmAMAKqcxM3jYsGEZNmxYrevK5XLOO++8/OQnP8nee++d5JM/cdetW7dce+21OfLII2sdd95552XHHXfMyJEjkyQjR47MuHHjct555+W6666rV30AyxP5DlA8sh2gmOQ7QDHJd4Dike0AlVPvO6Vvt912efHFF7PXXnvlvffey7vvvpu99947L7zwQr70pS81RY0AtBDmAIDKaa4MnjRpUqZOnZqddtqpellVVVW22267PPjgg3WOmzBhQo0xSbLzzjsvdgwA8h2giGQ7QDHJd4Biku8AxSPbASqn3ndKT5IePXrk9NNPb+xaAFgGmAMAKqc5Mnjq1KlJkm7dutVY3q1bt7z22muLHVfbmIX7q83cuXMzd+7c6uczZ85sSMkAyzz5DlA8Rcv2RL4DJMXLd9kO8An5DlA8sh2gMup9p3QAAKD4SqVSjeflcnmRZUs7ZvTo0encuXP1o1evXg0vGIAlIt8Biqc5sj2R7wDNzWd3gGKS7wDFI9sB/kNTOgAAUK179+5Jssg38adNm7bIN/Y/O66+Y0aOHJkZM2ZUP6ZMmbIUlQOwOPIdoHiaM9sT+Q7QXHx2Bygm+Q5QPLIdYFGa0gEAgGprrrlmunfvnrvuuqt62bx58zJu3LhsvfXWdY7baqutaoxJkjvvvHOxY6qqqtKpU6caDwCahnwHKJ7mzPZEvgM0F5/dAYpJvgMUj2wHWFSb+mxcLpczefLkrLbaamnfvn1T1QRAC2QOAKicxs7gWbNm5V//+lf180mTJmXixInp2rVrevfuneOPPz5nnHFG+vXrl379+uWMM85Ihw4d8vWvf716zPDhw9OzZ8+MHj06SXLcccdl2223zVlnnZU999wzt9xyS+6+++7cf//9S10vQFHJd4Dike0AxSTfAYpJvgMUj2wHqKx6N6X369cvzzzzTPr169dUNQHQApkDACqnsTP40Ucfzfbbb1/9fMSIEUmSQw89NFdccUVOPPHEfPDBBznqqKMyffr0bLHFFrnzzjvTsWPH6jGTJ09Oq1b/+cNLW2+9da6//vqcfPLJOeWUU7L22mvnhhtuyBZbbLHU9QIUlXwHKB7ZDlBM8h2gmOQ7QPHIdoDKqldTeqtWrdKvX7+88847GhIBljPmAIDKaewMHjp0aMrlcp3rS6VSRo0alVGjRtW5zdixYxdZtu+++2bfffdd6voAlhfyHaB4ZDtAMcl3gGKS7wDFI9sBKqvV529S09lnn50f/vCHefrpp5uiHgBaMHMAQOXIYIBiku8AxSPbAYpJvgMUk3wHKB7ZDlA59bpTepIcfPDBmTNnTgYOHJi2bdumffv2Nda/++67jVZckrz++us56aST8re//S0ffPBB1llnnVx22WXZdNNNG/V1APh8zT0HAPAfMhigmOQ7QPHIdoBiku8AxSTfAYpHtgNUTr2b0s8777wmKKN206dPz5AhQ7L99tvnb3/7W1ZbbbW8/PLL6dKlS7PVAMB/NOccAEBNMhigmOQ7QPHIdoBiku8AxSTfAYpHtgNUTr2b0g899NCmqKNWZ511Vnr16pXLL7+8elnfvn2b7fUBqKk55wAAapLBAMUk3wGKR7YDFJN8Bygm+Q5QPLIdoHJaNWTQyy+/nJNPPjkHHnhgpk2bliS5/fbb88wzzzRqcbfeems222yzfO1rX8tqq62WQYMG5be//W2jvgYA9dNccwAAi5LBAMUk3wGKR7YDFJN8Bygm+Q5QPLIdoDLq3ZQ+bty4DBgwIA8//HBuuummzJo1K0ny5JNP5mc/+1mjFvfKK6/k4osvTr9+/XLHHXfkO9/5Tr73ve/lqquuqnPM3LlzM3PmzBoPABpHc84BANQkgwGKSb4DFI9sBygm+Q5QTPIdoHhkO0Dl1Lsp/Uc/+lFOO+203HXXXWnbtm318u233z4TJkxo1OIWLFiQTTbZJGeccUYGDRqUI488Mt/61rdy8cUX1zlm9OjR6dy5c/WjV69ejVoTwPKsOecAAGqSwQDFJN8Bike2AxSTfAcoJvkOUDyyHaBy6t2U/tRTT2WvvfZaZPmqq66ad955p1GKWmj11VdP//79ayxbf/31M3ny5DrHjBw5MjNmzKh+TJkypVFrAlieNeccAEBNMhigmOQ7QPHIdoBiku8AxSTfAYpHtgNUTr2b0rt06ZI333xzkeWPP/54evbs2ShFLTRkyJC88MILNZa9+OKL6dOnT51jqqqq0qlTpxoPABpHc84BANQkgwGKSb4DFI9sBygm+Q5QTPIdoHhkO0Dl1Lsp/etf/3pOOumkTJ06NaVSKQsWLMgDDzyQH/zgBxk+fHijFnfCCSfkoYceyhlnnJF//etfufbaa3PppZfm6KOPbtTXAWDJNOccAEBNMhigmOQ7QPHIdoBiku8AxSTfAYpHtgNUTr2b0k8//fT07t07PXv2zKxZs9K/f/9su+222XrrrXPyySc3anGDBw/OzTffnOuuuy4bbrhhfvGLX+S8887LQQcd1KivA8CSac45AICaZDBAMcl3gOKR7QDFJN8Bikm+AxSPbAeonDb1HbDCCivkmmuuyc9//vM8/vjjWbBgQQYNGpR+/fo1RX35yle+kq985StNsm8A6qe55wAA/kMGAxSTfAcoHtkOUEzyHaCY5DtA8ch2gMqpd1P6QmuvvXbWWmutJEmpVGq0ggBo+cwBAJUjgwGKSb4DFI9sBygm+Q5QTPIdoHhkO0Dza9WQQZdddlk23HDDtGvXLu3atcuGG26Y3/3ud41dGwAtkDkAoHJkMEAxyXeA4pHtAMUk3wGKSb4DFI9sB6iMet8p/ZRTTsm5556bY489NltttVWSZMKECTnhhBPy6quv5rTTTmv0IgFoGcwBAJUjgwGKSb4DFI9sBygm+Q5QTPIdoHhkO0Dl1Lsp/eKLL85vf/vbHHjggdXL9thjj2y00UY59thjhTZAgZkDACpHBgMUk3wHKB7ZDlBM8h2gmOQ7QPHIdoDKaVXfAfPnz89mm222yPJNN900H3/8caMUBUDLZA4AqBwZDFBM8h2geGQ7QDHJd4Biku8AxSPbASqn3k3pBx98cC6++OJFll966aU56KCDGqUoAFomcwBA5chggGKS7wDFI9sBikm+AxSTfAcoHtkOUDltlmSjESNGVP+7VCrld7/7Xe68885sueWWSZKHHnooU6ZMyfDhw5umSgAqxhwAUDkyGKCY5DtA8ch2gGKS7wDFJN8Bike2A7QMS9SU/vjjj9d4vummmyZJXn755STJqquumlVXXTXPPPNMI5cHQKWZAwAqRwYDFJN8Byge2Q5QTPIdoJjkO0DxyHaAlmGJmtLvvffepq4DgBbKHABQOTIYoJjkO0DxyHaAYpLvAMUk3wGKR7YDtAytKl0AAAAAAAAAAAAAAADLriW6U/qnffjhh7ngggty7733Ztq0aVmwYEGN9Y899lijFQdAy2IOAKgcGQxQTPIdoHhkO0AxyXeAYpLvAMUj2wEqp95N6UcccUTuuuuu7Lvvvtl8881TKpWaoi4AWiBzAEDlyGCAYpLvAMUj2wGKSb4DFJN8Byge2Q5QOfVuSr/tttvy17/+NUOGDGmKegBowcwBAJUjgwGKSb4DFI9sBygm+Q5QTPIdoHhkO0DltKrvgJ49e6Zjx45NUQsALZw5AKByZDBAMcl3gOKR7QDFJN8Bikm+AxSPbAeonHo3pZ9zzjk56aST8tprrzVFPQC0YOYAgMqRwQDFJN8Bike2AxSTfAcoJvkOUDyyHaBy2tR3wGabbZYPP/wwa621Vjp06JAVVlihxvp333230YoDoGUxBwBUjgwGKCb5DlA8sh2gmOQ7QDHJd4Dike0AlVPvpvQDDzwwr7/+es4444x069YtpVKpKeoCoAUyBwBUjgwGKCb5DlA8sh2gmOQ7QDHJd4Dike0AlVPvpvQHH3wwEyZMyMCBA5uiHgBaMHMAQOXIYIBiku8AxSPbAYpJvgMUk3wHKB7ZDlA5reo7YL311ssHH3zQFLUA0MKZAwAqRwYDFJN8Byge2Q5QTPIdoJjkO0DxyHaAyql3U/qZZ56Z73//+xk7dmzeeeedzJw5s8YDgOIyBwBUjgwGKCb5DlA8sh2gmOQ7QDHJd4Dike0AldOmvgN22WWXJMmXv/zlGsvL5XJKpVLmz5/fOJUB0OKYAwAqRwYDFJN8Byge2Q5QTPIdoJjkO0DxyHaAyql3U/q9997bFHUAsAwwBwBUjgwGKCb5DlA8sh2gmOQ7QDHJd4Dike0AlVPvpvTtttuuKeoAYBlgDgCoHBkMUEzyHaB4ZDtAMcl3gGKS7wDFI9sBKqfeTenjx49f7Pptt922wcVAU3r9vQ8yffa8z92uNGd2NmiGemBZZA4AqJzmzuC+ffvmtddeW2T5UUcdld/85jeLLB87dmy23377RZY/99xzWW+99Rq1NoAiac58l+0AzcNnd4Bi8tkdoJjkO0Dx+N0MQOXUuyl96NChiywrlUrV/54/f/5SFQRN4fX3PsgO54zLBx99/vnZft6Hee7///sLHdo2bWGwjDEHAFROc2fwI488UmOfTz/9dHbcccd87WtfW+y4F154IZ06dap+vuqqqzZqXQBF05z5LtsBmofP7gDF5LM7QDHJd4Di8bsZgMqpd1P69OnTazz/6KOP8vjjj+eUU07J6aef3miFQWOaPntePvhofs7bf+N8cbWVFrttac7s5NxP/t2jS/tmqA6WHeYAgMpp7gz+7C89zjzzzKy99tqf++fuVltttXTp0qXR6wEoqubMd9kO0Dx8dgcoJp/dAYpJvgMUj9/NAFROvZvSO3fuvMiyHXfcMVVVVTnhhBPyz3/+s1EKg6bwxdVWyoY9Fz2Ha5hd78sClhvmAIDKqWQGz5s3L1dffXVGjBhR4y4CtRk0aFA+/PDD9O/fPyeffHKtf3oOgP+oVL7LdoCm47M7QDH57A5QTPIdoHj8bgagchqt+3bVVVfNCy+80Fi7A2AZUt85YPTo0bnpppvy/PPPp3379tl6661z1llnZd11123CKgGKqTk+h48ZMybvvfdeDjvssDq3WX311XPppZdm0003zdy5c/OHP/whX/7ylzN27Nhsu+22tY6ZO3du5s6dW/185syZjV06wDKrqfO9qbI9ke8AdVmWP7sn8h2gLj67AxSTfAconmX5dzOyHVhW1Lsp/cknn6zxvFwu580338yZZ56ZgQMHNlphALQ8jTUHjBs3LkcffXQGDx6cjz/+OD/5yU+y00475dlnn82KK67Y2GUDFEIlP4dfdtllGTZsWHr06FHnNuuuu26NLxdttdVWmTJlSn75y1/W+cuT0aNH59RTT230egGWJZXK96bK9kS+AxTxs3si3wF8dgcoJvkOUDxF/N2MbAeWFfVuSt94441TKpVSLpdrLN9yyy3z+9//vtEKA6Dlaaw54Pbbb6/x/PLLL89qq62Wf/7zn4v9BQrA8qxSn8Nfe+213H333bnpppvqPXbLLbfM1VdfXef6kSNHZsSIEdXPZ86cmV69ejWoToBlVSXyvSmzPZHvAEX87J7IdwCf3QGKSb4DFE8Rfzcj24FlRb2b0idNmlTjeatWrbLqqqumXbt2jVYUAC1TU80BM2bMSJJ07dp1qfYDUGSV+hy+8ItDu+22W73HPv7441l99dXrXF9VVZWqqqqlKQ9gmVeJfG/KbE/kO0ARP7sn8h3AZ3eAYpLvAMVTxN/NyHZgWVHvpvQ+ffo0RR0ALAOaYg4ol8sZMWJEttlmm2y44YZ1bjd37tzMnTu3+vnMmTMbvRaAlqwSn8MXLFiQyy+/PIceemjatKn5nw4jR47M66+/nquuuipJct5556Vv377ZYIMNMm/evFx99dW58cYbc+ONNzZ73QDLkubOd9kO0PR8dgcoJp/dAYpJvgMUj9/NAFROvZvSk+Tvf/97/v73v2fatGlZsGBBjXVN+ScuAKi8xp4DjjnmmDz55JO5//77F7vd6NGjc+qpp9Z7/wBF0tyfw+++++5Mnjw5RxxxxCLr3nzzzUyePLn6+bx58/KDH/wgr7/+etq3b58NNtggt912W3bddddGrwugaJoz32U7QPPw2R2gmHx2Bygm+Q5QPH43A1AZ9W5KP/XUU/Pzn/88m222WVZfffWUSqWmqAuAFqix54Bjjz02t956a8aPH5811lhjsduOHDkyI0aMqH4+c+bM9OrVa6leH2BZUonP4TvttFPK5XKt66644ooaz0888cSceOKJTV4TQNE0d77LdoCm57M7QDH57A5QTPIdoHj8bgagcurdlH7JJZfkiiuuyCGHHNIU9QDQgjXWHFAul3Psscfm5ptvztixY7Pmmmt+7piqqqpUVVUt1esCLMt8DgcoJvkOUDyyHaCY5DtAMcl3gOKR7QCVU++m9Hnz5mXrrbduiloAaOEaaw44+uijc+211+aWW25Jx44dM3Xq1CRJ586d0759+6XeP0AR+RwOUEzyHaB4ZDtAMcl3gGKS7wDFI9sBKqdVfQd885vfzLXXXtsUtQDQwjXWHHDxxRdnxowZGTp0aFZfffXqxw033NAIVQIUk8/hAMUk3wGKR7YDFJN8Bygm+Q5QPLIdoHLqfaf0Dz/8MJdeemnuvvvubLTRRllhhRVqrP/Vr37VaMUB0LI01hxQLpebojyAQvM5HKCY5DtA8ch2gGKS7wDFJN8Bike2A1ROvZvSn3zyyWy88cZJkqeffrrGulKp1ChFAdAymQMAKkcGAxSTfAcoHtkOUEzyHaCY5DtA8ch2gMqpd1P6vffe2xR1ALAMMAcAVI4MBigm+Q5QPLIdoJjkO0AxyXeA4pHtAJXTqtIFAAAAAAAAAAAAAACw7NKUDgAAAAAAAAAAAABAg2lKBwAAAAAAAAAAAACgwTSlAwAAAAAAAAAAAADQYJrSAQAAAAAAAAAAAABoME3pAAAAAAAAAAAAAAA0mKZ0AAAAAAAAAAAAAAAabJlqSh89enRKpVKOP/74SpcCAAAAAAAAAAAAAECWoab0Rx55JJdeemk22mijSpcCAAAAAAAAAAAAAMD/t0w0pc+aNSsHHXRQfvvb3+YLX/hCpcsBAAAAAAAAAAAAAOD/Wyaa0o8++ujstttu2WGHHT5327lz52bmzJk1HgAAAAAAAAAAAAAANI02lS7g81x//fV57LHH8sgjjyzR9qNHj86pp57axFUBAAAAAAAAAAAAAJC08DulT5kyJccdd1yuvvrqtGvXbonGjBw5MjNmzKh+TJkypYmrBAAAAAAAAAAAAABYfrXoO6X/85//zLRp07LppptWL5s/f37Gjx+fCy+8MHPnzk3r1v+vvfsOk6ow9wf+Lm3pCEgVRBQVsCCRqAgIimLQeOVqomiCWGKCgPkJsWDXqxFujF6sGLFfaxQ1JPYG2FBBMAQIQUCxQOyUFann94eX1ZW2O7A7M2c/n+eZ52HOnDP7vrsz3zlzeOdM1RLbFBYWRmFhYUWXCgAAAAAAAAAAAABQKeX0UHrv3r1jxowZJZadcsop0b59+zjvvPM2GEgHAAAAAAAAAAAAAKBi5fRQer169WLPPfcssaxOnTrRuHHjDZYDAAAAAAAAAAAAAFDxqmS7AAAAAAAAAAAAAAAA8ldOnyl9YyZMmJDtEgAAAAAAAAAAAAAA+D/OlA4AAAAAAAAAAAAAQMYMpQMAAAAAAAAAAAAAkDFD6QAAAAAAAAAAAAAAZMxQOgAAAAAAAAAAAAAAGTOUDgAAAAAAAAAAAABAxgylAwAAAAAAAAAAAACQMUPpAAAAAAAAAAAAAABkzFA6AAAAAAAAAAAAAAAZM5QOAAAAAAAAAAAAAEDGDKUDAAAAAAAAAAAAAJAxQ+kAAAAAAAAAAAAAAGTMUDoAAAAAAAAAAAAAABkzlA4AAGzgsssui4KCghKX5s2bb3abiRMnxr777hs1a9aMnXfeOW655ZYKqhaA0pDtAOkk3wHSR7YDpJN8B0gn+Q7wnWrZLgAAAMhNe+yxRzz//PPF16tWrbrJdRcsWBBHHHFEnH766XHvvffGq6++GoMHD44mTZrEscceWxHlAlAKsh0gneQ7QPrIdoB0ku8A6STfAb5lKB0AANioatWqbfFT/OvdcsstseOOO8bo0aMjIqJDhw4xZcqU+OMf/+jgCUAOke0A6STfAdJHtgOkk3wHSCf5DvCtKtkuAAAAyE1z586Nli1bRtu2baN///4xf/78Ta77+uuvR58+fUosO/zww2PKlCmxevXqjW6zcuXKWLp0aYkLAOWrvLM9Qr4DZIN8B0gf2Q6QTvIdIJ38vyrAtwylAwAAG9h///3jnnvuiWeeeSbGjh0bixcvjgMPPDA+//zzja6/ePHiaNasWYllzZo1izVr1sRnn3220W1GjhwZDRo0KL60bt16m/cBwHcqItsj5DtARZPvAOkj2wHSSb4DpJP/VwX4jqF0AABgA3379o1jjz029tprrzj00EPjiSeeiIiIu+++e5PbFBQUlLieJMlGl693/vnnx5IlS4ovH3zwwTaqHoCNqYhsj5DvABVNvgOkj2wHSCf5DpBO/l8V4DvVsl0AAACQ++rUqRN77bVXzJ07d6O3N2/ePBYvXlxi2SeffBLVqlWLxo0bb3SbwsLCKCws3Oa1AlA65ZHtEfIdINvkO0D6yHaAdJLvAOnk/1WBysyZ0gEAgC1auXJlzJ49O1q0aLHR27t27RrPPfdciWXPPvtsdOnSJapXr14RJQJQRrIdIJ3kO0D6yHaAdJLvAOkk34HKzFA6AACwgbPPPjsmTpwYCxYsiDfeeCN+9rOfxdKlS2PgwIER8e1XxJ100knF6w8aNCjef//9GD58eMyePTvuuOOOuP322+Pss8/OVgsA/IBsB0gn+Q6QPrIdIJ3kO0A6yXeA71TLdgEAAEDu+fDDD+OEE06Izz77LJo0aRIHHHBATJ48Odq0aRMREYsWLYqFCxcWr9+2bdt48sknY9iwYXHTTTdFy5Yt4/rrr49jjz02Wy0A8AOyHSCd5DtA+sh2gHSS7wDpJN8BvmMoHQAA2MCDDz642dvvuuuuDZb17Nkz3n777XKqCICtJdsB0km+A6SPbAdIJ/kOkE7yHeA7VbJdAAAAAAAAAAAAAAAA+ctQOgAAAAAAAAAAAAAAGTOUDgAAAAAAAAAAAABAxgylAwAAAAAAAAAAAACQMUPpAAAAAAAAAAAAAABkzFA6AAAAAAAAAAAAAAAZM5QOAAAAAAAAAAAAAEDGDKUDAAAAAAAAAAAAAJAxQ+kAAAAAAAAAAAAAAGTMUDoAAAAAAAAAAAAAABkzlA4AAAAAAAAAAAAAQMYMpQMAAAAAAAAAAAAAkDFD6QAAAAAAAAAAAAAAZMxQOgAAAAAAAAAAAAAAGTOUDkBWTJo0KY466qho2bJlFBQUxOOPP57tkgAAAAAAAAAAAIAMGEoHICuKioqiU6dOceONN2a7FAAAAAAAAAAAAGArVMt2AQBUTn379o2+fftmuwwAAAAAAAAAAABgKxlKByAvrFy5MlauXFl8fenSpVmsBgAAAAAAAAAAAFivSrYLAIDSGDlyZDRo0KD40rp162yXBAAAAAAAAAAAAIShdADyxPnnnx9LliwpvnzwwQfZLgkAAAAAAAAAAACIiGrZLgAASqOwsDAKCwuzXQYAAAAAAAAAAADwA86UDgAAAAAAAAAAAABAxpwpHYCsWL58ebz77rvF1xcsWBDTp0+PRo0axY477pjFygAAAAAAAAAAAICyMJQOQFZMmTIlDj744OLrw4cPj4iIgQMHxl133ZWlqgAAAAAAAAAAAICyqpLtAjZn5MiR8eMf/zjq1asXTZs2jX79+sWcOXOyXRYA20CvXr0iSZINLgbSAQAAAAAAAAAAIL/k9FD6xIkTY8iQITF58uR47rnnYs2aNdGnT58oKirKdmkAAAAAAAAAAAAAAEREtWwXsDlPP/10iet33nlnNG3aNKZOnRoHHXRQlqoCAAAAAAAAAAAAAGC9nB5K/6ElS5ZERESjRo02uc7KlStj5cqVxdeXLl1a7nUBAAAAAAAAAAAAAFRWVbJdQGklSRLDhw+P7t27x5577rnJ9UaOHBkNGjQovrRu3boCqwQAAAAAAAAAAAAAqFzyZih96NCh8fe//z0eeOCBza53/vnnx5IlS4ovH3zwQQVVCAAAAAAAAAAAAABQ+VTLdgGlceaZZ8b48eNj0qRJ0apVq82uW1hYGIWFhRVUGQAAAAAAAAAAAABA5ZbTQ+lJksSZZ54Zjz32WEyYMCHatm2b7ZIAAAAAAAAAAAAAAPienB5KHzJkSNx///3xl7/8JerVqxeLFy+OiIgGDRpErVq1slwdAAAAAAAAAAAAAABVsl3A5owZMyaWLFkSvXr1ihYtWhRfHnrooWyXBgAAAAAAAAAAAABA5PiZ0pMkyXYJAAAAAAAAAAAAAABsRk6fKR0AAKh4I0eOjB//+MdRr169aNq0afTr1y/mzJmz2W0mTJgQBQUFG1z++c9/VlDVAGyJfAdIJ/kOkD6yHSCd5DtAOsl3gO8YSgcAAEqYOHFiDBkyJCZPnhzPPfdcrFmzJvr06RNFRUVb3HbOnDmxaNGi4suuu+5aARUDUBryHSCd5DtA+sh2gHSS7wDpJN8BvlMt2wUAAAC55emnny5x/c4774ymTZvG1KlT46CDDtrstk2bNo3tttuuHKsDIFPyHSCd5DtA+sh2gHSS7wDpJN8BvuNM6QAAwGYtWbIkIiIaNWq0xXU7d+4cLVq0iN69e8dLL7202XVXrlwZS5cuLXEBoOLId4B0ku8A6SPbAdJJvgOkU3nku2wH8oWhdAAAYJOSJInhw4dH9+7dY88999zkei1atIhbb701xo0bF48++mjsvvvu0bt375g0adImtxk5cmQ0aNCg+NK6devyaAGAjZDvAOkk3wHSR7YDpJN8B0in8sp32Q7ki2rZLgAAAMhdQ4cOjb///e/xyiuvbHa93XffPXbffffi6127do0PPvgg/vjHP27ya+nOP//8GD58ePH1pUuXOoACUEHkO0A6yXeA9JHtAOkk3wHSqbzyXbYD+cKZ0gEAgI0688wzY/z48fHSSy9Fq1atyrz9AQccEHPnzt3k7YWFhVG/fv0SFwDKn3wHSCf5DpA+sh0gneQ7QDqVZ77LdiBfOFM6AABQQpIkceaZZ8Zjjz0WEyZMiLZt22Z0P9OmTYsWLVps4+oAyJR8B0gn+Q6QPrIdIJ3kO0A6yXeA7xhKBwAAShgyZEjcf//98Ze//CXq1asXixcvjoiIBg0aRK1atSLi26+I++ijj+Kee+6JiIjRo0fHTjvtFHvssUesWrUq7r333hg3blyMGzcua30AUJJ8B0gn+Q6QPrIdIJ3kO0A6yXeA7xhKBwAAShgzZkxERPTq1avE8jvvvDNOPvnkiIhYtGhRLFy4sPi2VatWxdlnnx0fffRR1KpVK/bYY4944okn4ogjjqiosgHYAvkOkE7yHSB9ZDtAOsl3gHSS7wDfMZQOAACUkCTJFte56667Slw/99xz49xzzy2nigDYFuQ7QDrJd4D0ke0A6STfAdJJvgN8p0q2CwAAAAAAAAAAAAAAIH8ZSgcAAAAAAAAAAAAAIGOG0gEAAAAAAAAAAAAAyJihdAAAAAAAAAAAAAAAMmYoHQAAAAAAAAAAAACAjBlKBwAAAAAAAAAAAAAgY4bSAQAAAAAAAAAAAADImKF0AAAAAAAAAAAAAAAyZigdAAAAAAAAAAAAAICMVct2AQBAZj76akV8WbSq1Osv+LQo9ijHegAAAAAAAAAAAKicDKUDQB766KsVceg1E2PF6rWl3qbWqm/iqHKsCQAAAAAAAAAAgMrJUDoA5KEvi1bFitVrY/Tx+0S7pnVLtU3B10UR/1POhQEAAAAAAAAAAFDpGEoHgDzWrmnd2HOHBqVbucjLPgAAAAAAAAAAANtelWwXAAAAAAAAAAAAAABA/jKUDgAAAAAAAAAAAABAxgylAwAAAAAAAAAAAACQMUPpAAAAAAAAAAAAAABkzFA6AAAAAAAAAAAAAAAZM5QOAAAAAAAAAAAAAEDGDKUDAAAAAAAAAAAAAJAxQ+kAAAAAAAAAAAAAAGTMUDoAAAAAAAAAAAAAABkzlA4AAAAAAAAAAAAAQMYMpQMAAAAAAAAAAAAAkDFD6QAAAAAAAAAAAAAAZMxQOgAAAAAAAAAAAAAAGTOUDgAAAAAAAAAAAABAxqpluwBIm3c/WV6m9RvWqRE7bFernKoBAAAAAAAAAAAAgPJlKB22kYZ1akSt6lXjrIeml2m7WtWrxvO/62kwHQAAAAAAAAAAAIC8ZCidvPTRVyviy6JVpV6/rGcvz8QO29WK53/Xs8x1nfXQ9PiyaJWhdCCn+RYIAAAAAAAAAAAANsVQOnnno69WxKHXTIwVq9eWabta1atGwzo1yqmqb+2wXS1DmECq+BYIAAAAAAAAAAAAtsRQOnnny6JVsWL12hh9/D7RrmndUm/nrL0AZedbIAAAAAAAAAAAANgSQ+nkrXZN68aeOzTIdhkA28RHX60o8+B3RfEtEAAAAAAAAAAAAGyOoXSyLpcHMQEqwkdfrYhDr5kYK1avLdN2tapXjYZ1apRTVVuvrHntGy0AAAAAAAAAAADyk6F0siqtg5hlZXATKrcvi1bFitVrY/Tx+0S7pnVLvV2uZkHDOjWiVvWqcdZD08u0Xa3qVeP53/XMyZ4AAAAAAAAAAADYtLwYSr/55pvj6quvjkWLFsUee+wRo0ePjh49emS7rNQr6xnMM/HuJ8tTNYhZVlszuHnLgH2jcRkG89PyOyNd5HtJ7ZrWjT13aJDtMrbaDtvViud/17PM34Jx1kPT460FX8SXlfD1AHJVWXN64sSJMXz48Jg5c2a0bNkyzj333Bg0aFAFVgxAach3gHSS7wDpI9sB0km+A6SPbAf4Vs4PpT/00ENx1llnxc033xzdunWLP/3pT9G3b9+YNWtW7LjjjtkuL2+UdcD886JVMeh/p5b5DOaZqFW9avy4baNKOVSYyeDm+r/NwDveLNPPymSQvaIYKq2c0p7vZcndsn5bQj7YYbtaZXpeV+SHdDIhp6iMyprTCxYsiCOOOCJOP/30uPfee+PVV1+NwYMHR5MmTeLYY4/NQgcAbIx8B0gn+Q6QPrIdIJ3kO0D6yHaA7xQkSZJku4jN2X///eNHP/pRjBkzpnhZhw4dol+/fjFy5Mgtbr906dJo0KBBLFmyJOrXr1/qn1sRZwmvKJkOmFfaQb+iooi6/3eW3uXLI+rUyW49G5HLHzLIRJoea5lmTmWUrXzPREU852pVrxrP/65n+T5GczzfcjnbfLCHypjvZc3p8847L8aPHx+zZ88uXjZo0KB455134vXXXy/Vz6yMv2e2ge+9vs3810exx64tS7XZPz5aEj+94ZX425ndU/FNJZRdZc0c+U5ekO1shcqaOfKdvCDfyVBlzRvZTt6Q72SosmaOfCcvyHa2QmXMHNlO3pDvZKgsmZPTZ0pftWpVTJ06NUaMGFFieZ8+feK1114rt5/70Vcr4tBrJubsAG8malWvGneful+ZhukMueWusp6BOCLKfEb2ipLpmd8zkclQqedB+chWvkdU3OBzWXPXYy13s60icyoTuTwwnwnPhdyQSU6//vrr0adPnxLLDj/88Lj99ttj9erVUb169XKrF4DSke8A6STfAdJHtgOkk3wHSB/ZDlBSTg+lf/bZZ7F27dpo1qxZieXNmjWLxYsXb3SblStXxsqVK4uvL1myJCK+ndQvrQ8WL4mi5cti1DF7xc5NcussspnarnaNaLldYRm3Wh1Ll64ul3pyWlHRd/9eujRibTo+nFCvSkS9egXZLmMDO9YrjMdO7xxffV2+Q6VffL06znpwWgwYM6FM29WsXiXGD+0eLUs5ILk+a3L8SyiyLlv5/vFXK+I/bnwlvlm9rkz11qxeJW7u3zka1S79jn/Zc7cCMjeF+VYR2VZROZWJTLMtl9WsXiVGl/H5lokmdQujSf2apV6/suV7Jjm9ePHija6/Zs2a+Oyzz6JFixYbbLMtsj0i4tOl38Sny1dueUVSqWDF19Hh//69fPmyUj9+li9bGutWfh1/n78oli8r22OO3CTbt0y+ky9kO98n37dMvpMv5DvryfYtk+3kE/nOevJ9y+Q7+UK2833yffNkO/lEvrNeeWZ7Tg+lr1dQUHLYLEmSDZatN3LkyLj88ss3WN66desy/9xfjC7zJqRNy9J9RQXp1uHqsm+zbNmyaNDA15VsSbbyPRNHZvA4yGnyjRyWy8+3ypbvZcnpTa2/seXrZTvbSaEftS/zJt73UdmyPUK+k2dkOxmS7/KdHCffyYBsl+3kAflOBuS7fCfHyXYyVNnyXbaTd+Q7GShNtuf0UPr2228fVatW3eBTQ5988skGnxZa7/zzz4/hw4cXX1+3bl188cUX0bhx480G/fctXbo0WrduHR988EHUr18/8wZyhH5yV5p6idBPxLc7icuWLYuWBn43q6LzPW2PzVzgd7rt+Z1ue9vyd1rZ8j2TnG7evPlG169WrVo0btx4o9tsi3339dL2HNJP7kpTLxGVu5/Klu0R+Znv5a0yPwdyXZp6idBPRZLv36nM+f59ufx4zYR+clua+smlXmT7d2T7d3LpMbotpKmfNPUSoZ/yJN+/U5nzPZcek1srTb1E6CfX5XI/lS3fZfuGcvnxmQn95LY09ZPLvZQl23N6KL1GjRqx7777xnPPPRf/+Z//Wbz8ueeei6OPPnqj2xQWFkZhYWGJZdttt11GP79+/fo598fdGvrJXWnqJUI/lemTnpnKVr6n7bGZC/xOtz2/021vW/1OK1O+Z5LTXbt2jb/+9a8llj377LPRpUuXqF69+ka32Zb77uul7Tmkn9yVpl4iKm8/lSnbI/I738tbZX0O5IM09RKhn4oi378l30vK1cdrpvST29LUT670Itu/Jds3lCuP0W0lTf2kqZcI/ZQX+f4t+Z47j8ltIU29ROgn1+VqP5Up32X7puXq4zNT+sltaeonV3spbbZXKec6ttrw4cPjtttuizvuuCNmz54dw4YNi4ULF8agQYOyXRoAW0G+A+S2LeX0+eefHyeddFLx+oMGDYr3338/hg8fHrNnz4477rgjbr/99jj77LOz1QIAGyHfAdJJvgOkj2wHSCf5DpA+sh3gOzl9pvSIiOOPPz4+//zz+K//+q9YtGhR7LnnnvHkk09GmzZtsl0aAFtBvgPkti3l9KJFi2LhwoXF67dt2zaefPLJGDZsWNx0003RsmXLuP766+PYY4/NVgsAbIR8B0gn+Q6QPrIdIJ3kO0D6yHaA7+T8UHpExODBg2Pw4MEV9vMKCwvj0ksv3eArL/KVfnJXmnqJ0A9lV1H57m+57fmdbnt+p9ue3+nW21xO33XXXRss69mzZ7z99tvlXNXGpe3vrZ/claZeIvRTWeVTvpe3tD1m0tRPmnqJ0A8VQ75vXNoer/rJbWnqJ0295DPZvmlpe4ymqZ809RKhH8qHfP9Omh6TaeolQj+5Lm39pIFs/07aHp/6yW1p6ictvRQkSZJkuwgAAAAAAAAAAAAAAPJTlWwXAAAAAAAAAAAAAABA/jKUDgAAAAAAAAAAAABAxgylAwAAAAAAAAAAAACQsUo7lH7zzTdH27Zto2bNmrHvvvvGyy+/vNn1J06cGPvuu2/UrFkzdt5557jlllsqqNLSKUs/jz76aBx22GHRpEmTqF+/fnTt2jWeeeaZCqx288r6t1nv1VdfjWrVqsU+++xTvgWWUVn7WblyZVx44YXRpk2bKCwsjF122SXuuOOOCqp2y8raz3333RedOnWK2rVrR4sWLeKUU06Jzz//vIKq3bxJkybFUUcdFS1btoyCgoJ4/PHHt7hNrmcBG8o0U9i4TJ43bN7IkSPjxz/+cdSrVy+aNm0a/fr1izlz5mS7rLw2ZsyY2HvvvaN+/frF+zpPPfVUtstiK2SSPbm8D5Lp8z5X90My6SdX35NsbSbn2nuSTPvJ1fckmfaTy3lAxUjTcQbHGHL3+ZymYwxl7SVXX9epXNKU9RHyPlfzPk1ZHyHvyQ9pynfZnpvZHpGufJft5CJ5kZt5EZG+zNia/8fOtX2HTHrJ9X0H8tdll10WBQUFJS7Nmzcvvv3f//53nHzyydGyZcuoXbt2/OQnP4m5c+du8X6/+uqrGDJkSLRo0SJq1qwZHTp0iCeffLI8W4mI8utn9OjRsfvuu0etWrWidevWMWzYsPjmm2/Ks5WIiPjoo4/il7/8ZTRu3Dhq164d++yzT0ydOrX49iRJ4rLLLouWLVtGrVq1olevXjFz5swt3u+4ceOiY8eOUVhYGB07dozHHnusPNsoVh79jB07Nnr06BENGzaMhg0bxqGHHhpvvvlmebcSEeX391nvwQcfjIKCgujXr185VJ+5SjmU/tBDD8VZZ50VF154YUybNi169OgRffv2jYULF250/QULFsQRRxwRPXr0iGnTpsUFF1wQv/3tb2PcuHEVXPnGlbWfSZMmxWGHHRZPPvlkTJ06NQ4++OA46qijYtq0aRVc+YbK2st6S5YsiZNOOil69+5dQZWWTib9HHfccfHCCy/E7bffHnPmzIkHHngg2rdvX4FVb1pZ+3nllVfipJNOitNOOy1mzpwZDz/8cLz11lvxq1/9qoIr37iioqLo1KlT3HjjjaVaP9ezgA1lmilsWlmfN2zZxIkTY8iQITF58uR47rnnYs2aNdGnT58oKirKdml5q1WrVjFq1KiYMmVKTJkyJQ455JA4+uijy7TzTm4pa/bk+j5IJs/7XN4PyaSfXH1PsjWZnIvvSTLtJ1ffk2TST67nAeUvTccZHGPI7edzmo4xlLWXXH1dp/JIU9ZHyPtczvs0ZX2EvCf3pSnfZXvuZntEuvJdtpNr5EXu5kVE+jIj0//HzsV9h0x6yeV9B/LfHnvsEYsWLSq+zJgxIyK+Hajt169fzJ8/P/7yl7/EtGnTok2bNnHooYdu9v9uVq1aFYcddli899578cgjj8ScOXNi7NixscMOO+RlP/fdd1+MGDEiLr300pg9e3bcfvvt8dBDD8X5559frn18+eWX0a1bt6hevXo89dRTMWvWrLjmmmtiu+22K17nD3/4Q1x77bVx4403xltvvRXNmzePww47LJYtW7bJ+3399dfj+OOPjwEDBsQ777wTAwYMiOOOOy7eeOONvOxnwoQJccIJJ8RLL70Ur7/+euy4447Rp0+f+Oijj/Kyn/Xef//9OPvss6NHjx7l2EWGkkpov/32SwYNGlRiWfv27ZMRI0ZsdP1zzz03ad++fYllv/nNb5IDDjig3Gosi7L2szEdO3ZMLr/88m1dWpll2svxxx+fXHTRRcmll16adOrUqRwrLJuy9vPUU08lDRo0SD7//POKKK/MytrP1Vdfney8884lll1//fVJq1atyq3GTEVE8thjj212nVzPAja0LfKRTSvN84ay++STT5KISCZOnJjtUlKlYcOGyW233ZbtMtgGSpM9+bQPkiSle97n035IpjmWK+9Jvq8sveTqe5LvK00/uf6e5PtK00++5QHbXpqOMzjGkD/P5zQdY8j0fV8uvq6TXmnK+iSR9/mS92nK+iSR9+SmNOW7bM+PbE+SdOW7bCcXyIv8yIskSV9mlKWfXNx3+L7S9JLr+w7kt809N+bMmZNERPKPf/yjeNmaNWuSRo0aJWPHjt3kfY4ZMybZeeedk1WrVm3rcreoPPoZMmRIcsghh5RYNnz48KR79+7bpOZNOe+88zb7M9atW5c0b948GTVqVPGyb775JmnQoEFyyy23bHK74447LvnJT35SYtnhhx+e9O/ff+uL3ozy6ueH1qxZk9SrVy+5++67t6reLSnPftasWZN069Ytue2225KBAwcmRx999LYqe5uodGdKX7VqVUydOjX69OlTYnmfPn3itdde2+g2r7/++gbrH3744TFlypRYvXp1udVaGpn080Pr1q2LZcuWRaNGjcqjxFLLtJc777wz5s2bF5deeml5l1gmmfQzfvz46NKlS/zhD3+IHXbYIXbbbbc4++yzY8WKFRVR8mZl0s+BBx4YH374YTz55JORJEn8+9//jkceeSSOPPLIiih5m8vlLGBD2yIfIRuWLFkSEZH11+W0WLt2bTz44INRVFQUXbt2zXY5VJB82wcpzfM+n/ZDMsmxXHlP8kOl7SVX35P8UGn6yeX3JD9Umn7yLQ/YttJ0nMExhvQ9n/Pptb2scvV1nXRKU9ZHyPuIdOV9mrM+Qt5TvtKU77I9Xdkeke58l+1sS/Ii3XkRkY7MyMV9h0zk8r4D6TB37txo2bJltG3bNvr37x/z58+PiIiVK1dGRETNmjWL161atWrUqFEjXnnllU3e3/jx46Nr164xZMiQaNasWey5555x1VVXxdq1a8u3kf+zrfvp3r17TJ06Nd58882IiJg/f348+eST5f76tf65//Of/zyaNm0anTt3jrFjxxbfvmDBgli8eHGJ16LCwsLo2bPnZt9Xber1q7xnrcqrnx/6+uuvY/Xq1eX++lWe/fzXf/1XNGnSJE477bRyq39rVMt2ARXts88+i7Vr10azZs1KLG/WrFksXrx4o9ssXrx4o+uvWbMmPvvss2jRokW51bslmfTzQ9dcc00UFRXFcccdVx4lllomvcydOzdGjBgRL7/8clSrllsP50z6mT9/frzyyitRs2bNeOyxx+Kzzz6LwYMHxxdffBF33HFHRZS9SZn0c+CBB8Z9990Xxx9/fHzzzTexZs2a+I//+I+44YYbKqLkbS6Xs4ANbYt8hIqWJEkMHz48unfvHnvuuWe2y8lrM2bMiK5du8Y333wTdevWjcceeyw6duyY7bKoIPm0D1La532+7IdkmmO58p7k+0rbSy6/J/m+0vaTy+9Jvq+0/eRTHrDtpek4g2MM6Xs+58treyZy8XWd9EpT1kfI+4h05X2asz5C3lO+0pTvsj1d2R6R7nyX7WxL8iLdeRGR/5mRq/sOmcjlfQfy3/777x/33HNP7LbbbvHvf/87rrzyyjjwwANj5syZ0b59+2jTpk2cf/758ac//Snq1KkT1157bSxevDgWLVq0yfucP39+vPjii/GLX/winnzyyZg7d24MGTIk1qxZE5dcckne9dO/f//49NNPo3v37pEkSaxZsybOOOOMGDFiRLn2Mn/+/BgzZkwMHz48LrjggnjzzTfjt7/9bRQWFsZJJ51U/Hq7sdei999/f5P3u6nXr/KetSqvfn5oxIgRscMOO8Shhx66Tev/ofLq59VXX43bb789pk+fXp7lb5VKd6b09QoKCkpcT5Jkg2VbWn9jy7OlrP2s98ADD8Rll10WDz30UDRt2rS8yiuT0vaydu3aOPHEE+Pyyy+P3XbbraLKK7Oy/G3WrVsXBQUFcd9998V+++0XRxxxRFx77bVx11135cwnGMvSz6xZs+K3v/1tXHLJJTF16tR4+umnY8GCBTFo0KCKKLVc5HoWsKFM8xGyYejQofH3v/89HnjggWyXkvd23333mD59ekyePDnOOOOMGDhwYMyaNSvbZVFB8mkfpCzP+3zYD8kkx3LxPUlE6XrJl/ckEaX/2+TDe5KI0veTT3lA+UnTcQbHGNL1fM6H1/ayytXXddIvTVkfIe/TlPdpzPoIeU/FSVO+y/b0ZHtEOvNdtlNe5EX68iIi/zMjH/YdyiIf9h3IX3379o1jjz029tprrzj00EPjiSeeiIiIu+++O6pXrx7jxo2Lf/3rX9GoUaOoXbt2TJgwIfr27RtVq1bd5H2uW7cumjZtGrfeemvsu+++0b9//7jwwgtjzJgxednPhAkT4ve//33cfPPN8fbbb8ejjz4af/vb3+KKK64o117WrVsXP/rRj+Kqq66Kzp07x29+85s4/fTTN/g9ZjI3lY1Zq/LsZ70//OEP8cADD8Sjjz5a4oz45aE8+lm2bFn88pe/jLFjx8b2229fbrVvrfz+qFcGtt9++6hateoGn9z45JNPNvjUwXrNmzff6PrVqlWLxo0bl1utpZFJP+s99NBDcdppp8XDDz9c7p/8KI2y9rJs2bKYMmVKTJs2LYYOHRoR3z6ZkySJatWqxbPPPhuHHHJIhdS+MZn8bVq0aBE77LBDNGjQoHhZhw4dIkmS+PDDD2PXXXct15o3J5N+Ro4cGd26dYtzzjknIiL23nvvqFOnTvTo0SOuvPLKvPvkby5nARvamnyEbDjzzDNj/PjxMWnSpGjVqlW2y8l7NWrUiHbt2kVERJcuXeKtt96K6667Lv70pz9luTIqQr7sg5TleZ8P+yGZ5FiuvSdZr7S95Pp7kvXK8rfJ5fck65Wln3zJA8pHmo4zOMaQvudzPry2l1Wuvq6TbmnK+gh5H5GuvE9j1kfIeypGmvJdtqcr2yPSme+ynfIgL9KZFxHpyIxc3nfIRC7vO5A+derUib322ivmzp0bERH77rtvTJ8+PZYsWRKrVq2KJk2axP777x9dunTZ5H20aNEiqlevXmLQu0OHDrF48eJYtWpV1KhRo9z7WG9b9HPxxRfHgAED4le/+lVEROy1115RVFQUv/71r+PCCy+MKlXK57zRLVq02ODb4jt06BDjxo2LiG9fhyK+PfP5919DtzQ3tanXr/KetSqvftb74x//GFdddVU8//zzsffee2/DyjeuPPqZN29evPfee3HUUUcVL1u3bl1ERFSrVi3mzJkTu+yyyzbtIxOV7kzpNWrUiH333Teee+65Esufe+65OPDAAze6TdeuXTdY/9lnn40uXbpE9erVy63W0sikn4hvP7V48sknx/333x9HHnlkeZdZKmXtpX79+jFjxoyYPn168WXQoEHFZ0bdf//9K6r0jcrkb9OtW7f4+OOPY/ny5cXL/vWvf0WVKlWyPqCYST9ff/31Bi+s63co1n8COJ/kchawoUzzESpakiQxdOjQePTRR+PFF1+Mtm3bZrukVEqSJFauXJntMqggub4PksnzPpf3QzLNsVx8T1LWXnL9PUkmf5tcfk+SST+5ngeUrzQdZ3CMIX3P51x+bc9ELr6uUzmkKesj5H1EuvI+bVkfIe+pOGnKd9mermyPSF++y3bKi7xIX15EpCczcnnfIRO5vO9A+qxcuTJmz569wQeFGjRoEE2aNIm5c+fGlClT4uijj97kfXTr1i3efffd4mHaiG8fsy1atKjQgfSIbdPPpl6/kiQp19evbt26xZw5c0os+9e//hVt2rSJiIi2bdtG8+bNS7wWrVq1KiZOnLjZualNvX6V96xVefUTEXH11VfHFVdcEU8//fRmP2CwLZVHP+3bt9/g9es//uM/4uCDD47p06dH69aty6+hskgqoQcffDCpXr16cvvttyezZs1KzjrrrKROnTrJe++9lyRJkowYMSIZMGBA8frz589PateunQwbNiyZNWtWcvvttyfVq1dPHnnkkWy1UEJZ+7n//vuTatWqJTfddFOyaNGi4stXX32VrRaKlbWXH7r00kuTTp06VVC1W1bWfpYtW5a0atUq+dnPfpbMnDkzmThxYrLrrrsmv/rVr7LVQgll7efOO+9MqlWrltx8883JvHnzkldeeSXp0qVLst9++2WrhRKWLVuWTJs2LZk2bVoSEcm1116bTJs2LXn//feTJMm/LGBDW3rMUnZbet5QdmeccUbSoEGDZMKECSVel7/++utsl5a3zj///GTSpEnJggULkr///e/JBRdckFSpUiV59tlns10aGSrra3au74OU5nmfT/shmfSTq+9JMunlh3LpPUkm/eTye5JM+sn1PKD8pek4g2MMuf18TtMxhrL2kquv61Qeacr6JJH3uZz3acr6JJH35L405btsz91sT5J05btsJ9fIi9zNiyRJX2aUtZ8fyqV9h7L2kuv7DuS33/3ud8mECROS+fPnJ5MnT05++tOfJvXq1SvO8j//+c/JSy+9lMybNy95/PHHkzZt2iTHHHNMifsYMGBAMmLEiOLrCxcuTOrWrZsMHTo0mTNnTvK3v/0tadq0aXLllVfmZT+XXnppUq9eveSBBx5I5s+fnzz77LPJLrvskhx33HHl2subb76ZVKtWLfn973+fzJ07N7nvvvuS2rVrJ/fee2/xOqNGjUoaNGiQPProo8mMGTOSE044IWnRokWydOnSTfbz6quvJlWrVk1GjRqVzJ49Oxk1alRSrVq1ZPLkyXnZz3//938nNWrUSB555JESr1/Lli3Ly35+aODAgcnRRx9dnq2UWaUcSk+SJLnpppuSNm3aJDVq1Eh+9KMfJRMnTiy+beDAgUnPnj1LrD9hwoSkc+fOSY0aNZKddtopGTNmTAVXvHll6adnz55JRGxwGThwYMUXvhFl/dt8Xy7tJK5X1n5mz56dHHrooUmtWrWSVq1aJcOHD8+pwcSy9nP99dcnHTt2TGrVqpW0aNEi+cUvfpF8+OGHFVz1xr300kubfS7kYxawoc09Zim7LT1vKLuN/T4jIrnzzjuzXVreOvXUU4uf902aNEl69+5tID3PZfKancv7IKV53ufTfkgm/eTqe5JM/zbfl0vvSTLtJ1ffk2TaTy7nARUjTccZHGPI3edzmo4xlLWXXH1dp3JJU9YnibzP1bxPU9YnibwnP6Qp32V7bmZ7kqQr32U7uUhe5GZeJEn6MiOTv8/35dK+Qya95Pq+A/nr+OOPT1q0aJFUr149admyZXLMMcckM2fOLL79uuuuS1q1apVUr1492XHHHZOLLrooWblyZYn76Nmz5wZZ8dprryX7779/UlhYmOy8887J73//+2TNmjV52c/q1auTyy67LNlll12SmjVrJq1bt04GDx6cfPnll+Xez1//+tdkzz33TAoLC5P27dsnt956a4nb161bl1x66aVJ8+bNk8LCwuSggw5KZsyYsdl+kiRJHn744WT33XdPqlevnrRv3z4ZN25cebeSJEn59NOmTZuNZuqll16al/38UC4OpRckSR5+xw0AAAAAAAAAAAAAADmhSrYLAAAAAAAAAAAAAAAgfxlKBwAAAAAAAAAAAAAgY4bSAQAAAAAAAAAAAADImKF0AAAAAAAAAAAAAAAyZigdAAAAAAAAAAAAAICMGUoHAAAAAAAAAAAAACBjhtIBAAAAAAAAAAAAAMiYoXQAAAAAAAAAAAAAADJmKB0AyFu9evWKs846K9tlAAAAAAAAAABATttpp51i9OjR2S6DFDOUDlvh1ltvjV69ekX9+vWjoKAgvvrqq2yXBFCpPProo3HFFVdsk/tatGhRnHjiibH77rtHlSpVDLsDpJS8B6h45fVh0kcffTQOP/zw2H777aOgoCCmT5++zX8GANkxadKkOOqoo6Jly5ZRUFAQjz/+eLZLAmArOSYDkB8mTJiQ0fxLptsBAOliKB22wtdffx0/+clP4oILLsh2KQCVUqNGjaJevXrb5L5WrlwZTZo0iQsvvDA6deq0Te4TgNwj7wHSo6ioKLp16xajRo3KdikAbGNFRUXRqVOnuPHGG7NdCgDbiGMyAAAA6WcondRLkiT+8Ic/xM477xy1atWKTp06xSOPPFJ8+8yZM+PII4+M+vXrR7169aJHjx4xb9684tvvuOOO2GOPPaKwsDBatGgRQ4cOLb7trLPOihEjRsQBBxxQoT0B8K3vn3Fxp512iquuuipOPfXUqFevXuy4445x6623llj/ww8/jP79+0ejRo2iTp060aVLl3jjjTeKt7/uuuvipJNOigYNGlR0KwCV1iOPPBJ77bVX1KpVKxo3bhyHHnpoFBUVxcknnxz9+vWLq666Kpo1axbbbbddXH755bFmzZo455xzolGjRtGqVau44447StzfeeedF7vttlvUrl07dt5557j44otj9erVxbfLe4CKdfLJJ8fEiRPjuuuui4KCgigoKIh58+bFaaedFm3bto1atWrF7rvvHtddd12J7TZ2dvV+/frFySefXHx9wIABcckll8Shhx5aAZ0AsN7TTz8d3bt3j+222y4aN24cP/3pT0scU58xY0Yccsghxfv4v/71r2P58uXFt6/f17/88sujadOmUb9+/fjNb34Tq1atKl6nb9++ceWVV8YxxxxTob0BVFa9evWKoUOHxtChQ4vz/aKLLookSSIi4t57740uXbpEvXr1onnz5nHiiSfGJ598Urz9+rPjvvDCC9GlS5eoXbt2HHjggTFnzpzidRyTAciOXr16xZlnnhlnnXVWNGzYMJo1axa33nprFBUVxSmnnBL16tWLXXbZJZ566ql477334uCDD46IiIYNG0ZBQUHxsZiVK1fGb3/722jatGnUrFkzunfvHm+99VZExGa3A6BibWnf/vsmTJgQNWrUiJdffrl42TXXXBPbb799LFq0qCLLJkUMpZN6F110Udx5550xZsyYmDlzZgwbNix++ctfxsSJE+Ojjz6Kgw46KGrWrBkvvvhiTJ06NU499dRYs2ZNRESMGTMmhgwZEr/+9a9jxowZMX78+GjXrl2WOwJgU6655pro0qVLTJs2LQYPHhxnnHFG/POf/4yIiOXLl0fPnj3j448/jvHjx8c777wT5557bqxbty7LVQNUXosWLYoTTjghTj311Jg9e3ZMmDAhjjnmmOKDIi+++GJ8/PHHMWnSpLj22mvjsssui5/+9KfRsGHDeOONN2LQoEExaNCg+OCDD4rvs169enHXXXfFrFmz4rrrrouxY8fG//zP/2SrRYBK77rrrouuXbvG6aefHosWLYpFixZFq1atolWrVvHnP/85Zs2aFZdccklccMEF8ec//znb5QJQCkVFRTF8+PB466234oUXXogqVarEf/7nf8a6deuKv120YcOG8dZbb8XDDz8czz//fImTvUREvPDCCzF79ux46aWX4oEHHojHHnssLr/88ix1BEBExN133x3VqlWLN954I66//vr4n//5n7jtttsiImLVqlVxxRVXxDvvvBOPP/54LFiwYKPDhhdeeGFcc801MWXKlKhWrVqceuqpFdwFABtz9913x/bbbx9vvvlmnHnmmXHGGWfEz3/+8zjwwAPj7bffjsMPPzwGDBgQTZo0iXHjxkVExJw5c2LRokXFJxI499xzY9y4cXH33XfH22+/He3atYvDDz88vvjii2jduvUmtwOg4m1u3/771p8cZsCAAbFkyZJ455134sILL4yxY8dGixYtslA5aVCQbOwjEJASRUVFsf3228eLL74YXbt2LV7+q1/9Kr7++uvYaaed4sEHH4w5c+ZE9erVN9h+hx12iFNOOSWuvPLKzf6cCRMmxMEHHxxffvllbLfddtu6DQA2oVevXrHPPvvE6NGjY6eddooePXrE//7v/0bEt9+U0bx587j88stj0KBBceutt8bZZ58d7733XjRq1KjU9wtA+Xn77bdj3333jffeey/atGlT4raTTz45JkyYEPPnz48qVb79PHX79u2jadOmMWnSpIiIWLt2bTRo0CBuu+226N+//0Z/xtVXXx0PPfRQTJkyZYPb5D1AxShN3g4ZMiT+/e9/F3+73ca26devX2y33XZx1113ldj2vffei7Zt28a0adNin3322fYNALBZn376aTRt2jRmzJgRr7/+epx33nnxwQcfRJ06dSIi4sknn4yjjjoqPv7442jWrFmcfPLJ8de//jU++OCDqF27dkRE3HLLLXHOOefEkiVLivf/1ysoKIjHHnss+vXrV9GtAVQavXr1ik8++SRmzpwZBQUFERExYsSIGD9+fMyaNWuD9d96663Yb7/9YtmyZVG3bt3i/yt9/vnno3fv3hHxbf4feeSRsWLFiqhZs+YGP88xGYCK0atXr1i7dm3xWXDXH1c/5phj4p577omIiMWLF0eLFi3i9ddfj2+++WaD+ZeioqJo2LBh3HXXXXHiiSdGRMTq1atjp512irPOOivOOeccczMAOWJL+/brs3v9N5WuWrUqDjjggNh1111j5syZ0bVr1xg7dmwWOyDfOVM6qTZr1qz45ptv4rDDDou6desWX+65556YN29eTJ8+PXr06LHRgfRPPvkkPv744+IDJwDkvr333rv43wUFBdG8efPirxCdPn16dO7ceYsD6QBUnE6dOkXv3r1jr732ip///OcxduzY+PLLL4tv32OPPUoMpDRr1iz22muv4utVq1aNxo0bl/i66EceeSS6d+8ezZs3j7p168bFF18cCxcurJiGACi1W265Jbp06RJNmjSJunXrxtixY+U1QJ6YN29enHjiibHzzjtH/fr1o23bthERsXDhwpg9e3Z06tSpeCA9IqJbt26xbt26mDNnTvGyTp06FQ+kR0R07do1li9fXuJbkACoWAcccEDx0ErEt9k8d+7cWLt2bUybNi2OPvroaNOmTdSrVy969eoVEbHBPvz3j9GvP7Pi94/bAJAd38/n9cfVv3+svVmzZhGx6cyeN29erF69Orp161a8rHr16rHffvvF7Nmzy6lqADK1uX37H6pRo0bce++9MW7cuFixYoUPjrLVDKWTauvWrYuIiCeeeCKmT59efJk1a1Y88sgjUatWrU1uu7nbAMhNP/yQUUFBQfFrgVwHyD1Vq1aN5557Lp566qno2LFj3HDDDbH77rvHggULImLjub65rJ88eXL0798/+vbtG3/7299i2rRpceGFF8aqVasqpiEASuXPf/5zDBs2LE499dR49tlnY/r06XHKKaeUyOsqVarED7/gcfXq1RVdKgAbcdRRR8Xnn38eY8eOjTfeeCPeeOONiPj2zFpJkpT4T8/v29Tysq4DQMX65ptvok+fPlG3bt24995746233orHHnssImKDYy7fP26zPtPXH7cBIHu2dKx9S5m9/hjND/fXN7f/D0D+eO211yIi4osvvogvvvgiy9WQ7wylk2odO3aMwsLCWLhwYbRr167EpXXr1rH33nvHyy+/vNH/1KxXr17stNNO8cILL2ShcgC2tb333jumT59uBxogxxQUFES3bt3i8ssvj2nTpkWNGjWK/2OzrF599dVo06ZNXHjhhdGlS5fYdddd4/3339/GFQNQVjVq1ChxBpaXX345DjzwwBg8eHB07tw52rVrF/PmzSuxTZMmTWLRokXF19euXRv/+Mc/KqxmADbu888/j9mzZ8dFF10UvXv3jg4dOpT4tqOOHTvG9OnTo6ioqHjZq6++GlWqVInddtuteNk777wTK1asKL4+efLkqFu3brRq1apiGgFgA5MnT97g+q677hr//Oc/47PPPotRo0ZFjx49on379s5+DpBiNWrUiIgocSynXbt2UaNGjXjllVeKl61evTqmTJkSHTp02OR2AGTHpvbtq1atusG68+bNi2HDhsXYsWPjgAMOiJNOOskHS9kqhtJJtXr16sXZZ58dw4YNi7vvvjvmzZsX06ZNi5tuuinuvvvuGDp0aCxdujT69+8fU6ZMiblz58b//u//Fn+N6GWXXRbXXHNNXH/99TF37tx4++2344Ybbii+/8WLF8f06dPj3XffjYiIGTNmGHgEyFEnnHBCNG/ePPr16xevvvpqzJ8/P8aNGxevv/568Trrv1Fj+fLl8emnnxZ/uwYA5eONN96Iq666KqZMmRILFy6MRx99ND799NPig9hl1a5du1i4cGE8+OCDMW/evLj++us3OuAu7wEq1k477RRvvPFGvPfee/HZZ59Fu3btYsqUKfHMM8/Ev/71r7j44ovjrbfeKrHNIYccEk888UQ88cQT8c9//jMGDx4cX331VYl1vvjiixIZPmfOnJg+fXosXry4oloDqHQaNmwYjRs3jltvvTXefffdePHFF2P48OHFt//iF7+ImjVrxsCBA+Mf//hHvPTSS3HmmWfGgAEDolmzZsXrrVq1Kk477bSYNWtWPPXUU3HppZfG0KFDo0qVb//bavny5cX77RERCxYsiOnTp8fChQsrtF+AyuSDDz6I4cOHx5w5c+KBBx6IG264If7f//t/seOOO0aNGjXihhtuiPnz58f48ePjiiuuyOhnOCYDkPvatGkTBQUF8be//S0+/fTTWL58edSpUyfOOOOMOOecc+Lpp5+OWbNmxemnnx5ff/11nHbaaZvcDoDs2NS+/Q+tXbs2BgwYEH369IlTTjkl7rzzzvjHP/4R11xzTRaqJi0MpZN6V1xxRVxyySUxcuTI6NChQxx++OHx17/+Ndq2bRuNGzeOF198MZYvXx49e/aMfffdN8aOHVv8NUUDBw6M0aNHx8033xx77LFH/PSnP425c+cW3/ctt9wSnTt3jtNPPz0iIg466KDo3LlzjB8/Piu9ArBpNWrUiGeffTaaNm0aRxxxROy1114xatSoEp8E7dy5c3Tu3DmmTp0a999/f3Tu3DmOOOKILFYNkG7169ePSZMmxRFHHBG77bZbXHTRRXHNNddE3759M7q/o48+OoYNGxZDhw6NffbZJ1577bW4+OKLN1hP3gNUrLPPPjuqVq0aHTt2jCZNmsRPfvKTOOaYY+L444+P/fffPz7//PMYPHhwiW1OPfXUGDhwYJx00knRs2fPaNu2bRx88MEl1hk/fnx07tw5jjzyyIiI6N+/f3Tu3DluueWWCusNoLKpUqVKPPjggzF16tTYc889Y9iwYXH11VcX3167du145pln4osvvogf//jH8bOf/Sx69+4dN954Y4n76d27d+y6665x0EEHxXHHHRdHHXVUXHbZZcW3T5kypXi/PSJi+PDh0blz57jkkksqpE+Ayuikk06KFStWxH777RdDhgyJM888M379619HkyZN4q677oqHH344OnbsGKNGjYo//vGPGf0Mx2QAct8OO+wQl19+eYwYMSKaNWsWQ4cOjYiIUaNGxbHHHhsDBgyIH/3oR/Huu+/GM888Ew0bNtzsdgBUvE3t2//Q73//+3jvvffi1ltvjYiI5s2bx2233RYXXXRR8YkCoKwKkiRJsl0EAAAAAAAA6XfyySfHV199FY8//ni2SwHg//Tq1Sv22WefGD16dLZLAQAAtoJ9e7LNmdIBAAAAAAAAAAAAAMiYoXQAAAAAAAAAAAAAADJWkCRJku0iAAAAAAAAAAAAAADIT86UDgAAAAAAAAAAAABAxgylAwAAAAAAAAAAAACQMUPpAAAAAAAAAAAAAABkzFA6AAAAAAAAAAAAAAAZM5QOAAAAAAAAAAAAAEDGDKUDAAAAAAAAAAAAAJAxQ+kAAAAAAAAAAAAAAGTMUDoAAAAAAAAAAAAAABkzlA4AAAAAAAAAAAAAQMb+P4y/tzFQTMnpAAAAAElFTkSuQmCC", + "text/plain": [ + "
" + ] + }, + "metadata": {}, + "output_type": "display_data" + } + ], "source": [ - "# code" + "# plot posterior samples of each parameter with the groound truth value for comparison\n", + "fig, ax = plt.subplots(1,8, figsize=(30,5))\n", + "labels = ['ecc1', 'inc1', 'sma1', 'tau1', 'aop1', 'pan1', 'mtot', 'plx']\n", + "truth = [ecc, inc, sma, tau, aop, pan, mtot, plx]\n", + "for i in range(0,8):\n", + " ax[i].hist(nested_sampler.results.post[:, lab[labels[i]]], bins = 20, density = True, histtype='step')\n", + " ax[i].set_ylabel('number of samples')\n", + " ax[i].set_xlabel(labels[i])\n", + " ax[i].axvline(truth[i], color='red')\n", + "fig.suptitle('Comparison of Posterior Samples and Ground Truth Values', fontsize = 16)\n", + "plt.tight_layout()" ] }, { @@ -142,7 +175,9 @@ "metadata": {}, "outputs": [], "source": [ - "# code" + "# save results in .hdf5 format\n", + "filename = 'myposterior.hdf5'\n", + "nested_sampler.results.save_results(filename)" ] } ], @@ -165,7 +200,7 @@ "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", - "version": "3.9.12" + "version": "3.10.4" }, "orig_nbformat": 4 },