Skip to content

Instantly share code, notes, and snippets.

@stharrold

stharrold/20150123_binstarsolver_examples

Last active August 29, 2015 14:13
Show Gist options
  • Star 0 You must be signed in to star a gist
  • Fork 0 You must be signed in to fork a gist
  • Save stharrold/e481d1239fd0e47e77bd to your computer and use it in GitHub Desktop.
Save stharrold/e481d1239fd0e47e77bd to your computer and use it in GitHub Desktop.
Download ZIP
Example use of binstarsolver from 2015-01-23.
Raw
20150123_binstarsolver_examples
{
"metadata": {
"name": "",
"signature": "sha256:2fde707fc28198eef57c33a4f9f18d0bc7bd1d9d0a726833b7938eaea4c1f120"
},
"nbformat": 3,
"nbformat_minor": 0,
"worksheets": [
{
"cells": [
{
"cell_type": "heading",
"level": 1,
"metadata": {},
"source": [
"Estimate physical quantities of a binary system using observed quantities."
]
},
{
"cell_type": "heading",
"level": 2,
"metadata": {},
"source": [
"Imports"
]
},
{
"cell_type": "code",
"collapsed": false,
"input": [
"from __future__ import absolute_import, division, print_function\n",
"import sys\n",
"import numpy as np\n",
"import scipy.constants as sci_con\n",
"import astropy.constants as ast_con\n",
"import binstarsolver as bss\n",
"%matplotlib inline"
],
"language": "python",
"metadata": {},
"outputs": [],
"prompt_number": 1
},
{
"cell_type": "heading",
"level": 2,
"metadata": {},
"source": [
"Example from Budding, 2007, Introduction to Astronomical Photometry"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Example from section 7.3 of [Budding, 2007, Introduction to Astronomical Photometry](https://books.google.com/books?id=g_K3-bQ8lTUC) (Object: YZ Cas, Source: Kron, 1942, Astrophysical Journal, 96, 173).\n",
"\n",
"Assumptions: Spherical stars. No limb darkening. No orbital eccentricity."
]
},
{
"cell_type": "code",
"collapsed": false,
"input": [
"print(80*'=')\n",
"print(\"OBSERVED VALUES\")\n",
"print()\n",
"print(\"Light levels:\")\n",
"light_ref = 1.0 # Between minima.\n",
"light_oc = 0.898 # During occultation minima.\n",
"light_tr = (0.745 + 0.733) / 2.0 # During transit minima.\n",
"print(\"Between minima: light_ref = {lr}\".format(lr=light_ref))\n",
"print(\"During occultation minima, primary eclipse: light_oc = {lo}\".format(lo=light_oc))\n",
"print(\"During transit minima, secondary eclipse: light_tr = {lt}\".format(lt=light_tr))\n",
"print()\n",
"print(\"Phases of eclipse events:\")\n",
"p0 = 0.0 # Mid-eclipse.\n",
"p1 = p0 - np.deg2rad(12.3) # Begin ingress.\n",
"p2 = p0 - np.deg2rad(3.5) # End ingress.\n",
"p3 = p0 + np.deg2rad(3.5) # Begin egress.\n",
"p4 = p0 + np.deg2rad(12.3) # End egress.\n",
"print(\"Mid-eclipse: p0 = {p0} degrees phase\".format(p0=np.rad2deg(p0)))\n",
"print(\"Begin ingress: p1 = {p1} degrees phase\".format(p1=np.rad2deg(p1)))\n",
"print(\"End ingress: p2 = {p2} degrees phase\".format(p2=np.rad2deg(p2)))\n",
"print(\"Begin egress: p3 = {p3} degrees phase\".format(p3=np.rad2deg(p3)))\n",
"print(\"End egress: p4 = {p4} degrees phase\".format(p4=np.rad2deg(p4)))"
],
"language": "python",
"metadata": {},
"outputs": [
{
"output_type": "stream",
"stream": "stdout",
"text": [
"================================================================================\n",
"OBSERVED VALUES\n",
"\n",
"Light levels:\n",
"Between minima: light_ref = 1.0\n",
"During occultation minima, primary eclipse: light_oc = 0.898\n",
"During transit minima, secondary eclipse: light_tr = 0.739\n",
"\n",
"Phases of eclipse events:\n",
"Mid-eclipse: p0 = 0.0 degrees phase\n",
"Begin ingress: p1 = -12.3 degrees phase\n",
"End ingress: p2 = -3.5 degrees phase\n",
"Begin egress: p3 = 3.5 degrees phase\n",
"End egress: p4 = 12.3 degrees phase\n"
]
}
],
"prompt_number": 2
},
{
"cell_type": "code",
"collapsed": false,
"input": [
"print(80*'=')\n",
"print(\"CALCULATED VALUES\")\n",
"print()\n",
"print(\"Smaller-sized star (_s) is brighter, primary.\")\n",
"print(\"Greater-sized star (_g) is dimmer, secondary.\")\n",
"print()\n",
"print(\"Integrated flux of stars relative to integrated flux of total system:\")\n",
"(flux_intg_rel_s, flux_intg_rel_g) = bss.utils.calc_fluxes_intg_rel_from_light(light_oc=light_oc, light_ref=light_ref)\n",
"print(\"flux_intg_rel_s = {firs}\".format(firs=flux_intg_rel_s))\n",
"print(\"flux_intg_rel_g = {firg}\".format(firg=flux_intg_rel_g))\n",
"print()\n",
"print(\"Ratio of stellar radii from relative light levels:\")\n",
"radii_ratio_lt = bss.utils.calc_radii_ratio_from_light(light_oc=light_oc, light_tr=light_tr, light_ref=light_ref)\n",
"print(\"radius_s / radius_g = radii_ratio_lt = {rrl}\".format(rrl=radii_ratio_lt))\n",
"print()\n",
"print(\"Numerical solution for orbital inclination:\")\n",
"phase_orb_ext = p4\n",
"phase_orb_int = p3\n",
"incl_rad = bss.utils.calc_incl_from_radii_ratios_phase_incl(radii_ratio_lt=radii_ratio_lt, phase_orb_ext=phase_orb_ext,\n",
" phase_orb_int=phase_orb_int, show_plots=True)\n",
"print(\"incl = {incl} degrees\".format(incl=np.rad2deg(incl_rad)))\n",
"print()\n",
"print(\"Stellar radii in units of the star-star separation distance\\n\"\n",
" \"(recalculated using timings from eclipse events):\")\n",
"sep_proj_ext = bss.utils.calc_sep_proj_from_incl_phase(incl=incl_rad, phase_orb=phase_orb_ext)\n",
"sep_proj_int = bss.utils.calc_sep_proj_from_incl_phase(incl=incl_rad, phase_orb=phase_orb_int)\n",
"(radius_sep_s, radius_sep_g) = bss.utils.calc_radii_sep_from_seps(sep_proj_ext=sep_proj_ext, sep_proj_int=sep_proj_int)\n",
"print(\"radius_sep_g = {rg}\".format(rg=radius_sep_g))\n",
"print(\"radius_sep_s = {rs}\".format(rs=radius_sep_s))"
],
"language": "python",
"metadata": {},
"outputs": [
{
"output_type": "stream",
"stream": "stdout",
"text": [
"================================================================================\n",
"CALCULATED VALUES\n",
"\n",
"Smaller-sized star (_s) is brighter, primary.\n",
"Greater-sized star (_g) is dimmer, secondary.\n",
"\n",
"Integrated flux of stars relative to integrated flux of total system:\n",
"flux_intg_rel_s = 0.102\n",
"flux_intg_rel_g = 0.898\n",
"\n",
"Ratio of stellar radii from relative light levels:\n",
"radius_s / radius_g = radii_ratio_lt = 0.539115831462\n",
"\n",
"Numerical solution for orbital inclination:\n"
]
},
{
"metadata": {},
"output_type": "display_data",
"png": "iVBORw0KGgoAAAANSUhEUgAAAYUAAAEoCAYAAAC3oe14AAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJztnXeYVdX5tu+HIiBFJBYUUYwdOwqCdazBbowFkmgwRo2J\nUWNi1JSfGBPLF5MYE2OMPSbGEkuwtzhWFAuCDSwIioggoiI2yvv9sdaZ2RxOn3Nmn5l57+va19ll\nlefs9u7V3iUzw3Ecx3EAOqUtwHEcx6kf3Cg4juM4TbhRcBzHcZpwo+A4juM04UbBcRzHacKNguM4\njtOEG4UsJF0i6ZeJ7eMlvSfpY0krS9pB0muSFkg6IE2tlSJpuqTd09ZRK+K1GVRh3EZJR1dXUcuQ\nNEjSUklt4nmVNEbSo4ntpuuR/XzVIO+1Y36qVR558l0q6autmWetaBM3WbWIL8NP4wt+vqTHJR2X\nvIHM7Hgz+00M3xX4PbC7mfUxs/nAr4GLzKy3mY1L55+0GItL2cRzuFuV9VSVeG2mVxqdCs9NWyAN\no5e8Hsnnqxpk349m9lbMr91ew1rToYwC4WHfz8z6AGsD5wGnAVfkCd8f6A68kti3NvByJZlL6lxJ\nvDrDgFb9CnOqSlkvS0ldaiWkhLxLeV78fqw2ZtZhFuBNYLesfUOBJcDguH01cDawAbAQWAosAB4E\nXo9hPwU+BroCKxGMyixgZozbKaY1Bngc+APwPqGUsQJwATADmA1cAnSP4RtiGqcA78U0xyS09iCU\nXKYDHwKPJuIOB54A5gPPA7sUOQ+nAy8BHwBXAt0Sx/eLacyP+jeP+69N/P8FwKnxfJ0Sjw+I5+sH\ncXs9YF6xdOOxNYGbgTnANOBHiWNjgRuBa+J5fxHYpsD/Wwp8NXE9LwbuiHGfzByLx/cEpsTz+Weg\nETg6cfy7hI+AD4B7gLWz8vkR8AYwF/h/gMqIexzwajwff0kc6xTvkbkx7R/G8Jn7qtg99xjwu5jv\nNGBkPPZbYDHwWbx+F+U4d4NiXt8l3KONcf9NwLvxPD1MfF7isa8A44CPgKeinkcLXI+z81y3MSz/\nvHwV+F/cngv8E1gpz/3404T+zPlYM2qbB7wGfC9P3tvF/5e8fl8HJsX1YcD4eK1mxXula57/2Miy\n99CYrPOxMXB/1DQFODRxbB/Cc/lxvLY/afX3ZGtnmOZCDqMQ988AjovrVwG/juvrJG+wXGkAtxJe\n7D2AVeNDcWziZlhEeKg7EUodfwRuA/oCveINe04M3xDDjwU6A3sTDFPmIbg4PiBrxPSGE4zMgPjQ\nZB7+PeL2KnnOw3Rgcoy3MuElcnY8tjXBIA0lfIEdGf9z1zz//yhgXFz/JsFwXh+3vwvcWizd+F+e\nBX4JdAHWJbwM94pxxxJeZCNj3HOA8QWuc/ZL6H1g23hO/wn8Ox5bJT58B8djJ8fz/914/EDCi2Sj\nqPEXwONZ+TwYr+VAYCrxZVBi3HFAnxh3DvC1eOz7hNJp5vo8RHj5ZV50xe65L4Gj47n6PvBOIt+H\nMv8vz7kbFLVdHdPvlki3Z7xefwQmJuJcH5cewKaEl9kjea5H0/OVI+8xLP+8rAfsHvNdhWCQ/ljg\neczoz5yrR4C/EJ6TLeN53jVP/q8DeyS2bwJ+FteHEAxDJ8J74WXgpDz/cZlzTMIoxHP4NvCdmNZW\nBGO3cTz+LrBDXF8J2LrV35OtnWGaS/YNlNg/HjgjcdNmXpDL3GDZaQCrA58Tv9bjvtHA/xI3w4zE\nMQGfsOyX6ghgWlxvIHz1JPN7L3Ezfkri6zoR5jTgH1n77gGOLHAejk1s7w28HtcvyX5oCV8zO+U6\nh/Gh/SD+t0uAY4G347FrgJOLpLsz4SttRtaxM4Ar4/pY4L7EscHApwWuc/ZL6O9Z//WVuH4k8ERW\n3LdpNgp3s+zD3YlgpAcm8tkrcfx44IEy4m6fOH4DzS+g/2Vdnz1j+E4l3nOvJY6tGOOuFrcfIvEV\nm+PcDYrhBxUI0zeG6U0wpl8CGyaO/5b8JYWm5ytHumOy74McYQ4Cnsv3TCf0dyIY28VAz8Txc4Cr\n8qR9NnBFXO9NeFYH5gl7MnBLnv9YyCgcTsJgxn2XAv8X12cQnqE+hc5DLZeO1qaQj7UIL7ZyWYfw\nBfNubLieD/yN8PWW4e3E+qqEh/TZRPi7CV9AGeaZ2dLE9qeEEsUqhC+nN/LoODSTZkx3B0KbSD6S\nut4iFLMzaf0kK621EseXwczeILzstgJ2IlTTzJK0IeGF/3CRdNeIx9bMOnYGsFoiq/eyzkn3Mnrj\nJON+RjifxP80Myts8rysA/wpoWle3D8gT/js81gs7uzEeuY6Qzgn2ekmNRW755rSNbNP42qvxHGj\nOE35S+ok6TxJr0v6iPAiNsI9uSqhdJdPb7kk00HS6pKulzQz5n0tobqqFNYEPjCzhVnaBuQJfx1w\nsKQVCKXHZ83s7ahjQ0l3SHo36vhtGTqSrANsl3Wvf5Ng7AG+QahCmh47BQyvII8WkVojUr0gaSjh\n5nmsguhvA18AX8l6kSdJPoDvE15Kg83s3TLzep/whbg+oeonyVvAtWZ2bBnprZ21/k4ird+a2Tl5\n4uV6oTwMHEqoYpol6WHC19HKhDaEgunGG/9NM9uwjDyrwSxCNU9GhwhflxneInzV/rtAGmvT3BEh\n+zwWi5uPd1n++mQo5Z4rRKnnMhnuW8ABhF54MyT1pbl0OJfwNb42ofosW29L9Z1DqDrbzMw+lHQQ\noT4/X/gks4B+knqZ2ScJbdkfAiEhs1ckzSCUJr9JMBIZLiFUcR5uZgslnUx4gediIaGaKEPy4+wt\n4GEz2yuPhmeAg2Ij+48IbWktOZ9l0xFLCgKQ1EfSfsC/CS/Ul5LHSyG+2O8D/iCpd/yiWk/SznnC\nLwUuAy6UtGrUMUBSzhskR9wrY15rSOosaUT8qvknsL+kveL+7pIaJOX7IhLww5h3P0J99w3x2GXA\n9yUNU6CnpH0lZb403yNUGSV5GDiBUH8LoaHtBEKROfPQFkp3ArBA0s8k9Yj/YTNJ2yb0VkqhuHcB\nm0r6euxlcyLLPsB/A34uaTCApJUkHZqVxk8l9ZU0MMa/oYy42TozWm8ETozXZ2VCpwCg/HsuB7mu\nXzF6EQzRB5J6El7UGT1LgFuAsfHaDSbUl+ej3GvZi/CS/Tjez6dmHc/7f+JX/hPAuZK6SdqC0M71\nzwL5XUeoGtqJ0KaQ1LEA+FTSxoSqwnw8Tyhx9JC0PqF9J8OdwIaSvi2pa1yGSto4rn9L0krxvC4g\nGMRWpSMahdslfUyw2GcQevMclThuLPv1UezL6khCI1aml8lNNL9YstOCUP//OvBkLIbeDyS/kAvl\n91PgBeBpQnXEuYT2h5mEL96fExrS3gJ+Qv7ra8C/CC+XNwgNor8BMLNngWMIjXMfxGNHJuKeC/wy\nFn1PifseITw0GaPwOKHRMbNdMN1o8PYjVEFNI3x9/p3QCJvRm31eCp2n7OuXM66ZvU8o4ZxHKImt\nT6LEaGa3AecD18dr9QLwtay0/kv4gpxIqDq7ssS4uTQlDei9wCTgGUKvrGT4cu+55PafgEMkfSDp\nQnKTHf8fhLrudwg9v8ZnhTmBcP1nE/7/leR/hnLpK3TsLEIj70fA7Sx/LnLdj8njowntDLMIxuv/\nzOx/efKH8JG4M/CgmSWrlH9KKD18TLg3r8/xvzL8kdDO8h6hDeWfNN9zC4C9gFGE8/lu/A8rxLjf\nBt6M98yxhFJaq6LmDznHccpB0lJgfTOblrYWx6kWHbGk4DiO4+TBjYLjVI4Xs512h1cfOY7jOE14\nScFxHMdpwo1CG0bSGZIuq0I6TW5/VUPXxrG73b21SLtWSBor6doWxD9X0knV1NRaqE5ddkvaSdKU\nKqd5gaTvVzPNtopXHzk16UWj4D9/GtClwkFWdYGkMwnn5ogK4q5K6Kq6npl9UXVxNaaW17Deem5J\n6k8YL7OemS1KW0+a1NUXgNMuaetujVuifwxwZ1s0CK1E3dwbZjab4IurTU6cVU3cKKSApNMk3ZS1\n70+S/hTXx0h6Q2EyoGmSvpknnaaqjURR/0hJMyTNlfTzRNhOkn6u4L/mY0nP5BrxLOlqSWfH9QYF\nnzOnKMw+N0vSmETYfSVNlPSRpLfiV3WGzMC1D2N+w7X8jFzbS3pa0oeSJkgakTjWKOnXkh6L8e+V\nlNPXTBxRfIekOXFQ1u3J/1YsrcQ5e1/SL1VgIqH4P56Ig6Wel7RLrnCRkTT7fiLqWpBYlkg6soRz\nsaakcZLmKcz6973EsbGSbpJ0bfxvkyVtEKsW34v/a89E+JUkXRGv5UxJZ2eqh+I9ckG8d94A9i3w\n35C0STy38yW9KGn/rHN+dGK76dpLytwbk+J5ODQr3W7xPGya2LeqwgRZq8T7MumbaU1JN8frP03S\nj+L+7pI+Uxi1j6RfSFqkODo//vc/JrJuLPafOwRpeeLryAvBl8lCoFfc7kwYcTmM4DPlI2CDeGx1\nEr7rs9I5k+CiA5q9Q14KdAO2IPhK2igeP5XgMymT7hZAv7ie07UxxV157wJsGtc3J4xoPTBur8Py\nHmbH0Owtsh/BN/23CB8nowijc1eOxxsJo57XJzgCfAg4N8956Efwfd+dMLL2RqLL7mJpETyuLgC2\nJzia+x1hNGrGE+7YxDku10X5HPLM+xDP5cyYZrFzkdf9M81uxfeM1+gagmv0M+L294heeGP4Qm63\nC7rsztLflTAy/3SCD7VdCaN9M/dXXk+h2fdcnvNzBfCbxPYPgbsS92XGE28xt+sPAwfH9fvifTAy\ncV4PTOSRcYKX+jsizcVLCilgZm8BzxFeZAC7EVxBT4jbS4HNJfUws/fMLN9Mb7mK32eZ2RdmNpng\nJmHLuP97wC/M7LWoYbItO4w/X7qLCEZiiZndTXAnvFFM42GLPqPM7AXC0P9dcqSRi32BqWb2LzNb\nambXs2zx3Qgujl83s88JL/qtciVkZh+Y2a1m9rkFx2fnJHQUS+sQwnwQT1ioS/4/8o8/+DbhxXRP\nzPcBghuKffKE70swOMug4EH2auAwM3un0LlQ8Km0PXCamX1pZpOAy1nW9cgjZna/BX85/yF47zwv\nbt8ADFLw9bU6wRj92Mw+M7O5wIUEIwRwGGGugncsTD17Dvmv43CCS+rzzGyxmT1EcPORs1RbAdcl\ndMHyDuoyDCUY5d9EHW8Szk8m7sPALgoO5jYHLorb3QlzbDySSGsB4Zp1aNwopMd1BL8sEG74fwFY\ncPN7OOGrbVasFtmojHTzuWNei9xut4uRz5U3kraT9FAstn9ImEmsHLfG2S6WZ7Csi+7kf0m6vF4G\nSStKujRW+3xEeBGsJC0zeXu+tJZxn21mn9Hs5jqbcl2Uzyf45U9qXYngL+kXZvZEQkO+c7EGxd0/\nz8n6b+9b/PSN2xD+bzG324VcdmezZlbYpOZq0AisqOBAcRDh4+bWHOGKuV1/mFCyGELwP/UA4YNh\nO8IcIvMTafUmzCzXoXGjkB7/ATKeTA8i8RVkZvdZcK3bn/DFmK/baTldx94mVJ+UQqnpXkeYRW4t\nM+tLeMFk7qliabxDeKCTrEOz6+ly+AnBqeAwM1uJ8NAnvY4WYhbBYAIgqQf5DVvGRfnKiaW3mf2/\nPOEnE0tVMe1OhHP2oJldnghX6Fw0uX9OHMvr/rkISbfbGf0rmdnm8Xghl93ZzAIGZhne5PUr5D66\nKLGUcyPhw2k0cHuWYczwNsHtevKa9DGz/eLx8YRr8HXC1KKvxP+1D8HwJNmEZlfvHRY3CikRi+6N\nhGqEaWY2FUDSapIOVHBRvIjwcOVzn1tO743LgbMlra/AFpkGuBxplppuL2C+mX0paRihxJMxBnMJ\n1WD53DTfTXAhPFpSF0mHE+auvSNLS6k6PgM+iv/pzBxh8qV1M8HteMYN+dgCYct1UX4Xy1Zj/ZYw\nydLJOcLlPBcWPOCW6/45J1bc7XZel905eJJQavyZgsvnBoKn2+vj8ULuo6E0F96ZKqR8VUdQxO26\nhUmGniW0SWQa/Z8glMQfzkprF8J92aFxo5Au1xHmn03e8J2AHxO+uOYR/Lrn891ejpvvPxAe+vsI\nDdmXERpds+OVk+YPgF8ruCL/Fc1zCWQext8Cjyv0CNoumbaZzSO8RH5CaKz9KbBfVjtHqS6XLyQ0\nnL5PeODvzhE2Z1qxTeRHhJfZLEK98hzCF3V22HJdlP8D2CfWX0N4wW0HzFdzD6TR8T8XOheF3D+X\n4la8VLfbxVx2NycY2l/2J7RRzCU0hB9hZq/GIHndR0fGAtfEKp9D8uQxgdCGtQbLv6wz12QJhd2u\nQ3j5dyEYkMx20tU7ktYglBRuy6WlI1HzwWuSRhIe2s7A5WZ2fo4wDYSbqCuhPrShpqIcJw+xmmY+\nYWDVjCqk91tgjpn9qcXinJoh6QJCG8Pf0taSNjU1CrHFfyqh6947hMlhRsd6vUyYvoRJWb5mZjMl\nrWJh8hPHaRVi//oHCdVGvweGmtk26apynHSodfXRMIL1nR6Lm9eTmBM38k3g5lg0xw2CkwIHED5a\n3iHUc48qHNxx2i+1NgoDWLbbWmawTpINCL0rHlIYZVu2jxnHaQlmdkzstdLXzPbMjOVwnI5Ilxqn\nX0rdVFdCH+LdCT0zxkt60h9Mx3Gc1qfWRuEdYGBieyDL969+m9C4/BnwWfSLsiVhOHoTktydq+M4\nTgWYWcnd12tdffQMsIGCs7YVCCN1x2WF+S+wY+xfvCKhy15Otw75fHWktZx55pmpa3BN7UuXa3JN\n1V7KpaYlBTNbLOkEQt/nzsAVZvaKpOPi8UvNbIqkewijP5cCl1l+Xz+O4zhODal19REWnKjdnbXv\n0qztC4ALaq3FcRzHKYyPaG4BDQ0NaUtYDtdUOvWoyzWVhmuqHW1mOk5J1la0Oo7j1AuSsDpqaHYc\nx3HaEG4UHMdxnCbcKDiO4zhNuFFwHMdxmnCj4DiO4zThRsFxHMdpwo2C4ziO04QbBcdxHKcJNwqO\n4zhOE24UHMdxnCaKOsSTNAQYDewMDCJMnDMDeAS4zswm1lKg4ziO03oU9H0k6S5gPmEOhAnAu4TJ\nzdcgzL+8P9DXzPatuVD3feQ4jlM25fo+KmYUVjez94pkuJqZzSlDY0W4UXAcxymfqhqFHIl3B8zM\nvqhEXEtwo+A4jlM+5RqFgm0KkjoBBxHaFLYnNExL0hJgPPAv4DZ/WzuO47QPilUfPQI8SmhTeD5T\nQpDUDdgaOADY0cx2rrlQyR56yJBoWsJ+Cu4rJUwl8Tp1gs6dy/9NpuU4jlNrqt2m0K1YVVEpYaqB\nJNt5Z8MMMpIz6/n2lRKm0nhLl4ZlyZLyfs2CgajUqCR/u3SBrl0L/1Z6LFeYFVaAbt2al+zt7H1d\nu7oBdJy0qbZR6Fcospl9UIa2FtFe2hQyBqVcY5Lrd/HisCxalPu30LFKwnz5JXzxRfNSbHvx4uKG\nJLm94orQo0f4Ta6Xu69bNzdGjpOh2kZhOmFcQk7MbN2y1LWA9mIUOhJLl5ZuSD7/PCyffhqWzz4r\n/Fvo2KJFwUD07g29eoXffEuh4yutBH37BiPjOG2VmvY+ShM3Ck6pLFkSDMSCBZUtn3wSfj/6CObP\nD9VgK6/cvPTtu+x29v6vfAVWXRX69QtVfY6TJjUxCpK+DjxkZh/G7b5Ag5ndVrHSMnGj4KSBGSxc\nGIzDhx+G38ySb3vePJg7N2z36xcMRGZZbbXlt1dfHdZcE/r08Wovp/rUyihMMrMts/Y9b2ZbVaCx\nItwoOG2NxYubDcScOeE3syS3Z8+Gd94JcdZcMywDBuT/XWGFdP+X07ao6jiFZLo59nnB2HEK0KVL\nKAWsvnpp4T/+GGbNCgYi8zttGjz6aPP27NkhvXXXDcugQcv+DhgQ8nWcSim1pHAVwQfSxQQD8UNg\nZTMbU1N1y2rwkoLT4Vm0CGbOhOnTw/Lmm8v+zpkDa60FG20Ulo03bl7v39+rpzoitao+6gX8Ctg9\n7rof+I2ZLaxIZQW4UXCc4nzxRTAOr74KU6bA1KnNyxdfwIYbBkOx2Waw5ZZhcWPRvvHeR47j5OSD\nD4JxmDIFXngBJk0KS6dOzQZiiy1gyBDYZBPvOdVeqFVJYSPgp4T5FDI1lmZmu5UQdyRwIaEN4nIz\nOz/reAPwX2Ba3HWzmf0mRzpuFBynypiF9oqMgZg0CZ59NlRDbbstbLdd89K/f9pqnUqolVGYDFwC\nPAcsibvNzJ4tEq8zMBXYA3gHeBoYbWavJMI0AKeY2QFF0nKj4DitxLx5MGECPPkkPPVUWO/TB3ba\nCXbdFRoaQsO2VzvVP7UyCs+a2TYViBkBnGlmI+P26QBmdl4iTAPwEzPbv0habhQcJyXMQjvFww9D\nYyM89FDoGpsxEHvtFbrMOvVHuUah1Dmab5f0Q0lrSOqXWUqINwB4O7E9M+5LYsD2kiZJukvS4BI1\nOY7TSkihB9Oxx8J114Uqp3vvheHD4a67QsP1kCHwq1+F0sWSJcXTdOqTUksK08nhA6mY7yNJ3wBG\nmtkxcfvbwHZm9qNEmN7AEjP7VNLewJ/MbMMcaXlJwXHqlMWL4Ykn4M47wzJnDuy9Nxx6aChF+IC7\n9Kir3keShgNjE9VHZwBLsxubs+K8CWyT7YFVkp155plN2w0NDTQ0NNREt+M4LWP6dLj9drjxRnjp\nJTjwQDj8cNh99+BLyqkdjY2NNDY2Nm2fddZZVfWSuruZPRi/+HOVFG4pmLjUhdDQvDswC5jA8g3N\nqwNzzMwkDQNuNLNBOdLykoLjtEFmzoSbboIbboA33oBvfQuOPho23zxtZR2DarvOPsvMzpR0NbmN\nwlElCNqb5i6pV5jZuZKOi/EvlfRD4HhgMfApoSfSkznScaPgOG2cN96Aq66Cq68OLjmOPhpGjQo9\nm5zaUFfVR9XEjYLjtB+WLAkN1VdcEXoyHXkknHgifPWraStrf1S195GkMbEKKN/xFSQVLS04juMk\n6dwZ9tkHbr4ZJk+G7t1h2DA4+GB47LG01XVsilUfnQAcDUwhDDybTXCI1x/YFtgYuMzM/lpzoV5S\ncJx2zcKFcM018Ic/BKd+Y8eGMRBOy6h69ZEkATsAOwJrx90zgMeAJ1rrTe1GwXE6BosXh7EQZ58d\n2h3cOLQMb1NwHKddsHgx/Otf8OtfB8+uF1wQHPU55VHt3kd/TmwazZPtGICZnViJyEpwo+A4HZMv\nvoCLL4Zzz4XDDgslh1VXTVtV26Habi6ejUs3YAjwKvAasDXgYxQdx6k53brBKacEl9+dO8Omm8KV\nVwZ/TE71KdXNxVPAjma2KG53BR4zs+1qrC+pwUsKjuPw/PNwzDHQqxdcemmYOMjJT60c4vUFksNL\nesd9juM4rcpWWwWnewcdBNtvD+ef7w74qkmpJYWjgLFAY9y1C8Gn0dW1EpZDg5cUHMdZhhkz4Igj\nQrXSP/4BAwemraj+qFnvI0lrANsRGpmfMrPZlUmsDDcKjuPkYsmSUFq48EL461/hkEPSVlRf1NIo\nrAxsCHSnuffRI5WIrAQ3Co7jFOLpp4Mn1oMPhvPOgy55fTF0LGo189oxwInAWsDzwHBgfClzNFcL\nNwqO4xTjgw9g9OgwxuGGG2CVVdJWlD61amg+CRgGzDCzXQldUj+qQJ/jOE7N6NcvzAQ3bBhsuy28\n8ELaitoepRqFz83sMwBJ3c1sCrBR7WQ5juNURufOYaDbeefBHnsEL6xO6ZRqFGbGNoXbgPsljQOm\n10yV4zhOCxk1KlQhjRoF//532mraDmX7PpLUQBizcI+ZfVkLUXny9TYFx3HK5sUXg5vu00+HH/wg\nbTWtT7ltCkXb5+N8Ci+a2cYAZtZYuTzHcZzWZbPN4OGHYbfdQvfVH/0obUX1TVGjYGaLJU2VtI6Z\nzWgNUY7jONVk3XVD28Juu4WeST/+cdqK6pdSe/L2A16SNAFYGPeZmR1QG1mO4zjVZdAgaGyEXXcN\nM70df3zaiuqTUo3Cr3Ls8wp+x3HaFGuvDfffDzvtFMYwHHpo2orqj6pMsiNpvJmNqIKeQnl4Q7Pj\nOFVh0iTYc8/QK2n33dNWU1tqNXitGN2rlI7jOE7N2XJLuOmm0F110qS01dQX1TIKjuM4bYpddoE/\n/xkOPBDmzElbTf3gRsFxnA7LqFHwrW8Fz6pfttqoq/rGjYLjOB2as88OPpN8/EKgZKMgqb+k/SXt\nJ2m1rMNHVlmX4zhOq9CpE1x7bRjgdu21aatJn1JdZx8G/A54OO7aGTjVzG6qobZsDd77yHGcmjF5\ncuiJ9NhjsFE7cvdZq/kUJgN7mNmcuL0q8KCZbVGx0jJxo+A4Tq259NIwe9uTT0KPHmmrqQ616pIq\nYG5ie17c5ziO02449ljYeGP42c/SVpIepRqFe4B7JY2RdBRwF3B3KREljZQ0RdJrkk4rEG6opMWS\nDi5Rk+M4TlWR4G9/g1tv7bjzMJRafSTgYGBHgnuLR83s1hLidQamAnsA7wBPA6PN7JUc4e4HPgWu\nMrObc6Tl1UeO47QKd94JJ5wQ2hl6905bTcuoSZtCC8SMAM40s5Fx+3QAMzsvK9zJwJfAUOAONwqO\n46TNd78L3brBJZekraRlVLVNQdLj8fcTSQuylo9LSH8A8HZie2bcl8xjAHAgkDn1/uZ3HCd1/vAH\nuOOO4Fm1I1HQS6qZ7RB/e1WYfikv+AuB083MYjVVXos2duzYpvWGhgYaGhoqlOU4jlOYvn3hoovC\nbG3PPw8rrJC2otJobGyksQWWrNQ2hWvN7Ihi+3LEGw6MTVQfnQEsNbPzE2Gm0WwIViG0KxxjZuOy\n0vLqI8dxWhUz2HdfaGhouz2SajVOYaKZbZ3Y7gJMNrPBReJ1ITQ07w7MAiaQo6E5Ef4q4HYzuyXH\nMTcKjuO0Om+8AdttBxMnwsCBaaspn2q3Kfxc0gJg82R7AjAHGFcoLoSpPIETgHuBl4EbzOwVScdJ\nOq5UkY7jOGmx3nqhJ1JHmcKz1JLCeWZ2eivoKaTBSwqO46TC55+HQW3XXhtmbWtL1KxLqqSVgQ1I\nTKhjZo9oOAA2AAAcnElEQVSUrbBC3Cg4jpMm//wnXHwxPPFEGOTWVqiJmwtJxwCPAPcBZxGqg8ZW\nItBxHKct8s1vhhLDLcu1eLYvSnVzcRIwDJhuZrsCWwMf1UyV4zhOndGpE5x/Ppx+OixalLaa2lGq\nUfjczD4DkNTdzKYA7ci5rOM4TnH22gvWXRf+/ve0ldSOUo3CzNimcBtwv6RxwPSaqXIcx6lTzj0X\nzjknVCW1R8r2fSSpAegD3GNmrTarqTc0O45TL+y3H+yzTxjtXO9UvfdRHID2oplt3FJxLcGNguM4\n9cKECXDIIfDaa8FpXj1T9d5HcQDaVEnrtEiZ4zhOO2HYMNh0U7j66rSVVJ9SB689SuhxNAFYGHeb\nmR1QQ23ZGryk4DhO3fDkkzBqVCgtdO2atpr8lFtSKOglNcGvcuzzN7TjOB2W4cNh0CD4z39g9Oi0\n1VSPqkyyI2m8mY2ogp5CeXhJwXGcuuL22+Gss+Dpp+t3lHNNRjSXQPfiQRzHcdoX++4Ln3wCj7Sa\nw5/aUy2j4DiO0+Ho1AlOOQUuuCBtJdXDjYLjOE4LOOKI0EV1ypS0lVQHNwqO4zgtoEcPOOYYuOSS\n4mHbAtVqaN7czF6ogp5CeXhDs+M4dclbb8HWW4ffnj3TVrMs1Z557fH4+0ly5rW4fJwJV2uD4DiO\nU8+svTbssANcf33aSlpOVUoKrYGXFBzHqWfuvht+9St45pm0lSxLtUsKfeJvv1xLS8U6juO0F772\nNfjggzBmoS1TsKQg6U4z21fSdHKMYDazdWuoLVuLlxQcx6lrzj8fXn0VrrgibSXN1GyO5rRxo+A4\nTr3z7rsweDDMnFk/Dc5V9X0kaUih42b2XKkZOY7jtHfWWCM0ON9ySxi/0BYpVn3USKg26gFsA0yO\nh7YAnqm1v6MsLV5ScByn7vnPf8KYhQcfTFtJoKoNzWbWYGa7ArOAIWa2jZltQ3CjPatlUh3Hcdof\n++8PkybBjBlpK6mMUkc0b5wci2BmLwKb1EaS4zhO26VbNzj8cLj22rSVVEapk+xcD3wC/BMQ8E2g\nl5m1mhdxrz5yHKet8PTTYY6F115L36V2rVxnHwW8DJwEnBjXjypfnuM4Tvtn222hc+e2OWbBu6Q6\njuPUgP/7P1i4EH7/+3R11KSkIGlDSf+R9LKkN+MyrcS4IyVNkfSapNNyHD9Q0iRJEyU9K2m3UsU7\njuPUK4cfDjfeCEuXpq2kPEqtProK+BuwGNgVuAb4V7FIkjoDfwFGAoOB0ZKyG6gfMLMtzWxrYAzw\n9xI1OY7j1C2bbgorrQRPPJG2kvIo1Sj0MLMHCNVN081sLLBvCfGGAa/HOIuA64EDkwHMbGFisxfw\nfomaHMdx6ppRo+CGG9JWUR6lGoXP41f/65JOkHQwUMog7gHA24ntmXHfMkg6SNIrwN2EhmzHcZw2\nz+GHw003wZIlaSspnYJuLhKcBKxIeGGfDfQBvlNCvJJahs3sNuA2STsB1wIb5Qo3duzYpvWGhgYa\nGhpKSd5xHCcVNtgA1lwTHn4Ydmul1tLGxkYaGxsrjl+091EsIZxvZj8tO3FpODDWzEbG7TOApWZ2\nfoE4bwDDzGxe1n7vfeQ4TpvjnHNg9my46KJ08q967yMzWwLsKFU0BOMZYANJgyStABwOjEsGkLRe\nJu2MA75sg+A4jtNWOegguO02aCvftKVWHz0P/FfSTcCncZ+Z2S2FIpnZYkknAPcCnYErzOwVScfF\n45cC3wCOlLSIMGp6VAX/w3Ecpy7ZZJPg+mLiRBhS0O90fVCqm4uryT3JTquNavbqI8dx2iqnngo9\nesCvf936efskO47jOHXG44/DD34QvKe2NrXyfeQ4juNUyPDhobF5Wkl+INLFjYLjOE6N6dwZDjgA\n/vvftJUUx42C4zhOK3DggTBuXPFwaVN2m4KkO8xsvxrpKZSvtyk4jtNmWbgQ+veHd96BPn1aL9/W\naFNYzk2F4ziOU5iePWHEiPqZuzkflRiFiVVX4TiO0wHYZx+4++60VRTGu6Q6juO0ElOnwh57wFtv\ntd40nd4l1XEcp07ZcEPo2hVeeiltJflxo+A4jtNKSLD33vVdhVTqdJyHlrLPcRzHKUy9G4VSfR9N\njNNlFtxXS7xNwXGc9kCma+qsWdC7d+3zK7dNoaCXVEl7A/sAAyRdBGQS7g0sqlil4zhOB6VnT9h2\nW3j00dAbqd4oVn00C3gW+Dz+ZpZxwNdqK81xHKd9svvu9TteodTqo65mlmrJwKuPHMdpL4wfD8cf\nD88/X/u8quo6W9ILBeKamW1RjriW4EbBcZz2wuLFsMoq8NprsOqqtc2rqm0KwP4t1OM4juNk0aUL\n7LQTPPQQHHZY2mqWpaBRMLPppSQiabyZjaiKIsdxnA7AbruFdoV6MwrVGrzWvUrpOI7jdAjqtbHZ\nRzQ7juOkwGabwccfw4wZaStZFjcKjuM4KdCpU3MVUj3hRsFxHCcl6rEKqWTX2ZL6A0MBAyaY2ZzE\nsc3NrFD31RbjXVIdx2lvvPpqMAy1dKVdE9fZkg4DngIOBQ4DJiQd4tXaIDiO47RHNtgAFi2qr3aF\nYuMUMvwSGJopHUhaFXgQuKlWwhzHcdo7Uhiv8OijMGhQ2moCpbYpCJib2J5Hs3M8x3Ecp0IyRqFe\nKLWkcA9wr6TrCMbgcKCOPYI7juO0DXbeGf7617RVNFOqQzwBBwM7EhqaHzWzW2usLVuDNzQ7jtPu\nWLIEvvKV0Oi82mrVT78mDc0WuNnMfmxmp5RjECSNlDRF0muSTstx/FuSJkmaLOlxSa3mZM9xHCdt\nOneG7beHxx5LW0mgoFGQ9Hj8/UTSgqzl42KJS+oM/AUYCQwGRkvaJCvYNGDn6HH1bODvlfwRx3Gc\ntsrOO9dPu0JBo2BmO8TfXmbWO2vpU0L6w4DXzWx6nI/heuDArDzGm9lHcfMpYK3y/4bjOE7bZaed\n4JFH0lYRKHWcwrWl7MvBAODtxPbMuC8fRwN3laLJcRynvbDttjB1avCFlDal9j7aLLkhqQuwTQnx\nSm4ZlrQr8F1gh3xhxo4d27Te0NBAQ0NDqck7juPULd26BcMwfjx8rYUTHTc2NtLY2Fhx/GIzr/0c\nOAPoAXyWOLQI+LuZnV4wcWk4MNbMRsbtM4ClZnZ+VrgtgFuAkWb2ep60vPeR4zjtll/8Iky+c9ZZ\n1U23qr2PzOwcM+sNXJDVntCvmEGIPANsIGmQpBUI4xvGZQlem2AQvp3PIDiO47R3RowIJYW0Kcch\n3srABiQm1DGzok0jkvYGLgQ6A1eY2bmSjovxL5V0OfB14K0YZZGZDcuRjpcUHMdpt8ydG3whffBB\ncKtdLcotKZQ6eO0Y4ERgIDARGA6MN7PdKhVaLm4UHMdp72ywAdx6a5iAp1rUZPAacBKhe+l0M9sV\n2Br4qHAUx3EcpxxGjIAnn0xXQ6lG4XMz+wxAUnczmwJsVDtZjuM4HY/hw9NvVyjVKMyMbQq3AfdL\nGgdMr5kqx3GcDkg9NDaX3NDcFEFqAPoA95jZl7UQlSdfb1NwHKdds3gxrLwyvP029O1bnTSr3qYg\nqYukKZltM2s0s3GtaRAcx3E6Al26wDbbwFNPpaehqFEws8XAVEnrtIIex3GcDk3aVUilurnoB7wk\naQKwMO4zMzugNrIcx3E6JiNGpDvpTqnjFBpy7DYze7jqivJr8DYFx3HaPXPmwEYbwbx51RnEVm6b\nQkklBTNrLJLpeDMbUWqmjuM4Tm5WWw369YMpU2Dw4NbPv1qDqbsXD+I4juOUQprtClX0sOE4juNU\ng6FD4emn08nbjYLjOE6dMWyYGwXHcRwnstVW8Mor8PnnrZ93tYzCkVVKx3Ecp8PTo0fogTRpUuvn\nXdAoSHo8/n4iaUHW0jSbqJm9UGuhjuM4HYm02hUKdkk1sx3ib6/WkeM4juNAMAqPPdb6+RYrKfSJ\nv/1yLa0j0XEcp+ORVkmh4IhmSXea2b6SpgPLBTSzdWuoLVuLj2h2HKfDsGhR8Jj67rvQu3fl6VR1\nRLOZ7Rt/B1UuyXEcxymXrl1hiy3g2WehoaH18i1oFCQNKXTczJ6rrhzHcRwnQ6YKqW6MAvAHQrVR\nD2AbYHLcvwXwDOD+jhzHcWrE0KEwblzr5lmwodnMGsxsV2AWMMTMtjGzbYCt4z7HcRynRqTR2Fzq\n4LWNk2MRzOxFYJPaSHIcx3EANtgA5s+HuXNbL89SjcJkSZdLapC0q6TLgBTG2jmO43QcOnUK03M+\n80wr5lliuKOAl4GTgBPj+lG1EuU4juMEhg6FCRNaL7+SZl6rB3ycguM4HZGbb4arroI77qgsfrnj\nFEqdjnND4BxgMKEnEoTpOL9akcoKcKPgOE5H5K23Qmlh9mxQya/2Zso1CqVWH10F/A1YDOwKXAP8\nq3x5juM4TjkMHBh+3367dfIr1Sj0MLMHCCWL6WY2Fti3lIiSRkqaIuk1SaflOL6xpPGSPpf0k9Kl\nO47jtH+k1u2aWqpR+FxSZ+B1SSdIOhjoWSxSjPMXYCSh6mm0pOyurPOAHwEXlC7bcRyn41CPRuEk\nYEVCz6NtgW8D3ykh3jDg9Vi6WARcDxyYDGBmc83sGWBRyaodx3E6EK05PWcxNxeZr/3DzeynwAJg\nTBnpDwCSNWEzge3KEeg4jtPRGTo0OMZbujSMXaglRZM3syXAjlIl7d7Lu9t2HMdxymOVVYIb7dde\nq31eRUsKkeeB/0q6Cfg07jMzu6VIvHeAgYntgYTSQkWMHTu2ab2hoYGG1nQd6DiOkyKZQWwbbVQ4\nXGNjI42NjRXnU+o4havJPclOwVHNkroAU4HdCQ70JgCjzeyVHGHHAgvM7Pd50vJxCo7jdFh+97vQ\nLfWii8qLV9VJdjKY2ZjyZDTFWyzpBOBeoDNwhZm9Ium4ePxSSf2Bp4E+wFJJJwGDzeyTSvJ0HMdp\njwwbBrcUq5upAsWm4xwLXGJm7+U5vgbwfTM7szbylsnLSwqO43RYFiyA/v3hww/DrGylUu2SwjPA\n9ZJWAJ4D3gUE9AeGAF/g4wscx3FqTu/eMGgQvPgibL117fIptU1hILADsHbcNQN43MwqbjQuFy8p\nOI7T0RkzBkaMgOOOKz1OTRziZWXQGehpZh+XFbGFuFFwHKejc/HFMHEiXH556XFq4hBP0r8l9ZHU\nE3gBeEXSz0qX5TiO47SU1hjZXOrYuMGxZHAQcDcwCDiiVqIcx3Gc5dliC3j9dfj00+JhK6VUo9BF\nUleCUbg9+jHyuhzHcZxWpFs3GDw4VCHVilKNwqXAdKAX8IikQcBHtZHkOI7j5KPW03OWZBTM7CIz\nG2Bme5vZUkLvo91qJ8txHMfJRa3daJfa0LyKpD9LmijpOeBCwghkx3EcpxWpdWNzqdVH1wNzgIOB\nQ4C5wA21EuU4juPkZuON4b33YP782qRfqlHob2Znm9mbZjbNzH4DrF4bSY7jOE4+OncOI5qfeaY2\n6ZdqFO6TNFpSp7gcDtxXG0mO4zhOIWrZ2FzMId4nNHc97QksjeudgIVm1rs2snJq8RHNjuM4wI03\nwnXXwW23FQ9bczcXaeFGwXEcJzBjRmhwnj0bis2JWZP5FGLCKwMbAN0z+8zskVLjO47jONVh7bWD\n++w33oD1169u2iUZBUnHACcSptOcCAwHxuNjFRzHcVodCXbYAR5/vPpGodSG5pOAYcB0M9sV2Bof\n0ew4jpMaGaNQbUo1Cp+b2WcAkrqb2RSgyPTRjuM4Tq2olVEotU3h7dimcBtwv6T5BF9IjuM4Tgps\nuSW89RZ88AH061e9dCuZZKeB4OLiHjP7snpSiubrvY8cx3ES7L47nHIK7Ltv/jA1mWQniZk1mtm4\n1jQIjuM4zvJsv331q5DKNgqO4zhOfbDjjvDYY9VN0wevOY7jtFE++QT69w8O8nr2zB2m5tVHjuM4\nTn3QqxcMGVLd0oIbBcdxnDbMHnvAAw9ULz03Co7jOG2YahuFkn0fOY7jOPXH0KEwbRrMnQurrtq8\nf9Ys+O9/y0/PSwqO4zhtmK5dYZdd4MEHm/c99RRstRWMH19+ejU3CpJGSpoi6TVJp+UJc1E8PknS\n1rXW5DiO057Ye2+4446w/vHH8I1vwJVXwj/+UX5aNTUKkjoDfwFGAoOB0ZI2yQqzD7C+mW0AHAtc\nUktN1aSxsTFtCcvhmkqnHnW5ptJwTcty0EFw553w5Zdw9tmw556w336VpVXrksIw4HUzm25mi4Dr\ngQOzwhwAXANgZk8BfSW1ifmf/cYsjXrUBPWpyzWVhmtaljXWgMGD4YYb4PLL4dxzK0+r1kZhAPB2\nYntm3FcszFo11uU4jtOuOPhgOPLI8Nu/f+Xp1NoolDoEOXu0nQ9ddhzHKYNRo8LvySe3LJ2aurmQ\nNBwYa2Yj4/YZwFIzOz8R5m9Ao5ldH7enALuY2XtZabmhcBzHqYCazNFcIc8AG0gaBMwCDgdGZ4UZ\nB5wAXB+NyIfZBgHK+1OO4zhOZdTUKJjZYkknAPcCnYErzOwVScfF45ea2V2S9pH0OrAQOKqWmhzH\ncZz8tBkvqY7jOE7tqfsRzaUMfmslHVdKek/SC4l9/STdL+lVSfdJ6tvKmgZKekjSS5JelHRi2rok\ndZf0lKTnJb0s6dy0NSW0dZY0UdLt9aBJ0nRJk6OmCXWiqa+k/0h6JV6/7epA00bxHGWWjySdWAe6\nzojP3guSrpPUrQ40nRT1vCjppLivLE11bRRKGfzWilwVdSQ5HbjfzDYEHozbrcki4MdmtikwHPhh\nPD+p6TKzz4FdzWwrYAtgV0k7pqkpwUnAyzT3bktbkwENZra1mQ2rE01/Au4ys00I129K2prMbGo8\nR1sD2wCfAremqSu2kx4DDDGzzQnV46NS1rQZ8D1gKLAlsJ+k9crWZGZ1uwAjCHNBZ7ZPB05PUc8g\n4IXE9hRg9bjeH5iS8vm6DdijXnQBKwJPA5umrYkw9uUBYFfg9nq4fsCbwFey9qWmCVgJmJZjf13c\nTzH/vYBH09YF9AOmAisT2mZvB/ZMWdMhwOWJ7V8CPytXU12XFCht8FuarG7NPaXeA1IbiR2/XLYG\nniJlXZI6SXo+5v2Qmb2Utibgj8CpwNLEvrQ1GfCApGckHVMHmtYF5kq6StJzki6T1DNlTdmMAv4d\n11PTZWYfAL8H3iL0rPzQzO5PUxPwIrBTrC5aEdiH8DFUlqZ6NwptphXcghlORa+kXsDNwElmtiBt\nXWa21EL10VrAzpJ2TVOTpP2AOWY2keUHSqaiKbKDhSqRvQlVfzulrKkLMAT4q5kNIfQGXKaqIeX7\nfAVgf+Cm7GMp3FPrAScTag/WBHpJ+naamsxsCnA+cB9wN/A8sKRcTfVuFN4BBia2BxJKC/XCe5L6\nA0haA5jT2gIkdSUYhGvN7LZ60QVgZh8BdxLqgdPUtD1wgKQ3CV+Zu0m6NmVNmNm78XcuoY58WMqa\nZgIzzezpuP0fgpGYXQ/3E8F4PhvPF6R7rrYFnjCzeWa2GLiFUN2d6rkysyvNbFsz2wWYD7xKmeep\n3o1C0+C3+JVwOGGwW70wDvhOXP8OoU6/1ZAk4ArgZTO7sB50SVol07tBUg9CPevENDWZ2c/NbKCZ\nrUuofvifmR2RpiZJK0rqHdd7EurKX0hTk5nNBt6WtGHctQfwEqG+PLX7PMFomquOIN3nbwowXFKP\n+BzuQejEkOq5krRa/F0bOBi4jnLPU2s1grSg8WRvQoPO68AZKer4N6Hu8EtCO8dRhMamBwjW+D6g\nbytr2pFQR/484cU7kdBDKjVdwObAc1HTZODUuD/Vc5XQtwswLm1NhPr75+PyYubeTvs8EXqtPA1M\nInz9rpS2pqirJ/A+0DuxL+1z9TOC0XyB4Om5ax1oeiRqep7QC7Ds8+SD1xzHcZwm6r36yHEcx2lF\n3Cg4juM4TbhRcBzHcZpwo+A4juM04UbBcRzHacKNguM4jtOEGwWnRUh6vMJ4DQkX1vurQrfoklaS\ndHxie01Jy7lBSJvoJrtfmXFuiO4UsvePkfTn6qlrOZL+kO2mw2mbuFFwWoSZ7VCFNG63xLzdZbIy\n8INEWrPM7NCWaqoBZQ0IkrQ+0NPM3qiRHhSpUnKXEBwOOm0cNwpOi5D0SfxtkNQo6SaFCVr+mQgz\nVNLjChPvPBUd+CXTaPrylXS1pD/F8G9I+kbc30vSA5KeVZiY5oAY/TxgPYXJV86XtI6kF2Oc7tHj\n5+To9bMhkd8tku6OE4/kNEiSfiVpgsKkJZcm9jdKOi/+l6lxvoiM64obFSZeuUXSk5KG5Ej32zHu\nREl/k5TrORxFwqWLpKNiXk8RfDll9q+qMCnOhLhsn9h/v8JkK5dlSirRZcxUSdcQRuIOlHRqjDtJ\n0thCOhUmKro6npPJkk4GMLPXgEFKYfIkp8q09nB1X9rXAiyIvw3AhwSPkQKeILy8VgDeALaJ4XoR\nJiRpoHlegzHAn+P61cANcX0T4LW43pno4gBYJbF/HZad42JQZhv4CdG/PLARMAPoFvN7A+gdt6cD\nA3L8t5UT6/8A9ovrDwG/i+t7EyYwAfgpcElc35QwCdKQuP0mwd3AJoSXfee4/6/AETnyvjsRd42o\n/SsEVwqPARfFY9cRvK0CrE3wgwVhcqrT4vrXCO5Q+sXzswQYFo/tBVwa1zsRfPfslEPnxcARBAd5\n9yV0rpRYvwbYO+170peWLV1wnOoxwcxmASjMp7AusAB418yeBTCzTMkiXxpGdNhlZq9Iyvh+7wSc\nG+utlwJrRudfhao/dgAuimlNlTQD2DDm8aBFN+OSXia8LN/Jir+bpFMJkwX1I/gouiMeuyX+Phfj\nZvK7MOb3kqTJWekJ2J3gNfaZeA56ALNzaF8HeDeub0eYl2Je1HtD/B8QHLFtkjifvRUc7O0AHBS1\n3CtpfiLtGWY2Ia7vBewlaWLc7gmsT/CBlK3zPYLR+KqkiwgecO9LpDsrcS6cNoobBaeafJFYX0K4\nvypxrvVlYj3ztvsWoYQwxMyWKLjB7l5CWvmMRrbWzstEkroTvo63MbN3JJ2Zld8XibjJ56iUOvpr\nzOznJYTLpGVZ6Yrm8ypgOzP7cpmI4UWeT8vCrO1zzezvWfFPyKdT0hYEx4vfBw4Djs6hy2mjeJuC\nU0uM4OF2DUnbAkjqrTD3drn0IUyUs0Rh0p514v4FhGqgXDxKMCYouINem+DyONfLMntfxgDMi20g\npTReP054SSJpMMFjbBIjzJF7iKRVY7h+Cm6Os5lBqDYCmADsEsN2zdJyH3Bi05+QtsyhZS9Cg3wu\n7gW+G0sXSBoQteXUKekrQBczuwX4FaE6KcMahKo4pw3jJQWnpVie9bDDbJGkw4E/K8yv8ClhjoXk\nDFDZs0HlWv8XcHusknkGeCWmPy82Sr8A3EWoo8/E+StwSYyzGPhO1JNr9qllts3sQ0mXEaqMZhOm\nOS12Dv4KXCPpJYLxeQn4KCvdVyT9ErgvNjAvIvSeeisrzccIE7k8a2bvxgbg8YR2m4mJcCcCF0ua\nRHieH47pnQX8W9IRMd5sggHtk/yvZna/pE2A8bF0sQD4dgGdnwNXJRrHkzOzbU3CQDltE3ed7ThV\nIr4ou5rZFwrjC+4HNrQwM1e5aX2V0Pi+b4VaVgCWxJLVCOBiC1Ns1oRYErvAzA4oGtipa7yk4DjV\noyfwv1jFI+D4SgwCgJlNk7RA0npW2ViFtYEbo6H6EjimEh1l8H3g/9U4D6cV8JKC4ziO04Q3NDuO\n4zhNuFFwHMdxmnCj4DiO4zThRsFxHMdpwo2C4ziO04QbBcdxHKeJ/w9iZRqLvXmZjgAAAABJRU5E\nrkJggg==\n",
"text": [
"<matplotlib.figure.Figure at 0x109702a50>"
]
},
{
"metadata": {},
"output_type": "display_data",
"png": "iVBORw0KGgoAAAANSUhEUgAAAZYAAAEoCAYAAAB7ONeTAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJzsnXm8XPP5x98fWQgSEZSIkJbYqVgiFcu1NqJiaYtQ+5Ii\n9iWoX3tV1b7TRO17UKqhtjRcS2gIiSALQRAkVBOJrbI8vz++30lOxszcuffOzJl77/N+vc5rzjnf\n7Tnfc+Y857s9j8wMx3EcxykVS6UtgOM4jtOycMXiOI7jlBRXLI7jOE5JccXiOI7jlBRXLI7jOE5J\nccXiOI7jlBRXLE1E0lBJ5yaOj5U0U9IcSStK6ivpHUlzJQ1IU9bGImmapJ3TlqNcxHvTo5Fp6yQd\nWVqJmoakHpIWSmoW/29Jh0l6PnG86H5k/7/KUPaasTyVq4w85S6U9JNKlllJmsWDlxbxhfpNVBKz\nJI2WNCj5EJrZsWb2pxi/HXA5sLOZdTKzWcAfgWvMrKOZjUjnSpqMxa3BxDrcqcTylJR4b6Y1NjmN\nrJvmQBqKM3k/kv+vUpD9PJrZh7G8FnsP08AVS2EM+IWZdQLWBC4ChgA354m/GrAMMClxbk1gYmMK\nl9SmMemqDAMq+jXolJQGvXAltS2XIEWUXcz/xZ/HSmBmvuXZgPeBnbLObQUsADaMx7cB5wM9ga+B\nhcBcYBQwNcb9BpgDtANWICimT4DpMe1SMa/DgNHAFcB/CK2d9sBlwAfADGAosEyMXxPzOBWYGfM8\nLCFrB0ILahowG3g+kbYP8CIwCxgP7FBPPZwFvAX8F7gFWDoR/ouYx6wo/ybx/J2J658LnBHr69QY\n3i3W13HxeG3gi/ryjWGrAw8CnwHvASckwmqB+4HbY72/CWxR4PoWAj9J3M/rgUdj2n9nwmL4rsDk\nWJ/XAnXAkYnwIwgfEv8FngDWzCrnBOBd4HPgEkANSDsIeDvWx3WJsKXiM/J5zPv4GD/zXNX3zL0A\nXBrLfQ/oF8MuAOYD38b7d02OuusRyzqC8IzWxfMPAJ/GenqW+H+JYSsBI4AvgTFRnucL3I/z89y3\nw/jh/+UnwNPx+HPgLmCFPM/j6Qn5M/WxepTtC+Ad4Kg8ZW8dry95//YBXo/7vYGX4r36JD4r7fJc\nYx1LPkOHZdXH+sDIKNNk4NeJsP6E/+WceG9PS/u9aWauWApWTg7FEs9/AAyK+7cCf4z7ayUf0lx5\nAH8nKIcOwCrxj3VM4oGaR3gxLEVo/VwJPAx0BpaPD/2fY/yaGL8WaAPsTlBumT/S9fFP1jXm14eg\nqLrFP17mBbJLPF45Tz1MAybEdCsSXkTnx7BeBKW2FeFL8JB4ze3yXP/hwIi4fyBB+Q6Px0cAf68v\n33gtrwLnAm2BHxNeqLvFtLWEl2G/mPbPwEsF7nP2i+w/wJaxTu8C7o1hK8c/8L4x7ORY/0fE8L0I\nL6P1ooy/A0ZnlTMq3svuwBTiC6XItCOATjHtZ8DPY9hvCa3kzP15hvACzbws63vmvgeOjHX1W+Dj\nRLnPZK4vT931iLLdFvNfOpHvcvF+XQmMS6QZHrcOwEaEF+Jzee7Hov9XjrIP44f/l7WBnWO5KxOU\n2pUF/o8Z+TN19RxwHeF/8tNYzzvmKX8qsEvi+AHgzLi/OUG5LEV4L0wETspzjUvUMQnFEuvwI+DQ\nmNdmBIW5fgz/FOgb91cAeqX93jRzxVK4cvIrlpeAs+P+rSx+yS7xkGbnAawKfEdsNcRzA4GnEw/U\nB4kwAV+x5Bfzz4D34n4N4esrWd7MxAP9DYmv/EScIcAdWeeeAA4pUA/HJI53B6bG/aHZf3zCV9V2\nueow/vH/G69tKHAM8FEMux04uZ58tyd8LX6QFXY2cEvcrwWeSoRtCHxT4D5nv8j+mnWtk+L+IcCL\nWWk/YrFieZwlXxBLERR990Q5uyXCjwX+1YC02yTC72PxS+zprPuza4y/VJHP3DuJsGVj2h/F42dI\nfE3nqLseMX6PAnE6xzgdCQr5e2DdRPgF5G+xLPp/5cj3sOznIEecvYHX8v2nE/IvRVDY84HlEuF/\nBm7Nk/f5wM1xvyPhv9o9T9yTgYfyXGMhxbI/CaUbz90A/D7uf0D4D3UqVA+V3nyMpXGsQXg5NpS1\nCF9Sn8bJALOAYYSvyAwfJfZXIfzRX03Ef5zwJZbhCzNbmDj+htCyWZnwBfduHjl+nckz5tuXMEaU\nj6RcHxK6DDJ5nZaV1xqJ8CUws3cJL8zNgO0IXU6fSFqXoDSerSffrjFs9ayws4EfJYqamVUnyzRg\nllQy7beE+iRe0/SsuMl6WQu4OiHTF/F8tzzxs+uxvrQzEvuZ+wyhTrLzTcpU3zO3KF8z+ybuLp8I\nN+pnUfmSlpJ0kaSpkr4kvMyN8EyuQmhl5pO3oSTzQdKqkoZLmh7LvpPQ9VYMqwP/NbOvs2Trlif+\nPcC+ktoTWrGvmtlHUY51JT0q6dMoxwUNkCPJWsDWWc/6gYQPBoBfErrDpsWJFn0aUUbJSW2grbki\naSvCA/hCI5J/BPwPWClLGSRJ/on/Q3ixbWhmnzawrP8QvlTXIXRjJfkQuNPMjmlAfmtm7X+cyOsC\nM/tznnS5XkrPAr8mdJd9IulZwlfaioQxlYL5xj/P+2a2bgPKLAWfELqsMnKI8JWb4UPC1/W9BfJY\nk8WTO7Lrsb60+fiUH96fDMU8c4Uoti6T8Q4CBhBmR34gqTOLW6mfE1oFaxK6ArPlbap8fyZ0A25s\nZrMl7U0Y38gXP8knQBdJy5vZVwnZsj8mQkZmkyR9QGjVHkhQNBmGErpr9zezryWdTFACufia0OWV\nIfmB9yHwrJntlkeGscDeceLCCYSxxabUZ0nwFkv9CEBSJ0m/AO4lvJTfSoYXQ1QOTwFXSOoYv+zW\nlrR9nvgLgRuBqyStEuXoJinnQ5Yj7S2xrK6S2kj6Wfy6ugvYU9Ju8fwykmok5fsyE3B8LLsLof//\nvhh2I/BbSb0VWE7SHpIyX7wzCd1fSZ4FBhP6syEMXg4mNP8zf/xC+b4MzJV0pqQO8Ro2lrRlQt7G\nUijtY8BGkvaJs59OZMmXwDDgHEkbAkhaQdKvs/I4XVJnSd1j+vsakDZbzoys9wMnxvuzImGiBdDw\nZy4Hue5ffSxPUGb/lbQc4WWfkWcB8BBQG+/dhoTxg3w09F4uT3hRz4nP8xlZ4XmvJ7Y2XgQulLS0\npE0J4353FSjvHkI313aEMZakHHOBbyStT+j2zMd4Qsung6R1CONdGf4JrCvpN5LaxW0rSevH/YMk\nrRDrdS5BqaaOK5b6eUTSHMKXw9mEWVaHJ8KNJb+C6vvCO4QwMJiZ/fMAi19O2XlBGA+ZCvw7NqlH\nAskv9ULlnQ68AbxC6Fq5kDAeM53w5X0OYXDyQ+A08j8PBtxNeEG9Sxhk/hOAmb0KHE0Y8PxvDDsk\nkfZC4NzYjD81nnuO8MfLKJbRhIHczHHBfKPS/AWhO+09wlfwXwkD2xl5s+ulUD1l37+cac3sP4SW\n1kWEFuE6JFquZvYwcDEwPN6rN4CfZ+X1D8KX7DhCN+AtRabNJVNSCT8JvA6MJcyWS8Zv6DOXPL4a\n+JWk/0q6itxkp7+D0Pf/MWFG3ktZcQYT7v8MwvXfQv7/UC75CoWdRxg4/xJ4hB/WRa7nMRk+kDDu\n8glBAf7ezJ7OUz6ED83tgVFmluweP53QiplDeDaH57iuDFcSxp1mEsaU7mLxMzcX2A04gFCfn8Zr\naB/T/gZ4Pz4zxxBai6mjxR+IjuOUE0kLgXXM7L20ZXGccuItFsdxHKekuGJxnMrh3QNOq8C7whzH\ncZyS4i0Wx3Ecp6S4YmnlSDpb0o0lyGeRGXCV0dR5nF75ZDnyLheSaiXd2YT0F0o6qZQyVQrVY8Jf\n0psNmPrckHJLmq+kTSWNLlV+LR3vCnNKQjlmPCn45HgPaNvIxX1VgaQ/EOrm4EakXYUwNXltM/tf\nyYUrMy3lHgJI+icw1MweTVuWasdbLE5zoLmbOW+K/IcB/2yOSqUFcjfBwrRTD65YmimShkh6IOvc\n1ZKujvuHSXpXwUnZe5IOzJPPom6aRLfFIZI+kPS5pHMScZeSdI6CDag5ksbmWq0v6TZJ58f9GgW7\nTacqeNb8RNJhibh7SBon6UtJH8av+wyZBZOzY3l99ENvg9tIekXSbEkvS/pZIqxO0h8lvRDTPykp\np72muBr+UUmfxcWAjySvrb68EnX2H0nnqoCDs3gdL8ZFeuMl7ZArXqQfi+2nEeWam9gWSDqkiLpY\nXdIISV8oeDQ9KhFWK+kBSXfGa5sgqWfsJp0Zr2vXRPwVJN0c7+V0SednurriM3JZfHbeBfYocG1L\nON6Kctwv6fYox5uStsiTbqikS7PO/UPBdMoSXk8VOCs+t/+RdJ+ChQJiWafG/W4Kz/9x8XhtSV8k\ningW2FnBoZ9TiLStYPrWuI1gD+hrYPl43IawWrg3we7Ql0DPGLYqCX8YWfn8gWCiBhZber0BWBrY\nlGBvbL0YfgbB7lgm302BLnE/p6lz6jftvwOwUdzfhLAae694vBY/tBZ9GIstv3Yh+Ls4iPCRdABh\nZfmKMbyOsGJ/HYJBzmeAC/PUQxeCP41lCKvC7yea8K8vL4L15LnANgSDj5cSVlJnrFrXJuq4oS4L\nPiOPL5lYl9NjnvXVRV5z8Cx2M7BrvEe3E1wlnB2PjyJa1I7xC5nhL2jCP8c1LLI2TAPcHRBMqHyY\nOF6RYJhztRz5nkQw1bJ6vD/DgHtiWFFuHBLlfEmwQ5b6O6Cat9QF8K0JNy847jo47u/KYlP2y8WX\nzL5Ah3rySL70ehBe5KsnwscA+8X9KcCeefLJaeqcAqb98+RzFXBFljz5FMvBwL+z0r8IHBr3nwHO\nSYQdCzxeZN1uRrB0S315Ab8H7k6EdSDYysqlWBrqsmAJE/OJ8+vGetymvrqgHnPwUb4nE2F7EhRl\nZgy2Y7wPnajfDH9eE/55ri9bsRTl7oCgeD5gsXuGo4kuCHLkO5ElTeV3jfW6FEW6cUiknQ5sW47/\nc0vavCuseXMP4U8N4WvrbgALZr/3J3w9fhK7eNZrQL75zLOvQW4z/PWRz7Q/kraW9EzsgppN6MNu\niJnzbJPrH7Ckyf7ktSRN4C+BpGUl3RC7UL4kdHusICk5PpIvryXM6ZvZtyw2e5/NWjTMZcEswos9\nKesKBJtjvzOzFxMy5KuLrtRvDv6zrGv7j8U3aTyGcL1rUdgMfyET/sVQlLuDKNtwcjz/OegB/D0h\n70SCol3VinfjkKEjwSumUwBXLM2bvwEZq8R7kzDbbWZPWTC1vRrBQVa+KcUNmRb4EaErqBiKzfce\ngofMNcysM+EllXku68vjY8KLLslaLDZF3xBOI7QCepvZCoQuuqQF4UJ8QlC6AEjqQH7lmHFZsGJi\n62hml+SJP4HgVTKT91KEOhtlZjcl4hWqi0Xm4BNhec3B10PSDH9G/hXMbJMYXsiEf6m5l2Agcy1C\nF/CDeeJ9SOh6TNb5srbYFcUSbhzi8WEs6caB+D9rz2Jz/04eXLE0Y8zsc0Lf/22EPvApAJJ+JGkv\nBZPl8whfZPnMaTdkxtJNwPmS1okDopsqmNHPlWex+S4PzDKz7yX1Jnx5ZhTK54RulHxm2x8nmBQf\nKKmtpP0J/sGT00EbIse3wJfxmv6QI06+vB4kuCHIuCWoLRC3oS4LHiMouQwXEJy/nZwjXs66sGDN\nuqHm4HNi9Zvhz2vCv9SY2XjC+NRNwBNmNidP1GHAnyWtCWEKt6QBifBi3DhAuA+jzGxe6a6iZeKK\npflzD8HHd9LJ0FLAKYSv1S8ITfx8/iCM4s3+X0F4cTxFGMS8kTCQnZ2uIXkeB/xRwTXB/7HYPwkW\nvBleAIxWmKm1dTJvM/uCYD7/NMIL5nTgF7ak+fJCciW5ijA28h/CS/jxHHFz5mXBN88JhK6ZTwjj\nE58Rvuyz4zbUZcEdQH9JmXo+gOCaeZYWzwwbGK+5UF0UMgefq14KHRcyw1+fCf9CFCNHNvcAO7Hk\n85/N1cAI4Kn4nL1EaOFkqNeNQ+QggpJy6qHsCyQl9SP8adsAN5nZxTniXEOY4fINcJiZjUuEtSE8\noNPNbM94rgvhBbQWYfbKfmbm/Z5OVRC7nGYRFkV+UIL8LgA+M7Ormyyc0yhiK2+omfVNW5bmQFlb\nLFEpXEeYPrghMFDSBllx+hP+gD0JMzKGZmVzEuHLKKkBzwJGWnBNO4oyNrcdpxgk7RknACwHXAZM\nKIVSATCz37lSSRczm+BKpXjK3RXWmzAFdlrslxxOwmd4ZABhWh9mNgboLGlVAElrAP0JfajKlSb+\n7l22K3Cc4hhA6Hr8mDAmdEC64jhOepRbsXRjyamHmcVcxca5krAoL9vG0KpmlpmWOJMwt95xUsPM\njo6zjTqb2a5m9k7aMjlOWpRbsRQ7gJM9g0aSfkHoVx6XI3xxAWGQqLwDRY7jOE7RtC1z/h8TVv1m\n6M4P585nx1kjnvslMCCOwSwDdJJ0h5kdAsyUtJqZzZDUlSUXdwEgyZWN4zhOIzCzJhl+LXeLZSzQ\nU8G4YXvCavARWXFGEKYvIqkPMNvMZpjZOWbW3cx+TOivfjoqlUyaQ+P+oYQFdj8gbbMG1bL94Q9/\nSF2Gatm8LrwuvC4Kb6WgrC0WM5svaTBhXnsb4GYzmyRpUAy/wcwek9Rf0lTCQr7D82WX2L8IuF/S\nkcTpxmW7CMdxHKdBlLsrDDN7nLDYLHnuhqzjwfXk8SwJmz0WFn3tUkIxHcdxnBLhK+9bATU1NWmL\nUDV4XSzG62IxXhelpcW6JpZkLfXaHMdxyoUkrMoH7x3HcZxWhisWx3Ecp6S4YnEcx3FKiisWx3Ec\np6S4YnEcx3FKiisWx3Ecp6S4YnEcx3FKiisWx3Ecp6S4YnEcx3FKiisWx3Ecp6TUa4RS0ubAQGB7\noAfByvAHwHPAPRYccTmO4zgOUI+tMEmPAbMI/k9eBj4leHPsSvBnvyfQ2cz2KL+oDcNthTmO4zSc\nUtgKq0+xJH3L54vzIzP7gQfHtHHF4jiO03DKrlhyFLgMwc38/5pSaCVwxeI4jtNwSqFYCo6xSFoK\n2JswxrINYbBfkhYALwF3Aw/7G9xxHMfJUF9X2HPA84QxlvGZloqkpYFewABgWzPbvgKyNghvsTiO\n4zScSoyxLF1ft1cxcdLAFYvjOE7DqYSjr+Ukdcm3ARShePpJmizpHUlD8sS5Joa/LqlXPLeMpDGS\nxkuaKOnCRPxaSdMljYtbvwZet+M4jlMm6lvH8hph3Uo+flwosaQ2wHXALsDHwCuSRpjZpESc/sA6\nZtZT0tbAUKCPmX0naUcz+0ZSW+AFSX3NbHSU6Qozu6LeK3Qcx3EqSkHFYmY9mph/b2CqmU0DkDQc\n2AuYlIgzALg9ljdGUufMNGcz+ybGaQ+0IaypydCkplpz49NPoWvXtKVwHMepn6JMukjaR1LnxHFn\nSXsXkbQb8FHieHo8V1+cNWI5bSSNB2YCz5jZxES8E2LX2c1J2Voi8+ZBnz7w8stpS+I4jlM/xdoK\nqzWz2ZmDuF9bRLpiR8+zWx8Wy1lgZpsRFM32kmpi+FBCN9xmBGsAlxdZTrOkXTuorYWTTgKfj+A4\nTrVTr62wSK5upzZFpPsY6J447k5okRSKs0Y8twgz+1LSP4EtgbrkSn9JNwGP5Cq8trZ20X5NTQ01\nNTVFiFydHHooXH893HMPHHRQ2tI4jtNSqKuro66urqR5FrXyXtKthPGN6wlK5nhgRTM7rJ50bYEp\nwM7AJwR7YwNzDN4PNrP+kvoAV5lZH0krA/PNbLakDsCTwHlmNkpSVzP7NKY/BdjKzA7MKrvFTTce\nPRoOOAAmT4bllktbGsdxWiKVmG6c4QRgHnAfMBz4jqBcCmJm84HBBKUwEbjPzCZJGiRpUIzzGPCe\npKnADcBxMXlX4Ok4xjIGeMTMRsWwiyVNkPQ6sANwSpHX0azp2xe23RYuvjhtSRzHcfLTIFthzYmW\n2GIB+PBD6NULXnsN1lorbWkcx2lpVMwIpaT1gNMJ/lgy4zJmZjs1pfBy0lIVC4SB/EmT4L770pbE\ncZyWRiUVywTCTKzXgAXxtJnZq00pvJy0ZMXyzTew/vpw992w3XZpS+M4TkuikorlVTPboikFVZqW\nrFgA7r0XLr0UXnkF2hQzP89xHKcIKjl4/4ik4yV1zbYV5qTDAQfAssvCbbelLYnjOM6SFNtimUaO\nxY5mVtBWWJq09BYLwNixsOeeMGUKdOqUtjSO47QEKu5BsjnRGhQLwBFHwMorwyWXpC2J4zgtgUr4\nY9k5Lkj8JblbLA81pfBy0loUy4wZsPHG8NJL0LNn2tI4jtPcKbtrYmB7YBSwJ7ntflWtYmktrLYa\nnHEGnH46/OMfaUvjOI7jXWEtgv/9DzbcEIYNg113TVsax3GaM2WfFSbpsGjvK194e0mHN0UAp+ks\nvTRcfjmcfDLMn5+2NI7jtHbqm268PMHr472STpV0oKSDJJ0m6V6CDa8O5RfTqY+99grdYsOGpS2J\n4zitnXq7wiQJ6AtsC6wZT38AvAC8WK39Ta2pKyzDG2/AzjsHcy8rrZS2NI7jNEd8unEBWqNiATju\nuLAS/9pr05bEcZzmSCWmGydfT8Zih18ZD48nNqXwctJaFct//gMbbAB1dbDRRmlL4zhOc6MSJl1e\njdvSwObA28A7QC+gfVMKdsrDyivDuee6G2PHcdKjWJMuY4BtzWxePG4HvGBmW5dZvkbTWlssAPPm\nwWabwfnnw777pi2N4zjNiUoaoewMJK1RdYznnCqkXbswxnLqqcHEvuM4TiUpVrFcBLwm6XZJtxP8\nslxYPrGcprLTTtC7t7sxdhyn8hQ9K0xSV2BrwsD9GDObUU7Bmkpr7grLkHFjPHYs/Lhq7VA7jlNN\nVLIrDOA74FNgNrCupO2LSSSpn6TJkt6RNCRPnGti+OuSesVzy0gaI2m8pImSLkzE7yJppKS3JT0l\nybvlcrDmmnDKKaFLzHEcp1IUpVgkHQ08BzwB1AJPxt/60rUBrgP6ARsCAyVtkBWnP7COmfUEjiG4\nQMbMvgN2NLPNgE2BHSX1jcnOAkaa2boEI5lnFXMdrZHTT4cJE+Cpp9KWxHGc1kKxLZaTgN7AB2a2\nI2G68ZdFpOsNTDWzaXFG2XBgr6w4A4DbAcxsDNBZ0qrxODP03B5oA8zKThN/9y7yOlodyywDV14J\nJ54I33+ftjSO47QGilUs35nZtxC6qMxsMrBeEem6AR8ljqfHc/XFWSOW1UbSeGAm8IyZTYxxVjWz\nmXF/JrBqkdfRKtlzzzDGcs01aUviOE5roD5/LBmmS1oReBgYKWkWMK2IdMWOnmcPFGVW9i8ANpO0\nAvCkpBozq1sioplJyllObW3tov2amhpqamqKFKdlIcHVV8M228BBB0HXrmlL5DhOtVBXV0ddXV1J\n82ywrTBJNYQ1LU+YWcHOFUl9gFoz6xePzwYWmtnFiTjDgDozGx6PJwM7JFokmXj/B3xjZpfHODVm\nNiPOVnvGzNbPit/qZ4VlM2QIfPop3HFH2pI4jlOtVGRWmKS28UUOgJnVmdmI+pRKZCzQU1IPSe2B\n/YERWXFGAIfEsvoAs81spqSVM7O9JHUAdgXGJ9IcGvcPJbSknHo491x4+ml48cW0JXEcpyVTr2Ix\ns/nAFElrNTTzmHYwYRbZROA+M5skaZCkQTHOY8B7kqYCNwDHxeRdgafjGMsY4BEzGxXDLgJ2lfQ2\nsFM8duqhY8ewYHLwYFiwIG1pHMdpqRRrK+x5wkywl4Gv42kzswFllK1JeFdYbsxg++3hN7+BQYPS\nlsZxnGqjYv5Y4rhKNmZmzzal8HLiiiU/48fDz38eHIJ16ZK2NI7jVBNV4+hL0ktm9rMmZ1RCXLEU\n5rjjwmyx669PWxLHcaqJalIs48ysV5MzKiGuWArzxRew4YZhRf5Pf5q2NI7jVAuVthXmtCBWWgnO\nOw9OOMEdgjmOU1pcsbRijj4avvoK7r03bUkcx2lJuGJpxbRpA9ddB2eeCXPmpC2N4zgthYb4Y1kN\n2IpgbuVlM/ssEbaJmb1RHhEbh4+xFM8RR0DnznDFFWlL4jhO2lRyuvF+wKVAZnrx9sAZZvZAUwov\nJ65Yiufzz2GjjWDUKNhkk7SlcRwnTSqpWCYAu2RaKZJWAUaZ2aZNKbycuGJpGEOHhrGWZ58N05Ad\nx2mdVHJWmIDPE8df8EOLxE4z5phj4Ntv4c4705bEcZzmTrEtlkuBnwL3EBTK/sAEMzuzvOI1Hm+x\nNJxXXoEBA8KK/M7u7NlxWiWV7AoTsC+wLWHw/nkz+3tTCi43rlgax6BB0L49XHtt2pI4jpMGVbPy\nvhpxxdI4vvgiDOQ/9hhsvnna0jiOU2nKPsYiaXT8/UrS3KzNVz60QFZaCS64INgSW7gwbWkcx2mO\neIvF+QELF0LfvnDUUXDkkWlL4zhOJanYrDBJP5grlOuc0zJYain4y1/gnHNC15jjOE5DKHa68cbJ\nA0ltgS1KL45TLfTqBfvtF5SL4zhOQ6hvjOUcSXOBTZLjK8Bn/NB3vdPCOP98eOQRePnltCVxHKc5\nUex044vM7KwKyFMyfIylNNx5J1x9NYwZE4xWOo7TsqnYGIuZnSVpRUm9JW2f2YoUsp+kyZLekTQk\nT5xrYvjrknrFc90lPSPpLUlvSjoxEb9W0nRJ4+LWrxhZnIbzm9/AssvCX/+atiSO4zQXim2xHA2c\nCHQHxgF9gJfMbKd60rUBpgC7AB8DrwADzWxSIk5/YLCZ9Ze0NXC1mfWJ1pRXM7PxkpYHXgX2MrPJ\nkv4AzDWzvPZ4vcVSOt54A3beGd58E370o7SlcRynnFTSVthJQG9gmpntCPQCviwiXW9gqplNM7N5\nwHBgr6w4A4DbAcxsDNBZ0qpmNsPMxsfzXwGTgG6JdG6rrEJssgkcfDCc1aw6Qx3HSYtiFct3ZvYt\ngKRlzGwysF4R6boBHyWOp7OkcsgXZ41kBEk9CMpsTOL0CbHr7GZJbtmqzNTWwlNPwfPPpy2J4zjV\nTtsi401MltMoAAAgAElEQVSXtCLwMDBS0ixgWhHpiu2Lym59LEoXu8H+BpwUWy4AQ4E/xv3zgcuB\nHyzlq62tXbRfU1NDTU1NkeI42XTsCFddBcceC6+9FuyJOY7T/Kmrq6Ourq6keTZ45b2kGqAT8ISZ\nfV9P3D5ArZn1i8dnAwvN7OJEnGFAnZkNj8eTgR3MbKakdsCjwONmdlWeMnoAj5jZJlnnfYylxJjB\nHnvA9tt7t5jjtFQqMsYiqW182QNgZnVmNqI+pRIZC/SU1ENSe4K5/ez1LyOAQ2JZfYDZUakIuBmY\nmK1UJHVNHO4DVJVb5JaKBNdfD5ddBu+/n7Y0juNUK/UqFjObD0yRtFZDM49pBwNPAhOB+8xskqRB\nkgbFOI8B70maCtwAHBeT9wV+A+yYY1rxxZImSHod2AE4paGyOY3jxz+G00+H448PLRjHcZxsip1u\n/Dxh8Pxl4Ot42sxsQBllaxLeFVY+5s0LJl9qa+FXv0pbGsdxSkklHX3V5DhtZvZsUwovJ65YyssL\nL8ABB8DEidCpU9rSOI5TKqrG0Zekl8zsZ03OqIS4Yik/Rx0VVuVfc03akjiOUyqqSbGMM7NeTc6o\nhLhiKT8Zb5OPPgpbbpm2NI7jlIJKrrx3nB+w0kpwySXw29/CggVpS+M4TrXgisVpEgcfHBZPXn99\n2pI4jlMteFeY02QmT4bttoPx46FbtsEex3GaFdXUFXZIifJxmiHrrx9MvZx8ctqSOI5TDRRssUga\nbWZ9JX3FD+1+mZlV7URTb7FUlm+/DVaQr7kG+vdPWxrHcRpL1cwKq0ZcsVSekSPhmGPgrbfCNGTH\ncZofZVcskjqZ2RxJXXKFm9l/m1J4OXHFkg4HHghrrQUXXpi2JI7jNIZKKJZ/mtkekqaRwwS+mf24\nKYWXE1cs6TBjRugSe+YZ2HjjtKVxHKeheFdYAVyxpMcNN8BttwWzL23apC2N4zgNoRItls0LJTaz\n15pSeDlxxZIeCxdCTQ3stx8MHpy2NI7jNIRKKJY6QhdYB2ALYEIM2hQYW232wZK4YkmXzNqW116D\n7t3TlsZxnGIp+zoWM6sxsx2BT4DNzWwLM9uCYEL/k6YU7LRs1l8fTjjB/bY4Tmuk2AWS65vZIi+N\nZvYmsEF5RHJaCmedBe++C3/7W9qSOI5TSYr1xzIc+Aq4CxBwILC8mQ0sr3iNx7vCqoMXXwzOwN56\nC1ZcMW1pHMepj0o6+uoAHAtsF089Bww1s++aUng5ccVSPRx/PHz/Pdx4Y9qSOI5THz7duACuWKqH\nOXOC35a77oIddkhbGsdxClExI5SS1pX0N0kTJb0ft/eKTNtP0mRJ70gakifONTH8dUm94rnukp6R\n9JakNyWdmIjfRdJISW9LekpS52JkcdKhUye47rpg7uW7qm3jOo5TKoodvL8VGAbMB3YEbgfuri+R\npDbAdUA/YENgoKQNsuL0B9Yxs57AMcDQGDQPOMXMNgL6AMdLWj+GnQWMNLN1gVHx2Kli9torrMj/\n05/SlsRxnHJTrGLpYGb/InSdTTOzWmCPItL1BqbGNPOA4cBeWXEGEBQVZjYG6CxpVTObYWbj4/mv\ngElAt+w08XfvIq/DSZFrr4W//hXeeKP+uI7jNF+KVSzfxdbHVEmDJe0LLFdEum7AR4nj6SxWDoXi\nrJGMIKkHYe3MmHhqVTObGfdnAqsWIYuTMl27hhbL0Ue7K2PHacm0LTLeScCywInA+UAn4NAi0hU7\nep49ULQonaTlgb8BJ8WWy5IRzUxSznJqa2sX7dfU1FBTU1OkOE65OOqoMIh//fVw4on1x3ccp7zU\n1dVRV1dX0jzrnRUWWyoXm9npDc5c6gPUmlm/eHw2sNDMLk7EGQbUmdnweDwZ2MHMZkpqBzwKPG5m\nVyXSTAZqzGyGpK7AM2a2Pgl8Vlj1MnkybLttMPey5pppS+M4TpKKzAozswXAtpIaU9BYoKekHpLa\nA/sDI7LijCC6No6KaHZUKgJuBiYmlUoiTabFdCjwcCNkc1Ji/fXhpJPguOPc3IvjtESKXSA5DFgd\neAD4Jp42M3uoiLS7A1cBbYCbzexCSYNiBjfEOJmZY18Dh5vZa5K2JSzEnMDirrGzzeyJ6HjsfmBN\nYBqwn5nNzirXWyxVzPffw+abw+9+BwOr1n6D47Q+Krny/jZyO/o6vCmFlxNXLNXPyy/DgAFhltgq\nq6QtjeM44CvvC+KKpXlwxhnw0UcwfHjakjiOAxVcee845eKPfwyD+A/7KJnjtBi8xeKkznPPhXGW\nN990C8iOkzbeFVYAVyzNi8GD4Ztv4JZb0pbEcVo3qXSFSXq0KQU6Ti4uvBCefhqefDJtSRzHaSqN\nGWPJNsniOE2mY8dgR+yYY2Du3LSlcRynKTRGsYwruRSOA+y2G+yyS3Bp7DhO88XHWJyqYvZs2Hhj\nuPtudwrmOGng042dFkfnzvCXv8CRR4bBfMdxmh/eYnGqkgMPhNVXh8suS1sSx2ldVNI18a+LOec4\npeLqq4N5/TFj6o/rOE51UWxX2DlFnnOckrDKKkG5HHEE/O9/aUvjOE5DKNgVFi0T9yeYux/OYodc\nHYENzax32SVsJN4V1vwxg333hQ03hAsuSFsax2kdlH3lvaSfElwC/xH4PxYrljkE51qzmlJ4OXHF\n0jKYMQN++lN49FHYaqu0pXGclk8lzea3M7N5TSmo0rhiaTkMH77YWOUyy6QtjeO0bCrRYnmjQFoz\ns02bUng5ccXScjCDX/8afvITuOSStKVxnJZNJRRLj0KJzWxaUwovJ65YWhaffw6bbgoPPgjbbJO2\nNI7Tcqka68aSXjKznzU5oxLiiqXl8eCDcPbZMH48LLts2tI4Tsukmlbee8+3U3Z++UvYckv43e/S\nlsRxnEKU3aSLpH6SJkt6R9KQPHGuieGvS+qVOH+LpJnZYz2SaiVNlzQubv3KfR1OdXDttXD//cE5\nmOM41UlZFYukNsB1QD9gQ2CgpA2y4vQH1jGznsAxwNBE8K0xbTYGXGFmveL2RFkuwKk6VloJhg2D\nww+Hr75KWxrHcXLRtsz59wamZgb5JQ0H9gImJeIMAG4HMLMxkjpLWs3MZpjZ8wUmEDSpD9Bpvuy5\nZxhvGTIErr8+bWkcJ13mzIF33oEpU+Dtt8P+t9/CQw+lJ1PRikXSasBWhNbCy2b2WSL4kDzJugEf\nJY6nA1sXEacbMKMekU6QdAgwFjjNzGbXE99pQVx1FWyySViZv/POaUvjOOVnwYKgNMaNC9trr8HE\nifDll7DuumFbbz34+c9h7bXTlbUoxSJpP+BS4Nl46jpJZ5jZAwBmlm+9S7HTsrJbH/WlG0qwBgBw\nPnA5cGR2pNra2kX7NTU11NTUFCmOU+107gw33hjM60+YAJ06pS2R45SWWbPgxRfh+edh9OgwG3KV\nVWDzzaFXLzjttOC7qFs3WKoJgxp1dXXU1dWVTG4ofuX9BGCXTCtF0irAqPoWSErqA9SaWb94fDaw\n0MwuTsQZBtSZ2fB4PBnYwcxmxuMewCNmtkmeMnKG+3Tj1sExx4QFlDfemLYkjtM0vvwSRo2Cf/0L\nXngB3n8ftt4att0W+vYNMyJXXLH8cpRiunGxXWECPk8cf0FxYxxjgZ7x5f8JwZjlwKw4I4DBwPCo\niGZnlEpeYaSuZvZpPNwHKGQhwGnBXHZZWDj5+OOw++5pS+M4xbNwYejOeuIJePLJ0CLp2xd23TVM\nTtlsM2jXLm0pG0exiuUJ4ElJ9xAUyv7A4/UlMrP5kgYDTwJtgJvNbJKkQTH8BjN7TFJ/SVOBr4HD\nM+kl3QvsAKwk6SPg92Z2K3CxpM0IXWbvA4OKvA6nhdGpE9xyCxxyCLz+epg15jjVyrx5UFcXJp88\n/DB06QL9+sG558L220OHDmlLWBqK7QoTsC+wLeFl/ryZ/b3MsjUJ7wprXZxyCnz8Mdx3H8jnCzpV\nxPffh1bJgw8GK909e4bFvvvsA+usk7Z0P6RqTLpUI65YWhfffhv6oM85Bw46KG1pnNaOWfB+escd\nYUHvhhvCfvvB3nvDGmukLV1hKmGEcrSZ9ZX0FT+cqWVmVrVzcVyxtD7GjQtTLceOhTXXTFsapzXy\n4Ydw++1w552h5XzIIeFDp0ePtCUrHm+xFMAVS+vkwgth5Mgws6YpUzAdp1gWLAiD78OGhWnBBxwA\nhx4aHNM1x27ZSjr6utPMDq7vXDXhiqV1smBBGAT91a/CuIvjlIuZM8PEkb/+FVZeGX7726BUllsu\nbcmaRiWnG2+cVXBbYIumFOw45aBNm9ANsfXWsNtusNFGaUvktDTefBMuvzzM6vrlL+GBB8L4nrOY\ngp0Fks6RNBfYRNLczAZ8Rlh/4jhVx09+ErrEfvObMCPHcZqKWehe7dcvrDPp2RPefRduusmVSi6K\n7Qq7yMzOqoA8JcO7wlo3ZrDXXsHkxZ//nLY0TnNlwYIwhf2SS8IalNNOC4PxSy+dtmTlo6KD95JW\nBHqScOplZlXrFcMVizNzZli9/MADwSyG4xTL/Plw773wpz8F+1znnBMsOzTHwfiGUrExFklHAycC\n3YFxQB/gJWCnphTuOOVk1VXDTJ3MqvyOHdOWyKl25s+Hu+8OCmX11WHoUNhxx9ahUEpJsV1hbxJM\n5r9kZptJWh+40Mz2KbeAjcVbLE6Go44KvzfdlK4cTvWycGFoofzhD9C9e/htrcbQKzkr7Dsz+1YS\nkpYxs8mS1mtKwY5TKa68MnSJ/eMfYdzFcTKYwVNPBadxHTqEj4/WqlBKSbGKZXocY3kYGClpFjCt\nbFI5Tgnp2DFMQd5337BobfXV05bIqQbGjg0KZfr0MItwn328y6tUNHjlvaQaoBPwhJlV7WRO7wpz\nsjnvvODn4sknfVV+a2batKBQXnghdHkdcQS0LbeT9mZEKbrC6v17SWobnW8BYGZ1ZjaimpWK4+Ti\nd78LxiqvuCJtSZw0+OaboEi22CIsnH377eAozpVK6alXsZjZfGCKpLUqII/jlI22beGuu8KahNde\nS1sap1KYBQvDG2wQlMn48fD73zd/0yvVTLG6ugvwlqSXCc64IFg3HlAesRynPPToAVdfDQMHBuXi\nL5eWzYQJcOKJMHt2GGfbfvu0JWodFDvduCbHaTOzZ0suUYnwMRanEIceCu3bw403pi2JUw6++ip0\ne911VxhbO/roYEfOqZ+qMZsv6SUz+1mTMyohrlicQsydC716wcUXB0OCTsvhn/+E448PrZPLLw8r\n553iqeQ6lvpYpv4ojlM9dOwYVlgPGAC9e4dFcU7z5tNP4aSTQhfnTTfBLrukLVHrpeyTLiX1kzRZ\n0juShuSJc00Mf11Sr8T5WyTNlPRGVvwukkZKelvSU5I6l/s6nJbH1luHF9EhhwRjg07zZOHC4BNl\n002DD/k33nClkjZlVSyS2gDXAf2ADYGBkjbIitMfWMfMegLHAEMTwbfGtNmcBYw0s3WBUfHYcRrM\nkCFh1tAll6QtidMYPvww+N256SZ4+ulgybpDh7SlcsrdYukNTDWzaWY2DxgOZBvVGADcDmBmY4DO\nklaLx88Ds3LkuyhN/N27DLI7rYCMY7CrroKXX05bGqdYzODmm8OalJ13hhdfhE02SVsqJ0OpxlgO\nyXO+G/BR4ng6sHURcboBMwqUt6qZzYz7M4FVixfVcZake3e4/no48MDQP9+pU9oSOYX4+OMwy2vG\njNBKcYVSfRRULJJGm1lfSV8B2VOszMw6xZ03fpg6xClSjuwZCEVP5zIzk5Qzfm1t7aL9mpoaaty6\nnJOHX/0qGCP87W/DoL7bjKo+zMK9OfXUMOvrnHOgXbu0pWr+1NXVUVdXV9I8SzLdOG/mUh+g1sz6\nxeOzgYVmdnEizjCgzsyGx+PJwA6ZFomkHsAjZrZJIs1koMbMZkjqCjxjZutnle3TjZ0G8c03YUD/\n5JPhyCPTlsZJMns2DBoEb70Fd9wBm2+etkQtl7LbCpPUKf52ybUVkf9YoKekHpLaA/sDI7LijCB2\npUVFNDvRzZWPEcChcf9QgtVlx2kSyy4b3NAOGRJeYE518MILwe3Bj34Er7ziSqU5ULDFIumfZraH\npGnk6J4ysx/XW4C0O3AV0Aa42cwulDQopr8hxsnMHPsaONzMXovn7wV2AFYCPgN+b2a3RqV2P7Am\nwXz/fmY2O6tcb7E4jeLWW+Gyy8JLbNll05am9TJ/fvDkOGxYsJCw555pS9Q6qJqV99WIKxansZiF\ntS1LL+1eJ9Pigw/goIPC1OHbb3cfOpWk7IpFUsFGZ6ZlUY24YnGawty5YSprbW2YLeZUjoceCpMo\nzjgDTjvNfedUmkooljpCF1gHYAtgQgzaFBhbbfbBkrhicZrK+PGw665hjUTPnmlL0/KZNy+Mbz30\nUDBz37t32hK1Tso+eG9mNWa2I/AJsLmZbWFmWwC94jnHabFstlmwjLvffvDdd2lL07L5+OPga37K\nFHj1VVcqzZ1iG5nrJ9eqmNmbwAYF4jtOi+DYY2HttUO3jFMe/vUv2HJL+MUv4JFHYKWV0pbIaSrF\n+mMZDnwF3EVYzHggsLyZDSyveI3Hu8KcUjF7djCxf/nlsO++aUvTcli4EC64AIYODQsfd9wxbYkc\nqOCsMEkdgGOB7eKp54ChZla1HQSuWJxSMmZMmO768svBC6XTNL78Msz6mjMHhg/3WV/VhE83LoAr\nFqfUXHFFeAk+/3yYiuw0jilTYK+9wsSIK65wsyzVRiVbLOsCfyaYvs8YpTYz+0lTCi8nrlicUmMW\nvE127RqMVjoN57HH4LDDgnn7o45KWxonF2WfFZbgVmAYMB/YkWCq/u6mFOw4zQ0prMp/8km45560\npWlemMGFFwarxA8/7EqlpVNsi+U1M9tc0hsZY5CZc2WXsJF4i8UpF5n1Lc8+CxtumLY01c/XXwej\nnu+9F9aorLFG2hI5hahki+W76A1yqqTBkvYFlmtKwY7TXNlsM7j44mBq/6uv0pamupk+HbbbDtq3\nh+eec6XSWii2xbIVMBnoDJwPdAIuMbN/l1e8xuMtFqfcHHkkfPut+2/Jx7hxMGAADB4MZ57pddRc\nqMjgfWypXGxmpzeloErjisUpN99+C336BD8hxx2XtjTVxSOPwBFHhDUqv/pV2tI4DaEUiqVe18Rm\ntkDStvI3teMsQYcO8OCDsM02YeW4myEJXHMNXHQRPPpocJzmtD6K7QobBqwOPAB8E0+bmT1URtma\nhOtBp1L8/e9wyinBxlVrNkcyf36oh6efhn/+0xeSNlcquY7lNnI7+jq8KYWXE1csTiU5/XSYNCl0\nAbVGM+9z58IBB8D338MDD0DnzmlL5DQWX3lfAFcsTiWZNw922gl22w3+7//SlqayfPYZ9O8fZssN\nHeor6Zs7lfB5Xytp1QLhXSWd1xQBHKcl0K4d3HdfcKP72GNpS1M53n0X+vYNiuXGG12pOIH6Bu/H\nAsMltQdeAz4lWDdeDdgc+B9wWVkldJxmwuqrBwdV++4Lo0fDOuukLVF5GTcO9tgjtNCOPTZtaZxq\noj5HX49GR18HAKMJJl3mAS8A+5vZTmZW8PtMUj9JkyW9I2lInjjXxPDXJfWqL21sSU2XNC5u/Yq/\nZMcpH337BnfG++zTshdPjhoFP/85XHutKxXnhzR4jCWua1nOzOYUGXcKsAvwMfAKMNDMJiXi9AcG\nm1l/SVsDV5tZn0JpJf0BmGtmVxQo28dYnFQwC4snv/46WENuaQsD77sPTjghDNLvsEPa0jilpmIm\nXSTdK6mTpOWAN4BJks4sImlvYKqZTTOzecBwYK+sOAMIRi0xszFAZ0mrFZG2hf1dnZaCBH/5S7CN\ndVkL6yi+9lo47bTg9dGVipOPYidGbhhbKHsDjwM9gIOLSNcN+ChxPD2eKybO6vWkPSF2nd0sySc3\nOlXFMssEg4tXXAEjR6YtTdMxg/POC4rlhRdg003TlsipZupdeZ+JJ6kdQbFcb2bzJBXTz1RsX1RD\nWx9DgT/G/fOBy4EjsyPV1tYu2q+pqaGmpqaBxThO4+neHe69N6zv+Pe/m++CQTM44wx46qng5GzV\nvPNEneZIXV0ddXV1Jc2zWMVyAzANmAA8J6kH8GUR6T4GuieOuxNaHoXirBHjtMuX1sw+y5yUdBPw\nSK7Ck4rFcdKgpgaGDAmD+aNHw7LLpi1Rw1iwINhBGz8e6uqgS5e0JXJKTfZH93nnNX0FSVFdYWZ2\njZl1M7PdzWwh8AGwUxFJxwI9JfWIU5b3B0ZkxRkBHAIgqQ8w28xmFkorqWsi/T6EcR/HqUpOPjn4\nbRk0KHz9NxfmzYNDDoG33w5jKq5UnGIpdvB+ZUnXxqm9rwFXEUznF8TM5gODgSeBicB9cVbXIEmD\nYpzHgPckTSW0jI4rlDZmfbGkCZJeB3YATin+kh2nskhh8eAbb8DVV6ctTXH873/w61/DrFlhwWfH\njmlL5DQnirUV9i/gWeAuwnjIgUCNme1SXvEaj083dqqNadPgZz+D228Ppl+qla+/Dl13K6wQfM20\nb5+2RE4lqaQRyjfNbOOsc4vcFFcjrlicauS554J/kuefh/XWS1uaHzJ3bjDPss46oZXVtthRWKfF\nUEnXxE9JGihpqbjtDzzVlIIdpzWy/fZwwQXBs+KsWWlLsyRz58Luu4fxoJtvdqXiNJ6CLRZJX7F4\nyvBywMK4vxTwtZlVbc+rt1icauakk2Dy5OC3pBpe4HPmBKWyySZhcWdrNP3vBNxsfgFcsTjVzPz5\noctpo43gyivTlWXOHOjXL5i9v+46VyqtnYq4Jk4UtiLQE1gmc87MnmtK4Y7TWmnbNtjc2nrr0Eo4\n4oh05Pjyy6BUNt88KJWWZtfMSYdiB++PBk4kLFIcB/QBXjKzYtaypIK3WJzmwOTJYdzloYdg220r\nW/aXXwYLxVttFfzUu1JxoLKD9ycRjEJOi2b0e1HcynvHcQqw/vpwxx1hzcgHH1Su3Nmzw5Tn3r1d\nqTilp1jF8p2ZfQsgaRkzmwxU4WRJx2l+9OsHZ54ZZopVwodLZkylT5+wYNOVilNqilUsH8UxloeB\nkZJGEGyHOY5TAk4+GbbcEg48MNjnKhfffAO/+AX06gVXXeVKxSkPjXH0VUMw5/KEmX1fDqFKgY+x\nOM2N779fPOX3qqtKn/9338Gee0K3bnDLLT77y8mNTzcugCsWpzkyaxZssw0cfzwMHly6fL//Hvbd\nF5ZfPphpadOmdHk7LYuKTjd2HKf8rLhiWDTZty/8+Mewxx5Nz3P+fDjooDDF+c47Xak45cdbLI5T\nhbz0UhjMHzkyLFxsLAsXwqGHwmefwYgRsPTSpZPRaZlUcrqx4zgV5Gc/g+uvD8rl448bl4cZ/Pa3\n8NFH8Pe/u1JxKod3hTlOlbLffvDuu2HA/bnnwvhIsZjBqacGHzBPPdX8PFc6zRvvCnOcKsYMjjoK\nPv88tDqKHR+54AK4//7gTnjFFcsqotPC8K4wx2nhSDBsWHC+deqpxaUZNixMJ37iCVcqTjq4YnGc\nKqddO3jwQRg1Cq64onDc+++H888P3V9du1ZGPsfJpuyKRVI/SZMlvSNpSJ4418Tw1yX1qi+tpC6S\nRkp6W9JTkjqX+zocJ006d4bHHw8m9ocPzx1n5Eg44YTgo37ttSsrn+MkKatikdQGuA7oB2wIDJS0\nQVac/sA6ZtYTOAYYWkTas4CRZrYuMCoeO3moq6tLW4SqoTnXRffuQWmceGIYO0kyZkwwB/O3v8FP\nf1pcfs25LkqN10VpKXeLpTcw1cymmdk8YDiwV1acAcDtAGY2BugsabV60i5KE3/3Lu9lNG/8T7OY\n5l4Xm2wS/Ljst1+Y8QUwcSLstRfceitst13xeTX3uiglXhelpdyKpRvwUeJ4ejxXTJzVC6Rd1cxm\nxv2ZwKqlEthxqp0ddwxWiffYA0aPDpaKL700GJd0nGqg3OtYip3vW8zUNuXKz8xMks8rdloVAwfC\n9OnBOdiVV8LBB6ctkeMkMLOybQRPk08kjs8GhmTFGQYckDieTGiB5E0b46wW97sCk3OUbb755ptv\nvjV8a+q7v9wtlrFAT0k9gE+A/YGBWXFGAIOB4ZL6ALPNbKakLwqkHQEcClwcfx/OLripC3wcx3Gc\nxlFWxWJm8yUNBp4E2gA3m9kkSYNi+A1m9pik/pKmAl8DhxdKG7O+CLhf0pEEh2P7lfM6HMdxnOJp\nsSZdHMdxnHRolivvJZ0t6S1Jb0i6R9LSku6TNC5u70salyftNEkTYryXKy17qclTF70lvRyv8RVJ\nW+VJW+/i1eZEE+uiNTwXP5X0UrzOEZI65knbGp6LYuuipT0XJ8V6eFPSSfFcUQvOG/RclHPwvkwT\nAnoA7wFLx+P7gEOz4lwGnJsn/ftAl7Svo5x1ATwD/Dye2x14JkfaNsDUmEc7YDywQdrXlEZdtKLn\n4mVgu3jucOCPrfi5qLcuWuBzsTHwBrBMvM8jgbWBS4AzY5whwEVNfS6aY4tlDjAPWFZSW2BZYJHH\nCkkijLncWyCPljKwn6suPgFmACvEOJ1J1E+CYhavNieaUhcZWvpzsa6ZPR/j/Av4ZY60reW5KKYu\nMrSU52J9YIyZfWdmC4BnCdddzILzBj0XzU6xmNl/gcuBDwkPyGwz+1ciynbATDN7N18WwL8kjZV0\ndHmlLS956mIkwcTN5ZI+BC4lTNXOppjFq82GJtYFtI7n4i1JmZfBr4HuOZK3lueimLqAFvRcAG8C\n28Wur2WB/sAaFLfgvEHPRbNTLJLWBk4mNMlWB5aXdFAiykDgngJZ9DWzXoRukeMlNcAIRnVRoC5u\nBk40szWBU4BbciRvUbM2mlgX0DqeiyOA4ySNBZYHvs+RvLU8F8XUBbSg58LMJhOWaDwFPE7ozlqQ\nFSezluUHyRtSVrNTLMCWwItm9oWZzQceArYBiE3dfQj9qDkxs0/j7+fA3wlNvOZKrrroC/Q2s7/H\nOPwW2AwAAAcKSURBVH8j9zV+zJJfad0JXyHNlabURWt4LrYxsylm9nMz25LQlZGrVd8anoti66Kl\nPReY2S1mtqWZ7QDMAt4GZirYZ0RSV+CzHEkb9Fw0R8UyGegjqUMcT9kFmBjDdgEmmdknuRJKWjYz\n+0PScsBuhMGs5kp2XexMqIt3JO0Q4+xEeHiyWbR4VVJ7wgLUEZUQukw0ui5awXOxCzBR0ioAkpYC\nziVaEs+ipT8XRddFC3wukPSj+LsmsC+hdyez4BzyLDinoc9F2jMVGjm74UzgLcJNvh1oF8/fChyT\nFXd14J9x/yeE5t94Qn/j2WlfSznqgvCVNiZe50tAr+y6iMe7A1MIsz1abV20kueiPXBSvN9TgD8n\n4ra256Koumihz8VzsS7GAzvGc10IExjeJnSTdW7qc+ELJB3HcZyS0hy7whzHcZwqxhWL4ziOU1Jc\nsTiO4zglxRWL4ziOU1JcsTiO4zglxRWL4ziOU1JcsTgVRdLoRqarkfRI3N+zsebcJa0g6djE8eqS\nHmhMXuUkmmvv0sA090UTJtnnD5N0bemkazqSrmjO5lGcwrhicSqKmfUtQR6PmNnFjUy+InBcIq9P\nzOzXTZWpDDRogZmkdYDlLL/x1SajSImyGwqcUaK8nCrDFYtTUSR9FX9rJNVJekDSJEl3JeJsJWm0\npPGSxkhaPiuPRV/gkm6TdHWM/66kX8bzy0v6l6RXFRw1DYjJLwLWVnDcdLGktSS9GdMsI+nWGP81\nSTWJ8h6S9LiCM6ScSk3S/yk4FXtD0g2J83WSLorXMkXStvH8spLuV3BC9ZCkf0vaPEe+v4lpx0ka\nFs2QZHMACRMbkg6PZY0h2tKL51eR9Lco58uStkmcH6ngAOrGTIspmvCYIul2wsr17pLOiGlfl1Rb\nSE5JbeI9eiPW68kAZvYO0EN5nEo5zZy0TQz41ro2YG78rQFmE8xGCHiR8AJsTzAIuEWMtzzByVAN\n8Eg8dxhwbdy/Dbgv7m8AvBP32wAd4/7KifNrAW8k5OmROQZOA26K++sBHwBLx/LeBTrG42lAtxzX\ntmJi/w7gF3H/GeDSuL87MDLunw4MjfsbEfyGbB6P3yeY2tiAoDDaxPN/AQ7OUfbjibRdo+wrEcza\nvABcE8PuIVjsBVgTmBj3rwOGxP2fAwtj+T0IFnB7x7DdgBvi/lLAIwRXFdlyXg8cDGwOPJWQc4XE\n/u3A7mk/k76VfmuL46THyxYNhkoaD/wYmAt8amavAphZpoWTLw8jGs0zs0mSMr4klgIujP34C4HV\nowG+Ql05fYFrYl5TJH0ArBvLGGVmc6MsEwkv3GynYTtJOoPgTKoLwb7UozHsofj7WkybKe+qWN5b\nkiZk5ZcxprkFMDbWQQeC87Js1gI+jftbEzxlfhHlvS9eBwQjjBsk6rOjgoHFvkQHT2b2pKRZibw/\nMLOMW97dgP9v735CbArDOI5/fzGiaWbBalIWQ1MUkz8lK1mYlLJSFkiRGqVZSyQrJSuaUSw0CykL\nG1FmooRmQ5KmsSIKQ03RRMaYHovnvebMce78uR2G2/NZnXvnvO997mnuee77575vhya3/m4EVgHt\nBXF+wBNPq6TzwC18LaqKd5lrEepIJJYwn8YyxxP4/2Mti9dl99Ko3DH34i2VDWY2IekVviXrTKol\nnnysC6YUkhbj39I3mtlbSadyrzeWKZv93M1mzKLXzI7P4rxKXZarV0xeVwGbzWzK/iMpGVSL5Uvu\n8Rkzu5Qrf7RanJLWATuATnx310MFcYU6EmMs4V9i+OqpLZI2AUhqkrRg+mKFmoGPKalsw7/Rg7eI\nmqqUeYAnJCS14V1FLyi+4eafqySRkTQmNJsJAY/wGy2S1gBrc3834C6wW5PLvC+VL3me9xrvAgPf\nz31rOrchF0sf0PXrTUjtBbF04JMcitwBDqZWDpKWp9gK45S0DFhoZjeAk3jXWEUL3q0Y6ky0WMLf\nZlWO/QmzcUl7gAuSlgBfge3pXMuUq1ZP5fgqcDN1Lz0GhlL9I2mg/zlwGx+zqJTpAS6mMj+AAyme\nol31pjw2s0+SLuPdX8P4Uv0zXYMeoFfSIJ7ABoHPuXqHJJ0A+tKg/Tg+q+1Nrs6H+BYBT8zsfRpU\nH8DHsZ5mzusCuiU9wz//91N9p4FrkvancsN4Em7Ovlcz65e0GhhIrZxRYN80cX4DrmQmHBzLxLKe\nTJIL9SOWzQ9hnqSbbYOZjcl/f9IPtJnvdDjXulrxCQ07a4xlETCRWnhbgG4z+22GWllSi/Ccme2a\n8eTw34kWSwjzpxG4l7qrBBypJakAmNlLSaOSVlptv2VZAVxPye47cLiWOOagEzj7h18jzJNosYQQ\nQihVDN6HEEIoVSSWEEIIpYrEEkIIoVSRWEIIIZQqEksIIYRSRWIJIYRQqp8KV88o2OGHrQAAAABJ\nRU5ErkJggg==\n",
"text": [
"<matplotlib.figure.Figure at 0x10b2b7590>"
]
},
{
"output_type": "stream",
"stream": "stdout",
"text": [
"INFO: Inclination yields self-consistent solution for model.\n",
"incl = 88.8889180428 degrees\n",
"\n",
"Stellar radii in units of the star-star separation distance\n",
"(recalculated using timings from eclipse events):\n",
"radius_sep_g = 0.138957176768\n",
"radius_sep_s = 0.0749140127382\n"
]
}
],
"prompt_number": 3
},
{
"cell_type": "code",
"collapsed": false,
"input": [
"# Define functions for evaluation:\n",
"def is_calc_close_to_ref(calc, ref, tol=0.001, name=None, verbose=False):\n",
" \"\"\"Evaluate whether a calulated value is within a tolerance of the reference value.\n",
" \n",
" Parameters\n",
" ----------\n",
" calc : float\n",
" Calculated value to be tested.\n",
" ref : float\n",
" Reference value to be tested against.\n",
" tol : {0.001}, float, optional\n",
" Tolerance value to determine if `calc` is close to `ref`.\n",
" name : {None}, string, optional\n",
" Name of quantity being tested for printed output.\n",
" verbose : {False, True}, bool, optional\n",
" Print outcome of test if `True`.\n",
" \n",
" Returns\n",
" -------\n",
" is_close : bool\n",
" Boolean for whether or not calculated value was\n",
" within the tolerance of the reference value.\n",
" \n",
" \"\"\"\n",
" is_close = None\n",
" abs_diff = abs(ref - calc)\n",
" if (abs_diff < tol) or (abs_diff == tol):\n",
" is_close = True\n",
" if verbose:\n",
" print((\"INFO: Reference and calculated values match within tolerance.\\n\" +\n",
" \" name = {name}\\n\" +\n",
" \" ref = {ref}\\n\" +\n",
" \" calc = {calc}\\n\" +\n",
" \" abs(ref-calc) = {arc}\\n\" +\n",
" \" tol = {tol}\").format(name=name, ref=ref, calc=calc, arc=abs_diff, tol=tol)) \n",
" else:\n",
" is_close = False\n",
" if verbose:\n",
" print((\"ERROR: Reference and calculated values do not match within tolerance.\\n\" +\n",
" \" name = {name}\\n\" +\n",
" \" ref = {ref}\\n\" +\n",
" \" calc = {calc}\\n\" +\n",
" \" abs(ref-calc) = {arc}\\n\" +\n",
" \" tol = {tol}\").format(name=name, ref=ref, calc=calc, arc=abs_diff, tol=tol),\n",
" file=sys.stderr)\n",
" if is_close is None:\n",
" raise AssertionError(\"Program error. `is_close` was not set to `True` or `False`.\")\n",
" return is_close\n",
"\n",
"\n",
"def summarize_failures(tests_failed):\n",
" \"\"\"Print summary of tests that have failed.\n",
" \n",
" Parameters\n",
" ----------\n",
" tests_failed : array\n",
" Array of string names of tests that have failed.\n",
" \n",
" Returns\n",
" -------\n",
" None\n",
" \n",
" Notes\n",
" -----\n",
" The names from `tests_failed` will be formatted with a message and printed to stdout.\n",
" If `tests_failed` == [], then the message will state that all tests passed.\n",
" \n",
" \"\"\"\n",
" if tests_failed == []:\n",
" print(\"INFO: Tests complete. All tests passed.\")\n",
" else:\n",
" print((\"ERROR: Tests complete. Some tests failed:\\n\" +\n",
" \"tests_failed = {tf}\").format(tf=tests_failed), file=sys.stderr)\n",
" return None\n"
],
"language": "python",
"metadata": {},
"outputs": [],
"prompt_number": 4
},
{
"cell_type": "code",
"collapsed": false,
"input": [
"print(80*'=')\n",
"print(\"COMPARE TO PUBLISHED VALUES\")\n",
"radii_ratio_lt_ref = 0.539\n",
"radius_sep_s_ref = 0.075\n",
"radius_sep_g_ref = 0.139\n",
"incl_deg_ref = 88.9\n",
"tests_failed = []\n",
"for (name, ref, calc) in zip(['radii_ratio_lt', 'radius_sep_s', 'radius_sep_g', 'incl_deg'],\n",
" [radii_ratio_lt_ref, radius_sep_s_ref, radius_sep_g_ref, incl_deg_ref],\n",
" [radii_ratio_lt, radius_sep_s, radius_sep_g, np.rad2deg(incl_rad)]):\n",
" if not is_calc_close_to_ref(calc=calc, ref=ref, tol=0.002*ref, name=name, verbose=True):\n",
" tests_failed.append(name)\n",
"summarize_failures(tests_failed)"
],
"language": "python",
"metadata": {},
"outputs": [
{
"output_type": "stream",
"stream": "stdout",
"text": [
"================================================================================\n",
"COMPARE TO PUBLISHED VALUES\n",
"INFO: Reference and calculated values match within tolerance.\n",
" name = radii_ratio_lt\n",
" ref = 0.539\n",
" calc = 0.539115831462\n",
" abs(ref-calc) = 0.000115831461792\n",
" tol = 0.001078\n",
"INFO: Reference and calculated values match within tolerance.\n",
" name = radius_sep_s\n",
" ref = 0.075\n",
" calc = 0.0749140127382\n",
" abs(ref-calc) = 8.59872617724e-05\n",
" tol = 0.00015\n",
"INFO: Reference and calculated values match within tolerance.\n",
" name = radius_sep_g\n",
" ref = 0.139\n",
" calc = 0.138957176768\n",
" abs(ref-calc) = 4.28232323877e-05\n",
" tol = 0.000278\n",
"INFO: Reference and calculated values match within tolerance.\n",
" name = incl_deg\n",
" ref = 88.9\n",
" calc = 88.8889180428\n",
" abs(ref-calc) = 0.0110819572264\n",
" tol = 0.1778\n",
"INFO: Tests complete. All tests passed.\n"
]
}
],
"prompt_number": 5
},
{
"cell_type": "heading",
"level": 2,
"metadata": {},
"source": [
"Example from Carroll and Ostlie, 2007, An Introduction to Modern Astrophysics"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Examples 7.3.1, 7.3.2 of [Carroll and Ostlie, 2007, An Introduction to Modern Astrophysics](http://books.google.com/books/about/An_Introduction_to_Modern_Astrophysics.html?id=M8wPAQAAMAAJ).\n",
"\n",
"Assumptions: Spherical stars. No limb darkening. No orbital eccentricity. Orbital inclination = 90 deg."
]
},
{
"cell_type": "code",
"collapsed": false,
"input": [
"# Assumptions:\n",
"incl_rad = np.deg2rad(90.0)\n",
"print()\n",
"print(\"Assumed inclination: {ideg} deg\".format(ideg=np.rad2deg(incl_rad)))\n",
"print()\n",
"print(80*'=')\n",
"print(\"OBSERVED VALUES\")\n",
"print()\n",
"print(\"Smaller-sized star (_s) is brighter, primary.\")\n",
"print(\"Greater-sized star (_g) is dimmer, secondary.\")\n",
"print()\n",
"print(\"Stellar radial velocities:\")\n",
"velr_s = 33.0*sci_con.kilo # in m/s, smaller star\n",
"velr_g = 3.1*sci_con.kilo # in m/s, greater (larger) star\n",
"print(\"velr_s = {vs} m/s = {vskms} km/s\".format(vs=velr_s, vskms=velr_s/sci_con.kilo))\n",
"print(\"velr_g = {vg} m/s = {vgkms} km/s\".format(vg=velr_g, vgkms=velr_g/sci_con.kilo))\n",
"print()\n",
"print(\"Durations of eclipse events:\")\n",
"period = 8.6*sci_con.year # in seconds\n",
"ingress = 11.7*sci_con.hour # in seconds\n",
"totality = 164.0*sci_con.day # in seconds\n",
"egress = 11.7*sci_con.hour # in seconds\n",
"print(\"period = {per} sec = {peryr} yr\".format(per=period, peryr=period/sci_con.year))\n",
"print(\"ingress = {ing} sec = {inghr} hr\".format(ing=ingress, inghr=ingress/sci_con.hour))\n",
"print(\"totality = {tot} sec = {totdy} days\".format(tot=totality, totdy=totality/sci_con.day))\n",
"print(\"egress = {egr} sec = {egrhr} hr\".format(egr=egress, egrhr=egress/sci_con.hour))\n",
"print()\n",
"print(\"Times of eclipse events:\")\n",
"t0 = 0.0 # Mid-eclipse.\n",
"t1 = t0 - (totality/2.0) - ingress # Begin ingress.\n",
"t2 = t0 - (totality/2.0) # End ingress.\n",
"t3 = t0 + (totality/2.0) # Begin egress.\n",
"t4 = t0 + (totality/2.0) + ingress # End egress.\n",
"print(\"Mid-eclipse: t0 = {t0} sec = {t0d} days\".format(t0=t0, t0d=t0/sci_con.day))\n",
"print(\"Begin ingress: t1 = {t1} sec = {t1d} days\".format(t1=t1, t1d=t1/sci_con.day))\n",
"print(\"End ingress: t2 = {t2} sec = {t2d} days\".format(t2=t2, t2d=t2/sci_con.day))\n",
"print(\"Begin egress: t3 = {t3} sec = {t3d} days\".format(t3=t3, t3d=t3/sci_con.day))\n",
"print(\"End egress: t4 = {t4} sec = {t4d} days\".format(t4=t4, t4d=t4/sci_con.day))\n",
"print()\n",
"print(\"Light levels as magnitudes:\")\n",
"mag_ref = 6.3 # Brightness in magnitudes between minima.\n",
"mag_pri = 9.6 # Brightness in magnitudes during primary minimum, occultation.\n",
"mag_sec = 6.6 # Brightness in magnitudes during secondary minimum, transit.\n",
"print(\"Between minima: mag_ref = {mr} mag\".format(mr=mag_ref))\n",
"print(\"During occultation minima, primary eclipse: mag_pri = {mp} mag\".format(mp=mag_pri))\n",
"print(\"During transit maxima, secondary eclipse: mag_sec = {ms} mag\".format(ms=mag_sec))"
],
"language": "python",
"metadata": {},
"outputs": [
{
"output_type": "stream",
"stream": "stdout",
"text": [
"\n",
"Assumed inclination: 90.0 deg\n",
"\n",
"================================================================================\n",
"OBSERVED VALUES\n",
"\n",
"Smaller-sized star (_s) is brighter, primary.\n",
"Greater-sized star (_g) is dimmer, secondary.\n",
"\n",
"Stellar radial velocities:\n",
"velr_s = 33000.0 m/s = 33.0 km/s\n",
"velr_g = 3100.0 m/s = 3.1 km/s\n",
"\n",
"Durations of eclipse events:\n",
"period = 271209600.0 sec = 8.6 yr\n",
"ingress = 42120.0 sec = 11.7 hr\n",
"totality = 14169600.0 sec = 164.0 days\n",
"egress = 42120.0 sec = 11.7 hr\n",
"\n",
"Times of eclipse events:\n",
"Mid-eclipse: t0 = 0.0 sec = 0.0 days\n",
"Begin ingress: t1 = -7126920.0 sec = -82.4875 days\n",
"End ingress: t2 = -7084800.0 sec = -82.0 days\n",
"Begin egress: t3 = 7084800.0 sec = 82.0 days\n",
"End egress: t4 = 7126920.0 sec = 82.4875 days\n",
"\n",
"Light levels as magnitudes:\n",
"Between minima: mag_ref = 6.3 mag\n",
"During occultation minima, primary eclipse: mag_pri = 9.6 mag\n",
"During transit maxima, secondary eclipse: mag_sec = 6.6 mag\n"
]
}
],
"prompt_number": 6
},
{
"cell_type": "code",
"collapsed": false,
"input": [
"print(80*'=')\n",
"print(\"CALCULATED VALUES\")\n",
"print()\n",
"print(\"Smaller-sized star (_s) is brighter, primary.\")\n",
"print(\"Greater-sized star (_g) is dimmer, secondary.\")\n",
"print()\n",
"print(\"Semimajor axes of stellar orbits:\")\n",
"axis_s = bss.utils.calc_semimaj_axis_from_period_velr_incl(period=period, velr=velr_s, incl=incl_rad)\n",
"axis_g = bss.utils.calc_semimaj_axis_from_period_velr_incl(period=period, velr=velr_g, incl=incl_rad)\n",
"print(\"axis_s = {axs} m = {axsau} AU\".format(axs=axis_s, axsau=axis_s/ast_con.au.value))\n",
"print(\"axis_g = {axg} m = {axgau} AU\".format(axg=axis_g, axgau=axis_g/ast_con.au.value))\n",
"print()\n",
"print(\"Stellar masses:\")\n",
"mass_ratio = bss.utils.calc_mass_ratio_from_velrs(velr_1=velr_s, velr_2=velr_g)\n",
"mass_sum = bss.utils.calc_mass_sum_from_period_velrs_incl(period=period, velr_1=velr_s, velr_2=velr_g, incl=incl_rad)\n",
"(mass_s, mass_g) = bss.utils.calc_masses_from_ratio_sum(mass_ratio=mass_ratio, mass_sum=mass_sum)\n",
"print(\"mass_s = {ms} kg = {mssun} M_sun\".format(ms=mass_s, mssun=mass_s/ast_con.M_sun.value))\n",
"print(\"mass_g = {mg} kg = {mgsun} M_sun\".format(mg=mass_g, mgsun=mass_g/ast_con.M_sun.value))\n",
"print()\n",
"print(\"Stellar radii:\")\n",
"radius_s = bss.utils.calc_radius_from_velrs_times(velr_1=velr_s, velr_2=velr_g, time_1=t1, time_2=t2)\n",
"radius_g = bss.utils.calc_radius_from_velrs_times(velr_1=velr_s, velr_2=velr_g, time_1=t1, time_2=t3)\n",
"print(\"radius_s = {rs} m = {rssun} Rsun\".format(rs=radius_s, rssun=radius_s/ast_con.R_sun.value))\n",
"print(\"radius_g = {rg} m = {rgsun} Rsun\".format(rg=radius_g, rgsun=radius_g/ast_con.R_sun.value))\n",
"print()\n",
"print(\"Light levels, relative:\")\n",
"depth_pri = bss.utils.calc_flux_intg_ratio_from_mags(mag_1=mag_pri, mag_2=mag_ref)\n",
"depth_sec = bss.utils.calc_flux_intg_ratio_from_mags(mag_1=mag_sec, mag_2=mag_ref)\n",
"light_ref = 1.0 # Integrated flux between minima.\n",
"light_oc = depth_pri # Integrated flux during occultation minima.\n",
"light_tr = depth_sec # Integrated flux during transit minima.\n",
"print(\"Between minima: light_ref = {lr}\".format(lr=light_ref))\n",
"print(\"During occultation minima, primary eclipse: light_oc = depth_pri = {lo}\".format(lo=light_oc))\n",
"print(\"During transit minima, secondary eclipse: light_tr = depth_sec = {lt}\".format(lt=light_tr))\n",
"print()\n",
"print(\"Stellar radiative flux ratio:\")\n",
"flux_rad_ratio = bss.utils.calc_flux_rad_ratio_from_light(light_oc=light_oc, light_tr=light_tr, light_ref=light_ref)\n",
"print(\"flux_rad_ratio = {frr}\".format(frr=flux_rad_ratio))\n",
"print()\n",
"print(\"Stellar effective temperature ratio:\")\n",
"teff_ratio = bss.utils.calc_teff_ratio_from_flux_rad_ratio(flux_rad_ratio=flux_rad_ratio)\n",
"print(\"teff_ratio = {tr}\".format(tr=teff_ratio))"
],
"language": "python",
"metadata": {},
"outputs": [
{
"output_type": "stream",
"stream": "stdout",
"text": [
"================================================================================\n",
"CALCULATED VALUES\n",
"\n",
"Smaller-sized star (_s) is brighter, primary.\n",
"Greater-sized star (_g) is dimmer, secondary.\n",
"\n",
"Semimajor axes of stellar orbits:\n",
"axis_s = 1.42442349898e+12 m = 9.52168297795 AU\n",
"axis_g = 1.33809480207e+11 m = 0.894461128231 AU\n",
"\n",
"Stellar masses:\n",
"mass_s = 2.61291629258e+30 kg = 1.31361736091 M_sun\n",
"mass_g = 2.78149153727e+31 kg = 13.9836686806 M_sun\n",
"\n",
"Stellar radii:\n",
"radius_s = 760266000.0 m = 1.09310892182 Rsun\n",
"radius_g = 2.56521546e+11 m = 368.826161597 Rsun\n",
"\n",
"Light levels, relative:\n",
"Between minima: light_ref = 1.0\n",
"During occultation minima, primary eclipse: light_oc = depth_pri = 0.0478630092323\n",
"During transit minima, secondary eclipse: light_tr = depth_sec = 0.758577575029\n",
"\n",
"Stellar radiative flux ratio:\n",
"flux_rad_ratio = 3.94386308928\n",
"\n",
"Stellar effective temperature ratio:\n",
"teff_ratio = 1.40922538433\n"
]
}
],
"prompt_number": 7
},
{
"cell_type": "code",
"collapsed": false,
"input": [
"print(80*'=')\n",
"print(\"COMPARE TO PUBLISHED VALUES\")\n",
"axis_s_ref = 9.5*ast_con.au.value\n",
"axis_g_ref = 0.90*ast_con.au.value\n",
"mass_s_ref = 1.3*ast_con.M_sun.value\n",
"mass_g_ref = 13.9*ast_con.M_sun.value\n",
"radius_s_ref = 1.1*ast_con.R_sun.value\n",
"radius_g_ref = 369.0*ast_con.R_sun.value\n",
"depth_pri_ref = 0.048\n",
"depth_sec_ref = 0.76\n",
"flux_rad_ratio_ref = 3.97\n",
"teff_ratio_ref = 1.41\n",
"tests_failed = []\n",
"for (name, ref, calc) in zip(['axis_s', 'axis_g', 'mass_s', 'mass_g', 'radius_s', 'radius_g',\n",
" 'depth_pri', 'depth_sec','flux_rad_ratio', 'teff_ratio'],\n",
" [axis_s_ref, axis_g_ref, mass_s_ref, mass_g_ref, radius_s_ref, radius_g_ref,\n",
" depth_pri_ref, depth_sec_ref, flux_rad_ratio_ref, teff_ratio_ref],\n",
" [axis_s, axis_g, mass_s, mass_g, radius_s, radius_g,\n",
" depth_pri, depth_sec, flux_rad_ratio, teff_ratio]):\n",
" if not is_calc_close_to_ref(calc=calc, ref=ref, tol=0.02*ref, name=name, verbose=True):\n",
" tests_failed.append(name)\n",
"summarize_failures(tests_failed)"
],
"language": "python",
"metadata": {},
"outputs": [
{
"output_type": "stream",
"stream": "stdout",
"text": [
"================================================================================\n",
"COMPARE TO PUBLISHED VALUES\n",
"INFO: Reference and calculated values match within tolerance.\n",
" name = axis_s\n",
" ref = 1.42117977165e+12\n",
" calc = 1.42442349898e+12\n",
" abs(ref-calc) = 3243727331.2\n",
" tol = 28423595433.0\n",
"INFO: Reference and calculated values match within tolerance.\n",
" name = axis_g\n",
" ref = 1.3463808363e+11\n",
" calc = 1.33809480207e+11\n",
" abs(ref-calc) = 828603422.675\n",
" tol = 2692761672.6\n",
"INFO: Reference and calculated values match within tolerance.\n",
" name = mass_s\n",
" ref = 2.58583e+30\n",
" calc = 2.61291629258e+30\n",
" abs(ref-calc) = 2.7086292584e+28\n",
" tol = 5.17166e+28\n",
"INFO: Reference and calculated values match within tolerance.\n",
" name = mass_g\n",
" ref = 2.764849e+31\n",
" calc = 2.78149153727e+31\n",
" abs(ref-calc) = 1.66425372668e+29\n",
" tol = 5.529698e+29\n",
"INFO: Reference and calculated values match within tolerance.\n",
" name = radius_s\n",
" ref = 765058800.0\n",
" calc = 760266000.0\n",
" abs(ref-calc) = 4792800.0\n",
" tol = 15301176.0\n",
"INFO: Reference and calculated values match within tolerance.\n",
" name = radius_g\n",
" ref = 2.56642452e+11\n",
" calc = 2.56521546e+11\n",
" abs(ref-calc) = 120906000.0\n",
" tol = 5132849040.0\n",
"INFO: Reference and calculated values match within tolerance.\n",
" name = depth_pri\n",
" ref = 0.048\n",
" calc = 0.0478630092323\n",
" abs(ref-calc) = 0.000136990767736\n",
" tol = 0.00096\n",
"INFO: Reference and calculated values match within tolerance.\n",
" name = depth_sec\n",
" ref = 0.76\n",
" calc = 0.758577575029\n",
" abs(ref-calc) = 0.00142242497082\n",
" tol = 0.0152\n",
"INFO: Reference and calculated values match within tolerance.\n",
" name = flux_rad_ratio\n",
" ref = 3.97\n",
" calc = 3.94386308928\n",
" abs(ref-calc) = 0.0261369107164\n",
" tol = 0.0794\n",
"INFO: Reference and calculated values match within tolerance.\n",
" name = teff_ratio\n",
" ref = 1.41\n",
" calc = 1.40922538433\n",
" abs(ref-calc) = 0.000774615665908\n",
" tol = 0.0282\n",
"INFO: Tests complete. All tests passed.\n"
]
}
],
"prompt_number": 8
},
{
"cell_type": "heading",
"level": 3,
"metadata": {},
"source": [
"Including methods from Budding, 2007 with example from Carroll and Ostlie, 2007"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"**Note:** Calculating the ratio of stellar radii by using light levels is not valid for stars that differ greatly in radius, such as those in examples 7.3.1, 7.3.2 from Carroll and Ostlie, 2007."
]
},
{
"cell_type": "code",
"collapsed": false,
"input": [
"print(80*'=')\n",
"print(\"CALCULATED VALUES\")\n",
"print()\n",
"print(\"Smaller-sized star (_s) is brighter, primary.\")\n",
"print(\"Greater-sized star (_g) is dimmer, secondary.\")\n",
"print()\n",
"print(\"Integrated flux of stars relative to integrated flux of total system:\")\n",
"(flux_intg_rel_s, flux_intg_rel_g) = bss.utils.calc_fluxes_intg_rel_from_light(light_oc=light_oc, light_ref=light_ref)\n",
"print(\"flux_intg_rel_s = {firs}\".format(firs=flux_intg_rel_s))\n",
"print(\"flux_intg_rel_g = {firg}\".format(firg=flux_intg_rel_g))\n",
"print()\n",
"print(\"Ratio of stellar radii from relative light levels (invalid):\")\n",
"radii_ratio_lt = bss.utils.calc_radii_ratio_from_light(light_oc=light_oc, light_tr=light_tr, light_ref=light_ref)\n",
"print(\"radius_s / radius_g = radii_ratio_lt = {rrl}\".format(rrl=radii_ratio_lt))\n",
"print()\n",
"print(\"Numerical solution for orbital inclination (no solution):\")\n",
"phase_orb_ext = bss.utils.calc_phase_orb_from_time_period(time_event=t4, period=period, time_mideclipse=t0)\n",
"phase_orb_int = bss.utils.calc_phase_orb_from_time_period(time_event=t3, period=period, time_mideclipse=t0)\n",
"incl_rad = bss.utils.calc_incl_from_radii_ratios_phase_incl(radii_ratio_lt=radii_ratio_lt, phase_orb_ext=phase_orb_ext,\n",
" phase_orb_int=phase_orb_int, show_plots=True)\n",
"print(\"incl = {incl} degrees\".format(incl=np.rad2deg(incl_rad)))\n",
"print()\n",
"print(\"Stellar radii in units of the star-star separation distance\\n\"\n",
" \"(recalculated using timings from eclipse events, valid):\")\n",
"sep_proj_ext = bss.utils.calc_sep_proj_from_incl_phase(incl=incl_rad, phase_orb=phase_orb_ext)\n",
"sep_proj_int = bss.utils.calc_sep_proj_from_incl_phase(incl=incl_rad, phase_orb=phase_orb_int)\n",
"(radius_sep_s, radius_sep_g) = bss.utils.calc_radii_sep_from_seps(sep_proj_ext=sep_proj_ext, sep_proj_int=sep_proj_int)\n",
"print(\"radius_sep_s = {rs}\".format(rs=radius_sep_s))\n",
"print(\"radius_sep_g = {rg}\".format(rg=radius_sep_g))\n",
"print()\n",
"print(\"Stellar radii:\")\n",
"sep = bss.utils.calc_sep_from_semimaj_axes(axis_1=axis_s, axis_2=axis_g)\n",
"radius_s = bss.utils.calc_radius_from_radius_sep(radius_sep=radius_sep_s, sep=sep)\n",
"radius_g = bss.utils.calc_radius_from_radius_sep(radius_sep=radius_sep_g, sep=sep)\n",
"print(\"radius_s = {rs} m = {rssun} Rsun\".format(rs=radius_s, rssun=radius_s/ast_con.R_sun.value))\n",
"print(\"radius_g = {rg} m = {rgsun} Rsun\".format(rg=radius_g, rgsun=radius_g/ast_con.R_sun.value))"
],
"language": "python",
"metadata": {},
"outputs": [
{
"output_type": "stream",
"stream": "stdout",
"text": [
"================================================================================\n",
"CALCULATED VALUES\n",
"\n",
"Smaller-sized star (_s) is brighter, primary.\n",
"Greater-sized star (_g) is dimmer, secondary.\n",
"\n",
"Integrated flux of stars relative to integrated flux of total system:\n",
"flux_intg_rel_s = 0.952136990768\n",
"flux_intg_rel_g = 0.0478630092323\n",
"\n",
"Ratio of stellar radii from relative light levels (invalid):\n",
"radius_s / radius_g = radii_ratio_lt = 2.2458916679\n",
"\n",
"Numerical solution for orbital inclination (no solution):\n"
]
},
{
"metadata": {},
"output_type": "display_data",
"png": "iVBORw0KGgoAAAANSUhEUgAAAZgAAAEoCAYAAABl8ecgAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJztnXm8VVX5/98fEAQBwakQREEFp1LBFEzL64yoYJapfdUw\nLXNI+2WmZiVqpWWDmTl8+5b6zQG1zC/OU141B0wFZxBElMERRXFken5/rHW4m8MZ9rmcffe53Of9\neq3X2cManr32PvvZa61nPUtmhuM4juPUm055C+A4juOsmriCcRzHcTLBFYzjOI6TCa5gHMdxnExw\nBeM4juNkgisYx3EcJxNcwWSIpEsk/SSxf6ykNyS9L2ktSTtJmiZpgaTRecraWiTNlLR73nJkRbw3\nA1uZtlnSUfWVaOWQNFDSUknt4r8vaaykBxP7y+5H8f8rg7I3jOUpqzLKlLtU0sZtWWZWtIuHrBGJ\nL9aPorJ4V9JDko5JPoxmdqyZ/TzG7wL8FtjdzNY0s3eBs4ELzayXmU3I50pWGouhZmId7lZneepK\nvDczW5ucVtZNeyAPBZq8H8n/Vz0ofh7N7NVY3ip7D7PGFUzrMWA/M1sT2BA4DzgV+EuZ+H2BbsAL\niWMbAs+3pnBJnVuTrsEwoE2/Dp26UtOLV9JqWQmSouw0/xd/HuuNmXloRQBeBnYrOrY9sATYMu5f\nAZwDDAY+BJYCC4B7gekx7kfA+0AXoDdBQc0FZse0nWJeY4GHgN8BbxNaP12B3wCvAK8DlwDdYvym\nmMcPgDdinmMTsnYntKhmAvOBBxNpRwAPA+8Ck4FdqtTDacBzwDvAX4HVE+f3i3m8G+X/fDz+t8T1\nLwBOifX1g3i+f6yv4+L+JsC8avnGc/2AfwBvAjOA7yXOjQOuB66M9f4ssF2F61sKbJy4n38Cbolp\nHy2ci+f3BKbE+vwj0AwclTj/LcIHxTvAHcCGReV8D3gJeAv4NaAa0h4DvBjr46LEuU7xGXkr5n18\njF94rqo9c/8Gzo/lzgBGxnO/ABYDH8f7d2GJuhsYy/oW4RltjsdvAF6L9XQ/8f8Sz60DTADeAyZG\neR6scD/OKXPfxrLi/2Vj4F9x/y3gKqB3mefxhwn5C/XRL8o2D5gGHF2m7OHx+pL37yvAU3F7B+CR\neK/mxmelS5lrbGb5Z2hsUX1sDtwdZZoCHJQ4N4rwv3w/3tuT2/w92dYFriqBEgomHn8FOCZuXw6c\nHbc3Sj6spfIA/klQEt2B9eIf7DuJB2sR4QXRidAa+j1wE9AH6Bkf/l/G+E0x/jigM7APQckV/lB/\nin+29WN+IwgKq3/8AxZeJHvE/XXL1MNM4OmYbi3CC+mceG4oQbltT/gyPCJec5cy138kMCFuf4Og\nhMfH/W8B/6yWb7yWJ4CfAKsBgwgv1r1i2nGEl+LImPaXwCMV7nPxC+1t4AuxTq8Cro3n1o1/5APj\nue/H+v9WPD+G8FLaLMp4BvBQUTn3xns5AJhKfLGkTDsBWDOmfRPYO577LqHVXLg/9xFepIWXZrVn\nbiFwVKyr7wJzEuXeV7i+MnU3MMp2Rcx/9US+PeL9+j0wKZFmfAzdga0IL8YHytyPZf+vEmWPZcX/\nyybA7rHcdQnK7fcV/o8F+Qt19QBwEeF/sk2s513LlD8d2COxfwPwo7g9jKBkOhHeC88DJ5W5xuXq\nmISCiXU4C/hmzGtbguLcPJ5/DdgpbvcGhrb5e7KtC1xVQvHDmDj+CHB63L6clpftcg9rcR7AZ4FP\niK2IeOxQ4F+JB+uVxDkBH7D8F/SOwIy43UT4GkuW90biwf6IxFd/Is6pwP8WHbsDOKJCPXwnsb8P\nMD1uX1L8AiB8ZX2pVB3GF8A78douAb4DzIrnrgS+XyXfLxO+Hl8pOnc68Ne4PQ64K3FuS+CjCve5\n+IX230XX+kLcPgJ4uCjtLFoUzO0s/6LoRFD4AxLl7JU4fyxwTw1pv5g4fx0tL7N/Fd2fPWP8Timf\nuWmJc2vEtJ+J+/eR+LouUXcDY/yBFeL0iXF6ERTzQmBI4vwvKN+CWfb/KpHv2OLnoEScA4Any/2n\nE/J3IijuxUCPxPlfApeXyfsc4C9xuxfhvzqgTNzvAzeWucZKCuZgEso3HrsM+FncfoXwH1qzUj1k\nGXwMpv5sQHhJ1spGhC+r16LRwLvApYSvygKzEtvrEf7wTyTi3074Miswz8yWJvY/IrR01iV80b1U\nRo6DCnnGfHcijCGVIynXq4SuhEJeJxfltUHi/HKY2UuEF+e2wJcIXVFzJQ0hKI/7q+S7fjzXr+jc\n6cBnEkW9UVQn3Wqwqkqm/ZhQn8Rrml0UN1kvGwF/SMg0Lx7vXyZ+cT1WS/t6YrtwnyHUSXG+SZmq\nPXPL8jWzj+Jmz8R5ozrLypfUSdJ5kqZLeo/wUjfCM7keodVZTt5aSeaDpM9KGi9pdiz7b4QuuTT0\nA94xsw+LZOtfJv41wIGSuhJatU+Y2awoxxBJt0h6LcrxixrkSLIRMLzoWf8G4cMB4KuEbrKZ0SBj\nRCvKWClyG3RbFZG0PeFB/Hcrks8CPgXWKVIKSZJ/5rcJL7gtzey1Gst6m/DluimheyvJq8DfzOw7\nNeS3YdH2nERevzCzX5ZJV+rldD9wEKEbba6k+wlfbWsRxlwq5hv/RC+b2ZAayqwHcwldWQU5RPjq\nLfAq4Wv72gp5bEiLEUhxPVZLW47XWPH+FEjzzFUibV0m4/0XMJpgTfmKpD60tFrfIrQSNiR0ERbL\nu7Ly/ZLQPfg5M5sv6QDC+Ee5+EnmAmtL6mlmHyRkK/6oCBmZvSDpFUIr9xsEhVPgEkI37sFm9qGk\n7xOUQSk+JHSFFUh+6L0K3G9me5WR4XHggGjg8D3C2OPK1GfNeAtm5RCApDUl7QdcS3g5P5c8n4ao\nJO4CfiepV/zS20TSl8vEXwr8GbhA0npRjv6SSj5sJdL+NZa1vqTOknaMX1tXAftL2ise7yapSVK5\nLzUBx8ey1yaMD1wXz/0Z+K6kHRToIWlfSYUv4DcI3WJJ7gdOIPR3QxjkPIHQLVB4AVTK9zFggaQf\nSeoer+Fzkr6QkLe1VEp7G7CVpK9Ea6kTWf5lcCnwY0lbAkjqLemgojx+KKmPpAEx/XU1pC2WsyDr\n9cCJ8f6sRTDIAGp/5kpQ6v5VoydBqb0jqQfhpV+QZwlwIzAu3rstCeML5aj1XvYkvLDfj8/zKUXn\ny15PbH08DJwraXVJWxPGBa+qUN41hO6vLxHGYJJyLAA+krQ5oTu0HJMJLaHukjYljIcVuBUYIukw\nSV1i2F7S5nH7vyT1jvW6gKBc2xRXMCvHzZLeJ3xJnE6wyjoycd5Y/quo2hffEYQBxIK10A20vKSK\n84IwXjIdeDQ2te8Gkl/ulcr7IfAM8B9Cl8u5hPGa2YQv8R8TBjFfBU6m/LNiwNWEF9VLhMHonwOY\n2RPAtwkDo+/Ec0ck0p4L/CQ2738Qjz1A+AMWFMxDhAHfwn7FfKPy3I/QzTaD8FX834QB8IK8xfVS\nqZ6K71/JtGb2NqHldR6hhbgpiZasmd0E/AoYH+/VM8DeRXn9H+HLdhKhe/CvKdOWkimpjO8EngIe\nJ1jXJePX+swl9/8AfE3SO5IuoDTF6f+XMDYwh2DB90hRnBMI9/91wvX/lfL/oVLyVTp3FmGA/T3g\nZlasi1LPY/L8oYRxmbkERfgzM/tXmfIhfHB+GbjXzJLd5j8ktGreJzyb40tcV4HfE8al3iCMOV1F\nyzO3ANgLOIRQn6/Fa+ga0x4GvByfme8QWo9tilo+Ch3HyQtJS4FNzWxG3rI4Tr3wFozjOI6TCa5g\nHKcx8K4EZ5XDu8gcx3GcTPAWjOM4jpMJrmAcACSdLunPdchnmatxZehOPZpg3plF3lkhaZykv61E\n+nMlnVRPmdoKNegyAZK+JGlKnfP8jaTv1jPP9op3kTl1JQtrKIX1P2YAq7VyQmBDIOlMQt0c3oq0\n6xHMlzcxs0/rLlzGZHkPG80CT1JfwnysTcxsUd7y5ElDfU04ThXauyv1lZF/LHBre1QubUTDPBtm\n9jrBN167XESwnriCaedIOlXSDUXH/iDpD3F7rKSXFBZGmyHpG2XyWdZ9k+jOOELSK5LekvTjRNxO\nkn6s4E/qfUmPl5rpL+kKSefE7SYFH1A/UFjVc66ksYm4+0qaJOk9Sa/Gr/0ChUmW82N5I7TiSodf\nlPQfSfMlPSZpx8S5ZklnS/p3TH+npJK+n+JM+lskvRknEN6cvLZqeSXq7G1JP1GFRdXidTwcJ/ZN\nlrRLqXiRkbT4YiPKtSARlkg6IkVd9JM0QdI8hdVUj06cGyfpBkl/i9f2tKTBsfv0jXhdeybi95b0\nl3gvZ0s6p9AFFp+R38Rn5yVg3wrXhqQtYt2+K+lZSfsX1flRif1l915S4dl4KtbDQUX5rh7rYavE\nsfUUFgtcNz6XSV9p/ST9I97/GZK+F493k/SxgrcKJJ0haZGiV4p47b9PFN1c7Zo7BHl52WzUQFj7\n4gXCzOcbie7ti+IMIHg5fY4wG/nEEnFOJnhFXTvu70mYSf10/C3p5rsV8m5IcH/RM+53Jsw03oHg\nw+g9YHA891kSa28U5XMmwc0NtHiRvQxYHdia4Ltss3j+lHgdhXy3TlxnSXfqVF8+YBdgq7j9ecJM\n7jFxfyNW9EQ9lhavsmsT1tb4L8JH0yGEWelrxfPNhNn+mxKcfN4HnFumHtYmrN3RjTCj/HriMgHV\n8iJ4Zl4AfJHgRPJ8wizsgsfscYk6rnVZhDcps25NrMvZMc9qdVHW5TwtSxnsGe/RlYTlGE6P+0cT\nvXXH+JVc/VdcJqBI/i4EjxSnEfwj7kqY5V54vsp6FC5+5srUz1+Anyf2jwduSzyXBY/d1ZZ6uB84\nMG7fFZ+DkYl6HZMoo+DgMvd3Wp4hdwFyvfjwcF1edGxPWtZ/OA84r0S6vsC2cbsnwTHfFonzAwgu\n7l+m5cW7LdA3bm8FzK7jdTwIHJ6Qv+Auv0d82RwIdK+SxzhWVDD9EucnAl+P21OB/cvkU9KdOhWW\nDyiTzwXA74rkKadgDgceLUr/MPDNuH0f8OPEuWOB21PW7bYEL7pUywv4GXB14lx3gt+tUgqm1mUR\nlnNjnzg+JNbjF6vVBVVczkf57kyc25+gMAtjtb3ifViT6q7+yy4TUOIavgS8VnTsGuDMRJ2vjILZ\nnfifiPsPAYclnsuCgqm21MPZBPc4nQluWb5HcM3SjfBsr1V0vS+tzP96VQgdvYtsBQsHM7vbWgYh\nJxLcwBfHed3MJsftDwhfakkX9L8DflSUZrKFvlkIfp+6S+qy8pcAhD/joXH7GwTfYFhwLX4w4Wty\nbuz62ayGfMu5gN+A0q7+q1Fu+QAkDZd0X+yamE9YobEWV+rFbt1fYfl7kryWpJv95ZC0hqTLYtfW\ne4Sv1t6Skn385fJazmW/mX1Mi2v9YjaitmUR3iW84JOy9ib4LzvDzB5OyFCuLtanusv5N4uu7W2L\nb8y4D+F6N6Kyq/9KywQU068oblLmetAMrKHgHHUgoeX2zxLxNqLyUg/3ExTSMII/uHsILe/hBAX2\nbiKvXoQVOzs0HV3BVBsY/BbBS275DMIDO5SgjJA0htA6KXaDn+SrhOZzvSxM/g4UPB4fQMI1uJnd\nZcGdd1/CwGM5U+RazAlnEbqI0pA232sIq3NuYGZ9CC+rwvNZLY85hJdDko1ocXdfCycTWgU7mFlv\nwgsk6Z24EnNJfJBI6k55JVlYFmGtROhlZr8uE/9pwoqWhbw7EersXjP7n0S8SnWxzOV84lxZl/NV\nSLr6L8jf28w+H89XWiagmLnAgCIlnrx/lVzWV8WCN+HrCR9hhwI3FynZArMISz0k78maZrZfPP8I\n4R58hbD88wvxukYRlFiSLWhZXqLD0iEVjKRHJU0ivGxHKwwuT1LC1b2kM4CFZnZNhXx6El7uJ5nZ\nB5LWIHghPjMZrSjNVoSut2PqdT1m9hbhAb+C0Ec+NZb1GUljFNyiLyL8Ucu57K7FCud/gHMkbarA\n1oXBzxJ5ps23J/CumS2UtAOhJVZQLG8RukHKuYa/neC2/FBJq0k6mLBW+S1FsqSV42PgvXhNZ5aI\nUy6vfxCWOigsfTCuQtxal0W4jaDsCvyCsODc90vEK1kXFjxl1+pyviRW3dV/2WUCSvAooTX7IwU3\n800Ej9jj4/lKLush3bIB1xDGo4rXZklScakHCwuuPUEYwykYXDxM6CG4vyivXQjPZYemQyoYMxth\nZkMJg5YTzGxoDHdBsFIhfJWUdW8du7f+AVxlwZ06hId8IMGi5WXC1+wTkj4T02xAMBw43MxervNl\nXUPoa07+eToB/4/wJTiP0Nddbu0JI/3SAr8jvEDuIhgR/JnQD12crpY8jwPOVlj+4Ke0rIVS+GP/\nAnhIwbJreDJvM5tHeCGdTBgo/yGwny3vIr2SXEkuIIydvE14edxeIm7JvCysA/Q9wotxLmH84k3C\nl35x3FqXRfhfYJSkQj0fQuiaeVctlmSHxmuuVBeVXM6XqpdK+5Vc/VdbJqAlw9CS359grPAWwQjh\ncDN7MUYp67I+Mg64MnZrfa1MGY8Rli1enxVf/IV7soTKSz1AUCSrEZRRYT+5vASS1ie0YG6io5P1\nIA/BvHIKweLi1DJxLoznnwKGpk1LkaWWtQzKTYvp9qoiWxMrDvKPJFiHlbTmiXFE+MP/vkr+yUH+\nPvH6Dsi6zj00RiC8eBYBG9Upv18QWsu5X5uHivfpN8B385ajEULWFd2ZYH44kDAgOJmEtVWMM4oW\nk8HhRAuYamkpbam1ZYzXJaabTgmrlUQeuxAtRBLHphEGGCfFcHE83o8w0Q1g56jYJifijSyRf1K2\nnxC+oCYlQlkl5qF9BsKX+BqEMYNLcVNVDx04rEa27ECwrpgJIGk8oVvghUSc0QR7e8xsosJEt74E\nG/RKaQuWWv+XyGsMcK2FJvdMSdOjDI+WEs7M7qeo79TMBpeJO5c4ccrM/k2K7kUzG5TY/jlxpUdn\nlWY0oXUrwmqhh+QrjuPkR9ZjMP1Z3vywMBksTZxi08VlaStYai1nJlqmPMfJDDP7tgXroz5mtqeZ\nTctbJsfJi6xbMGlNVFNbMEXTzx8TJjKlSV+L+a3jOI5TJ7JWMHMIYyUFBrCizX1xnA1inC5l0iYt\ntQrxn4iWRaXyWmEuhCRXOo7jOK3AzNJPachygIegwF4iKISuVB/kH0HLIH/VtDFeqUH+rrT4EVKJ\nNNZonHnmmXmLsAIuU3oaUS6XKR0uU3riu7MxBvnNbLGkEwj28J2Bv5jZC5KOiecvM7PbJI2KA/If\nAkdWSluqmER5z0u6nmCXvxg4LlaK4ziO08Zk3UWGmd1O0cQmM7usaP+EtGlLxNm4aP+XBAd+juM4\nTo50yJn8jUhTU1PeIqyAy5SeRpTLZUqHy5QdHXLJZEnec+Y4jlMjkmoa5M+8i6xR6dQJpNoCrHya\nTp1C6Ny5JRTvlzq2snFWWw26dMk+dO0Kq68eynQcp2PTYRXMokVglj5AbfHLpVm6NIQlS0JIbpc7\nVo84ixeHay6Ejz5afr9eYeFC+PTToGBWXx26dVv+t9SxNOe6dYM11gihR4/y2926tSh2x3HyxbvI\nnLpjFhTaJ58EZVP4TW7XeuyTT4JS/Ogj+PDD8tuffgrdu1dXRL16wZprhpDcLrW/xhqutBwHau8i\ncwXjrFIsXdqicMopow8/hAUL4P33W34LIblf2P7kkxWVTu/esPbasNZa4Te5XXxs9dXzrhXHqQ+u\nYFLgCsaphcWLg7JJKp/58+Hdd+Gdd0JIbhfvd+3aonTWXRc+85kVw2c/27Lds6e3mJzGxBVMClzB\nOG2FWWgxFZTN22/Dm29WDkuWLK94+vcPYYMNlt/u3dsVkdO2uIJJgSsYp5H58MMWZfP66zB3LsyZ\nA7Nnh9/C9pIlyyudAQNg0KCWsOGGwbLPceqFK5gUuIJxVgXef79F4cyZA6++Ci+/3BJeew3WX395\npTNoEAweDJttBn365H0FTnvDFUwKXME4HYFFi1ZUOjNmwLRpMHVqGOvZfPOWsNlm4XfDDX0ek1Ma\nVzApcAXjdHTMQqtnypSWMHVq+J03D7baCrbZZvnQu3feUjt5U3cFI2kYcCjwZYLrfCOsWf8AcI2Z\nTWq1tDnhCsZxyrNgATzzDDz1FEyeHH6ffRbWWy8omqFDYfhw2GGHYBnndBzqqmAk3Qa8C0wAHgNe\nI6weuT5hrfv9gT5mtu/KCN3WuIJxnNpYsgReeikomyefhIkT4fHHwxjP8OEwYkT43XprNyxYlam3\ngvmsmb1RpcDPmNmbNciYO65gHGflWbIEnn8+KJtHHw2/M2cGZbPrrtDUBNtv7wpnVSLTMRhJ3Qgr\nmn3aGuEaBVcwjpMN8+fDAw9AczPcd19o9ey4Y1A4e+8N227rc3faM/VuwXQCDiCMwXyRsH6MgCXA\nI8DVwE3t7W3tCsZx2oZ33gkK51//gjvugA8+gFGjQthjj+B2x2k/1FvBPAA8SBiDmVxouUhaHRgK\njAZ2NrMvr5TUbYwrGMfJh2nT4LbbQnj44TBu89WvwoEHBq8FTmNTbwWzerXusDRxGg1XMI6TPx98\nAHfdBX//e1A4w4bBQQe5smlk6q1gKhohmtk7NcjWMLiCcZzG4uOPQxfaDTcEZbPjjnDkkTBmjHuj\nbiTqrWBmEua9lMTMBtUkXYPgCsZxGpePPoJ//hMuvzzMwznkEPjWt8L8GzcQyBefyZ8CVzCO0z6Y\nOROuvBKuuCJM9DzxRPj618MSCE7bk4mCkfQV4D4zmx/3+wBNZnZTqyXNEVcwjtO+WLIkdJ1deGHw\nKvDd78Ixx0DfvnlL1rGoVcF0ShlvXEG5AMTtcSkFGilpiqRpkk4tE+fCeP4pSUOrpZV0Tow7WdK9\nkgbE4wMlfSxpUgwXp7w+x3EamM6dYf/94e674Z57gqfoLbaAY4+FV17JWzqnHGkVTCmNVdXfqqTO\nwEXASGBL4FBJWxTFGQVsamaDge8Al6RI+2sz28bMtgVuAs5MZDndzIbGcFzK63Mcp52w1VZw6aXw\n4othyYFhw+Coo2D69Lwlc4pJq2CekPQ7SZtI2lTS74EnUqTbgfDCn2lmi4DxwJiiOKOBKwHMbCLQ\nR1LfSmnNbEEifU/g7ZTX4TjOKsJ668G554a5NQMGBMuzo44KXqKdxiCtgvkesAi4jvCi/wQ4PkW6\n/sCsxP7seCxNnH6V0kr6haRXgW8C5yXiDYrdY82Sdk4ho+M47Zi114Zx44Ki+cxngsPNn/wkLMjm\n5MtqaSKZ2QdAyfGTaklTxqvZ+NDMzgDOkHQa8HvgSGAuMMDM3o3LDNwkaauiFg8A48aNW7bd1NRE\nU1NTrSI4jtNA9OkTWjTHHgs/+xkMGQJnnw1HHw2d0n5KO8vR3NxMc3Nzq9OntSLbDPghYT2YglIy\nM9utSroRBAOBkXH/dGCpmf0qEedSoNnMxsf9KcAuwKBqaePxDYHbzOxzJcq/DzjZzJ4sOu5WZI6z\nijNpEhwXR2EvuSQ42nRWjqysyG4AngR+ApySCNV4HBgcrbu6AgcT/JolmQAcAcsU0vy4REDZtJIG\nJ9KPASbF4+tG4wAkbQwMBmakvEbHcVYhhg6Fhx4KkzT32gv+3/+DDz/MW6qORdoWzBNmtl2rCpD2\nAS4gWJ39xczOlXQMgJldFuMUrMU+BI4stDhKpY3H/w5sRvDq/BJwrJm9KelA4GzCeNFS4GdmdmsJ\nmbwF4zgdiLfegh/8IKxbc+WV8MUv5i1R+ySriZbjgLeAG4Flji3dF5njOO2JG2+E44+Hb34TzjrL\n/ZzVSlYKZiYlBuzdF5njOO2NN98MngBeeik41xwyJG+J2g/uiywFrmAcp2NjBpddBj/9KVx8cVgm\nwKlOvb0p725m90r6KqVbMDe2Tsx8cQXjOA7Ak08G5TJqFPz2t+5Esxr1VjBnmdmZkq6gtII5slVS\n5owrGMdxCsyfD0ccAQsWwD/+ESZuOqXxLrIUuIJxHCfJkiVw2mlw001wyy2w2WZ5S9SY1HUejKSx\nksrO9pfUVVK7bMU4juMU6NwZzj8/KJkvfxlWYvK6k6Caq5iewH/i7Pr/AK8T3Lr0Bb4AbA78OVMJ\nHcdx2oijjoKNNw6Lmv3P/8Do0XlL1L6p2kUmScBOwM7AhvHwK8C/gYfbY1+Td5E5jlOJxx+H/faD\n3/wGDjssb2kaBx+DSYErGMdxqvH887D33qHb7Pg0vuM7ALUqmIpdZJL+mNg1WrweG4CZnVizhI7j\nOO2ALbeEBx6A3XYLYzTf/W7eErU/qo3BFBYV+yJhVcnrCErmIOC5DOVyHMfJnUGDwhLNTU3QrRuM\nHZu3RO2LtK5iJgI7x5UlkdQF+LeZDc9YvkzwLjLHcWph6tTQkvntb+GQQ/KWJj/q2kWWoA+wJjAv\n7veKxxzHcVZ5NtsM7rwT9tgjLNW8++55S9Q+SLsezHnAk5KulHQlYW2Yc7MTy3Ecp7H43Ofg+uvh\n0EPhmWfylqZ9kNqKTNL6wHDCAP9EM3s9S8GyxLvIHMdpLePHw49+BA8/DBtskLc0bUtmZsqS1gKG\nAN1osSJ7oDVC5o0rGMdxVoZf/xquuSasmNmjR97StB1ZrQfzbeBEYANgMjACeMTMdmutoHniCsZx\nnJXBLFiULVoEV18NSv3Kbd/U1RdZgpOAHYBXzGxXYCjwXivkcxzHafdIcOmlMGUKXHBB3tI0Lmmt\nyD4xs48lIambmU2R5P5GHcfpsHTvHpZgHjECtt0Wdt01b4kaj7QtmNlxDOYm4G5JE4CZmUnlOI7T\nDhg4EK66Cr7xDXi93Zo9ZUfNvsgkNRHmxNxhZguzECprfAzGcZx68tOfwmOPwe23Q6e0n+3tkLoP\n8sf1YJ41s81XVrhGwRWM4zj1ZPHisI7M174GP/hB3tJkR90H+c1sMTBV0kYrJZnjOM4qymqrBWuy\nc8+FJ5/MW5rGIW1jbm3gOUn/knRzDBPSJJQ0UtIUSdMknVomzoXx/FOShlZLK+mcGHeypHslDUic\nOz3GnyLCvQZYAAAgAElEQVRpr5TX5ziOs1IMGhQsyo44Aj79NG9pGoO082CaShw2M7u/SrrOwFRg\nD2AOYVXMQ83shUScUcAJZjZK0nDgD2Y2olJaSb3MbEFM/z1gGzM7WtKWwDXA9kB/4B5giJktLZLL\nu8gcx6k7ZnDAAbD11nDOOXlLU38ycXZpZs1VCn3EzHYscWoHYLqZzYzxxgNjgBcScUYDV8ZyJkrq\nI6kvMKhc2oJyifQE3o7bY4Bro9fnmZKmRxkeTXOdjuM4K0Nhfsw228CBB8LQodXTrMrUy96hW5nj\n/YFZif3Z8ViaOP0qpZX0C0mvAmNpcbzZL8arVJ7jOE5mrL8+nH8+HHkkLGyXdrb1I+1Ey9aSth+q\nZkcLZnYGcIak04ALgCNrkWHcuHHLtpuammhqaqpVBMdxnJIccURwinnBBcExZnulubmZ5ubmVqfP\nWsHMAQYk9gewfAujVJwNYpwuKdJCGHO5rUJec0oJllQwjuM49USCiy6C4cODe/8BA6qnaUSKP77P\nOuusmtJnPSXocWCwpIGSugIHA8XWZxOAIwAkjQDmm9kbldJKGpxIPwaYlMjrEEldJQ0CBgOPZXNp\njuM45dlkEzj+eDj55LwlyY/ULZg48L49ocvpMTN7M3H6iFJpzGyxpBOAO4HOwF+iFdgx8fxlZnab\npFFxQP5DYldXubQx63OjL7QlwEvAsTHN85KuB54HFgPHubmY4zh5cdppsNVWcNddsFcHnDSR1kz5\n68D5QMEs+cvAKWZ2Q4ayZYabKTuO01bccktoxTz7LHTpkrc0K0dW68E8DexRaLVIWg+418y2brWk\nOeIKxnGctmSvvWD0aDjhhLwlWTmyWg9GwFuJ/Xm0wvLLcRynI3L++fDzn8N7HWwVrbQK5g7gTklj\nJR1JsNq6PTuxHMdxVh222Qb22Qd+9au8JWlb0naRCTgQ2JkwyP+gmf0zY9kyw7vIHMdpa2bPDopm\n8uT2a7acyRjMqoYrGMdx8uCMM+C11+Cvf81bktZRVwUj6SEz20nSB6w4I97MbM1WypkrrmAcx8mD\nd9+FwYPh0Udh003zlqZ2vAWTAlcwjuPkxVlnwcsvwxVX5C1J7WRiRSbpb2mOOY7jOJX5/vfh1lvh\nxRfzliR70lqRfS65E5dR3q7+4jiO46za9O4NJ520aq4XU0xFBSPpx5IWAJ+XtKAQgDdZ0aeY4ziO\nk4ITT4Q774SpU/OWJFvSmimfZ2antYE8bYKPwTiOkzdnnw2zZsGf/5y3JOnJbJBf0loE78TLFhcz\nswdqlrABcAXjOE7ezJsXLMqeey4sUtYeyMoX2beBEwlrrUwCRgCPmNlurRU0T1zBOI7TCJx4Iqyx\nBpx3Xt6SpCMrBfMswVX/I2a2raTNgXPN7CutFzU/XME4jtMIzJwJX/gCzJgBa7aDWYVZObv8xMw+\njgV0M7MpwGatEdBxHMcJDBwIe+8Nl12WtyTZkFbBzI5jMDcBd0uaAMzMTCrHcZwOwimnwAUXwMKF\neUtSf2qeyS+pCVgTuMPM2mWVeBeZ4ziNRFMTHHssHHxw3pJUpu5jMHFS5bNmtvnKCtcouIJxHKeR\n+Mc/4Pe/h3//O29JKlP3MRgzWwxMlbTRSknmOI7jlGTMGHj1VZg0KW9J6ktaK7IHgaHAY8CH8bCZ\n2egMZcsMb8E4jtNonHsuTJvW2K78szJTbipx2Mzs/hpkaxhcwTiO02i89RYMGQLTp8M66+QtTWly\ncdcv6REz23GlM2ojXME4jtOIjB0LW2wBp56atySlyWoeTDW6VY/iOI7jVOK444JvsqVL85akPtRL\nwZRF0khJUyRNk1RSL0u6MJ5/StLQamklnS/phRj/Rkm94/GBkj6WNCmGi7O+PsdxnHqx/fbBdcz9\n7XLwYUUyVTCSOgMXASOBLYFDJW1RFGcUsKmZDQa+A1ySIu1dwFZmtg3wInB6IsvpZjY0huOyuzrH\ncZz6IsHRR7cvD8uVyLoFswPhhT/TzBYB44ExRXFGA1cCmNlEoI+kvpXSmtndZlZoRE4ENsj4OhzH\ncdqEww6D224L3pbbO/VSMEeUOd4fmJXYnx2PpYnTL0VagG8BtyX2B8XusWZJO6eQ3XEcp2FYe23Y\nd1+4+uq8JVl5Vqt0UtJDZraTpA+AYrMrM7M148YzZbJIa6qV2iqhSL4zgIVmdk08NBcYYGbvShoG\n3CRpKzNbUJx23Lhxy7abmppoampqjQiO4zh15+ijgyv/730vdJvlRXNzM83Nza1OXxcz5bKZSyOA\ncWY2Mu6fDiw1s18l4lwKNJvZ+Lg/BdgFGFQpraSxwLeB3c3skzLl3wecbGZPFh13M2XHcRqWpUvD\nnJirr4bhw/OWpoW6milLWjP+rl0qpMj/cWBwtO7qChwMTCiKM4HYxRYV0nwze6NSWkkjgVOAMUnl\nImndaByApI0JK3DOSCGn4zhOw9CpE3zzm/C3v+UtycpRsQUj6VYz21fSTEp0d5nZoKoFSPsAFwCd\ngb+Y2bmSjonpL4txCtZiHwJHFlocpdLG49OArsA7sZhHzOw4SV8FzgIWAUuBn5nZrSVk8haM4zgN\nzYwZofUydy506ZK3NIFcZvK3N1zBOI7THth55zCrf//985YkUKuCqTbIP6zS+eKxDcdxHKd+HHYY\nXHVV4yiYWqnWRdZM6BrrDmwHPB1PbQ083p78jyXxFozjOO2BefNg442DK//evfOWps6D/GbWZGa7\nEsx/h5nZdma2HcF1/9yVE9VxHMepxDrrwG67wY035i1J60g70XLz5FwXM3sW2KJCfMdxHKcOHHZY\n+7UmS7sezHjgA+AqwqTIbwA9zezQbMXLBu8icxynvfDJJ9C/P0yeDAMG5CtLVu76jwSeB04CTozb\nR9YunuM4jlML3brBV74C11+ftyS142bKjuM4Dc6dd8KZZ8Kjj+YrR1ZLJg8Bfklwm989HjYz27hV\nUuaMKxjHcdoTixbB+uvDE0/ARhvlJ0dWXWSXA5cCi4FdCe71VwFfn47jOI1Ply5wwAHw97/nLUlt\npFUw3c3sHkKLZ6aZjQP2zU4sx3EcJ8lBB8ENN+QtRW2kVTCfRCeS0yWdIOlAoEeGcjmO4zgJdtsN\npk+HV17JW5L0pFUwJwFrECzIvgAcBnwzK6Ecx3Gc5WmP3WRVFUxsuRxsZgvMbJaZjTWzA80sZ3sG\nx3GcjkV76yarqmDMbAmws5TnumqO4zhOe+smS9tFNhn4P0mHS/pqDAdmKZjjOI6zPF26wJgx7cc3\nWVoF0w2YB+wG7BdDO3Ug7TiO03454AD4v//LW4p0+Ex+x3GcdsTHH0PfvmHFy3XWaduys5po6TiO\n4zQA3bvD7rvDLbfkLUl1XME4juO0M9pLN5l3kTmO47QzCitdvv56aNG0FZl3kUlqBw0zx3GcVZd1\n1oFhw+Cee/KWpDKt6SLrX3cpHMdxnJoYMwZuuilvKSrTGgUzqe5SOI7jODUxZgzcfDMsWZK3JOWp\nWcGY2bdqiS9ppKQpkqZJOrVMnAvj+ackDa2WVtL5kl6I8W+U1Dtx7vQYf4qkvWq9PsdxnPbAoEFh\njZhHHslbkvJkakUW/ZhdBIwkLFZ2qKQtiuKMAjY1s8HAd4BLUqS9C9jKzLYBXgROj2m2BA6O8UcC\nF0tySznHcVZJDjigsbvJsn757gBMj2vILALGA2OK4owmLGCGmU0E+kjqWymtmd1tZktj+onABnF7\nDHCtmS0ys5nA9JiP4zjOKsfo0Y09HyaVgpF0UJpjJegPzErsz2ZFI4FycfqlSAvwLeC2uN0vxquW\nxnEcp90zdCi89x689FLekpRmtZTxfgwUO4kudayYtJNNWuWpWdIZwEIzu6ZWGcaNG7dsu6mpiaam\nptaI4DiOkxudOsGoUXDrrXDiifXPv7m5mebm5lanr6hgJO0DjAL6S7qQFkXQC1iUIv85wIDE/gCW\nb2GUirNBjNOlUlpJY6Nsu1fJa04pwZIKxnEcp72y775w2WXZKJjij++zzjqrpvTVusjmAk8An8Tf\nQpgA7J0i/8eBwZIGSupKGICfUBRnAnAEgKQRwHwze6NSWkkjgVOAMWb2SVFeh0jqKmkQMBh4LIWc\njuM47ZI99wyWZB98kLckK1KxBWNmTwFPSbo6DrTXhJktlnQCcCfQGfiLmb0g6Zh4/jIzu03SKEnT\ngQ+BIyuljVn/EegK3B3XQXvEzI4zs+clXQ88DywGjnOfMI7jrMr06gU77AD33hvmxjQSFX2RSXqm\nQlozs63rL1L2uC8yx3FWJX73O5gyBf77v7Mtp1ZfZNUUzMBKiaMpcLvDFYzjOKsSU6cGF/6zZkGW\ni9vXqmCqdZHNTFnoI2a2Y9pCHcdxnPoxZAh06wZPPQXbbpu3NC3Ua6Jltzrl4ziO49SIFKzJbr01\nb0mWx92oOI7jrAK4gnEcx3EyYZdd4Lnn4O2385akBVcwjuM4qwCrrw677gp33JG3JC2kVjCS+kra\nX9J+kj5TdPqIOsvlOI7j1EijdZNVNFNeFkn6OnA+cH889GXgFDOr5ousIXEzZcdxVkXmzIGtt4Y3\n3oDV0nqarIG6mikn+AmwvZm9GQtZD7iX6s4uHcdxnDaif38YMAAmToSddspbmvRdZALeSuzPo5Ue\nkB3HcZzsGDUKbrutery2IK2CuQO4U9JYSUcS1l+5PTuxHMdxnNbQSOMwacdgBBwI7ExYX+VBM/tn\nxrJlho/BOI6zqrJ4MXz2s/D006HLrJ7U1RfZqoorGMdxVmUOPTT4Jjv66PrmW6uCqdhFJumh+PuB\npAVF4f2VFdZxHMepP40yDuMtGMdxnFWMt96CwYPhzTeha9f65VvXFkwi07+lOeY4juPkz3rrweab\nw4MP5itHWiuyzyV3JK0GbFd/cRzHcZx60AjdZNXGYH4saQHw+eT4C/AmMKFNJHQcx3FqphEUTFoz\n5fPM7LQ2kKdN8DEYx3FWdZYuhfXXh0cegY03rk+emZkpS1oLGExicTEze6BmCRsAVzCO43QExo6F\n7beH44+vT35ZDfJ/G3gAuAs4C7gTGNcaAR3HcZy2Ie9Z/WkH+U8CdgBmmtmuwFDgvcykchzHcVaa\nPfcMlmQffZRP+WkVzCdm9jGApG5mNgXYLDuxHMdxnJWlTx8YNgyam/MpP62CmR3HYG4C7pY0AZiZ\nJqGkkZKmSJom6dQycS6M55+SNLRaWkkHSXpO0hJJwxLHB0r6WNKkGC5OeX2O4zirJHlak9U8k19S\nE7AmcIeZLawStzMwFdgDmAP8BzjUzF5IxBkFnGBmoyQNB/5gZiMqpZW0ObAUuAw42cyejHkNBG42\ns89XkcsH+R3H6RA88wyMHg0zZoBWcpGVug/yS1pN0pTCvpk1m9mEasolsgMw3cxmmtkiYDwwpijO\naODKmPdEoI+kvpXSmtkUM3sxRfmO4zgdms99DpYsgSlTqsetN1UVjJktBqZK2qgV+fcHZiX2Z8dj\naeL0S5G2FINi91izpJ1rF9lxHGfVQcqvmyztkslrA89Jegz4MB4zMxtdJV3afqh6rY45FxhgZu/G\nsZmbJG1lZguKI44bN27ZdlNTE01NTXUSwXEcp7EYNQr+8Ac4+eTa0jU3N9O8EhYCaWfyN5U4bGZ2\nf5V0I4BxZjYy7p8OLDWzXyXiXAo0m9n4uD8F2AUYlCLtfSTGYEqUX/K8j8E4jtOR+OCDMKt/zhxY\nc83W55PJRMs47lIclikXSY+USfo4MDhad3UFDmZFH2YTgCNiPiOA+Wb2Rsq0kGj9SFo3GgcgaWOC\n54EZaa7RcRxnVaVnT9hxR7j33rYtN62ZcjW6lToYx29OIMz8fx64LlqBHSPpmBjnNmCGpOkEq7Dj\nKqUFkPQVSbOAEcCtkm6PRe4CPCVpEnADcIyZza/TNTqO47Rb8pjVX5cFxyRNMrOh1WM2Bt5F5jhO\nR2PaNNhll9BN1lpz5Uy6yBzHcZz2zeDB0KMHPPVU25XpCsZxHKeD0NbmyvVSMEfUKR/HcRwnI/bd\nF26+ue3KqzgGI+khM9tJ0gesOKfFzGwlDN7yw8dgHMfpiCxcCH37Bvcx/dNMWy+irmMwZrZT/O1p\nZr2KQrtULo7jOB2Vrl1hv/3gppvapryKCkbSmvF37VKhbUR0HMdx6sWBB8I//tE2ZVXrIrvVzPaV\nNJMSbl/MbFCGsmWGd5E5jtNR+fjj0E02fTqst15taWvtIqvLPJj2hisYx3E6Ml//Ouy9Nxx1VG3p\n6qpgkot5laKcD7BGxxWM4zgdmeuugyuvrN1kud4KppnQNdYd2A54Op7aGnjczHasTbzGwBWM4zgd\nmQULghXZq6+GZZXTUm8rsiYz25XgBn+YmW1nZtsBQ+Mxx3Ecp53Rqxc0NWU/JybtRMvNzeyZwo6Z\nPQtskY1IjuM4TtYccghcc022ZaRdD2Y88AFwFcE9/jeAnmZ2aLbiZYN3kTmO09H56KPQTfbCC8Gq\nLA1ZObs8kuAy/yTgxLh9ZNpCHMdxnMZijTVgzBi49trsynAzZcdxnA7KvffCKafAkyntgTNpwUga\nIunvkp6X9HIMvlKk4zhOO6apCd58E557Lpv803aRXQ5cCiwGdgWuBK7ORiTHcRynLejcGQ4/HP76\n12zyTzvI/6SZDZP0jJl9PnksG7GyxbvIHMdxAi+/DNtvH+bErLFG5bhZDfJ/IqkzMF3SCZIOBHqk\nLcRxHMdpTAYNghEjwuz+epNWwZwErEGwIPsCcBjwzfqL4ziO47Q1xx0Hf/oT1Ltjp2oXWWy5/MrM\nfljfovPDu8gcx3FaWLIENtsMLr8cvvSl8vHq3kVmZkuAnSWlztRxHMdpP3TuDKefDuecU9980w7y\nXwr0A24APoqHzcxurK84bYO3YBzHcZZn4UIYMgTGjw9jMqXIapC/GzAP2A3YL4b90ySUNFLSFEnT\nJJ1aJs6F8fxTkoZWSyvpIEnPSVpSvKSApNNj/CmS9kp5fY7jOB2arl3htNPgzDPrNxaT6Uz+OH4z\nFdgDmAP8BzjUzF5IxBkFnGBmoyQNB/5gZiMqpZW0ObAUuAw4ubAujaQtgWuA7YH+wD3AEDNbWiSX\nt2Acx3GKWLgQttkGzjsvuJEppq4tGEnjJH22wvn1JZ1VIYsdgOlmNtPMFgHjgWKxRxMmbmJmE4E+\nkvpWSmtmU8zsxRLljQGuNbNFZjYTmB7zcRzHcarQtStcdBGcdBJ88MHK57dalfOPA+MldQWeBF4j\neFPuCwwDPgV+UyF9f2BWYn82MDxFnP6EMZ9qaYvpBzxaIi/HcRwnBbvvHsJ3vgNXXw0rY95VUcGY\n2S3ALZIGADsBG8ZT/yaYLs+ukn/afqgsLdRKyjBu3Lhl201NTTQ1NWUoguM4Tvvhootgxx3hmGOa\n6devudX5VGvBAGBmswhdVIVxlR5m9n6KpHOAAYn9AYRWRaU4G8Q4XVKkrVbeBvHYCiQVjOM4jtNC\n9+5wyy2w++5N7LNPE7/4BfToAWedVWlEZEXSelO+VtKaknoAzwAvSPpRiqSPA4MlDYzdbAcDE4ri\nTACOiOWMAOab2Rsp08LyrZ8JwCGSukoaBAwGHktzjY7jOE4LG2wADz0Eb70FAwbA3nvXnkdaM+Ut\nY4vlAOB2YCBweLVEZrYYOAG4k7BI2XXRCuwYScfEOLcBMyRNJ1iFHVcpLYCkr0iaBYwAbpV0e0zz\nPHB9jH87cJybizmO47SOddcN4zDPPhsG/msl7UTL54BtCSbAfzKzZklPm9nWtReZP26m7DiOUztZ\nTbS8DJgJ9AQekDQQeK9W4RzHcZyOQ6smWka/ZKvF+SntDm/BOI7j1E5WSyavK+mPkiZJehK4AFiz\ntUI6juM4qz5pu8jGA28CBwJfA94CMliexnEcx1lVSDvI/6yZfa7o2LLlk9sb3kXmOI5TO1kN8t8l\n6VBJnWI4GLirdSI6juM4HYGKLRhJH9DiaqUHwYMxBMX0oZn1yla8bPAWjOM4Tu3U2oKp5ous58qL\n5DiO43REUvkiA5C0FsH1SrfCMTN7IAuhHMdxnPZPKgUj6dvAiQRHkpMILloeIaxw6TiO4zgrkHaQ\n/yTCwl0zzWxXYCg+k99xHMepQFoF84mZfQwgqZuZTQE2y04sx3Ecp72TdgxmVhyDuQm4W9K7BN9k\njuM4jlOSmn2RSWoiuIm5w8wWZiFU1riZsuM4Tu3UaqbcKmeX7R1XMI7jOLWT1Ux+x3Ecx6kJVzCO\n4zhOJriCcRzHcTLBFYzjOI6TCa5gHMdxnExwBeM4juNkgisYx3EcJxMyVzCSRkqaImmapFPLxLkw\nnn9K0tBqaSWtLeluSS9KuktSn3h8oKSPJU2K4eKsr89xHMcpTaYKRlJn4CJgJLAlcKikLYrijAI2\nNbPBwHeAS1KkPQ2428yGAPfG/QLTzWxoDMdld3X1pbm5OW8RVsBlSk8jyuUypcNlyo6sWzA7EF74\nM81sETAeGFMUZzRwJYCZTQT6SOpbJe2yNPH3gGwvI3sa8YFymdLTiHK5TOlwmbIjawXTH5iV2J8d\nj6WJ069C2s+a2Rtx+w3gs4l4g2L3WLOknVdSfsdxHKeVpF7RspWkdfiVxreNSuVnZiapcHwuMMDM\n3pU0DLhJ0lZmtiClHI7jOE69MLPMAmHlyzsS+6cDpxbFuRQ4JLE/hdAiKZs2xukbt9cHppQp/z5g\nWInj5sGDBw8eag+16ICsWzCPA4MlDSS0Lg4GDi2KMwE4ARgvaQQw38zekDSvQtoJwDeBX8XfmwAk\nrQu8a2ZLJG0MDAZmFAtVizdQx3Ecp3VkqmDMbLGkE4A7gc7AX8zsBUnHxPOXmdltkkZJmg58CBxZ\nKW3M+jzgeklHERY++3o8/mXgbEmLgKXAMWY2P8trdBzHcUrTIdeDcRzHcbKnQ83kTzPps43k+Kuk\nNyQ9kzhWcvJoG8o0QNJ9kp6T9KykE/OWS1I3SRMlTZb0vKRz85YpIVvnaK14cyPIJGmmpKejTI81\niEx9JP1d0gvx/g1vAJk2S0zEniTpPUknNoBcp8f/3jOSrpG0egPIdFKU51lJJ8VjNcnUYRRMmkmf\nbcjlUY4klSaPtgWLgP9nZlsRDCyOj/WTm1xm9gmwq5ltC2wN7BpNz/OuK4CTgOcJA580gEwGNMUJ\nxjs0iEx/AG4zsy0I929K3jKZ2dTCRGxgO+Aj4J95yhXHmb9NMEj6PGFI4JCcZfoccDSwPbANsJ+k\nTWqWKUsrskYKwI4sb5V2GnBajvIMBJ5J7E8hzO8B6EsZy7g2lO8mYI9GkQtYA/gPsFXeMgEbAPcA\nuwI3N8L9A14G1ik6lptMQG9gRonjDfE8xfL3Ah7MWy5gbWAqsBZhXPxmYM+cZfoa8D+J/Z8AP6pV\npg7TgiHdpM88qTR5tE2JX1RDgYnkLJekTpImx7LvM7Pn8pYJ+D1wCsGQpEDeMhlwj6THJX27AWQa\nBLwl6XJJT0r6s6QeOctUzCHAtXE7N7nM7B3gt8CrBIvZ+WZ2d54yAc8CX4pdYmsAowgfVjXJ1JEU\nTLuxZrDweZCLvJJ6Av8ATrKiCap5yGVmSy10kW0AfFnSrnnKJGk/4E0zm0SZCcI53b+dLHT77EPo\n3vxSzjKtBgwDLjazYQQL0eW6U3J+zrsC+wM3FJ/L4ZnaBPg+oVejH9BT0mF5ymRmUwjTQO4Cbgcm\nA0tqlakjKZg5wIDE/gBCK6ZReEPBBxuS1gfebGsBJHUhKJe/mdlNjSIXgJm9B9xK6DfPU6YvAqMl\nvUz4+t1N0t9ylgkzey3+vkUYU9ghZ5lmA7PN7D9x/+8EhfN6IzxPBEX8RKwvyLeuvgA8bGbzzGwx\ncCOhSz/XujKzv5rZF8xsF+Bd4EVqrKeOpGCWTfqMXy8HEyZsNgqFyaOQmDzaVkgS8BfgeTO7oBHk\nkrSuWpZi6E7ol56Up0xm9mMzG2BmgwhdLP8ys8PzlEnSGpJ6xe0ehLGFZ/KUycxeB2ZJGhIP7QE8\nRxhfyO05T3AoLd1jkO//bwowQlL3+D/cg2BAkmtdSfpM/N0QOBC4hlrrqa0GjRohEL5apgLTgdNz\nlONaQl/rQsK40JGEgb57CF8JdwF92limnQljCpMJL/FJBEu33OQCPg88GWV6GjglHs+1rhLy7QJM\nyFsmwnjH5BieLTzbedcTwfroP8BThK/y3nnLFOXqAbwN9Eocy7uufkRQwM8QPMR3aQCZHogyTSZY\nc9ZcTz7R0nEcx8mEjtRF5jiO47QhrmAcx3GcTHAF4ziO42SCKxjHcRwnE1zBOI7jOJngCsZxHMfJ\nBFcwTsMg6aFWpmtKuM3fX61cikFSb0nHJvb7SVrBlUjeRNf8a9eY5rrokqT4+FhJf6yfdCuPpN8V\nu7px2ieuYJyGwcx2qkMeN5vZr1qZfC3guERec83soJWVKQNqmrwmaVOgh5m9lJE8KFKn7C4hOBN1\n2jmuYJyGQdIH8bdJUrOkGxQWq7oqEWd7SQ8pLEI2MTrnTOax7Itc0hWS/hDjvyTpq/F4T0n3SHpC\nYZGu0TH5ecAmCgtR/UrSRpKejWm6Rc/AT0fvwE2J8m6UdHtchKmkcpP0U0mPKSzgdFnieLOk8+K1\nTI3r3RTcv1yvsAjVjZIelTSsRL6HxbSTJF0qqdR/+hASbpEkHRnLmkjwrVY4vp7CAmGPxfDFxPG7\nFRae+nOhBRXdLk2VdCVhBvoASafEtE9JGldJToVF266IdfK0pO8DmNk0YKByWEjOqTNt7abBg4dy\nAVgQf5uA+QTPsgIeJrwIuwIvAdvFeD0JizM10bIuy1jgj3H7CuC6uL0FMC1udya6CQHWTRzfiOXX\n6BlY2AdOJq6PAWwGvAKsHst7CegV92cC/Utc21qJ7f8F9ovb9wHnx+19CIs5AfwQuCRub0VYEG5Y\n3H+Z4LJjC4Li6ByPXwwcXqLs2xNp14+yr0NwR/Jv4MJ47hqCV2aADQl+6SAs1Hdq3N6b4FJo7Vg/\nS1SGT6QAAANoSURBVIAd4rm9gMvidieCL60vlZDzT8DhBOeXdyXk7J3YvhLYJ+9n0sPKhdVwnMbk\nMTObC6CwHswgYAHwmpk9AWBmhRZPuTyM6IzPzF6QVFi7ohNwbuznXwr0i479KnXx7ARcGPOaKukV\nYEgs416LSxtIep7w4p1TlH43SacQFk5bm+Az7JZ47sb4+2RMWyjvgljec5KeLspPwO4E79KPxzro\nDrxeQvaNgNfi9nDCujrzorzXxeuA4GRxi0R99lJwnrkTcECU5U5J7ybyfsXMHovbewF7SZoU93sA\nmxJ8khXL+QZBAW0s6UKCp+y7EvnOTdSF005xBeM0Kp8mtpcQntXWOM5bmNguvDn/i9ByGWZmSxRc\n73dLkVc5BVQsa+flEkndCF/t25nZHElnFpX3aSJt8j+ZZkzjSjP7cYp4hbysKF/RUq8ChpvZwuUS\nBqVQTpYPi/bPNbP/Lkp/Qjk5JW1NcKr6XeDrwFEl5HLaKT4G47QXjOAJe31JXwCQ1EtS58rJSrIm\nYdGwJQoLmG0Ujy8gdHWV4kGCYkLBBf2GBDfrpV68xccKymReHDNKYzjwEOGFi6QtCZ6lkxhhTfSv\nSVovxltbwbV6Ma8QusYAHgN2iXG7FMlyF3DisouQtikhy14EY4hS3Al8K7Z6kNQ/ylZSTknrAKuZ\n2Y3ATwldZgXWJ3Q3Ou0Yb8E4jYSV2Q4HzBZJOhj4o8L6MB8R1ohJrqxXvMpeqe2rgZtjt9PjwAsx\n/3nRIOAZ4DbCmEYhzcXAJTHNYuCbUZ5Sq/ott29m8yX9mdAt9jphKepqdXAxcKWk5wiK7DngvaJ8\nX5D0E+CuOLi/iGAF92pRnv8mLGr1hJm9FgffHyGMc01KxDsR+JOkpwjvhvtjfmcB10o6PKZ7naCM\n10xeq5ndLWkL4JHY6lkAHFZBzk+AyxOGCckVL4eSUHZO+8Td9TtOAxJful3M7FOF+St3A0MsrHhY\na14bEwwf9m2lLF2BJbHFtyPwJwvLIGdCbCH+xsxGV43sNDTegnGcxqQH8K/YjSXg2NYoFwAzmyFp\ngaRNrHVzYTYEro9KbyHw7dbIUQPfBX6dcRlOG+AtGMdxHCcTfJDfcRzHyQRXMI7jOE4muIJxHMdx\nMsEVjOM4jpMJrmAcx3GcTHAF4ziO42TC/wfg8LffnuyBZwAAAABJRU5ErkJggg==\n",
"text": [
"<matplotlib.figure.Figure at 0x109545a90>"
]
},
{
"metadata": {},
"output_type": "display_data",
"png": "iVBORw0KGgoAAAANSUhEUgAAAaoAAAEoCAYAAAAJ97UlAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJzsnXm4ndP1xz/fjFJBxBQZCE0MUUMMkRJ1UZqmRNAitATF\nj6qphhhaQWusmYa2hhQRQaupMaFuzEPMFSEhIaaEmENIZP3+2Pskb07OdHPvOe85967P87zPeYc9\nrHc473r33muvJTPDcRzHcaqVVmkL4DiO4ziFcEXlOI7jVDWuqBzHcZyqxhWV4ziOU9W4onIcx3Gq\nGldUjuM4TlXjiqrKkDRS0umJ7SMkzZL0uaSVJW0raaqkLyQNTlPWZUXSDEk7pS1HuYj3pucy5q2X\ndEjTStQ4JPWUtFBSTbwvJA2T9Ehie9H9yP5/laHutWJ9KlcdeepdKGndStZZSWriwWsuxBf0V1Hp\nfCLpMUmHJx9qMzvCzP4Y07cFLgJ2MrMVzewT4CzgcjNbwczGpXMmjcbi0mDiNdyxieVpUuK9mbGs\n2VnGa1MLpKGIk/cj+f9qCrKfRzN7O9bXbO9hGriiqiwG7GpmKwJrAecBJwPX5knfBVgOeDWxby1g\n8rJULqn1suSrMgyo6Neq06Q06AUuqU25BCmh7lL+L/48VgIz86VCCzAd2DFr31bAd0CfuH0DcDbQ\nG5gLLAS+AB4EpsW0XwGfA22BlQiK7j3gnZi3VSxrGPAYcDHwEaE11g74M/AW8AEwElgupq+LZRwP\nzIplDkvI2oHQwpsBfAo8ksjbH3gc+AR4Adi+yHUYDrwCfAxcB7RPHN81lvFJlH/juP/GxPl/AZwY\nr9fx8Xi3eL2OjNvfB+YUKzce6wrcAcwG3gR+mzg2AhgLjIrX/X/AFgXObyGwbuJ+XgXcFfM+mTkW\nj+8MTInX8wqgHjgkcfxgwofJx8B9wFpZ9fwWeAP4ELgAUAPyHg68Hq/HlYljreIz8mEs+zcxfea5\nKvbMPQpcGOt9ExgYj/0JWAB8He/f5TmuXc9Y18GEZ7Q+7r8NeD9ep4nE/0s8tgowDvgMeCrK80iB\n+3F2nvs2jKX/L+sC/43bHwI3ASvleR5PSMifuR5do2xzgKnAr/PUvXU8v+T92wN4Ma73A56I9+q9\n+Ky0zXOO9Sz5DA3Luh4bABOiTFOAXySODSL8Lz+P9/Z3ab83zcwVVUUvdg5FFfe/BRwe168Hzorr\naycf+lxlAP8iKJsOwGrxj3pYPDYMmE940bQitM4uAe4EOgEd45/onJi+LqYfAbQGfkpQlpk/5lXx\nT7tmLK8/QfF1i3/kzAvpx3F71TzXYQbwUsy3MuHFdnY81pegJLcifKkeEM+5bZ7zPwgYF9f3Iyjz\nMXH7YOBfxcqN5/IscDrQBliH8ILeJeYdQXi5Dox5zwGeKHCfs1+MHwFbxmt6E3BLPLZqfCHsGY8d\nG6//wfH47oSX2/pRxtOAx7LqeTDeyx7Aa8QXVIl5xwErxryzgZ/EY/9HaMVn7s9DhBdy5uVb7Jn7\nFjgkXqv/A95N1PtQ5vzyXLueUbYbYvntE+UuH+/XJcDziTxj4tIB2Ijwgn04z/1Y9P/KUfcwlv6/\nfB/YKda7KkFJXlLg/5iRP3OtHgauJPxPNo3XeYc89U8DfpzYvg04Ka5vTlBWrQjvhcnAMXnOcYlr\nTEJRxWs4EzgwlrUZQQFvEI+/D2wb11cC+qb93jRzRVXZi51fUT0BnBLXr2fxS3uJhz67DGANYB6x\nVRP3DQX+G9eHAW8ljgn4kiW/6H8IvBnX6whfh8n6ZiX+IF+RaIUk0pwM/CNr333AAQWuw2GJ7Z8C\n0+L6yOwXCeGrb7tc1zC+SD6O5zYSOAyYGY+NAo4tUu6PCF+zb2UdOwW4Lq6PAMYnjvUBvipwn7Nf\njH/NOtdX4/oBwONZeWeyWFHdy5IvnFaED4ceiXp2SRw/AnigAXm3SRy/lcUvxf9m3Z+dY/pWJT5z\nUxPHvhfzrh63HyLxtZ/j2vWM6XsWSNMpplmBoOC/BdZLHP8T+VtUi/5fOcodlv0c5EgzBHgu3386\nIX8rwgfAAmD5xPFzgOvzlH02cG1cX4HwX+2RJ+2xwD/znGMhRbUPCSUe910D/CGuv0X4D61Y6DpU\nevExquqgO+Fl21DWJnzpvR+NMz4BriZ85WaYmVhfjfDieDaR/l7Cl2KGOWa2MLH9FaHltSrhC/ON\nPHL8IlNmLHdbwhhbPpJyvU3oIsmU9bussronji+Bmb1BeAFvBmxH6GJ7T9J6BCU0sUi5a8ZjXbOO\nnQKsnqhqVtY1Wa4BVnDJvF8TrifxnN7JSpu8LmsDlyVkmhP3d8uTPvs6Fsv7QWI9c58hXJPscpMy\nFXvmFpVrZl/F1Y6J40ZxFtUvqZWk8yRNk/QZQTkY4ZlcjdAKzidvQ0mWg6Q1JI2R9E6s+0ZCV2Mp\ndAU+NrO5WbJ1y5N+NLCnpHaEVvazZjYzyrGepLskvR/l+FMD5EiyNrB11rO+H+EDBGAvQvffjGj4\n0n8Z6mhyUhuodAKStiI80I8uQ/aZwDfAKlnKJUnypfAR4UXZx8zeb2BdHxG+pHsRuu2SvA3caGaH\nNaC8tbLW302U9SczOydPvlwvuYnALwjdg+9Jmkj4ilyZMCZVsNz4Z5xuZus1oM6m4D1CF11GDhG+\nwjO8Tfj6v6VAGWux2Ngm+zoWy5uP91n6/mQo5ZkrRKnXMpluf2Awwfr1LUmdWNyK/pDQalmL0PWZ\nLW9j5TuH0O35AzP7VNIQwvhQvvRJ3gM6S+poZl8mZMv+OAkFmb0q6S1Cq3s/guLKMJLQPb2Pmc2V\ndCxBqeRiLqGLL0Pyg/FtYKKZ7ZJHhknAkGhI8lvC2GxjrmeT4C2qyiMASStK2hW4hfCSfyV5vBSi\nshkPXCxphfjl+X1JP8qTfiHwN+BSSatFObpJyvnQ5sh7XaxrTUmtJf0wfv3dBOwmaZe4fzlJdZLy\nfTkK+E2suzNh/OTWeOxvwP9J6qfA8pJ+JinzRT6L0N2XZCJwFGE8AMJg8lGE7o7Mi6RQuU8DX0g6\nSVKHeA4/kLRlQt5lpVDee4CNJO0RrduOZsmXytXAqZL6AEhaSdIvsso4QVInST1i/lsbkDdbzoys\nY4Gj4/1ZmWD4AjT8mctBrvtXjI4E5fixpOUJyiMjz3fAP4ER8d71IYy/5KOh97Ij4cX/eXyeT8w6\nnvd8YmvoceBcSe0lbUIYN72pQH2jCd162xHGqJJyfAF8JWkDQjdvPl4gtMw6SOpFGC/McDewnqRf\nSmobl60kbRDX95e0UryuXxCUdOq4oqo8/5H0OeHL5hSCFd1BiePGkl9pxb5ADyAM1Gasu25j8csu\nuywI40nTgCdjF8IEINmSKFTfCcDLwDOErqRzCeNZ7xBaBqcSBovfBn5H/ufLgJsJL7w3CIP+fwQw\ns2eBQwkD0B/HYwck8p4LnB67LY6P+x4m/JEziuoxwsB6ZrtguVEJ70roPnyT8JX+V4KhQUbe7OtS\n6Dpl37+cec3sI0JL8DxCi7UXiZa1md0JnA+MiffqZeAnWWX9m/Cl/Tyh2/O6EvPmkimp1O8HXgQm\nEawhk+kb+swlty8Dfi7pY0mXkpvs/P8gjJ28S7C4fCIrzVGE+/8B4fyvI/9/KJd8hY6dSTBk+Az4\nD0tfi1zPY/L4UMK41XsEhfoHM/tvnvohfLj+CHjQzJLDAScQWlmfE57NMTnOK8MlhHG7WYQxuZtY\n/Mx9AewC7Eu4nu/Hc2gX8/4SmB6fmcMIrdnU0eIPTsdxaglJC4FeZvZm2rI4TjnxFpXjOI5T1bii\ncpzaxbtDnBaBd/05juM4VY23qBzHcZyqxhWV02RIOkXS35qgnEUhC1TGsAzRFPf+cpRdLiSNkHRj\nI/KfK+mYppSpUqhIuBFJ/2uAmXxD6m3SciVtIumxpiqvJeBdf07VUQ5rNoV4RG8CbZZxompVIOkM\nwrX51TLkXY1gxv59M/umyYUrM83lHgJIuhsYaWZ3pS1LLeAtKqelUeshGRoj/zDg7lpUUs2Qmwne\n650ScEXlIOlkSbdl7btM0mVxfZikNxQCPr4pab885Szqlkp00xwg6S1JH0o6NZG2laRTFfy3fS5p\nUi5PFpJukHR2XK9T8Ll2vELU4/ckDUuk/Zmk5yV9Junt2PrIkJn8+2msr7+WjgS7jaRnJH0q6WlJ\nP0wcq5d0lqRHY/77JeX0tRY9RdwlaXac2Pqf5LkVKytxzT6SdLoKBIuM5/F4nHD6gqTtc6WLDGSx\n70OiXF8klu8kHVDCtegqaZykOQrRpn+dODZC0m2Sbozn9pKk3rFbeFY8r50T6VeSdG28l+9IOjvT\ntRefkT/HZ+cN4GcFzm2JIIZRjrGSRkU5/idpizz5Rkq6MGvfvxXcFC0RkVqB4fG5/UjSrQreO4h1\nHR/Xuyk8/0fG7e9LmpOoYiKwk0JwVKcYaXvFbQ4LIfbOq4SZ/P8khsXIStOD4NX4FcLs+qNzpPkd\nwQty57i9M8EzwEvxd4dE2n1iff8DzstR1l6xrM1LkH8tgpuYjnG7NWEmfT+Cz7DPgN7x2BokYgFl\nlXMGwR0ULPYifQ3QHtiE4Ctw/Xj8xHhemXI3SZx3zrAMFA9Dsj2wUVzfmOCpYPe4vTZLe6IfxmKv\n0p0JsX72J3zA7UvwurByPF5P8GbRi+Cc9yHg3DzXoTMhltByBI8JY4nhRoqVRfDM/gWwDcH564UE\nLwMZj/kjEte4oeFVZpMnjla8lu/EMotdi7yhK1gcEmXneI9GEcK6nBK3f0301h/TFwoZUjDcSI5z\nWOTJnAaEZiG4K3o7sb0ywUlvlxzlHkNwi9Q13p+rgdHxWEkhZxL1fEbwIZj6O6zal9QFqLWF8LK8\nPmvfziyOP3MeuRVHF2CzuN6R4EBzw8TxHoTQGNNZ/MLeLPFn2Qh4J66vQnAps0rcvoElQw2sEF8m\nj1OCoop5HgF+lTifTNiN5QkvrT2BDkXKGMHSiqpr4vhTwN5x/TVgtzzl5AzLQIEwJHnKuRS4OEue\nfIrqV8CTWfkfBw6M6w8BpyaOHQHcW+K13YzgRZtiZQF/AG5OHOtA8HOXS1E1NLzKEuEwEvvXi9dx\nm2LXgiKhK6J89yeO7UZQvJnx8BXifViR4iFD8oYbyXN+2YqqpNAsBEX2FotDyRxKDJeSo9zJLPlf\nWzNe11aUGHImkfcdYEApz1BLX7zrr+EsZX1iZhNs8eDuU4TwEdlpPjCzF+L6l4QvxWToiouBk7Ly\nvGBmmZAJk4EOsatgXULMn0xXwoMs6Un5bILC/IbSxzRGE14SEL4Gb44yzCW03v6PED7jLknrl1gm\n5A8l0Z3cIUOKkS8MCZK2lvRQ7HL7lDAG0JCQDNnhId5iyXuUPJdkuI4lkPQ9SdfELqPPCN08K0lK\n3ot8ZS0R+sPMvmZxiI5s1qZh4VU+ISiKpKwrEfwFnmZmjydkyHct1qR46IrZWef2kcU3c9yGcL5r\nUzhkSKFwI6VQUmiWKNsYcjz/OegJ/Csh72SC4l7DSg85k2EFQsRipwiuqBpOsRf/wQSv2PkLCNZL\nfQlKDUm7E1pL2eEzkuxFiE8zn9ClsL6ktRW8bg8hKkdJmwPdzCwjw1KKNQ+3AxmP50NIhBgws/EW\nwgJ0IQQbzGeCXmpdEF5AvUpMW2q5ownRi7ubWSfCSy/zjBcr413CizPJ2iwOm9EQfkdopfQzs5UI\nXZJJ7+SFeI/Eh46kDuRXtpnwKisnlhXM7II86V8iRPzNlN2KcM0eNLO/J9IVuhaLQlckjuUNXVGE\nZMiQjPwrmdnG8XihcCNNzS0EZ7lrE7q878iT7m1CV2vymn/PFofNWSLkTNwexpIhZ4j/s3YsDk3i\nFMAVVYlIelLS84SX9GCFQfvnlQiRIek04FszG12gnI4EpXCMmX0p6XsEr+NnJJNl5dmI0EI6HMDM\nPiF0F91K6OKbDnwXv9gvJnhazllWPszsQ8LYyQ2EMYTXYt2rS9pdIbzCfMIXYz7X/w2xSPs7cLak\nXnGAehOFkB+5yiy13I7AJ2b2raR+hC/jjIL6kNBtlC/ExL2E8AdDJbWRtA+wAeGrOClLqXJ8DXwW\nz+mMHGnylXUHIWRKJoTKiAJpGxpe5R6C0szwJ0IgzWNzpMt5LSx4ym9o6IqcWPGQIXnDjTQ1sbfj\nI8JzeZ+ZfZ4n6dXAOZLWgmDyL2lw4ngpIWcg3IcH44enUwRXVCViZv3NrC9hMHicmfWNy3gIlnGE\nyJh53eLHbrs7gJsshGGA8OLsCbwoaTrha/pZSavHPN0JBhq/MrPpCXnuijJtA7welxUIY1n1saz+\nwLjYyiqF0cBOLBmwrRVwHOFreg6hSyNfLByj9BAlFxNeROMJg8p/IxgWZOdrSJlHAmcphFH5PYtj\nM2Eh0uyfgMcULPG2TpYdu1F3JbSGPiIo+11tyVALheRKcilhbOkjwkv93hxpc5ZlIS7ZbwldUe8R\nxndmE1oe2WkbGl7lH8AgSZnrvC+wNfCJFlv+DY3nXOhaFApdkeu6FNouFDKkWLiRQpQiRzajgR1Z\n8vnP5jJgHDA+PmdPEFpgGYqGnInsT1B6TimUexCMYHUzhWDldHKeNJfH4y8CfYvlJVglTSC8nMcD\nnRLHTonppwC7JPbXx33Px2W1uL894YU2FXgSWLvI+dSxtDHFQII1X05rq5hGhBfFJUXKTxpTdIrX\nZEiOdKvH35Xj+fTKkeYhSjSm8KU6F8JLb36x57IB5f2J0JpP/dxa6kKwcH0sbTlqaSn3DWlNGE/p\nSRg0fYGEpVtMMwi4J65vTbQ2KpQXuAA4Ka6fTLSyI1j2vBDT94z5M9ZGOV/ahK/wv8T1fYjmpAXO\naXvguqx9UwmDzRklmCmvK2GCJcAAQtfTC4l0A3OUn1RUpwNfJtI/n1GGhK++V+Kydx5ZXVHV4EKw\nlPseweLyasLYZOpy+eJLWktZXSgpTBI8w8wGxu3hAGZ2XiLN1cBDZnZr3J5CaLWsky9vTLO9mc2S\n1AWoN7MNJJ0CLDSz82Oe+4ARZvakpIeAEyxEek3KeF+s56lomPC+ma2G46SEgr/EnxNa4c8AR5rZ\n1HSlcpz0KPcYVTeWNC/NTCgsJU3XAnnXMLOM6ekswnwMyDLtjetJ8+JR0QAi6eR0Uf1mtoDFA+CO\nkwpmdqgFa7JOZrazKymnpVNuRVVqc60UayrlKs9Ck7CUevY3sx8QjAG2k9Rgp56O4zhO5WlT5vLf\nJcxkz9CDpedbZKfpHtO0zbE/M6dllqQuZvaBpDVZPMEwV1nvAliY04AFk/DRBEudG+PxtQgT89oQ\n3PEkLb0AkFS+PlLHcZxmipk12hF0uVtUk4DeCg5K2xGMFcZlpRlHMFFFUn/g09itVyjvOII7F+Lv\nnYn9+0pqJ2kdoDfwdJxjsmqsoy1hsPrlHGX9nODlISdpDyhWy3LGGWekLkO1LH4t/Dr4tci/NBVl\nbVGZ2QJJRxHmQrQGrjWzVyVlJq5eY2b3SBokaRphMulBhfLGos8Dxko6hODwcu+YZ7KksSx2a3Kk\nmVmcN3JfVFKtCabtGe8K1wI3SppKmCe0bzmvieM4jtMwyt31h5ndS5jwmNx3Tdb2UaXmjfs/JniK\nzpXnHIKTzOS+ucCWedJ/Q1R0juM4TvXhnimcBlNXV5e2CFWDX4uAX4fF+LVoejwUfYlIMr9WjuM4\npSMJqwFjCsdxHMdpFK6oHMdxnKrGFZXjOI5T1biichzHcaoaV1SO4zhOVeOKynEcpwmZPbt4Gqdh\nuKJyHMdpIt5/HzbYAL75pnhap3RcUTmO4zQRN98Me+wB7dunLUnzwhVVA1i4MG0JHMepVszghhvg\nwAOLJnUaiCuqBjBxYtoSOI5TrTz3HHz1FQwYkLYkzY+iikrS5pIulPSUpFmSPojrF0rqWwkhq4Xr\nr09bAsdxqpVRo0JrqpV//jc5BX39SboH+IQQs+lp4H1CpN01CYEHdwM6mdnPyi9qukiylVYy3n4b\nVlwxbWkcx6kmvv0WunWDp5+GddZJW5rqoal8/RUL83GQhSCG2bwZlzGSVm+sELVCXR3cfjscfHDa\nkjiOU03cfTdstJErqXJRsJGaraQkLSepfVaaFjNrYNiwMFjqOI6TxI0oykuxrr9WwBBgKLANQbEJ\n+A54ArgZuLMlxL+QZN9+a3TvDo89Br16pS2R4zjVwIcfQu/eMHMmrLBC2tJUF5UK81EPbAH8GVjX\nzNY0sy7AunHfVkCLsYVr2xb22y8MmjqO4wCMHg277eZKqpwUa1G1j6HaG5WmOZAJnPjii+GhnDHD\nrXscp6VjBptuCpddBjvskLY01UelWlTLS+qcbwEoQZENlDRF0lRJJ+dJc3k8/mLS5D1f3lj/BEmv\nSxovqVPi2Ckx/RRJu+Soa5yklxPbwyR9KOn5uBQ0ldh0U1hlFXjooUKpHMdpCUyaBHPnwvbbpy1J\n86aY1d9zQKHxp4I2LpJaA1cCPwbeBZ6RNM7MXk2kGQT0MrPekrYGRgL9i+QdDkwwswuiAhsODJfU\nB9gH6AN0Ax6QtJ6ZLYx17Ql8kXVOBtxiZkcXuRaLOOigMHi6006l5nAcpzly3XXhfeC9K+WlmNVf\nTzNbJ99SQvn9gGlmNsPM5gNjgN2z0gwGRsX6ngI6SepSJO+iPPF3SFzfnaB05pvZDGBaLAdJHYHj\ngD8SDEIyKGu7KPvtB//5D3z+eUNyOY7TnPjqK7j11mAN7JSXkr4DJO2R1b3WSdKQQnki3YCZie13\n4r5S0nQtkHeNhOn8LGCNuN41pkvm6RrXzyYYgHyVVb8Be0l6SdJtkroXO6lVV4Udd4SxY4uldByn\nuXLHHdC/P3Qv+sZwGkupDdYRZvZpZiOujyghX6lm66W0aJSrvGgaX6geSdqMYLX47xx1/QdY28w2\nASawuKVWkGHD3KWS47RkrrsODjkkbSlaBsXGqDLkUiStS8j3LtAjsd2DJVs8udJ0j2na5tj/blyf\nJamLmX0gaU0gM+k4X1n9gS0lTSec8+qS/mtmO5rZx4n01wIX5DuZESNGLFrfbrs6pk+vY/Jk6NMn\nXw7HcZojb7wBr7wSLICdxdTX11NfX9/k5RY0T1+USLqe4PPvKoLS+g2wspkNK5KvDfAasBPwHsFf\n4NAcxhRHmdkgSf2BS82sf6G8ki4A5pjZ+ZKGE/wNZowpRhPGpboBDxAMNSxR39rAXWa2cdzuYmYf\nxPU9gBPNbJsc57LUvOZTTw0B0i66qOgldBynGXHaafD113DxxWlLUt00lXl6qYqqI/B7gtKA0EX2\nRzObW0LenwKXElpg15rZuZIOBzCza2KaK4GBwFyCf8Hn8uWN+zsDY4G1gBnA3pmuSUmnAgcDC4Bj\nzOz+LHl6AuNiVx+SziEYZywA5gBHmNnrOc5jKUX1xhvwwx+GGekeKM1xWgYLFsDaa8P998MPfpC2\nNNVNRRWVk1tRQTBRP/xw2HvvFIRyHKfi3HMPnHkmPPVU2pJUP5Wa8JupbH1Jf4uTbB+Ky38bW3lz\n4Ne/hr//PW0pHMepFNde60YUlabUrr+XCBNxnyM4pIVgcPdsGWWrKvK1qObNgx49PA6N47QEZs+G\n9dbD49KVSEVbVMB8MxtpZk+Z2aS4tBglVYjlloP99w+mqo7jNG9uugl2392VVKUptUU1AvgQ+Cew\nyLdflml3syZfiwrgf/+DgQODo9o2pRr8O45TU5gF44mRI+FHP0pbmtqgUhF+MwwjTKo9IWu/d3YR\nHt4ePeC++2DXXdOWxnGccvDkkzB/Pmy3XdqStDxKUlRm1rPMctQ8v/41/O1vrqgcp7lyzTVw2GGg\nRrcPnIZSLB7VTmb2oKS9yO2+6J/lFK6aKNT1B/Dll6FV9cor0LVr3mSO49Qgn3wSjKWmTQu+Pp3S\nqJQxRaYndrc8ixPp2BF+8QuP/us4zZEbb4Sf/tSVVFr4hN8SKdaigmCiPnQoTJ3q8Wkcp7mQMaK4\n6iqoq0tbmtqiIi2qGP027ziWpHaSDmqsEM2FrbYKLSuP/us4zYfHHoPvvvMovmlSzJiiIyGy7hTg\nGeADglPaLsCWwAbA38oqYQ0hwaGHhkFXj/7rOM0DN6JIn6Jdf5IEbAsMIDiBBXgLeBR4vGh/WDOh\nlK4/gM8+g549YfJkWHPN8svlOE75mDMHvv/94IB6lVXSlqb2cKe0FaZURQXh62utteD008sslOM4\nZeWSS+DZZ4NHCqfhVERRSboisWksDqBoAGZ2dGMFqBUaoqiefz64WZk+HVqXEl7ScZyqwww23DDM\nj/RJvstGpczTn41Le2Bz4HVgKtAXaNfYypsrffuGuVR33522JI7jLCsPPxysdwcMSFsSp1Rff08B\nA8xsftxuCzxqZluXWb6qoSEtKgjzqcaMgXvvLaNQjuOUjf32g623hmOOSVuS2qXSEX5fA7Yxszlx\nuzPwhJmt31gBaoWGKqqvvw7jVE89BeuuW0bBHMdpcj76CHr1Ct33K6+ctjS1S6XDfJwHPCdplKRR\nhLhU5za28uZMhw5wwAHBtNVxnNrihhvCOLMrqeqgZKs/SWsCWxMMKZ4ysw/KKVi10dAWFQQPFdtu\nCzNnQvv2ZRLMcZwm5bvvQnDE0aND15+z7FS6RQUwD3gf+BRYT1JJEVkkDZQ0RdJUSSfnSXN5PP6i\npL7F8krqLGmCpNcljZfUKXHslJh+iqRdctQ1TtLLie32km6NeZ6UtHZJV6MEeveGTTeF229vqhId\nxyk3990HnTtDv35pS+JkKElRSToUeBi4DxgB3B9/i+VrDVwJDAT6AEMlbZiVZhDQy8x6A4cRQt4X\nyzscmGBm6wEPxm0k9QH2iekHAn+R1CpR157AFyzpCf4QYE6s/xLg/FKuSakccUQItOY4Tm1w1VXw\nm9+4J4pqotQW1TFAP+AtM9uBYJ7+WQn5+gHTzGxGtBgcA+yelWYwMArAzJ4COknqUiTvojzxd0hc\n3x24xczmm9kMYFosB0kdgeOAP7J4Plh2WXcATer8aPDgMCD78svF0zqOky7TpsGkSbDPPmlL4iQp\nVVHNM7MMwjukAAAgAElEQVSvASQtZ2ZTgFIs/roBMxPb78R9paTpWiDvGmY2K67PAtaI611jumSe\nTHSos4E/A1/lq9/MFgCfRavGJqFNm+D/z1tVjlP9jBwJBx0UjKGc6qHUUPTvSFoZuBOYIOkTYEYJ\n+Uq1Piilka1c5ZmZSSpUjyRtBqxrZsdJ6lmiTEsxYsSIRet1dXXUlejz/9BDYeON4fzzYYUVlrV2\nx3HKydy5Yf7jpElpS1K71NfXU19f3+TllhqKPtO1NkJSPbAiYbyqGO8CPRLbPViyxZMrTfeYpm2O\n/e/G9VmSupjZB9EacXaRsvoDW0qaTjjn1SX918x2jHnWAt6LIU1WMrOPc51MUlE1hG7d4Mc/Dn+C\no45apiIcxykzo0fDNtsEp9LOspH9AX/mmWc2SblFu/4ktYlhPgAws3ozG2dm35ZQ/iSgt6SektoR\nDB3GZaUZBxwQ6+oPfBq79QrlHQccGNcPJLT0Mvv3jXGy1gF6A0+b2dVm1s3M1iF4gX89Kqnssn5O\nMM5ocn77W7jiCli4sBylO47TGMyCEYV/SFYnRVtUZrZA0muS1jaztxpSeMx7FMFKsDVwrZm9Kunw\nePwaM7tH0iBJ04C5wEGF8saizwPGSjqE0AW5d8wzWdJYYDKwADgyx+Sn7C7Ea4EbJU0F5gD7NuQc\nS2XAAPje92D8eBg4sBw1OI6zrDz2WPAm8+Mfpy2Jk4tSXSg9QrD0e5qgTCAMDw0uo2xVxbJM+M3m\n+uvhttvgnnuaSCjHcZqEoUOhf3/369fUVNrXX12O3WZmExsrQK3QFIrq669h7bXD11vv3k0kmOM4\njeL996FPnzCNpFOn4umd0qmqwImSnjCzHza6oCqmKRQVwKmnBuuiyy5rAqEcx2k0Z50F770HV1+d\ntiTNj2pTVM+bWd/iKWuXplJUM2cGt0ozZsCKKzZeLsdxlp1vv4V11glukzbeOG1pmh9p+PpzmoAe\nPRabqjuOky5jx4Yovq6kqhtXVCngpuqOkz5mcMklcOyxaUviFMMVVQoMGADLLx9M1R3HSYdHH4Uv\nvoBBg9KWxClGyYpKUhdJu0naVdLqWYcPaGK5mjUSHH00XH552pI4Tsvl0kuDOXor/1yveko1T98b\nuBDImKP/CDjRzG4ro2xVRVMZU2SYNy+Eqn/00RCkzXGcyjF9Omy1VTBq6tgxbWmaL5WeR/US8GMz\nmx23VwMeNLNNGitArdDUigrgtNPg00+D6xbHcSrH8ceHyAYXXJC2JM2bSiuql4FNMm/qGIzwRTNr\nMbYy5VBUmYmG06bBKqs0adGO4+Th88+DSfrzz4deDad8VNo8/T7gfknDJB0E3APc29jKWzprrglD\nhvhEQ8epJNddF6aIuJKqHUptUQnYk+B53IBHzOxfZZatqihHiwpC5N9ddgl95e3bN3nxjuMk+O67\n4L5s9Ojg288pL03Voio1HpURwrTf0dgKnSXZeGPYZJPwxznooLSlcZzmzbhxsPrqrqRqjYJdf5Ie\ni79fSvoia/m8MiI2f373O7j44jAB0XGc8nHppXDccWlL4TSUgorKzLaNvx3NbIWsxT3VNRE77xzm\nVvkEYMcpH88+C2++CXvumbYkTkMpyZhC0o2l7HOWDSmYy150UdqSOE7z5cILg7uktm3TlsRpKKUa\nUyzhHV1SG+AlM+tTTuGqiXIZU2T45pvFXpw3aTGz0xynMrz5ZpjgO326Ry2oJBUxT5d0qqQvgI2T\n41PAbGBcYyt3FtO+PRx1VBirchynabn4YjjsMFdStUqpLarzzGx4BeSpWsrdogL4+GP4/vfhlVeg\na9eyVuU4LYaPPgpuyl55JcxddCpHRSf8mtlwSStL6ifpR5mllLySBkqaImmqpJPzpLk8Hn9RUt9i\neSV1ljRB0uuSxkvqlDh2Skw/RdIuif33SXpB0iuSrpXUNu4fJulDSc/H5eBSzqscdO4Mv/wlXHll\nWhI4TvPjqqtgr71cSdUypbaoDgWOBnoAzwP9gSfMbMci+VoDrwE/Bt4FngGGmtmriTSDgKPMbJCk\nrYHLzKx/obySLgA+MrMLogJbOSrTPsBoYCugG/AA0NvMTFJHM/sy1nk7cKeZ3STpQGALMzu6yLmU\nvUUFoS+9X7/w690UjtM4vvoKevaERx6B9ddPW5qWR6VdKB0D9ANmmNkOQF/gsxLy9QOmmdkMM5sP\njAF2z0ozGBgFYGZPAZ0kdSmSd1Ge+Dskru8O3GJm881sBjAN2DqWnVFSbYF2wEcxj+JSFay7bjBX\n/+tf05bEcWqf66+Hbbd1JVXrlKqo5pnZ1wCSljOzKUApt74bMDOx/U7cV0qargXyrmFms+L6LGCN\nuN41pstZn6T7Y/qvzey+uNuAvSS9JOk2Sd1LOK+yMnx4iDz6zTdpS+I4tcuCBWHKx0knpS2J01hK\ncqEEvCNpZeBOYIKkT4AZJeQrta+slBaNcpUXu/UK1WOJtD+R1B64VdKBZjYK+A8w2szmSzqM0ELb\nKVdBI0aMWLReV1dHXV1dCWI3nE03hc02g1GjgqWS4zgN5447glHSD3+YtiQth/r6eurr65u83JLG\nqJbIINUBKwL3mdm3RdL2B0aY2cC4fQqw0MzOT6S5Gqg3szFxewqwPbBOvrwxTZ2ZfSBpTeAhM9tA\n0nAAMzsv5rkPOCN2KSbl+hWwtZkdlbW/NTDHzDqRRaXGqDI88kjw/ffaa9C6dcWqdZxmgRlsuSWc\ncQYMHpy2NC2Xio1RSWoTFQMAZlZvZuOKKanIJKC3pJ6S2gH7sPT8q3HEUPZRsX0au/UK5R0HHBjX\nDyS09DL795XUTtI6QG/gaUnLR4WWmay8K8EohDgelmEwMLmE8yo7220HXbrA7benLYnj1B7//W8w\npNh117QlcZqCol1/ZrZA0muS1jaztxpSeMx7FHA/0Bq4NlrtHR6PX2Nm90gaJGkaMBc4qFDeWPR5\nwFhJhxC6IPeOeSZLGktQNguAI2PX4PLAv2O3n2KZ18WyjpY0OKafAwxryDmWk1NOCVGA9947uFly\nHKc0/vSnMNbbqtRReKeqKdU8/RGCpd/TBGUCYXioxTSqK931B6H7YtNNQ7jsgQMrWrXj1CyPPRbm\nI77+uvv1S5um6vorVVHV5dhtZjaxsQLUCmkoKghxqq65Bia2mCvtOI1j0CDYfXc4/PC0JXEqqqhK\nEOYJM2vWtjVpKaoFC4L7l5tugm22qXj1jlNTPPtsUFJvvOERs6uBSk/4LcZyTVSOk0WbNnDiiXDu\nuWlL4jjVzznnhP+LK6nmRVO1qJYIA9IcSatFBTBvXnBWe9dd0LdZX2XHWXZeeQV22im4H/ve99KW\nxoHqa1E5ZWS55cJX4tlnpy2J41Qv55wTAiO6kmp+eIuqRNJsUQF8/XVoVXlgRcdZmmnTggeKN95w\nZ87VRLW1qA5oonKcPHToACecAGedlbYkjlN9nHce/OY3rqSaKwVbVJIeM7NtJX3J0n72zMxazGOR\ndosKwkz7ddeFCRNg441TFcVxqoa33w5jt1OnhphuTvVQVebpLYFqUFQAF14IkybBrbemLYnjVAdH\nHgkrrADnn188rVNZKqKoJK1oZp9LyvmdYmYfN1aAWqFaFNXcuaFV9dBD0KdP2tI4Trq89RZsvnlw\n3rzqqmlL42RTKUV1t5n9TNIMcofYWKexAtQK1aKoIHw5vvAC3HJL2pI4Trocdhistlrw7edUH971\nV2GqSVF9+WVoVU2cCBtumLY0jpMOb74J/foFn34+NlWdVKpFtXmhzGb2XGMFqBWqSVFBmDPyyitw\n881pS+I46XDwwdCjB5x5ZtqSOPmolKKqJ3T5dQC2AF6KhzYBJjV3/35Jqk1Rff459OoF9fU+VuW0\nPKZNg/79w2+npcKcOtVCReZRmVmdme0AvAdsbmZbmNkWhJAf7zW2cmfZWXHF4K3iD39IWxLHqTxn\nnw1HH+1KqqVQapiPyWbWp9i+5ky1taggzKvq3RvGjYMttkhbGsepDK+9BgMGhNbUSiulLY1TiErH\noxoDfAncRIiQux/Q0cyGNlaAWqEaFRXAVVcFZ7X33pu2JI5TGfbfHzbaCE49NW1JnGJUWlF1AI4A\ntou7HgZGmtm8xgpQK1Srovr2W1h/fbjxxvCV6TjNmVdegR12CD79VlghbWmcYrh5eoWpVkUFcMMN\ncN11wVxdjX4kHKd6GTIEttsOfve7tCVxSqGiTmklrSfpdkmTJU2Py5sl5h0oaYqkqZJOzpPm8nj8\nRUl9i+WV1FnSBEmvSxovqVPi2Ckx/RRJuyT23yfpBUmvSLpWUtu4v72kW2OeJyWtXcp5VRO//CXM\nng3jx6ctieOUjyeegOeeC85nnZZFqd7TrweuBhYAOwCjgKIzeCS1Bq4EBgJ9gKGSNsxKMwjoZWa9\ngcOAkSXkHQ5MMLP1gAfjNpL6APvE9AOBv0iL2hg/N7PNzGwjYKWYDuAQYE6s/xKg5jyGtWkTrKBO\nOw2qtNHnOI3CDIYPhxEjQnw2p2VRqqLqYGYPELoKZ5jZCOBnJeTrB0yLeeYDY4Dds9IMJig+zOwp\noJOkLkXyLsoTf4fE9d2BW8xsvpnNAKYBW8eyvwSILal2wEc5yroD2KmE86o69toLvvsO/vWvtCVx\nnKbn3nvhww/hAA8o1CIpVVHNiy2caZKOkrQnsHwJ+boBMxPb78R9paTpWiDvGmY2K67PAtaI611j\nupz1Sbo/pv/azO7Lrt/MFgCf5XPCW820ahX8nZ1+OixYkLY0jtN0LFwIp5wSvLG0aZO2NE4alHrb\njwG+BxwNnA2sCBxYQr5SO6JKGWxTrvLMzCQVqscSaX8iqT1wq6QDzWxUgXxLMWLEiEXrdXV11NXV\nNSR72fnpT+GCC+D66+HQQ9OWxnGahltuCeHld8/ui3Gqjvr6eurr65u83KKKKrak9jGzE4AvgGEN\nKP9doEdiuwdLtnhypeke07TNsf/duD5LUhcz+0DSmsDsAmW9m9jGzL6RdAehS3BUPL4W8J6kNsBK\n+cKXJBVVNSKFeFVDhsB++8HypbR5HaeK+fZb+P3vw8eXW7RWP9kf8Gc2kSPGol1/ZvYdMCBhlNAQ\nJgG9JfWU1I5gwDAuK804Yih7Sf2BT2O3XqG841jcojsQuDOxf19J7SStA/QGnpa0fFRoRGW0K/B8\njrJ+TjDOqFm22gp+9CO4+OK0JXGcxvPXv8IGG8D226ctiZMmpU74vZow/nMb8FXcbWb2zxLy/hS4\nFGgNXGtm50o6PBZwTUyTse6bCxyU8cqeK2/c3xkYS2gJzQD2NrNP47FTgYMJForHmNn9klYH7gLa\nE7oQ7wdOit2G7YEbCf4L5wD7RkOM7POo2nlU2bz5ZlBYkyfDGmsUT+841cgXX8B66wVDis02S1sa\nZ1motGeKG8g9PnRQYwWoFWpJUQEcfzzMmwd/+UvakjjOsnH66fD22/CPf6QtibOsuGeKClNrimrO\nnNBl8sgj4ddxaom334a+feHFF6F797SlcZaVinqmcGqPVVaBk04KZr2OU2uceiocdZQrKSfgLaoS\nqbUWFYSuv/XXD1GA3WGtUys8/TTssUcI59GxY9rSOI3BW1ROUZZbLkySPP74MGnScaods/C8/vGP\nrqScxTRYUUm6qxyCOOVh6FBo3TqEAXGcaueOO2DuXHeV5CxJg7v+JD1vZn2Lp2xe1GLXX4ZMV8qU\nKR7Dx6levvkGNtwQ/v532HHHtKVxmoI0u/6eL57EqSb69YOddw7dgI5TrVxxBWy8sSspZ2ncmKJE\narlFBfD+++El8OST0KtX2tI4zpJ88EF4Ph99NBgAOc0Dn0dVYWpdUQGcf34IPnfnncXTOk4lGTYM\nVl89OFV2mg+uqCpMc1BU33wDG20EI0eGrkDHqQYefxz23htefdXHUJsblQ5F/4tS9jnVTfv2cNFF\ncOyxHrPKqQ6++y5M7L3wQldSTn5KNaY4tcR9TpUzeDB07QpXXZW2JI4TvKOvsALsu2/akjjVTMGu\nv+i9fBAhxMYYFgc4XAHoY2b9yi5hldAcuv4yTJkC220X/Kh17Zq2NE5L5aOPoE8feOAB2GSTtKVx\nykFFxqgkbUoIf3EW8HsWK6rPgYfM7JPGClArNCdFBXDaafDGGzBmTNqSOC2Vww6DDh3gssvSlsQp\nF5UO89HWzOY3trJaprkpqq++gh/8AK65xg0rnMrzzDOhG/rVV6FTp7SlccpFpVpULxfIa2bWYhrs\nzU1RAdx9dzCsePnl4BfQcSrBggVhEvoxx8CBBxZP79QulVJUPQtlzhUJt7nSHBUVwJ57wqabwhln\npC2J01K4+OLwkfTAA6BGv8Kcaqaq5lFJesLMftjogqqY5qqoZs4MAercY4VTCd56C7bYIkw87907\nbWmcclNtYT6846hG6dEDhg+H3/wmhFhwnHJhFp6zY491JeU0jLLHo5I0UNIUSVMlnZwnzeXx+IuS\n+hbLK6mzpAmSXpc0XlKnxLFTYvopknaJ+zpIulvSq5L+J+ncRPphkj6U9HxcDi7Plahejjkm+Fq7\n6aa0JXGaM7ffDtOnh8jTjtMQyqqoJLUGrgQGAn2AoZI2zEozCOhlZr2Bw4CRJeQdDkwws/WAB+M2\nkvoQ5nz1ifn+Ii3qBb/AzDYkmNtvK2lg3G/ALWbWNy7XNfV1qHbatoXrroMTToBZs9KWxmmOfPpp\naEldcw20a5e2NE6tUe4WVT9gmpnNiObtY4Dds9IMBkYBmNlTQCdJXYrkXZQn/g6J67sTlM78aOgx\nDdjazL42s4mxjvnAc0C3mEcsnh/WYtlii+AY9Oij05bEaY6ccgrsuisMGJC2JE4tUrKiktRF0m6S\ndpW0etbhfPE4uwEzE9vvsFhBFEvTtUDeNcws8+0/C1gjrneN6fLWF7sJdyO0xCC0qPaS9JKk2yR1\nz3MuzZ4RI+D55927utO0TJwI48bBeeelLYlTq7QpJZGkvYELgYlx15WSTjSz2wDMLN98q1KH50tp\n0ShXeWZmkgrVs+iYpDbALcBlCdP6/wCjzWy+pMMILbSdchU0YsSIRet1dXXU1dWVIHbt0KFDiK66\n335QV+cTMZ3GM3cuHHxw8Ni/8sppS+OUm/r6eurr65u83FI9U7wE/NjMZsft1YAHi034ldQfGGFm\nA+P2KcBCMzs/keZqoN7MxsTtKcD2wDr58sY0dWb2gaQ1Ce6cNpA0HMDMzot57gPOiF2KSLoO+NzM\njs0jb2tgjpkt9YpurubpuTjySJg/H/72t7QlcWqdo4+GTz6BG29MWxInDSptni7gw8T2HEprBU0C\nekvqKakdwdBhXFaaccSuw6jYPo3deoXyjgMyc9oPBO5M7N9XUjtJ6wC9gadj2X8EVgSOW+LEwnhY\nhsHA5BLOq1lz3nkwfjw8+GDxtI6Tj4kT4Y473Jef03hK6voD7gPulzSaoKD2Ae4tlsnMFkg6Crgf\naA1ca2avSjo8Hr/GzO6RNEjSNGAucFChvLHo84Cxkg4BZgB7xzyTJY0lKJsFwJGxa7A7ISzJq8Bz\n0RDwimjhd7SkwTH9HGBYidek2bLiisE66+CD4aWXYKWV0pbIqTXmzoVDDgldfp07py2NU+uU2vUn\nYE9gAGHM5xEz+1eZZasqWlLXX4YjjgjOa0eNKp7WcZJ4l58DVeZCqSXQEhXV3Lmw2WZw/vnBJ6Dj\nlEJ9Pey/f3B27K2plk2lnNI+ZmbbSvqSpS3uzMxWbKwAtUJLVFQQfAAOGQIvvABduhRP77RsPvkk\nfNxcfTX89KdpS+OkjbeoKkxLVVQAp58eFNV//uPerp38mMHQobDaanDFFWlL41QDFbX6k7RUT3Ou\nfU7z5A9/gPfeC3OsHCcfN90UuvsuuCBtSZzmRqnGFM+bWdJZbBvgJTPrU07hqomW3KICmDwZtt8e\nHn/cPV87SzN9egiGOGFC6PpzHKhQi0rSqZK+ADaW9EVmAWaz9HwopxnTpw+ceSbssw98803a0jjV\nxIIF8Ktfwcknu5JyykOpLarzzGx4BeSpWlp6iwrCGMRee4UYVj6J08lw1llhcu+ECdCq7IGDnFqi\n4sYUklYmeHpYFCTRzB5urAC1giuqwCefhIjAl10Gu2f7wXdaHA89FHxDTpoE3bLdTTstnooqKkmH\nAkcDPYDngf7AE2a2Y2MFqBVcUS3m8cdhjz3gmWdgrbXSlsZJiw8+COFhbrgBdt45bWmcaqTSvv6O\nIcSHmmFmOxCCD37W2Mqd2mSbbeC448KX9IIFaUvjpMF334X7/+tfu5Jyyk+pimqemX0NIGk5M5sC\nrF8+sZxq56STYPnl4bTT0pbESYOzzgpz6v7wh7QlcVoCpTqlfSeOUd0JTJD0CcEZrNNCadUqzJvZ\ncstglrzXXmlL5FSK8ePDnLpnn4XWrdOWxmkJNNgzhaQ6QriM+8zs23IIVY34GFVuJk0KrnIefhg2\n3DBtaZxy8/bbsPXWMHo07LBD2tI41U7FjCni5N7/mdkGja2slnFFlZ/rrw+Oa59+OoQIcZonX30F\nAwYEh7O/+13a0ji1QKWt/v4NHG1mbzW2wlrFFVVhjjgCZs0KgfLcH2Dzwwx++ctwb2+80e+xUxqV\nVlSPECz9niYEN4TgPX1wYwWoFVxRFeabb4KLpcGD4dRT05bGaWr+/Ge45RZ49FHo0CFtaZxaoakU\nVanGFL/Psc/f2s4i2rcPran+/WGDDTx+VXNi/Hi46CJ46ilXUk46NEmYD0lPmNkPm0CeqsVbVKXx\n7LMwcCDcd1+YDOrUNlOmhJby2LHh13EaQqUn/BZjueJJnJbAFlvAX/8agi2++27a0jiNYfZs+NnP\ngqGMKyknTcruQlLSQElTJE2VdHKeNJfH4y9K6lssr6TOkiZIel3SeEmdEsdOiemnSNol7usg6W5J\nr0r6n6RzE+nbS7o15nlS0trluRIthz32gKOOgt12C+Hsndrj66+DL8f99oNhw9KWxmnplFVRSWoN\nXAkMBPoAQyVtmJVmENDLzHoDhwEjS8g7HJhgZusBD8ZtJPUB9onpBwJ/kRbZJ11gZhsSjEK2lTQw\n7j8EmBPrvwQ4v2mvQsvkpJNCyAd3s1R7LFwYwnass07wQOE4aVPuFlU/YJqZzTCz+cAYINvn9mBg\nFICZPQV0ktSlSN5FeeLvkLi+O3CLmc03sxnANGBrM/vazCbGOuYDzwHdcpR1B7BTk5x5C0eCq6+G\nefOC6boP79UOw4eHqQbXX+9m6E510FSK6oA8+7sBMxPb77BYQRRL07VA3jXMbFZcnwWsEde7xnR5\n64vdhLsRWmJL1G9mC4DPJHXOcz5OA2jXLlgCvvgi/D6X3ahTdVx0EYwbB3feGSw5HacaKGieLukx\nM9tW0pcsbY5uZrZiXHk5TxGlfkeX8t2mXOWZmUkqVM+iY9HLxi3AZbHF1SBGjBixaL2uro66urqG\nFtHi6NgR7r4btt0W1lgDfvvbtCVy8nHddXDFFfDII7DKKmlL49Qi9fX11NfXN3m5BRWVmW0bfzsu\nY/nvEmJYZejBki2eXGm6xzRtc+zP2JHNktTFzD6QtCYwu0BZSduzvwKvmdnlWfWvBbwXFdlKZvZx\nrpNJKiqndFZbLczFGTAgrO+7b9oSOdn8859w+ulQXx8iODvOspD9AX/mmWc2SbkFu/4krRh/O+da\nSih/EtBbUk9J7QiGDuOy0owjdh1K6g98Grv1CuUdBxwY1w8keHXP7N9XUjtJ6xAiEj8dy/4jwZnu\ncTnqz5T1cxZ3CTpNSM+ecO+9cOyxoVvJqR4eeAD+7/9Cy3e99dKWxnGWpphniluAnxGMD3J1r61T\nKLOZLZB0FHA/0Bq41sxelXR4PH6Nmd0jaZCkaQT3TAcVyhuLPg8YK+kQQriRvWOeyZLGApOBBcCR\nsWuwO3Aq8CrwXDQEvMLMrgOuBW6UNBWYA/j3fpnYeOPwMhw0CNq0gV13TVsiZ+LEYJl5xx3Qt2/x\n9I6TBk3imaIl4J4pmo6nnw5KatSoECLESYf6evjFL+DWW2HHHdOWxmmOVMQpraTNC2U2s+caK0Ct\n4IqqaXniiTCh9KabYJdd0pam5eFKyqkElVJU9YQuvw7AFsBL8dAmwKTm7t8viSuqpufRR4Pz2r//\nPXhddyrDf/8bDFrGjgU3XHXKSUV8/ZlZnZntALwHbG5mW5jZFgTvDu81tnKnZTNgQBizOuwwuPnm\ntKVpGfz730FJ3XabKymndig1zMcGyblSZva/bFdIjrMsbLUVPPhg8Lj++efBi4VTHq6/PsQKu+ce\n2HLLtKVxnNIpVVG9JOnvwE2Eibf7AS+WTSqnRbHRRsH6bOed4aOPwnwed93TtFx4IVx1VbjOboLu\n1BqlRvjtABwBbBd3PQyMNLN5ZZStqvAxqvLz/vvBGnDTTYOfwHbt0pao9vnuOzjhhDDh+v77oXv3\ntCVyWhIVDUXvuKKqFHPnwtCh4feOO6BTp+J5nNx8/nmYIzVvXhiTWnnltCVyWhoVDZwoaT1Jt0ua\nLGl6XN5sbOWOk83yy8O//gU/+AFssw288UbaEtUmM2YE/4rduwePIK6knFqmVO/p1wNXE7w97EAI\ni+F2Wk5ZaN0aLrssBF/84Q/hrrvSlqi2mDgxKPlf/xpGjoS2bdOWyHEaR6ljVM+Z2eaSXjazjZP7\nyi5hleBdf+nw+OOw995wyCFwxhnQquwxqWuXhQuD0cQllwSvHz/5SdoSOS2dpur6K9Xqb16MuDst\n+t97D1i+sZU7TjG22QYmTQrK6umnwwt49dXTlqr6+OQTOOAAmDMHnnnGPaA7zYtSv0+PAb4HHA1s\nCfySxR7HHaesdOkS5lptumkIb3/vvWlLVF08/DBsvjn06uVhOpzmSdGuv9iSOt/MTqiMSNWJd/1V\nB/X1cOCBweXSBRdAhw5pS5Qe8+aFyMk33wx//at7o3eqj4pZ/ZnZd8AAyadgOulTVwcvvBAmBm++\neQI+aCQAABGJSURBVIhG2xJ59tng1ePNN+HFF11JOc2bUo0prga6ArcBX8XdZmb/LKNsVYW3qKqP\nf/0rhLYfNAjOP79lmGB/9lloRY0dC3/+M+y/v3vxcKqXis6jApYjBBXcEdg1Lrs1tnLHaQx77AGv\nvBLMrzfaKBhaLFyYtlTlwSyE5NhoI/j663Dev/ylKymnZeCeKUrEW1TVzVNPwXHHhXGbiy6CHXZI\nW6Km49FH4cQTg4K66qowkddxaoGKtKgkjZC0RoHja0o6s7FCOE5j2XpreOwxGD48zLnadddgzl7L\nvPwyDBkSuvd+8xt47jlXUk7LpFjX3yRgjKTHJF0h6VRJp8X1xwjeKZ4qVICkgZKmSJoq6eQ8aS6P\nx1+U1LdYXkmdJU2Q9Lqk8ZI6JY6dEtNPkbRLYv+fJL0t6YusuodJ+lDS83E5uMg1caoUKcy3evXV\nEOL+5z8Pk14ffTRtyRrGE08Eq8ZddoHttoPXXgvdfD7Z2WmplGpM0QPYFlgr7noLeMzM3imSrzXw\nGvBj4F3gGWComb2aSDMIOMrMBknaGrjMzPoXyivpAuAjM7sgKrCVzWy4pD7AaGAroBvwANDbzExS\nP+BtYKqZrZCo/0BgCzM7usi5eNdfjfHtt/CPf8C558Kqq8KRRwZFVo0m7d98E4xDRo6EmTNDV9+w\nYdUpq+OUSkWNKcxsppmNMbMLgIuAe4spqUg/YJqZzTCz+cAYYPesNIMJvgMxs6eATpK6FMm7KE/8\nHRLXdwduMbP5ZjYDmAZsHct+2sw+yCGj4uI0M9q1C/7uXn89WMrdeiustRYcf3zwdpH2d4dZ6N47\n8cQwSffvfw/+DV9/PQSQdCXlOIFSvaffImlFScsDLwOvSjqphKzdgJmJ7XfivlLSdC2Qdw0zmxXX\nZwGZcbSuMV2h+rIxYC9JL0m6TZJH7GlmtG4dxqzuuQeefDJ4aB86FHr3htNOC11tCxZURpaFC8Mc\nqN//HjbcMMjVqlUYX3vgAfjFL6BNqY7NHKeFUGqvdx8z+5zQcrkX6An8qoR8pX6zltKiUa7yYn9c\noXqKyfAfYG0z2wSYwOKWmtMM+f734eyzQ6tl7NigoI44AlZbDfbcE664Iiiur74qXlYpzJsXjDr+\n8peghFZfHX71q2DBN2pUCMdx/vlBaTqOk5tSv93aSGpLUFRXmdl8SaUooXeBpOexHizZ4smVpntM\n0zbH/nfj+ixJXczsA0lrArMLlPUuBTCzjxOb1wIX5Es7YsSIRet1dXXU1dUVKtqpYqTg2WLzzYOi\n+OCD0KJ5+GG44YZgkLHuukGBrLtuWFZbDTp3DhOL27dfXM68efDpp8Ex7OzZMH16WKZODct668EW\nW8Buu8Gll0K3Ym18x6lR6uvrqa+vb/JySzWmOBo4GXgJ+BnBqOJGM9uuSL42BIOInQge15+msDFF\nf+DSaEyRN280pphjZudLGg50yjKm6MdiY4peSSsISV9kGVN0yYxdSdoDONHMtslxLm5M0YKYNy8o\nqzfeCG6Kpk8Pbps+/jgopG+/XTzG1a5dUF4rrxyMNtZdF9ZZJ7TeNtoIllsu3XNxnLRINRR99PvX\nJho5FEv7U+BSoDVwrZmdK+lwADO7Jqa5EhgIzAUOMrPn8uWN+zsDYwkKcwawt5l9Go+dChxMCPJ4\njJndH/dfAAwF1gTeB/5mZmdJOodgnLGA4H3jCDN7Pcd5uKJyHMdpABVVVJJWBc4ABhDGfB4BzjKz\nOY0VoFZwReU4jtMwKu3rbwxhHGhP4OfAh8Ctja3ccRzHcYpRaovqf2b2g6x9i8LStwS8ReU4jtMw\nKt2iGi9pqKRWcdkHGN/Yyh3HcRynGAVbVJK+ZPE8pOWBTBCFVsDcpPVcc8dbVI7jOA2jqVpUBedR\nmVnHxlbw/+2deZBU1RWHv5+jiCAuY4yCUYgIilkskQiRUGJMEDXRuGup5VYmalliVeKW0oqxKhFM\nYhJU0BijZEHRhCQQN0Yj7kKpoAgMiwoaQbSIGGJcAE/+uKfh0XbP9PRMOz3d56t61ffdvue+8w53\n+nDfu/ecIAiCIGgPJQdrkbQjMICURBEAM3usEkoFQRAEQY6SHJWkc4GLSFEf5gDDgKdJGX+DIAiC\noGKUuphiDCnawzIzOwTYH3i3YloFQRAEgVOqo/rAzN4HkNTdzJqBvSunVhAEQRAkSn1H9bq/o/ob\n0CTpHVLooiAIgiCoKG2O9SdpJLAd8ICZfVQJpaqRWJ4eBEHQNjo1KG09Eo4qCIKgbXzakSmCIAiC\noFMIRxUEQRBUNeGogiAIgqomHFUQBEFQ1YSjCoIgCKqacFRBEARBVROOKgiCIKhqKu6oJI2W1Cxp\niaTLirQZ79+/IGn/1mQlNUpqkrRY0gxJO2S+u8LbN0salan/iaTXJK3Nu/bWkqa4zDOS+nasBYIg\nCIL2UFFHJakBuBEYDewLnCJpUF6bI4C9zGwA8F1gYgmylwNNZjYQeNjPkbQvcJK3Hw1MkJTbbPZ3\nUmDdfM4BVvv1fwmM64Bbr2lmzpzZ2SpUDWGLRNhhE2GLjqfSM6oDgaVmtszM1gF3AUfntTkKmARg\nZrOAHSTt2orsRhn//I6XjwbuNLN1ZrYMWAoM9b5nm9mbBXTM9vUX4NB23G9dEH+ImwhbJMIOmwhb\ndDyVdlS7Aa9nzv/ldaW06dOC7C5mtsrLq4BdvNzH27V0vaI6mtl64F1Jja3IBEEQBJ8SlXZUpQbH\nKyUWlAr15wH4WrpOBOgLgiDoyphZxQ5SJuAHMudXAJfltbkZODlz3kyaIRWV9Ta7erk30Ozly4HL\nMzIPAEPzrrc27/wBYJiXtwTeLnIvFkccccQRR9uOjvAlpeajKpdngQGS+gErSAsdTslrMw24ELhL\n0jBgjZmtkrS6BdlpwBmkhQ9nkPJk5eonS7qe9EhvADC7FR1zfT0DHE9anPEJOiICcBAEQdB2Kuqo\nzGy9pAuBB4EG4DYzWyjpe/79LWZ2n6QjJC0F3gPOaknWux4L3C3pHFICxxNdZoGku4EFwHrgglxu\nDknXkRzdNpJeB241s2uA24A/SFoCrAZOrqRNgiAIgrYR+aiCIAiCqqbuI1P4BuH5kuZJmuwbgPeT\n9LSkFyVNk9SriGyrm5m7Eu20xTJvM0dSa49bqx5JY9wOL0ka43VFN5rnydbauGiPLephXJzgfzcb\nJA1uQbYexkWptmjbuKjkYopqP4B+wCvA1n4+hfS+ajYwwuvOAq4pINtA2qfVD9gKmAsM6ux76gxb\n+HevAo2dfR8dZIsvAvOA7v7v3AT0B64DLvU2lwFj62BclG2LOhoX+wADgUeAwUVk62VctGqLcsZF\nvc+o/gOsA3pI2hLoQVq4MdDMHvc2DwHHFZAtZTNzV6I9tshRKwtO9gFmmdkHZrYBeJR038U2mmep\ntXHRHlvkqOVxcayZNZvZ4lZk62FclGqLHCWPi7p2VGb2b+AXwGukH+U1ZtYEzJeUG0QnALsXEC9l\nM3OXoZ22gLQU9SFJz0o6t+IKV5aXgBH+eKsHcATwOYpvNM9SU+OC9tkCantcHEmyRSnU+rhoiy2g\njeOi0svTqxpJ/YGLSdPxd4F7JJ0KnA2Ml3QVafn6RwXEa2oVSjttATDczFZK2hloktScmYl1Kcys\nWdI4YAZpJepcYENeG5NUaAzU1Lhopy2gtsfFHODjUsUrplgn0E5bQBvHRV3PqIAhwFNmttpS+KSp\nwEFmtsjMDjOzIaQp+ssFZN9g89nF7mwevqmr0R5bYGYr/fNt4K8UDgDcZTCz35nZEDM7GHgHWAys\nUopDiaTewFsFRGttXLTHFrU+LtYAi0oUrfVx0RZbtHlc1LujagaGSdpGkoBvAAvcyyNpC+BKPKJ7\nHhs3M0vqRtqQPO1T0rsSlG0LST1yqwEl9QRGkV60dlkkfdY/9wCOBSazaXM4bL7RPEutjYuybVEH\n4+IYki02a1JEtNbHRcm2KGtcdPbqkc4+gEuB+W6oSUA3YAzpfweLgJ9m2vYB7s2cH+5tlgJXdPa9\ndJYtgD1Jj4Tmkp5d14ItHnNbzAUO8bpG0oKSxaRHHjvUybgoyxZ1NC6OIb1/eh94E7i/jsdFq7Yo\nZ1zEht8gCIKgqqn3R39BEARBlROOKgiCIKhqwlEFQRAEVU04qiAIgqCqCUcVBEEQVDXhqIIgCIKq\nJhxV0GWR9GSZciMlTffyt8tNuSBpe0nnZ877SLqnnL4qiadUaGyjzBQPq5Vff6akGzpOu/Yj6XpJ\nIzpbj6ByhKMKuixmNrwD+phuZuPKFN8RuCDT1wozO6G9OlWANm2WlLQX0NPMCobL6gjkdFB3E4FL\nOqivoAoJRxV0WST91z9HSpop6R5JCyX9MdPmK5KelDRX0ixJ2+b1sXGGIOkOSb/29i9LOs7rt5X0\nkKTnPNnbUS4+Fujvyd/GSeor6SWX6S7pdm//vKSRmetNlXS/UtLBgk5S0lWSZislprslUz9T0li/\nl0WSvub1PSTdrZS0bqqkZ1QgcZ2k01x2jqSbPTRWPieTCe8j6Sy/1izgoEz9zpL+7HrOlnRQpr5J\nKaHerbkZnYcPWiRpEin6ye6SLnHZFyRd3ZKekhr832ie2/ViADNbAvRTkeSNQQ3Q2WE44oij3ANY\n658jSUEx+5Diiz1F+kHtRgqie4C325aU5G0kMN3rzgRu8PIdwBQvDwKWeLkB6OXlz2Tq+wLzMvr0\ny50D3wd+6+W9geXA1n69l4Fefr4M2K3Ave2YKf8e+JaXHwF+5uXDgSYv/wCY6OUvkHKLDfbzV0kh\njwaRHFCD108ATi9w7fszsr1d951ICf+eAMb7d5NJUbAB9gAWePlG4DIvH0aKqt3o9tkAHOjfjQJu\n8fIWwHRgRAE9bwJOBwYDMzJ6bp8pTwIO7+wxGUdljrpO8xHUFLPNbAWApLnA54G1wEozew7AzHIz\nsGJ9GB5c1cwWSsrlWNoCuNbfg3wM9PGAnC09uhoOjPe+FklaTsp8asDDZrbWdVlA+gF/I0/+65Iu\nISWwbCTFRPuHfzfVP5932dz1fuXXmy/pxbz+BBwKHAA86zbYhhSPLZ++wEovDwUeMbPVru8Uvw9I\ngYsHZezZy4OMDscTKZrZg5LeyfS93MxyqcdHAaMkzfHznsBewH4F9FxFcmR7ShoP3EuKMZhjRcYW\nQY0RjiqoFT7MlDeQxnY5gSyz+bZyv8CnkmZSg81sg6RXSSm4W6OYI8vXtWEzIak7aRZxgJm9IelH\nedf7MCOb/Rsu5Z3PJDP7YQntcn1ZXr9ik10FDDWzzXKUuXMppst7eefXmtlv8uQvLKanpC8Do4Hz\ngBOBcwroFdQY8Y4qqFWMFKm6t6QhAJJ6SWpoWawg2wFvuZM6hDTjgDRj61VE5nGSg0PSQNKjsWYK\n/4Dn1+Wc0mp/p1bKAo0nST/cSNoX+FLe9wY8DByvTalbGpVSNOSznPTID2A2cLC33SpPlxnARRtv\nQtqvgC6jSItOCvEgcLbPwpC0m+tWUE9JOwFbmtlU4CrSo8AcvUmPUYMaJGZUQVfGipRThdk6SScB\nN0jaBvgf8E1vaxm5Yv3kyn8CpvvjtGeBhd7/al94MQ+4j/TOJyczAZjoMuuBM1yf/Ot9QnczWyPp\nVtLjvjeBWSXYYAIwSdJ8kkOcT8rUnO13oaQrgRm+iGIdadXia3l9PkFKpPmcpSysVwNPk94Dzsm0\nuwi4SdILpN+SR72/HwN3Sjrd5d4kOfXtsvdqZk2SBgFP+yxsLXBaC3p+ANyeWQByeUaX/ck4zaC2\niDQfQVAD+I/3Vmb2odL+pyZgoKVszW3ta0/SApMjy9SlG7DBZ6BfBW4ys0+sQOwofMb6czM7qtXG\nQZckZlRBUBv0BP7pj+cEnF+OkwIws1ckrZXU38rbS7UHcLc7z4+Ac8vRow2cB1xX4WsEnUjMqIIg\nCIKqJhZTBEEQBFVNOKogCIKgqglHFQRBEFQ14aiCIAiCqiYcVRAEQVDVhKMKgiAIqpr/A4Fe+ip6\nq5C0AAAAAElFTkSuQmCC\n",
"text": [
"<matplotlib.figure.Figure at 0x10b25d050>"
]
},
{
"output_type": "stream",
"stream": "stdout",
"text": [
"incl = 90.0003314284 degrees\n",
"\n",
"Stellar radii in units of the star-star separation distance\n",
"(recalculated using timings from eclipse events, valid):\n",
"radius_sep_s = 0.000481306259875\n",
"radius_sep_g = 0.163880773626\n",
"\n",
"Stellar radii:\n",
"radius_s = 749987287.227 m = 1.07833020932 Rsun\n",
"radius_g = 2.55364426119e+11 m = 367.162456965 Rsun\n"
]
},
{
"output_type": "stream",
"stream": "stderr",
"text": [
"WARNING: From eclipse timing events, ratio of smaller star's radius\n",
" to greater star's radius is < 0.1. The radii ratio as calculated\n",
" from light levels may not be valid (e.g. for a binary system with\n",
" a main sequence star and a red giant).\n",
" VALID:\n",
" radii_ratio_rad = radius_s / radius_g from eclipse timings = 0.00293692938608\n",
" MAYBE INVALID:\n",
" radii_ratio_lt = radius_s / radius_g from light levels = 2.2458916679\n",
"WARNING: Inclination does not yield self-consistent solution for model.\n",
" Input parameters cannot be fit by model:\n",
" radii_ratio_lt = 2.2458916679\n",
" phase_orb_ext = 0.165111260919\n",
" phase_orb_int = 0.164135455619\n",
" incl_init = 1.4835298642\n"
]
}
],
"prompt_number": 9
},
{
"cell_type": "code",
"collapsed": false,
"input": [
"print(80*'=')\n",
"print(\"COMPARE TO PUBLISHED VALUES\")\n",
"radius_s_ref = 1.1*ast_con.R_sun.value\n",
"radius_g_ref = 369.0*ast_con.R_sun.value\n",
"tests_failed = []\n",
"for (name, ref, calc) in zip(['radius_s', 'radius_g'],\n",
" [radius_s_ref, radius_g_ref],\n",
" [radius_s, radius_g]):\n",
" if not is_calc_close_to_ref(calc=calc, ref=ref, tol=0.02*ref, name=name, verbose=True):\n",
" tests_failed.append(name)\n",
"summarize_failures(tests_failed)"
],
"language": "python",
"metadata": {},
"outputs": [
{
"output_type": "stream",
"stream": "stdout",
"text": [
"================================================================================\n",
"COMPARE TO PUBLISHED VALUES\n",
"INFO: Reference and calculated values match within tolerance.\n",
" name = radius_s\n",
" ref = 765058800.0\n",
" calc = 749987287.227\n",
" abs(ref-calc) = 15071512.7732\n",
" tol = 15301176.0\n",
"INFO: Reference and calculated values match within tolerance.\n",
" name = radius_g\n",
" ref = 2.56642452e+11\n",
" calc = 2.55364426119e+11\n",
" abs(ref-calc) = 1278025881.19\n",
" tol = 5132849040.0\n",
"INFO: Tests complete. All tests passed.\n"
]
}
],
"prompt_number": 10
}
],
"metadata": {}
}
]
}
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment