{
"cells": [
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Distances from parallaxes\n",
"\n",
"In principle, if we have a very good measurement of the parallax, we can simply inverse it to get the distance (in kpc, i.e in 1000 pc, where 1 pc=1 parsec = 3.26 light-years). This would work typically if the relative error on the parallax is below 20% at most. If the relative error on the parallax is larger than, then the distribution of distances is no more Gaussian and we cannot make such an approximation. You can see this below, by playing with the relative error on the parallax (relative_error_plx) and see how this changes the distribution of distance."
]
},
{
"cell_type": "code",
"execution_count": 7,
"metadata": {},
"outputs": [],
"source": [
"import pandas as pd\n",
"import numpy as np\n",
"from scipy import stats\n",
"import pylab as plt\n",
"%matplotlib inline\n"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"First, we define a function to compute and plot the distribution of the distance, assuming a normal (i.e. Gaussian) distribution of parallaxes, based on the relative error. The parallax is fixed to 1 mas (i.e. a nominal distance of 1,000 pc), as everything can be scaled to it (this is particularly true, because the error does not depend on the parallax, but on the brightness of the object, its colours, the number of observations and the position in the sky). \n",
"The mean distance should thus be this nominal distance, and the standard deviation should be given by the error_plx (as the parallax is 1). Note, however, how quickly the distribution of distances becoms non-Gaussian and the standard deviation should be asymmetric, and how the mode becomes different from the mean or median."
]
},