{ "cells": [ { "cell_type": "raw", "id": "eaf29ef7", "metadata": { "raw_mimetype": "text/restructuredtext" }, "source": [ ".. _cookbook-spb:\n", "\n", ".. include:: \n", "\n", "Light Curves of an SPB Star\n", "===========================\n", "\n", "This :ref:`cookbook` recipe demonstrates how to synthesize light curves of a slowly pulsating B-type (SPB) star, using the MSG Python interace. The light curves are based on surface perturbation data for a :math:`5\\,\\Msun` model, produced by the `GYRE `__ stellar oscillation code. \n", "\n", "Preparation\n", "-----------\n", "\n", "Download :download:`spb-light-curve.tar.bz2 <../downloads/spb-light-curve.tar.bz2>` and unpack it into your working directory using the :command:`tar` utility. This archive contains the following items:\n", "\n", "* :file:`summary.h5` --- surface perturbation data produced by GYRE\n", "* :file:`pygyre` --- Python module for reading :file:`summary.h5`\n", "* :file:`pg-*.h5` --- :f-schema:`photgrid` files for each passband\\ [#photgrid-files]_\n", "* :file:`gyre.in` --- GYRE input file used to create :file:`summary.h5`\n", "* :file:`spb.mesa` --- stellar model file used to create :file:`summary.h5`\\ [#model-file]_\n", "\n", "Note that :file:`gyre.in` and :file:`spb.mesa` aren't needed for the recipe; they're included for those who wish to recreate :file:`summary.h5`, or are curious about its provenance.\n", "\n", "Initialization\n", "--------------\n", "\n", "To start, import modules and set other parameters:" ] }, { "cell_type": "code", "execution_count": 1, "id": "2eb4244c", "metadata": {}, "outputs": [], "source": [ "# Import standard modules\n", "\n", "import numpy as np\n", "import matplotlib.pyplot as plt\n", "import scipy.special as ss\n", "import astropy.constants as con\n", "\n", "# Import pymsg and pygyre\n", "\n", "import pymsg\n", "import pygyre\n", "\n", "# Set constants\n", "\n", "SIGMA = con.sigma_sb.cgs.value\n", "\n", "# Set plot parameters\n", "\n", "plt.rcParams.update({'font.size': 12})" ] }, { "cell_type": "raw", "id": "1a08b07c", "metadata": { "raw_mimetype": "text/restructuredtext" }, "source": [ "Define a Synthesis Function\n", "---------------------------\n", "\n", "Now let's create a function that synthesizes a light curve using the semi-analytical approach described in Section 2 of :ads_citet:`townsend:2003`. The :py:func:`synthesize_light_curve` function, defined below, takes the following arguments:\n", "\n", "* ``photgrid_file`` --- name of :f-schema:`photgrid` file\n", "* ``phi`` --- oscillation phase array (radians)\n", "* ``i`` --- stellar inclination angle (radians)" ] }, { "cell_type": "code", "execution_count": 2, "id": "02d6f490", "metadata": { "editable": true, "slideshow": { "slide_type": "" }, "tags": [] }, "outputs": [], "source": [ "def synthesize_light_curve(photgrid_file, phi, i):\n", " \n", " # Create the PhotGrid object\n", "\n", " photgrid = pymsg.PhotGrid(photgrid_file)\n", " \n", " # Read the GYRE output file into a table\n", "\n", " tbl = pygyre.read_output('summary.h5')\n", " \n", " # Select the table row corresponding to the l=m=2, g_30 mode\n", " \n", " mask = (tbl['l'] == 2) & (tbl['m'] == 2) & (tbl['n_pg'] == -30)\n", "\n", " row = tbl[mask][0]\n", " \n", " # Extract stellar parameters\n", "\n", " M_star = row['M_star']\n", " R_star = row['R_star']\n", " L_star = row['L_star']\n", "\n", " Teff = (L_star/(4*np.pi*SIGMA*R_star**2))**0.25\n", " g = con.G.cgs.value*M_star/R_star**2\n", " \n", " # Extract surface perturbation parameters (the T03\n", " # annotations indicate which equation from Townsend 2003\n", " # is being implemented)\n", "\n", " l = row['l']\n", " m = row['m']\n", " \n", " omega = row['omega']\n", "\n", " dR_R = row['xi_r_ref'] # Delta_R (T03, eqn. 4)\n", " dL_L = row['lag_L_ref']\n", " dT_T = 0.25*(dL_L - 2*dR_R) # Delta_T (T03, eqn. 5)\n", " dg_g = -(2 + omega**2)*dR_R # Delta_g (T03, eqn. 6)\n", " \n", " # Evaluate photometric P-moments & their partial derivatives \n", " # (T03, eqn. 15; the I_{l;x} in that equation is equivalent to\n", " # P_l here)\n", "\n", " x = {'Teff': Teff, 'log(g)': np.log10(g)}\n", "\n", " P_0 = photgrid.P_moment(x, 0)\n", " P_l = photgrid.P_moment(x, l)\n", "\n", " dP_l_dlnT = photgrid.P_moment(x, l, deriv={'Teff': True}) * Teff\n", " dP_l_dlng = photgrid.P_moment(x, l, deriv={'log(g)': True}) * np.log(np.exp(1.))\n", "\n", " # Set up differential flux functions (T03, eqns. 12-14)\n", "\n", " R_lm = (2+l)*(1-l) * P_l/P_0 * ss.sph_harm(m, l, 0., i)\n", " T_lm = dP_l_dlnT/P_0 * ss.sph_harm(m, l, 0., i)\n", " G_lm = dP_l_dlng/P_0 * ss.sph_harm(m, l, 0., i)\n", " \n", " # Construct the light curve (T03, eqn. 11)\n", "\n", " dF_F = ((dR_R*R_lm + dT_T*T_lm + dg_g*G_lm) * np.exp(1j*phi)).real\n", "\n", " # Return it\n", " \n", " return dF_F" ] }, { "cell_type": "raw", "id": "5b4d9d5f", "metadata": { "editable": true, "raw_mimetype": "text/restructuredtext", "slideshow": { "slide_type": "" }, "tags": [] }, "source": [ "Evaluate Light Curves\n", "---------------------\n", "\n", "Using the :py:func:`synthesize_light_curve` function, evaluate light curves in each passband:" ] }, { "cell_type": "code", "execution_count": 3, "id": "affede99", "metadata": {}, "outputs": [], "source": [ "# Set up phase array\n", "\n", "phi = np.linspace(0, 4*np.pi, 101)\n", "\n", "# Set stellar inclination\n", "\n", "i = np.radians(75)\n", "\n", "# Evaluate light curves\n", "\n", "dF_F_Kepler = synthesize_light_curve('pg-Kepler.h5', phi, i)\n", "dF_F_TESS = synthesize_light_curve('pg-TESS.h5', phi, i)\n", "dF_F_Gaia_B = synthesize_light_curve('pg-Gaia-B.h5', phi, i)\n", "dF_F_Gaia_R = synthesize_light_curve('pg-Gaia-R.h5', phi, i)\n", "dF_F_JWST = synthesize_light_curve('pg-JWST.h5', phi, i)" ] }, { "cell_type": "raw", "id": "93b3427b", "metadata": { "raw_mimetype": "text/restructuredtext" }, "source": [ "Note that the :file:`pg-JWST.h5` :f-schema:`passband` file is based on the JWST NIRCam F150W filter |---| chosen to demonstrate what the light curve will look like in the near infrared, but with the recognition that it's very unlikely JWST will ever observe an SPB star!\n", "\n", "Plot the Light Curves\n", "---------------------\n", "\n", "As a final step, plot the light curve in each filter:" ] }, { "cell_type": "code", "execution_count": 4, "id": "8ac1f847", "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAkEAAAGyCAYAAADj6hCHAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjkuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8hTgPZAAAACXBIWXMAAA9hAAAPYQGoP6dpAADkSklEQVR4nOydd3gU1deA3930hBRCQkkhoffee6hSLIiggApWFGzYP2w/FBQrCopgp6gIYkVApIXee++BkAoJJKQnu/P9cXazCaRtSLIp932efeYyc2f3kGTunnPuKTpN0zQUCoVCoVAoqhh6WwugUCgUCoVCYQuUEqRQKBQKhaJKopQghUKhUCgUVRKlBCkUCoVCoaiSKCVIoVAoFApFlUQpQQqFQqFQKKokSglSKBQKhUJRJbG3tQDlFaPRSGRkJO7u7uh0OluLo1AoFAqFoghomsb169fx8/NDry/Y16OUoHyIjIwkMDDQ1mIoFAqFQqEoBuHh4QQEBBQ4RylB+eDu7g7ID9HDw8PG0igUCoVCoSgKiYmJBAYGZn+PF4RSgvLBvAXm4eGhlCCFQqFQKCoYRQllUYHRCoVCoVAoqiRKCVIoFAqFQlElUUqQQqFQKBSKKolSghQKhUKhUFRJlBKkUCgUCoWiSqKUIIVCoVAoFFUSpQQpFAqFQqGokiglSKFQKBQKRZVEKUEKhUKhUCiqJEoJUigUCoVCUSVRSpBCoVAoFIoqiVKCFAqFQqFQVEmUEqRQKEoUTYP4eDAabS2JQlF1SUmB5GRbS1H+UUqQQqEoEVJS4LvvoH17qFED3N2hQwd48EGYMQP277e1hApF5WffPnjkEfD2hmrVoF49GDoUXnoJfvwR0tNtLWH5QqdpmmZrIcojiYmJeHp6kpCQgIeHh63FUSjKLRcvwhdfiAIUH5//PL0e3n4bpkwBO7uyk0+hqOxkZcGvv8pzuG1bwXM7dIAlS6BBg7KRzRZY8/2tlKB8UEqQQlE4GzeKlZmSIv8ODoannoJx4+DaNTh2TF5bt8LKlTJnwACxSGvVspXUCkXlIS0Nhg+H1avl3w4OMGoUPP00NGoEx4/LM3j0KPz8M8TFgYeHGC0jR9pU9FJDKUElgFKCFIqCyakAde0Kr70m/87Py7NgAUyaJPNr15YFuW/fspVZoahM5FSAXF3h1VdhwgR5vvIiPBzGjBGjBMRg+fhjcHYuM5HLBGu+vytdTNCBAwcYNmwYdevWxcXFBW9vb7p168aPP/5oa9EUikrDpk0WBWjwYNiwAe64o+BtrvHjYfduaN4coqPFI7R4cdnJrFBUJtLS4O67LQrQqlXw1lv5K0AAgYHyrP7f/8m/58yR57AqxwlVOiXo2rVrBAYG8t5777Fy5UoWLlxIcHAwDz74INOnT7e1eApFhWfzZosCNGgQ/PFH0S3J5s1FERo3TrLHHntMXPUKhaLopKfDiBHw77+iAK1cCb17F+1eBwdJVFi5Ery8xCv00kulKm65pspsh3Xt2pXIyEguXrxYpPlqO0yhuJmdO6F/f0m9HTgQ/voLXFysfx+DAYYMgTVroFkzUYzc3EpeXoWismEwwF13wYoV8uytXAkhIcV7r5UrYdgwGS9dKrFElYEqvR2WHz4+Ptjb29taDIWiwpKSAvfff+sKEMi22Y8/gp+fBG5OnCj1hRQKRcHMnm1RgFasKL4CBOLRNW+NPfoonD5dIiJWKCqtEmQ0GsnKyuLy5ct8+eWXrF69mldffTXf+enp6SQmJuZ6KRQKC2+/DWfPgr8/LFtWfAXITM2aEhNkZweLFkm2ikKhyJ+wMHjjDRnPmlUyiQXTpkGvXnD9uniCUlNv/T0rEpVWCZo0aRIODg7UrFmT559/ntmzZ/PEE0/kO3/GjBl4enpmvwIDA8tQWoWifLN/P3zyiYy//FJSbEuC3r3h3Xdl/MwzcPBgybyvQlHZ0DTxmKakQJ8+4rkpCezt4ZdfwNdXnr/nniuZ960oVNqYoIsXLxIbG0tsbCzLly/n66+/5oMPPuClfCLA0tPTSc8RIp+YmEhgYKCKCVJUebKyJAV+716xFJcuLdn3Nxols2zlSmjcGA4fBkfHkv0MhaKi8/PPsh3t5ASHDsmzUpKsXSuJDpomnt577inZ9y9LVJ2gPJg4cSLffvstkZGR+Pr6FjpfBUYrFMLMmfDii5JJcvx4wSm4xSUuDlq2lNT5L76Q+iUKhUK4ckUSCK5cgenT4fXXS+dzXntNMscaNJBn3cGhdD6ntFGB0XnQuXNnsrKyOHfunK1FUSgqDOfPw5tvyvijj0pHAQLpNfa//8n4nXcgKal0PkehqIi8+KIoQC1bwssvl97nvPaaxOqdPQvfflt6n1OeqDJK0IYNG9Dr9dSvX9/WoigUFYLSikHIj0cfhYYNITYWPv20dD9LoagorF0LCxeCTieKSWluFVerZjF63nmnanShr3RK0IQJE3jppZdYunQpGzdu5LfffmP06NEsWrSIF198sUhbYQqFAtatk2q0Tk7w9deyCJcmDg7i6gfxOl2+XLqfp1CUdzTN4vl5+mno0qX0P3PCBKhfX7amZ80q/c+zNZVOCerWrRu7du3iqaeeYsCAATz22GNER0ezaNEiPvzwQ1uLp1BUGGbMkOOECSUfhJkfo0ZBu3aSrmv+fIWiqrJ6NRw4IIVEp04tm890dJS0eYAPPpB4vcpMlQmMthYVGK2oyuzcKRlh9vYSH1C3btl99n//wW23yWJ86hQEBZXdZysU5YnevaVNzYsvSqPTssJohPbtJWX+pZfEM1uRUIHRCoXiljB7YR54oGwVIJBq1H37QkaGJVhaoahqbN0qCpCjI7zwQtl+tl5vWQM+/1y6z1dWlBKkUChycfSotMTQ6aCAIuulhk4H778v44UL4ciRspdBobA1ZiVk/HhpL1PWDB4sCRHp6WW3FWcLlBKkUChyYVZARoyApk1tI0PnzvL5mgYqlE9R1Th4UPqC6fXwyiu2kSGnMbJgAURE2EaO0kYpQQqFIpvz56WfF8CUKbaVxeyFWrJEZYopqhZm5WPUKCkbYSu6dpW4JINBMkQrI0oJUigU2Xz8sSx4AwdChw62laVzZ+jYUWKDVHNVRVXhzBlLaxpbGyJgqd7+9dfyLFY2lBKkUCgAqQtiVjZee822spiZNEmO8+aJcqZQVHY+/FCys4YOhTZtbC0NDB8uleKjo+HPP20tTcmjlCCFQgHA7NkSBNm1qwRElgdGjwZvb7hwQRqsKhSVmehoib+B8uEFAslOmzBBxnPm2FaW0sDe1gIoFArbk5kJ338v45desq469NaLW9kavpWz8Wc5e1Ve8anxDKw/kHFtxjG44WAc7YpX69/FBR55RLbp5syRbvMKRWVl4ULZcurSBXr2LPp9cSlx/HHiD07Hnc5+Bs9fPY+/hz8PtHqAB1o/QKBnYLHlmjAB3n0XNm2Cw4ehVativ1W5QxVLzAdVLFFRlVi+HO68E3x94dKlovUnik2OZfK/k1l8ZHGB82q41GBMyzE82fFJWtRsYbVs585JcKimSfHERo2sfguFotyjadIp/uRJ+OYbeOyxotyj8dPhn3h+9fNcSbmS7zwdOvrW68v4NuMZ22os9nrr/R+jRsGyZfDkkzB3rtW3lynWfH8rJSgflBKkqEqMGAF//CFF2T75pOC5mqax4OACXvzvReJT49Hr9Nzd9G6a+TSjoXdDGng3wNHOkaVHl/LT4Z+ITooGwNHOkfl3zWdMqzFWyzdsmGyHFUU+haIism0b9OgBrq4QFQWFfe2cv3qeiSsmsvrsagCa+TSjf73+2c9gsFcwuyJ2sejQIkLDQrPvu63BbSwdtRQPJ+u+1zZsgH79pIVHZGTh8tkSpQSVAEoJUlQVYmPB3x+yssTV3bJl/nMvJV7ioT8fYt35dQC0rd2Wb+74ho5+HfOcn2XMYu25tczcPpM159YA8G6/d5nScwo6K/bcVq4URcjLS+qVuLoW+VaFokLw+OPSJX78eJg/P/95mqYxe+dsXlv/GimZKTjZOfFWn7d4ufvLONg55HlP2LUwFh1cxIwtM0jNSqV1rdb8M+Yfq7bINE3WhmPHpIr0009b+R8sQ1TbDIVCUWR+/FEUoM6dC1aArqRcof/C/qw7vw5ne2fe7/8+ux7bla8CBGCvt2dww8Gsun8Vz3d9HoDX17/OY38/RqYhs8gyDh4M9erBtWuWOkYKRWUhORl++UXGjzxS8NwPtn7A5NWTSclMoU9QHw5NPMRrvV7LVwECCPYK5s0+b7Lp4U3UcqvFoZhDdP2uK/uj9hdZRp3Okq355ZeiFFUGlBKkUFRhNM0SEP3ww/nPS81M5c7Fd3Iq7hR1Pety6MlDvNrz1QIX3pzY6e2YedtMvhjyBXqdnu8PfM/Qn4eSmJ5YpPv1epg4UcZz5lSeBVihAIm1SUqS2LdevfKf9+OhH5myTtLG3uv3HuvHr6dxjcZF/pyOfh3Z+dhOWvi2IPJ6JL1+6MXK00VPu3zwQahWDY4fl+2xyoBSghSKKsyePdIrzNlZ0tHzwmA0MPb3sWy/tB0vZy9W3b+KRjWKF538VOen+Gv0X7g5uLH23FrG/TGOou7IP/KIBGzv36/6iSkqFzkNkfx2idedW8cjf4mb6MVuLzKl1xT0Ouu/woO8gtjyyBYG1B9AcmYy9yy9h0Mxh4p0r4eHNFUGSyp/RUcpQQpFFca8+N5zj8Tb3IimaTz373P8eeJPHO0c+Xv03zT3bX5Ln3l749tZN24djnaO/HXyLz7f9XmR7qtRQwrIgdoSU1QeTp+W1HO9HsaNy3vOweiD3L3kbjKNmdzX4j4+HHhrDfW8nL1YOXYlQxoOIS0rjfuW3UdSRlKR7jUrQX/8AamptyRGuUApQQpFFSUlBX7+Wcb5xSF8tO0j5uyegw4dP979I72CCvDVW0GXgC58PPBjAF767yX2RO4p0n1jTIllv/yitsQUlQNzEPRtt0FAwM3XLyZcZOjPQ7mecZ0+QX1YMHxBsTxAN+Jg58DCuxfi5+7HiSsneHpl0SKdu3WDunXh+vXKUcBUKUEKRRXljz8gMRGCgyEk5ObrWy9u5dW10sV05m0zGdViVIl+/tOdn2Z40+Fi3S67j4S0hELvuf12iUk4fx527ixRcRSKMsdgsGwr5RWTZ9SMjP1tLJHXI2nh24I/R/+Jk71TiX2+j6sPi+9ZjF6nZ8HBBSw4UPgel15v2TqvDB5ZpQQpFFUU81bYQw/JwpaTDEMGT/zzBADj24xnctfJJf75Op2O7+/8niDPIM5dPceEfyYUGh/k6gp33SXjyrAAK6o2a9ZIyQdvbylWeiPf7fuOreFbcXNw45+x/+Dl7FXiMvQO6s3bIW8DMGnlJI5fPl7oPWaP7D//iCFVkVFKkEJRBbl4EdavlyDMhx66+fon2z7h6OWj+Lj68Mmg0qtOWN2lOr+M/AV7vT1Ljy7lm33fFHqP2QpdulQ1VVVUbMxbYQ88AE43OHhikmJ4Ze0rAEzrO41gr+BSk2NKzyn0r9eflMwU7l12L6mZBQf7tGkDTZtKr8GK3lRVKUEKRRXk99/l2KsXBAXlvnbu6jne2fQOADMHzaSGa41SlaVrQFfe6/ceAJP/ncylxEsFzh80CKpXl2aTGzeWqmgKRamRmiqeFLAEG+fkxf9e5FraNdrVbsczXZ4pVVns9Hb8OOJHarrV5EjsEaZtmlbgfJ0ud3xeRUYpQQpFFWTZMjnec0/u85qmMWnFJNKy0uhXrx8PtM5jdS4FXuz+Ij0Ce5CalcqbG94scK6jI4wcKWO1JaaoqKxeLUUS69aFjjfUG11zdg0/Hf4JHTq+uv2rYvX6spba1Wrz1e1fAfDpjk8JTwgvcL7ZI7tmDVzJv21ZuUcpQQpFFSMqSvoUgfQMy8mSo0tYfXY1TnZOzB0216rWFreCXqfP3nZbcGABB6MPFjjfbIX+9pt03VYoKhq//SbHESNy1wZKzUxl4gqpDPp056fp5N+pzGS6q8ld9A7qTVpWGm9seKPAuY0bQ/v2Um3ebFRVRJQSpFBUMf74Q9LLu3TJnZJ7Le0ak/+dDMDrvV63qhJtSdAloAv3tbgPDY2X17xc4NzevaFOHbh6VSxqhaIikZEBy5fL+EZv7Hub3+Ps1bP4ufsxvd/0MpVLp9Nll65YdHBRoW01zMZIRfbIKiVIoahimC3QGxffqaFTiUmOoalPU17p8UrZCwa81/89HPQOrDm3htVn8tdu7Ozg3ntlXJEXYEXVZN06SEiA2rWhe3fL+XNXz/HB1g8A+HzI51Z3ei8JOvl3YkzLMdnGSEEZm/fdJ8fNm+FSwaF85RalBCkUVYgrVyzBxDmVoKjrUczbMw+A2YNnl2gtEmuoX70+z3SWINCX1ryEwZh/+pfZCv3rL4mtUCgqCmZD5O67c5enmLF5BpnGTAbWH8jdTe+2jXCIMeJo58i68+v498y/+c4LDJTkCk2DJUvKUMASRClBCkUV4q+/JK28bVuoX99y/uNtH5NuSKd7YHcG1B9gM/kAXu/9Ol7OXhyJPcKCg/kXb+vcWf4PKSmWrQWForyTlWVJK89piFxMuJj99z41ZGqZxePlRbBXMM92fhaAl9e8TJYxK9+5FX1LTClBCkUVIq+tsMvJl5m3V7xAb/R6w4rFNxNYATwADAbGAJOAN4DPgZhiyejt4s2bvSVD7M0Nb5KckbebR6ezuOPNKf8KRXln0yaIi5NeeH36WM5/sOUDMo2Z9A3uS/fA7vm/wU2cB94GBgF3A48ALwHvAuuB4vWXea3Xa1R3rs7Ry0eZf2B+vvNGjhRv1t69cOFCsT7KpiglSKGoIly7BmvXyjinEvTpjk9JyUyho19HBjccXIR3OgC8AAQAtwM/AauBX4C5yOL7LFAPWYxjrZb1qU5PUc+rHpHXIwtssGqusvvvvypLTFExMBsid90F9qbM98jrkXy3/zuAbAOgYBKB74DeQH1gKrAG+BP4AfgEMUb6Az1M16xThqq7VOetPm8B8NaGt0jPSs9znq+vJa7JXPeoIqGUIIWiivDPP5CZCc2ayQsgPjWeL3Z9ARTFC5SCeH3agfYpxMXCXhf40x2WOsPfNWC1P4QGwy5fiEtFFuN6wCvA5SLL6mTvxNSQqQDM3jmbDEPeGk7nzlCzpjRz3LSpyG+vUNgEo1GyMyG3IfLR1o9IN6TTI7AHIcEhhbzLCiAIeAzSN8MZYF11WOIAv3vAylqwri5sDYKT9mDYjniJegFrrZJ3UqdJ+Lv7E5UUxeIj+e933XGHHCvitrRSghSKKoLZAjUXGgRRMK5nXKd1rdbc2SSP5kXZhAE9IO0nOKiD5W7i/DmZCinXISsNkuIgLgIiw+DMZbm+0QuupgAfAS2BglNuczK65Wj83P2ISopiyZG8oy71emmqChVzAVZULbZvlzpdHh7Qv7+ci02O5au9UqTwzd5vFmCIGIHpoN0O4dcg1A1+18MuIOYqGDIhLRGuxUDMRdmb2psFf7vBKXswbAUGAi9SVK+Qo50jz3aR2KCZ22fmmylmVoI2bBCDpCKhlCCFogqQlCRbRmCxQBPTE5m1cxZQmBdoHdARYg/AKj0c1SApGexcIWg09P4Lbj8JA7dCrz+g89cQ/CDo9BBxDVYBmz0gIRYIATYXSWZHO8fsTLGZOwpfgJcvlywVhaK8YjZE7rjD0ivsk22fkJqVSie/TgxqMCifO68DIyHrTdiOPEKRyaAZoXpbaPsBDD0Cg/dB39XQ7Udo8x64BkJqMuzJguXu4jUyzgQeBfIPds7J4+0fx83BjcOxh1l3fl2ec5o2hQYNZEt6zZqi/jTKB0oJUiiqAKtWQVqaLFStW8u5L3Z9wbW0azTzacY9ze/J587PQBsIx+NEF0o1gmdz6L4Y7omFHosh4E7waAy+3SFwODR8HLovhGHHIGgMoIPwRPhXD5GJiGt+RZHkntBhAq4OrhyIPkBoWGiecwYOlC+U8+fh2DFrfioKRdmhaZYAfrMhEpcSx5zdc4CCvEBnga6Q+Aes1olTVmcHzV6BYcdhyH5o/gp4tQDvdlBnENS7H1pMgTtOQ8cvwKWOeGx3Adt0YPgBGAWkFSp3dZfqPNLuEUC8QXmh01XcLTGlBCkUVYCccQg6HaRkpmQvaK/1eg29Lq+l4FvIeB62aLKLpQFBY2HQTggeDfZuBX+oRxPo8TMMPQS1+oHBCJt0cDENGA78XKjc3i7ePNz2YUC8QXnh5gb9+sm4oi3AiqrD/v2yQ+XqCrfdJudm75xNcmYybWq14fbGt+dxVxwwCC4eg391kKCBc23ovx7afQCeTQv+UDsnaPwU3HEW2n0Cege4qMFGHWT9CQxDvEwF81yX59ChY9WZVRy7nLelYU5SWLFCynBUFJQSpFBUcgwGS2sJs7X2y5FfiEuNo55XPUa3HJ3HXVsgbSL8B4Qji2enL6H7j+BQzToBvFpCyCqoex8YNdiqgzNZSJB14cVFzAvwP6f+4eSVk3nOqahWqKLqsHKlHAcMEEUow5CRXZritV6v5eEFygTuhaPnYAuQpUHN3uL5qdnbug+3d4FmL0CfFbKNHa3Bej1krEfKWxScWtnAuwHDmw4H4LMdn+U5p2dP8PSEy5dh1y7rxLMlSglSKCo5u3ZBfDx4eUHXrnLOXB36yY5P5tGh+gIY7obNWZKJ6xoAA7ZAo4m5Oz1ag50jdP8JGk6QfYFdwHENeBw4UeCtjWo0yg7azm8BNgdHb98ui7BCUd5YtUqOQ4fK8c8TfxKbHEvtarXzqQ79PIStB3Mv4WYvQ7914FK7+ELUGQj91oKDF1wxwlo9pG4DXi301he6vQDAwoMLuZx880Pm4ABDhsi4IhkjSglSKCo55sV30CCpS7Ivah+7I3fjoHfI3mqykAzanbDzimS0O3hB3zXg0/nWBdHbQad50Ny04O4HTiYD9wKpBd5qXoDnH5zPlZQrN10PDJQq2JpmsbgVivJCfDzs2CFjs6JgNkQea/cYDnYON9zxFVyZA6Z7aPoitPsQbjJYioFvNxiwEZxrwTWjxPplfgb8VeBtPQJ70MmvE+mGdL7c/WWecyqiR1YpQQpFJcesBJkX36/2SDruyOYj8XXzzTHTCIyHo4cswZe9lhUed2ANOh20fR9aT5N/7weuHAaeL/C2XnV70aFOB9Ky0rK/PG6kIi7AiqrBmjVSI6h5c6hbF05eOcmGsA3odXoe7/D4DbM3QvJTsAl5JP3vlOyvkqR6axi4BVz8xdu7G9DGA/mXfNbpdNnGyJzdc0jLujmoesgQaW585AiEhZWsyKWFUoIUikpMbCzs2SPjwYMlLf6nwz8B8ESHJ26Y/RFc+A0Omf7Z6Uuo3b90BGvxOgSOlEV+K5DxFZB/B8acC/AXu77Is3iiWQlavRrS8y5uq1DYhBsNka/3fg3A0EZDqetZN8fMaMgcAaEGSdyq3la2kfV2JS+Ue0Po8YsYO2HAuQRgNBKLlDf3NLuHQI9ALqdc5ufDNyc2VK8usUFQcYwRpQQpFJUYc0B0u3ZQuzb8dOgnkjOTaerTlN5BOYMrz0Dcmznc7y9I/E5podNBl2+hWn1IRj5Xe0zkyIdRzUdRu1ptYpJjWHn65j2vDh2gTh2pibRxY6lJrlBYhdFoqdE1ZAikZqYy/+B8AJ7s8OQNkyfD1nhIQGJ/+iy3PhHBGmr2hNbTZbwHuLYDeD3f6Q52DkzqNAkgu83HjVQ0j6xSghSKSkxOC1TTNObumQvI4mvJRtHA8CRsywQD4Hc7tP2w9IVz9IQeSyTz7BJwKgm4j/yKuDnYOfBg6wcB+OHADzddz1k9+u+/S0dkhcJaDhyAmBgp5dCzJyw7toz41Hjqeta9oVffGjizBCKR1PbeyyUpobRp/grUGSzP/hYg8yMkLTRvxrUZh16nZ1v4tjyzNc1KUGgoJCaWhsAlS6VTgtavX88jjzxC06ZNcXNzw9/fn7vuuou9e/faWjSFokzJmRo/ZAjsuLSDw7GHcbZ3ZlybcTlmLoGj66RciEtN6L6odNzveVGjI7T7WMb7gbh9wDf5TjcHcq84tYKYpJu71Kvq0YryhjlQv39/KeppTot/vP3j2GU/Z2mQMkF6E4PU9KnRsWwE1Omh20Jw8csRH/QM+aXN+7n7ZStveXWXb9xYXpmZlvWnPFPplKC5c+cSFhbGc889x8qVK5k1axaxsbF07dqV9evX21o8haLMuDE13tyfaHTL0VR3qW6adQ2uPQ3m+mcdvgRHr7IVtPEzEHC3JT4o6w0gPs+pzXyb0cW/CwbNkB3blBPzF83Fi3Ci4Mx7haJMyOmNPRxzmG3h27DT2fFou0ctk7R3YXeYOEF9Oks5irLE2dcUH6SX+KALp4A5+U43GyMLDy3EYLy5MqLZI2veBizPVDolaM6cOaxfv56JEyfSp08fRo4cyZo1a6hRowbvvfeercVTKMqMnKnxiZnxLDkqgce5AqK112BXnCkL5XYIHFH2gup00PV7yVRJAk7EA2/nO928AP9w4Ieb+om5ukKvXjKuaD2MFJWPG1PjzYbI8KbDqeNexzTrBFyaARGIB7bzd6KMlDU1e0HL/8l4P5A5FamTcTN3NL4DbxdvIq9HsubczQ/aIFMLtDVryr9HttIpQTVr1rzpXLVq1WjevDnh4eE2kEihsA05LdCFBxeSlpVGm1pt6OLfxTRjF5yZC1eQirKd5ha/GOKt4ugF7T6S8VEg5Qss7qncjG45Gmd7Z47EHmFv1M3b3AMHylEpQQpbkzM13qdOMosOLQKkSKmgQcZjsMfkTWn2f1Jh3VY0f0WSFVKBY4nAm3lOc7J34v5W9wN5x+f16gWOjhAeDqdOlaK8JUClU4LyIiEhgX379tGiRYt856Snp5OYmJjrpVBUVG5MjV9wcAEgDUklINoAKQ+LxQfQ5sOyCcIsiKDR4NtTAjT3G5HaQTebkZ7OnoxoJh6rH/bfvACblaDQUIlLUChsRU5D5M8Tf5KYnkj96vXpV8/U7I4f4cBWUTrcg6HlGzaS1ISdM7Q39eg7DiR9hSVQKTdmj+yfJ/4kPjX39rWrqyVVvrwbI1VCCXrqqadITk7m9dfzT/2bMWMGnp6e2a/AwMAylFChKFlypsbH2x3jQPQB7PX23NfiPtOMxbDnmMQg1Ghf9jEIeaHTQYfZgE5qtsX+R37d5s0L8M9Hfr6paFubNuDrK6ny5q0IhaKsuTE1/ucjUlfngVYPmBoWp8PllyxVITr/IEqIrfG/E2oPlC3yfQCTycsYaVenHW1qtSHDkMHiwzf3AKwoHtlKrwS9+eab/PTTT3z66ad06NAh33lTpkwhISEh+6W2zhQVmZwWqHmBGtxwMDVcawBZEP1/kpau00OX+WWXDVYY3u2gwWMy3ovUTckjS6VvcF8CPQK5lnaNv07kLvev10uANJT/BVhRecmZGt+k/WVWnxHLZGyrsTJB+w72xsq4wUNQK8QWYt6MTgcdPpMiipeA6I3Ab3lOfajtQ0DeW2JmJWjDhvLtka3UStDbb7/N9OnTeffdd3n66acLnOvk5ISHh0eul0JREcmZGj94sJZtgZr38NF+gkMRMm70GHi1soGUBdBmOjh4wlXg/Fng5j5Fdno7xrcZDxS8ACslSGErzIZI//6w/MwyDJqBDnU60MSnCZAOl96SJEh7R9mOLk94NodGT8l4L2B8ibyMkftb3Y+93p69UXs5HHM417V27aBGDbh+vXx3la+0StDbb7/N1KlTmTp1Kq+99pqtxVEoyoy9eyUrxdMT9HV3ce7qOdwc3Lij8R1AFkRNkWBoO3toMdXG0uaBc01oZcpSOQhkvA/c3AfDbIX+d/Y/LiVeynXNrATt2gXXrpWapApFvuSs0WU2RCxeoG/gUJyMmzwvKerljdZTwamGVK8+fQG4uU2Gr5uvaV252RipKB7ZSqkETZs2jalTp/LGG2/wv//9z9biKBRlirkcVt++sPSYLFzDmw7HzdENtEVwKEomNJ4ELnXyeRcb0/hp8Ggi/ZOOxQALbprSwLsBvYN6o6Gx6OCiXNcCA6FJE4nLUOXBFGVNcrIlHq1Z1wtsubgFHTpTTF4aXHhLlAsHF2j2qi1FzR/H6paWGoeBzPeQrIXcmOPzfjz0I1nG3NXezcbIf/kXoLY5lU4J+uSTT3jrrbcYPHgww4YNY8eOHbleCkVlZ906OYb0y8quDSQWaBZces3igm9m40yUgtA7WDpnnwYy3iOvdhrjWkvl66XHlt50TW2JKWzFli0SB1O3LmxP+gWAkOAQ/D38wfgVHLkqE5u+LMpGeaXB4+DRSHbCTp8G/rhpypBGQ/Bx9eFyymVCw0JzXcvpkU1IKHVpi0WlU4KWm7q2/fvvv3Tr1u2ml0JRmUlPlwUYwK3lBmKSY6jhUoOB9QeCtgAORcvFJs+VTxd8TvzvkNiETODUBeBmRWd40+HY6ew4EH2A03Gnc11TSpDCVpi9j/36weJcW2FpEDZV2lM4uUHTF20lYtHQ20FzUzjJCSBrOjdmitnr7bm76d0A/Hr011zXgoKgUSOJU9ywoQzkLQaVTgkKDQ1F07R8XwpFZWbHDkhLk47xWxKkrcS9Le7FwQ64MMXkgneGZlNsKmeR0OmhhclbdQLInI7k7Vqo4VqD/vUl8ODXY7kX4JAQsLODs2fh/PnSF1ehMGNWghr3PMKhmEM46B24p9k9YJwHh6/JxWavg0MFSMAJvh/cAmVr+uxB4OaGYPe2uBeA30/8nu+WWHk1RiqdEqRQVGXMi2/v/qn8fvx3wGSBGn+AI6YS+E1fKt8u+JzUvVcq2GYAZ44D/9w05d7msgDfqAR5eEjPNCi/C7Ci8nH1KuzbJ+PIGuIFGtpoKNVdXODcVEgGnD2g8XM2k9Eq9A7Q3GQ0HQcM02+aEhIcgo+rD1dSruS7JVZen0GlBCkUlQhzPJB3lxVcz7hOXc+6dA/slsMF7ypxCBUFvR20MLnjjwOGadzoji9oSyxnDyOFoizYuFEC8hs30VgelmMrzLAAjpgCY1r8D+xdbSilldR/GFxqQQoQthXYnOtyQVtiffuKR/b0abhwoYzktQKlBCkUlYSkJNi5U8bnqsniO6blGPSsg2OmjLBmr1QMF3xOgh8EV39xx5/bA+QOLihoS8xsha5bJ3EJCkVpY/bGtrhtOxcSLlDNsRq3Nx4GYe+KEuHiCQ0n2VRGq7FzthhPRwHjuzdNyW9LzNMTupjaFZZHY0QpQQpFJWHLFsjKgqDGCYRGSLuJsa3GQvRb4gWyd4BGz9tWyOJg5yiNJUF6qhqn3TQlvy2xTp1kEc65RaFQlCZmJSi1oRgidze9G1f7PXDS1IWgyfPloz2GtTR8QhodJwEXV2PqqZFNRd0SU0qQQlFJMG+FBd/2DxmGDJr5NKNVTVc4YXIPNRhT8bxAZho8KtlsyUBYKHAo1+X8tsTs7cUdD+VzAVZULqKj4ehRQGdkf5rE5I1uORpi34RrSIHShs/aUsTi41ANmrwg46OA9mmuy/Z6e0Y0lcbGN26J5fTIGnPnNtgcpQQpFJUEswWaEiS1PEY0G4EucTqYdsJo/JZtBCsJ7F0koBvEG6TlbqVR0JaYuWptaGgpy6io8pjTwBuF7CYmJQp3R3f612sEJ0wxNPXvqThJCXnR5BlwcJMs04hfkNLzFka1GAXcvCXWubP0UIuLgyNHylDeIqCUIIWiEhAfD/v3A/apHE2T1tUjmg2CU6ZS9wHdwL2B7QQsCRpNlAU4EYiajwwsmLfElh7NXU8oJESOW7eW70aOioqP2RDx6iqGyNBGQ3FK/RhMrfpo8rZtBCspHL0sPcVOZgG5W2XktyXm4AA9e8p448YykbTIKCVIoagEbNwImgaBfdaSkpVMoEcg7arvgHOmb/0mN6e1Vjgc3KH+ozI+nQ78mOuyeUvsYMxBTsWdyj7fvLk0ckxJgT17ylBeRZVDlCCNSE9RgkY2GwonTS1f/DpIK5iKTqNJ0mk+BkiYTc7aXQVtifXpI8fy5pFVSpBCUQkwxwNV6ySL791N70J37mNp9ePlDzX72k64kqSRKasmAkj6jJzp8rm2xHIswHo99O4t4/JmhSoqD2FhcO4c6GudICLtFI52jgwLjoJzpua/lcEQAXALAr+hMj59CcjdGCy/LTGzR3bTpvIVF6SUIIWiErB+PaDP4pLb3wA83MYPTpmLI74hlltlwKMJ1A6R8ZnTwJZcl/PLEjMvwOXNClVUHsxbYf79xRDpX68fLhc/k5Z3nrWg9m02k63EaWwK7j4PZH6e61J+W2IdO4KrK1y5AseOlZ2ohaGUIIWighMVBcePA3W3ct0Qh7eLN60zlkpNEmcXCHrY1iKWLI1MlXbPAobcC3DOLbFzV89lnze74rdulTICCkVJY1aCDI3+BOCpDs3glLlX35TKY4gA1B4A1epKX78Lq4CL2ZdyFk40V60HiQvq3l3G5ckjq5QghaKCY158a/YWC/SRtiHoTx6Qk40eBzsn2whWWvjfDq61IB24+DsSnCDUcK1Br6BeACw/uTz7fKtWUL26FJRU9YIUJY2mmbakPS4RqduNDh0DnTdKSQcnZwieYGsRSxadHhpPlvEpDbSvcl0e3nQ4AMtPLc/Vs9PskVVKkEKhKDHMwZhpwX8CMKlhIsQBeh00fM2GkpUSento+IyMTxuA73JdvqPxHQD8fepvyy054oLUlpiipDlxQmoE2bf8E4AhDTvheHa/XGx4v5R4qGzUf0gKmV4DrnyJNPgT+tXrh6uDK5cSL3Eg+kD2ebNH1pzIUR5QSpBCUcHZtAmofYBE/QVc7Z0JvrJdLgR2k34/xSE1Ff74A8aPl2qDXbtC27bQpIm8xo2DhQshIqLQtyoVGjwmfcWuAPGfIxHgglkJ2nRhEwlpCdnncy7ACkVJsmmTHKt1+hOAN9p7QIzpW77hLdTnOnYM3n4bbrtNcsw7dIAWLaBBA3GrTJsG27bZpvaDY3UIGivjU9eAP7IvOds7M7C+VEhcfsrike3UCVxcIDZWFMfygFKCFIoKTFQUnDkDNJMF6LVubdFdSJaLDd6w7s0yMmDZMhg9Gnx9YcQIUXRCQ6Up2cGDcOqUvBYtEgUpIACaNoU334SEhEI/osRwqQWB98j4dDSwMvtSoxqNaOrTlCxjFv+e+Tf7vNkVv3mzigtSlCybNwMu8SR4hQLQKc3kBarTDNzqWvdmx4/DW29JbYcWLWDqVPjvPwlo27dPFKNz50Sbf+st6NFDakDcdZfMKUsamzyy4UDqZ7ku3dnkTiC3EuTkBN26ybi8eGSVEqRQVGA2mwrROrf9E4DHa16RYEU3D6hlRTbKxo3Qpg2MGgVLlkByMgQFwYsvwi+/wF9/werVYvKuWgX/939i1ul0cPIkTJ8ODRvCnDllZ5U2Mi3AYUDGvFyX7mwsC3DOLbHWraWP2PXrcOBA2YioqBps2gQ0/gdNZ+C+Zg2wD4uTCw1eLvqbJCTA00+L4jNtmihDjo5w++0wbx789husWCH735s3w9y5cM894O0tf9R//y3eopEjTZZRGeDdHmq0k1JBZ3cAls8d1mgYOnTsidxD5PXI7PPlLi5IU+RJQkKCBmgJCQm2FkWhyJenn9Y0qp/RmIrm/b5eM67RadpPaNqhJ4r2BrGxmjZ+vKbJFr2m+fpq2ssva9quXZpmNBZ+f3y8pv38s6Y1aWJ5j8aNNe3PP4t2/61gNGraisby/z2u1zTtcvalzRc2a0xF83rfS8vIysg+f8cdIuLHH5euaIqqQ1iY/E3pRt+tMRXtxIFG8je5zFHTstIKfwOjUdMWL9a02rUtz9DQoZq2aJGmXbtW+P0Gg6bt3atpjz2maXq93O/goGmTJ2valSu3/h8sjHML5f/7B5pmeCPXpa7fdtWYivbVnq+yz23cKCLWrl16S4Q139/KE6RQVGA2bSJ7K+zj7kHoYk1xCPWLEBC9YIFsZS1YIB6dJ58Ur86HH1q8PIVRvTqMGQOHD4sXyMdHtsuGD4dHHoH09GL/3wpFp4NGk2V8zggsyb7ULaAbNVxqcC3tGlvDLVsE5bVqraLismkT4JCCrtG/ONpBo4TzciF4aOGZmRcuwODB8gxFR0PjxpJmtmIFPPCAuC4LQ6+H9u3hm2/ExTl4sHhjP/tM0iJLu0x63VHg6CYlOWK+IWcF6ewkhZMWj2znzrItFh0tS4WtUUqQQlFBuXpVdA8ay577KGeTC75O04LjEIxGeOkleOghaTrWurUEV86dK0pNcXBwgEmTxA3/f/8nC/P8+dI++sqVQm8vNnXvkwDpa8DVudmn7fR2DGs8DMidKm9WgjZvBoMBheKW2bwZqLceo10qk9t4o48wBZw1mFrwjfv3Q5cuEu/j5CQB0IcOQb9+xRemVSvZrl69WhIYoqIkLXLZsuK/Z2HYOVsCpM/HkLOAqVkJWnd+HSmZKQA4O1vigsrDlphSghSKCsrWraA5XYO6W2npA9UumRqKNngp/5syMiSz65NP5N/vvAN790r2V0ng6QkzZogl6+Eh3xCdO5deiVgnb/A3lfA/dxQ4nn0pZ6q8ZsrHbdtWxEpIkO8bheJW2bQJaCSB+c/XNkonlxp+4NUm/5vWrhWNPCZGjJDDhyXI2amEanoNGgS7dsGQIZLpOWqUxO2VVl56PVNPv3Ag89vs0y1rtiTYK5i0rDTWnlubfb48eWSVEqRQVFA2bwbqrwG9gZmtqkEq4OQI/g/mfcP163DHHfDTT2BvL5lfb74p45Jm8GDYvh3q14fz58X0W7265D8HoL6pEF0YYJyfffq2BrfhaOfImfgznIw7Cch/1dzNujwswIqKTWwsnDypQaOV+LlDrZhrcqHBxPxv+uUXGDpUnseQENGiGjUqeeE8PCRYevJk+febb8oWW1payX9Wjc7gUVcqVVxciixGoNPp8twSyxkcbet6QUoJsgGaJga5QnErmC1Qez2EGGXRIXiwFDC7kdhYqffz33/SwGf5cngwH2WppGjeXFLre/WCxERRwMydXkuSOreBs6dUkI78DnNMgruTOyHBIUDuLbFyl52iqLBs3gz4HgevC/yvuR26RMBOD0HP5n3DrFkS/5OZKd6Zf/8tWtxPcbG3h08/lewye3v4+We4//6S3wvW6aDekzI+lw78lX3JrAT9c+ofjJo8m126SOJbZCScPVuyoliLUoLKmFmzwN8fPvjA1pIoKjIpKbB7jxEarWRMA3CIMi1qDd6+eXJysrjF9+6VwOUNG8RTUxb4+Ijr/557ZOG/+26JhShJ9A4QPF7G5+OA0OxL5gU4Z60Ssyu+vHWzVlQ8Nm8GGq0A4H4308mgbuDgcfPkb76xeGWefhoWLy657a/CeOIJ2aJ2dITff5fPL2kXTD2TUXUZSLKUrOgT3Ad3R3dikmPYEylB2i4ulh14W3tklRJUxtjbS6zali2Fz1Uo8mPHDjDU3AfVYnnVX2eKQ6gDXm1zTzQYYOxYKbLm6yt/eJ07l62wjo7w44/igrl+XRSyc+cKvc0qzDEJEUD619mnzUrQ1vCtxKVI4Hj79uDmJoHlR4+WrBiKqoXZG9s3ANzMhkj9PCpEr10LE01bZK+9BrNng51dmckJSJzQTz+J12bePIkRKklcA6B2Dxmf2whI81hHO0cGNxSjK+eWmNkYsfV3oVKCypgepr+R7dtVdoqi+IgFuhJ3R2iWZLLoGjx688SXXpK4ACcnKXjYpEmZypmNszP8+acUZIyJkTYAsbEl9/7VW0P1RrITduEPpHMlBHkF0bpWa4yakZWnJXjV3t6SnVLWBXYVlYeEBDhwIgHqbuHVICALcHcH34G5Jx47JgUMDQbZipo+3XYd5UeOhM8/l/Fbb4l3qiSpZ1L0zgPaj9mn8/LIPvywPH9ff41NUUpQGdOqlTwn16+b0psVimJgtkCfbgz6BORJrvt87klffim1QkBqAZm/+W2Fp6ek7wYFSSr9sGHyIJQU5gX4XAY5+xiZF+AVp1dknzMbI7a2QhUVl23bQAtei94+ixDzyeA7cys4MTHyd56QIBH5331nOwXIzFNPweuvy/jJJ8U4KikC7wZ7Z7FBLlu2xIY2Gopep+dQzCEuJV4CoF496N5dHMW2RClBZYydnbJCFbdGRgZsO3AZ/HfxiLfppF9TcPS2TFq5Ep4xtZV49124774ylzNP6tSRLLEaNaSI2+OPl1xsQvD9oNNDPJDwZfbpIQ2HALDm3BoMRnG/mpUg9QwqiovZG3tXEDhdNp0MmmKZkJoq/bzCwqSlzJ9/ll0MUGFMmybFTI1G2S4vqW6m9q4QNFLG584CBwGo4VqDzv6yDb/6TClliRYTpQTZAHOKrrJCFcVh3z5IC1hNDVeNBuaepUFPWCacPi1Kj9EoBRGnTMnrbWxHkyayRWdnJ33Kvv++ZN7XuSb49ZXx+e1AFABdArrg6eRJfGp8dmBm165SzzEsDCIiSubjFVWLjZskNf4FfyQmz9sXPFpYJkyYINmR3t4SlFyjhq1EvRmdDr76SgozpqRI0+SSSp2vZypZcRHIsjzbgxtIXNC/Z//N4ybboZQgG6CsUMWtYLZAX24MuhTAXg/+poUnM1NqgSQlSWr6V1/Z3v2eF927i4cKxGNVUsUUc8YkGH8FwF5vz8AGEqdh7irv7i7hSaCeQ4X1pKbCrosHcPCMpoupQDTBoywTfv5ZkgHs7CQbq3Fjm8hZIPb2sGiRJEwcPAgvW9HstSB8e0K1WhIjFf4ToiGSHRy95uwasoxZ+d9fxiglyAZ06SLPRng4XLxoa2kUFY2NmwzQ8F/uN2fhBrQVNzSIm3vXLvDykkwQW2+4F8TLL0tbjdRU8Vylpt76e/rfDo6uUqvtsiXo87YGtwGw+qzFFa88sorismsXZAWv5IEG4BBvOmneCrtwQVrIALzxhiUNqjzi5yfxggBffCFbdreKTgfBpiSNC3HALgA6+nXE28WbhPQEdl7aeeufU0IoJcgGuLlBu3YyVlaowhqMRth4dieBNa/if9V0Mvg5OW7davGuzJsHgYE2kbHI6PViidaqBUeOwAsv3Pp72jlBgARCc+EIkjNvUYJ2RuwkPlW+tZRHVlFczN7YZ2qaTtQKBJcAyQAbN04Cobt0ESWovDNkiGSRgsQJlYRlHmSqGRQFpIuSZae3Y2B98cjmNEZsjVKCbISyQhXF4ehRSKqzktcagi4dcHKA2mOkIvODD4qW9MAD5ScQujBq1RJFCERxK4lGj0GPyDEcMEpn+UDPQFr4tsCoGbN7GJmVoAMHSjZJTVH5Wbc9Dpd622ltDqMJNn3pf/yxpG66ucl2WGm0pCkN3n0XOnWS4lljx0LWLW5XeTYFr2DZCbv0C+Yq7uYtMfO2dHlAKUE2QlmhiuKwbRvQaCX3mHa/qNtVKiY/+6z06AoKErd2RWLgQHj1VRk//jhER9/a+9XqB05u0kYj5rvs0zcuwAEB8uMyGiV+VaEoCgYD7Lz8H0820rBLRL5FA1+SjIU335RJs2dLRlhFwdFRepp5eMiXUkm0NAh6WI4XrgI7AItHdk/kHi4nX87nxrJFKUE2wqwEHToknlOFoiis2RFF02b78TXHIQS/It6TBQtke+nHH0u3F1FpMW2alHK+dg2ef77Q6QWit4fAETK+cAxxCeWOCzJ3lVceWYW1HD0KqQErmWBO9vJrDAZnKYSYmQkjRkglwIpG/fowZ46Mp02TWl63Qt375RgDpM0HoI57HdrUaoOGxppza27t/UsIpQTZiDp15G9O06QFgkJRFDZF/MebDZHMCzdncOgpfYAA/u//LN/qFQ0HB6leq9eLRfrvLbrLc26JGRYD0CuoFy72LkRej+RI7BFAeWQV1rNlqxGP5qtokmQ6EfwYzJghtXbq1JESyOUxI7Mo3H+/tNdIT5dCirdSw8u9AXg3lC2x8KXcuCVWXuKClBJkQ5QVqrCGmBi47L6GYeaEr6B+Uvo+JkZq77yVR8+iikT79vCcKch74kSpX1JcfHuBiwdkAtGyJeZs70zfelJHyLwlZn4Gt2+/9TAIRdVg5d7DvNg2Dl0KaPZ6SBpk2T764ovyVQ/IWnQ6qTTv7Azr1kmG6a0QZM4SSwC2ATk8smdWZ3eVtyVKCbIhygpVWMPWbUY6dVuBZ5zpRMI9Fvf1nDnlpxrtrfDOO5LVFhYGb79d/PfR20GgqW7LhVNAGHBzwbYWLWT3MDlZtqYVisLYHvsf40w7zrqAlvDcS1LGfcgQuPtu2wpXEjRoYDGonn8e4uIKnl8QdcfKMRZIEWOkR90euDm4EZMcw8Hog7cmawmglCAbYrZCd+yQrWSFoiD+3nmYV5pdAyNo7m7wylcS1TtmDPTvb2vxSoZq1SyK3SefSBG34mLeErtE9paY2RW/5eIWkjKS0OulbiMoj6yicGJiwFhnBUHmbMJDnaVDvJOTNCatqNtgN/LSS9CyJVy5cmtFFN3qgk9TGYf/BhhwtHOkf31Zr8pDlphSgmxI06ZQvbrUiDtwwNbSKMo7oeH/McCUcavb3UAqtrm7i7JQmbjjDrjnHknDeeIJORYHn67gWl3ipyKlfH9D74bU86pHhiGD0LBQQHlkFUVnw5ZUJvfbgi4NtEwdvGNqyjtlinhQKgsODpb27j/8AKGhxX+voMfleOE6IJZGeWqhoZQgG6LXqwVYUTTS08HD73e8rgKJwOfn5cL06RKMWdmYPVvSdXfulIDp4qDTQ9AYGV84A5xDp9PdlCqfMzavpHq5Kionv+3ewtjqJqV8lTdERYnyYy7xUJno1k2Co0GOxd2uCBwtxytA8rcA3NZQ4oK2hW8jMT3xFgW9NSqlEnT9+nVeeeUVBg0ahK+vLzqdjqlTp9paLAD+OfUPI5eOZMEBqaJpVoKUK15RENt3p/Jyz52yFbbMDq5dh7ZtLeX5Kxt+fpbq12++Wfw6EnVNqcoRQJYEed6oBHXqJDXtIiOl44FCkR/74/+hQTJwAXTLTCXbv/hCAokrIzNmQM2acPKk9CEsDq5+ULOVjC/+ARioX70+jbwbkWXMYv359SUmbnGolEpQXFwcX3/9Nenp6QwfPtzW4uTiUMwhfjv+G3+d/AuwWKFbtyorVJE/i7duYYizBmdBt85kic6dW3Eq0haHJ56QrLcrV4pfvM27A1TzAQMQIYZH3+C+2OvtOXv1LGfjz+LqCh06yHTlkVXkR3o6jOy0FH0KaIsAgxFGjoTBg20tWunh5WVJUJg6tfjGSJCpwfOFZMxZYuWlenSlVIKCgoK4evUqGzduZMaMGbYWJxcD6g8AYP359RiMBjp2lGKd0dFw7pyNhVOUW84kLsHnKrDYdGL8eOja1ZYilT4ODvDhhzL+9NPi9TTS6aCuaUvs4lkgAncnd7oFdANg3fl1gPLIKgpn7fYYHqwbDftBdxwJhq5s8Xh58dhj0KyZZIm9917x3iPwXnkW44Gk+QCMbTWW9/u/zzOdnykxUYtDpVSCdDodunIapd+hTgc8nTxJSE9gb9RenJ2lPApIrRKF4kY0Dca0+gP2AscR1/u0abYWq2y44w4ICYG0NHj99eK9R93xcowEsn4FLMaIuY+YOUNMPYOK/Phl1yqaJgFLTCeeew7q1rWlSGWDvT189JGMZ82S8hXW4lwTfFvI+NIfgEbXgK682vNVWtRsUVKSFotKqQSVZ+z0dvSr1w+wLMBmg14twIq82H0shts94uEX04lnny3/HeJLCp1OmlKCtATZu9f696jeHly9ZEss+gfAogStO78Oo2bMfgYPH1bNVBV5U9dzLnarkZIL1b2kQntVYehQ6NdP9gRfe6147xH4kBzDrwK2rw9kRilBJtLT00lMTMz1Ki1utEK7iWdetc9Q5MlvO5dSazmy+HpVq1qLL0jAzgMPyPill6wPntPpIHC4jMMPA1fp5NcJd0d34lPjORB9AH9/0SuNRtizpySFV1QGjEaNUb57YZnpxOtvSH2TqoJOJ1t/Oh0sXizlOawl4F45XgZSfyxR8W4FpQSZmDFjBp6entmvwFK0tM1K0NbwraRkpmQrQQcPSuVahSInrR2+QGd2wb/+v6q1+Jp5913ZBgwNheXLrb8/0JwlpoHxLxzsHAgJDgFuNkaUR1ZxI6v3HaTVnwa4CpqfBzz1lK1FKnvatoVx42RcHGPELRBq1JdxxC8Fzy1DlBJkYsqUKSQkJGS/wsPDS+2zGnk3ItAjkAxDBlsvbiUwEPz9pSacskIVOdE0jWHrT8NVyKppb2mWWtWoW9fSXf6VV6xv9OXTA5xcIQOIlfL9ObfEQClBivyJPvQOdr/JWDf9/cqbEl8Y06eDiwts3gx//WX9/QGmzvLhEUD5yARSSpAJJycnPDw8cr1KC51Ol++WmFqAFTnZv3cFnkvF4kp5flzVXXxBtgF9fKRmibWNHfV2ECAF2gjfAaRmP4ObL2wmLSst17a0KlehyMnQv1dACmTW08P4J2wtju0ICIAXXpDxm2/K/rE1BJq2taOBjMUFTi0rlBJkI7KVoPNKCVLkT/UZk6RbdRB4vDTX1uLYFg8PSx+jadOsr2Ab+JgcL2WBtppmPs2oU60OqVmpbA/fTrt2kvV85QqcOVOyoisqLulnT1JzRQYAYY+ESKn/qsyLL0rX4SNH4LffrLvXozF41gENiFhYKuJZS6X9ba5atYply5ax3BQ/cOzYMZYtW8ayZctISUmxsXTQv540kNsftZ8rKVdyKUHKClUAEBND0ArZlg2/vwbYO9pYoHLAU0+Bry+cPQuLFll3b63+4OAIqcCVb27yyDo6WoomKmNEYSZhyv3oskBrDtUetPJvrjJSvTpMnizjqVOt7+0XMEqOl04BMSUoWPGotErQxIkTGTVqFI88Ip2kf/31V0aNGsWoUaOIjY21sXRQq1otWtZsiYbGhvMbaN9eiiZevqyKJioEwwfT0acDDeBAy0dsLU75wM3N0qfJWm+QnRP49ZLxpfVAVrYxojyyijy5eBGfP6Qsw9X7HakT5GdjgcoJkydLNeljx+DXX62715wqHwlkLStoZplQaZWgsLAwNE3L8xUcHGxr8QAYUM9ihTo5qaKJihzExqKbOw8A473g3/htGwtUjpg4EWrVkqJt8+dbd2+gqXx/eBpom+hfX5SgPZF7uJp6VSlBity8/z76LKAFbAvoZGtpyg9eXpbYoLffts4bVL0tuHlJ3a6o70teNiuptEpQ+ceo4oIU+fPRR+jTsqA+nOvgSps2LraWqPzg6mqplTR9OmRkFP1ev2FgZwdJwLXvCPAIoKlPU4yakdCw0OxnUBVNVBAejvbtNzK+Gy6kvmlbecobzz0nW2MnTsAvVqS863QQeJeMww8Aqot8FWM+0Br4nN5BvbHX23Pu6jnOXT2nKkcrhNhY+PJLGY+AdVE9KnWf1GLxxBNQp470E/veCmvS3g3qdJRx+D+Alssj6+dnKZq4e3fJi62oQLz/PrrMLGgGiR11dGp3m60lKl94eEi9IIB33rGubEXAo3KMMILhn5KXzQqUElTWxOyC7Yfh4re4O7nTNUA0n3Xn1mVboYcOqaKJVZpPPoGUFKgPWns4dWWqrSUqf7i4WLxB774r5fyLSoA5SywROKA8soqbuXQJvv1WxiNgS4I/bdvaVKLyyTPPQI0acOqUVJIuKj7dwdkNMoHYb0tNvKKglKCyJsYezgMXjgIJFiv0/FpVNFEh+dlz5sj4boitpqdnm262lam8MmEC+PnJF9Z33xX9Pv8R4pK/BlxfQEhwCHqdnlNxp7iYcFG1sVHABx9ARgZaU6A5/HFoLI4qOfNm3N0tZSus8Qbp7SDgdhmH1ysd2YqIUoLKmoCH5BilgeFvS9Xac9LIUVmhVZxPPoHkZIwNgHawOjqYrl11tpaqfOLsDFOmyPjDD4ueKebkDTWbyDjidzydPens3xnI7ZFVRROrKBER8I3EAulGQIYr2POCjYUqxzz1lHiDzpyBZVZke9V7Bpq/Cg1t24JEKUFlTfV24FINsoCY7+ns35lqjtWIS43jYPRBpQRVZeLi4PPPAdDfBejg240PUaeObcUq1zzyiNQNunABliwpfL4Z/zFyjAgHInN5ZFXRxCrOhx9CejqG1g7QHA5oLvTrUsvWUpVfqlWDZ5+V8fvvF91y8O0Bbd8H7/alJ1sRsFoJio+Pzy42aDAYWLhwIT/99BOaMpmKhk4HAbLgcmkbDnY6egf1BmBD2AZVNLEq88UX4gVq5gXt4bo7ePGwraUq37i6Wgq3vf9+0cv4+5t6GMUCGUvpV68fABvOb8DBQVNFE6sqly9ne4H0d2aCDuYf6JKdtKLIh6eflhpeBw/Cv//aWhqrsFoJGjRoEKdOnQKk6ehHH33EzJkzed7c3FBROP6mwMyIDNC20je4LyBKkCqaWEVJTs72AhnuTgcdbEt3pV/HABsLVgGYNEliE44ehRUrinaPewPwrCnl+yMX0i2wG052TkQlRXEq7pTyyFZVvvgCUlOhdQC6pmB0hN/WPU2AegwLxttbMjZBjJEKhNVK0OnTp2lrCpP/8ccfWbVqFevWrWOJNa7oqk6tAWBvL+X7r36bbYVuDNuInUOWKppYFfn+e9kOq1cX++apAMzZ0UtZoEXBy0sKKALMmFF0F6q/KTDz0iGc7TW6B3YHYP359UoJqookJ4sSBGTeB+jgrAt0qTMAnQrLK5wXXgAHB9i0CbZts7U0RcZqJcjR0ZHk5GR2795NrVq1CAgIwN3dnWSV01107JygjsnffmkVbWq1wcvZi+sZ19kXtU8twFWNzEwJiAZ4pBE6nQRjrlrzCO3a2Va0CsPkyRLIs307bN5ctHvMqfJRBjD8m8sjm7NoYlJSyYurKId8+y3Ex0PDhtg1iQBgWUQdenX2tLFgFQR/fxg3TsYVyBtktRI0ZswYQkJCGDduHOPHjwdg//795aYVRYXB/yE5RsRhpz9Hn6A+gFihZut/507biKYoY5YulcBeX1+yuslW8yEHaO3RDycnG8tWUahTBx56SMYzZhTtnhpdwNlFapVc/pa+9SxKUO06xuyiiapcRRUgMxNmzpTxpP7o0zU0O5i5YpTyxlrDK69I3Ovy5dJlvgJgtRI0a9Yspk+fzpw5c5hsCkjU6XR89tlnJSxaJcdvJOiAq0DyQktgZtgGunSRKQcPyva0ohKjaZKNAvDMJOzipWv8wvN+9OrgY0PBKiAvvwx6vQRmHjhQ+HydHvwlKYFLG+ns3wlXB1eupFzhaOzR7OdQ1QuqAixZItXHa9VC63kSgBh3iD96Z3aQvKIING4M99wj4w8+sK0sRcRqJUin0+Hk5MTSpUu56667ePLJJ0lISKBfv36lIV/lxdkHfOrLOGJJtit+y8Ut1PbPoHZtqTu1b58NZVSUPqtXS4lwNzcY6YwuU4Ixv1p1j7JAraVBA7j3XhkX1R2fnaSQjKP+MD3r9gSUR7ZKkdMQee45tMvyC/8vVU/r6t1xdbWhbBURcyX3xYvh/HnbylIErFaCPv74Y8aOHYuPjw/Dhg3D19eXBx98kI8++qg05KvcBJgW7EunaVHTD19XX1IyU9gduSvbClULcCXHbC1NmAApfwAQ5gYZpwdnx6UorMC8AP/6a9HSK2sPBTs9JAPXvqNf8M0eWVU0sZKzapUEf7m7w/je6K+J+33mzvb06KIaF1tNhw4wcKC0PjBvMZZjrFaCZs6cyYYNG5g+fToTJkxg2rRprF+/npkV4D9b7vA31YCJBX3Wn4QEhwBSq0S54qsAu3ZBaKhkCk6ejBaxH4Dl13XUTOtF3bq2Fa9C0qYN3HabBPPMmlX4fHtXqN1axhF/ZccFbbywkTZtDdjZQXQ0hIeXoswK22I2RJ54AlKk/UqyBxzceZfyxhaXV16R4w8/wNWrtpWlEKxWghwcHAi4oWiCn58fjqqxivV4NAYPbzACUT9kxwWtD1Ou+CqB2QV///3gdR5dUiaaHj7d0ZbuHdxVWm5xecHU4uC77+DatcLn+z8ox0tRtK9TCw8nD66lXeNU4gHatJFL6jmspOzcKSndDg4mQ2Q1ALt1wPl+SgkqLv37Q+vWUnbg669tLU2BFFkJioyMJDIykpdffpmRI0eyefNmzp8/z6ZNmxg9ejSvmDU/hXX43ybHS7voGyzxCNvDt9OiTSo6ncTqRUXZUD5F6XD+PPwh21+89BJESKHEBA+4cHCIWnxvhYEDoWVLWYBN1X8LxH+sHOPBPu2XXBXc1bZ0Jce8g3H//VDLA6IjAfgu0gHvtE40aGBD2SoyOp3FGJk9GzIybCtPARRZCQoICCAwMJBnn32W1atX06dPHxo0aEBISAirVq3iWXPvEIV1+JuqbEZm0rh6DH7ufqQb0jmSsJ0WLeSSWoArIbNny5bNoEHQsiVaRCgAmw1AWIhSgm6FGxfgwhqrutSGGoEyjvgpOy5o/fn1alu6MnPhAvz2m4yffx5i5qIzgMEFft7ch26dHZQ39lYYPRpq14bISCkDUk4pshJkNBoxGAwYjcY8XwaDoTTlrLz49JQibxmgi/vOUrDt/Aa1JVZZSUiQwmwgi29GFFyOA2BuhB26iO507GhD+SoDY8dCrVpw6ZIESRdGwF1yjDhKv3pSOXrzxc106CwK1N69RW9Sr6ggfP65BO+at24ifgbgjAsYz/dXhsit4uQEzzwj45kzy212geoib2v0dlCnk4wj/ssVF6Rc8ZWU776TMsTNm0sQb+QsdBpkusOq3V1o3dQNNzdbC1nBcXKSpo4g1bgLW4D9TKnyMUZa1YjE28WbpIwkEtz24OUFaWlSyUBRSbh+3bJV+vzzoGloEUcB+OM6EBaSvf4qboEnngAXF9i/X5JAyiHFUoISEhJYsmQJn376KUuWLCEhIaGk5apamLPEIuMYWL8xALsidtGqg9Tr371bDBZFJSArS7ZoQFo96HQQ8TsARxyA8/3U4ltSPPmkLMD79knwa0F4tQZXNzCAPvb77EzNjRc20LmzTFHGSCXihx8gMRGaNIEhQ+DqSnSpWWj2MOukK0R2oFMnWwtZCahRAx42fb+V0wxyq5Wgbdu2Ub9+fWbPns2uXbuYPXs29evXZ1sFaphW7qgzQqpHJ0Cg3XqCvYLJMmZxxW0L1aqJ0+DYMVsLqSgR/vhDYhF8fOCBB8CYiRZ1FoBf4oGwvkoJKil8fMDU2qfQBVinA79eMo7cSD/TtrQqmlgJMRjA3OFg8mSpMh4xD4A4D4g+0YemjR3w8rKVgJUMs7H3zz9w8qStpbkJq5Wg5557jrlz57J161YWL17M1q1bmTdvngqMvhUcvcC3noxzVI/eeGF9tjWiFuBKwqefynHiRPFSXFmMLsOI5ghfnrGH8G5KCSpJnn9ejsuXw6lTBc/1f1SOEde5rUEQAFvDt9K+UzqggqMrDX//LdmZ3t6Whp8R0nQ3NBO1FVbSNGoEd94pY/P6V46wWgk6c+YM95h7g5gYMWIEZ8+eLTGhqiT+d8sx4gQD6vUAIDQsVMUFVSZ27JAu546OMGmSnIuQ4myR7pB0vjvuLi40bWpDGSsbjRvDHXdITFBh/Q1rDZPq0SnQQLeBWm61SMtKQ1dXHr5Tp8p93TdFUTB7BZ98ElxdIfUcxEtIx+eXgPPKG1vimLM1FyyAuDjbynIDVitBrVq1Yu7cubnOffXVV7Qw53Mriof/43KMNTIwKAuAvVF7ad0pEVBWaKXAbAWNHSupowARuwH4NxUI60unTmBnZxvxKi1mb9DChQUXT7R3gVrNANBF/ZUdF3TgaigNG8qUXbtKT0xFGbBnD2zZIsURn3pKzkV+BkC6J2wKd4fodkoJKml69YL27SXDoCi1u8oQq5QgOzs75s2bx6effkr9+vXp06cP9evX55NPPmHevHmlJWPVwL0JVPMAI/gm/kmD6g0wakaMAVsAOHpUEhoUFZSLFy01SSZPluP1vZCYCjqYdRHlhi8tQkIsxRO//77guf73yTHiHIMayC8jp0dWGSMVHLMhct994Ocn44jlABywAy70xtnRnlatbCNepUWnA3PIzJw55arehFVKkKZpNG/enOPHjxMcHMyjjz7KwoULOXHiBC1btiwtGasGOh3495Fx5BZCgmV8MCGUunXFm79njw3lU9waX34pAZn9+pHdiyFCelslecGRK05wqatSgkoDnQ6ee07G5tow+eH3kByvwJBAaaS5/dJ22ndOA9S2dIUmZ9E+s3fQkApRFwBYFAec70uHDuIoUpQwo0dDzZpSu+vPP20tTTZWKUE6U/lMR0dHNm7cSNOmTenZsycOpr+Yw4cPl7yEVQn/CXKMSOKuxuJ/Dw0Lzc5OUVZoBSUlxeICNn8ZA0SsAWAXoF3sBlnOSgkqLe6/XwJhw8IkSDo/3ALByxeA2on/UrtabdKy0nBrLNrPzp3ltuabojC++kpKVPTsKVszADHzwaBhdIHvwlDe2NLEyUnqBkHRmhuXEVZ7gvIaAxgMBtq2bcu+fftKRrKqiO8gsLeDNAjxPA1IXFCbzhIXpKzQCsrPP0N8PAQHw7Bhci7zGlyOBuDbGCCsL3XrWkKFFCWMiwtMMBkZ5jpN+eE/EABd5C76mjyyEQ6hODnJr/HMmdIUVFEqpKeDOWQjZyZz5CIAwt0gM8ULotsqJag0mThR3Gxbt0oZ9nJAiVaMvlExUliJnSPUkW1F9yvraejdEKNmxKGhpG8qK7QCommWL92nn7ZEPUfPAyMYqsGv4aiMlLJg0iT5+W/YUHD5Z/8n5RiVzsimkiq/+VIo7drJaWWMVECWLoXYWAgIgOHD5ZymQcR+AFamgPF8b9Ds1HNYmtSpA/feK+PCjJEyQrXNKG/4jZFjxAWGNOgGQKRDKPb2EB0t8bWKCsTGjXD4sKTiPvKI5XzEEgDOOoMxyxkiOqvFt7QJDIQRI2T8+ef5z/PuDk6OkAn9q0ltoe3h2+nYVeKC1LZ0BUPTLNsvkyZZAn4Sd0JyGuhhdhho50OoVQvq1rWZpFUDsyful18gJsa2snCLStDly5dLSg6FGXNgZjyMbSQP65aI0OxYWpWiW8EwWzvjxkH16jLWjBApfYr+SgL7yB5gcFJKUFlgjsn68Ue4ciXvOXo78OsAgEfcDmpXq026IZ3qrcQFpJ7BCsaOHbL14uQEjz9uOR85B4Ck6nAigexq7apzfCnTuTN07QoZGRKnZWNuSQm68847qV69OiEhIUyePBmdTkdERERJyVZ5SUiQXhh54VILakhgSDvjQQD2Re2jTRcp5qVc8RWIsDD46y8Zm7spA1xdC2mZYA+fn4eMk32xs7PEaipKke7dLfVKvv02/3n+DwGgi4zmziaSmZDoHQrAgQNyu6KCYDZE7r9fWqmYiVgLwA4NHA3VIaa1MkTKCrM3aO5cUYZsSLGywwBWrFjBO++8Q79+/QgLC2POnDlomsbw4cPx8PCgd+/ePP/88yxatIijR4+WuOAVlo8/Bn//gjVgv9sAcIo5QiNvqRfk0VLqBSklqALx5ZdgNMKAAdIx3kyEFBu9Vl1HeBJwoQ+tW8uOmaKUyZkuX1C9ktr3ST+/RHioviyTB65twMdHbjlwoEykVdwqERGwbJmMcxoiGfHZiQk/xIDdpd6g6ZUSVFaMHCl1mqKjLWULbESxPUFDhgzh9ddf57fffiMsLIwrV67w33//8d577zFkyBCioqKYPXs248ePp3Xr1iUpc8XGy0uKts2Zk3+9En9TGmF0Og+3li/PxOqhgHh1y1GdKUV+JCdb0uJv7KsXIR3NN2dp2GsuENE5u1O5ogy47z5LvRKzp+5GHD3BV4Ki22knANhxaUd2XJDaEqsgzJsnafG9e0PbtpbzUV+DBkZ3+PUipB7ri06H6hxfVjg4SKYYSKKCDSl2ivyNeHt7M2DAAF599VWWLFnC6dOnuXr1KuvWreOjjz66ZUErDWPHSmzI+fOwcmXec6p3AWcnyIL7al4D4FBiKJ6ekJoq1aMV5ZyffpIWDfXrw9ChlvOpEZJnDXwVCdXie4DBUVmgZYmTkyVdvqAAaf/bZXrMKQLcJS6oTieJilYe2QpAWprF436jIRIpiQnnXHVkGoGwEJo1Aw+PshWxSvPkk5IqX9C2dBlglRJkNBqtenMPDw/69u3LC+bmaQrZ83jU1K36iy/ynqPTg19HAIKSTgKwL3of7bqquKAKgaZZvlyfeip3M7AoCcY0eOlYFQHJR0MAlBJU1jz5pPxeNm3KP13eT1LldbFZTGwrXW0z/UMB9QxWCJYuhcuXJSvwrrss5zUjRB4B4I8EDWetOsS2Us9gWePjIzF6No5EVynytmDSJPnF//cfnDyZ9xz/8QDYRcfS2a8+Rs1IzU6WekGKcszGjXDkyM1p8QARfwBwykmPUYPM0yF4eKA6x5c1/v6WdPn8jBGPFuDmBka4z0cKlp7XQgE4ezb/5DJFOSCnITJpEtjbW67FrYH0LHCAL8LAPb6PigeqwiglyBbUqwd33CHj/Bbg2veBXgfX4cVWUsY/wy8UUEpQucf8O33wQYkBM2PMhGipO7P4qgEHXCCyE506gV49iWWPOVD2xx+ztyhzodOBv9TqCko5B8Du6O00bCY9xXbvLhMpFcVh1y5ptujkBI89lvtapFSOTvK25+J1SDwUAihvbFWl2EvvrFmzuHTpUknKUrV4+mk5zp8PiYk3X3fwAF+p2tXfMRaAcyYr9PjxvG9RlAPCwy3NAc2/YzOX/4BMI5oTzD0HPqndVTyQLenZU5rZpqbCDz/kPcdPPHn2UddoW8uXDEMG9XqquKByj9kLNGZM7rR4gIiNgCQmAKSf7IOrK6ge4FWTYitB7733HkFBQXTp0oWPPvqIs2fPlqRclZ8BA2QPJCkJFi7Me46f9JnyvhqOXgdHruwnsNE11VG+PDNvnmT9hYTcvKpGfAdAYg1HrqRC1pkQQFmgNkOnsyiq+WVr1hwOdnpIhddb1QLArqFksyglqJwSE2NJu77REEm9BFevAvBlhAE3vdQH6tAh946ZoupQbCUoKiqKtWvX0qlTJz777DMaN25M27ZtmT59OseOHStJGa0mKSmJyZMn4+fnh7OzM23btuWXX36xqUw3kXMB/uILqSdzI/6SQqi7nMXYpv4YNSNBvVS9oHJLWhp8/bWMc9YkMRO5HYDQdPnn5d0hgFKCbEph2Zr2LlCrIQD9nGTLLMZFPAm7dqlefuWSr7+WOiLdukGHDrmvRX4JgNFbx6pL4Jui6gNVdYqtBOn1evr27csXX3xBREQEmzZtom/fvnz77be0atWKZs2a8cYbb7B///6SlLdIjBgxggULFvC///2PVatW0alTJ8aMGcPPP/9c5rIUyLhx4O4uwdHr1t183b0FVJPAzGeDnACwaxAKKCWoXLJ0qUTLBgTAnXfmvpZ0AhKvgw4+vZCBk94FIjoRFAS1atlGXAUSvG6OGckvPs//bgCqx0djr4ejCTtwdEslPl4CpBXliMxMS7f4G71AABG/A3DRzQ2DpryxihIMjO7RoweffvopYWFhbN++nTvvvJMlS5bQsWNH6tevzyuvvFJSH1UgK1euZM2aNXz55Zc88cQT9O3bl2+++YaBAwfy8ssvY8ivQKEtcHeHhx6ScV71SnQ68JPAzFZZ0qct2jkUUB3lyyXmL9GJE2/2rUeYUuN97NgUBQFadzA4qSKJ5YGJEwvO1jSnyscZGdu4OhmGDBr0UXFB5ZI//oDISLEsRo7Mfc2QDtGnAfg5Xlo1RG0PAVDPYRWmVHJSOnfuzAcffMDp06fZu3cv999/PytWrCiNj7qJP/74g2rVqjFq1Khc5x9++GEiIyPZWd5WraeekuM//8C5czdf95fATOfY69R0g9NJ+9G7XiM6WgreKsoJO3dKupCjY+4mjWYi/wEgzNUdDXCICAGUBVouKCxb0y1Ysvw0eLauCwDurUIBVTm63GH+/T3xhDyLObn8J2QZwRlmncrA3b46hsjW1K4tpYQUVZNST8xt27Yt06ZNK7P+YUeOHKFZs2bY32CJm1t3HDlyJM/70tPTSUxMzPUqE5o0gUGDxK0zd+7N12venR2YOaW5N0bNSHBvFRdU7jAvvqNHg69v7mtZyRBzAYD5l8UCjd0VAiglqNxQWLamXy8AWmZKUO11b4kLUs9gOeLgQdi8WbywTzxx8/WI7wFI8XUjNgWCsMQDqc7xVZdKV50kLi4Ob2/vm86bz8XFxeV534wZM/D09Mx+BZalaWAOov3uO0hJyX3NzhlqS2DmSC/Z/6rWMhRQC3C5IWc2Sl4B0dG/gFFDc4MvTqbgbOdC/OFOqnN8eWLAADFI8svW9JM2G06xqdT1gLPpO8A+jf37IT29jGVV5I3ZELnnHmnOeSOR2wDYq3MHwP5SCKC2wqo6lU4Jgtzd7ot6bcqUKSQkJGS/wsPDS0u8mxkyRFzyV69CXsHb/sMB8LuegJ0OErxDAeWKLzd88w1kZIhbp2PHm69HLgAgxdeDa2nQ0LkbGJxo1Up1ji833JiteWPAnc9gcLSHDHilsRsZxnQ8mu8gI0McEAobEx8v/fog74Do68fhehLo4MOz0n4oZlcfQHljqzpWK0Hx8fGkmLwVBoOBhQsX8tNPPxXYXLUsqVGjRp7ennhTRdi8vEQATk5OeHh45HqVGXZ2ltigvBZgv0kA6OOMDA6A8Iz94JTAnj3SIFlhQzIzLduYeXmBNA0ipbTwbuRvz/NqCKAW33LH+PH5Z2vq7aFOMwBGeIohVbNzKKCMkXLBd99J0cu2baFHj5uvR8gzavS1Z+WFVDwdvYja3xqdLm+7RVF1sFoJGjRoEKdOSen/KVOm8NFHHzFz5kyef/75EheuOLRq1Yrjx4+TdYN2cPjwYQBalteyoA8/DC4uYlZu2ZL7mmtQdvuF54NcMWLEpckWUlJUR3mb8+eflmyUG4LxAbi2B1LSwA5mnL4GQNLhEEC54csd7u6iCEHe2Zp+9wJQKyEZRzvIDAgF1La0zTEY4Eup/8Mzz+Qd4BO5HIBYz1oYNWjq0hs0O5o2BU/PMpRVUe6wWgk6ffo0bdu2BeDHH39k1apVrFu3jiVLlpS0bMXi7rvvJikpid9++y3X+QULFuDn50eX8mp+e3vD/ffLOK8MFf/eAHTVi3Ln2ykUUAuwzTF/WU6YcHM2CmT3KTLWdOS/sGs42ztzKlS0n/L6p1ilMXtkly+HsLDc1+pIXJA+QePuuhCpk7gg9QzamBUr5Hfl7S1tMm4kMwliJTFhRbLEA3nEhwDqGVQUQwlydHQkOTmZ3bt3U6tWLQICAnB3dyc5Obk05LOaIUOGMHDgQCZOnMg333zDhg0bmDBhAv/++y8ffvghdnZ2thYxf8x72b/9BhERua+ZAjPdrmQQ6A4Z/qGAUoJsSs5slCefzHtO5L8AxHgGANDKqxup151wd1ed48slTZvCwIGyjWn2Lphxrgk+kvk30c+RTC0dAnZw+nTe/VcVZYTZEHnsMfGm30jMEklMqAbTDkis5/UjKh5IIVitBI0ZM4aQkBDGjRvHeJPreP/+/QQHB5e0bMXm999/58EHH+Stt95i8ODB7Ny5k8WLF3O/2dNSXmnTBnr1EvfuV1/lvlZjMDjZQyY8VQ9i9fvAKUHFI9iSwrJR0uPhSiQAfyd5AVAnvS8AnTpJKJiiHGI2Rr799uZsTb/+AHTWSZub6u1DAdVR3mYcPw5r14JeL0Uv8yJCEhMya9XgQkIyXs5eHA9tAyglSFEMJWjWrFlMnz6dOXPmMHnyZEAyrj777LMSFq34VKtWjVmzZhEVFUV6ejoHDx5k9OjRtharaJiDa7/6Knfurd4O6jQH4F4vB4wYoe4Wjh6F69dtIGdVp7BsFICoBaCB5gkf7Bd3vMFUpl/FA5Vjhg2D4GDJ1ly8OPc1f0lScLmSRRMvcGyo6gXZlDlSiZ077pDf2Y3kSEw45RQEQPsavUm4aoezs+ocryiGEqTT6XBycmLp0qXcddddPPnkkyQkJNCvX7/SkK/qMXy4eBViY+HXX3NfMwVmBiVn4qCXqrWqo7yNKCwbBSBSyh1k1q7N+WtxuNi7cH6rigcq99jZwSRRdvj889zZml49wcURDPBMMMS5bldxQbYiMREWiJcnz8xMkMSEVElM+CFS+i/WThVvbIcO4OBQFoIqyjNWK0Eff/wxY8eOxcfHh2HDhuHr68uDDz7IRx99VBryVT0cHCzxJTcGSNd5AnSgT4QRdcG+YSigrNAypyjZKEYDREkBmZNOUuyyi193jh+WhVgpQeWcRx4BZ2eJ+9q61XJepwO/tgDcWU1PFungv1P18rMFCxZIcctmzSA/I9yUmKDVduTbg5IhrJqmKnJitRI0c+ZMNmzYwPTp05kwYQLTpk1j/fr1zJw5szTkq5qYM4127syt4Tj5gE9NAB6tBddcJS5IKUFlTGHZKABxmyA9ExxgwSVnAOrrQ9A0aTJfp07ZiasoBjVqwAMPyHj27NzX/CW2MCDJiLM92DUIJS4u79Z/ilLCaLQYiU8/nX/fiwhJTEjyaUxiRhLVnatzdpu0UFJb0goohhLk4OBAQEBArnN+fn445pUerCgetWpJDyq4uV6JKTCzu16HZooLUlZoGVNYNgpApAS2a3VcWHj4AACOkSGAskArDOYtlt9/h5wV5Gs9DHrQJcODweDaIhRQHtkyZc0aOHUKPDxg3Li856THQZwkJmw3NACgR2BvDh2Urz31HCrACiUoMjKSyMhIXn75ZUaOHMnmzZs5f/48mzZtYvTo0bzyyiulKWfV49ln5bh0KURFWc77Sx0Tt6sajT1BV38DUVGqo3yZcfRo4dkoAJFScTjZtxWXU67gYu9C1N5OgFp8KwytW0NIiGx/5mxu7OAONcUQfNAHUrxVXFCZY/bOPfwwVKuW95zIRaABXvDtSSnh0kAfQmYm1KwJQUFlIqminFNkJSggIIDAwECeffZZVq9eTZ8+fWjQoAEhISGsWrWKZ81f2oqSoUMH6N5d2jLkTJf37A6uThKYGQTOTUMBZYWWGWYv0PDheWejAKSEw9UrAGzPkmJAPer2YM8OiQdSbvgKhHld+/prCYQ34z8YgI4aGHRSL0g9g2XEqVOwcmXufm95EfkjAJpfAKvO7ADAMVKCojt3Vp3jFUKRlSCj0YjBYMBoNOb5MhgMpSln1cS8AM+bZ0mX1+nAT2pcDKsGaV77wfmaWoDLgvh4S4fxgpT+yO/k6KNj4clEANpXDyEiQhxIHTqUrpiKEuSOO8RlEBeXO13eT7bKXK5CuxpAcKjqKF9WmA2RYcOgYcO85xizIFISE2I8O5FkigcK39MKUN5YhYVidZFPSEhgyZIlfPrppyxdupTExMSSlksBMGKEpMvHxOROl/eXgM2gZHC0N0LdzUoJKgu+/Va8AW3aQO/e+c+LkN+V5teAf89KHzhz09SWLfP33ivKIfb2llYas2dbgu/cW4OHK2jwTF2wb7RBdZQvCxISYP58GT/3XP7zrmyAzCxwgpXxgQD0Ce7Drp0qHkiRG6uVoG3btlG/fn1mzZrFrl27mDVrFvXq1WPbtm2lIV/VxsHBUq9k1izLAlzrEbADfQqMCwKCQ9m7V3WUL1WysiyF2Z59Nn9felYqRJ8A4Er13lwxxQPFH1HxQBWWxx4DV1dLmxQzfvI7HeQMxjo7wD5VGSOlzQ8/SFp88+bQv3/+8yJMIQR1XPn1xEkAOtYIyc7g69SplOVUVBisVoKee+455s6dy7Zt21i8eDFbt25l3rx5KiaotJgwAZycpCKieYW1d8sOzBznC/oGG0hJgSNHbChnZeevv+DiRfDxgbFj858X8w8YjOAKqy/L76hH3R7s3iHZk127loWwihKlenV48EEZ50yX938YAL/r4OqUAYHblRJUmhgMlq2wggwRgIj1ABj9urHlotR5MjdNbdYMvLxKUU5FhcJqJejMmTPcc889uc6NGDGCs2fPlphQihz4+lpq0eRagIcA0EEDY8394BKvFuDSxPyznzBBiujlR+S3cvT35M9TRwHoU7dvdlVv5QmqoJjT5f/4Ay5ICxR8x4CDDl06PFEPCN6gnsHSZOVKKcbk5WWp4ZQXSWch8Sro4KxTH5IykvB28SbmsIoHUtyM1UpQq1atmJszXRT46quvaNGiRYkJpbgB8wL8668QKXUv8Jf9cJdr0NkXFRdUmhw4AJs2SXyIeXsyLzQNIsTq1OqEsPGC9JUKNISQkoLqHF+RadFCtl+MRsu2qN4R6kj9mbHeQHAoZ85IDLWiFDAbIo8/Dm5u+c+L+F6Ovjr+Pi8hBH2CVDyQIm+sVoLmzZvHp59+Sv369enTpw/169fnk08+Yd68eaUhnwKgfXvo2TN3XIpbC/CqJoGZgUA9ZYWWGrNmyXHkSPD3z39ewiFISQY7uOA8iCspV3B1cCXpVEdA0nJV5/gKjHnL/5tvJC4FwH84AC0zgYAd4JDCrl02ka5yk7M+lzlQPT8ilsnRvzFrzkmsap+gkOz1UW1JK3JilRJkZ2dH8+bNOXHiBAsXLuSpp55i4cKFnDhxgpaqHW/p8sILcpw3D1JSZOwnT/NAJyB4A8ePS09BRQkSG2tJjS4s7i3iaznWsmPl+QwAegT2YO8uiQdSFmgF5/bbJSX72jVLhlKd59B04JgItwVkQeBWZYyUBmYv0PDhBVc5zLwOsacBMNQZwZaLkp0ZpIWQmCjx7eqrSpETq5QgzZSd5ODgwOeff0737t3p2bMnDqoVb+lz551Qv37uWjX+EwComQjeQYfQnONUR/mSZu5cKf7SqVPhJmTEcjn6t2LdeckiCgm2WKBKCarg6PUwebKMZ82SrTHnAHQ+3gBM8gOCQ5UnqKS5csWy5hWUFg8Q/RcYNagGB5PbkZyZjLeLN3HHRfPp0EF2tRUKM1YpQboc0fi//vorkeb4FBOHDx8uGakUN2NnZ1kAPv1UFuAaI8DJDl0mPB0MBG9UVmhJkpZm2X584YWCs1HSrsAV6S9lrDOa0LBQADr6hHD8uExRSlAl4KGHJFvszBn45x8559cHgF72QL0N7NqlevmVKPPmybPYoQP06lXwXHM8kH91/j0nHqE+QX3YvUu+6tRWmOJGiuUJyguDwUDbtm3Zt2/fLQulyIeHH5aGgadOwapVoLcDvyYAjPICgkOVElSS/PQTXL4MgYFwQ0bkTUQtkaMXHE/pQHxqPG4ObhjDO6Fp0mGjVq3SFlhR6ri5SYYgwMyZcvSX1g1e1yCo8U7iriehkmVLiLQ0S7f4wgwRzQiR0h4D//5sCNsAQL96/ZQ3VpEvxaoYnR8FKUmKEsDdXTIjIMcCPBKApumgq7dedZQvKTTN8jN+7jkpXFkQEeYtytr8d1Y8or2DerN3t9ynLNBKxNNPy57Kxo2wfz949gU3R3RGmFzfqOKCSpLFi6VifkAAjBpV8Nz4PZCWCvaQ4T0+Ox6oa61+HDokU9RzqLiRElWCFGXAM8/I1tj69VLBts4zaHqwT4YR7Y4SnXiZixdtLWQlYPVqOHZMelw89ljBc42ZELVfxn5Dsi3QvsF92WEyTJUFWokICIB775Xxp5+Kd8Jf+vndUQ2otyH79664BXIaIs8+WwRD5Bs51rFnZ5QraVlp1HKrRXJYM4xGSewsKLlTUTW5JSXo8uXLJSWHoqgEBVm2Zj77DBx80NWsCcCE2qi4oJLik0/k+Nhj4OlZ8NzLGyEzE5wgq/r47PpAfYOVG77S8vzzcly8WGp3mfr51UsBhwbr1DNYEqxZI2Xwq1WzeMALInKFHP3asO78JgD61uvLzp2yhaaeQUVe3JISdOedd1K9enVCQkKYPHkyOp2OiIiIkpJNkR/mdPmff4boaPAfAEA3PRCsrNBb5tAhS02SwrJRAC6Z+hT5ObE/xonE9ES8nL3wSmvL5ctiwLZrV7oiK8qYjh0lSNdcu6vm42j2oE+Dh7vuZd/R66Sm2lrICo7ZEHn00cL7XKRcgvgoGfs9kMsbq+oDKQqi2NlhK1as4J133qFfv36EhYUxZ84cNE1j+PDheHh40Lt3b55//nkWLVrE0aNHS1zwKk2XLtCtG2RkyALsJ/Vr3BOgSbtVSgm6VT79VI733CMRzQWhaRCxRsYBXVl/XrxAkpEilRHbti2404aigmL2Bs2bB2lGdHXqAjDeV8Pgv5n9+20oW0XnyBH477+iGyIRP8rRB1Ls72LHJVkEcypByhOkyItiV0wYMmQIQ4YMyf53fHw8+/btY+/evezbt499+/Yxe/ZsNE1Dp9NhMBhKRGCFiRdekEDBOXPg1VfRPJzRJabxfLvzPDs3mvT02jg52VrICkhUlGSFAbz4YuHzE49DUoKYE7UfZcNGubdvcF92/iJT1OJbSbnzTmjQAM6ehe+/h6HDIHwu7YyYPLJD6d7d1kJWUMyxQCNGQL16hc+/ZHpm/QPZFn6WDEMGAR4BOKc0JCJCwig7dCg9cRUVlxJLkff29mbAgAG8+uqrLFmyhNOnT3P16lXWr1/P+++/f8uCKm7g7ruleu3Vq/Dtt+j85Qkf6goZfqEcPGhj+Soqc+ZIfE/37kXTXi59J8daOjL0Q9h8UYok5kzLVW74SoqdnUVR/uQT8H0aDXC5Dl27LVce2eISHW0xRMxb/wWRmQQxx2QccDcbzlu2wnbtkt2LVq0KbjemqLpYpQQZjcYiz9U0jS1btvD555/zxhtvWC2YohDs7OCll2Q8cybUlMDMgCRwafyvWoCLw/XruYsjFoWI3+QY0ITdESdJyUzBx9WHhp4tMJfMUp6gSsxDD0HNmtJZfvl+NO9qADzb7iRb98XbVraKyqxZstXftats+xdG9HIpHlsN8HiE9WHrAVUfSFE0SjxF/uzZs7z++usEBgZyxx138Mcff5CVlVXSH6MAGD9eKvCFh8MaRzRHHbpMmDhIWaHF4uuvpS9U48bSo6gwUmPgygUZ+49h/XlZfPsG9+XwIT0ZGVCjhuyYKCopLi6WnnIffojevwcAAxwh0mEjKk/EShIS4MsvZfx//1e0ey59K0d/d65n1GN3xG5AlahQFI0SUYLS0tJYtGgRISEhNG7cmBkzZhAYGMjcuXP54YcfSuIjFHnh7GwJGvx4Jsbasnf+YEA8mw+H2U6uikh6uiUO4ZVXitbuPXKpHL0B13F5ZqR06VJwkVtFJWDSJEnjPnQITsjep891qN3mL5Uqby1z50oX6ObN4Y47Cp9vNEDkVhkH9GfzxS0YNAP1q9fHv1pQdi9FtSWtyI9bUoJ2797Nk08+Se3atRk/fjynTp3ipZde4tixY2zfvp0JEybgVVhqo+LWmDhRKkkfPYrdmY4AtMiESw5riYmxsWwViYULTfVe/OGBB4p2T8QCOfrXJC2rNtvCtwHihlcWaBWienVLK40vN2J0s0NngBduW8727bYVrUKRmiq1zwBefVUywwrjylYxYBwB38dyxQMdPixv6ekJTZqUmtSKCk6xlaDWrVvTtWtX5s+fz4ABA1i+fDmXLl3igw8+oGnTpiUpo6IgvLzgiSdk/N1FND04pMA9QxaoLbGiYjDAhx/K+IUXKFJaXVYKRB2QccDtbA/fTrohnTrV6tC4RuPsn72yQKsIkydLK43QUIxX6wMwolY8G/dFFnyfwsL8+dIio25dGDOmaPdEmLbC/OxBPyDPau2dOxdNn1JUTYr9p3HkyBF0Oh3PP/88X375JcOGDUOv/tJsw+TJUpFvyw4yL7sDMLHjLrbvUE3EisRvv0lXcG9vi0VfGNGrRXlyBbwezY4H6levH5cv6zh3TrbBlCeoihAYCPffD4D9MlcA6iXD8eR/ycy0pWAVhKws+OgjGb/0UuEtMsxErJSjf3uupqawL0qyEfrW65vthStKbLWi6lJsreWzzz6jVatWfPDBBwQEBDB8+HD++usvFQRtC/z94cEHAbD/S1o8dLXPYP2hE7aUqmKgaWAu4fDMMxLbURTMFqi/M+i65LJAzYtv8+aFd9xQVCJeeUWOKw+hXQZ9Jjxy51fZzTsVBbB0KZw/Dz4+UiG6KCSehMQ4+RarM55NFzahodGkRhP83P2UEqQoEsVWgp599ln279/Prl27ePTRR9m4cSMjRozA39+fF198kcOHD5eknIrCePll0OnQr7uEdgncroO913conbQQ/vtPOoG7uooSVBQ0I0SEyjigF0kZqeyMkAjYfvX6qcW3qmIO5tU0sv5zAWBci/1s3648sgWS0xB57jl5FotCxCI51gQcR97gjRXnLihvrKJgbnn/qmPHjsydO5eoqCjmz59Ps2bN+Oyzz2jbti0dO3Zk2bJlJSGnojCaNpXqqoBxhfxanxm4BNWxpBDMi++ECZLPXhTidkFaitRbr/koWy5uIcuYRZBnEPWq18tWglQ8UBXktdcAsP8vHWKhtTGTtXvP2Fiocs6qVXD4sHhhn3qq6PddWiJH/2CgZnZ9oJzxQE2bSty6QpEfJRbE4+zszIMPPkhoaCinTp3i1VdfJTo6mp/MlT8Vpc+bbwKg32yESzCgRgTbtqt2JfmybRuEhkr8QVGLIwJEmMo++OnBbhhrz60FYED9AWRmwm4pU6I8QVWRrl1h0CB0WUa05eCQCh7VP7W1VOUXTYP33pPxk08WXWNJuwJXTMplwCiik6I5EnsEHTr61rMoQeoZVBRGqUQyN2jQgPfee4+LFy/y999/c9ddd5XGxyhupE0buPtudBrwJ/gka+w4+Z+tpSq/TJ0qx3HjJLC1qFz6U47+LYBquZSgQ4ckLdfLS6xQRRXkf/+T40bgMjzQ4XeuXLGpROWXtWth61bJyDQ3pC0Kkb+BBngBbuNYd24dAO3qtMPH1UdtSSuKTKmmc+n1em6//XZ+//330vwYRU7eegsAbQfoLkHXejNsLFA5ZcsWWLNG0pqtaeuSeBISYkEH+D/I5eTLHIyRRm0544G6dFFpuVWW7t1hwAB0BuAv6F4thh0qU/NmNC17vWLiRPDzK/q94aaefYFeQAvWnjcZIvUGkJUFu3bJZaUEKQpDLdOVjbZt4a67sr1Bt9Xdw9WrNpapPGK21h95BIKDi35fuCkYsxbgODo7GLNNrTbUdKup3PAKwfz3tQmqhcH+4yo28iZWr4YdO6T1yKuvFv2+zCSI2ivjwGFoGrm8sUeOQHIyeHhIrLpCURBKCaqMmK2rbRAclsq27Qm2lae8sXEjrF8vsUCvv27dveGmGLfAICAw1+ILKDe8QujZE/r1AwPwN7T2/cDWEpUvcnqBJk2C2rWLfm/UP5aGqZ6PcSruFJcSL+Fk50TPuj2VN1ZhFepPpDLSvj3a7UNBA/0ySIi18ou+smO20h97TKrTFpXkcIgPk3HAGDRNY825NQD0r9ef2FhUkUSFhRyxQT0yD2FQOQoWVqyQDAJXV0t9paIS/o0cA11A1zPbEOke2B0XBxdliCisQilBlRTd1HdksA06XV9iW2HKExs2iCfI0TE7nbnIXPpFjr6AyzjOXT3HhYQLOOgd6BXUSxVJVOSmd2+SuzYAA9RYksnRQ+G2lqh8oGkWBfGpp6BmzaLfa0iHiM0yDuwH2FvigW7wxqoSFYqiUKmUoOvXr/PKK68waNAgfH190el0TDVnAFU1OnQgrqsvaNDg1ysYspQZmssFP2ECBARYd3/4fDkG+AJNsy3QboHdqOZYTVmgiptwfu9rAHShkLprkm2FKS/8/Tfs2wdublLk1Rqi10JWJrgANR4ny5iV3TR1QP0BXLliKZKolCBFUahUSlBcXBxff/016enpDB8+3Nbi2JyUaVKfRL8FLv3yPxtLUw5Yt06ywpycYMoU6+5NuwKXj8s4cASgy5WRAioeSHEzdn37cb2NPRig9eLVthbH9hiNFi/QM8+Ar691918ybYUF2INuEPui9pGQnoCnkycd6nRQRRIVVlOplKCgoCCuXr3Kxo0bmTFDpYYHDriflD460MB3ZhUv2KZpllT4J5+0Lh0XIOI3eY/qQLWHMRgN2ZlhqkiioiA2jxsEgPOmTNi90cbS2JjffoODB6U69EsvWXev0QCXTHXPArsBLtne2H71+mGnt1OGiMJqKpUSpNPp0Ol0thajXLFuZGOwA9f9KVIXp6qydCns3Cku+P/7P+vvN9clCXAHOnEg+gDxqfG4O7rTyb9TriKJTZqUpOCKio7/iM8w9gSdBpmTx9taHNuRnm559l54oehtasxc3gLpqeAI1JQmqyo7U3GrVCol6FZIT08nMTEx16sycKLWy2iDTP948QlxR1c1ci6+r75qXTouQOZ1iDbXJRkK6LMX3771+mKvt1dpuYp8aRPciBPjHcAOHLZdkCrJVZEvv5T0ydq1rY8FArj0vRwDdKC/i5TMFLaGbwVECVJFEhXFQS3XJmbMmIGnp2f2K9CaNgrlmHt7jyR2HOAMHD4Pv/xia5HKni++gLAw2QKzpkeYmUhTXRJ3wPMRABUPpLCK33VdoL/pHy+/VPWMkfh4mDZNxtOmyXaYNWgahP8l44DWgBdbLm4hw5BBoEcgjbwbqSKJimJhb2sB8iM0NJS+ffsWae7+/ftp27btLX3elClTeCHHF2RiYmKlUISCannyvVaTR+6IhV+R4oD33CPBwVWBuDiYPl3G06fLdpi1hH8rx0An0PUlLSuNLRe3ABY3vKoUrSiIuMynMYzdgt0m4MBBMUbGjrW1WGXHu+/C1avQsiU8/LD198fvhZQE+caqLffn3ArT6XTZhkjnzoV7Yw0GA5mZmdbLobAp9vb22NnZlWjYS7lVgpo0acI333xTpLl1rSl4lw9OTk44VVLF4K9T9/LwsC/QrUE8InPnwuTJNpaqjJg2Da5dg9atpVGqtRjSINJUlySgL+DAtvDNpGWl4efuR1OfpqpIoqJQHuh1GyfCocXtwDKqljFy9ix8/rmMP/4Y7Oysf4/wBXL0A+xHATfHAxXFENE0jejoaK5du2a9DIpygZ2dHTVr1sTT07NElKFyqwTVqVOHxx57zNZiVApa1R5DXO0v8LkH+A7xiDz0kETxVmZOn4Y5c2Rc3MU3anWOuiSmrbAbLNBt22SqKpKoyI/2zb14bWtdZgy5iLYOdGFhEiNjTef0isqUKZCZCYMGwW23WX+/pkH4UhkHNAD8uJJyhf3R+wGp1g5kP4cFKUFmBahmzZq4urqqRJoKhKZpZGVlkZiYSFRUFKmpqdSpU+eW37fcKkGKkuP+kM78e8yeB/pkof2nQxceJ+niX3xha9FKlylTICsLhgyBgQOL9x4XpNgdde1ANxSA/85Kmq558d0iO2P07HlL0ioqMTodrDp8P+/2nYF+BGKMTJ0Ko0dDCSzk5Zbt2+HXX+UH8NFHxXuPawfheizYAf7izTUbIi1rtqRWtVrExEiRRJ0ufyXIYDBkK0A1rM1MU5Qb3N3dcXJy4sqVK9SsWRO74hi3Oah0gdGrVq1i2bJlLF++HIBjx46xbNkyli1bRkpKio2lsw1NG9szc2svcADdg5qc/PJL2LvXtoKVJhs3Sk0SvR4+/LB475GVAhGmsgJB3QE3YpNj2WvqYD2ogaTdbZUEFXr0uEWZFZWa7n63c7YaEAJaMxdITLS+Vk5Fwmi0eLoeeUS2pIvDBVNWmB/gMAaA1Wel8OTgBoMByzPYsmX+Dm5zDJCrq2vx5FCUG9zc3NA0rUTiuiqdEjRx4kRGjRrFI4/I1sWvv/7KqFGjGDVqFLGxsTaWzjbodJB1dSTXqgMtgGH+4mKeOJFK2dUxPR2eeELGEybIylgcIlfIVpgbUONJANacFaWobe221K5Wm9RUiy6pPEGKgri3R2cWxznKqjs+TR7Mn3+G9ettLVrpMG+e1OZyd4d33inee2gaXPhZxkHBQCM0TePfM/8CcFtD2V6zxhurtsAqPiX5O6x0SlBYWBiapuX5Cg4OtrV4NuO2hoNYa1KatdFXJY90924oYvB5heL99+HkSahVC957r/jvk2sr7A4A/j0ri6/ZAt29W8Id6tSBKvznpSgCXTrZ8+naAWhOoAvUYHwfuTBpEmRk2Fa4kiYy0tKa5r33rK/QbiZuNyTHSeCGn2SFHYo5RHRSNK4OrvSsK1qP8sYqikulU4IUeXNnzwZ8fLQ26ECnT4E3JsiFKVOgMnnITpywKD6zZhW/gVBmEkSGyjioN+COUTOy+ozJDd9QlKCcFqgyMBUF4eICXteHcdHDdOL+q9JB/eRJ+OQTm8pW4jz3nGz3dekiHuficvEHOfoD9g8AZHuB+gb3xdnemZQU6ccKyhursB6lBFUROnbUsXvbnVz3Np0YcBratZP08VdesaVoJYemSV+wjAwJhr733uK/V8TfYMiCakB1URgPRB/gcspl3B3d6RYo0ZfKAlVYQ9+6g1hyXcbatSPw0fvyj2nTpHxFZeCff2DZMsnG/Oqr4mVlAmhGuGAq7hrUAKgP5IgHMhkiu3ZJ/oO/P5RAtZQKzaFDh3j00Udp0KABLi4uuLi40KhRI5544gn27NlTrPcMDQ1Fp9MRGhpaYnKaW1yZX25ubjRr1oy3336b5OTkEvucoqCUoCqCiwvU1w1kdZbpxKU1Ehyt08GCBbBpk03lKxHmz5eAaBcXSY2/FdeMeSssyN6yFWayQPvV64ejnSNGoyUtV1mgiqIwtGsDPtgdgOYCugwD9E2CPn2k8dyzz9pavFsnKQmeekrGL7wAbdoU/70ub4PUa+AA1JFeYdfTr2cXKlXe2Nx89dVXdOjQgZ07d/Lcc8/xzz//sGLFCiZPnszRo0fp1KkTZ8+etfp927dvz/bt22nfvn2Jyjty5Ei2b9/O9u3b+euvvxg5ciTvvPMO44pTz+0WUCnyVYiBDfrx7jkYWRO4ngLdk+Dxx+Hrr+HRR2H/fuvL2ZcXLl+2ZNq8/TbUq1f898q4BlGmlTWoHxIZbVGCzIvvsWPiSHNzu7W1XlF16N5dR/yCoZzr/zUNUoEL8+DLX+UPaPlyCZSuyJWk//c/uHhRAuT+979bey/zVlgAYCdZYRvCNpBpzKRB9QY09G4IKG8swNatW5k0aRLDhg1j2bJlODo6Zl/r168fTz31FL/++isuLi5Wv7eHhwddu3YtSXEBqFWrVq73HTBgABcuXOCnn34iLS0NZ2fnEv/MvFCeoCpEv27eHDjYlXjzlljY+xJEHBAgRTaK01erPKBpEoMQHy9fJrdaDfvSn2A0gAfgKQU7E9IS2BYubp/bGuTOSOnaFeyVOaEoAn5+4Jt4G99dNZ2IOA6NAuDNN+XfEyfChQs2k++W2L0bPvtMxl9+WbwWNWaMBrj4q4zrNgGCAYshYn4GDQbljQV47733sLOz46uvvsqlAOVk1KhR+OUIUN+zZw+jR48mODgYFxcXgoODGTNmDBdu+PvLazusqPdai7kK9K3W/rEGtXRXIbp3Bz4fyl/pO3gY4OJmaOcBCxdC//6SKTZ0KAwfbltBrWXhQli8WGIPvv4aHBxu7f0umDLmghxANwyA9efXY9AMNKnRhHrVxcukLFBFcegbNIDZJ/VMb21En6TBpc/gtdfg33+luOCDD8KGDcWPpbEFiYlS+NFohDFjJCbvVojdCGnXwRGoLeVOcqbGm72xR4/KR1erBq1aWf8xmgblrXycq6t123oGg4ENGzbQsWNHqyooh4WF0aRJE0aPHo23tzdRUVHMnTuXTp06cezYMXx8fErlXjPmCtAASUlJbNy4kQULFjB69GgcbnUNtwZNkScJCQkaoCUkJNhalBKlToc9Wv2ZaNqvaNpPaFrUErnw8suaBppWo4amRUbaVkhrOHFC01xdRfbp02/9/dKuaNrPevnZXBuSfXrC3xM0pqI9u/LZ7HP16snH/vffrX+souowZ46mMb6vdmit6RlcX08unD2radWqyR/Ve+/ZVkhrMBo1bfRokTsoSNOuXr3199x5v/xsdqBp2gVN0zTt1JVTGlPRHN5x0K6nX9c0zfSzRNMGDiz8LVNTU7Vjx45pqamp2eeSkuT+8vRKSrLuRxUdHa0B2ujRo2+6lpWVpWVmZma/jEZjvu+TlZWlJSUlaW5ubtqsWbOyz2/YsEEDtA0bNlh9b34Aeb6GDBmiJRXhB5DX7zIn1nx/q+2wKkb/5u04F1mLaHPm+IWZcpw2Ddq2la7rDz8sFl15Jy1NrM+UFOjbF/7v/279PcN/k6wUL8DzccBkgZ7NbYFGRsL581KQuhS2yxWVmF69gNNDmWuuTBF9HtJioH59S6PRt96CYmbzlDk//AC//CKeq8WLb70noTETwv+UcVAzQFK+zF6gXkG9qOYosYvmLWnljc2bDh064ODgkP36JEcphqSkJF599VUaNmyIvb099vb2VKtWjeTkZI4fP17g+97KvWbuvfdedu/eze7du9m0aROzZ89mz549DB48mPT09Fv6f1uD2g6rYvTprefHf4awpOd8ngO4uAc6pks3659/hvbtYfVqya565hlbi1swr7wCBw6Ajw/8+GPJbB9kZ4U5AqLwnLhygosJF3Gyc6JPsBS4M2+FtWkjBXEViqLSogV4xAzl+7MvM7s92F8DLn4AjWfC+PGSYv7bb3D//VIA51Zia0qb48fh6adlPH16wd1Li0r0WkhPBiegpqWJ9o2FSsHyHBY3HsjVVRLayhPWdvXw8fHBxcUlz3icn3/+mZSUFKKiorjzzjtzXRs7dizr1q3jzTffpFOnTnh4eKDT6Rg6dCipqakFfuat3GvG19eXjh07Zv+7V69e+Pr6MmbMGObPn88T5qr/pYxSgqoYvXsDnw7l/dPzebYB6FINEPkDBD4JzZpJt/Wnn4aXXxYXR6dOthY5b/76y2I1L1hQ/Iq0OUm+CDGmHhhBg5HW8RYLtE9wH1wdZIVSFqiiuOj10Lt5M/6JC2Kv4wW6AFxYLEqQTie1dbZvh1OnpJr0/PnlM/c7NRXuu0+OAwaUXL2xc3PlGATo7wMgLSuNDec3AJZWGeHhkohmZyc1GYuDTle+dcyiYGdnR79+/fjvv/+IiorKFRfUvHlzQGJ4cpKQkMA///zD//73P/4vhwc9PT2d+Pj4Aj/vVu4tjNam/nIHDx68pfexBrUdVsVo1Ah8EgcSfd2OC56mkxfmWCZMmgR33CH9t+64o3xmqpw/L1t2AC++KMHcJfK+C+RYE6j2SPZpc3E2c0YK3LoFqqja9Omtg9ND+SzadOJyNCSdkXGNGrBokXy7L1x4a61fSgtzRubhw1L1etEi0e5ulYxrcGmVjOu1QUpFw5aLW0jNSqVOtTq0qikR0OZnsG3bilvZo6SYMmUKBoOBJ598skhNRXU6HZqm4eTklOv8t99+i6GQfpK3cm9hHDhwAICaNWve0vtYg/IEVTF0Ogjp6sWyiz2Zf20jUwEijkJmIjh4yISffpJv90OH4Pbbxe3h6VnIO5cRly/DbbfB1avipSqpLwhNg/Nfybi+GyDZLamZqWy8sBGwxAMlJckuHChPkKJ49OoFfDOUX8/NZWFXcLgCXJgBLb6TCf36iadz0iR44w1o2FC8LuWF996TbFKdThSg2rVL5n0vLgFjFngC3pZ2GzmzwszNM5U31kKPHj2YM2cOzzzzDO3bt2fChAm0aNECvV5PVFQUv/32GyA1f8zH3r1789FHH+Hj40NwcDAbN27ku+++w6uQmK5buTcnMTEx7NixA4C0tDQOHDjA9OnT8fLy4mGzkVsWFBo6XUWprNlhmqZps2ZpGj0+0LzeRzP+bcpQOftu7kkXL2panTqSrnDbbZqWmWkbYXNy/bqmdepkyUKJiCi5947dJj+HX9C0jAnZp1ecWqExFS1wZmB2ZsXatRYRFIrikJGhaS6eSRpvOGnr15iewRXVb574/PPyx+bkpGlbt5a9oHnx7beWVKbZs0v2vVe3kZ/FUTtN065mn272RTONqWi/HP4l+1zbtiLC0qVFe+vCMooqAwcOHNAefvhhrV69epqTk5Pm7OysNWzYUBs3bpy2bt26XHMvXbqk3XPPPVr16tU1d3d3bfDgwdqRI0f+v73zDoviaAP47+i9CChW7AVssQU1FDX2FhONJdFYkqhRo2kmtohGE6PxS7XE2BIj1qhRYxewRLC32BV7wY6AgMDN98dwwEmHo8n8nmefW3Zmdt9hb/benXmLcHNzE++8805SvbS8w7LaNj14zivM1NRUVK5cWQwYMEBcvHgx0/aG9A5TSlA6vMhK0NGjQlDypMAPcTZQIx86OyulrnjoULL7+ZAh0hW2oIiNlcqYzo3/7FnDnn//APl/+BchxIGkw0M2DBH4IYZsGJJ0zM9PitGnj2FFUBQvWrUSgrfbivYLEWJZoiL08IB+pfh4Ibp0kV84Z2fpRl+QrF8vhJGRlGfMGMOeO/y8/B/4I8TTzkmHLz64KPBDmEw2EY+iH8mq4cliZPVdqDgoQcUF5SKvyBV16oBtjAeEl2fWPSEP3rkMkaH6FRs2lB5jGg3MnQtTp8r3v/xGq5VpPbZula4T//wDNWoY7vwJMXB1hdyvXB6QHgtCCDac3wBAlxrJnhVqGl5hCLy9gQsd2HINYlwTH8WXvtSvZGyc7LV5/760f7t1K99lBWRo5jfflOOxf3/5PDAklxPTZLgClsmeQbox6FXBCwcLB0DajWu1MjuHIXwiFMUXpQQVQ4yNwesVaZj523lIKJnoeRLql7py166giy0xYQJ89FH+xhBKSJAGmH/+KXNTrF6dc1eQ9LixHuKeghVQajAg/x9H7xzlZsRNrE2taVGpBQBxcclh+r29DSuGonihixckgL/jEsfUlZ1SKU+JtbXMK1auHJw7J0O/nzuXv8IeOiTtA2NipCI2b55hPdaEFi7Pl/uV7IA2SUXrz60H9F9EdkkzPTUGFblGKUHFFN0DOCYeQswTc82E/iVz9jzPRx/B99/L/R9/lAke8yOYVUSETOHxyy/y7wULch+OPy0uz5WflQBN36TDG87JN9DWVVpjYSKT+R06JGMzOjtDovepQpEjXn4ZTCOqwoNqTDiDVMKfxcH1hakrlykDu3dL986rV+U05P79+SPoqlXygfHokRR65crcp6Z5nru74Ok9mTG+XD/kDjyKfsTuq7sB6Fy9c1J1nRLk62tYMRTFD6UEFVO8vYHLLSHBjC/OxspnztOncOfvtBuMGiWjwZqawooV8m3wyZO8E1D3oN+4ESws5LX79TP8daJvw+3EJ2olT3TRaQHWn098A62e/AaqyyHo7W0Yj2BF8cXKCho1Ai504MIjeOCaGLAm9Pu0G1SqJP3CGzeWkd1btoRNm/JOQCFkJPk330yeAdq2LW8C64Tq8vUBJsnhKbZc3EKCSMDdxZ0qJaoAEBUFBw7Ich8fw4uiKF6ox3gxpVEjsDCygSs+7L0F0eVkYEBCM3A579VLPnRtbCAgQGoCJ08aXrh9+6BJE3luV1f52terl+GvA3DlTzkV7wzYJbvk3nxykyO3j6BBQ8fqHZOO695A1cNXYQh0dkEAsx8kPo7vXITIy2k3cHGRY69dO/nS0qWLXK5OTERpMGJi4O23ZfoOkLPB69dDoou1QYmLhOtr5H6lykD9pCKdPVDKWaDgYNnd8uWlXqhQ5AalBBVTzMwSc14lPoA3JpSQBTeOQMz99Bu++qrUBEqWhOPH4aWXZKTYqKjcCxUeLvN/tWgBd+/KKGgHDkiFKC8QAkJ1S2GmwOtJRRvPbwTAs5wnJa1l4K64uGSjaDUNrzAEXl7AVW808VZ8ezwCUSrxkRz6VfqNbGykQtKvn7SZ+/RTOUZ00yO5ZetW+Zbk7y/t8ObNg//9L++y2l//C+JjwRZwTrbJi0uIY/NFGTgxraUwH5/CGUhbUbRQSlAxxssLONsVgFH7byMcAa2AKzMzbEeDBjKn0euvy4fwjBkyIdLGjTkTJC5O2v1UrQrffgvPnslz790rX/fyikdHIDxUjgK314HksLO6pbCUD98jR6SuV6IE1K6dd2Ipig/Nm4MmwQJxvh1RcXDFSUZIJnRl2vZ5OkxNZTqNefPA0RGOHpVvNcOGyZeJnHDihAxE2q4dnDolI1dv2wbvvZez82UVXZqMShrQvJ10eO+1vTyOeYyzlTOe5ZKzFOuWpNWLiMIQKCWoGOPtDTyuhOmDetyK1PKoTDlZcGl+5q7wZcvKJI8bNoCbm7Th6dxZ2itMnw6hoRm3FwIuXoQ5c6QCNWKEdAGuWVOec/XqvE/qcyHR4Lo8YJacqDHqWRQ7Q3cC+h4pyh5IYWgcHKBuXeBsNwCmhGrADHgaBXcyeanQaKSCcvYs9O0rx9Ts2XKN6N13YcsW+YKREVFRsl7//nLmdds2qWB9/LHMXdaiRe47mRHhZ+Dufjn5U8kLSPZ31y2FdazWEWMjOQv19KmyB1IYFpU2oxjTtKmc7Y47+Rr4HmfxYxc+Nr4B4ffhQTA4N8v8JJ06yQflV19J24RDh+T2+edyxqhNG2lHYGEBlpZSe9i/H3bu1M9L5uICkybJh7pJPnwtnz2CK/5yv7ozkPyw3x66ndiEWCo5VMLdJdkFTNkDKfICb284/ltHjIQJC09eY157W4wvRkDoVCjTNfMTlCwpc4wNGABDh0r3+QUL5OboKMNcVK4sx59uHN68KcdgcLC+ovTmm/DNN7J+fnBhtvwsC1i/n3RYCJHkGp9yNjYkRE4UlykDVarkj4iKFxulBBVjrK2lnnLgTDfwnYTfvtN85G2G5uozCJ0MzluyfqJp06Tx5Nq1chYnMFCuHx05kn47U1OpibVvL3Mk5YXRZXqELoaEZ+AAOA8Fku0ddK7xnat3TspTFB8Pe/bIcqUEKQyJlxf8/LMjVvd8iSy5g+PWtWlAMNw4JO3zLJyzdqIWLeC//6Qr/erVsGYNhIXJZbOMcHODVq3kC4inZ8Z1DUlcJFxODAdQzYaUNnln75/l0qNLmBmb0aZKcsyglK7xyh5IYQiUElTM8faGA9/VxTa+EhFc5qaLJ+WuhsigbQ2iwCQbS1KlSsGQIXK7dw/WrZPG0zExEB0tt9hYGWCnVSv59M/rJa+0EFo4/4Pcr64BTfIbqFZo2XhBLkOkXAo7elQmTk1avlAoDISXl/yMPNgNOu5g8slY1jkCjwRc+QFqTsn6yUxMpOu8LgHr3r2webOM8aMbhzExcty1aCHHYeXKBaNRXFkqg5TaAq7vAZZJRbqlsBYVW2Brbpt0XM3GKgyNUoKKOT4+8N13GowvvAa1vmfuNSem2ACR8XB5DlT7NGcndnHJe4PKnHJnB0Rek7GR3DoB5ZKKDtw8wN2ou9iZ2+Hl5pV0XPfw9fLKOycZRfHE1VVmgTl3tit0HMb6c0eJ61kN00cX4PxsqDEZNDkwQjM2lgO8MGoMQsCFxHhI1QDNUL3itJbCYmLkchgoo+i0WLx4MQMGDODgwYM0atQo6fj9+/dp164dZ8+eZe3atbRu3dqg1/Xz82PSpEmIgkipZACUeWcxx9tbPisfB0vDzFmH/0VbLdEj69x0OWvyonH+R/lZCTAdoVekWwprX7U9ZsZmSceVR4oiL2nZEogoS6n4JggEe4yaSyU98hHcSieAaVHm3r/w+Jxcha7si9SEEoui7hF8IxiAzjWSlaD9++VEsqurDJytyJwbN27g5eVFaGgoO3bsMLgC9CKglKBijp1dYtTa682wNXLhccxjztm/IecIn9yD25sLWkTDEnUVbiX2qZob0EqveO3ZtYD+G2hCgrIHUuQtLVsm7pyRLyP/O3QLqiaGbDg7oWCEyksu/Cw/KwJmI/WKNpzfgFZoqVeqHhXskyO4p3wRUfZAmXPhwgWaN29OeHg4u3btwjM/7b1ywdOnT/P1ekoJUsgHsDCm1GNpA/PriRiokrg+f258wQmWF1yYK6fiSwH2H5JyCJy6e4oz989gZmxGp+qdko4fOyYzhNjZSS9ihcLQ6GYYw3a9BsC2S4HEVBogXcfDTsHjPIjMXlBEh8kAiQDVSwKd9IpXnV4FwBu13tA7nhf2QEIIop5FFarNEMtKx44d45VXXsHExIS9e/dSp06dpLILFy7Qp08fSpYsibm5ObVq1WLWrFl67YOCgtBoNPz55598/PHHuLq6YmlpiY+PD0ePHs2SDCtWrKBp06ZYW1tjY2ND27ZtU7Xt378/NjY2nDx5kjZt2mBra0urVq3SOWPeoGyCFLRsKb1iHwV3g7YLWH16A9+/NxDNuVlw+5iM5WFfq6DFzD0JsXApMUJ0dVOgv16x7uHbpkob7C3sk44reyBFXuPsLBXsY8dqUsasJreenWXTLQ9eL2cE17Vwbhy8vL6gxTQMl+bLQJDOgONwUv4MPYx+yI7QHQD08OiRdDw2Vnrzg2GXpJ/GPcXmG5vMK+YjkWMisTbLucPI3r178fPzo3z58mzbto3SpUsnlZ0+fZpmzZpRoUIFZs6ciaurK1u3buXDDz/k/v37TJw4Ue9cY8eOpUGDBsyfP5/w8HD8/Pzw9fXl6NGjVM4gjMLXX3/N+PHjGTBgAOPHj+fZs2fMmDEDLy8vDhw4gHuK7NPPnj2jS5cuDB48mC+++IJ4Q6eAyQSlBClo3lym0XhwsBVWHW24GXGT45EdqF9uNtwQcG4CNFld0GLmnmurIfaxdEIp2xsooVesU4LedH9T77huGl4thSnykpYt5axjyQfduGX7DStOB/C6T1u4vhku/wP17oGFS0GLmTu08XAxcSmsmhHwrl7x32f/Jl4bT52SdajpXDPp+IED0jC6VClpRK5In48++gh7e3sCAgJwcdH/vnz88cfY2tqyd+9e7BJDkrRu3ZrY2FimTZvGhx9+iKOjY1J9FxcX1q5dmxQq5JVXXqFatWp88803/Pbbb2le//r160ycOJHhw4fz008/JR1v3bo11apVY9KkSaxYsSLpeFxcHF9++SUDBgww2P8gOyglSIGlJTRrBkFBFtQ0bs+R+FWsPL2X+rVbwo2dcHkd1HsA5k4FLWrOEQLO/0/uVwWMhusVn7p7itP3TmNmbKbnGq/sgRT5RcuWMkXXvT2vQYdv2HRhE7Fdt2JeYjM81MLFGVB7ekGLmTturIOnYWAOVHgNKK1XvPL0SgB6uPfQO66bjfX2Nqw9kJWpFZFjIg13QgNgZWqVq/ZdunRh/fr1jBo1ij/++APjxOnrmJgYdu7cydChQ7GystKbcenQoQO//PILISEhtG/fPul4nz59khQgADc3N5o1a0ZgYGC619+6dSvx8fH069dP7xoWFhb4+Pik2faNN95IdSy/UEqQApAP4KAgML3YDcqvYu3ZtXzdciE47oRHCXDxe/DIRrySwkbYTnhwRHqjVK0LNNYrTm8p7MQJePxY5qxs0CD/xFUUP3SemjcPNsK1e1nuPL3JjsuP6FijOgSfh/NzoNYUSOG1WKQQAk4lJoatBhjrv4iktxQGeeedqdFocrX0VBiZMGEC9evXZ/LkyWi1Wv7880+MjY158OAB8fHx/Pzzz/z8889ptr1/Xz95tqura6o6rq6uHD9+PN3rh4WFAdC4ceM0y42eyzlkZWWVNCtVECglSAFIJejLL+HC5g6YDjHl7P2znLlvT60alSDkMpz/GWpNBCPTghY1Z/yXqMBVASw/TFWc3lLYDvlMxts7f7J5KIovtrYyGXxwsBHuxq9xh1msObOGjp384FgfiI6Ea0uhUsEsG+SaW5vh0Qn5q1OjGuCrV5zeUlh0tIz5CHmfyuxFYdKkSWg0GiZNmoRWq2Xp0qU4OjpibGxM3759GTZsWJrtKlWqpPf3nTt3UtW5c+cOTk7prwo4O8sI56tXr8bNzS1TWTUF7OqnHusKQD58ra3h4S17vFzasCfsH5b/t4JJ3hPg2ECIfgLXVkDFtzM/WWHj3j64u0s6gtVyBt7SK05vKQxg+3b5qcJrKPKDli2lAbDR2e5QdhZrzq5hdscfMa9uB8efwNnJULF/0fMRf34WyHw00vUtmfSWwvbulYbRZcvK/MqKrOHn54eRkRETJ05ECIG/vz8tWrTg6NGj1K1bFzOzzGcUly1bxscff5ykqFy9epV9+/bRr1+/dNu0bdsWExMTLl26VKDLXFlFKUEKQKbx8vaWEfYrPOkD/IP/f/74+Z5EU20EnIyC0+PBrU/OotcWJKcSZ4EqAdajAQu94vSWwmJiku2BlBKkyA9atoSpU+HUJi/KflyWmxE32XxxB69VGQH/TYVHV+DOdijdJtNzFSru7oL7IfJFpKYr0FevOKOlsJQvIkVN9ytovvzyS4yMjJgwYQJCCGbOnImvry9eXl4MHTqUihUrEhERwcWLF9mwYQMBAQF67e/evUu3bt147733CA8PZ+LEiVhYWDBmzJh0r1mxYkUmT57MuHHjCA0NpV27djg6OhIWFsaBAwewtrZm0qRJed31LKOUIEUSLVtKJejBv12wetmKiw8vcujWSRpX/wDOzIDHV2V8jwo9Mj9ZYeHhUTkNrwHcbYEhqaqktxS2d69UhMqUkenOFIq8plkzMDeH2zeNGVCmF4vOzcT/pD+v1ZwNVafBuQQ4MQpcTxUtjeBUyuXoz5GW0cmktxQGajY2t4wfPx4jIyPGjRuHVqtl//79TJs2jfHjx3P37l0cHByoVq0aHTp0SNX266+/5uDBgwwYMIAnT57QpEkTli9fTpUqVTK85pgxY3B3d+fHH39k2bJlxMbG4urqSuPGjRkyJPUzuEARijQJDw8XgAgPDy9oUfKNw4eFACFsbYXotaq3wA8xavMoIcRdIY6bCrEUITaWFyIhvqBFzTq7u0u59yKE+DJV8X9h/wn8EGZfmYnH0Y/1ykaPlv+Pfv3ySVaFQgjRsqX83n3+w2GBH8JiioUIjwkX4ml/IZYjv883NhS0mFnn3n4psz9CRDoKISJTVWn3ZzuBH2Jy0GS942Fh8n8Bcj83REdHi9OnT4vo6OjcnagYEBgYKACxatWqghYlTTK7l9n5/S5i6xqKvKR+fXB0hIgIaGQu7WaWn1pOgrYE1PwAzIDw63B1WYHKmWXCzyRHpvWwBNI3iH5+KQzUG6iiYNCl0Liw5yVqOtckJj6GdWfXgeUUqJ4YrfPER0Unr5/ecvQngL43VkZLYTt3ys969aBkyTyWU1EsUUqQIgkjo2Tvi6cn2uBk6cSdyDsEXQkCs/FQK9Ez7ORoGfSssHN6GiBkkniHD4DUHg3pLYXduwe6CO+vvpq3YioUKdFlDQgK1NDbow8A/if9gbLgPkAaMTy6CNfXFJiMWebRCbi5IXE52hpI7ZWklsIUBYlSghR66N5CdwWaJnlpyAewM1QfIZfyI2/D5d8LTMYsEXEJriyV+x4mwMepqmTkFaZ7A61TR2atVijyi0aNpLv8w4dQ36Q3ADtCdxAWGQbmk6Bm4mzQyc9k+onCzKmp8rM8YDcCcEhVJT2vMCGUElRQ+Pr6IoSge/fuBS1KnvNCKUEBAQEMHDiQmjVrYm1tTdmyZenatSuHDx8uaNGKDDol6N9/4Y0a8i30rzN/ERMfA6ZjwT3RrfLkGJmLq7By9DMQCTIgrdMgoEyqKr8fl4pc+6rt1VKYotBgYiI9NQEu7K9Kk7JNSBAJibOWZaDmu4lL01fg2vIClDQT7u2Da1LBwcMMGJWqyq2IW2y7tA2AnrV76pWdOwc3bkhDcS+vPJZVUWx5oZSgOXPmcOXKFUaOHMmmTZv48ccfuXv3Lp6enqlc/xRpU7MmlC4tvaK0V5pT3q484bHhbL6wGXCCaiNl7q2n92QixMLInZ1wY62cgn9JA4xOVSVeG8+SE0sA6F+/v16ZegNVFDS6JbHt26FP7ZRLYoDZRKiV6Nh78vPCuTQttHA40QavCuD4PlAqVbUlx5egFVqalW9GdafqemXbpG7EK6/I1D4KRV7wQilBs2bNIiAggKFDh+Lj40P37t3Zvn07Tk5OfP311wUtXpFAo4G2beX+1i1G9K4tp+P9/0t8AJt8Dh6J7q2nvoT4pwUgZQZo4+HwSLlfDXDoB6TOdrzl4hbuRN7BxcqFjtU66pWdPw/Xr8uksro3coUiP9GNwaAg6FT5TYw0RgTfCCb0UShQGqq/L5emI24WzqXp0EXw8DCYAvXMgM9SVRFCsOjYIgAG1E8dBVu9iCjygxdKCSqZhvuAjY0N7u7uXL9+vQAkKpro8udt3gx96si30A3nNvAk9gngBFVGSQeP6Ifw3+SCEjNtLsyF8FPyB6KOFTA1zWq6h+/bdd/G1Fg/FYju4du8OVjlLpehQpEjatWCChVkpORzh0rTspJcp17+X+Lyl+l4cE+cDTr+KcQ+LCBJ0+BZOBwfK/drAxajgQqpqoXcCOHcg3NYmljypoe+Y0JcXHK+MKUEKfKSF0oJSovw8HCOHDmCh4dHhvViY2N58uSJ3lZcad1aeoqdOQP2MXVxd3EnNiGWtWfWygrGo6Fh4vz0mRnSA6QwEPsATn4p9+sC5uOAsqmq3X96nw3nNgDqDVRRONFonnsZSVwSW3pyKUII5GzQB2APxDyGo58UlKip+e8riLkLdkD1MsAXaVbTvYh0d++Onbl+As2QEIiMBGdnGbpDocgrXnglaNiwYURFRTFu3LgM633zzTfY29snbeXLl88nCQsfjo4yci3Ali2apAewzpAYSkC5ydLjQ2jhwIDC4aVy4kt49kg6oFRxIy2PMIClJ5YSp42jYemG1ClVR68sLg4CA+W+UoIUBYlOCdq0CbrVfB1zY3NO3zvNkdtHZIHxBGiSqDyELoawwAKRU48n5+Dcj3K/AWA8nefjAgE8jXuaNKuV0YtIq1byhUyhyCsK7dcrKCgIjUaTpe3YsWNpnmPChAksXbqU77//noYNG2Z4vTFjxhAeHp60Fffls5QP4H71+mGkMSLwSiDn7p9LrDESGtaSa/4PjsCFWQUlquTRCbg4V+43BIxm8nyOMB2Ljy8G0n74Hjggg0WWKAEvvZQ3oioUWaFlS5nTLzQU7l635/VarwPw6+FfE2s4g8sP0vYN4MAgiI8uCFGTOfwRiHjpjFmmKdAnzWprz6wl4lkEFR0q4lPRJ1W5mo1V5BeFVgmqUaMGv/32W5a2ChVSrzdPmjSJKVOmMHXqVIYPH57p9czNzbGzs9PbijM6JSggAEpalE8yHp53eF5iDVOwWgD1E/88/jlEXctvMSXaeDg4VM5KlQdK+QKvp1n12J1jHLtzDDNjM3rX6Z2qPOUbqLFxnkmsUGSKrW2yYf7mzTCkkcy55H/SP9E+D6A/1GsmPTYjLidHZy4Irq2C25vlr0oDgB95PlO8Dt1SWP96/TF6LiHz48fyZQSUEpQVsjpZEBQUxJUrVzKs4+fnl3ReIQTLly/Hy8uLkiVLYmFhQbly5Wjbti3z5+t7Bj948CApX5i1tTX29vbUrFmTvn37cuJEITGXSIdCm0C1dOnSvPvuuzlqO2nSJPz8/PDz82Ps2LEGlqx4UL++DBJ4547MpD6k0RA2nN/A4uOLmdpqKhYmFkBTqPo+XJkH92KkIuKzMf8TO574Eu7vk9/mlzTAD6T78D0qH76v1XyNEpYlUpXr3HLVw1dRGGjfXgbu3LwZPvzQi1rOtThz/wxLTyxlaOOhgAbMfoNGdWFPApz+Ftx6g0Pt/BU0MhT2Jz6vawF27wCN06x69fFVAi7LkCXv1H8nVXlAAGi1UL26NA5XZExwcLDe31999RWBgYGpwsK4u7vz8KE0oB8xYgR9+qSepStXrlzS/pgxY/j222957733+Oyzz7C1teXq1asEBATw999/J/0+R0ZG4unpSWRkJJ999hn16tUjOjqa8+fPs2bNGo4dO0bdunUN3W3DYcikZoWByZMnC0CMHz8+V+cpjglUn6d/f5m48OOPhYhPiBdu37sJ/BBLji9JUeuhEI9LCLEsMbHjlRX5K+TNzfK6SxHiCkKIwelWjY2PFU7fOgn8EJvOb0pVHhYmhEYj+3z9eh7KrFBkkVOn5PfR3FyIqCghfgz5UeCHqDO7jtBqtSlqjhViV+I42NJYiIS4/BMyPkaIzY3ktbciRIK1EOJWutUnBU0S+CFaLG6RZrnuuTNqlGHFLC4JVN955x1hbW2dZtnly5cFIGbMmJHhOZ4+fSrMzc1Fv3SyRyckJCTtL1y4UAAiICAg07qGQiVQTYeZM2fy5Zdf0q5dOzp27EhISIjepsgeKb1TjI2Meb/h+wDMPTQ3RS1HsP8J3BP/3D8QHv+XPwI+vQnBfeV+NcDNFfgq3eobzm3gQfQDytiWoU2VNqnKN22SgRJfeglSvBApFAVGSlf5oCDoW7cvliaWnLx7kpAbKZ9p46BR+UQbvYNwJB+9xY6OhoeHZBTr5oDRVGSo9tRohZbFxxYDadvkJSTAP//I/c6d80Ta5xBAVCHbRJ72OCtERUURGxtL6dJp30ejFNbqDx48AMhS3cJI4ZYum2zYIN2et2zZQtOmTVNtiuzRurW0izlzBq5ehYEvDcTEyIR/r//LybCTKWr2AY8WUBKIj4JdnaSLbF6ijYd9fSD2PjgCDTSAP+CSbpMFRxcA0K9uP4yNUhv8rF8vP/Pn4atQZM7zrvKOlo5J6SXmHk75MmIFVvPAM/HP8z/B+dl5L+D1tfJaAE0B6y7Ah+lWD7oSxOXHl7E1s+UN9zdSlR84IJMX29vnV6qMp4BNIdvyPgCtVqslPj4+1abD2dmZqlWrMnv2bP73v/9x9uzZxNAMqdH9tvbr149169YlKUVFhRdKCQoKCkIIke6myB6OjqDTHTdvBlcbV16r+RqQ0kMFQAPG88DLBmyBqKuw+zVIiMk74U76wd3d0g6oOWA8EWiRbvWz98+y+eJmNGgY1GBQqvKYmGR7oC5dUhUrFAVGSk9NIWBIQ2kgvfLUSh5GpwyS2A7K94N6iX8e/hBubc07wSKvQMhAuV8LKFseWER69ngAP4T8AMBbdd7CyjR1JNLE91jat5eecYq84fPPP8fU1DTVtnfv3qQ6/v7+ODo68sknn1CrVi3s7e3p3LkzS5Ys0fs9bd68OZMnT+b48eN069YNZ2dnKleuzNChQwu9UTS8YEqQwvCkfABD8gP4j+N/EPksMkXNqmC+AnyQ0+L3g+UDMi+Uz9Df4VRiGpQmgF1LYHyGTXQP3y41ulC1RNVU5UFBEBUFZcpAgwYGlVahyBUpXeUvXIAmZZtQ37U+MfEx/HH8j+dqzwX3RlAJmUD43zfh8SnDCxUdBru7QNxjcALqGQHLgdTOBjrOPzjPhvNSyxnlOSrNOjolKP9mY62AyEK25X2Y+pEjR3Lw4MFUW/0UkSkbN27MxYsX2bJlC2PHjqVp06bs3LmTfv360aVLFz1FaMKECVy7do2FCxcyePBgbGxsmDt3Lg0bNmTZsmV53p9cYRgzpRcPZRgtOXJEGilaWwsREyNEgjZBVPupmsAP8dvh39JoMUOI2wjhn2ikecLPsAKd+SHZEPogQoiSQojbGTa5F3VPWEyxEPghdl3ZlWadDz6Q/Xz/fcOKq1AYgpYt5ffzhx/k33MPzhX4IWr8XOM5A2khhLghRLyrENsSx8m6ikJEhxlOmMgrQvxdVZ77L4SIRAjxbabNhm4cKvBDdPLvlGb55cuyj8bGQjx8aDhxdSjD6KwbRqfH/fv3ha+vrwDEP//8k2HdXbt2CSsrK+Hi4pKja2WEMoxW5Bs6V/moKOkqb6QxSsdAWscn4No32Tv2pB8c/hi0cbkTRAg4MRGOjJJ/10QGRWQp4Jph07mH5hITH0PD0g3xqpDa0ECIgngDVSiyTkq7IJA5/WzMbDj34By7ru56rnZZMP4bvMykiUnUFdj+Cjw+Sa4JPwPbmkPkRbA2htaAdXvg0wybPXj6IMkg+mPPtCO568bgK6/IpXhF4cPJyYlRo0YB8N9/GTvAeHt706ZNG+7du8fdu3lsI5oLlBKkyBCNBtq1k/u6B3D/+v0xMzbj8O3DBF8Pfr4FMA+qvgy6jBTnvoedLaQ3V04QWpkZXpestS7wEqD5Gng1w6ax8bH8cuAXAD5u+jGaNGIYnTghs8ZbWsogiQpFYaNDB/kZFARPn4KtuS1v1XkLgJ8P/JxGiyZgsQB8kasrERdg68tweUnOhXhwCHZ4QfRNsDeG1glgWx34g8x+Sn49/CvR8dHUd62Pb0XfNOuoF5HCQ1xcXLoGzmfOnAGgTJkyAISFhaHValPVS0hI4MKFC1hZWeHg4JBnsuYWpQQpMqWjDBbN33/LWRNnK2fervM2AFP2pBWh1gJYC3XKghdgqoF7/8KWBtnPb/ToOAR1hPOJD/pGyMzUmu+BMZk29z/pT1hUGOXsytHDvUeadXQP31dflYqQQlHYqFULKlaUrvI6A/4PX5ZeWGvOrOHU3bTsft4Gu9HQDjlZmhANwf3gwJDsOS3ER8HpGfJFJvYBlDCCVxPA6iVgD+CcYfPY+NgkRe2Tpp+k+SLy5Ely1nilBOU9165dSxVCJiQkhEuXLgEy8bibmxsDBgxg6dKl7N69m02bNjF69GgmTpxIrVq1eP11GZV/yZIl1KhRg4kTJ7Jx40b27NnDsmXLaNOmDadOneLTTz/FzMysILubMYZeq3tRUDZByURECGFhIdfrjx2Txy48uCCMJhkJ/BCHbx1Op+UFIYSHEE8Q4p9E+wR/IyFCBglxJ0CIhPj0L/roPyF2d0+2//HXCBGKEMJICLEoS3JrtVpRe3ZtgR9i+t7p6dZr0kT2bd68LJ1WoSgQPvpIfk/ffjv52Bsr3hD4IXqv7p1OqwQhxAQhEhDiOMnj6Z86QpyfLUT03fQvGB8txJnvhfirVHK7HUZCPEMI4S2EeJwluX8/9rvAD1FmZhkRGx+bZp2VK2XfqlfP0ilzhLIJSrYJSm976623hBBCxMbGiu+++060b99eVKhQQZibmwsLCwtRq1YtMXr0aPHgwYOkc54+fVp88sknolGjRsLFxUWYmJgIR0dH4ePjI5YsWZKmHLnFkDZBGiGU73haPHnyBHt7e8LDw4t9HjGAbt1g3TqYMAEmJ65K9V3blz9P/Em3mt1Y03NNOi0jgQEQvxoOAaEpiixLQ4U3oVRLiAuXsYVi78lM1Df+Ro5LDbjZQe1wsDcDVgCvZUnm7Ze20+bPNlibWnPj4xs4WDikqnPnDuhifN26lbyvUBQ2/v1X2svY28Pdu2BmJnPhvfTrSxhpjDgz7AzVnaqn03ot0A9uRcI+I3iWuHyhMQbXV6FCTzC1leMv5p4cizfWQvQtWc/aFupEQUUtGHUCViITlmWMEIKXfn2J42HH+abVN3zxyhdp1uvXD5YsgU8+ge++y+5/JmvExMRw+fJlKlWqhIVF2smVFUWDzO5ldn6/1XKYIku8kRjX7K+/ko+NfWUsGjSsPbv2ueCJKbEBVoLJN/Ay0AqoYgGm5hB9G879CLu7ymn6o5/K3Ec31gECyttDBwHNw+WTn81kVQEC+F/I/wAY9NKgNBUgSI5O26iRUoAUhZumTeV3NDxc5hMDqO9an07VO6EVWr7e83UGrbsBB6BMdeiolTn2HG2lG/3trTLS+94ecPADODkRLsySCpCVDTQxhs4RUFkLRv2ANWRFAQIIvBLI8bDjWJlaJTlUPE9CQnIIDrUUpshvlBKkyBKdOslYJadPw9mz8lgtl1p0d+8OwNQ9UzNorQG+AM1mKOUKL8fA67HgDbjZg2NJcC0Pbm5Qww3qOUF7wCscHGyAccAloGWW5T0ZdpItF7egQcNIz5Hp1lNRohVFBSMjOSML+i8jE7wnAPDniT8JfRSaRksdtYADYPkG1BLQPgI6AXXMwMkZXMpAuQpQtSJ4lAdPY+gcCVUTwKgFEAT8jszNkTVm7JsByBQZaSUsBggOhgcPpEdY8+ZZPrVCYRCUEqTIEg4O0nAY9B/A471lkMKVp1Zy9v7ZTM7SDrgGbADj7lDOTM7ytL8LLa9D86vQ8Cp4PABHG6Th8xVgCjIiW9b5fMfnAHR3705lx8pp1omOhu3b5b5SghRFAd2M7Lp1oMty0KRsE9pUaUOCSGDa3mmZnMEeWA1cBCaCXSWo8wza3ofWt8D7GjS5AvWuQ+UEMPYBAoEAZCTUrBN4OZAtF7dgYmSSbnBE0I8SbWKSrUsoFLlGKUGKLKN7AK9enXysbqm6dK3RFYHIZDpehyny9XMVcAeYC4xARnyeCSwE1iGVn6/JrvIDsDN0J5svbsbEyISpLdOfodq+XSpC5crJeEgKRWHH2xucnOTMya4U4YF0s0GLjy3mevj1LJypCuCHnGHdg4zz8ykwFZgNLAMOImd/fLMtp1Zo+XS7jB00uOHgNKO0g/Q2XbdO7qsXEUVBoJQgRZbp2lUmVD12TIbw16F7APuf9OfSw0vZOKMjMBj4CZn9/WNgANCVnCg/oP/w/aDRB1RzqpZu3eXL5ecbb8h4SApFYcfEBF57Te6nnJF9pcIr+Fb0JU4bx/R/p2fjjBrgFWBG4jYWGAr0QsajyBnLTi7jyO0j2JrZMtFnYrr1jh2D8+fBwiI5FIdCkZ8oJUiRZZydwSdxRjzlA7hhmYa0r9qeBJHAV7u/KhjhEvnzxJ8cu3MMO3M7JvhMSLdeVJSMewTQu3c+CadQGADdjOzatZAyRp3uZeS3I79xLfxaAUgmiYmPYWzAWADGvDIGF2uXdOvq0kp16gS2tvkhnUKhj1KCFNkiLS8xIOlt7/fjv7Pn6p58lkoSHRfNuIBxgPRcc7ZKP4jbhg0y8m7lytCkSX5JqFDknlatpLPknTuwb1/y8RYVW+Dj5kNsQiwjNo8oMPl+2v8T18KvUc6uXIa2QFpt8mysehFRFBRKCVJki27d5NLR/v1w40by8ZfLvcz7DaQL7OCNg3mW8CzfZfsh5AduPLlBBfsKSdF000P3Btqrl1oKUxQtzMyS7WdSvoxoNBpmd5yNqZEp68+tZ93Zdfku2/2n95M8Rae2nIqlafqu9Pv2yXQ1dnbJaUEUivxGKUGKbFG6NDRrJvfXPBcfcdqr0yhpXZIz988w498Z+SrXvah7fLP3GyDzh++jR8l50NQbqKIoopuRXbNGGhfrcHdx57NmnwEwfNNwImIj8lWur3Z9xZPYJ9R3rc/bdd/OsK7uRaRbN2kTpFAUBEoJUmSb9JbEHC0d+aHtDwB8tfsrLj68mG8yfRn4JRHPImhQugF96vTJsO7atRAXBx4eULt2PgmoUBiQtm3B2hquXYNDh/TLxnuPp7JjZW5G3GRCYPp2cYbm7P2zzD40G4DvWn+HkSb9n5f4eFi1Su736pUf0ikUaaOUIEW2Scybx5490i4hJb1q96J15dbEJsQy9J+h5EdWlvXn1jP38FwAZraZmeHDF5LfQNUskKKoYmmZvISUMmQFgKWpJXM6zgFkhvnDtw7nuTwx8TH0Wt2LeG08Hat1pFXlVhnWDwiAe/eks0WrjKsqssGJEycYNGgQVapUwdLSEktLS6pVq8bgwYM59Ly2nEWCgoLQaDQE6TLcGgCNRqO32dnZ0axZM5bpHs75iFKCFNnGzQ1efllOw/v765dpNBrmdJyDhYkFO0J3sOy/vP1SXwu/Rv91/QH4yPMjfCv6Zlj/zh35AAb1Bqoo2rz5pvxculSmnkhJmypt6F27N1qh5f2N7xOvjc9TWT7Z+gnHw47jYuXCb51/y7S+7reuRw8ZiV6Re3799VcaNmzI/v37GTlyJBs3buSff/5h1KhRnDp1isaNGydlic8ODRo0IDg4mAYNGhhU3u7duxMcHMy+ffuYO3cuT548oU+fPvg//6OS1xgys+uLhMoinzFz58qszx4eQmi1qcun7p4q8EO4THcRtyNu54kMz+KfiWYLmgn8EI3mNUo3Q3VKfvpJyt2kSZ6IpFDkGzExQpQoIb/PmzenLr8dcVs4THMQ+CG+3fttnsmx6tQqgR8CP8TWi1szrR8dLYSdnZR79+48EyuN6764WeT37t0rjIyMROfOnUVsbNrPwZUrV4qbN2/ms2RpA4hhw4bpHbty5YoAhLe3d6btDZlFXs0EKXJEr17SmPHUqdQ2CQCfNvsUDxcP7j29Ryf/TkQ+izS4DBODJrLv+j7szO1Y0X0FZsZmmbZRS2GKFwVzc3g70fZ40aLU5a42rkx/VQZO/GLHF3niLRb6KJRB6wcBMiZQmyptMm2zeTM8eSIjtatcYYbh66+/xtjYmF9//RUzs7Sfgz169KBMmTJJfx86dIhevXpRsWJFLC0tqVixIr179+bq1at67dJaDstq2+zg5uaGi4sLYWFhOT5HTlBKkCJH2NsnG0gvXJi63MzYjHW91uFi5cLh24fpsaoHcQlxBrv+1otbk7zB5neen25+sJRcuSKTNWo0yUsJCkVRZuBA+blunUyl8TzvNniX9xu8j0DQ+6/e7Lu+L3WlHPIs4Rk9V/fkSewTmpdvzuQWk7PULmV4CqOC/gUSAuKjCteWTTvKhIQEAgMDadSoEaVLl85yuytXrlCjRg1++OEHtm7dyrfffsvt27dp3Lgx9+/fz7O26REeHs7Dhw+pXr16jtrnFJWuTpFjBg6U9gj+/jBzJlhZ6ZdXLVGVjX020uL3Fmy5uIXBGwezoMsCNLkMzHPx4UX6ru0LwNBGQ+nh0SNL7VaskJ++vpDihUihKLLUqwcvvQRHj8pxOOK5GIkajYZZHWdxK/IWG89vpPOyzuwbuI8azjVydV2t0PLh5g85dOsQJSxLsOyNZZgYZf5zEhEBGzfK/UIxG5vwFFbaFLQU+rwZCSbWWa5+//59oqOjcXNzS1WWkJCg55xibGyc9Pzt3r073bt316vbqVMnSpUqhb+/Px9+mH6stdy01SGEID4+HiEEV65c4dNPP8XKyoqJE9NPs5IXFLQerijC+PpCxYpyanvt2rTrNCnbhBXdV2CkMWLRsUX4Bfnl6pp7r+3Fc74n957eo16pevyv7f+y1E4I+OMPua8MohUvErrZoLRmZAFMjExY/sZympRtwsPoh7Rb2o47kXfSrpwFouOi6bm6J78e/hWARV0XUd6+fJba/vWXTFpcrZpU3hR5S8OGDTE1NU3aZs6cmVQWGRnJ559/TtWqVTExMcHExAQbGxuioqI4c+ZMhufNTVsds2fPxtTUFDMzM6pXr87mzZtZtmwZDRs2zFWfs4uaCVLkGCMjGDAAJk6UD+C33kq7XqfqnZjTcQ6DNw5m8u7JWJhYMLr5aIyNjLN1veX/Leedde/wLOEZjcs0Zn3v9ViYZC3K2q5dcPq0jK3Ss2e2LqtQFGr69IFPPpHJSI8eTVu5sDazZkPvDTRb0IxLjy7Rfml71ry5hkqOlbJ1rbtRd+m6vCshN0IwNTJlYdeFdKnRJcvtZ82Sn4MGFZJI7cZWcualMGFslXmdFDg7O2NpaZmmPY6/vz9Pnz7l9u3bdOmif5/69OnDzp07mTBhAo0bN8bOzg6NRkOHDh2Ijo7O8Jq5aavjzTff5LPPPiMuLo6TJ08yZswYevXqxZEjR6hWLf3E1wYnC4bcxRLlHZY1rlwRQqORnh6hoRnXHb9zfJIXSZPfmojjd45n6RparTbJ2ww/xGvLXxNRz6KyJWf37lLGwYOz1UyhKBL07Cm/38OHZ1zvwoMLwnm6s8APYTXVSszcN1PEJcRl6Rpn7p0RlX6oJPBDOE5zFEGXg7Il4/79UkZzcyHu3ctWU4PwInuHdezYUZiamopbt26lWX758mUBiBkzZgghhHj8+LHQaDTCz89Pr15MTIwwNjYW77zzTtKxwMBAAYjAwMBst00P0vAO2717t9BoNKJjx46ZtlfeYYpCg5tbcrCz33/PuO7kFpP5tdOv2JnbceDmARrOa8i4neOIiY9Js350XDQr/ltBmz/bJCVG/cjzI1b3WI2Vadbflm7eTF6uGzYsy80UiiKDbkls6VKISXs4AdJOL3hQML4VfXka95RPtn2C53xPjt05lm6bk2EnGb19NJ7zPbn8+DKVHSsTPCgYn4o+2ZJRNwvUs6cMkqgwHGPGjCEhIYEhQ4YQF5e5A4pGo0EIgbm5ud7x+fPnk/B80CkDts0ILy8v+vXrxz///ENwcHCOz5Nd1HKYItcMHAg7dkg33S+/TN/jQ6PR8H7D9+lYrSMjNo9g7dm1fL33axYfX0x91/pUcaxCFccqlLEtw7ZL21h5eiVPYp8AYKQx4sd2PzK8yfBsy/fbbzKYnJcX1KmTm54qFIWTVq2gfHmZkHT9+oy9H6uWqEpAvwAWHl3Ip9s/5fDtwzSa14iXy72cNAarlKjCvah7/HHiDz0FybOcJ+t7rcfF2iVb8t2/n+yY8MEHOeigIkOaN2/OrFmzGDFiBA0aNOD999/Hw8MDIyMjbt++zV+JOY7s7OySPr29vZkxYwbOzs5UrFiRXbt2sWDBAhwcHDK8Vm7aZsZXX33FihUrmDBhAjt27MjVubKKUoIUuea118DBQeYxCgiAV1/NuH5Zu7Ks6bmGNWfWMHzTcG5F3OJWxK0067rZu9G3bl/61etHNafsrxPHxcG8eXJfPXwVLyrGxtC/P3z1lbTPyywEhEajYVCDQXSs3pEPN3/IqtOr2Hd9X5ou9KZGpnSq3ol+9frRsVpHTI2zH+J5wQKIjYWGDaFJk2w3V2SBIUOG0LRpU3788Ue+//57bt26hUajoVy5cjRr1oydO3fSsmXLpPr+/v6MHDmS0aNHEx8fT/Pmzdm+fTsdO3bM9Fq5aZsR5cuXZ8SIEcyYMYPdu3fj7e2dq/NlBU3i+pziOZ48eYK9vT3h4eFJ2rMifT74AObMkWHwV67MeruoZ1GE3Ajh0qNLXHp4iUuPLnEt/Bp1StahX71+eLl5ZZoLLCNWrpTT76VKSSUtnThiCkWRJzQUqlSRBseXLkGlbNg8n7p7ilP3TiWNwUuPZHqFHu496OnREycrpxzLlZAg5bp6VSpoAwbk+FS5IiYmhsuXL1OpUiUsVNr6Ik1m9zI7v99qJkhhEAYPlkrQX3/BxYtQtWrW2lmbWdOqcitakTdZFHV2CO+/rxQgxYtN5crQpg1s2ybjdv3yS9bbepT0wKOkR57ItXmzVIAcHVV4CkXhQxlGKwxCvXrQvj1otTBjRkFLI/nvP9i9Wy4VDB5c0NIoFHnPF1/IzwULIJ+zD6SL7kVk4ECwtCxYWRSK51FKkMJgjB0rPxcvhltpm/jkK7Nny8/XXoOyZQtUFIUiX/D1hZdflh5iP/xQ0NLIWeEtW+QS3dChBS2NQpEapQQpDMYrr8jt2TP4X9YCOecZT57AkiVyX7nFK4oLGg2MGSP3Z8+G8PCClWfOHPnZrp20C1IoChtKCVIYFN0DeO5cePiw4OT46SeIjIRateTbsUJRXOjcGTw85IuAbja0ILh7V3lmKgo/SglSGJT27aV9UFQU/Pxzwchw/36yXdKECYUkPL9CkU8YGSXbBn3/PTx9WjByfP21fBFp2BA6dCgYGRSKzFBKkMKgpJyO183G5DfffCPfguvXV3nCFMWTXr1kcuN799JPrJqXXLmSvBQ2bVr6AVQVioJGfTUVBqd7d+ki//ChjNacn1y7luwa/M036uGrKJ6YmMBnn8n9GTNk0ND85MsvpW1gq1aZB09VKAoS9ROhMDjGxvD553L/u+9kpNj8ws9PPnx9faFt2/y7rkJR2BgwIDlIqL9//l335En480+5/803+XddhSInKCVIkSf07Svd0m/dgunT8+eap08nJ3GdNk3ZAimKN5aW8NFHcn/cOLlEnB+MHQtCyBnhxo3z55oKRU55oZSgY8eO0bFjRypUqIClpSUlSpSgadOm/Kl7LVHkG+bmycbJU6bA2bN5f81x42Swxm7dZKwUhaK48+GH0jX95s3kOF55yd69sHGjnA2eMiXvr6dQ5JYXSgl6/Pgx5cuX5+uvv2bTpk388ccfVKxYkb59+zJFjch8p1cv6S327JlMW6HV5t21goNh3TppAzR1at5dR6EoSlhawq+/yv3Zs2Ff6vyoBkOIZK+0gQOhRo28u5ZCn8WLF6PRaDh06BDfffcdGo2G/fv369XRarWUKFECjUbDuXPn9MqePXuGlZUVr7/+etKxBw8eMGbMGNzd3bG2tsbe3p6aNWvSt29fTpw4AchEvFnZgoKC8vx/kFNeqNxhvr6++D4XFKZTp05cvnyZefPmMX78+IIRrJii0UgPEQ8P2LNHGknnRfqKZ89g5Ei5P2CAjA2kUCgkrVrJDPOLF8N778GRI3Km1tD4+8O//4KFBUycaPjzK7JGixYtAAgMDOTlFFPix48f59GjR1hbWxMYGEiNFFrq/v37iY6OTmobGRmJp6cnkZGRfPbZZ9SrV4/o6GjOnz/PmjVrOHbsGHXr1iU4OFjv2l999RWBgYEEBAToHXd3d8+r7uaaF0oJSg9nZ2fu3r1b0GIUS9zc5LT4Rx/B6NEykFuZMoa9xhdfwMGD4OAgDaMVCoU+330HmzZJu7lvv5XeW4bk7NnkF5wvvlBpagqSl156CQcHB4KCgvhCNzUHBAUFUaZMGXx8fAgMDGTIkCF6ZZCsQK1atYqLFy8SEBCQdEzHxx9/jDZxWt/T01OvzMXFBSMjo1THCzMvpBKk1WrRarU8evSIVatWsXXrVn7JJKVybGwssSncmJ7klxVhMWDECPmWePCg3P/rL8Ode+1aGRAOpFF0uXKGO7dC8aLg5AQ//gi9e8vl4h49DDdj+vSpPF9UlPTKLFIT7kIUXDTJ9LCyypVXh5GREd7e3gQEBBAfH4+JifyZDwoKwtfXF29vb758TgsOCgrCxcUFDw8PQC6FAZQuXTrda7wovDg9ScEHH3yAqakpJUuW5KOPPuKnn35icCbrMN988w329vZJW/ny5fNJ2hcfY2O5FGZsDGvWyM0QhIbK5S+ATz6BLl0Mc16F4kWkZ08ZuVlno5eQYJjzjhgB//0n3fH9/eU4LzI8fQo2NoVrM4BS1qJFCyIjIzl48CAgJwZ2796Nj48PPj4+3L17l9OnTwPSHig4OBhfX180icpX06ZNAejXrx/r1q1LUopeRAqtEhQUFJRlo6tjx47ptR07diwHDx7kn3/+YeDAgQwfPpzvvvsuw+uNGTOG8PDwpO369et52LviR716ycHb+vWTNkK5ISZGvn2Gh0PTpioeiUKRGTobPWtr6cXVv3/uFaHFi2VEaiMjWLYM0pk4UOQzuiUs3TLXsWPHePz4MT4+PtSsWZNSpUoRGBgIQEhIiJ49EEDz5s2ZPHkyx48fp1u3bjg7O1O5cmWGDh2aZBT9olBol8Nq1KjBb1kMN1yhQoVUf+uOdUhMWjNmzBjeeecdXFxc0jyHubk55nlhLahIYuJEOHwYtm+XXmNbtsis8znhk0+kgaeTE6xYAaamhpVVoXgRqVABliyRLxB//ikVo0WLcjZ7899/yYlR/fzgOdORooGVVcHk9skIK6tcn6Ju3bo4OTkRFBTEmDFjCAoKwtXVNckY2tvbm8DAQIYNG5bKHkjHhAkTeP/999m0aRP79+9n3759zJ07l/nz5/PHH3/Qu3fvXMtZKBDFgIULFwpAhISEZLlNeHi4AER4eHgeSlb8ePpUiNathQAhbGyE2Ls3e+1jYoQYMUK2ByE2b84bORWKF5nVq4UwNpZj6J13hIiPz177oCAhSpeW7du0ESIhIU/ENCjR0dHi9OnTIjo6uqBFMTiLFi0SgDh48GDSsddff11YW1uLZ8+eiU6dOomePXsmlf3yyy/C2dlZaLVa4evrK1xdXbN0nV27dgkrKyvh4uKSZvk777wjrK2tc9eZLJDZvczO73ehXQ4zJIGBgRgZGVG5cuWCFqXYY2kJf/8t8wlFRkK7dtKtNiuEhkLz5snZ6adNk+0VCkX2eOMNWL5czgD9/ju8+27Wlsa0WmlY3bIl3L4N7u5yZukFspN9YWjRogVRUVHs37+fPXv24OPjk1Tm4+PD/fv3OXz4MCEhIalmgdLD29ubNm3acO/evRfG47rQLoflhPfffx87OzuaNGlCqVKluH//PqtWrWLFihV89tln6S6FKfIXnSLUuTMEBEDr1jBoEAwbBjVrpt1m9WpZ58kTuQT2xx/SyFOhUOSM7t2lHU/v3tK25+xZGW/r9dfBzCx1/bAwmQ5n+3b59zvvwKxZ0sZIUfjQKTbff/894eHhejH0PDw8cHJy4ptvviEmJiaVEhQWFpbk7p6ShIQELly4gJWVFQ4ODnndhXzhhVKCmjZtyqJFi/j99995/PgxNjY21KtXjyVLlvD2228XtHiKFFhZwYYNMsXFtm0y8/svv0iF6IMPpJPE6dNyO3kyOdJts2byDVY57ykUuadHD7mw3K8fhITIrXRpGfOnbVu4dCl5HO7ZAw8eyJeY2bOlYbWicKFJ4Vrv4eFByZIlWbt2LS4uLtRKERNBo9Hg7e3N2rVrgdT2QEuWLOHXX3+lT58+NG7cGHt7e27cuMH8+fM5deoUX375JWZpacpFkBdKCRowYAADdD7TikKPlZU0jt65UypA69fLt0zdm+bzfP45fPWVMoJWKAzJm2+ClxfMmwdz58plLj+/tAOP1qoFq1bJKPCKwsPTRLf65517fH19Wblypd5SmA4fHx/Wrl1L2bJlqVq1ql5Zx44duXPnDps2bWLOnDk8evQIW1tb6tat+8JNKmiEEKKghSiMPHnyBHt7e8LDw7GzsytocYoFV67IN0x/fznF7u6evDVqpHIRKRR5zbNnMo7XrFlw4QJUr64/Dr288iblRn4QExPD5cuXqVSpEhYWFgUtjkEZOXIkv/zyC48fP8bW1ragxclzMruX2fn9fqFmghRFm4oVYfp0uSkUivzHzEwmPu7Vq6AlUWSFw4cPc/DgQRYuXEiXLl2KhQJkaJQSpFAoFApFEaR79+6Eh4fTpUsXfvrpp4IWp0iilCCFQqFQKIogly9fLmgRijwquoNCoVAoFIpiiVKCFAqFQlFsUL5ARR9D3kOlBCkUCoXihcc0MbbGUwNkaVcULFFRUWg0mqR7mhuUTZBCoVAoXniMjY1xcHBISvdgZWWlF1xQUbgRQhAfH8+TJ0948uQJDg4OGOck8+9zKCVIoVAoFMUCV1dXgBcm71VxxNjYmNKlS2Nvb2+Q8yklSKFQKBTFAo1GQ+nSpSlZsiRxcXEFLY4im5iYmGBsbGzQGTylBCkUCoWiWGFsbGyQpRRF0UcZRisUCoVCoSiWKCVIoVAoFApFsUQpQQqFQqFQKIolSglSKBQKhUJRLFFKkEKhUCgUimKJ8g5LB11Y7idPnhSwJAqFQqFQKLKK7nc7K+k1lBKUDhEREQCUL1++gCVRKBQKhUKRXSIiIjINqqgRKptcmmi1Wm7duoWtra3BQ6s/efKE8uXLc/36dezs7Ax67sKA6l/R50Xvo+pf0edF76PqX84RQhAREUGZMmUwMsrY6kfNBKWDkZER5cqVy9Nr2NnZvZBfbh2qf0WfF72Pqn9Fnxe9j6p/OSOraTWUYbRCoVAoFIpiiVKCFAqFQqFQFEuUElQAmJubM3HiRMzNzQtalDxB9a/o86L3UfWv6POi91H1L39QhtEKhUKhUCiKJWomSKFQKBQKRbFEKUEKhUKhUCiKJUoJUigUCoVCUSxRSlAOiIyMZNSoUZQpUwYLCwvq16/P8uXLs9T27t279O/fH2dnZ6ysrGjatCk7d+5Ms+6OHTto2rQpVlZWODs7079/f+7evWvIrqRJTvu3Zs0aevfuTdWqVbG0tKRixYq89dZbXLhwIVVdX19fNBpNqq1du3Z50aVU5LSPixcvTlNujUbDnTt3UtUvavcwvfuSVh8L8h5GREQwevRo2rRpg4uLCxqNBj8/vyy3LwrjMDd9LApjMTf9KwrjMDf9KwrjMCAggIEDB1KzZk2sra0pW7YsXbt25fDhw1lqX1jGoAqWmANef/11Dh48yLRp06hevTr+/v707t0brVZLnz590m0XGxtLq1atePz4MT/++CMlS5Zk1qxZtGvXjh07duDj45NUd9euXbRv356OHTvy999/c/fuXT7//HNatWrFoUOH8tSiPqf9+/bbb3F1dWXcuHFUrlyZ69ev8/XXX9OgQQNCQkLw8PDQq1+5cmWWLl2qd8zBwSEvupSKnPZRx6JFi6hZs6beMScnJ72/i+I9nD17dqp8eU+fPqVdu3Y0bNgQV1dXvbKCuocPHjxg3rx51KtXj9dee4358+dnuW1RGYe56WNRGIu56Z+OwjwOc9O/ojAO58yZw4MHDxg5ciTu7u7cu3ePmTNn4unpydatW2nZsmW6bQvVGBSKbPHPP/8IQPj7++sdb926tShTpoyIj49Pt+2sWbMEIPbt25d0LC4uTri7u4smTZro1W3cuLFwd3cXcXFxScf+/fdfAYjZs2cbqDepyU3/wsLCUh27efOmMDU1FYMGDdI77uPjIzw8PAwjdDbJTR8XLVokAHHw4MFMr1MU72FaLF68WABi/vz5escL8h5qtVqh1WqFEELcu3dPAGLixIlZalsUxqEQuetjURiLuelfURiHuelfWhS2cZjWdywiIkKUKlVKtGrVKsO2hWkMquWwbLJ27VpsbGzo0aOH3vEBAwZw69Yt9u/fn2HbGjVq0LRp06RjJiYmvP322xw4cICbN28CcPPmTQ4ePEjfvn0xMUmerGvWrBnVq1dn7dq1Bu6Vvow57V/JkiVTHStTpgzlypXj+vXrBpc1p+Smj1mlqN7DtFiwYAE2Njb07NnTkGLmCt2Uf04oCuMQctfHojAWc9O/rFKQ99DQ/Sts4zCt75iNjQ3u7u6ZfscK0xhUSlA2+e+//6hVq5beDQGoW7duUnlGbXX10mp76tQpvXOkVzeja+SW3PQvLUJDQ7l69Wqq6XeAS5cuUaJECUxMTKhSpQrjxo0jOjo658JnEUP0sVOnThgbG1OiRAlef/31VG1elHt44cIF9uzZQ69evbCxsUlVXlD3MDcUhXGYFxTGsZhbCvM4NCRFZRyGh4dz5MiRNL9jKSlMY1DZBGWTBw8eULly5VTHS5QokVSeUVtdvYza6j7Tq5vRNXJLbvr3PPHx8QwaNAgbGxs++ugjvbJXXnmFnj17UrNmTaKjo9m8eTPTp09n7969BAYGZpr5Nzfkpo86OwtPT0/s7Ow4efIk06ZNw9PTk3///Zd69erpnaOo38MFCxYAMGjQoFRlBXkPc0NRGIeGprCOxZxSFMahISkq43DYsGFERUUxbty4DOsVpjGolKAckNEUZ2bTn9lpm17dvJ5Czk3/dAghGDRoEHv27OGvv/6ifPnyeuVTpkzR+7tDhw5UrFiRTz/9lL///ptu3bplX/BskNM+tmvXTs/rwtvbm44dO1KnTh2+/PJL/v777yydqyjcw/j4eH7//Xc8PDzw9PRMVV7Q9zA3FIVxaCgK+1jMCUVlHBqCojIOJ0yYwNKlS/n5559p2LBhpvULyxgsfCp+IcfJySlN7fPhw4dA2hprdtvqvBvSq5vRNXJLbvqnQwjBu+++y59//snixYvp2rVrlq799ttvAxASEpINibOPIfqYkooVK/LKK6/oyV3U7yHApk2buHPnDu+++26Wr51f9zA3FIVxaCgK+1g0JIVtHBqKojAOJ02axJQpU5g6dSrDhw/PtH5hGoNKCcomderU4cyZM8THx+sdP3nyJAC1a9fOsK2uXkZtdZ/p1c3oGrklN/2D5IfuokWLmD9/ftJgzA55PX2b2z6mhRBCT+6ifA91LFiwADMzM/r27ZttGQrjEoqOojAODUFRGIuGpjCNQ0NR2MfhpEmT8PPzw8/Pj7Fjx2apTaEag7n2LytmbNq0SQBi+fLlesfbtWuXqfvx7NmzBSBCQkKSjsXFxQkPDw/x8ssv69Vt0qSJqF27tt75goODBSDmzJljoN6kJjf902q1YtCgQUKj0Yh58+Zl+9rffvutAMS6deuy3TY75KaPaREaGipsbGzEa6+9pne8KN5DHbdv3xYmJibizTffzNa18+sepiS77sdFYRw+T3b7WFTGog5DuJAXtnGYkpz2r7CPw8mTJwtAjB8/PlvtCtMYVEpQDmjdurVwdHQU8+bNEwEBAeK9994TgPjzzz+T6gwcOFAYGxuLK1euJB2LiYkRHh4eonz58mLp0qVi+/btolu3bsLExEQEBQXpXSMwMFCYmJiIbt26ie3bt4ulS5eK8uXLi9q1a4uYmJhC2b/hw4cLQAwcOFAEBwfrbUeOHEmqt3v3btG2bVsxd+5csW3bNrF+/XoxdOhQYWxsLFq2bCkSEhLytH+56WOrVq3EpEmTxNq1a8XOnTvFDz/8IMqUKSNsbW3FyZMn9a5RFO+hjmnTpglAbNu2Lc3zF4Z7uGnTJrFq1SqxcOFCAYgePXqIVatWiVWrVomoqCghRNEeh7npY1EZizntX1EZhzntn47CPA6/++47AYh27dql+o4FBwcn1SvsY1ApQTkgIiJCfPjhh8LV1VWYmZmJunXrimXLlunVeeeddwQgLl++rHf8zp07ol+/fqJEiRLCwsJCeHp6iu3bt6d5nW3btglPT09hYWEhSpQoIfr165dmgCpDk9P+ubm5CSDNzc3NLanehQsXRIcOHUTZsmWFubm5sLCwEHXq1BFTp07Nlx+W3PRx1KhRwt3dXdja2goTExNRpkwZ8fbbb4tz586leZ2idg91VK9eXVSsWDEp2NvzFIZ7mNH3TdenojwOhch5H4vKWMxp/4rKOMzNd1SIwj0OfXx80u1bykWmwj4GNUIIkZ3lM4VCoVAoFIoXgaJl9aZQKBQKhUJhIJQSpFAoFAqFoliilCCFQqFQKBTFEqUEKRQKhUKhKJYoJUihUCgUCkWxRClBCoVCoVAoiiVKCVIoFAqFQlEsUUqQQqFQKBSKYolSghQKRRKLFy9Go9EkbSYmJpQrV44BAwZw8+bNVPUOHTpUgNLmjD/++AMXFxciIiIKTAZfX198fX3z5NyPHj3CwcGBdevW5cn5FYoXCaUEKRSKVCxatIjg4GC2b9/Oe++9x7Jly/Dy8iIqKqqgRcsVT58+ZezYsXz++efY2toWtDh5gqOjIx999BGfffYZz549K2hxFIpCjVKCFApFKmrXro2npyctWrRg4sSJjB49msuXLxf52YXff/+dBw8e8O677xa0KHnKkCFDuHLlCqtXry5oURSKQo1SghQKRaZ4enoCcPXqVb3jERERDB06FGdnZ5ycnHj99de5deuWXp0VK1bQpk0bSpcujaWlJbVq1eKLL75INasUGhpKr169KFOmDObm5pQqVYpWrVpx7NixVOdr2rQp1tbW2NjY0LZtW44ePZqlfsyZM4fOnTvj4OCgd1yr1fLzzz9Tv359LC0tcXBwwNPTk/Xr1wMwaNAgSpQowdOnT1Ods2XLlnh4eGT5XOnx7NkzpkyZQs2aNTE3N8fFxYUBAwZw7949vXoBAQH4+vri5OSEpaUlFSpU4I033tCTrVSpUrRu3Zq5c+dm6f+iUBRXlBKkUCgy5eLFiwC4uLjoHX/33XcxNTXF39+f6dOnExQUxNtvv61X58KFC3To0IEFCxawZcsWRo0axcqVK+ncubNevQ4dOnD48GGmT5/O9u3bmTNnDi+99BKPHz9OqvP111/Tu3dv3N3dWblyJUuWLCEiIgIvLy9Onz6dYR9u3LjByZMnadGiRaqy/v37M3LkSBo3bsyKFStYvnw5Xbp04cqVKwCMHDmSR48e4e/vr9fu9OnTBAYGMmzYsCyfKy20Wi1du3Zl2rRp9OnTh3/++Ydp06axfft2fH19iY6OBuDKlSt07NgRMzMzFi5cyJYtW5g2bRrW1taplr58fX35999/9f5/CoXiOQySi16hULwQLFq0SAAiJCRExMXFiYiICLFx40bh4uIibG1txZ07d/TqffDBB3rtp0+fLgBx+/btNM+v1WpFXFyc2LVrlwDE8ePHhRBC3L9/XwDihx9+SFe2a9euCRMTEzFixAi94xEREcLV1VW8+eabGfZtxYoVSX1Lye7duwUgxo0bl2F7Hx8fUb9+fb1jQ4cOFXZ2diIiIiLb5/Lx8Un6e9myZQIQf/31l169gwcPCkDMnj1bCCHE6tWrBSCOHTuW4fmFEGL79u0CEJs3b860rkJRXFEzQQqFIhWenp6Ymppia2tLp06dcHV1ZfPmzZQqVUqvXpcuXfT+rlu3LqC/bBYaGkqfPn1wdXXF2NgYU1NTfHx8ADhz5gwAJUqUoEqVKsyYMYP//e9/HD16FK1Wq3furVu3Eh8fT79+/YiPj0/aLCws8PHxISgoKMM+6ZbpSpYsqXd88+bNAHqzOWkxcuRIjh07xr///gvAkydPWLJkCe+88w42NjbZOtfzbNy4EQcHBzp37qzXt/r16+Pq6prUt/r162NmZsb777/P77//TmhoaLrn1PUzpVefQqHQRylBCoUiFX/88QcHDx7k6NGj3Lp1ixMnTtC8efNU9ZycnPT+Njc3B0havomMjMTLy4v9+/czZcoUgoKCOHjwIGvWrNGrp9Fo2LlzJ23btmX69Ok0aNAAFxcXPvzwwyRX9rCwMAAaN26Mqamp3rZixQru37+fYZ9017KwsNA7fu/ePYyNjXF1dc2wfdeuXalYsSKzZs0CZJiAqKgoPYUnq+d6nrCwMB4/foyZmVmqvt25cyepb1WqVGHHjh2ULFmSYcOGUaVKFapUqcKPP/6Y6py6fur6rVAoUmNS0AIoFIrCR61atWjUqFGuzxMQEMCtW7cICgpKmv0B0rRTcXNzY8GCBQCcP3+elStX4ufnx7Nnz5g7dy7Ozs4ArF69Gjc3t2zLomv/8OFDSpcunXTcxcWFhIQE7ty5o3f8eYyMjBg2bBhjx45l5syZzJ49m1atWlGjRo1snyst2ZycnNiyZUua5Snd+b28vPDy8iIhIYFDhw7x888/M2rUKEqVKkWvXr2S6j18+FCv3wqFIjVqJkihUOQZGo0GSJ4h0vHrr79m2K569eqMHz+eOnXqcOTIEQDatm2LiYkJly5dolGjRmluGVGzZk0ALl26pHe8ffv2gPQcy4x3330XMzMz3nrrLc6dO8fw4cNzfK6UdOrUiQcPHpCQkJBmv1IqWjqMjY15+eWXk2amdP8nHbqlMnd392zJolAUJ9RMkEKhyDOaNWuGo6MjQ4YMYeLEiZiamrJ06VKOHz+uV+/EiRMMHz6cHj16UK1aNczMzAgICODEiRN88cUXAFSsWJHJkyczbtw4QkNDadeuHY6OjoSFhXHgwAGsra2ZNGlSurK8/PLLWFpaEhISomfL5OXlRd++fZkyZQphYWF06tQJc3Nzjh49ipWVFSNGjEiq6+DgQL9+/ZgzZw5ubm6pPNyyc66U9OrVi6VLl9KhQwdGjhxJkyZNMDU15caNGwQGBtK1a1e6devG3LlzCQgIoGPHjlSoUIGYmBgWLlwIwKuvvqp3zpCQEJycnKhTp04W7pRCUUwpaMtshUJReNB5fR08eDBH9QIDAwUgAgMDk47t27dPNG3aVFhZWQkXFxfx7rvviiNHjghALFq0SAghRFhYmOjfv7+oWbOmsLa2FjY2NqJu3bri+++/F/Hx8XrXWLdunWjRooWws7MT5ubmws3NTXTv3l3s2LEj0/717dtXuLu7pzqekJAgvv/+e1G7dm1hZmYm7O3tRdOmTcWGDRtS1Q0KChKAmDZtWprXyMq5nvcOE0KIuLg48d1334l69eoJCwsLYWNjI2rWrCkGDx4sLly4IIQQIjg4WHTr1k24ubkJc3Nz4eTkJHx8fMT69ev1zqXVaoWbm1sqTzqFQqGPRgghClQLUygUinzi0KFDNG7cmJCQEF5++eUcneOTTz5hzpw5XL9+PZVheGFh586dtGnThlOnTiUtAyoUitQoJUihUBQrevbsSVRUFBs3bsxWu5CQEM6fP8/gwYMZPHgwP/zwQ94IaABatGhB1apV+e233wpaFIWiUKNsghQKRbFi5syZLFiwgIiIiGwlUW3atClWVlZ06tSJKVOm5KGEuePRo0f4+PjwwQcfFLQoCkWhR80EKRQKhUKhKJYoF3mFQqFQKBTFEqUEKRQKhUKhKJYoJUihUCgUCkWxRClBCoVCoVAoiiVKCVIoFAqFQlEsUUqQQqFQKBSKYolSghQKhUKhUBRLlBKkUCgUCoWiWPJ/SmbQfzfPwTcAAAAASUVORK5CYII=", "text/plain": [ "
" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "# Plot light curves\n", "\n", "plt.figure()\n", "\n", "plt.plot(phi/(2*np.pi), dF_F_Gaia_B, label='Gaia B', color='blue')\n", "plt.plot(phi/(2*np.pi), dF_F_Kepler, label='Kepler', color='green')\n", "plt.plot(phi/(2*np.pi), dF_F_TESS, label='TESS', color='yellow')\n", "plt.plot(phi/(2*np.pi), dF_F_Gaia_R, label='Gaia R', color='orange')\n", "plt.plot(phi/(2*np.pi), dF_F_JWST, label='JWST', color='red')\n", "\n", "plt.xlabel('Phase (cycles)')\n", "plt.ylabel(r'$\\Delta \\overline{F}^{\\mathrm{obs}}\\,/\\,\\overline{F}^{\\mathrm{obs}}$')\n", "plt.legend();" ] }, { "cell_type": "raw", "id": "a4e8ca55", "metadata": { "editable": true, "raw_mimetype": "text/restructuredtext", "slideshow": { "slide_type": "" }, "tags": [] }, "source": [ "A clear trend can be seen in the light curves: from red wavelengths to blue (i.e., JWST |rarr| Gaia R |rarr| TESS |rarr| Kepler |rarr| Gaia B), the amplitude becomes larger and maximum/minimum light occur at later phases. The amplitude trend arises because the visible and NIR fall on the Rayleigh-Jeans tail of the star's spectral energy distribution, and therefore the flux variations due to temperature perturbations vary as :math:`\\lambda^{-4}`. The phase trend is a consequence of non-adiabatic effects that cause (for this particular mode) the radius perturbations to lag the temperature perturbations by :math:`\\sim 45^{\\circ}`.\n", "\n", "Note that the overall scaling of the light curves is not physically meaningful; as with any linear oscillation calculation, the normalization of perturbations is arbitrary. In the case of GYRE, this normalization is chosen to yield a mode inertia of :math:`M R^{2}`.\n", "\n", ".. rubric:: Footnotes\n", "\n", ".. [#photgrid-files] Generated by applying the :command:`specgrid_to_photgrid` tool to an extended-wavelength version of the demo grid.\n", ".. [#model-file] Generated using the `MESA `__ stellar evolution code." ] } ], "metadata": { "celltoolbar": "Raw Cell Format", "kernelspec": { "display_name": "Python 3 (ipykernel)", "language": "python", "name": "python3" }, "language_info": { "codemirror_mode": { "name": "ipython", "version": 3 }, "file_extension": ".py", "mimetype": "text/x-python", "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", "version": "3.12.2" } }, "nbformat": 4, "nbformat_minor": 5 }