{ "cells": [ { "cell_type": "code", "execution_count": 1, "id": "107e4ef2", "metadata": { "execution": { "iopub.execute_input": "2026-08-11T01:27:06.865023Z", "iopub.status.busy": "2026-08-11T01:27:06.864851Z", "iopub.status.idle": "2026-08-11T01:27:07.290359Z", "shell.execute_reply": "2026-08-11T01:27:07.289463Z" }, "nbsphinx": "hidden" }, "outputs": [], "source": [ "%matplotlib inline" ] }, { "cell_type": "code", "execution_count": 2, "id": "78b6ded6", "metadata": { "execution": { "iopub.execute_input": "2026-08-11T01:27:07.292593Z", "iopub.status.busy": "2026-08-11T01:27:07.292333Z", "iopub.status.idle": "2026-08-11T01:27:07.300232Z", "shell.execute_reply": "2026-08-11T01:27:07.299315Z" }, "nbsphinx": "hidden" }, "outputs": [], "source": [ "%run ../tutorials/nb_setup" ] }, { "cell_type": "markdown", "id": "16be264b", "metadata": {}, "source": [ "# Defining the MilkyWayPotential model\n", "\n", "## Introduction\n", "\n", "`gala` provides simplified mass models for the Milky Way to use in orbit integration or dynamical calculations. Some of these mass models come from other publications or packages (e.g., the Law and Majewski 2010 model `LM10Potential`). Some of the potential models are defined and provided by Gala. This document describes how we determined the parameters of the Gala Milky Way models.\n", "\n", "We determine parameters of the Gala Milky Way models using compilations of enclosed mass measurements of the Milky Way and measurements of the mass structure of the Galactic disk. We then fit for the parameters of a multi-component model (e.g., disk, bulge, halo, etc.) using these measurements." ] }, { "cell_type": "code", "execution_count": 3, "id": "6d004f7f", "metadata": { "execution": { "iopub.execute_input": "2026-08-11T01:27:07.302202Z", "iopub.status.busy": "2026-08-11T01:27:07.302034Z", "iopub.status.idle": "2026-08-11T01:27:08.162613Z", "shell.execute_reply": "2026-08-11T01:27:08.161697Z" } }, "outputs": [], "source": [ "import astropy.units as u\n", "import matplotlib.pyplot as plt\n", "import numpy as np\n", "from astropy.constants import G\n", "from astropy.io import ascii\n", "from scipy.optimize import leastsq\n", "\n", "import gala.potential as gp\n", "from gala.units import galactic" ] }, { "cell_type": "markdown", "id": "361db562", "metadata": {}, "source": [ "## `MilkyWayPotential` version 1 (circa 2017)\n", "\n", "This model was previously just known as `MilkyWayPotential` in Gala, now known as \"version 1,\" and represents an older model based on measurements that are now out of date. We still describe the process of fitting for this model, for completeness.\n", "\n", "The source data for this model was compiled from published values and is included with Gala:" ] }, { "cell_type": "code", "execution_count": 4, "id": "3ad90814", "metadata": { "execution": { "iopub.execute_input": "2026-08-11T01:27:08.165165Z", "iopub.status.busy": "2026-08-11T01:27:08.164829Z", "iopub.status.idle": "2026-08-11T01:27:08.177593Z", "shell.execute_reply": "2026-08-11T01:27:08.176604Z" } }, "outputs": [ { "data": { "text/html": [ "
Table length=16\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "
rMencMenc_err_negMenc_err_posref
float64float64float64float64str27
0.0130000000.010000000.010000000.0Feldmeier et al. (2014)
0.12800000000.0200000000.0200000000.0Launhardt et al. (2002)
8.189502860861.524294994562473.7977144858963492.608627Bovy et al. (2012)
8.3110417867208.190554475949382.6968844387023236.020782McMillan (2011)
8.4102421035406.9035616733918715.62994415468328224.531876Koposov et al. (2010)
19.0208023299175.3043844317988008.3810134833267089.920685Kuepper et al. (2015)
50.0539884832748.4897519995734543.31433268490735257.4718Wilkinson & Evans (1999)
50.0529886965173.187269997867269.65930238536752776.21277Sakamoto et al. (2003)
50.0399914690706.92847109976539940.5371172696676468.2511Smith et al. (2007)
50.0419910425325.726839991469076.96850638172735113.64386Deason et al. (2012)
60.0399914690957.518869985070910.6335464344945146.92987Xue et al. (2008)
80.0689852841359.0248299936018002.3314110361048549.1029Gnedin et al. (2010)
100.01399701417307.7747899808054059.3811831336271726.1903Watkins et al. (2010)
120.0539884832260.29584199957345314.72906123854764645.32489Battaglia et al. (2005)
150.0750000000000.0250000000000.0250000000000.0Deason et al. (2012)
200.0679854974257.2006409912558030.05396313652012195.06256Bhattacherjee et al. (2014)
" ], "text/plain": [ "\n", " r Menc ... Menc_err_pos ref \n", "float64 float64 ... float64 str27 \n", "------- ------------------ ... ------------------ ---------------------------\n", " 0.01 30000000.0 ... 10000000.0 Feldmeier et al. (2014)\n", " 0.12 800000000.0 ... 200000000.0 Launhardt et al. (2002)\n", " 8.1 89502860861.52429 ... 4858963492.608627 Bovy et al. (2012)\n", " 8.3 110417867208.19055 ... 4387023236.020782 McMillan (2011)\n", " 8.4 102421035406.90356 ... 15468328224.531876 Koposov et al. (2010)\n", " 19.0 208023299175.30438 ... 34833267089.920685 Kuepper et al. (2015)\n", " 50.0 539884832748.48975 ... 268490735257.4718 Wilkinson & Evans (1999)\n", " 50.0 529886965173.18726 ... 38536752776.21277 Sakamoto et al. (2003)\n", " 50.0 399914690706.92847 ... 72696676468.2511 Smith et al. (2007)\n", " 50.0 419910425325.7268 ... 38172735113.64386 Deason et al. (2012)\n", " 60.0 399914690957.5188 ... 64344945146.92987 Xue et al. (2008)\n", " 80.0 689852841359.0248 ... 110361048549.1029 Gnedin et al. (2010)\n", " 100.0 1399701417307.7747 ... 831336271726.1903 Watkins et al. (2010)\n", " 120.0 539884832260.29584 ... 123854764645.32489 Battaglia et al. (2005)\n", " 150.0 750000000000.0 ... 250000000000.0 Deason et al. (2012)\n", " 200.0 679854974257.2006 ... 313652012195.06256 Bhattacherjee et al. (2014)" ] }, "execution_count": 4, "metadata": {}, "output_type": "execute_result" } ], "source": [ "mwdata1 = ascii.read(\"data/MW_mass_enclosed.csv\")\n", "mwdata1" ] }, { "cell_type": "markdown", "id": "53f9d71d", "metadata": {}, "source": [ "We can now plot the above data and uncertainties:" ] }, { "cell_type": "code", "execution_count": 5, "id": "deff9ff9", "metadata": { "execution": { "iopub.execute_input": "2026-08-11T01:27:08.179406Z", "iopub.status.busy": "2026-08-11T01:27:08.179178Z", "iopub.status.idle": "2026-08-11T01:27:08.906079Z", "shell.execute_reply": "2026-08-11T01:27:08.905065Z" } }, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAABJwAAAMMCAYAAAAW770RAAAAOnRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjExLjEsIGh0dHBzOi8vbWF0cGxvdGxpYi5vcmcvctoD+AAAAAlwSFlzAAAewgAAHsIBbtB1PgAAeZNJREFUeJzs3Xl4VdW9P/5PSAJEURARAYlINaDECVGccG7FeWiLVi21gKhVqx20tdr7tcPV1utca6+KgjjW6hXnuVULOIITBCpBEQmDIFAGCRDC+f3hj1NiBgjZOSchr9fz5Hn23mvtdT6H7h5P3ll77ZxUKpUKAAAAAEhIq2wXAAAAAMDmReAEAAAAQKIETgAAAAAkSuAEAAAAQKIETgAAAAAkSuAEAAAAQKIETgAAAAAkSuAEAAAAQKIETgAAAAAkSuAEAAAAQKIETgAAAAAkSuAEAAAAQKIETgAAAAAkSuAEAAAAQKIETgAAAAAkSuAEAAAAQKIETgAAAAAkKi/bBdB8rFy5MiZNmhQREdttt13k5bl8AAAAoDlbs2ZNLFiwICIi9thjj2jbtm0i40oM2GiTJk2K/v37Z7sMAAAAoBG8/fbbsd9++yUyllvqAAAAAEiUGU5stO222y69/fbbb0fXrl2zWA0AAABN0csvvxwrV66Mtm3bxje/+c1sl8MGzJ07N3030/q/9zeUwImNtv6aTV27do3u3btnsRoAAACaok6dOkV5eXkUFBT4vbGZSXKtZrfUAQAAAJAogRMAAAAAiRI4AQAAAJAogRMAAAAAiRI4AQAAAJAoT6kDAAAAElFaWhqjR4+OOXPmRLdu3aJ3795RVFSU7bLIAoETAAAA0GCjRo2K4cOHR2VlZfrYmDFjYsSIETFkyJAsVkY2uKUOAAAAaJDS0tJqYVNERGVlZQwfPjxKS0uzVBnZInACAAAAGmTkyJHVwqZ1KisrY+TIkRmuiGwTOAEAAAANMmPGjAa1s/kROAEAAAAN0rNnzwa1s/kROAEAAAANMnTo0MjNza2xLTc3N4YOHZrhisg2gRMAAADQIEVFRTFixIhqoVNubm6MGDEiioqKslQZ2ZKX7QIAAACA5m/IkCExYMCAuOKKK2LOnDnRrVu3uOaaa4RNLZTAqQn58MMP4+mnn445c+bEXnvtFcOHD6+179KlS+OFF16Ijz76KBYvXhzdu3ePY489NnbdddcMVgwAAAD/UVRUFGeffXaUl5dHQUGBsKkFEzg1Af/3f/8XP//5z2PmzJnpYyeffHKtgdPll18et9xyS6xcubLK8Z/97Gdx0UUXxS233BKtWrlbEgAAAMgOgVMTMGnSpJg5c2bsscce0a9fv7jnnnvq7P/mm29G69at48QTT4xevXpFu3bt4vXXX4+nnnoq/vznP8fOO+8cP/nJTzJSOwAAAMDXCZyagO9+97sxZMiQ6NGjR0yePHmDgdO1114be+21V7Rt27bK8Ztuuil+9rOfxUMPPSRwAgAAALLGfVdNwO677x49evTY6P77779/tbApItKPmVy6dGlitQEAAADUlxlO61m6dGk6rNl6661j66233qRxVq9eHUuXLo2OHTtmdC2ldWtA7b///hl7TQAAAICva7EznNasWRNjx46Nm266Kc4888zo1atXdOjQIQoLC6OwsDBuvPHGeo23fPny+H//7/9FUVFRtGnTJrbbbrto3bp1DBgwIO67775Gehf/sXbt2vjZz34WBQUFccUVVzT66wEAAADUpsXOcPr000/j0EMPrXKsTZs2sWrVqnqPNWPGjDjqqKNixowZERGRl5cX7du3j4ULF8b48eNj/PjxMWbMmPjb3/4WeXnJ/5OnUqm44IIL4h//+Efce++90atXr8RfAwAAAGBjtdgZTvn5+TFgwID4yU9+Evfff3989NFH8dJLL9V7nJUrV8axxx4bM2bMiNatW8ctt9wSy5Ytiy+++CLmzJkTP/jBDyIiYsyYMXHJJZck/TaisrIyzjnnnLjzzjvj9ttvj+9///uJvwYAAABAfbTYGU49evSIsWPHVjk2f/78eo9z8803x0cffZTe/tGPfpRu69q1a4wePToWLlwYzzzzTNx+++1x7rnnxl577dWw4v9/q1atirPOOisee+yxGDFiRAwbNiyRcQEAAAAaosXOcEpCKpWKv/zlLxHxVYB13nnn1djv6quvjoiv1lm6/fbbE3nt5cuXx/HHHx9jxoyJkSNHCpsAAACAJqPFznBKwsSJE2PWrFkREfHtb3+71ifS7bXXXtGrV6+YNm1ajBkzJv73f/+3Qa+7aNGiOOaYY+Ldd9+Ne++9N84666wGjbdOWVlZne1z585N5HUAAACAzZvAqQEmTJiQ3j7ggAPq7HvggQfGtGnT4vPPP4/Zs2fHDjvskG775JNP0k/FW7RoUUREfPjhh3HRRRdFRMSuu+6a3o6I+M53vhPvvPNO7LrrrvHGG2/EG2+8UeW12rRpEzfccEO9309hYWG9zwEAAAD4OoFTA/zrX/9Kb3/jG9+os+/67VOnTq0SOM2ZMyduu+22Kv1nzJiRPnbUUUdVCZzWzTT617/+VaWGdbbccstNCpwAAACguZs2bVpUVFREfn5+4k9xb8yxNzcCpwZYuHBhertTp0519t1uu+3S2+tmMa2z8847x6233lrruV+fefRf//VfsXjx4lr7t27dus5aarPu9sDazJ07N/r3779JYwMAAEAmTJs2LcrLy6OgoKBRAqfGGntzI3BqgOXLl6e3CwoK6uy7fvv650V89TS79WcwbUhSazZ9Xffu3RtlXAAAAKBl8ZS6BsjNzU1vr127ts6+67evfx4AAADA5sYMpwZo165denvFihV19l2/ff3zmrLi4uIq+xUVFVmqBAAAAGhOzHBqgC5duqS3582bV2ffdQt9R0Rsv/32jVYTAAAAQLaZ4dQAffr0SW+XlpbGwQcfXGvf0tLS9PbXZw41VSUlJVX2y8rKqi1gDgAAAPB1Zjg1wIEHHpjeHjt2bK39UqlUur13796xzTbbNHptAAAAANkicGqAoqKi2HPPPSMi4rHHHosvv/yyxn4vvvhi+pa77373uxmrDwAAACAbBE4N9Mtf/jIiIv7973/H5ZdfXq19xYoVcemll0ZExJZbbhkXXnhhRusDAAAAyLQWvYbTggULYtWqVVX211m6dGmUlZVV6d+9e/dqY5xxxhlx//33x3PPPRd//vOfY+HChXHBBRdEly5dYvLkyfHb3/42Jk+eHBER1157bXTt2rWR3g0AAABA09CiA6dBgwbFa6+9VmPbTTfdFDfddFOVYxUVFZGXV/WfLCcnJx555JE444wz4qmnnoqHHnooHnrooSp98vPz4/e//32zm9309cXNKyoqslQJAADA5mXatGlRUVER+fn50atXr2yXs9koLS2N0aNHx5w5c6Jbt27Ru3fvKCoqynZZLVKLDpy222672GGHHTa6f05OTo3Ht9xyy3jyySfj2WefjQcffDCmTJkSS5Ysic6dO8dBBx0Uw4YNq/JEOwAAAFq2adOmRXl5eRQUFDSZwKm5h2CjRo2K4cOHR2VlZfrYmDFjYsSIETFkyJAsVtYytejA6ZFHHkl0vOOOOy6OO+64RMfMppKSkir7ZWVlUVhYmKVqAAAAaExNMQTbWKWlpdXCpoiIysrKGD58eAwYMMBMpwyzaDgAAADQrI0cObJa2LROZWVljBw5MsMV0aJnOAEAAADNf+2jGTNmNKid5AmcAAAAoAXbHNY+6tmzZ4PaN0ZzD+UyzS11AAAA0EJtaO2j0tLSLFVWP0OHDo3c3Nwa23Jzc2Po0KENGn/UqFGx2267xaOPPhqvv/56PProo7HbbrvFqFGjGjTu5kzgRK2Ki4ur/Bx55JHZLgkAAIAEbS5rHxUVFcWIESOqhU65ubkxYsSIBs1E2lxCuUwTOAEAAEALtTmtfTRkyJCYOnVqfPe7342DDjoovvvd78bUqVMbfFvg5hLKZZo1nKhVSUlJlf2ysrIoLCzMUjUAAAAkLRNrH2VSUVFRnH322VFeXh4FBQWJrLG0OYVymWSGEwAAALRQjb32UX1MmzYtSkpKYtq0aRl7zY2xuYVymSJwAgAAgBaqMdc+qq9p06bFlClTmlzg1JRCueZE4AQAAAAtWGOtfbS5aEqhXHNiDScAAABo4Rpj7aPNyZAhQ2LAgAFxxRVXxJw5c6Jbt25xzTXX+Heqg8AJAAAAYAOEcvUjcKJWxcXFVfYrKiqyVAkAAADQnFjDCQAAAIBEmeFErUpKSqrsl5WVRWFhYZaqAQAA2DyUlpbG6NGj02sB9e7du8XfnuXfZPMjcAIAAIAMGTVqVJxzzjmxdu3a9LHHHnss7rrrrhb7VLhRo0bF8OHDo7KyMn1szJgxMWLEiBb7b7I5cEsdAAAAZEBpaWkMGzasStgUEbF27doYNmxYlJaWZqmy7CktLa0WNkVEVFZWxvDhw1vkv8nmQuAEAAAAGXD99ddHKpWqsS2VSsX111+f4Yqyb+TIkdXCpnUqKytj5MiRGa6IpAicAAAAIAPGjh3boPbN0YwZMxrUTtMlcAIAAACyomfPng1qp+kSOAEAAEAGHHLIIQ1q3xwNHTo0cnNza2zLzc2NoUOHZrgikiJwolbFxcVVfo488shslwQAANBsXXrppdGqVc2/hrdq1SouvfTSDFeUfUVFRTFixIhqoVNubm6MGDEiioqKslQZDZWX7QIAAACgJSgqKoq77rqr2lPZNrdwpVevXlFRURH5+fkb1X/IkCExYMCAuOKKK2LOnDnRrVu3uOaaazb536O+r0/jEDhRq5KSkir7ZWVlUVhYmKVqAAAAmr914crIkSNjxowZ0bNnzxg6dOhmEzZFfBX41FdRUVGcffbZUV5eHgUFBQ3699iU1yd5AicAAADIoKKiovjDH/6Q7TKgUQmcAAAAIMOmTZuWvu3LjBw2RwInAAAAyLBp06albx+rT+AkqKK5EDgBAABABpWWlsaoUaNi3rx50aVLl+jdu/dGr1m0qUEVZJrACQAAADJk1KhR1Z5S98QTT8SIESNiyJAhWawMktUq2wUAAABAS1BaWlotbIqIqKysjOHDh0dpaWmWKoPkCZwAAAAgA0aOHFktbFqnsrIyRo4cWef5paWlMXr06Lj55ptj9OjRAiqaNLfUAQAAQAbMmDFjk9truhVvzJgxbsWjyTLDCQAAADKgffv2m9TuVjyaIzOcqFVxcXGV/YqKiixVAgAA0HJtzK14f/jDHzJcFdTNDCcAAADIgCVLlmxSe0NuxYNsMcOJWpWUlFTZLysri8LCwixVAwAA0Lz17Nlzk9o39bz66tWrV1RUVER+fn4i49GymeEEAAAAGTB06NDIzc2tsS03NzeGDh2a6Hn11atXryguLo5evXolMh4tm8AJAAAAMqCoqChGjBhRLTzKzc2NESNGRFFRUaLnQTa5pQ4AAAAyZMiQITFgwIC44oorYs6cOdGtW7e45pprNhgabep5kC0CJwAAAMigoqKiOPvss6O8vDwKCgo2OjTa1PMgG9xSBwAAAECiBE4AAAAAJMotdQAAAEDW9erVKyoqKiI/Pz/bpZAAgRMAAACQdb169cp2CSRI4AQAAADNhFlANBcCJwAAAGgmzAKiubBoOAAAAACJEjgBAAAAkCi31FGr4uLiKvsVFRVZqgQAAABoTsxwAgAAACBRZjhRq5KSkir7ZWVlUVhYmKVqAAAANh+eNsfmTuAEAAAAGeZpc2zu3FIHAAAAQKIETgAAAAAkyi11AAAAABvB2lsbT+AEAAAAsBGsvbXx3FIHAAAAQKIETgAAAAAkSuAEAAAAQKIETgAAAAAkyqLhAAAANNi0adPST++ysDIgcAIAAKDBpk2bFuXl5VFQUCBwAtxSBwAAAECyBE4AAAAAJMotdQAAADRIaWlpjB49OubMmRPdunWL3r17R1FRUbbLArJI4AQAAMAmGzVqVAwfPjwqKyvTx8aMGRMjRoyIIUOGZLEyIJvcUgcAAMAmKS0trRY2RURUVlbG8OHDo7S0NEuVAdkmcAIAAGCTjBw5slrYtE5lZWWMHDkywxUBTYXACQAAgE0yY8aMBrUDmy9rOFGr4uLiKvsVFRVZqgQAAGiKevbs2aB2YPNlhhMAAACbZOjQoZGbm1tjW25ubgwdOjTDFQFNhRlO1KqkpKTKfllZWRQWFmapGgAAoKkpKiqKESNGVFs4PDc3N0aMGBFFRUVZrA7IJoETAAAAm2zIkCExYMCAuOKKK2LOnDnRrVu3uOaaa4RN0MIJnAAAAGiQoqKiOPvss6O8vDwKCgqETYA1nAAAAABIlsAJAAAAgEQJnAAAAABIlMAJAAAAgEQJnAAAAABIlMAJAAAAgETlZbsAAAAAmr9evXpFRUVF5OfnZ7sUoAkQOAEAANBgvXr1ynYJQBPiljoAAAAAEiVwAgAAACBRAicAAAAAEiVwAgAAACBRAicAAAAAEiVwAgAAACBRAicAAAAAEiVwAgAAACBRAicAAAAAEiVwAgAAACBRAicAAAAAEiVwagIqKytj3Lhx8Ytf/CL69OkTHTp0iLPOOqvOc2bMmBF/+tOf4lvf+lZ07NgxOnToEF9++WWGKgYAAACoXV62CyDi6quvjquuuqrKsbrCo08++SR23nnnasdTqVTitQEAAADUlxlOTUCrVq1iwIABce2118ajjz66wf5r166Nnj17xo9//ON48cUXo7i4OANVAgAAAGwcM5yagCuvvDJ+/etfR0TE5MmTN9h/5513jk8++SS9n5fnf0YAAACg6ZBUNAE5OTmN2h8AAGgcpaWlMXLkyPjwww+ja9eucdppp8XRRx+d7bIAsk7gFBGrVq2KDz74ICZOnBhLliyJiIgBAwbEgAED6jVOKpWKt99+O6ZMmRJLliyJzp07x4EHHhg9e/ZsjLIBAIAsGjVqVAwfPjwqKyurHLvrrrtiyJAhWawMIPtabOC0ePHiuPzyy2PChAkxadKkqKioqNJ+1VVX1Stwev755+PCCy+scqvbOkcffXTcfvvtgicAANhMlJaWVgubIr5ab3X48OExYMCAKCoqylJ1ANnXYhcNX7hwYdx5553x7rvvRqtWrWK//faLE044YZPGuv/+++O4446LTz75JLbYYos4/fTT45JLLomjjjoqcnJy4sUXX4z9998/pk+fnvC7AAAAsmHkyJHVwqZ1KisrY+TIkRmuCKBpabEznDp27Bi333577LvvvrHnnntGfn5+jBs3Lp5++ul6jTN9+vQYNmxYpFKp6N27d7z00ktRWFiYbn/hhRfi5JNPjgULFsSgQYNi4sSJ0apVi835AABgszBjxowGtQNs7lps8tGxY8c477zzol+/fpGfn7/J41x11VWxevXqiIh48MEHq4RNEREDBw6MK6+8MiIi3n///Xj44Yc3vWgAAKBJ2NByGcuWLctQJQBNU4sNnJLw5ZdfxpgxYyIi4vDDD4999tmnxn4XXnhhtG7dOiK+uv0OAABo3oYOHVrnnQsvvPBClJaWZrAigKZF4NQAr732WpSXl0fEVzOZatOxY8fo379/RES8/PLL1RYobyrKysrq/Jk7d262SwQAgCahqKgojjnmmFrbreMEtHQtdg2nJEyaNCm9vddee9XZt2/fvjFu3LhYvXp1lJaWRp8+fdJtb775Zvo/VmvXro2IiGeeeSY6dOgQEV/Nnnr88cerjLfnnnvGZ599FhH/ma7bvXv3dPu0adOic+fO9Xo/X78dEAAAqN1WW21VZ7t1nICWTODUAJ988kl6u0ePHnX23XHHHdPbH3/8cZXAac2aNbFkyZIq/dc/tnz58mrjLV26tNo56++vC64AAIDGsaF1nDbUDrA5Ezg1wNKlS9PbW2+9dZ1912//+gKCBx54YCxevLjWc2ta1PzDDz+sM1Rq3759nfXUZNasWXW2z507N31rIAAAtHRDhw6N6667LiorK6u15ebmxtChQ7NQFUDTIHBqgBUrVqS31y0KXpu2bdvWeF7EV/8xWnf73MbaUMC1Kda/JQ8AAKhdaWlpjBw5MvbZZ5+YMGFCpFKpdFtubm6MGDEiioqKslghQHYJnBqgoKAgvb169eo6+65cubLG8wAAgOZl1KhRMXz48Gozm3bZZZfYe++945prrhE2AS2ewKkB6rpN7uvWb9/Q4oJNRXFxcZX9pvp0PQAAyJTS0tIaw6aIr9Z4veyyy4RNABHRKtsFNGc77bRTenvdE+Nqs367xQMBAKB5GjlyZI1hU8RXD+556aWXMlwRQNNkhlMD7LHHHuntSZMmxcCBA2vt++GHH0ZERF5eXvTu3bvRa0tCSUlJlf2ysrIoLCzMUjUAAJB9M2bMqLP9888/z1AlAE2bGU4NcNhhh0WbNm0iIuLFF1+std/SpUvjzTffjIiII444YoMLjAMAAE3Thu5W2H777TNUCUDTJnBqgK233jpOOOGEiIj4+9//HlOmTKmx35133pleNPzMM8/MWH0AAECyhg4dGrm5uTW25ebmxvDhwzNcEUDTJHBqoN/97neRm5sba9eujbPOOisWLFhQpX38+PFx1VVXRUTEbrvtFoMHD85GmQAAQAKKiopixIgR1UKn3NzcGDFiRBx99NFZqgygaWnRazg9+OCDVRbznjlzZnp73Lhx8cc//rFK/1/+8peRk5NT5VifPn3iT3/6U1x44YXx/vvvx6677hqDBg2KLl26xOTJk+OJJ56INWvWRPv27ePhhx+u9a8hAABA8zBkyJAYMGBAjBw5MmbMmBE9e/aMoUOHejodwHpyUqlUKttFZMvhhx8er7322kb3r6ioiLy8mjO6hx9+OC6++OKYP39+tbb+/fvHyJEjo7i4eJNrzYav11tRURGlpaURETFr1qzo3r17NsoCAAAAErL+A8KS/F2/Rc9wOvPMM+OAAw7Y6P6tWtV+B+Lpp58e3/72t9NrOS1ZsiQ6d+4cBx10UPTt2zeJcgEAAACahRY9w4n6aazUEwAAAMiOxvpd36LhAAAAACRK4AQAAABAogROAAAAACSqRS8aTt1qekodAAAAwIaY4QQAAABAosxwolYlJSVV9tdfuR4AAACgNmY4AQAAAJAogRMAAAAAiRI4AQAAAJAogRMAAAAAiRI4AQAAAJAoT6mjVsXFxVX2KyoqslQJAAAA0JyY4QQAAABAosxwolYlJSVV9svKyqKwsDBL1QAAAADNhRlOAAAAACRK4AQAAABAogROAAAAACRK4AQAAABAogROAAAAACRK4AQAAABAogROAAAAACQqL9sF0HQVFxdX2a+oqMhSJQAAAEBzYoYTAAAAAIkyw4lalZSUVNkvKyuLwsLCLFUDAAAANBdmOAEAAACQKIETAAAAAIkSOAEAAACQKIETAAAAAIkSOAEAAACQKIETAAAAAIkSOAEAAACQKIETAAAAAIkSOAEAAACQqLxsF0DTVVxcXGW/oqIiS5UAAAAAzYkZTgAAAAAkygwnalVSUlJlv6ysLAoLC7NUDQAAANBcmOEEAAAAQKIETgAAAAAkSuAEAAAAQKIETgAAAAAkSuAEAAAAQKIETgAAAAAkSuAEAAAAQKIETgAAAAAkSuAEAAAAQKIETgAAAAAkSuAEAAAAQKIETgAAAAAkKi/bBdB0FRcXV9mvqKjIUiUAAABAc2KGEwAAAACJMsOJWpWUlFTZLysri8LCwixVAwAAADQXZjgBAAAAkCiBEwAAAACJEjgBAAAAkCiBEwAAAACJEjgBAAAAkCiBEwAAAACJEjgBAAAAkCiBEwAAAACJEjgBAAAAkCiBEwAAAACJEjgBAAAAkCiBEwAAAACJEjgBAAAAkCiBEwAAAACJEjgBAAAAkCiBEwAAAACJEjgBAAAAkCiBEwAAAACJyst2ATRdxcXFVfYrKiqyVAkAAADQnJjhBAAAAECizHCiViUlJVX2y8rKorCwMEvVAAAAAM2FGU4AAAAAJErgBAAAAECiBE4AAAAAJErgBAAAAECiBE4AAAAAJCrxp9Sdd955sXDhwqSHreLOO++Mjh07NuprAAAAALBpEg+cnnnmmZg9e3bSw1Zx8803C5wAAAAAmqjEA6d1dtppp8jNzU10zE8//TQqKysTHRMAAACAZDVa4DR27Njo3r17omN279690WdPAQAAANAwFg0HAAAAIFECJwAAAAASlfgtdW+++WasWbMmunbtmvTQjTo2AAAAAMlIPHBKet2mTI0NAAAAQDLcUgcAAABAogROAAAAACRK4AQAAABAogROAAAAACRK4AQAAABAogROAAAAACRK4AQAAABAovKyXcDXpVKpeOedd+Ltt9+OZcuWRbdu3eKb3/xm7LDDDtkuDQAAAICN0KQCp4kTJ8bw4cPjvffei7Zt20a7du1i0aJFkZOTE0OHDo2bbropttxyy2yX2SgWLlwYzzzzTDz99NMxZ86cOPTQQ+Oaa66p85wpU6bEHXfcEVOnTo2CgoI44ogj4vzzz4+2bdtmqGoAAACA6ppM4PTKK6/E8ccfH8XFxfHyyy/H4YcfHrm5ubF48eK4995747/+67/ivffei1dffXWzC51uu+22uOSSS6KysjJ9rFOnTnWe8/jjj8fpp58eq1evTh978sknY9SoUfHaa69Fhw4dGqtcAAAAgDo1iTWcFi1aFIMGDYpDDz00xo0bF0cddVTk5uZGRMQ222wTl1xySYwdOzamTZsWl1xySZarTd6CBQuiQ4cOMXjw4PjDH/6wwf4LFy6MH/7wh7F69eoYMmRIvPDCC/Hggw9Gnz594sMPP4xLL700A1UDAAAA1KxJBE5/+tOfYvXq1XHfffdFmzZtauyz1157xbXXXhujRo2KmTNnZrjCxnXRRRfF559/Hvfee2+ccMIJG+x/zz33xJIlS+KMM86IkSNHxtFHHx1nnHFG/P3vf4+tttoq7r333li8eHEGKgcAAACorkkETs8//3yceOKJsd1226WPvfXWW/Hqq6/Gp59+mj42ePDgyMvLi+effz4LVTaeTp06pWd0bYyXXnopIiIuvPDCKse7dOkSgwYNioqKivjHP/6RaI0AAAAAG6tJBE5lZWXRq1evKsfOOuusOOKII+Kee+5JH9tyyy2jW7duMWvWrMRrWLZsWYwdOzaefvrpePrpp2PatGmbNM6///3veOedd+Lll1+ODz/8MCoqKhKuNGLq1KmRk5MT/fr1q9bWv3//dB8AAACAbGgSi4ZvscUWsWzZsirHjjzyyNh1112rBVHLli2Ldu3aNfg1v/zyyxg5cmRMmDAh3nnnnfjoo49i7dq16farrroqfvOb32z0eKWlpfHzn/88nn/++Soh0zbbbBNDhgyJ3/3ud4ktdr5w4cLo0KFDjU+j69KlS0R8tS4WAAAAQDY0icCpX79+8eKLL1Y5duedd1brN3HixFi4cGHst99+DX7NuXPnxsUXX5zeb9euXey4444xZcqUeo81bty4GDhwYKxYsSIiIoqLi6NLly4xderUmDNnTtx4443xwgsvxNixY2ObbbZpcO1r1qyJLbbYosa2vLyv/idtjJlVAAAAABujSdxSd+6558akSZNi1KhRtfaprKyMyy67LIqKiuKwww5r8Gu2a9cuLrroorjnnnti8uTJsWTJkrjjjjvqPc4XX3wRJ598cqxYsSI6deoUr732WkyePDlefvnl+Oyzz+KGG26InJycKCkpibPOOqvBda+r/eszwtZZunRpug8AAABANjSJwOmII46IYcOGxY9+9KO47777qrUvW7YszjzzzBg7dmzceeed6Vk8DdGlS5e49dZb4+yzz47i4uJo1WrT/il+//vfp29fGz16dBx66KHpttzc3PjZz34W559/fkREPPfcc/Hyyy83uPYdd9wxVq5cWeNaVuvWnurRo0eDXwcAAABgUzSJwCki4n//93/jrLPOih/84Aex3377xZVXXhn/8z//E+eff3584xvfiGeffTYeffTROPzww7NdalpFRUU6INtrr73iuOOOq7Hf5Zdfng607r777ga/7rqFwZ966qlqbU8++WSVPgAAAACZ1iTWcIqIyM/Pj7vvvjvOPPPMuP322+OBBx6IZcuWRbdu3eLss8+On/70p7HDDjtku8wqxo8fH4sXL46IiBNPPLHWfjvuuGPstdde8d5778UzzzwTa9eu3eQZVRERZ5xxRtxxxx3x29/+No444ojYbbfdIiLiuuuui3fffTd23XXX2Geffeo9bllZWZ3tc+fO3aR6AQAAgJalyQRO6xx11FFx1FFHZbuMjfL++++ntze0kHn//v3jvffei2XLlsWMGTNi5513TrdNmjQpfvSjH0VEpBceHzt2bAwYMCA99k033ZTuf9hhh8Wpp54aY8aMiT322CP23HPPWLRoUcycOTNatWpVpW99FBYWbtJ5AAAAAOtrcoFTc1JaWpre3tCaSTvttFN6e9q0aVUCpyVLlsT48eOr9F+0aFH6WNu2bauNd99998Ull1wS9957b7z33nsREdG1a9e45ZZb4phjjqn3ewEAAABIisCpAf7973+nt7fZZps6+3bo0CG9vWTJkipte+65Z4wdO3ajzl1nyy23jLvuuiuuv/76+Pjjj6OgoCB69+4dubm5G1V7TWpahHx9c+fOtTYUAAAAsEECpwb48ssv09s1zUJa3xZbbJHeXr58eZW2rbfeOn37XH116NAh+vXrt0nnfl337t0TGQcAAABo2ZrMU+qao9atW6e316xZU2ff1atXp7fbtGnTaDUBAAAAZJsZTg2w1VZbpbe/Pmvp69afDbX+eU1ZcXFxlf2KioosVQIAAAA0J2Y4NcD6t6DNnj27zr5lZWU1ngcAAACwuTHDqQHWnwE0derUOOKII2rtO3Xq1IiIyMnJid12263Ra0tCSUlJlf2ysrIoLCzMUjUAAABAc9EogdPgwYNjwYIFiY13//33R6dOnRIbLymHHnpotGrVKtauXRuvvPJKXHDBBTX2W716dYwbNy4iIvbbb7/YcsstM1kmAAAAQEY1SuD0yiuvbPAWs/pYuXJlYmMlqUuXLnHooYfGq6++Gk899VTMmzcvunTpUq3fww8/HEuWLImIiNNPPz3TZQIAAABkVKOu4ZSbmxtt2rRp8E9OTk5jltkgv/nNbyIiYtWqVTFs2LAqT6OLiJg5c2ZcdtllEfFVQHXuuedmukQAAACAjGrUNZzy8vLi1FNPjWHDhsVRRx3V5IKjN954IxYuXJjenzJlSnp72rRp8fTTT1fpf/zxx1d7D4cddlhccsklccstt8Szzz4b/fv3j3POOSe6dOkSkydPjr/85S+xYMGCyMvLi5EjR0a7du0a900lyFPqAAAAgE2Rk0qlUkkP+pvf/Cb+93//N+bPn58+ttNOO8XQoUPjhz/8YZNZePrwww+P1157baP7V1RURF5e9Yxu7dq1cdVVV8Uf//jHWLNmTbX27bffPu6666444YQTGlRvptUUOJWWlkZExKxZszxtDwAAAJq59R8QluTv+o0SOEV8FU489dRTcffdd8cLL7wQlZWVERHRqlWrOProo+Occ86Jk046KfLz8xvj5TfKr371q5g0adJG93/yySejVava70KcOXNmPPLIIzFlypRYsmRJdO7cOQ466KA45ZRTYquttkqi5KxqrIsQAAAAyI5mFzitb/bs2TFq1KgYOXJkzJgxI318u+22i8GDB8ewYcOiT58+jV0GDSRwAgAAgM1LY/2u36iLhq+zww47xK9//ev4+OOP4+WXX44zzjgj2rZtGwsWLIgbb7wxiouL46CDDoq77747li9fnomSAAAAAGgkGQmc1snJyYmjjjoqHnzwwZgzZ0786U9/ir322isivlrA+5xzzomuXbvGggULMlkWAAAAAAnKaOC0vm222SZ+/OMfx/vvvx/jx4+PnXfeOSIili9fHqtWrcpWWQAAAAA0UPVHrmXQpEmT4u677477778/Fi5cmM1SqEFNT6kDAAAA2JCMB05Lly6Nv/71r3H33XfH22+/nT5eUFAQ3/72t2PYsGGxww47ZLosAAAAABKSscBp7Nixcffdd8cjjzwSK1asSB/v27dvnHPOOXHmmWdGhw4dMlUOG6GkpKTK/vor1wMAAADUplEDp3nz5sXo0aPj7rvvjtLS0vTxbbbZJs4888wYNmxY9O3btzFLAAAAACDDGiVwevbZZ+OOO+6IZ599NtasWRMRXz2h7vDDD49hw4bFd77znWjbtm1jvDQAAAAAWdYogdO5554bs2fPjoiIHXbYIc4+++wYNmxYfOMb32iMlwMAAACgCWnUW+q22GKL2G233WLixIkxceLETR7n/vvvj06dOiVYGQAAAACNpVEDpxUrVsTLL7/c4HFWrlyZQDUAAAAAZEKjBE5bbrllbLnllomN16pVq8TGYuMVFxdX2a+oqMhSJQAAAEBz0iiB00cffdQYwwIAAADQDDTqLXU0byUlJVX2y8rKorCwMEvVAAAAAM2Fe9UAAAAASJTACQAAAIBECZwAAAAASFTigdOFF14Y3/ve92LRokVJD92oYwMAAACQjMQDpyeeeCIefvjhWLFiRdJDN+rYAAAAACTDLXUAAAAAJCqvsQZ+5JFHYptttkl0TDObAAAAAJq+RgucfvaznzXW0GRIcXFxlf2KioosVQIAAAA0J26pAwAAACBRic9wGjduXKxZsybpYavo2rVro47PV0pKSqrsl5WVRWFhYZaqAQAAAJqLxAOnnXbaKekhAQAAAGhG3FIHAAAAQKIETgAAAAAkSuAEAAAAQKIETgAAAAAkSuAEAAAAQKIETgAAAAAkSuAEAAAAQKIETgAAAAAkKi/bBdB0FRcXV9mvqKjIUiUAAABAc2KGEwAAAACJMsOJWpWUlFTZLysri8LCwixVAwAAADQXWQmcysvLY+bMmfH5559HeXl5tGnTJjp27Bg77bRTtG/fPhslAQAAAJCQjAROa9asiRdffDGeffbZeO2112LKlCmxdu3aGvvutNNOcfDBB8exxx4bJ510Umy11VaZKBEAAACAhDRq4FRWVha33nprjBo1KhYsWLBR53z66afx6aefxgMPPBAFBQUxaNCg+MlPfhJ9+/ZtzFIBAAAASEijBE7z5s2L3//+9zFixIgqTzbLz8+PvfbaK/r16xfbbbddbLvtttGhQ4dYtmxZLFq0KL744ov48MMPY+LEifHll19GeXl53HvvvXHvvffGCSecEFdffXXsueeejVEyAAAAAAlJPHC69dZb48orr4xly5ZFRMR2220X3/ve92LQoEGx3377Rdu2bTc4RmVlZUyaNCmeeOKJuP/++2P69Onx9NNPx3PPPRfnn39+3HDDDdGmTZukSwcAAAAgATmpVCqV5IDdu3eP2bNnx4ABA+KXv/xlHHPMMZGX17Bc66233oqbb745/va3v8XatWtj1qxZ0b1794QqZmOt/5Q6/xsAAABA89dYv+u3SmSU9Rx00EHx4osvxtixY+OEE05ocNgUEbH//vvHQw89FJMmTYrTTjstWrVKvGwAAAAAEpL4LXV/+9vfajz+3nvvRXFxcbRu3XqTx+7Tp088/PDDm3w+AAAAAI0vY1OFTjzxxOjQoUN8/vnnmXpJAAAAALIgY4FT69ato7y8vMpT6wAAAADY/GQscNp+++0jIqK8vDxTLwkAAABAFiS+hlNt9t9//3jzzTfjww8/jKKioky9LA1QXFxcZd/sNAAAAGBjZGyG0/e+972IiPjLX/6SqZcEAAAAIAsyNsPpgAMOiJNOOimefPLJuPrqq+PKK6/M1EuziUpKSqrsl5WVRWFhYZaqAQAAAJqLjM1wGjRoUHTp0iVyc3Pj17/+dZx++umxaNGiTL08AAAAABmSscDpjTfeiDvvvDMqKysjIuJvf/tbdO3aNU499dQYPXp0TJkyJdauXZupcgAAAABoJBm7pe6ggw6Kl19+ORYvXpw+tnr16nj88cfj8ccfj4iILbfcMvbaa6/o169f9OvXL/bZZ5/o06dP5ObmZqpMAAAAABooY4HT3/72t4iI+Pjjj2PixInpn3fffTcdQn355Zfx+uuvx+uvv54+r6CgID755JPo0qVLpkoFAAAAoAEyFjits/POO8fOO+8cp512WvpYXSFUeXl5rFmzJtNlAgAAALCJMh441WRDIdQWW2yRxeoAAAAAqI8mETjVpKYQCgAAAICmL2NPqQMAAACgZRA4AQAAAJAogRMAAAAAiRI4AQAAAJAogRMAAAAAiRI4AQAAAJAogRMAAAAAicrLxIuUlpbGqlWrYptttokddthhg/0XL14cs2fPjoiI3XbbLXJzcxu7RAAAAAASkpEZTldddVXsscceMXDgwI3qf+mll8Yee+wRp5xyirAJAAAAoJnJSOA0ePDgiIgoKSmJ999/v86+q1ativ/7v/+LiIjvf//7jV0aAAAAAAnLSOB09NFHx/bbbx8REQ888ECdfZ9++ulYsmRJRAicAAAAAJqjjAROubm5ccYZZ0RExEMPPRRr166tte+DDz4YEREHHnhg7LLLLpkoDwAAAIAEZWTR8Iivbqu7+eabY/bs2fHqq6/GkUceWa3PkiVL4plnnkn3J7uKi4ur7FdUVGSpEgAAAKA5ycgMp4iIffbZJ/r06RMRtd9W9+ijj8aqVauidevWcdppp2WqNAAAAAASlLHAKeI/azL93//9X6xcubJa+7rb6Y477rjYdtttM1kaNSgpKany849//CPbJQEAAADNQEYDp7POOitycnKq3Dq3zpw5c+LVV1+NCIuFAwAAADRnGQ2cdtxxxzjssMMiIuL++++v0rZuMfEOHTrECSeckMmyAAAAAEhQRgOniP8sBv7ss8/G4sWL08fXret02mmnRZs2bTJdFgAAAAAJyXjg9N3vfjcKCgpi9erV8eijj0ZExL/+9a947733IsLT6QAAAACau4wHTltvvXWcdNJJEfGf2+rWzW7q2bNnHHzwwZkuCQAAAIAEZTxwivjPLKaxY8fGZ599ln463fe///3IycnJRkkAAAAAJCQrgdPAgQNju+22i1QqFZdcckl88sknEeHpdAAAAACbg6wETnl5efG9730vIiIef/zxiIjo379/9OrVKxvlAAAAAJCgrAROEdUXB7dYOAAAAMDmIWuB03777Re9e/eOiIj8/Pz0jCcAAAAAmre8bL74vffeG9OmTYuOHTtGp06dslkKAAAAAAnJauDUv3//6N+/fzZLAAAAACBhWbulDgAAAIDNk8AJAAAAgEQJnAAAAABIlMAJAAAAgEQJnAAAAABIlMAJAAAAgEQJnAAAAABIlMAJAAAAgEQlHjhVVlYmPWQ1a9eubfTXAAAAAGDTJB44HXLIIXHHHXdERUVF0kPHM888E/vuu2/MmTMn8bGbo+nTp8dFF10U++67b+yxxx7xne98J55//vlslwUAAAC0cIkHTp999lmcf/750atXr7j11lvjiy++aNB45eXl8dBDD8X+++8fJ5xwQkycODGhSpu3V155Jfbee++47bbbYuLEiTF58uR47LHH4thjj43f//732S4PAAAAaMESD5xuu+226NGjR3z66adx8cUXR7du3eLkk0+O++67Lz7++OONGmPu3Lnx+OOPx5AhQ2L77bePM888M95+++3Yeuut48Ybb4yuXbsmXXazsnr16vjBD34QX375ZQwaNCheeeWVeO+99+LWW2+NbbbZJq666qqYMGFCtssEAAAAWqicVCqVSnrQ8vLy+POf/xz/8z//U22G07bbbht9+/aNzp07R8eOHaN9+/axfPnyWLRoUSxcuDAmT54cn332WZVz2rRpE+edd15ceeWV0blz56TLbXZeeeWVOPLII+PAAw+M8ePHR05OTrrtr3/9a5xxxhlx9tlnxz333JPo65aVlUVhYWFERMyaNSu6d++e6PgAAABAZjXW7/p5iYzyNQUFBXHZZZfFhRdeGPfcc0/ccccd8eGHH0ZExMKFC+Pll1/eqHF23HHHGDJkSFxwwQWCpvXMmjUrIiIOPfTQKmFTRMThhx8eEREvvfRSpssCAAAAiIhGCpzW2WKLLeKCCy6ICy64IKZMmRJPP/10vPbaa/Huu+/GvHnzqvXv0KFD7LXXXnHwwQfHscceGwcffHC1QKUp+eSTT2LixImxZMmSiIjYZ599Yp999qn3ONOnT4+pU6fGkiVLonPnzrHvvvtGx44da+2/1VZbRUTEjBkzaqwpImLOnDmxdOnS2HrrretdDwAAAEBDNGrgtL4+ffpEnz594he/+EVERKxYsSLmz58f5eXl0bp169h2222jQ4cOmSpnkyxdujSuvfbamDBhQkyYMCEWLVpUpf2qq66qV+D09ttvx0UXXRTvvPNOleP5+fkxaNCguPnmm2O77bardt7BBx8c+fn58eijj8bdd98dQ4cOjZycnPRT69ZZvHixwAkAAADIuMQXDd9YW2yxRey0006x2267xc4779zkw6aIiPnz58c111wTL774YixatCh69uwZhxxyyCaN9dRTT8VBBx0U77zzTuTm5sZRRx0VZ511Vuy1115RUVERDz74YOy7774xe/bsaud27tw5Lr300li7dm2cc8450a5du+jSpUsUFRXF559/HnvssUdERKxdu7ZB7xcAAABgU2QtcGqO2rdvH//93/8dzz//fHzxxRfxySefxDXXXFPvccrKyuJ73/teVFZWxo477hgffvhhvPzyy3H//ffH+++/Hw888EDk5ubGZ599FqeffnqNY1x99dVx4403Rrdu3WLFihXx+eefR79+/eLvf/97euZVXbflAQAAADSWRgmcTj311Ljoooti1KhR8eGHH8aaNWsa42Uybrvttosrr7wyBg4cGNtuu+0mj/Ob3/wmVqxYERERDzzwQPTp06dK+5lnnhmXXnppRESMHz8+Hn/88Wpj5OTkxE9/+tOYPXt2zJs3L+bPnx8TJkyIVq1axezZs2OnnXaK9u3bb3KNAAAAAJuqUdZweuedd6qEJG3bto299tor+vXrl/4pLi6OvLyMLSHVZKxcuTIefvjhiIg44IADYsCAATX2++lPfxo33HBDrFmzJkaPHh2nnHJKrWNuv/326e3rrrsuIiJOPvnk5IoGAAAAqIeMJD4rV66Mt956K9566630sbZt28aee+5ZLYTKz8/PRElZ889//jOWL18eERHHHXdcrf2233776NevX7z11lvxwgsvRGVlZeTm5qbb77vvvmjbtm2cfPLJ0bp161i2bFn88Y9/jLvuuivatWsXP/3pT+tdW1lZWZ3tc+fOrfeYAAAAQMvTqIHTFltsESeddFL8+9//jnfffTfmz5+fblu5cmW8/fbb8fbbb6ePtWnTploItfvuu29WIdSHH36Y3t7QE+3WBU7l5eUxffr06N27d7rt448/jt/+9reRl5cXXbp0ic8//zwqKioiLy8vRo4cGT169Kh3bYWFhfU+BwAAAODrGmUNp5NPPjlycnJixYoVMWbMmOjXr1/MmDEjZs2aFY8//nj8v//3/+L444+Prl27Vjlv1apV8c4778Ttt98ew4cPj3322Se22mqrmDdvXmOUmRXTp09Pb++444519l0/NFr/vIiIH/zgB3HmmWdGXl5elJWVRWVlZRx++OHx2muvxaBBg5ItGoAWY9q0aVFSUhLTpk3LdikAADRjjTLD6bbbbovBgwfHRRddFBMnToyrr7467r333rjhhhti0KBBVdYXmjdvXkycODH98+6771a5tWvVqlWbzaLjERFLly5Nb3fo0KHOvusv+r3+eRER3/jGN+KBBx6IysrKmDt3bmy77bZRUFDQoNpmzZpVZ/vcuXOjf//+DXoNAJq2adOmRXl5eRQUFESvXr2yXQ4AAM1Uo91Sd8ABB8Tbb78dI0aMiCuvvDJmzZoVp512Whx55JHxpz/9KYqLiyMiokuXLnH88cfH8ccfnz53wYIFVQKo1q1bN1aZGffll1+mt9u0aVNn3/UDpPXPW19ubm507949kdqSGgcAAABo2Rrllrr04K1axXnnnRfTpk2LH/3oR9GqVav4xz/+EXvvvXf89Kc/jSVLltR43nbbbRfHHHNMXHnllfF///d/0blz58YsM6Patm2b3q6oqKiz76pVq2o8DwAAAKApa9TAaZ2OHTvGX/7yl5g4cWIcfPDBsWbNmrj55pujV69eMWrUqEilUpkoo0nYaqut0tvrnlZXm/Xb1z8vU4qLi6v8HHnkkRmvAQAAAGh+MhI4rbP33nvHuHHj4t57740uXbrE/PnzY+jQoXHggQfGhAkTMllK1qy/EPiG1kxav31TnjoHAAAAkA0ZDZzWGTx4cEybNi1+/vOfR35+frz11lux//77x/Dhw+OLL77IRkkZs/vuu6e3J0+eXGffde25ubmx6667NmpdNSkpKany849//CPjNQAAAADNT1YCp4ivbhG7/vrr44MPPohvfvObsXbt2rjrrruiV69esXDhwmyV1egOO+ywyMv7aq32v//977X2W7FiRbz++usREXHwwQdbwwmARldaWhqjR4+Om2++OUaPHh2lpaXZLgkAgGYqa4FTRUVFTJo0Kd57773o169fdOrUKSIiFi9eHOXl5dkqq9F17Ngxjj766IiIeO6552LGjBk19hs9enT6yXRnnHFGxuoDoGUaNWpU7LbbbvHoo4/G66+/Ho8++mjstttuMWrUqGyXBgBAM5SXiReZPXt2fPjhhzFp0qT48MMP48MPP4x//etfG3xK2+bqd7/7XTz33HNRWVkZgwcPjueee67KouCTJk2KK664IiIievbsGUOGDMlWqQC0AKWlpTF8+PCorKyscryysjKGDx8eAwYMiKKioixVBwBAc9QogdP9998fb731VkyaNCkmTZoUixYtqrVv586do2/fvtG3b9/Ye++9o2/fvrHDDjs0RlmJeOaZZ2Lu3Lnp/enTp6e333333bjrrruq9B82bFjk5ORUOdavX7/4/e9/H7/+9a9j/Pjxsfvuu8fgwYOjS5cuMXny5Lj//vvjyy+/jIKCgnjggQeiTZs2jfumalFcXFxlv6UGhACbu5EjR1YLm9aprKyMkSNHxh/+8IcMVwUAQHOWk0qlUkkP2r1795g9e3bVF8rJiZ49e1YLl7p165b0yzeqww8/PF577bWN7l9RUZFes+nrbrvttrj88stj+fLl1dqKiopi9OjRceCBB25yrQ1VU+C0bj2PWbNmRffu3bNRFgAJ+973vhcPP/xwre2nn356/PWvf81gRQAAZEpZWVkUFhZGRLK/6zfqLXXt2rWL008/Pc4444zYd999o3379o35chlx/PHHxy677LLR/Vu1qn2ZrAsvvDAGDx4cTz75ZEyZMiWWLFkSnTt3joMOOiiOPPLIyM3NTaLkTVZSUlJlf/2LEIDNR8+ePRvUDgAAX5eRGU4FBQWx5557xj777BP9+vWLffbZJ3bffffIz89P+qVpRI2VegKQXaWlpbHbbrvVeFtdbm5uTJ061RpOAACbqWY1w2n77bevEjiVl5fHW2+9FW+99Vb6WJs2bWL33XdPB1D77LNP7LnnnllbrwgAWqqioqIYMWJEtYXDc3NzY8SIEcImAADqrVFmOEVEfP755zFx4sQqP2VlZXWek5+fH3369EmHUP369Yu99torCgoKGqNE6skMJ4DNW2lpaVxxxRUxZ86c6NatW1xzzTXCJgCAzVxj/a7faIFTTebPn18thJo1a1ad5+Tm5sasWbOia9euGaqS2gicADZ/Tz/9dJSXl0dBQUGccMIJ2S4HAIBG1qxuqatN586d49hjj41jjz02fWzBggXVQqjPPvss3V5ZWVnro5ppXDU9pQ4AAABgQzIaONVku+22i2OOOSaOOeaY9LEvvviiSgDVunXrLFYIAAAAQH1kPXCqSadOnWLgwIExcODAbJfSopWUlFTZX3+aHQAAAEBtWmW7AAAAAAA2LwInAAAAABIlcAIAAAAgUQInAAAAABIlcAIAAAAgUU3yKXU0DcXFxVX2KyoqslQJAAAA0JwInACAtF69ekVFRUXk5+dnuxQAAJoxgRO1KikpqbJfVlYWhYWFWaoGgEzo1atXtksAAGAzYA0nAAAAABIlcAIAAAAgUQInAAAAABIlcAIAAAAgUQInAAAAABIlcAIAAAAgUQInAAAAABKVl+0CaLqKi4ur7FdUVGSpEgAAAKA5McMJAAAAgESZ4UStSkpKquyXlZVFYWFhlqoBAAAAmgsznAAAAABIlMAJAAAAgEQJnAAAAABIlMAJAAAAgEQJnAAAAABIlMAJAAAAgEQJnAAAAABIlMAJAAAAgETlZbsAmq7i4uIq+xUVFVmqBAAAAGhOzHACAAAAIFFmOFGrkpKSKvtlZWVRWFiYpWoAAACA5sIMJwAAAAASJXACAAAAIFECJwAAAAASJXACAAAAIFECJwAAAAASJXACAAAAIFECJwAAAAASJXACAAAAIFECJwAAAAASJXACAAAAIFECJwAAAAASJXACAAAAIFF52S6Apqu4uLjKfkVFRZYqAQAAAJoTM5wAAAAASJQZTtSqpKSkyn5ZWVkUFhZmqRoAAACguTDDCQAAAIBECZwAAAAASJTACQAAAIBECZwAAAAASJTACQAAAIBECZwAAAAASJTACQAAAIBECZwAAAAASJTACQAAAIBECZwAAAAASJTACQAAAIBECZwAAAAASJTACQAAAIBECZwAAAAASJTACQAAAIBECZwAAAAASJTACQAAAIBECZwAAAAASJTACQAAAIBE5WW7AJqu4uLiKvsVFRVZqgQAAABoTsxwAgAAACBRZjhRq5KSkir7ZWVlUVhYmKVqAAAAgObCDCcAAAAAEiVwAgAAACBRAicAAAAAEiVwAgAAACBRAicAAAAAEiVwAgAAACBRAicAAAAAEiVwAgAAACBRAicAAAAAEiVwAgAAACBRAicAAAAAEiVwAgAAACBRAicAAAAAEiVwAgAAACBRAicAAAAAEiVwAgAAACBRAicAAAAAEiVwAgAAACBRAicAAAAAEiVwAgAAACBRAicAAAAAEiVwasY++eSTOOecc2LPPfeMzp07R9++fePiiy+OefPmZbs0AAAAoAXLy3YBbJrx48fHN7/5zVi5cmX62IIFC+L999+P+++/P15//fXYdddds1ghAAAA0FKZ4dRM/eIXv4iVK1fGiSeeGO+9914sWLAg3njjjRgwYEAsXrw4rrrqqmyXCAAAALRQZjg1Ux988EG0bt06Hn744SgoKIiIiE6dOsVDDz0UhYWF8eGHH2a5QgAAAKClMsOpmercuXO0atUq8vKqZoZt2rRJtwMAAABkgxlOm2jVqlXxwQcfxMSJE2PJkiURETFgwIAYMGBAvcZJpVLx9ttvx5QpU2LJkiXRuXPnOPDAA6Nnz551nnfuuefGr371qzjvvPPimmuuiS5dusTMmTPjxz/+cUREDB8+fNPeGAAAAEADCZzqYfHixXH55ZfHhAkTYtKkSVFRUVGl/aqrrqpX4PT888/HhRdeGJ988km1tqOPPjpuv/32WoOnX/7ylxERceONN8aoUaMiLy8v1qxZE927d48777wzvv/979fjnQEAAAAkR+BUDwsXLow777wzIr66dW2//faL7bffPp5++ul6j3X//ffHD37wg0ilUrHFFlvEiSeeGF26dInJkyfHP/7xj3jxxRdj//33j9dffz122WWXaucvXbo0pk+fHitWrIiIiDVr1kRExPLly2PKlCmxatWq9O11AAAAAJkkcKqHjh07xu233x777rtv7LnnnpGfnx/jxo2rd+A0ffr0GDZsWKRSqejdu3e89NJLUVhYmG5/4YUX4uSTT44FCxbEoEGDYuLEidGq1X+W20qlUnHCCSfEuHHj4uc//3kMGzYsOnfuHLNnz45bbrklbr755pgzZ048/PDDib13AAAAgI1l0fB66NixY5x33nnRr1+/yM/P3+Rxrrrqqli9enVERDz44INVwqaIiIEDB8aVV14ZERHvv/9+teBo0qRJMW7cuDjmmGPi+uuvj9122y223Xbb2HPPPePuu++OvffeO/72t7/FggULNrlGAAAAgE0lcMqwL7/8MsaMGRMREYcffnjss88+Nfa78MILo3Xr1hHx1e1365s/f35ERLRt27bGc9fdSreuHwAAAEAmuaUuw1577bUoLy+PiK9mMtWmY8eO0b9//xg3bly8/PLLUVFRkZ5VVVxcHK1atYonn3wyrr/++hg2bFhss8028fnnn8fNN98cb731Vmy55Zax884716u2srKyOtvnzp1br/EAAACAlknglGGTJk1Kb++111519u3bt2+MGzcuVq9eHaWlpdGnT5+IiOjatWtcfPHFcfPNN8dll10Wl112WeTn51d5at5vfvObWmdA1ebrt/YBAAAAbAq31GXYJ598kt7u0aNHnX133HHH9PbHH39cpe3GG2+MO++8M72eVEVFRbRp0yYOPvjgePjhh+PSSy9NtnAAAACAjWSGU4YtXbo0vb311lvX2Xf99mXLllVpy8nJieHDh8fw4cMjlUrFqlWr6j2j6etmzZpVZ/vcuXOjf//+DXoNAAAAYPMncMqwFStWpLfXLQpem/UDpPXP+7qcnJwGh00REd27d2/wGAAAAABuqcuwgoKC9Pbq1avr7Lty5coazwMAAABoysxwyrC6bpP7uvXbt9pqq0arqTbFxcVV9tdflBwAAACgNmY4ZdhOO+2U3v7ss8/q7Lt+e8+ePRurJAAAAIBEmeGUYXvssUd6e9KkSTFw4MBa+3744YcREZGXlxe9e/du9Nq+rqSkpMp+WVlZFBYWZrwOAAAAoHkxwynDDjvssGjTpk1ERLz44ou19lu6dGm8+eabERFxxBFHbHCBcQAAAICmQuCUYVtvvXWccMIJERHx97//PaZMmVJjvzvvvDO9aPiZZ56ZsfoAAAAAGkrglAW/+93vIjc3N9auXRtnnXVWLFiwoEr7+PHj46qrroqIiN122y0GDx6cjTIBAAAANok1nOrpwQcfrLKY98yZM9Pb48aNiz/+8Y9V+v/yl7+MnJycKsf69OkTf/rTn+LCCy+M999/P3bdddcYNGhQdOnSJSZPnhxPPPFErFmzJtq3bx8PP/xw5ObmNu6bqoWn1AEAAACbIieVSqWyXURzcvjhh8drr7220f0rKioiL6/mXO/hhx+Oiy++OObPn1+trX///jFy5MhqoU8m1RQ4lZaWRkTErFmzonv37tkoCwAAAEjI+g8IS/J3fTOc6unMM8+MAw44YKP7t2pV+12Lp59+enz7299Or+W0ZMmS6Ny5cxx00EHRt2/fJMptEE+pAwAAADaFGU5stMZKPQEAAIDsaKzf9S0aDgAAAECiBE4AAAAAJErgBAAAAECiLBpOrWp6Sh0AAADAhpjhBAAAAECizHCiViUlJVX211+5HgAAAKA2ZjgBAAAAkCiBEwAAAACJEjgBAAAAkCiBEwAAAACJEjgBAAAAkChPqaNWxcXFVfYrKiqyVAkAAADQnJjhBAAAAECizHCiViUlJVX2y8rKorCwMEvVAAAAAM2FGU4AAAAAJErgBAAAAECiBE4AAAAAJErgBAAAAECiBE4AAAAAJErgBAAAAECiBE4AAAAAJCov2wXQdBUXF1fZr6ioyFIlAAAAQHNihhMAAAAAiTLDiVqVlJRU2S8rK4vCwsIsVQMAAAA0F2Y4AQAAAJAogRMAAAAAiRI4AQAAAJAogRMAAAAAiRI4AQAAAJAogRMAAAAAiRI4AQAAAJAogRMAAAAAicrLdgE0XcXFxVX2KyoqslQJAAAA0JwInIAWb9q0aVFRURH5+fnRq1evbJcDAADQ7AmcqFVJSUmV/bKysigsLMxSNdB4pk2bFuXl5VFQUCBwAgAASIA1nAAAAABIlMAJAAAAgEQJnAAAAABIlMAJAAAAgEQJnAAAAABIlMAJAAAAgEQJnAAAAABIlMAJAAAAgEQJnAAAAABIVF62CwDIptLS0hg9enTMmTMnunXrFr17946ioqJslwUAANCsCZyAFmvUqFExfPjwqKysTB8bM2ZMjBgxIoYMGZLFygAAAJo3t9QBLVJpaWm1sCkiorKyMoYPHx6lpaVZqgwAAKD5M8OJWhUXF1fZr6ioyFIlkLyRI0dWC5vWqaysjJEjR8Yf/vCHDFcFAACweTDDCWiRZsyY0aB2AAAAameGE7UqKSmpsl9WVhaFhYVZqgaS1bNnzwa1AwAAUDsznIAWaejQoZGbm1tjW25ubgwdOjTDFQEAAGw+BE5Ai1RUVBQjRoyoFjrl5ubGiBEjoqioKEuVAQAANH9uqQNarCFDhsSAAQPiiiuuiDlz5kS3bt3immuuETYBAAA0kMAJaNGKiori7LPPjvLy8igoKBA2AQAAJMAtdQAAAAAkSuAEAAAAQKIETgAAAAAkSuAEAAAAQKIETgAAAAAkSuAEAAAAQKIETgAAAAAkSuAEAAAAQKIETgAAAAAkSuAEAAAAQKIETgAAAAAkSuAEAAAAQKIETgAAAAAkSuAEAAAAQKLysl0ATVdxcXGV/YqKiixVAgAAADQnAiegxevVq1dUVFREfn5+tksBAADYLAicqFVJSUmV/bKysigsLMxSNdB4evXqle0SAAAANivWcAIAAAAgUQInAAAAABIlcAIAAAAgUQInAAAAABIlcAIAAAAgUQInAAAAABIlcAIAAAAgUQInAAAAABIlcAIAAAAgUQInAAAAABIlcAIAAAAgUQInAAAAABIlcAIAAAAgUQInAAAAABIlcAIAAAAgUQInAAAAABIlcAIAAAAgUQInAAAAABIlcAIAAAAgUQInAAAAABIlcAIAAAAgUXnZLoD6W7t2bTz77LN19tliiy3iyCOPzFBFAAAAAP8hcGqGVq9eHSeeeGKdfYqLi2Py5MkZqggAAADgPwROzVBubm4cf/zxNbZ98sknMXXq1Bg0aFCGqwIAAAD4isCpGcrPz4+nn366xrYTTjgh/vWvf8XZZ5+d4aoAAAAAvmLR8M1IWVlZPP/883HkkUfGTjvtlO1yAAAAgBZK4LQZGTlyZFRWVsbQoUOzXQoAAADQgrmlrp7Wrl0bH330UbzzzjsxYcKEmDhxYixZsiQiIi644IK44IIL6jXe66+/Hg899FBMmTIllixZEp07d46DDjooBg8eHD169KhXXSNHjoz27dvHqaeeWq8aAAAAAJIkcKqHzz77LIqLi2P58uU1ts+fP3+jx1q1alUMHz487rvvvmptzz33XFxzzTVx4403xvnnn79R47300ksxc+bMOP/886OgoGCj6wAAAABImsCpHlavXh0rVqyIXXfdNfbdd9/Yb7/9IicnJy6++OJ6jzV48OB45JFHIiLiuOOOix/96EfRpUuXmDx5cvzhD3+IadOmxY9+9KMoKCjYqAXAR4wYERHhdjoAAAAg6wRO9dCjR49YsmRJtGvXLn1s3Lhx9R5nzJgx6bDphz/8YYwaNSrdtu+++8Ypp5wS++23X0yfPj0uueSSOP7446NTp061jjd//vx48skno7i4OPbbb7961wMAAACQJIuG10N+fn6VsGlTXXPNNRER0a5du7jllluqtXfo0CFuuOGGiIhYsmRJ/OUvf6lzvNGjR0dFRYXZTQAAAECTIHDKsE8//TQmTJgQERGnnnpqbL311jX2O+GEE2K77baLiEjPhqrNXXfdFfn5+TF48OBkiwUAAADYBG6py7Dx48entw877LBa+7Vq1SoOOeSQeOyxx2Ly5MmxZMmSaN++fbV+//znP2PatGlxyimnpAOqTVVWVlZn+9y5cxs0PgAAANAyCJwyrKSkJL3dq1evOvuu3z5lypQ48MADq/VJcrHwwsLCBo8BAAAA4Ja6DFt/llC3bt3q7Lt+e02ziyorK2PlypXx3e9+N4499tjkigQAAABoADOcMmzZsmXp7S222KLOvuu3r3/eOrm5uRtc36k+Zs2aVWf73Llzo3///om9HgAAALB5EjhlWGVlZXo7Nze3zr55ef/5n2fNmjWNVtM63bt3b/TXAAAAADZ/bqnLsHbt2qW3y8vL6+y7YsWK9PZWW23VaDUBAAAAJMkMpwzbZptt0tsLFy6MHj161Np34cKFNZ6XKcXFxVX2KyoqMl4DAAAA0PyY4ZRhu+66a3p7xowZdfZdv3398wAAAACaMjOcMqxv377p7TfffDO+853v1Nr3jTfeiIiIbbfdNgoLCxu9tq8rKSmpsl9WVpaVOgAAAIDmxQynDNt///2ja9euERHx+OOPRyqVqrHf1KlTY+rUqRERccopp2SqPAAAAIAGEzhlWKtWrWL48OERETF9+vS47777aux31VVXRURETk5OnHvuuRmrDwAAAKChBE5ZcNlll6UXC//Rj34U9957b1RWVkZExJIlS+KSSy6JRx55JCIifvCDH0T//v2zVisAAABAfeWkarunixoNGTIk3nnnnfT+ihUr0ot7b7fddtG5c+cq/T/44IPIzc2tNs7UqVPjyCOPjHnz5kVERPv27aNTp05RVlYWq1atioiII488Mp555plo27ZtY72dOtX0lLrS0tKIiJg1a1Z07949G2UBAAAACVl/veYkf9e3aHg9zZgxo9pi2ussWLAgFixYUOVYbXnebrvtFpMnT47f/va38de//jUWLFgQS5YsiYiIPn36xPnnnx8XXnhhtGplEhoAAADQvJjhVE8zZsyIL7/8cqP777777hvsk0ql0oFT586do3379g0psdE0VuoJAAAAZIcZTk1Ez549Ex8zJycnOnfuXO12PAAAAIDmyP1aAAAAACRK4AQAAABAotxSx0Zbs2ZNenvu3LlZrAQAAABIwvq/36//e39DCZyoVXFxcZX9FStWpLf79++f6XIAAACARrRgwYLYaaedEhnLLXVstMrKymyXAAAAADQDOalUKpXtImgepk+fHkVFRRER8frrr6cfm7i5OPLIIyMi4h//+Mdm9/pJjN2QMep7bn36b0zfuvrMnTs3PWPv7bffjq5du25Ujc2F67rxxmis6zqJfq7r5vv6rmvXtes6+TFc19mzuV7XSY27qeNsynmu6+Rsjtf1mjVrYsGCBRERsccee0Tbtm0TGdctdWy09S+6wsLC6N69exarSV5+fn5ERNbeV2O+fhJjN2SM+p5bn/4b03djx+vatavruhm9/uZ6XSfdz3XdvF7fde26dl0nP4brOns21+s6qXE3dZxNOc91nZzN9bpO6ja69bmlDgAAAIBECZwAAAAASJTACQAAAIBEWTScjVZWVpZeKHzWrFmb3b24tEyuazZHrms2R65rNkeuazZHrmvWMcMJAAAAgEQJnAAAAABIlMAJAAAAgERZwwkAAACARJnhBAAAAECiBE4AAAAAJErgBAAAAECiBE4AAAAAJErgBAAAAECiBE4AAAAAJErgBAAAAECiBE4AAAAAJCov2wUAAM3Pl19+GW+//XasWrUq9ttvv9h2222zXRJAi1ReXh5vvfVWrFy5MvbZZ5/o3LlztksCiAgznGgiJk6cGD/+8Y+jT58+0bZt29hiiy1i//33j7vuuitSqVS2y4N6+/LLL+PBBx+Mk046KXr06BF5eXmx7bbbxkknnRTjxo3Ldnmwye6+++44/vjjo1OnTnHkkUfGscceG++99162y4IN+uCDD+L444+PrbfeOtq0aRP77LNPPPDAA9kuCzbZ/fffHyeffHJ06tQpjjjiiDj22GPj9ddfz3ZZsMmWLFkSI0eOjGOOOSZ22GGHyM3Nje233z5OO+20eP/997NdHpsgJ+W3eZqAnJycWtvOPffcuOOOOzJYDTTc9ddfH5dddlmNba1atYr7778/zjjjjAxXBQ3XqVOnWLhwYWyxxRbRrl27mD9/frz00kvxzW9+M9ulQa3eeeedOOyww6K8vLxa2/XXXx8///nPs1AVNMxOO+0UM2fOjLZt20aHDh1i3rx5MWbMmDjllFOyXRpskksvvTRuuOGGGttat24dTz/9dHzrW9/KcFU0hBlONAkHH3xwjBw5Mj766KNYvXp1fP755/GnP/0p8vLy4s4774zJkydnu0Sol3bt2sWZZ54ZTz31VHz66aexevXq+Ne//hXf/va3Y+3atXHxxRfHmjVrsl0m1NuwYcPi6aefjoULF8app56a7XJgo/zoRz+K8vLy+OEPfxjz5s2L8vLyGDFiROTl5cWVV14Zs2bNynaJUG+DBw+OJ554IhYuXBhnnXVWtsuBBuvQoUMMGzYsXnjhhSgrK4vVq1fH+++/H0cddVSsXr06Lr744myXSD2Z4USTNnz48Ljrrrvi/vvv9x9SNgsVFRWxyy67xGeffRbTp0+PnXfeOdslwSY7//zz44477jDDiSbtgw8+iL333jt69+4dJSUlkZubm2776U9/GjfffHNcc8018atf/SqLVULDrJsZYoYTm6Ply5dH165dY/ny5bF8+fLYcssts10SG8kMJ5q0dV8Kt99++yxXAsnIz8+Pnj17RkRE+/bts1wNwOZv7NixERFx+umnVwmbIiL9x6x1fQBoetq1axc77LBDtG3bNtq2bZvtcqgHT6nbDKxevTomTZoUEyZMiOnTp0cqlYo2bdrE1VdfXe+xpk2bFn//+99jzpw5UVBQELvssksce+yxsdVWWzVC5dWtWbMmUqlULFiwIB5//PEYNWpUFBcXx+GHH56R16fpWLlyZXzwwQcxYcKE+PTTTyOVSkXHjh3jiiuuqPdYkydPjldffTXmzZsX7dq1i969e8cxxxwTBQUFjVB53WbNmhVvvvlmfPOb34xOnTpl/PXJnvnz58eECRPinXfeiQkTJkRZWVmkUqkoLi6u98LFZWVl8cADD8T48eNj/vz5sfXWW0dxcXF873vfi/3337+R3gFsmlQqFdOnT48JEybEpEmTYvXq1RERcfHFF8eOO+5Yr7EWLlwYzz77bHz66adRWVkZ3bt3j4EDB0ZhYWGt50yfPj0iInbfffdqbbvvvnvk5OSk+8DGWLlyZbz//vvpz/OSkpL0bfJPPfVUndfj161atSoee+yxePrpp2PWrFmRSqWiR48eccIJJ8S3v/3taN26dWO9DaimsrIypkyZEhMmTIh//etfUVlZGRER//3f/13voOezzz6LF198McrKyiI3Nzd69uwZxx577CY91faDDz6IadOmxVlnnVXtDwc0cSmarVtvvTW13377pdq0aZOKiCo/W265Zb3GmjVrVurYY4+tNk5EpAoKClLXXnttqrKystbzKyoqNvpnzZo1tY6x/uu2atUq9f3vfz/1+eef1+u90LxdffXVqb59+6by8/OrXYs9evSo11ilpaWpQw89tMbreuutt0795S9/qfP8+lzXdf3/Y52VK1emDj744FS7du1SU6dOrdd7oXk75JBDarwOIyLVr1+/eo01YsSI1JZbblnreEOGDEmtWrWqxnPnz5+feu655zb6p6SkpM5azjvvvFREpF566aV6vQdahldffTV11FFHpTp06FDjtfrGG29s9FgrV65MXXrppam8vLxq4+Tk5KQGDx6cWrRoUY3nDhkypM7rdIsttkh17tx5k94jLc8tt9xS43W47qe0tHSjx5o8eXKqT58+tY61++67p6ZMmbJRY/385z9PRURqzJgxm/jOaMkefPDB1MEHH5zaYostarwWFy9evNFjLVy4MHXmmWemcnJyqo2Tn5+f+sUvfpFauXLlRo/373//O1VcXJzq3Llzavbs2Zvw7sgmM5yasZdffjneeeedyM/Pj759+8a+++4bzz//fL0Xvpw7d24ccsgh8emnn0ZOTk4MHDgwDjjggFi+fHk89thj8cknn8Qvf/nL+Oyzz+LPf/5ztfOnT58eRUVFG/16hx12WLz66qvVjufk5ERubm6kUqlYu3ZtrF27NsaPHx+vvfZaDBo0qF7viebr2Wefjffeey/atGkT++23X+y7777xt7/9LRYuXFivcT7++OM4+OCDY/78+ZGbmxsnnHBC7LPPPrF48eJ45JFHYvbs2XHBBRfEvHnz4re//W21899888048MADN/r1vvOd78Sjjz5aa/vKlStj0KBBMWHChHjiiSdi1113rdf7oXmbP39+dO7cOfr16xf77rtv7LvvvnHWWWfF8uXL6zXO/fffH8OHD4+IiG984xtx+eWXxx577BGff/553H777fH888/HqFGjory8PB566KFq57/11ltx4oknbvTrnXfeeXH77bfXq0ZYZ9KkSfH3v/89IiJ23nnn2HfffWPOnDn1vn2tsrIyTjvttHjyyScj4qtZSSeeeGLk5+fHK6+8EmPHjo377rsvpkyZEq+++mq0a9euXuOnLGdKPSxatCjy8vKqfJ4/8cQT8fjjj9drnFmzZsW3vvWtmDt3brRu3Tp+9rOfxcCBAyMi4vnnn48bb7wxJk+eHN/61rfirbfeih122KER3g18Zfz48TF+/PjIzc2N4uLi2HfffeOdd96JKVOm1GucpUuXxlFHHRXvv/9+REQccsghccQRR0RFRUU89dRTMXny5Pif//mfKC0tjUcffTRatap7hZ8lS5bE8ccfH7NmzYqXXnopunXrtqlvkWzJduLFpnv55ZdTb731VpWEeP/996/3DKfjjjsu/RfCv/71r1XaVq5cmTrppJPSqXRNfzWZPn16Kjc3d6N/jjzyyA3WtHDhwtRjjz2W6tGjRyoiUq+++upGvx+at2effTY1ceLE1OrVq9PHdt5553rPcFr3/4X8/PzUCy+8UKVt2bJlqcMPPzx9Xf/zn/+sdv6bb75Zr+v6tNNOq7WWpUuXpo444ohU27ZtU88999xGvwc2H/Pnz692rH379vWa4fTFF1+kZ4p84xvfSC1cuLBK+9q1a1ODBw9OX9dPPPFEtTHeeuut1MCBAzf6509/+lOdNZnhRF0mTZqUeumll6r8ZfzKK6+s9wynG2+8MX3O+eefX21G6XXXXZduv/DCC6ud/5Of/CQVEdW+46RSqdSKFStSEZEqKiqq35ujxVq0aFGqoqKiyrF1n4VRjxlO3/3ud9Pfv5999tlq7U899VR6zDPOOGOD45nhREOMHz8+NXbs2NSXX36ZPvad73yn3jOczj333PQ5N954Y5W2ysrK1PDhw9Ptt956a51jzZ8/P9W3b99U+/btU2+++Wa93xNNg8BpM1PfwGn8+PHp/9MPHz68xj6LFy9OdezYMRURqT59+iRZ7gY9//zzqYhInXPOORl9XZqW+gZOTzzxRPq6/tWvflVjn9mzZ6fatm2biojUIYcckmC1VS1YsCDVr1+/VEFBQerFF19stNeh+alv4PS73/2uzvA/lfrq83rd7Xb9+/dPsNqaCZyor/oGTuXl5envIL17967yx4j1fetb30pFRCovLy81c+bMKm233XZbKiJSv/71r6ud99Zbb6UiInXcccdt2huCVP0Dp48++ijd//TTT6+137pf+HNyclIff/xxnWMKnEhafQOnjz/+OJWbm5uKiNSxxx5bY59Vq1aldtlll1REpDp16lTrEgCffvppqlevXqkOHTqk3nnnnYa8DbLMU+pauPUXqr344otr7NOhQ4cYPHhwRERMmTIlPUUyE1asWBEREcuWLcvYa9L8rbuuc3Jy4qKLLqqxT7du3eI73/lORESMGzeu3reibozPPvssBgwYEFOnTo2nn346vvWtbyX+GrQcjzzySEREbLvttrXeFtehQ4f49re/HRERb7/9dsycOTNj9UFjeO6552LRokUREXH++edHfn5+jf3WfYdZs2ZN/O1vf6vSdsghh0RExMMPP5xe2Hmd++67LyK+ut0fMmX9W/CHDBlSa7+hQ4dGxFe3fdZ12z40BQ8//HB6kfHafq9s3bp1nHfeeRER8cUXX8SLL75Yrc+UKVPi4IMPji+++CL+/ve/x7777tt4RdPoBE4t3EsvvRQREV26dKnx6S3rHH300entmj4YGuKee+6JYcOGxauvvhqzZ8+OFStWxKeffhp33313+gPp0EMPTfQ12bytu66Li4vrvNd73XWdSqXS5yRl2rRpccghh8SsWbPiqaeeikMPPTTWrFlT5Sdl3RA20pIlS2LSpEkRETFgwIA6n9Cy/lM9x48fn3gtH3zwQTz//PNV1gycMGFC+tj8+fMTf01arvU/m+sK7Y844oh0GPX17yl77LFH9O/fP0pLS+Pss8+OmTNnxpIlS+LPf/5z/OUvf4mCgoI466yzGucNQA3WrWPWqlWrGDBgQK39DjnkkPQaN+PGjavWXlJSkv7s/fTTTyMi4r333ksfmzt3bvLFQy3WfV7n5+fXGeLX9XvlxIkT49BDD42VK1fGSy+9FHvuuWe17880LxYNb8HKy8vj448/joivvozVZc8990xvT548OdE6Vq1aFSNHjoyRI0fW2H744YenF8mFDZk9e3YsXrw4IrJ7Xd95553x2WefRUTEUUcdVWOfV155pUo4ALUpKSlJb/fq1avOvuu3r39eUv7whz/Eww8/XOXYr371q/T2mDFj4pRTTkn8dWmZ1gWteXl5dT5sYcstt4xvfOMb8dFHH9X4eX777bfHIYccEg8++GA8+OCDVdquvfZaCzKTUes+m3fYYYfYcssta+231VZbRdeuXWP27Nk1fp7fdNNNcffdd1c59rvf/S69fd9998X3v//9hKqGuq377N1ll12ioKCg1n7FxcWRm5sblZWV1T6vb7rppvSDgvr161fj+VOnTvXwnWZE4NSCffbZZ7F27dqIiCgsLKyzb7du3dIfDOv+gpKUIUOGRIcOHeK+++6LyZMnx+effx7t27eP3XffPU4//fQYMmRI5OW5VNk461+fG7qud9xxxxrPS0Jubm6ds1AivrrlDzbG+n+l3tATWtb/xbkx/rq99957x7///e9a27fffvvEX5OWa91n87rvIXXZcccd46OPPoq5c+fG6tWro3Xr1um2vn37xltvvRX/9V//FWPHjo1Vq1ZFcXFxXHrppenbqyFT5s2bFxEb/jyP+Oozffbs2TV+nu++++7pJ9vVxBO9yJQvv/wyFixYEBEb/v6dm5sb3bp1i1mzZlX7/u378+bHb/Et2PrrIm211VZ19m3VqlVsueWWsXTp0sTXU2rdunWcfvrpcfrppyc6Li1Tfa7r9duTvq6vvfbauPbaaxMdk5Zr+fLl6e26/hoeEbHFFlvUeF5SLr/88rj88ssTHxdqsu6zeUOf51/vs2zZsth2222rtBcXF8djjz2WbIFQT6tWrYqKioqI2PDnecR/PtPLy8ujsrKyyi/jP/nJT+InP/lJo9QJ9VGf79/r9/n69+/Ro0fH6NGjky2OrLKGUwtWXl6e3l7/r4C1adOmTbXzoKmpz3Wdn5+f/iuJ65qmbP01Czb0l7/1Z4Su+6UGmqt1n831+Z6y/nnQ1NTn8zyi6me69WtoqvxeSW0ETi3Y+n8FX7Vq1Qb7r1y5stp50NTU57pevXp1euFu1zVN2fp/Bd/Ql7N1T/eMiGjXrl2j1QSZsO6zuT7fU9Y/D5qaLbbYol5/7Fr3mZ6Xl1clVIWmxO+V1Ebg1IJtvfXW6e2lS5fW2beysjL9H7z1z4Ompj7X9frtrmuaso4dO6a31y2mWZv127fZZptGqwkyYd1n84Y+z7/eZ2Nu6YBsyMnJiQ4dOkTEhj/P1+/j85ymrD7fv9fv4/v35k/g1IL16NEjPU133dO0alNWVhaVlZUREbHzzjs3em2wqda/Pjd0Xc+cObPG86Cp6d27d3p7xowZdfZdfwFOT3GhuVv32Tx37twN3k607jN9xx13jPz8/EavDTbVus/0mTNnpmda12Tt2rXp7zI+z2nKCgoKomvXrhGx4e/fFRUV6UXwff/e/AmcWrDWrVunH5/9wQcf1PkfvPfeey+9vfvuuzd6bbCpOnfuHJ07d46IiPfff7/Ovq5rmosdd9wxPcvprbfeqrPvG2+8kd7u27dvo9YFjW2PPfaIiK9mWk+aNKnWfkuXLk2HsT7PaerWfTavWLGizuv6gw8+SN925/Ocpm7d5/XHH39c58N4Pvjgg/ST0n1eb/4ETi3cMcccExFfTdedMGFCrf2ef/759Paxxx7b6HVBQ6y7rqdNm1bnbJB113WrVq3i6KOPzkhtsClycnLi5JNPjoiIqVOnxr/+9a8a+6VSqXj88ccjIqJLly6x//77Z6pEaBTrPs8jqn4X+boXX3wxPRPb9xSaulNOOSW9PWbMmFr7rf9UxfXPgaZo3ed1ZWVlvPzyy7X283tlyyJwauG+//3vp7dvvPHGGvvMmzcvHnjggYiI2HfffU3ppcn7wQ9+kN6u7bqePn16PPHEExERcfTRR8f222+fkdpgU5177rnp7d/85jc19nnggQdi2rRpERExfPjwaNXKf+Zp3r75zW+mb9O4/fbbqyyKv04qlYqbbropIiLatm0b3/3udzNaI9TXkUceGUVFRRERcdttt8WCBQuq9fn888/jf//3fyPiq9vpDj300IzWCPV12mmnpRe2v+mmm2q8e+bLL7+MO+64IyIiunfvHocffngmSyQLfBNt4fr27Rvf+973IiLir3/9a9x6661V2pcsWRKnn356LF++PCIi/vjHP2a8Rqivo446Kr75zW9GRMRf/vKXuP/++6u0z58/P0477bRYs2ZNtGrVKq655ppslAn1csABB6T/SPDwww/HT3/60/SU9bVr18Zf//rXOP/88yMiorCwMH7xi19krVZISn5+fvz+97+PiK/WBfnhD39Y5QlIqVQqfvnLX8brr78eERE//elPo0uXLlmpFTZWXl5e3HjjjZGTkxMLFiyIY445psrM1alTp8axxx4bCxcujJycnLjpppsiNzc3ixXDhu2www7x4x//OCIixo4dG7/61a+qhE6rVq2Ks88+O8rKyiIi4uqrr06vJ8zmKydV18I9NGlTp06Nu+++u8qxBx98MObOnRv5+flx8cUXV2k7/vjj44gjjqg2zr///e845JBDYvLkyRERsffee8eBBx4Yy5cvj2eeeSYWLVoUEV/9Rf2qq65qpHcDX5kwYUL89a9/rXLs7rvvjn//+9+x9dZbx/Dhw6u0nXbaadG/f/9q48ybNy8OPPDA9ALK+++/f/Tr1y8WLVoUzzzzTPoX9T//+c9x4YUXNs6bgf/fvffeW2223aRJk2Lt2rVRUFCQXk9vnZtuuqnGz+vy8vI49thj47XXXouIr2Zz9OjRIxYsWJD+rN5uu+3ilVdeieLi4kZ6N7Bx1qxZE5dffnmVY+PHj48333wzIiLOPPPM9OyliIhevXpVmcm3viFDhsQ999wTERHdunWLgQMHRn5+frz66qvpWX1HHXVUPPvss9G6detGeDfwlfnz51e7Db+srCz9NLnddtutyjV4zDHH1PoH2xtuuCEuu+yySKVSkZOTEz169IiI/ywm3qpVq7jxxhvjkksuaaR3A1+ZM2dOte8pTz31VPrz9aKLLkrPXoqIGDBgQI23ea5atSoGDhyY/p7Su3fvOOyww6KioiKef/759GLh55xzTowYMaKR3g1NSopm66mnnkpFxEb//OEPf6h1rMWLF6eGDBmSysvLq3be9ttvn7rnnnsy+M5oyUaNGlWv63rEiBG1jjVv3rzUoEGDUq1atap23o477ph67LHHMvjOaMmuu+66el3XY8aMqXWsNWvWpG688cZU7969q5zTsWPH1Pnnn5+aP39+5t4Y1KG8vLxe1/1RRx1V61iVlZWp6667LtW+fftq57Vt2zZ16aWXpsrLyzP47mipZs2aVa/r+vTTT69zvNdeey31rW99K5Wbm5s+Jy8vL3X00Uenxo4dm6F3RUv33nvv1eu6vvDCC2sda8WKFamf/exnqbZt21Y7r0OHDqkbbrghVVlZmcF3RzaZ4dSMffLJJ1UWE9yQQw89tMaZIOtbuHBh/POf/4w5c+ZE27Zto6ioKA4++GDTeMmYkpKSeO655za6/8CBA9NPxajNvHnzYty4cTFv3rzYcsstY9ddd40DDjggcnJyGloubJQFCxbE7NmzN7p/z549o3379hvst2jRovj8889j6623ji5duvispkmprKxMr620MXr06BGDBg2qs8+qVavin//8Z8yYMSMqKyujsLAwDjvssNhqq60aWi5slIqKiigpKdno/ttss0165lJdysvLY/bs2ZGTkxPdunWLgoKChpQJ9bJgwYIYPXr0Rvfv27dvHHXUUXX2WbZsWbz66qtRVlYWubm50bNnzzjssMPMQm1hBE4AAAAAJMqi4QAAAAAkSuAEAAAAQKIETgAAAAAkSuAEAAAAQKIETgAAAAAkSuAEAAAAQKIETgAAAAAkSuAEAAAAQKIETgAAAAAkSuAEAAAAQKIETgAAAAAkSuAEAAAAQKIETgAAAAAkSuAEAAAAQKIETgAAAAAkSuAEAAAAQKLysl0AAACbbsaMGXHLLbdUO37YYYfFqaeeWuXY9OnT489//nNERJx//vmx6667ZqTGJH355Zdx5ZVXVju+2267xXnnnZeFigCAmgicAACasdmzZ9cYOEVEtcCprKws3feYY45ploFTeXl5je934MCBAicAaEIETgAAm4nzzjsvHSL17ds3y9U0jnbt2sVNN92U3r/uuutizpw5WawIAKiJwAkAYDNxyimnxDHHHJPtMhpV27Zt4yc/+Ul6/5577hE4AUATZNFwAAAAABIlcAIAAAAgUW6pAwDIgl/84hexevXqOPLII+Okk06K5cuXxyuvvBKTJk2KhQsXRiqViuuuuy5yc3MzXtuTTz4Z//jHPyIiYq+99oohQ4ZERMT8+fPjmmuuiYiIM888M/r37x+pVCpeeumlePfdd2PJkiWxww47xLHHHhs777zzRr1WKpWK9957L954442YN29erF27NnbYYYfo1atXHH744ZGX5+sqADRH/gsOAJBhM2fOjOuuuy4iIvbcc8/4r//6r7jlllti2bJl6T577LFHVsKm66+/Pn7xi19EKpWKY445Jn7/+9+n2xYtWpR+Qtzee+8d22yzTZx66qlRUlJSZYycnJz44Q9/GLfddlsUFBTU+lpPPPFE/OpXv4qpU6fW2N6pU6f4n//5n3TgBQA0HwInAIAMmzhxYnr797//fXz66aex/fbbx5lnnhk777xz5Ofnp582lymVlZVx0UUXxe233x4RXz3x7s9//nOtM4y++OKLOPTQQ6O8vDzOOeec2GWXXeKLL76Ixx9/PKZPnx6jRo2KmTNnxosvvlhjcPa73/1/7d1NSJRdGMbxq3RqNHFRmjYjlFL4EaWOaUFjCaXMopDCTAipiGJoU6tSKoIXihZZRAsLyjatXEgtktHMzKwkKMsI+lqojVFolMUUBDnvQhzGNz/Get551P6/1XHOuT332V6cOfOPjh8/LkmyWCxyuVxyOByyWq16+/atXr16pZaWFj169IjACQCAaYjACQAAIMyCA6euri653W5VVVUpOjralH6+fv2q0tJSeTwezZo1S6dOndKhQ4fGrTlx4oTsdrsaGxtls9kCn588eVJut1s1NTVqbm7W6dOndfjw4RG1dXV1gbBp5cqVqqurG/UreP39/erp6THghAAAINx4NBwAACDMggOnkpISVVdXmxY2eb1eOZ1OeTweWa1W1dbWThg2SUMhVW1t7YiwSRq6rXThwgWlp6dLks6cOaMfP34E5gcHB1VRUSFJmj9/vjwez5jvPcXFxcnhcPzu0QAAgIkInAAAAMJsOHBKTEzU5cuXTeujo6NDq1evVmdnp+Lj43X79m2VlJSEVLthwwZlZGSMOmexWOR2uyUNPTT+4MGDwNzjx4/1+vVrSdK+ffu0aNGiPzwFAACYigicAAAAwqi7u1v9/f2Shn7pLTY21pQ+6uvrtW7dOr17906pqalqb2/XmjVrQq5fu3ZtyPMdHR2B8cOHDwPjwsLCSXQMAACmE95wAgAACKPgr9OVlZWZ1sf58+clSREREaqrq1NKSsqk6hMTE0Oe//jxY2Dc19cXGNvt9kntCQAApg9uOAEAAITRcOBktVqVnZ1tWh979+7VnDlz9PPnTxUXF6u7u3tS9X6/f9z5wcHBwHj27NmjjoPXAACAmYXACQAAIIyGA6fMzExFRpp32Xzr1q26fv26oqKi9ObNGzmdTr18+TLkeq/XO+58b29vYBwfHx8YB9986urqCr1hAAAwrRA4AQAAhNFw4JSTk2NyJ5LL5VJDQ4NiY2Pl9XqVn5+vJ0+ehFTb0tIS8vyqVasC4+C3nerr6yfTLgAAmEYInAAAAMIk+MFwh8NhcjdD8vPz1dzcrLi4OPX19amgoED379+fsO7evXtqb28fde7bt2+qrq6WJC1evFi5ubmBuYyMDOXl5UmSrly5EvjFOgAAMLMQOAEAAIRJ8IPhU+GG07CcnBy1trbKZrNpYGBARUVFunnz5rg1NptNpaWlevbs2YjPP3/+rNLSUvX09EiSjh07poiIiBFrqqqqFBkZKZ/Pp40bN6qtrW3UPZ4+faqmpqY/OBkAADALgRMAAECYDAdOc+fO1fLly03uZqT09HS1tbUpJSVFPp9Pmzdv1rVr18ZcX1lZqejoaDkcDhUVFWn//v3avn27kpOTdePGDUnSjh07tGfPnl9qnU6nLl26JIvFop6eHuXn5yszM1O7d++W2+3Wpk2btGzZMmVlZY3bAwAAmLrMe6kSAADgLzMcOK1YsUIWi8Xkbn6VnJysu3fvqqioSM+fP9e2bdtUU1Oj8vLyX9bOmzdPra2tKi8vV2Nj44gbUVFRUaqoqNDRo0fH3Gvnzp3KyMhQZWWlbt26pc7OTnV2do5Yk5WVpcLCQuMOCAAAwobACQAAIEzKysrkcrmUmppqyv5Lly7V2bNnJQ3daBqNzWbTnTt3dPXqVfn9fn369ElfvnxRbGzsL2sXLlyohoYGvXjxQh0dHRoYGFBSUpIKCgoUExMzYT+5ublqamrShw8f1N7ervfv3ysyMlJ2u12pqalKTk7+swMDAADTEDgBAACEya5du0zdPykpSQcPHpxw3YIFC3TgwIGQ/29aWprS0tJ+u6+EhAQVFxf/dj0AAJh6CJwAAABmiIsXL8rj8UiS1q9fry1btpjckfF8Pp+OHDkS+Lu3t9fEbgAAwFgInAAAAGaI/z6wPRMDp+/fv+vcuXNmtwEAACZA4AQAADCNpaSkBN5lCpadnW1CN/+/mJiYUc+7ZMmS8DcDAADGROAEAAAwjdlstpDeZZoprFbrX3VeAACmKwInAAAAhCQhISFwuygvL8/kbgAAwFQ2y+/3+81uAgAAAAAAADPHbLMbAAAAAAAAwMxC4AQAAAAAAABDETgBAAAAAADAUAROAAAAAAAAMBSBEwAAAAAAAAxF4AQAAAAAAABDETgBAAAAAADAUAROAAAAAAAAMBSBEwAAAAAAAAxF4AQAAAAAAABDETgBAAAAAADAUAROAAAAAAAAMBSBEwAAAAAAAAxF4AQAAAAAAABDETgBAAAAAADAUAROAAAAAAAAMBSBEwAAAAAAAAxF4AQAAAAAAABDETgBAAAAAADAUAROAAAAAAAAMNS/rdm8Ze3r0fYAAAAASUVORK5CYII=", "text/plain": [ "
" ] }, "metadata": { "image/png": { "height": 390, "width": 590 } }, "output_type": "display_data" } ], "source": [ "fig, ax = plt.subplots(1, 1, figsize=(6, 4), layout=\"tight\")\n", "\n", "ax.errorbar(\n", " mwdata1[\"r\"],\n", " mwdata1[\"Menc\"],\n", " yerr=(mwdata1[\"Menc_err_neg\"], mwdata1[\"Menc_err_pos\"]),\n", " marker=\"o\",\n", " markersize=2,\n", " color=\"k\",\n", " alpha=1.0,\n", " ecolor=\"#aaaaaa\",\n", " capthick=0,\n", " linestyle=\"none\",\n", " elinewidth=1.0,\n", ")\n", "\n", "ax.set_xlim(1e-3, 10**2.6)\n", "ax.set_ylim(7e6, 10**12.25)\n", "\n", "ax.set_xlabel(\"$r$ [kpc]\")\n", "ax.set_ylabel(r\"$M(" ] }, "metadata": { "image/png": { "height": 390, "width": 590 } }, "output_type": "display_data" } ], "source": [ "r = np.logspace(-3, 3, 256) * u.kpc\n", "\n", "fig, ax = plt.subplots(1, 1, figsize=(6, 4), layout=\"tight\")\n", "\n", "ax.errorbar(\n", " mwdata1[\"r\"],\n", " mwdata1[\"Menc\"],\n", " yerr=(mwdata1[\"Menc_err_neg\"], mwdata1[\"Menc_err_pos\"]),\n", " marker=\"o\",\n", " markersize=2,\n", " color=\"k\",\n", " alpha=1.0,\n", " ecolor=\"#aaaaaa\",\n", " capthick=0,\n", " linestyle=\"none\",\n", " elinewidth=1.0,\n", ")\n", "\n", "# Use symmetry coordinates for spherical mass_enclosed calculation\n", "fit_menc = init_potential.mass_enclosed(R=r)\n", "ax.loglog(r.value, fit_menc.value, marker=\"\", color=\"#3182bd\", linewidth=2, alpha=0.7)\n", "\n", "ax.set_xlim(1e-3, 10**2.6)\n", "ax.set_ylim(7e6, 10**12.25)\n", "\n", "ax.set_xlabel(\"$r$ [kpc]\")\n", "ax.set_ylabel(r\"$M(" ] }, "metadata": { "image/png": { "height": 390, "width": 590 } }, "output_type": "display_data" } ], "source": [ "r = np.logspace(-3, 3, 256) * u.kpc\n", "\n", "fig, ax = plt.subplots(1, 1, figsize=(6, 4), layout=\"tight\")\n", "\n", "ax.errorbar(\n", " mwdata1[\"r\"],\n", " mwdata1[\"Menc\"],\n", " yerr=(mwdata1[\"Menc_err_neg\"], mwdata1[\"Menc_err_pos\"]),\n", " marker=\"o\",\n", " markersize=2,\n", " color=\"k\",\n", " alpha=1.0,\n", " ecolor=\"#aaaaaa\",\n", " capthick=0,\n", " linestyle=\"none\",\n", " elinewidth=1.0,\n", ")\n", "\n", "fit_menc = fit_potential.mass_enclosed(R=r)\n", "ax.loglog(r.value, fit_menc.value, marker=\"\", color=\"#3182bd\", linewidth=2, alpha=0.7)\n", "\n", "ax.set_xlim(1e-3, 10**2.6)\n", "ax.set_ylim(7e6, 10**12.25)\n", "\n", "ax.set_xlabel(\"$r$ [kpc]\")\n", "ax.set_ylabel(r\"$M(" ] }, "metadata": { "image/png": { "height": 390, "width": 590 } }, "output_type": "display_data" } ], "source": [ "r = np.logspace(-3, 3, 256) * u.kpc\n", "\n", "fig, ax = plt.subplots(1, 1, figsize=(6, 4), layout=\"tight\")\n", "\n", "ax.errorbar(\n", " mwdata2[\"r\"],\n", " mwdata2[\"Menc\"],\n", " yerr=(mwdata2[\"Menc_err_neg\"], mwdata2[\"Menc_err_pos\"]),\n", " marker=\"o\",\n", " markersize=2,\n", " color=\"k\",\n", " alpha=1.0,\n", " ecolor=\"#aaaaaa\",\n", " capthick=0,\n", " linestyle=\"none\",\n", " elinewidth=1.0,\n", ")\n", "\n", "ax.set_xlim(1e-3, 10**2.6)\n", "ax.set_ylim(7e6, 10**12.25)\n", "\n", "ax.set_xlabel(\"$r$ [kpc]\")\n", "ax.set_ylabel(r\"$M(" ] }, "metadata": { "image/png": { "height": 390, "width": 590 } }, "output_type": "display_data" } ], "source": [ "r = np.logspace(-3, 3, 256) * u.kpc\n", "\n", "fig, ax = plt.subplots(1, 1, figsize=(6, 4), layout=\"tight\")\n", "\n", "ax.errorbar(\n", " mwdata2[\"r\"],\n", " mwdata2[\"Menc\"],\n", " yerr=(mwdata2[\"Menc_err_neg\"], mwdata2[\"Menc_err_pos\"]),\n", " marker=\"o\",\n", " markersize=2,\n", " color=\"k\",\n", " alpha=1.0,\n", " ecolor=\"#aaaaaa\",\n", " capthick=0,\n", " linestyle=\"none\",\n", " elinewidth=1.0,\n", ")\n", "\n", "fit_menc = fit_potential_v2.mass_enclosed(R=r)\n", "ax.loglog(r.value, fit_menc.value, marker=\"\", color=\"#3182bd\", linewidth=2, alpha=0.7)\n", "\n", "ax.set_xlim(1e-3, 10**2.6)\n", "ax.set_ylim(7e6, 10**12.25)\n", "\n", "ax.set_xlabel(\"$r$ [kpc]\")\n", "ax.set_ylabel(r\"$M(" ] }, "metadata": { "image/png": { "height": 390, "width": 590 } }, "output_type": "display_data" } ], "source": [ "mwpot_v1 = gp.MilkyWayPotential(version=\"v1\")\n", "\n", "r_grid = np.linspace(0.5, 30, 128)\n", "\n", "fig, ax = plt.subplots(1, 1, figsize=(6, 4), layout=\"tight\")\n", "ax.plot(\n", " r_grid,\n", " fit_potential_v2.circular_velocity(R=r_grid),\n", " marker=\"\",\n", " lw=2,\n", " label=\"v2\",\n", " color=\"tab:blue\",\n", ")\n", "ax.plot(\n", " r_grid,\n", " mwpot_v1.circular_velocity(R=r_grid),\n", " marker=\"\",\n", " label=\"v1\",\n", " linestyle=\"--\",\n", " color=\"tab:orange\",\n", ")\n", "ax.scatter(\n", " eilers[\"R\"][eilers[\"R\"] < 15],\n", " eilers[\"v_c\"][eilers[\"R\"] < 15],\n", " marker=\"o\",\n", " s=8,\n", " color=\"k\",\n", ")\n", "ax.legend(loc=\"best\", fontsize=16)\n", "ax.set_xlabel(\"$R$ [kpc]\")\n", "ax.set_ylabel(\"$v_c$\")" ] }, { "cell_type": "markdown", "id": "ed85e8e6", "metadata": {}, "source": [ "Note that the parameters of this fit were further tweaked (as described in [Hunt et al. 2022](https://ui.adsabs.harvard.edu/abs/2022MNRAS.516L...7H/abstract)) to provide a better match to the vertical phase spiral morphology." ] } ], "metadata": { "jupytext": { "custom_cell_magics": "kql" }, "kernelspec": { "display_name": "gala", "language": "python", "name": "python3" }, "language_info": { "codemirror_mode": { "name": "ipython", "version": 3 }, "file_extension": ".py", "mimetype": "text/x-python", "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", "version": "3.12.13" } }, "nbformat": 4, "nbformat_minor": 5 }