
{
"cell_type": "code",
"execution_count": 8,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"Parallax= 1.0, with a relative error of 1.0%\n",
"Distances vary between 960.7 and 1045.9 pc\n",
"Formal distance: 1000.0 +/- 10.0 pc\n",
"Mean distance: 1000.0 +/- 10.0 pc\n",
"Median of distance distribution: 1000.0 pc\n",
"Mode of distance distribution: 999.0 pc\n"
]
},
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAYAAAAEICAYAAABWJCMKAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvnQurowAAIABJREFUeJzs3Xd4FFX3wPHv2d3QpFfpoYnSERBERGwIqKC8KioiCIi9vvrqT5SiYO9dQKQIWFARBQu9t4ABaVJCgBAgJJCQhIQku/f3x0xigCQkkGSym/N5njzZvTM7e2Z2d87MnTv3ijEGpZRSxY/L6QCUUko5QxOAUkoVU5oAlFKqmNIEoJRSxZQmAKWUKqY0ASilVDEVkAlARD4XkZfyaVn1RCRBRNz288UiMjQ/lm0v7zcRGZhfy8vD+44RkWgROZTL+Y2INLYf59v29Uf296FhDtPDReS6fHqvfP2++QsRCba/cx6nYwlkfpcA7B9XkojEi0isiKwUkQdFJGNdjDEPGmNeyeWycvyhGmP2GWPKGmO8+RD7KBH5+rTl9zTGTD7fZecxjrrAf4FmxpgL8/r6/Ny+/sj+PoQBiMgkERnjdEy5JSIlRGSm/dkYEel2lvkri8hPIpIoIntF5O5CCtVv5ec2FpHWIrLFPlh7KlN5kIissX/L58zvEoDtZmNMOaA+8DrwHPBlfr9JAB991AdijDFRTgeiHLEcuAfIzdnfJ0AKUAPoD3wmIs0LMLZ8l9XvOK+/7XPYF+TXNn4NeAZoDbwoIukHbE8DPxhj9ucxrlMZY/zqDwgHrjut7DLAB7Swn08CxtiPqwK/ArHAUWAZVuKbar8mCUgA/gcEAwYYAuwDlmYq89jLW2x/KGuBOOBnoLI9rRsQkVW8QA/7Q061329jpuUNtR+7gBeBvUAUMAWoYE9Lj2OgHVs0MDyH7VTBfv0Re3kv2su/zl5nnx3HpGxe/yxwEIgEBtvv3fh8tq89//dYP4o4e/s2z/Sek7B+DHOAeGAN0CjT9ObAPPt9DgMvZNpuzwO7gRjgu0yfSSnga7s8FlgH1Mhife8Dfsn0fBfwXabn+4E29mMDNAaG2Z9nir2Ov2T6zJ8BNtnr+S1QKpvtPAhYAXxkz7sduDbT9Mzfj8+AmZmmvQEsAOQcf0sRQLccpl9gr9tFmcqmAq9nM/9lwCp7Ox8EPgZKZJpugAeBncAx+7MWe5obeBvrex0GPEKm310W71UL+AHr+70HeDzTtFHATPtzPw4MzaasJPA+1nc80n5cMvNvGevg8hAw1YltDGzLFNNqexvXw9r/BJ3PvtQYExgJwC7fBzxkP57Evzuo14DPgSD778pMX7pTlsW/O9kp9gdTmqwTwAGghT3PD8DXmb802cVrfwm/Pm36Yv79gQ/G2vE0BMoCP6Z/8TLFMd6OqzVwErgkm+00BSs5lbNfuwMYkl2cp722B9YONn0dp5N9Asj19s20juX498cXmmnaJKyd+2WAB5gGfGNPK4e1U/kv1k69HNDRnvYk1o+jjr3cL4AZ9rQHgF+AMlg7mXZA+SzWuSHWjssF1MRKmgcyTTsGuOznWW6L0z7ztVg7qcpYP+IHs9nWg4A04Cl7+/XDSgTpCSzz96OM/TkOsrdzNFDHnlbPjj+7v7uzeO+z7ZzaAkmnlT1DpkR52rR2QCf7swu21/vJTNMN1sFCRTveI0APe9qDWMmvrr3NFpFNArA/o/XACKCE/fmEATdk+p2lArfY85bOpuxl+3tTHagGrAReyfQbScNKsiXt+Qt9G2MdMN2M9d0+BFQBZuW0zLz8+WsVUFYisb44p0vF+kHXN8akGmOWGXvL5mCUMSbRGJOUzfSpxpjNxphE4CXgjvSLxOepP/CuMSbMGJMA/B9w52mnn6ONMUnGmI3ARqxEcAo7ln7A/xlj4o0x4cA7wIBcxnEH8FWmdRyVw7x52r7GmIl2TCft5bYWkQqZZvnRGLPWGJOGlQDa2OU3AYeMMe8YY5LtZayxpz2AdTYUkWm5t9nbLRXrR9PYGOM1xqw3xhzPIq4wrLOONsBVwB/AARG52H6+zBjjy2E7nO5DY0ykMeYoVgJqk8O8UcD79vb7FvgHuDGLGE9gVSu8i3Uk+5gxJsKets8YUzGHv+l5iD1dWaxklFkcVvI9g71tVxtj0uzv3BdY2y6z140xscaYfVg7+fTtcoe9Dfbb2+y1HOLqAFQzxrxsjEmxP7vxwJ2Z5llljJlljPFl+h2fXtYfeNkYE2WMOQKM5tTfiA8YaYw5af/mnNjGzwAPAbOxDhKuwPqehonIzyKyRERuP4f3BaxMHShqYx09nu4trB3CnyICMM4Y8/pZlnW2erXM0/diHblVzV2YOaplLy/zsj1YdYPpMtcpnsD6Ap2uKtaR0enLqp2HONaf9trs5Hr72olpLHA71hFX+g61Kv/+CLJbv7pYVTxZqQ/8JCKZd9BerO021X7tNyJSEWvHOdwYk5rFcpZgHfk1th/HYu3ALref58Xp61Erh3kPnJY092Y3vzFmrYiEYR21fpfHmPIqASh/Wll5rB3QGUTkIqzk1B7rbMXDqd8jyP7zrcWZv6vs1AdqiUhspjI3VvVjuqx+w6eXZfV7y7zdjxhjknOIIz/kuI2NMXuBXgAiUgbrLOUGrCrDb7GqSzeLyAI7ceZJQJwBiEgHrJ3b8tOn2UeK/zXGNMQ6lXpaRK5Nn5zNIs92hpD5yns9rKPMaCAR64ufHpcba0eX2+VGYn25My87Das6Ji+i7ZhOX9aBXL7+IGeuY5byuH3vBvpgXYeogFVNACC5iGk/0CiHaT1POxorZYw5YB9VjzbGNAM6Y51J3JvNctITwJX24yVYCeAqsk8A+dGdbm2xs6etHtZ34Qwi8ghWlUQk1nWr9PL05srZ/fU/h7h2AB4RaZKprDWwJZv5P8OqxmlijCkPvEDuPlvIw3cO6/Pec9rnXc4Y0yvTPFl9LqeXZfV7y7zdT5m/CGzjEcAEY8xhoCUQYoyJw6pmanwO7+3fCUBEyovITcA3WHXrf2cxz00i0tj+gR3HOjJMb9J5GKv+MK/uEZFmdkZ+GevCnBfrwywlIjeKSBDWhdeSmV53GAjO3GT1NDOAp0SkgYiUBV4FvrWrQ3LNjuU7YKyIlBOR+litBr7O+ZUZvgMGZVrHkdnNmMftWw7rukUMVqJ8NQ+r9StwoYg8KSIl7fXqaE/7HGtd69sxVRORPvbjq0WkpZ2Mj2Mlxuya9C4BrgZK21Ury7Cuh1QB/srmNef6HcqsOvC43bTvduASYO7pM9lH2GOwqoEGAP8TkTZwSnPl7P6mZVpOSREpZT8tISKlTktA2MtMxLoO9bKIXCAiV2Al8KnZrEc5rG2cYFedPZSHbfCdvQ3qiEglrIv62VkLHBeR50SktIi4RaSFfSCYFzOwWtZUE5GqWDvYbH8jTm5jEWmGdXDymV20B7hGRGoATbCugeaZvyaAX0QkHutIYDjWaed92czbBJiPdaq1CvjUGLPYnvYa1hcgVkSeycP7T8W6+HcI64Lk4wB2Nn4YmIB1tJ2IlZ3TfW//jxGRDVksd6K97KVYH3Ay8Fge4srsMfv9w7DOjKbbyz8rY8xvWBdoF2JdlF6Yw+x52b5TsC+uAluxLsDlijEmHrge6yzjEFZLkqvtyR9g1ZH+aX8vVgPpyeFCrNYfx7EuSi4hmx+5MWaHvR7L7OfHsbbfCpP9fSBfAs3sdZyV2/U5zRqs7RiNVUV2mzEmJvMM9vWMr4E3jDEbjTE7sY6wp4pIydMXeBb/YLXOqo11rSMJ+0hYRF4Qkd8yzfsw1gXQKKwd5kPGmOzOAJ7BOsuLx6qT/zYPMY23Y9kIbMDaKWbJ/ixuxrp+sAdru03AOqvMizFACFZrrb/t982vezryext/AjyR6Xv4f1j7nS3Aq8aYXN3Qebr01hpKKQeIyCCsVj5dnI5FFT/+egaglFLqPGkCUEqpYkqrgJRSqpjSMwCllCqmivSNYFWrVjXBwcFOh6GUUn5l/fr10caYamebr0gngODgYEJCQpwOQyml/IqI5HQndQatAlJKqWJKE4BSShVTmgCUUqqYKtLXAJRS+SM1NZWIiAiSkwu6c0tVmEqVKkWdOnUICgo6p9drAlCqGIiIiKBcuXIEBweTRZ9kyg8ZY4iJiSEiIoIGDRqc0zK0CkipYiA5OZkqVarozj+AiAhVqlQ5r7M6TQBKFRO68w885/uZagJQSqliSq8BKGU7mpjCqt0xrAs/Ssjeo+yNPpExLcjjomXtCnQIrkT74Mq0r18Jj1uPn/LC7XbTsmXLjOezZs2isO70nzRpEiEhIXz88cdnTPv9998ZMWIEx48fp1SpUjRt2pS33nqLevVyGpTs/I0YMYKuXbty3XXXFej75OSsCUBE6mIN5HEh1hiu44wxH4jIKOB+4Ig96wvGmLn2a/4PGII18tLjxpg/7PIeWIN3uLGGNjvb2LxKFbi9MYmMWxrG9+sjSEnzUdqTTFvfLv4jEbgirEGtEm9OI3R/LG//uQOAupVLc/+VDbm9XV1Kl3A7Gb7fKF26NKGhoXl+XVpaGh5PwRyrbt68mccee4zZs2dzySWXADB79mzCw8MLPAG8/PLLBbr83MjNVk0D/muM2SAi5YD1IjLPnvaeMebtzDPbQ5fdCTTHGmB5vj2UHVij2lyPNUrWOhGZbYzZmh8rolReHYk/yZg5W/llYyQel4u+l9amX4fXaJH6I0FT7dH89n5k/b/N+hc7shYrfC2YcKwXI35O4oP5O3nsmsbce3kwLpfWsedVcnIyDz30ECEhIXg8Ht59912uvvpqJk2axJw5c0hOTiYxMZERI0YwcuRIatSoQWhoKH379qVly5Z88MEHJCUlMWvWLBo1asQvv/zCmDFjSElJoUqVKkybNo0aNWpk+/5vvPEGL7zwQsbOH6B3794Zj8ePH8+4ceNISUmhcePGTJ06lTJlyjBo0CBuuukmbrvN+mKULVuWhIQEDh48SL9+/Th+/DhpaWl89tlndO7cmSFDhhASEoKIMHjwYJ566qlTlvHyyy/zyy+/kJSUROfOnfniiy8QEbp160bHjh1ZtGgRsbGxfPnll1x55ZX5tv3PmgCMMQexBmzGGBMvItuwhjnLTh/gG2PMSWCPiOwCLrOn7TLGhAGIyDf2vJoAVKH7ffMhXvjpbxJOpnF/14YMuaIB1cuXAg5CeHajP0JFSeRG9xp6udawdmAYHy3cxahftvLHlsO8fUdralcsXXgrcY5G/7KFrZHH83WZzWqVZ+TNzXOcJykpiTZt2gDQoEEDfvrpJz755BMA/v77b7Zv30737t3ZscM6y1q1ahWbNm2icuXKLF68mI0bN7Jt2zYqV65Mw4YNGTp0KGvXruWDDz7go48+4v3336dLly6sXr0aEWHChAm8+eabvPPOO9nGtGXLFp55JvvRYPv27cv9998PwIsvvsiXX37JY49lP0rr9OnTueGGGxg+fDher5cTJ04QGhrKgQMH2Lx5MwCxsbFnvO7RRx9lxIgRAAwYMIBff/2Vm2++GbDOgNauXcvcuXMZPXo08+fPz/b98ypP51UiEgy0xRrD9ArgURG5F2tczf8aY45hJYfMY71G8G/C2H9aeUdOIyLDgGFAgZ+CqWJmVAVOmJKMSBvETO9VtKhdnvfuaEOTGuXyvCgR6NiwCpc1qMx3Ift5+Zet9HhvKa/c0oJb2uZ0fFR8ZVUFtHz58owd6sUXX0z9+vUzEsD1119P5cqVM+bt0KEDNWvWBKBRo0Z0794dgJYtW7Jo0SLAut+hX79+HDx4kJSUlDy1j4+JieHaa6/lxIkTDBs2jGeeeYbNmzfz4osvEhsbS0JCAjfccEOOy+jQoQODBw8mNTWVW265hTZt2tCwYUPCwsJ47LHHuPHGGzPizmzRokW8+eabnDhxgqNHj9K8efOMBNC3b18A2rVrR3h4eK7XJzdynQBEpCzwA/CkMea4iHwGvAIY+/87wGAgq/NgQ9Ytjs4YjcYYMw4YB9C+fXsdrUblm2OmLPelPMsm04jHrmnMY9c0oYTn/C7kigj9OtTj8oZV+e/3oTz5bSj7j57g0WsaF9lml2c7Ui9MOQ1IdcEFF5zyvGTJkhmPXS5XxnOXy0VaWhoAjz32GE8//TS9e/dm8eLFjBo1Ksf3b968ORs2bKB169ZUqVKF0NBQ3n77bRISEgAYNGgQs2bNonXr1kyaNInFixcD4PF48Pl8GeuQkpICQNeuXVm6dClz5sxhwIABPPvss9x7771s3LiRP/74g08++YTvvvuOiRMnZsSQnJzMww8/TEhICHXr1mXUqFGntO1PX0+3252xnvklV99+EQnC2vlPM8b8CGCMOWyM8RpjfMB4/q3miQDqZnp5HSAyh3KlCtyB2CRuSxnJVlOfz4Le57/dm/678x9V4d+/8OXntPx6H9Vk+sGb6OtaxjvzdjBq9hZ8Pj1+OZuuXbsybdo0AHbs2MG+ffto2rTpOS8vLi6O2rWtM7DJkyefdf7//e9/jB07lm3btmWUnTjxb+uv+Ph4atasSWpqakacYHVVv379egB+/vlnUlNTAdi7dy/Vq1fn/vvvZ8iQIWzYsIHo6Gh8Ph//+c9/eOWVV9iwYcMpMaTv7KtWrUpCQgIzZ848x7XPu9y0AhLgS2CbMebdTOU17esDALcCm+3Hs4HpIvIu1kXgJsBarDODJiLSADiAdaH47vxaEaWysysqngFfriXBVGJKidfp5Npu7exzI3yZ9X/UTWedNUi8vB30OZXTjjNh1Y3EJKbwXr82BGlz0Ww9/PDDPPjgg7Rs2RKPx8OkSZNOOdLPq1GjRnH77bdTu3ZtOnXqxJ49e3KcP/1C8r333kt8fDxVqlShXr16jB49GoBXXnmFjh07Ur9+fVq2bEl8fDwA999/P3369OGyyy7j2muvzThbWbx4MW+99RZBQUGULVuWKVOmcODAAe67776MM4bXXnvtlBgqVqzI/fffT8uWLQkODqZDhw7nvP55ddYxgUWkC7AM+BurGSjAC8BdQBusapxw4IH0hCAiw7Gqg9Kwqox+s8t7Ae9jNQOdaIwZm9N7t2/f3uiAMOp8RMYm8Z/PVpLqNUxJeYpmrn05v2BQovV/kl39MOlXuzyLBDAqLtPjfxOKMfD51Rt44/ft3NauDm/d1srx6qBt27ad0tJFBY6sPlsRWW+MaX+21+amFdBysq7Xn5vDa8YCZ+zc7fsEsn2dUvkpLimVQV+tJT45jW8f6ESzcWfZ+edVNmcRIvBQt0YkpXr5cMFOLixfimduOPdqDaUKit4JrAJScqqXYVNC2BOdyOT7LqN5rVxW+eSjp65rQtTxZD5etIsaFUoxoFP9Qo9BqZxo5aQKOMYY/jdzE2v2HOWdO9rQuXFVR+IQEcbc0oJrXRsYOWsTC1+6ypE4lMqOJgAVcKas2svsjZE8e0NTereu5WgsHreLj4M+5BLZy1OpD7P/6Imzv0ipQqIJQAWU0P2xjJmzlesuqc5DVzVyLpBMTUtLSwqfBn2ADxePTt9ASprv7K9XqhDoNQAVMGJPpPDItA1UL1eKt29vjevlik6HlKG+K4q3gr7gwYineHXuNkb1Ljo3Y6niSxOA8n+jKmAM/Df1GaKkHTMf7EzFMiWcjuoMPdzrGNKxAV8u30OH4Mrc2Kqmc8Hk9j6IXC8vLsfJ6V1Bp6am4vF4GDhwIE8++SQul4uQkBCmTJnChx9+mOVrw8PDWblyJXffnfVtQ5GRkTz++OPMnDkzx26fszNp0iS6d+9OrVpWdeHQoUN5+umnadasWa6X4a+0CkgFhG+8V7PAdykv9LqE1nWLzpH/6Z7veTFt6lZk+Ky/iYovPgO0p/cDtGXLFubNm5fRsRlA+/bts935g5UApk+fnuW0tLQ0atWqdV53z06aNInIyH87JZgwYUKx2PmDJgAVAA6YKoxN68/lri0MvDzY6XByFOR28c4drUlK8TL8p8059oUTqKpXr864ceP4+OOPMcawePFibrrJutFuyZIltGnThjZt2tC2bVvi4+N5/vnnWbZsGW3atOG9995j0qRJ3H777dx88810796d8PBwWrRokbH8/fv306NHD5o2bZqRZE6f5+2332bUqFHMnDmTkJAQ+vfvT5s2bUhKSqJbt26k34A6Y8YMWrZsSYsWLXjuuecyXl+2bFmGDx9O69at6dSpE4cPHy6MTZfvNAEov2aM4fnU+/Hh4k3POKveP/0CbBHVqFpZnunelHlbD/NzaPHsDqthw4b4fD6ioqJOKX/77bf55JNPCA0NZdmyZZQuXZrXX3+dK6+8ktDQUJ566inA6ip68uTJLFy48Ixlr127lmnTphEaGsr3339PTr0J3HbbbbRv3z5j/tKl/+3OOzIykueee46FCxcSGhrKunXrmDVrFgCJiYl06tSJjRs30rVrV8aPH58fm6XQaQJQfu2bdftZ5mvF/3mmU9d15OwvKCIGd2nApfUqMnL2FqKOF5+qoMyyOvu54oorePrpp/nwww+JjY3NdiSw07uKPn1alSpVKF26NH379mX58nPr4G/dunV069aNatWq4fF46N+/P0uXLgWgRIkSGWctBdFNc2HRBKD8VmRsEmPnbONy1xb6uxc4HU6euF3CW7e3JinVy/BZm8/+ggATFhaG2+2mevXqp5Q///zzTJgwgaSkJDp16sT27duzfP3pXUVndnq/SyJySvfNwCndLWcnp+q5oKCgjPcpiG6aC4smAOW3xszZSprPZ1X9iP/VpTeqVpanrruIeVsPM3+rf9Yhn4sjR47w4IMP8uijj56xs969ezctW7bkueeeo3379mzfvp1y5cpl9MKZG/PmzePo0aMZQ0VeccUV1KhRg6ioKGJiYjh58iS//vprxvzZLb9jx44sWbKE6OhovF4vM2bM4KqrAutubm0GqvzSsp1HmPv3IZ7pfhF1l/pP1c/phnRpwA8bIhj96xa6NKlKqaBCGmD+LM0281v6cJDpzUAHDBjA008/fcZ877//PosWLcLtdtOsWTN69uyJy+XC4/HQunVrBg0aRKVKlXJ8ry5dujBgwAB27drF3XffTfv2VqeYI0aMoGPHjjRo0ICLL744Y/5Bgwbx4IMPUrp0aVatWpVRXrNmTV577TWuvvpqjDH06tWLPn365NMWKRrO2h20k7Q7aJWVk2leer6/DJ8x/PFUV0qOybou+JzkpTvo/DAqjpW7orl7whqevK4JT153UYG8jXYHHbjOpztorQJSfufL5XsIi05kVO/mlPQU0hFzAercuCo3tarJZ4t3sy9G+wpShUcTgPIrB2KT+GjBLm5oXoNuTauf/QV+4sUbm+F2CS//usXpUFQxoglA+Qe7bf8bb76Czxheuimw7tS8sEIpnri2CfO3RbF0h/9e01D+RROA8hubfA2Y7buC+69sSJ1KZZwOJ39k6jV00BXB1K1cmlfnbsOrA8qrQqCtgJRfMAbGpvanCnE8sPIqWJXkdEj5rqTHzf9uuJjHZvzFjxsiuL19XadDUgFOzwCUX1joa8sa04wnPD9STgJv55/uplY1aV23Iu/8uYOkFK/T4agApwlAFXlpXh+vpd1FQ4nkLveZfb8EEhFheK9LOHQ8mYkr9jgdTr4REQYMGJDxPC0tjWrVqmV0p5BbwcHBREdH53d4xZYmAFXkfRcSwS5Th/95viFIAv+o+LIGlbm+WQ0+W7yb6ISTToeTLy644AI2b95MUpJ19jZv3jxq167tcFRKE4Aq0pJTvbw/fwft5B9ucBWfmwKf63ExSalePl202+lQ8k3Pnj2ZM2cOYHWzfNddd2VMO3r0KLfccgutWrWiU6dObNq0CYCYmBi6d+9O27ZteeCBB07pn+frr7/msssuo02bNjzwwAN4vYF/cJDfNAGoIu3r1XuJij/Js0Hfclq3MQGtcfWy9G1bm6/X7OVQXEH0Ftotn//O7s477+Sbb74hOTmZTZs20bFjx4xpI0eOpG3btmzatIlXX32Ve++9F4DRo0fTpUsX/vrrL3r37s2+ffsA6+7Xb7/9lhUrVhAaGorb7WbatGl5WH8FmgBUEXYiJY3Pl+zmisZV6OTKulfIQPb4tU3w+QyfLNrldCj5olWrVoSHhzNjxgx69ep1yrTly5dnXCO45ppriImJIS4ujqVLl3LPPfcAcOONN2b0A7RgwQLWr19Phw4daNOmDQsWLCAsLKxwVygAaDNQVWRNXrmX6IQUvri+KXzldDSFr27lMtzRoS7frNvHA1fl970Pi/NxWbnXu3dvnnnmGRYvXkxMTExGeVZ9kqX3FHp6j6Hp8w8cOJDXXnut4IItBvQMQBVJ8cmpfLF0N1c3rUa7+jn3/hiQ7JvDHg3tjXhT+GhBYJwFDB48mBEjRtCyZctTyrt27ZpRhbN48WKqVq1K+fLlTyn/7bffOHbsGADXXnstM2fOzBhR7OjRo+zdu7cQ1yQw6BmAKpImLg8n9kQqT1/f1OlQCk8Ww1jWkqPc7V7A1A0leKhbI4KrZj8Qij+oU6cOTzzxxBnlo0aN4r777qNVq1aUKVOGyZMnA9a1gbvuuotLL72Uq666inr16gHQrFkzxowZQ/fu3fH5fAQFBfHJJ59Qv379Ql0ff6fdQasiJy4plS5vLOTylFWMK/Fe4b55YXcHnQtRpiJXej/nxlY1efeONue0DO0OOnBpd9AqoExZGU58chpPeH50OpQiobrEck+n+vwcGqndRat8pQlAFSmJJ9P4csUerr24Os1dWqebbljXhrhF+GxJ4NwXoJynCUAVKdPW7CX2RCqPXNPY6VCKlBrlS3FHhzrMXL+fg3GB2xeSKlxnTQAiUldEFonINhHZIiJP2OWVRWSeiOy0/1eyy0VEPhSRXSKySUQuzbSsgfb8O0VkYMGtlvJHyalexi3dQ5fGVbm0XjFs+XMWD3RthDHwxRJt767yR27OANKA/xpjLgE6AY+ISDPgeWCBMaYJsMB+DtATaGL/DQM+AythACOBjsBlwMj0pKEUwLfr9hOdcJJH9eg/S3Url+HWtrWZsXYfR+IDo48g5ayzJgBjzEFjzAYQ/RW+AAAgAElEQVT7cTywDagN9AEm27NNBm6xH/cBphjLaqCiiNQEbgDmGWOOGmOOAfOAHvm6NspvpaT5+HzJbjoEV6Jjg3wc5D3APNStEaleHxOW61mAOn95ugYgIsFAW2ANUMMYcxCsJAGkD9BaG9if6WURdll25ae/xzARCRGRkCNHdGi84mLWXwc4GJfMI1c3zvLOT2VpWK0sN7aqxder9hJ3ItXpcHItv7qD7tatG+lNw3v16kVsbGy+xlnc5DoBiEhZ4AfgSWPM8ZxmzaLM5FB+aoEx44wx7Y0x7atVq5bb8JQf8/kMXyzdTbOa5bnqIv3Mz+bBqxqSmOLl6zX+00qqILqDnjt3LhUrVsyP8IqtXCUAEQnC2vlPM8akN84+bFftYP+PsssjgMxj2dUBInMoV8Xc/G2H2X0kkQeuaqhH/7nQvFYFul5Uja9WhJOc6j9dIOfUHXRiYiKDBw+mQ4cOtG3blp9//hmApKQk7rzzTlq1akW/fv0yEgicOjjMLbfcQrt27WjevDnjxo3LmKds2bIMHz6c1q1b06lTJw4fPlwYq+o3ctMKSIAvgW3GmHczTZoNpLfkGQj8nKn8Xrs1UCcgzq4i+gPoLiKV7Iu/3e0yVcx9sTSMOpVKc2PLmk6H4jce7NqQ6IST/LjhwLktoFs+/+VCTt1Bjx07lmuuuYZ169axaNEinn32WRITE/nss88oU6YMmzZtYvjw4axfvz7LZU+cOJH169cTEhLChx9+mNHRXGJiIp06dWLjxo107dqV8ePH5y7YYiI3fQFdAQwA/haRULvsBeB14DsRGQLsA263p80FegG7gBPAfQDGmKMi8gqwzp7vZWPM0XxZC+W31oUfZf3eY4zu3RyPW29LyVbmfoJGxXF5oyq0qlOBcUt3069DXdyuon/mlFN30H/++SezZ8/m7bffBiA5OZl9+/axdOlSHn/88YzXt2rVKstlf/jhh/z0008A7N+/n507d1KlShVKlCiRcZ2hXbt2zJs3r6BWzy+dNQEYY5aTdf09wLVZzG+AR7JZ1kRgYl4CVIHtiyW7qXxBCe5oXzfLztBU1kSEB7o24pHpG/hzyyF65vXsaXGBhHVWOXUH/cMPP9C06Zmd/52tWnDx4sXMnz+fVatWUaZMGbp160ZysjWITlBQUMbr3W43aWlp+bg2/k8PuZRjdhyOZ/62KAZeHkzpEm6nw/E7PVpcSP0qZfh8ye4s+9MvirLrDvqGG27go48+yliPv/76Czi1m+jNmzdnDBWZWVxcHJUqVaJMmTJs376d1atXF/BaBA5NAMox45eGUTrIzb2Xaxe+58LtEu6/siEbI+JYs8c/alOz6w76pZdeIjU1lVatWtGiRQteeuklAB566CESEhJo1aoVb775JpdddtkZr+3RowdpaWm0atWKl156iU6dOhX4egQK7Q5aOSIqPpkury/iTn7n5aBJTofzryLYHfQZRsVlPExO9XL5awtoV78yEwZm3/uvdgcduM6nO2gdEEY5YuqqvaT6fAwO+s3pUPxPpmslpUbFMaBTfT5atIuwIwk0rFbWwcCUv9EqIFXoklK8fL16L9dfUoNgl7bLPi+jKjBgxXUEmRS+fP8lp6NRfkYTgCp0P2yI4NiJVIZe2dDpUAJCNTnOre4V/OC9kqOJKdnOV5Sre9W5Od/PVBOAKlQ+n2Hi8j20rlOBDsHaGWx+GeKeSzIlmbY66+4hSpUqRUxMjCaBAGKMISYmhlKlSp3zMvQagCpUC7dHERadyId3tdVuH/LRRa4DXOUKZfKqkgy7qiElPac2q61Tpw4RERFoB4uBpVSpUtSpU+ecX68JQBWqCcvDqEU0vX5sDj/5nA4noAx1z2VAQht+Do20bqzLJCgoiAYNGjgUmSqqtApIFZotkXGsDjvKIM8feER3/vmti2szTWuU46sV4VrVo3JFE4AqNF+tCKd0kJt+7sVOhxKQROC+K4LZdvA4q8P848Yw5SxNAKpQRCecZHZoJLe1q0MFSXQ6nIB1S9vaVCoTxFcr9jgdivIDmgBUoZi2eh8pXh+Drgh2OpSAVmpsZfqf/I55Ww+yL+aE0+GoIk4TgCpwJ9O8TF29l6ubVqOR3qla4AZ45uHGx6SV4U6Hooo4TQCqwP06ug/RCScZHPaUdvlcCGpILDe5VvNdyH7ik/1n3GBV+DQBqAJljGFiWg+aSARdXJudDqfYGOz5jYSTaXwfEuF0KKoI0wSgClTI3mNsMQ0Y5P4dve+r8LRy7aFd/UpMXhWOz6dNQlXWNAGoAjVpRTgVSOBW9wqnQyl2BnUOZm/MCRbviHI6FFVEaQJQBSYyNonftxziTvciyshJp8Mpdnq0uJALy5fiqxXhToeiiihNAKrAfL16L8YYBnh0IG4nBLldDLi8Pst2RrMrKt7pcFQRpAlAFYjkVC8z1u6je7MLqSPRTodTbN3ZoS4lPC5tEqqypAlAFYifQw9w7ESq3vjlsCplS9KndS1+WH+AuCRtEqpOpQlA5TtjDF+tCOfiC8vRsUFlp8Mp9gZ2DiYp1cv3IfudDkUVMZoAVL5bs+co2w/FM6hzsPb5XwS0qF2By4IrM3lVOF5tEqoy0QSg8t3kleFULBPELW1rOx2Ksg3sHMz+o0ks2q5NQtW/NAGofHUgNok/thyi38kfKDW2snb9UER0b16DC8uXYvKqcKdDUUWIJgCVr762x6Qd4JnvcCQqsyC3i3s61dMmoeoUmgBUvklO9fLN2n1c36yGNv0sCkZVOOXvziXXUIJUJn/wotORqSJCE4DKN7M3RnLsRCoDOwc7HYrKQlU5zs2ulfzg7cpx7SVUoQlA5RNjDJNWhNO0Rjkub1jF6XBUNgZ5/uQEpZipvYQqNAGofBKy9xhbDx5nYMx7yOiKToejstHStYdLZQdTVmkvoSoXCUBEJopIlIhszlQ2SkQOiEio/dcr07T/E5FdIvKPiNyQqbyHXbZLRJ7P/1VRTpq0MpzyJHKL9vpZ5A30/EF4zAmW7DzidCjKYbk5A5gE9Mii/D1jTBv7by6AiDQD7gSa26/5VETcIuIGPgF6As2Au+x5VQA4GJfE75sP0U97/fQLPV1rqVauJJO1f6Bi76wJwBizFDiay+X1Ab4xxpw0xuwBdgGX2X+7jDFhxpgU4Bt7XhUApq3eh88YBri16ac/KCFe7kn6msX/HCFsRFOnw1EOOp9rAI+KyCa7iqiSXVYbyNzhSIRdll258nPpvX5ee3EN6rn0LlN/cZd7IUGkMcXb3elQlIPONQF8BjQC2gAHgXfs8qw6fjE5lJ9BRIaJSIiIhBw5onWURd2cTQeJSUxhkDb99CvVJY4bXauZ6e1Kwsk0p8NRDjmnBGCMOWyM8RpjfMB4rCoesI7s62aatQ4QmUN5VsseZ4xpb4xpX61atXMJTxUSYwyTV4XTuHpZrmisTT/9zSDPHyRQhh/Wa5PQ4uqcEoCI1Mz09FYgvYXQbOBOESkpIg2AJsBaYB3QREQaiEgJrAvFs889bFUU/LU/lk0RcQy8vL72+umH2rh201p26cDxxVhumoHOAFYBTUUkQkSGAG+KyN8isgm4GngKwBizBfgO2Ar8DjxinymkAY8CfwDbgO/seZUfm7wynHIlPfS9tI7ToahzNMjzB2FHElm2S7vuKI48Z5vBGHNXFsVf5jD/WGBsFuVzgbl5ik4VWYePJzNn00EGXF6fC0qe9WukiqgbXasZW/ZpJq8M56qLtMq1uNFfrjon09bsw+vzMnDdrbD+sNPhqHNUQrz071iPDxbsZE90Ig2qXuB0SKoQaVcQKs9OpnmZvmYv17hCCXbpzt/f9e9YjyC3MGVVuNOhqEKmCUDl2ZxNB4lOSGGQ+3enQ1H5oHr5UtzYsibfh0Rok9BiRhOAypP0Ad8bVy9LF9fms79A+YVBVzQg4WSaNgktZjQBqDzZsC+Wvw/EMbBzMNryM3C0qVuRNnUrMnmlNgktTjQBqDyZtDKccqU89NUB3wOHPWLYfYfGEBadyFLtJbTY0ASgcu1QXDK/bdxHv9SfueA1vfM30PR0raF6uZJ8tSLc6VBUIdEEoHJt6upwfLgY6P7T6VBUASghXu7pVJ8lO46wKyrB6XBUIdAEoHIlOdXL9DX7uN4VQl2XVhEEqrs71qOEx8WklXucDkUVAk0AKldm/XWAYydSuc+jTT8DWdWyJenTuhY/rD9A3AkdOD7QaQJQZ2WMYeKKPTSrWZ6Ost3pcFQBu++KBiSlevlm3T6nQ1EFTBOAOqsVu2LYcTiB+67Qpp/FQbNa5enUsDKTV4aT5vU5HY4qQJoA1Fl9tWIPVcuW4ObWtZwORRU0u0no4P3DiYxL5s+t2tVHINMEoHK0JzqRhf9E0b9jfUoFuZ0ORxWSa10bqCeHmbhcLwYHMk0AKkdfrdhDkMtF/071nA5FFSK3GAa5/yBk7zFC98c6HY4qIJoAVLZiT6TwfUgEfVhI9XdqWNUDqti4w72YciU9fKlnAQFLE4DK1vS1+0hK9TLE/ZvToSgHlJVk7upYj7l/H+RAbJLT4agCoAlAZSklzcfkleFc2aQqF7v2Ox2OcsjAzsGANfynCjyaAFSW5vwdyeHjJxnSpYHToSgH1a5Yml4tazJjzT4dKyAAaQJQZzDGMGHZHhpXL6vjxCqGdGlA/Mk0vlunZ4KBRhOAOsPqsKNsiTzOkC4NEL3zq9hrU7ci7etXYuKKPXh1rICAoglAnWH8sjAqX1CCW7XPf2XfGDY0cgQRx5L4ffMhpyNS+UgTgDrFzsPxLNwexcDLg/XGL5XhelcIwVXKMG7pbozRs4BAoQlAnWLc0jBKcZIBy67OOPpTyi2GoVc2ZGNEHGv3HHU6HJVPNAGoDIePJzMr9AB3uJdQWeKdDkcVMbe1q0PlC0owbmmY06GofKIJQGX4akU4Xp9hqHuu06GoIqjU2Mrce3I6C7ZHsfOwHiAEAk0ACoCEkTWYtmQTPWU19VxRToejiqh73fMoxUnGL9OzgECgCUAB8I33GuK5gGGeX50ORRVhlSWe291LmPVXJFHHk50OR50nTQCKlDQfX6b1pKNspbVLj+xUzoa655Lm8/HlCu0kzt9pAlDMCj3AQarwkGe206EoP1DfFUWvljWZtnofcUk6brA/0wRQzHl9hs+X7Ka57OEq1yanw1F+4qFujUg4mcbUVeFOh6LOgyaAYu7PLYcIO5LIQ57ZOt6vyrXmtSrQrWk1vloRTlKK1+lw1Dk6awIQkYkiEiUimzOVVRaReSKy0/5fyS4XEflQRHaJyCYRuTTTawba8+8UkYEFszoqL4wxfLZkN8FVytDTtdbpcJSfebhbY2ISU/guRDuJ81e5OQOYBPQ4rex5YIExpgmwwH4O0BNoYv8NAz4DK2EAI4GOwGXAyPSkoRwyqgIrRnRhU0QcD8R9gFv09n6VN5c1qEz7+pUYtzSMVK/P6XDUOThrAjDGLAVOv/e7DzDZfjwZuCVT+RRjWQ1UFJGawA3APGPMUWPMMWAeZyYVVcg+9famOsfo617mdCjK39jdhDwc+X8ciE1idmik0xGpc3Cu1wBqGGMOAtj/q9vltYHM54MRdll25WcQkWEiEiIiIUeOHDnH8NTZrPc1YaWvBfd75lBSdKAPdW6udoVysezlk8W7tKtoP5TfF4Gzuoxocig/s9CYccaY9saY9tWq6WAkBeWjtFupzHH6uxc4HYryYyLwmGcWYUcSmfv3QafDUXl0rgngsF21g/0/ve+ACKBupvnqAJE5lCsHbNwfy2JfG4Z65lBGTjodjvJzPV1raVK9LB8t3IlPzwL8yrkmgNlAekuegcDPmcrvtVsDdQLi7CqiP4DuIlLJvvjb3S5TDvho4S4qkMC97nlOh6ICgEsMj17TmB2HE/hzqw4Y409y0wx0BrAKaCoiESIyBHgduF5EdgLX288B5gJhwC5gPPAwgDHmKPAKsM7+e9kuU4VsS2Qc87cdZojnN8qK9uWi8sdNPzWnoUTy4bQfdcAYP+I52wzGmLuymXRtFvMa4JFsljMRmJin6FS++3jhLsqV9DBQT8BUPnKL4WHPzzyT+hALtkVxXbMaToekckHvBC5G/jkUz2+bD3HfFcFUkBNOh6MCTB/XSurJYT5cuFPPAvyEJoBi5L15OyhX0sPgLg2cDkUFoCDx8qh7Fpsi4pi/TceU8AeaAIqJzQfi+H3LIYZ4v6Him9q8VhWMvu5lBFcpwzt//qMtgvyAJoBi4t15O6hAAoPdvzkdigpgHvHx5HUXsd2ublRFmyaAYmDDvmMs3B7FMM+vlJckp8NRAe7m1rVoUr0s783foXcHF3GaAALdqAq89/nnVCGOQW5t+aMKntslPHndReyKSuCXjXq/Z1GmCSDArfU1ZZmvFQ96fuECvetXFZKeLS7k4gvL8f78HdpTaBGmCSCAGWN4I/VOqnOMe9zznQ5HFSMul/Df7k0JjznB9yERToejsqEJIIDN23qY9aYpT3lmUlpSnA5HFTPXXVKd9vUr8f78HTpqWBGlCSBApXl9vPnHPzSUSG53L3E6HFWc2GMFyOiKPHfwSaLiTzJxxR6no1JZ0AQQoH7ccIBdUQn8z/MtHtE6WOWMDq5/uO6S6ny+eDfHEvUstKjRBBCAklO9vDtvB23qVuQG1zqnw1HF3LM3XExiShqfLt7ldCjqNJoAAtCkleEcOp7M8z0vRrIaikepQtT0wnL0vbQOk1fuJeKY9kFVlGgCCDDRI+vyyW8buNr1F52mNHQ6HKUAePr6ixCBN3//x+lQVCaaAALMe2m3cYKSDPdMczoUpTLUqliaYV0bMntjJOv3HnM6HGXTBBBA/jkUzwzvNQxwz6OxS+/AVEWE3SrowRVXUr1cSV75dat2FFdEaAIIEMYYxszZSjlO8ITnR6fDUeoMF8hJnr2hKaH7Y/llkx6gFAWaAALBqAosGtGNZTujecLzI5UkwemIlMrSfy6tQ4va5Xn9t+16c1gRoAkgAKQYN2PS7qGhRDJAB3pXRZjLJbx0YzMOxiXz+ZLdTodT7GkCCAATvL0IM7V4yfM1QaJHVapo69iwCje1qsnnS3azL0abhTpJE4CfOxCbxEdpt9LdtY6r3aFOh6NUzuwLwi/+0xe3S3j51y1OR1SsaQLwc2N+3YpBGBE01elQlMq1C+UYT17XhPnbopi/9bDT4RRbmgD82JIdR/ht8yEe88yijkQ7HY5SeXLfFQ1oUr0so3/dQnKqVl06QROAnzqZ5mXU7C00qHoBQ91znA5HqTwLcrsY3ac5+48m8eki7SfICZoA/NTHC3exJzqR0b2bU1LSnA5HqbwbVYHOUxvRx7WCzxZuZ+fheKcjKnY0Afih7YeO89ni3fRtW5uuF1VzOhylzstLQVO5gCSe+2GTDiJfyDQB+Bmvz/DcD39TvnQQL97UzOlwlDpvVeU4I4KmsmFfLF+v3ut0OMWKJgA/M2llOBv3xzIy5R0qv1XNalanlJ+71bWcK5tU5c3ft3MgNsnpcIoNTQB+ZP/RE7z9xz9c3bQavV2rnA5HqXwjAq/e2hIDDP/pb4zRqqDCoAnAT/h8hme+34jbJYy5taUO9KICTt0Pa/KsbyKL/znC9yERTodTLGgC8BMTV+xhzZ6jjLi5GbUrlnY6HKUKxED3n3RybWH0L1vYf1S7iShomgD8wI7D8bz5xz9cd0kNbm9Xx+lwlCowLjG8HfQFIsJ/v9+orYIK2HklABEJF5G/RSRURELsssoiMk9Edtr/K9nlIiIfisguEdkkIpfmxwoEupQ0H099G0q5tGO8FnYrMrqiXvhVAa2ORDPS+zFr9xxl4oi7nQ4noOXHGcDVxpg2xpj29vPngQXGmCbAAvs5QE+gif03DPgsH9474H2wYAdbIo8zNuhLqslxp8NRqlDc5l5Kd9c63krrx7aD+r0vKAVRBdQHmGw/ngzckql8irGsBiqKSM0CeP+AsXxnNJ8u3s3t7erQwx3idDhKFRoReDXoSyqQwKPTN5B4Uu92LwjnmwAM8KeIrBeRYXZZDWPMQQD7f3W7vDawP9NrI+yyU4jIMBEJEZGQI0eOnGd4/isqPpknvw2lUbWyjO7T3OlwlCp0VeU4HwR9Qlh0IiN+1m6jC8L5JoArjDGXYlXvPCIiXXOYN6uGi2dc4THGjDPGtDfGtK9WrXh2c+D1GZ76NpT45FQ+uftSypTwOB2SUo7o7N7KY9c04YcNEcxcr01D89t5JQBjTKT9Pwr4CbgMOJxetWP/j7JnjwDqZnp5HUBHhs7Cp4t2sWJXDKN7N6fpheWcDkcpRz1xbRM6NqjMS7M2sytKO4zLT+ecAETkAhEpl/4Y6A5sBmYDA+3ZBgI/249nA/farYE6AXHpVUXqX4v/ieLd+Tvo06YW/ea0zBhBSaniyv1yRT6MvJMyqUd54L0ZxCenOh1SwDifM4AawHIR2QisBeYYY34HXgeuF5GdwPX2c4C5QBiwCxgPPHwe7x2QwqMTeXzGXzStUY7X+urdvkqlqyGxfBz0IeHmQp76diM+vT8gX5xz5bIxJgxonUV5DHBtFuUGeORc3y/QJZ5MY9jUEFwuYdyA9lrvr9RpLndv4yUzlVHbBvHBS0N4KugHGBXndFh+TfcyRUB6Pz+7ohKY4nmVeh9tdjokpYqkge4/2Wwa8IH3PzRz7eUGpwPyc9oVRBHw9p//8NvmQ/xfz0vo4tadv1LZEYExnom0ll08mfowmyJinQ7Jr2kCcNi0NXv5dPFu7u5Yj6FXNnA6HKWKvFKSyvgS71CZeAZPCtFO486DJgAHLdoexUuzNnN102q83Ls5old9lcqV6hLH5BJvkJLmZdBXa4k7oS2DzoUmAIdsiojlkekbaMYePg6/Cc8rlbS5p1J50NgVyTjfCPYfieX+Vz4gOdXrdEh+RxOAA7YfOs69E9dS+YISTCzxFhfISadDUsovdXJt5+2gz1lnmvLQ1+tJSfM5HZJf0QRQyMKOJHDPhDWU8riZPrQT1UUvYil1Pnq7VzHWM5FF/xzhiW/+Is2rSSC3NAEUov1HT9B/whqMga+HdqRelTJOh6RUQLjbs5AXb7yE3zYf4n8zN+mNYrmk9wEUkj3RifQfv5rEk2l843uGxp/uczokpQLK0AVtSfLcwjt/3YHZ+A1vvfIqHrce4+ZEt04h2HE4nju+WEVymo/p93eimUt3/koVhMc8s3jG8y0/+a7k0el/6TWBs9AEUMA2H4ij3xerEODbYZ1oUVtb+ihVkB71/MxLnin8vuUQw0a8RvLIqk6HVGRpFVBBGVWBJd5WPJz6BBVJYPqzt1O/ygVOR6VUsTDE8ztlOMkLaUO4O2U44xNOUqVsSafDKnL0DKCATE+7hsGpz1JPDjOz5Gjd+StVyO7yLOLToA/YYoLpO/ZrwkY01XttTqMJIJ/5fIbXf9vOC2lD6eL6m+9LvExNOep0WEoVSz3d65hRYgzxpgx9U0azxnex0yEVKZoA8lHsiRTum7SOz5fs5m73fL4MepuykmxNTB/YRY9AlCpUl7p28VOJEVSWePqnvMBXK/Zg9U6vNAHkk80H4rjpo+Ws3B3N2FtbMNYzEY9oCwSlioL6rih+KjGCbq6NjP5lK09+G8qJlDSnw3KcXgQ+T8YYvl6zjzG/bqXyBSX47oHLaVuvEvzmdGRKqcwqyAnGBb3LZ97evB16O9s2ruGDoI+5xLX/35mK2QAzegZwHqLikxk8aR0vzdpMp4ZV+DV5IG0nBms1j1JFlEsMj3h+ZkrQ6xwzZemTMobxab3wmeLZE68mgHNgjGHu3wfp8f4yVu6O4eU+zZl0XweqSLzToSmlcuFK92Z+L/k83VyhjE27h/6pL7DPV93psAqdVgHl0YHYJEb+vJn526JoUbs87/drQ+Pq5ZwOSymVR1Ukni+C3uM7bzdeThtA95Q3eOLFexjqnkuQeItFdZAmgFxKSfMxeWU4783fgTEwvNcl3HdFsPY1opQfE4F+nsV0dW9iZOpA3ki7i5+9nXklaBIdTq/KDcCEoAngLIwx/L75EK//vp29MSe4xrWB0Z5J1F0YDQudjk4plR9qylHGlXiPP7ztGZU6kNtTRtLDtZbnPTMIdh12OrwCowkgG8YYVoXF8N68HawLP8ZFNcoy6b4OdJtxt9OhKaUKyA3uELq6NjHeeyOfp93MgpRL6e+ez0Oe2dRwOrgCoAngNMYYVu2O4f35O1kbfpTq5Uoy1jOBfrGL8czQdv1KBbrSksLjnp+4072Id9NuY6r3eqZ7r+GunzfzULfGXFih1Kkt/fy4akgTgC0lzcecvyOZuDycvw/EcWH5Uozu3Zx+HepSaux/nA5PKVXIqkssrwdN4GH3bD7x9mHampJMX7uPm1vVYrAvmBaucKdDPG/FPgHsP3qC79dH8M2CtURRiUZygLGe3/jPszMpFeR2OjyllMPquaJ4wzWeR32z+NLbk+//uoofeZXLZBt3eRbSM9Xrt/sKKcp9YrRv396EhITk+3ITTqYxf+thZq6PYMXuaACulI0Mdv9GV9ffuKTobhNVwAYlWv8n2b23TvrVLr/JmXhUkXPclOY779VM8V7PPlODciTS272KW93LuVR24hrt/DjfIrLeGNP+bPMVmzOAuKRUlu44wpxNB1n0TxQn03zUrliaJ65twm3t6lDnA724q5Q6u/KSxFDPXAa7f2O17xK+83Zjprcr07zXUYtoev26lZ4tL6RN3Uq4XUX7DuOATQA+n2HboeMs3xnNwu1RhOw9htdnqMYx7nKv4cYSa2iXtAPXcgPLnY5WKeVvXGLo7N5KZ/dWXjFfMd93KXO8nZi8vCITlu+hEvFc1aYp3ZpWp3OjKlQvX8rpkM8QMAngZJqXrZHH2bAvlpDwo6wKiyH2RCoAF8s+HnD9xTUl/qKt7MStVTxKqXxUTpK41b2CW6crue4AAAeBSURBVN0riDNlWOprxSJvWxaFwqzQSAAaVy/L5Q2r0K5+JS6tV4m6lUsj4uwZgl8mgISTaew8HM/Wg8fZEmn9bTt4PGMA6Noc4Tr3VjoHbeFy11YdkEUpVWgqyAludq/mZvdqvEbYYoJZ5WvOyuhm/BjVlKmrSwNQtWwJmteqQPNa5WleqwKX1CxHvcplCrV3gSKdAFK9PhZsO0x4zAn2xiSyJzqR3VEJRMYlZ8xTrpSH5rXKM6hzMJfWq0jbepWo8W4g3rKhlPI3bjG0kj20cu3hAX7Fa4R/TF3+8jXmr6QmbN4ZzIodtUmzd8UlSKWBHKSRHKS+HCJYDlN/8ERa16lI6RL539Ko0BOAiPQAPgDcwARjzOvZzRuTmMKQyVYroHIlPdRP2UlHiaCx5wCN73qLZjXLU6dSaWR0RYgE1hbKKiil1Dlxi6GZ7KOZax/97b5kkk0QO00d/jF12OmrzS5Th62mPn/62lmJYdxq5pd4hsauyHy/6axQE4CIuIFPgOuBCGCdiMw2xmzNav7KCbuYVGIEwXKYSsQjJTNN/H52IUSslFIFq5Sk0lL20JI91mGxLc24iDRV2WuqU0/s/ojyeayRwj4DuAzYZYwJAxCRb4A+QJYJoASpXOraVYjhKaVU0eARH/UkinpEFdh7FOqNYCJyG9DDGDPUfj4A6GiMeTTTPMOAYfbTpkAMEF1oQRaeqgTeegXiOkFgrlcgrhME5nqdyzrVN8ZUO9tMhX0GkFWbp1MykDFmHDAu4wUiIbm5o83fBOJ6BeI6QWCuVyCuEwTmehXkOhX2aCYRQN1Mz+tgXb5VSilVyAo7AawDmohIAxEpwf+3d3YhWpRRHP/9MxSsrNWwIgrsg0xCQyskTOkDvy7aTAojSrKbokAFwRWjGwmyj5uuIlSym5LAyIgyE8pAV1TwYyV1dytQEoWyIgJt4XTxnLcdlnd293VXZ2bf84NhZp85s8yPM7xnZp6ZZ2AJEL25QRAEBXBZbwGZWY+kV4HtpP7uTWZ2dIDNPhhgfVUZiV4j0QlGptdIdIKR6XXJnEo9GmgQBEFw6YgvmgdBEDQpUQCCIAialMILgKTlkjokHZW0wtu2SDro0y+SDmbi10jqknRc0rzi9rx/crzuldTuXvslPeDtkvSeex2WNL3Yva9PjtM0SXskHZH0haRxmfhS5krSJklnJXVk2sZL2iGp0+ct3p6bG0lLPb5T0tIiXLI06DXZ83Ze0qo+/2e+56xLUtvl9uizL404Pes5Oixpt6RpmW1K4+T704hXqzvVfjdmZbYZ2jFoZoVNwD1ABzCW1CH9LXBnn5h3gdd9eQpwCBgDTAK6gVFFOjTiBXwDLPCYhcB3meWvSO9JzAT2Fu3QgNM+YI7HLAPWlT1XwGxgOtCRaXsLaPPlNmB9f7kBxgM/+bzFl1sq5DURuB94A1iViR/luboNGO05nFIRpwdrOQAWZHJVKqeL8Lqa3v7aqcCx4ToGi74CuBtoN7N/zKwH+B5YVFspScDTwMfe1Ap8YmbnzexnoIs0vETZyPMyoHaGfC2970C0Ah9Zoh24TtJNl3unByDP6S5gl8fsABb7cmlzZWa7gL5jhLcCm315M/BEpr1ebuYBO8zsdzM7R3Kff+n3Pp9GvMzsrJntA/7tE///cC1mdgGoDddSCA067fZcALST3jOCkjlBw15/m//iA1fR+/LskI/BogtABzBb0gRJY0lnW9kXxR4CzphZp/99M3Ays/6Ut5WNPK8VwNuSTgLvAGs8vgpeeU4dwOMe8xS9+auCU5YbzOw0gM8nenueR1X88rzyqILXYJxeJF25QTWcoB8vSYskHQO+JF1pwzB4FVoAzOxHYD2pcn1NujTryYQ8Q+/ZPwxiKIky0I/Xy8BKM7sFWAls9E1K79WP0zLgFUkHgGuAC75J6Z0GSZ7HSPHrS+W9JD1MKgCra011wirlZGafmdlk0lXBOm8eslfRVwCY2UYzm25ms0mXRJ0Akq4EngS2ZMIrM5REjtdSYKuHfErvLZFKeNVzMrNjZjbXzGaQinW3h1fCKcOZ2m03n9eGYMzzqIpfnlceVfDKdZI0FdgAtJrZb95cBScYRK781tHtkq5nGLwKLwCSJvr8VtIPfu2M/zFSZ8epTPg2YImkMZImkTohS/kZmByvX4E5HvIIXuxIXs/7EyczgT9rl4Jlop5Tpu0K4DXgfQ+vTK6cbaQCjc8/z7TXy812YK6kFn9aY663lY08rzyqMFxLXSc/LrcCz5nZiUx8FZwg3+sO7w/Fn0IbTRoleejHYJE94d6v8QPpewCHgEcz7R8CL9WJX0s6yzyOP1FTxqmeFzALOOBte4EZ3i7Sh3K6gSPAfUXvfwNOy4ETPr2JP61Q5lyRivFpUgfoKdLtggnATlJR3gmMHyg3pNtfXT69UDGvGz3mL+APXx7n6xZ6PruBtRVy2gCcAw76tD/zf0rjdBFeq4Gj7rQHmDVcx2AMBREEQdCkFH4LKAiCICiGKABBEARNShSAIAiCJiUKQBAEQZMSBSAIgqBJiQIQBEHQpEQBCIIgaFL+AyQei39WsxKfAAAAAElFTkSuQmCC\n",
"text/plain": [
"