{ "cells": [ { "cell_type": "markdown", "id": "94f8e434-a148-41b7-b0d8-e2bca39db509", "metadata": { "editable": true, "slideshow": { "slide_type": "" }, "tags": [] }, "source": [ "# SSP Models tutorial" ] }, { "cell_type": "code", "execution_count": 6, "id": "2fc9929a-bdf4-4a84-a713-265df77c99b3", "metadata": {}, "outputs": [], "source": [ "from milespy import SSPLibrary\n", "from astropy import units as u\n", "import matplotlib.pyplot as plt\n", "import numpy as np" ] }, { "cell_type": "markdown", "id": "3fb205ec-94c6-43bc-baa4-2e46f180d4a0", "metadata": {}, "source": [ "The level of verbosity of `milespy` can be changed using the `logging` module, setting the level to either `WARNING`, `INFO` or `DEBUG`." ] }, { "cell_type": "code", "execution_count": 7, "id": "21c37703-fa11-42cb-8ca4-70da5c3b0cff", "metadata": {}, "outputs": [], "source": [ "import logging\n", "logger = logging.getLogger(\"milespy\")\n", "# logger.setLevel(logging.DEBUG)" ] }, { "cell_type": "markdown", "id": "3cdffe87-3a4a-47c3-a700-b8e6c41b07fe", "metadata": {}, "source": [ "## Initialize" ] }, { "cell_type": "code", "execution_count": 8, "id": "09340c92-a1b8-457a-b56e-1e4a7b7cf080", "metadata": {}, "outputs": [], "source": [ "miles = SSPLibrary(\n", " source=\"MILES_SSP\",\n", " version=\"9.1\",\n", " imf_type=\"bi\",\n", " isochrone=\"P\",\n", ")" ] }, { "cell_type": "markdown", "id": "73e78767-0828-4cd2-bc1c-f08ee400d553", "metadata": {}, "source": [ "### Generate SSP\n", "\n", "1. Directly from parameters (interpolated)" ] }, { "cell_type": "code", "execution_count": 9, "id": "4bbd6fab-c2f4-40c8-92b5-60879e7833fa", "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAhwAAAGsCAYAAACW3H6UAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjcuNSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/xnp5ZAAAACXBIWXMAAA9hAAAPYQGoP6dpAABp50lEQVR4nO3deXwU5f0H8M/s5iKQBEISINzhCDdyCEQBUVBAvK0HUgtKPbHVelSxtWpbhdb6a621aq2KrQrFelcOEQVEkPsWkfu+IRfk3H1+f4TdzMzOzM7s7uwx+bxfL14ku7Mzz+xmd777PN/n+0hCCAEiIiIiG7li3QAiIiJyPgYcREREZDsGHERERGQ7BhxERERkOwYcREREZDsGHERERGQ7BhxERERkOwYcREREZDsGHERERGQ7BhxERERku5gFHEuWLMGVV16J/Px8SJKEjz76yNbjPfXUU5AkSfGvW7duth6TiIiI6sQs4Dhz5gz69u2Ll156KWrH7NmzJw4fPuz/t3Tp0qgdm4iIqCFLitWBx44di7Fjx+reX1VVhV/96leYOXMmiouL0atXL/zhD3/AiBEjQj5mUlISWrZsGfLjiYiIKDRxm8Nx3333Yfny5Zg1axY2btyIG264AWPGjMH27dtD3uf27duRn5+PgoICTJgwAfv27Ytgi4mIiEiPFA/L00uShA8//BDXXHMNAGDfvn0oKCjAvn37kJ+f799u1KhRGDRoEJ599lnLx5g7dy7Ky8tRWFiIw4cP4+mnn8bBgwexefNmZGRkROpUiIiISEPMhlSMbNq0CR6PB127dlXcXlVVhebNmwMAvv/+e3Tv3t1wP48++iimT58OAIrhmz59+mDw4MFo3749Zs+ejcmTJ0f4DIiIiEguLgOO8vJyuN1urFmzBm63W3FfkyZNAAAFBQXYunWr4X58wYmWpk2bomvXrtixY0f4DSYiIiJDcRlw9OvXDx6PB8eOHcOwYcM0t0lJSQlrWmt5eTl27tyJW2+9NeR9EBERkTkxCzjKy8sVvQu7d+/G+vXrkZ2dja5du2LChAn4yU9+gueffx79+vXD8ePHsXDhQvTp0wfjxo2zfLyHH34YV155Jdq3b49Dhw7hySefhNvtxvjx4yN5WkRERKQhZkmjixYtwsUXXxxw+8SJEzFjxgzU1NTg97//Pf71r3/h4MGDyMnJwZAhQ/D000+jd+/elo938803Y8mSJTh58iRyc3MxdOhQPPPMM+jUqVMkToeIiIgMxMUsFSIiInK2uK3DQURERM7BgIOIiIhsF/WkUa/Xi0OHDiEjIwOSJEX78ERERBQCIQTKysqQn58Pl8t6f0XUA45Dhw6hbdu20T4sERERRcD+/fvRpk0by4+LesDhKyO+f/9+ZGZmRvvwREREFILS0lK0bds25OVAoh5w+IZRMjMzGXAQERElmFDTIZg0SkRERLZjwEFERES2Y8BBREREtmPAQURERLZjwEFERES2Y8BBREREtmPAQURERLZjwEFERES2Y8BBREREtmPAQURERLZjwEFERES2Y8BBREREtmPAQUQKVbUevLZkF7YdKYt1U4jIQRhwEJHCx+sO4Zk5WzH6L0ti3RQichAGHERhEELA6xWxbkZEFVdU1/98ttpgSyIi8xhwEIXhrn+vwSXPL0JljSfWTYmYrEbJ/p9PnmHAQUSRwYCDKEQer8Dn3x3FnpNnsXL3KcV9Z6trMW/zYZytro1R60JX7anvsanxeGPYEiJyEgYcRCGS92qoL8xTP9iEu99ei199uDnazQrq9JlqPPXJFmw+WKJ5f63sXGo9+sNFh0sq0OvJ+Ziz6XDE20hEzsOAgyhEXlF/MS6trFHc9/H6QwCAD9cdjGqbzHjq0y2YsWwPrnhxKUY+vwgdHvsMR0oq/ffLgyePQX7KZf+3BOVVtbj3nbW2tpeInIEBB1GI5NfimtrESRyVT3fdefwMAOCxDzb6b6uR9Wqs3Xcae0+e0dxParLbphYSkRMx4CBSEcJk8CDbzNfbMXPlPjz/+Tb/7cluKZJNi4hkd+DbfveJ+qBC3sPx9Kff4aLnFgU8J3tPnkGSK/7OjYjil6WA46mnnoIkSYp/3bp1s6ttRFH3zoq9GDJtoamiV/Ihlc82Hcb+U2cx9YNNePHLHf7bJSk+LsoHTp/Fv5bvQUW1B0kaQVBNbX2QoZUoWl5Vn/xaUe3BRc8twpHSyoDt4lFljQfvrzmAY2WJ0V4ip0qy+oCePXviiy++qN9BkuVdEMUtX5Ln1A824oN7LzTcVv6d/+vtJzDsj1/Z2LLwXP7C1yitrMXS7Sewbl9xwP3KmSmBPTwVNR5kpNVNlz2dYLU5Xl+6G8/N34Y2zRph6aOX+G+vqvUgNYnDQkTRYnlIJSkpCS1btvT/y8nJsaNdRCE5cPosJr25Eku3nwhrP7Uminl5TQy9xEf/BlBaWddD8fl3RzXvr/Ua93BUVNfPyImTThvTvthad84HTlf4b3vpqx3o9eR8LNp2LFbNImpwLAcc27dvR35+PgoKCjBhwgTs27fPcPuqqiqUlpYq/hHZ5bH3N2HRtuP48esrwtqPmTQOUwFHglycgw2pvLui/n2+YX9xeMeKcW2PrYdL8dz8bajxCDw3f1vwBxBRRFgKOAYPHowZM2Zg3rx5ePnll7F7924MGzYMZWX6493Tpk1DVlaW/1/btm3DbjQ5T2WNB9/sOIHq2vAuRsHyClbtOYW3lu0xnxhqJEITU06UV2HqBxuxYX8xjpdVRWanFsmHUZbvPBlw/7+W7wVQN8Pl7rdDnwb720+/w4DfLcAPRyOzMFx5VS3+t/EQzlTpF1iTl54f8uxCjH3ha//vR0tj83wTNUSWAo6xY8fihhtuQJ8+fTB69GjMmTMHxcXFmD17tu5jpk6dipKSEv+//fv3h91ocp5H39+ICf9cgac/3WLrcW54ZTme/GQLFv1w3HA7SQL+MO97jPnLEpSpamz4mFlCRTIxqPLrDzdj5sr9uPqlb3D+M1/gkw2Hgu84wmpkQyq+qbJyFeeKnH13WLtYmJn1ZA6XVOCNb3ajtLIWz87ZGnB/ZY0Hn2w4hNMWyqlP/WAT7nt3Hf4w73u8ungn7vzXavxv4yFFUTaPLLhUB6TVtc4pSU8U78KaFtu0aVN07doVO3bs0N0mNTUVmZmZin9Ear5CWe+sMB6i8/F4BT5efzDg4mR2BGO3xkVV7eVFO/H9kTLM23xE835hoouj0sQFbZvq2/4zn30X9DGRJgTwyuKduGDaQt1tFv9wXDfJ0qPTY1RaWYNff7QJK3adxMjnF/tv1+qRePrTLfj5zHW4/a1Vptv96bng7F/L92La3O/x+XdHcd+769DtiXl4feluAEBVjX6vmcPW3SOKa2EFHOXl5di5cydatWoVqfYQmfLJhoO4f9Z6XPz8IsXtZnMmgl1n5GW/9UZfzFyshAAunP4lHpq9QXOBt+NlVYoaGIBxOfFQGFULlZs+93scKtEfklq1+xRcOk/wK4t2Yufx8oDbZ67Yh7e/3Yeb/vEtzsoTT1WhYWWNBzNX1vV+as2iCcXv/vcdzlbXYvuxwHbVt8N+a/aextvf7o3MMB5RArMUcDz88MNYvHgx9uzZg2XLluHaa6+F2+3G+PHj7WofkaZPN9St31F8VjncYTSEIS8/HuzDX3GN1tml2QvIweIKvL/2AP759a6A+85/5ouA2yJ5Wdp94gz6Pv15RPblEQJ6tb6eX/CDogfDZ8/Js9oPUO1nh0FQAADVtV787cvt2HRAe0hHT7izlcJVUlGD619ehl9/tBkLdGYIETUUlgKOAwcOYPz48SgsLMSNN96I5s2b49tvv0Vubq5d7SPSpPdN26iH49nP6vMGrHzZ1MtPsPqFda/exVclkhU8/7V8j6JoVzg8XgF3kLbVerx46pMt/mGo1CRzHzHBnst/Ld+DP33+A67821JT+SI+ycGOb1MXh8crUFXrwbp9p/23rVCtKEzU0Fiq2jVr1iy72kFkSSjX5PWy6Zy1XoE/zvsegwua46KuuSitrNEdypAPBciZmRYrpxckqfkCjqOllais8aBddjomvbkKHq/AvycPslS99M1v9lhqoxGPVwQdRvpw3UHMWLYHM5btwZ7p45CarH3BV5+B/LnMzUjFve+swYiuebjx/LpZbVsP1+e5FDw+B7cMbodnr+0dtM1r9pw2vN+uIZWfvrUKX207jn7tmvpvK6nQTj4maii4lgolJL1rrvpiXF3rxb3vrMHb3+5FlWzK7dzNh/H3RTsx8Y2V8HgF+jz1Ofr/boHmPvVyIKz2cJiNE5LOrXVy9d++wUXPLcLWw2VY/MNxLN1xAsfLzU/jPFxSEXwjCzxegXdW7DXc5qhqFohekqn6uZAnnR4vq8KcTUfwy/frF5RT9/q8u2KfqXoef/tKP6EdANKS3ajxePHOir3YcyJ4IrEZxWer8dW2ullQ8nyUClkOz96TZ7BmL3s8qGFhXXJKSGammwLAB2sPYM6mI5iz6Qjys9L8t58sr5/dsiTIFNmDxdoXbqs9HKYDjnMXV98UzjWybnmzvSRA5Few9XgFFm0LNp24vn3Hy6qQnqITcJx7/b7ZcQLtstN182G+2XECH647qKh06nO4OLS1US7qmouxvVrisQ82obLGg38t34vf/a9uZtCe6eOCPr6i2oMfjpahT5usgADX6xW6M61qar349/I9eGXxLv/f1KKHR6BDTuOQzoMo0TDgoITkkvXNHS2txHeHSzGia25AGFIs68aW93A0Sa3/0w/WazBj2R78ckwh3C7J/439WFklVgfprldTX1P1ek48QijyFORVQM3OONl6uDRixbV8zJR7l7fv/Ge+QG5GquZ2y3edxLIdJzDhn3UVYWffVaS5ne9+Lb8PcfqwVwgM61qXd1a3voxxEKX21Cdb8J/V+/HH6/sAUl2Pxp3DOwEAvvz+mG710hqPF098rKwz88PRMqSnuHHfzHWYMLgdrj6vdQhnRJQYGHBQ3Kmu9SIlSLKfvIdj8LN1tSP+dku/gF4E+UVenouR1SjZ/7Pet3C5295chU0HS/C/nw1FQW4TDHpGv16FHnXAoTckUFTQXHFxl29nZhhBCKGophlN/7fgB8XvRpVTb5EFE2YDKbmlO0KfgZIm+/v6StZrc6K8CjlNtIMkoC4p9j+r66bvPr9gm79SaVFBDnq3ycLKPfrDJOs0SsK7XRL+9Pk2rNx9Cit3n2LAQY7GHA6KO74PdLnis9U4Wy2bbaExsvDV98cDAg75sId8DF1+u1FhKJ8Vu0/hbLVHsfS8VepCYdU6wYMkKRdTk/fMmKnRcfqsPcmJM1eaK8oWip/NtF4uvXur0IoIeoVAI50g890gheee+Hiz/+czVfV/T77clX8sCZz67KOVNDr5rdWYvfqA4TGJnIIBB8WdU+XK6qFllTU477cLcN5v65M6tTIZhBABuR16uQGr99YPh1hJxLSat6F2sLgCVeeqj9borBsjhHL4okpWrdQXiBSfrcYlf1qEGd/sDni8XvJrPDtRbn3J+zV7rQ1p+Xi9QJpOMqve61tV68GyHSf8xckAKAJgM38VrPtFDR0DDooL8kkIbtVfpa9SZHWt1x9AaE0N9Qqh6OF46asdphbnmj73e9PtDOeiMXv1AVw4/Uvc/I9vASgXTJPzCqHoyajR+PnB2Ruw68QZPPVp9MugJzoBAZfFedVPfbJFMQQEKIvDLdx6FB+uS/yeijeW7kaHxz7DV98fwy6NyrFE4WAOB8WFJJfLP8SgvhikyCKQY2VVaJGZplmHQ339tmPp8Uh8SfVNldTLx/CKulwBH3nhLl8gsmjbsfo2CWGpNkdD5wsaB7ZvpujpAuqmycpV1Xqwfl+xomdDy6xV+zFrVeIvTPnbc7N1bptRt57NjNvOx4jCvFg2iRyEPRwUF+QVLN2qi6e8WuXgZxfiqU+2oLwysHrm4eIKw8myaTpFqKyI5HoYVTpDKl4hFEMqZbJz9a3qOqZXS/9tZmaPUKBerbMCbvP9XQkhUHy2Gle+uBQ3neuRsmLWnUNCatPGA8UhPc4uM5btiXUTyEEYcFDMfPn9Ufz0rVU4UV6lKOzkdknweAX2nzoLIQTGv6bsyp6xbA8Wfn9MvTus3nsaGwzW2qg0kRwaTCQv7Xo9HEIok0PLZGvA+G7vnJcRdD+kzRfPapWQ9z3X//x6N/r9bgF+OBrasELfNk1DetxVf/smpMeFwusVeP7zbVgsq0OjDsoXbTuOEpuSkKnh4ZAKxcztM1YDAK79+zeKC/kbS3dj5/FyzFy5H9f3b4MTFpI6bRfBiKPaoIejRjZLRV5s6+H3NuCCTs3RUlbETC8XxKyiguZYvutkWPtIJEMKmgPQXmel/NzMk2fmbA24zwqza8jE0rKdJ/2zrrY8PRrvrNirGZR/vOEgflLUIcqtIyeK/3cFOd7+UxWKPIVDJZX+MfP318ZXIl6wVU3NEkLgTLX2omofrz+EUp11N/adOotZq/YrpljWhtnDEXSBswh7727tIl/h+vDeCwzvf/eng/HQpV1x74jOALR7ON5feyAiPUYul4RxvVuFvR87bZAN3/z+s614do528vTpM+zhoMhgwEFkwbYIVe/81UebA9YdkXtUto6IFnklUl8Ox7p9p3HJ84sstyVFPS3IZnpTUkN1Zd98zL6rCHmZabrbFOQ0xgWdc/CzkV38ReXWaxTiAupmN4WjV+u6+iDP3dAnrP3YTZ5UbVRjpdqjvXghkVUMOIhi4N0V+zBN5xslgKC5A/LFznxDM5PfWo1dx60vQBbtIatIJO/KDWzfDIM6ZgckG8u9Pun8gNv0ejL+8sX2sNqTnlw3Up2eEr0R61NnqjFv82FFb1fJ2Ro8//k27JRNb91/6izeWbHXUi9OuEN2RD7M4SCKkWMGZb+DkS/M5uvh0Jq5Y4beN327BCtbb5VvGrXboLZGh+bpAbcZbR8OvSqmZpVV1mDjgRIcL6vCNf3MlTqf/NYqrNtXjPtHdkHHnMYYXJCNZ+d8j083HMKLX+7Ad78djfSUJNz8j29xsLgCS7ebLwuvl2tEZBUDDiKLvHEwDVX+DdX/rTYBSnH0bp2FdtmBF/9w+Ho2jAIIrTolVlbetcLM2jxGDpdU+het65mfiX2nziKrUTIGdsjWfYyvtssLC+t6Z5qmJyNJtsJhj9/Mx7bfj/GvUjt38xHT7ZmxbA9uv7Aj2mkEbURWcEiFyKLK2tiPacvXYfEXTEuAgOP1iQMjXqTMd95GQyraj7PnCUvSyInpYWHdlzOyBOpvdpzA5LdW40evLEdljfm/u+KzNahW/Z1Oecf6ejU+E17XrkVSFQfvhWjYsL8Y3Z+Yhz+rFieMhflbjigK/yUSBhxEFkWinke4qjUWdAv1AvrgpV0j0iYzstKTg29kkX9IxW3t/OU9IveO6BSx9iRrRH5eIVCQ29jU4w+X1CcTy0vXnzxjbb2ZUtUQ2xdbQ79I7T9VEXDbruPl6PPU5/htAyivf/VL36CixoMXFm5XBITRVnK2Bnf9ew0mvbkqpu0IFQMOIovO6kxnjSbFkMq5mh2hfl8vbJkRfKMI8fVC/DeCU2P9QyoWA668jPpl6JNcEn49rntE2pOsM+tn9l1F6GQi6DhUHHhxB5QF4GLhZHmVotLuja8uR1WtF29oLCDoZKfPWl9oMFJ2nqhPAD4ZwoKHscYcDooJTxzkQQB1Myas9licrY59N7J8SKUmzB6OZIs9A+Hw9SrkNEkNsqX1fbosfn16ZHShf/0TAe08DzO6tczA90fqp0snJwXuR4i6cx7eNRc7g8wkeu1r7SXuy0JMCo6Uq1/6BgdO1wVDvmrADZHekgTRsOVgfSVleXHARMEeDoqaPSfO+L+lbYjQmhFPXNEjrMdfazALQP4NWC4eujLls1RqPF7sOFaOshDbpfWN/K6LCkJumxHfRd3stX2QQaJk/T7r/k+yGHE0lwU9QljrIfrzTX39PzdOVX5vS3GHlzSqt8JxWWUN5m0+gmlzt8YkcdkXbADx84UhFmK5lID8y9FH6w76f46HRHYzGHBQVOw4VoYRf1qEC6Z/CSByxabC/XJu9K32/I7aF7tYf9MElD0ctR6By/68OOR9aQUcdw4LL+Bo06yR4f2ZacFzOWbdOQQPXNol6Hb+Hg7ZS9kz33ySJnBuyXqTf0u/urw7uraoH4ZS/y1r9RhFYgru2WoP7n57DV5dvAufbDgU9v7sVlHtSZgLoRW1MaxLIh9We/HLHdiwvxgdHvsMBY/PwcKtR2PWLrMYcFBUfH1u3r/vYh2pWgxul4Snr+oZ+uN1Ao7xg9rpPiYeejjkSaM1Hi/C+VzXCjjCncER7OHNGqfgb7f0Q+um+oGJ2yVBMtHv4Pb3mtRv+6cb+uptrsulExT8YlRXPDK60P/7HcMLFAGE+m+5aXpKwD6SzgUhZhYb1iuMdraqfihv+7EyLN1+AuP++jV++tbq4DuNom93ncTgZ79A99/Mw49eWYbf/+87zDKoZJpoqmPYw6FOBP7ZzHX+n3/3v/hP3mUOB0VFWrKymzlSq7y7XBKKOjYP/fE617Np1/XGlHe1pxGGOnQRSfJu3Tv/vSasfWn1NuldfM0yE7Bc0Scfa/cW6yYduqTA12fu/cMwc+U+/Gv5Xs223ndxZxwrq0Q3i4mwRkMq913SGRKAkooa9G/XzL+9jzxgG1KQjYkXtPf/3qZZIxw4XYGxvcyvq6KXUyRfe+fbXafw0lc7AQBbDpWa3rddhBD+gO/mf9RPoV27rxhrz9UI6dOmKXpY7HmKhpKzNfjj/O9xXf82GNC+mf92IQRW7z2N3q2zFNvHsoejVJU4vO/UWf/PpxNgVV8GHBQV6otapMaAXZIU1vCM4YVVp4nx0MMRyXHkJI0hgHBHAMz2kBgfRwpIBO3eKhO/vboX2jdv7P9GJz/Ww7KeCCuMkkZ9vRmPX14/i0VeEyMjrf5jdNadytk3H025EKt2n8KlPVoAMO75Gdi+GVbvPa17v/zvbo3BdrFQ4xFI0UiWlTtlcVpvtPxh/vd4d8U+vLNiH24Y0AY/H9kFbbPT8eX3xzD5rdXo27apYvtY5XCUVtbgg7UHde+3q3JuJHFIhaJCfeHwRqiL43hZleasADNm3TnE8MIodCKOb3aYLwttl0iWm7ZjSAUABuvkwMgZVwdVBiSDZPuT50lEIh1ICGvnnNWoPgfFKODNaZKKsb1b+YuBGf3ZtzBYfA6I3JomI7vlRWQ/cmaK4emtjhxrs8/NVAKA99YcwI9fr6vyOmPZHgB1Rb/koh1wVNV6IITA87LF9rS4pLop++v3FyumL8cTBhwUFeqx+Ej1cBwvq9KtexDMkILmmhe8nCaBY/ByVgsoqbtkI6E6gt26Wtf8sHM4ADx7Xe/g2xkcR1LdP3VsN//P8tctUhVDreymILcJ/vijPnjr9kGWHhfO8SP1nunbtqnlIadgfD0+t567WGupqvXGrHfA4xX4dMMhFKtqaNR4vP61iHz2njyLdftO6/5tRmsxO49X4KdvrUbhr+eh49Q5eEs2jKjFK4A7/rUa17z0TdwmFTPgIFuVVtagvKo24MM0Uj0cWY2SQw44AOWHvC/577r+bQK2626hNLXRMSIlkiu8al2wrda00JJkoovX6KWTJEnRNnlyZrKsgcGmwxolpvrozVIxWvflxoFtcVHXXNx8LsF4oGz8PxTBAqe/fbUjrP37pEZ48TwAqDqXd/K1waJwv/l4M7r8ai7mbzG/jkukvL/2AH42cx1+9Mpyxe1HZFVd5a79+zKs2HVS8z47g6byqlp8tvEwzlbX4uvtx/GFhZknHXMa45sddW1+fWl8FmNjDgfZZtfxclzyfN10TfkHtxAiYgHH7UM7hpXDIZ+l8vbkwTh1pgoXn+tynjy0I+ZsOoJR3VvgnxMHYsxfligKPJll15odkaIZcITZ5v7tm5naR7BZKPIgQN6rIc870ZvVMbxrLpb8cBy/v7ZX0Ha0y04P+Eb71cMjkN/UeJgDAM5r2xTfTh2J5kF6xoKJ1hC8HQHHpxsP4d4RnQ23KT6X1HjXv9dgz/RxqKr1wC1JmmvPRNq/z/UO7DhWjuc/34YHRnWF2yUpysir6RX4sjPg+PWHm/DR+kO4cWAbXGJx6KvGE7jcQbxhwEGWVNZ48M6KfbikWx465hiXaX53Rf1UOHk29aPvb8SPBrSNSHuyG6fUr5YaAvlFMbtxsiJPYED7bKz59Sg0OzfNMdSLcLzncmmdVrgBx2+u7BF2cq2kaof8Z3nwoZ4B5fPGxIE4XFKJtga9FO/8dDC+2XECNw1siw/XKRPyWmWlITXJXBGvllnBA5NgohWYpuo8X1Zd17+1P4nxj/O2BQ045A4WV+CqF5eirLIWm56+zPTzHCr53+KLX+7AruNn8NmmwyHty86iZx+trxsKmb36AEb3bGnpsfIKyJH6QhdpHFIhS178cjt+97/vcPGfFgXdtkJndcvZqw9E9E0bTna2K0guQPMmqf5tQr0exHsPh3bAEd4+M9OSw76IqNslfx7lw2h639iT3C7DYAMALuycg1+O6YYktyvgdQpnqC4UkV5FV49ej5BVVr+Byw37w5c4eaYa1R4v/rvmQETao6e0sga7TijLyYcabADRq7JqdcmFClnAEa+VYBlwkCWr9mhPx/N6BcpV32iNuh4jmUUdzge1/MIaLBcg9B6O+A44tNoXiYtfduMU/HJMaNNUAaBLXoaqh6P+PnmQGalvx/KX3yVFf5ph9IZUIvN8qS9qVt7T8oceLtYf1lDvf+vhUtNDGjUeL+ZtPoy+T39uul1mROti/uKX2y1tL5+q7WEPBzmB3kJfk99ahV5Pzsfek8YLU/nEyxtCnsMRLFHS6IJw65D2uvfFebxha0B074jOIa3C+sWDw9Eoxa14TeQBgPzvMDVC39jl+SSxCBLtPKY8oTVSORzqPIFQL8Q1Hi9KztZg5sp9KDEoXnX322sw9oWv8dj7mwLu++fXu3D/rHWKNkx+azXufnutqSKDuTrrJmlRz2qJhNNnqnGjKqHVar7YSVmdkzj5eA3AgIMscetclb/adhwA/KtvBhPpbwk/GhA4s8QM+ZBKsB4Oo2/9w7rk6B8jDiOOd3462P+z3d+sby1qj59d0hnv32N+SfomqXV1LvRyOOQidQGV7z4WL5mdx5TX+IhUD0crVd7KWZ0h1GC2HS3DPe+swdQPNuGB/6zT3W7+lroZG++vrR+C8XoFHnlvA37/2VZ8vP4Qlu44ASEEjpZWYskPx023wcrqxXb0cIz+yxKs3HMqYvtjDgc5QrCpjvI/dKMZCJU6H04Xdg6tTPmUi80nrOkJp4fDKKiIt3ijUbJbkWhpd+5AapIbD11WiAHtgxcCU5M/5/LgUN6rbkfCoZk1XMK14TeX4Zrz8m0/DgA0Tq1/joxyODJSg88jeGR0IX5zRQ8UdVK+V416J4ws2nYcy3bWTef0fXExa/Xe03hPlgNSXevFa1/vwuBnF1raT/PG5mcY2RFwHCsznuZ+8/n1Sfaje7ZA57wmhtszh4McIWjAIftDN7qO6U0523SgJKR2hfotPVNWMTJYD8d3h/XXrDA613j7siFJsf82r6WHRq0TSSeHw+Ot//uJXA+HJPs5IrsMIM9zyEpPRrrsAq/3d6I3jGn2fkD5fjMK0FY/MSrovsb1boXbh3YMCFRLKuxZy+PVxTvR4bHP/L9nN07Bsh0n8O/le7BfNvsNqCt89eyc7y0fIzuEgGPbkTLdL05WnA1SgXVU9zxMv76P//dGyW5c3qt+BsugDoGBfLyu0suAgyzRWndDzuzfuV5p7kdl1SSD+b8b61cEDXXYoovsm4LeyrE+Rlnj8Ths4jP7riJFtVMJyoXK1G0f1yf4QmMdmhvP/gjFzDuGBNymmAor+1neQxPuQnP1x6r/WS8gjjR5y/XymgYHWZzwsbHdcWXffMOKtvK1V4xyXsz0Ful9BkQq4Ph6+3EcKan072/aXGUA4ZIk3PLPFXji4y146L0NETmmlYCj1iuwbMcJjP7LEtz3rv4QkBn/XXMAPZ+cb7iN+otQRlqyIiJ+fdLAgBWL4zTeYMBB1gTrBTA7dqheZtnnPNVCSXp2PXu5oiKo0UUnM02/m1jeY+M28U1Rl8FDYx2LDOqYjSev7OH/Xf3NVP3UvXRL/6D7/PRnQ/H+PRdEpH1G5E2Tt/vCzjkY1DEbN4SYu6N9LPtfKKPhK733TrC/nwHtm+HF8f3w3A19dLdpml7fk6fubm+ZmYa/3HQevnhwuPGBztGbLhypgOPW11diyLSF6Pv051itmddg/Woa7G9aHnDI3ytaPF4v/vH1LgCwVAlUy18Xbg/aA5p8Lpj46/h+GNo5Bw+M6qJ4z2akJSM9RRkoxktSvhoDDrIk2JCK2W/60+duDbjN7TK/8qs6wNDrnZh2XW+0bqb/bVz+2WumFLdue4wWgZMdw6jM9tDO+omn4VIMoUD5kR1KDkdGWrJiKW+7yD845X8byW4XZt9VhOdu6Kv1sJBEIzBUTx2VHzPUa4Tv79aoumQT2dBNnsaMjGv6tUbnPO31VeY9MEzxu950YaOAQ73iqlnqUuSAcjaGWUMKslFUoN1TNLGovWJY7oo++dgzfZzuvjxe/ZLoVlTWeBQFEfX4hsyu6puPt386GM2bpAYEx+rXPl6LDTLgIEuC1SZQfiPV305vASQztQ+aaCS26XW8jB/UTvPN1+zcN76C3PpqqXoVK80w+/42Or8Mg56YCYPbWWyRkqSOOGTi5cNJ3tWffi7JUf5tPNRVgc2KxtOQ3Vh/NoReD0ewIN73NyV/rhY/MgLvymYiPXllTwDADQPaoHmTVFzZtz5ZNdgwabeWmWgvG0JL1nmzqSu12iWUwKx5k1T8fYJ2L8eDlxUqnuNgXzyqPZ6QljhQe+vcarTB6D3fRqI1JGgVS5uTJcE+nHw9D09+vBkzV5qbIqt4vImvmfLuYTOPk983rEsOLuycg3G96/IUcpqkYtHDI9DE4GJvhtmeHaMPs+v6t4EQwDyNxa3CrXopb1+ZajjLTA/H4I7ZWLE7ctP2tKQlu/HGpIHweOsqlQLKb27BhvPCFY0ejjuGd8T3R0ox5lzSn/ybqt64e7B2+f6meuZnolvLDOQ3bYT2zRujffPGeOKKHsjLSEX3VpnYPe1y/2t91/ACPHduuXMzf1vyYEbvM2Cl0d9HHHTxN2ucgs55TbDjWLnidvWsHXXvae/WWXj88u6YuXIfPtlwCCt2ReZ9sMrkNFitQLvao0xWVScORyKZ1Q7s4SBLgl2cJKmuiE2wpZT1hFrd0SjhU77LzEbJuPuiToqS1x1yGluah6/F7MXKKNck2S3hlVsHoEBjjZpwK7OqDyvfnZmn/CWdb4dyOSEsXiZU4/GXdGuBS3u08P+u6OEIJ8fGFPsjjvSUJLz84wG4+rzWAffp5nAE2adv8bMktwtz7x+G1ycO9N83eWhHf2+G/L0rXzDNzFCivGnygON31xgvjNcsPRnTruutGBp76NKuQY9nl0fH1CWl92pdPyMqxe1SvH/Vz8fI7nko6tQcLTLrPiNWyxJwwykTb/azTisgvObc349v7Sf1AniVNd64nKnCgIMsCTaTw6VaUtwqO9Yr6day/sPloi65oR0gCMNpsbKLqtHz57sgaH1MhPvRYbYHSE+wgCw3IxXXh5DAGaxqY4ec+sDQ7nohsR5aCjWolF8gJUky/Txdfy7p+oFR+gFA93NTleXBkLyn6dYh7Q1neMy5fxjGD2qHgedqsEgSMOnCDqbaZ4dLe7TA6l+Pwgf3XIiCnMYY2jkn4PnSCwS0Xp5wVqo+ZTIXRSvg6NIiAyt/NdI/bJas0eZ4HFbhkApZEuxD2aWec2mRmahf6/NUMnjfPz6uO0oqatAxt3HIFUmDMRtkqc8vPcXtX+XR6NTD7ZEObJ65einBNEp2o6LGg4cv64qdx82VtZcL9i0sIy0ZK3810vbVRIHoLZ6mR+81DvbSh9or+NyP+uDh0V3RKks7kfnCzs3x+sTzASgDDvXxjC6cvgD7kdGFyGmSgjG9WgVM4QTgzxHZezJ4EqUWK0N+vuB5wYMXab7n1D0cvlPXWsNFLxctmP2nzuquS6WmF4jmZdRXelX3cAB1i2c2SrH/fWMFezjIEq0hAfkb0SVJYXX/6/UAdAlSWc+o5yCrUTJeuXUAHh3TLWI1G9Tke33lxwPw6q0DNLdTj39PldUd8Y3naz1/Wt3t8qGHYNQBUefcDCS7JTRLD29V10WPjMA/bh2AGwa0NZwloaeRiUTdvIw0ZDUKzNuJtFjnzuoNqTQNcu6hzq5yuSTdYAOoWzzPl0htcr20AL7gpHFqEu67pAs65zXRTIKc/8BwzPn5sIDbjSx99GIMaN8M06/rjWv7BQ5R+ehNtXe76nuD5MGmOqDyvSrVGn/fZhaSE0Jg/6mzmLf5CL4/Ulc8cOTziw0fI69JY2ZkRB7ENW+cgpwmqbq1jmKJPRxkidY3+VWybxaSJIX1bVzvW6Y8UVSrXkKsC2/JA5lerTPRRmcqrrqd8kqnRqeg9aHTTCN5Vrd9qp1npSfjm8cuQXpKeB8BLTLTcFnPugTIQR2z8cY3uy09fvKwAizbeRJXRanEt5FY/AnJh9v03jfX9GuNtGS37jpFWt9uI0Fe7jvULxFavS9aQX+K2wXIYs8OzdOxx6C3454RndCmWbq/Fszs1YHPzZSLO6Fvm6YYrDMdVk/AZ9C5c9cKLmq9Al6v0DynutVtyzB/yxG8sLBu5dcueU1w48C2qA4SqAwpqK8eaqa2kfwL15onLg26faww4CBL5O8rIQQ8XqEo+S0h9IWDhBC63cPyC24HjaRKmycwBCVvtfrirki4U9cPUYy/6+//J0XtMXPlPtUxzV8htZ5WeZdsJIzu2QKv/WQgurfSruegJatRMv4bhQJiZsS6QJuetGQ3pl/fB4dLKrFYY0GyUIdU9Px1fD/M33wEk4d19N8W6nvabI+iervmTVLx0oT+GPfXpZrbq3MntL5wtMpq5A+GgzHTSr3ejBqvF6muwJ66uZuP4N531ipu236sHM/MUdYgKshpjF0nlMOR8qDHzFPfIz8T246GP1XXbhxSIUvkb2yvAB7570b8/rP6N9ALC7djyXZrCzDJGQ2NvH/PBbj6vHw896PAioqx7uFQrvlhkKCp+mBVJPz5hlQ0Hte9VSY2/OYyvDnpfNkxQ2ufXSRJwqU9Wuj27sS7aFQaVZNfTELtGAynYJ2Wq/rm46UJ/RW9X6EuBqbXtmnX9TZsd61XoGd+fan21k0bYWS3PP/v6vLsWru6xmCYRc3o7Jqm1/X06AUcvqHEpdtPYNG2Y/7b1cGGnvPaNcUTV/TwJ9/eP7KL4n71NHYtT1zRAzef3xYf3Bsfwbse9nCQJcoVO4VmsZ9f/GdDyPvXS/4UQmBA+2a61S1jH3Bo/6ym/mCUZ/wHO4Ws9GSMKMzFI6ML0TM/E/M16nWYPa4ZM247H7/6cDP+FMFqnvEs1j0cet9km6UbTzeOdA+HFqNv2S9P6I97dC6ueu/L8YPaobyyNuDbvo86mfjnIzvjpvPboc9T81FaWYtLZMEHoHztbhrYFveP6qJZIFCP1pDRn2/qi4Vbj+GWc0X3qmu1n4Raj0BJRQ0mvbkStV6BLx+6CAW5xjlncr8Y1RVts9MxeWhHCCH8Xw6Gds7B0h0ncP2A4IFTduMUxQJv8YoBB1ki/2yTLwgVKcGm3eqJ9ZRG+QerccCh6uFwy3s4Aj19VU/0biNbeE2SMOXizgCA+VvMr+MQSkA2ojAP3zx2ieXHJapYzFIRBr/9ckwhKqo9KGwZOER1y+B2eHdF3RBbuEXhzDBam8NXC0KLUTBkNEzj6014ZHQhlu084a9bsvChEThWVqmY6g4o/77/oNEDGoxWU67t1wbX9quf1abXw1FeXYt/Ldvjn+L9/ZEy0wFHq6w0RU0g+d/gPycOxIHTFUGXok8kDDjIEnlAMP61byO7b5cU8rc1SZIwqEM2TpRXBYyHRsqQgmx8q1Nl0DCHw+A+eQ+Hr/dI/uE38YIOuu2xNqRiftuGKp6eoj9e3wc3nt9W9/7usiAkGsG2UXCg/pt+feJATH5rNQDjLxBGNVh8h5tycWd/gA3U1XvJ1VgL5rIeLVGQ2xj922n3gAZjJkdleNdczRyaF774AbNXH/D/XlpRY7roltFQVVqy21HBBsAcDrLIzm+BLkkK68L4n7uGYP4vzK14GYrbL+yoe5+83Ua9Ceq7gvVwGLGyfSSGnHxVFe1Ymj4exGSWitD+WSvYUFwUFT1q9jfc6PopH2bt1jJDMZRhlDRqNG1TXYE2mEYpbix88KKQh//yMoNXGr5NJ/iXBxsAUFxRE3QWis+FNi7YGI/Yw0GW2Jkr4ZKkkIdUgLoPXjvLX5v9YDfaKlVV9EjeXt/uzX7YRruH48N7L8TfF+3EgzEsTW2nWCSNWkkVVda7saMt+oymxcp7Ja1MZmlvELiGkqMaTuB1Vd/W2HSg1HB4yOWS8Otx3RVJ8lqmz/0ery81Nz38qXOL6jUUDDgobhgNqTQPc60TO3RrmQFJkvDIaOUF2Cgoy0hLRu/WWdh0sERj28AhFSPWpsWGf4Xq3ioTL47vF/Z+4lUs8oCszFKR9whEOzgy6vpXTJWHMH3hv/q81jhcUql5kQ91Gm6o3C4Jv7myR9DtfjqsIGjAAQDHy6pMHTfLQi0dJ+CQCllitavTiiEF2bofVr+9OvbfBNQt694qE3PvH4ZLuikrfgbMtFE9ZbcP7VC/rWI6rcX2WNg+1rN4EkIMniJ5wmCw4lq/ONezNH5QWxR1qitm1SorsrVU9BgOqaimyvdv1xS9W2f5F43T43bVJUCf30GjVyH+1h3za5utX52VjLGHgyyx64tHUUFz/PFH2uOvPx7SzrAEc7QYLtCmWH3V+MrllQ3vyoMMq13C1nI4LO26QYrFkMrkoR1xvKwKo7q3wOtLdxluO6xLLtY9cSmapidDkiR8+dBFcdHzpww4BJLcLnxy34VhDXFEu4fDik/vG4rvj5ThJ2+sjMvy4fGMAQdZYtfHwMw7h2je/t7dRdrfgGJA/fmp93FqdHGXoPwwVUynPfe/6SEVCx/osV6YLBHE4ilKS3bjqavqeu/+GSTgAIBmsnLjVmo92Ek+DOqbnRHu31v8hht1hcCGFDRnsBECDqmQNVH+5pEWhVVCzTL7DTjYdnpPoZ0XPPZwBMenKDSR/NvyLdJ3QSdr659QYmDAQZqOlVbi5n8sx2cbDytuj/Y3j7j6Ym6yLQFrP6meNb0eDqt5FszhiKxY9wLF8SiCaeE+h//72VA8fnk3/Gpc8ATOWJtkUCOHtDHgcIiXvtqBPy/4IaTH1ni82HywRFGs5qlPt+DbXacw5V1lyeJ4Hlu1W8BHqSI7v57Wxd2XQDd5WEdFAl44i85Fe5aK07EXKDSRDNTaZqfjzuGdLJUljxV1GzvlBi4qSUoMOBygssaD5+ZvwwsLt+NYaaXlx//yvxtxxYtL8fLinf7b9p3SXho62vFGPF0nzX6wqi9ckiThrzefh01PXYae+Vn+GQbJbkkRNPh2f02/uuCkV2tl+ebA9phsOPTXqKF6sf5ba7ihfGKSV0p94ooemH1XUQxbkxjiP4ykoORBwJlqj+XH+xZge+mrHZgwuB22HyvXTYiK9odiPH8z1+th8AUmv726J/725Q48e20vSJKEjLS68emOOY3x1cMjkN04RREg+vZ3/8iu6NumKQZ3NB7HjnalUeeL9ZAKQ45E4pFNN5t0QQecrQ6+qmtDF9b3nunTp0OSJDzwwAMRag6FQr6wUq3Jkrpa3JKEkc8vxg2vLMcPR8s1t4n2Z2I8XSjNtsTXw/GTog5Y8fhIdM4LXHyrY05jZDVKVpZEP/duTEly4bKeLYMWBbLUw2F+0wYrjv7UKAG0lE3Vd7skpMZRgnu8CjngWLVqFV599VX06RP/S+I6nTyvwmwNfy2SBJw8U224jZ2Fv7TE07i6+oI0vKv2OgjyoZdgwzCKbS2GBVbGzpNsLPnuFLF+htjBkViu7NMKhS0y8KvLuwOo+6Jg5JdjCht8nkdIAUd5eTkmTJiA1157Dc2ahbY6H0WOPNkznA8to4WW6g8Q+v5DEeuZA3LygGBYlxxcJaukGGp3uNll7bXbY15KFJYwT3Tx1JtG8S8vMw3zfzEcdwwv8N/WwmARuHtHdMbCh0Zg6LkF23LioGhbtIX0KTRlyhSMGzcOo0aNCrptVVUVSktLFf8odF6vQIUqT0M+6yGsgEPnA7esskZ2LPZwAMDFhXkRCYYCV1IJ8cHnFOh8g4qnwC1exfopmjy0IwDgkm55sW0IhUxdAv68tk2Rm5GK134y0H/bn286D/eM6IT372l4SaaWk0ZnzZqFtWvXYtWqVaa2nzZtGp5++mnLDSNt17+yDOv2FWPtE5ci+1zVQfnCSuEEBHoBR++nPseWp0ejcWpSDGapxM+FUhEcqJrVvHFo31ZcFoZfAtsTuP2MSYPw2te7sGL3Sd08HNIWm9Vi613cLQ/LHrsELTKjsz4KRZ7887H+M1O5oF1uRioeHdMtBq2LPUs9HPv378f999+Pd955B2lp5t4UU6dORUlJif/f/v37Q2oo1Vm3rxgAsHDrUf9t8u78WlWIvflgCbYcKjG1b6Pr3eZzq5tGMt7wVRU0EleZ+5LmjwCAds3T8dyP+uD1iQNhhfw5tzykorF9u+bp+N01vdAuW3/pb9IWD7FtftNGuismJ4LEbXlkXNGnFYC6Ke+Nz9XpiKcvTbFmqYdjzZo1OHbsGPr37++/zePxYMmSJfjb3/6GqqoquN3KTN3U1FSkpja8sSq7yS/D8lkq8h6OimoPrnhxKQBg01OX4bWvd+Piwlz0a6edd5Ns8EHnS0aN5PW/TbNG+GjKhchOT9HdJo7CDVXNjMDn6oaBba3v0yCICd6e0O8lskM8vV9j4afDOqKi2oPbZCtCUz1LPRwjR47Epk2bsH79ev+/gQMHYsKECVi/fn1AsEE2kr2z5UMqtZ76n0sq6nMvXvpqJ/66cDuu/fsyAMDiH47j4fc2oLzK2tzxSM5S8XhF3fRQg+mfcdXBEUZvhJ6whlQYU0QUn0994we1AwDcdmGH2DYkzrVplo4//KgPurU0LtrXUFnq4cjIyECvXr0UtzVu3BjNmzcPuJ3sJb/wC53go0Y2RXbrYWWy7sQ3VgIAmstWnzxUol+l1HeMSAYAdiWg3juiE/6+aCceG2vfOGmkrk1aq8Wab4P+I3jxtI6zVPT99uqeuGFgG/RpnRXrplACY6XRBKUXZNTKqt/JA44anfocry4JviS28riRCxLU+SY6R7S830dGF+Lm89uhbXaj4BtboLgcRejiJN+N1bF7oyY8fnl3rNpzCpMv7Oi/bVDHbKzcfcpqExsMxhv6kt0u9NcZipXjU0hGwg44Fi1aFIFmkFXyy7C8p8DjFVjyw3HM23IEPxrQxn+7fKglUscNl0cn4PjFqK748xd1C9GFEt9IkoR2zSOfNCmF0Ruhv8/6ny0HHAb3dcxpjLW/vlRRW4UXA2OxnqVC5HTs4UhQ8gtxZU1970WtV2DyW3XDJfJ1OsKpQKp33HDV6KzX8rNLOtcHHJE7XNjsyOGQX+Qsz04I0gh1ITd+gzfG54fIXgw4Esi0OVv9P8tzOC7/69f+n+VVR+UrvuoNqZgl/P9HLgSo1ul1cbkkjCjMxemzNeiU2yRixwuXskhX5IdUrOYQ8PoYWXw+iezFgCNBnCivMpVvIc+LSHLVT0LSG76wKpI9HNW1+ivbvjnpfADxNYfdjh4OuQQuv+AI8fS3RuREXGAhQVSphh/0LvzywCJZtphQpGaERHKIo8Ygr0SSpDi8AEQ+H0L+ssRylsTvr6mbZfa3W/rFrA2xFnd/bkQOwx6OBOFV9VD8+qPNaJudjou65ipul/dwbNhfXP/4CEUKkZylEu4wT7TZ0cMhfz7tvuAZDQP9eEh7XN+/DRqlNNxaOow3iOzFHo4EoXWdn/jGyoCL9vGyKs3H7zgWmXU1vBGMEcxNi40fduRwKPZvMeKI9LPXkIMNgEMqRHZjD0eC8Oj0LPR56nPF73+Y933QfZ2ttlZdFKj/Jh7t1WLjVoSuTbkZqRjZLQ/JbpeptWXIPgw3iOzFgCNB6F3oK2r0Ey/13PvO2jDaEfJDA+Q00V9DJR7ZU4dDwuvnEmQtsxj88Qu8MT4/EcDnkAxwSCVBRDJ3YtG243HRjjdCvdDGiHJ5en6yOg1LmxPZiwFHgoh1uoNQ/R8Jfdo0jeDe7BfOyq624AWSiBIIA44EEak6GqH60/xtABp2DodyefoYNsTH4mvRgF86U+LiNU10/BsjAww4EkSsL/RbDpWea0dMmxFTdhf+slusg9Z4x2EyInsxaTRBRHI6ajgiEfikJLlwy6B2EWhN7CTiQl+RWk/HqRLvFY1DfBLJAAOOBFFSUWNqu5wmqThRrl2LIxLCTRrNz0rDV4+MQGpS4tV8SPQejtp4iVrjlDxp9I1JA2PYEiJn4pBKAqiq9eDHr68wva2dwr1mtcxKS8hgwwkKW2TGuglxTR5E9miVFbuGEDkUezgSwNES8z0W6jVXIumVxTtDXi32Lzedhw/XHcTvru4V4VZFjzJpNPG6OJ64ojuapifjuv6tY92UuKSc9hyzZhA5FgOOOPfJhkNYv6/Y9PbVNgYc0+d+j0t7tLD8uLcnD8bQLjm4pl9iX+jiblqsRU3TU/DEFT1i3Yz4leCvbzzg80ZGGHDEMSEEfj5zXayboRBKDke3Vhk2tCT64i2Hg3NOIktixBE2/k2SEeZwxLF4nMYYSpOcUsFRMaTCK5LjKHuw+PoSRRoDjjhW44nHgMN6m1wO+eyOtx4OiizmcBDZiwFHHIvHugmh9HAkYoKlFknnZ3IGp/ydEsUrBhxxzM4E0FCFksPBHg6ihoFvCzLCgCOO1cRhD0codb+cksMRb30cXBslsiK5EjIRBWLAEcfiMeAIpVqlcwKOeg48JZLhy0sUeQw44lhcBhwhJLI65eKc6HU4yDz2dRBFHgOOOFZdG38fe7UhZI06JuCQ/+yUkyIiihIGHHEsHhfbashDKvIgwxlnRHr4+hJFHgOOOBaPhb9CGVJxTMAh/9kZp0QyacluzZ+JKDJY2jyOhVJky26hBEGcFkuJoHFqEl75cX8IUfczEUUWezji2IHTFbFuQgCjHI43bztf83an5DvEW2nzUFfuJX1jerXC2N6tYt0MIkdiwBGnjpdV4f5Z62PdjABGM2cuLsyLYktiLPbxBlHcccqXC7IHA444tflgSViP/+OP+kSoJUrxmFcSLfE2LTYeelmI5Fg8jYww4IhT4XaXX9LNnt6GeFxQLhbi4Zsch1SIKJEw4IhT6i8KjVOsZc0naWRq9miVGU6TAACeOJyqGy3x1sNBRJRIGHDEKXXAYfUbtVsj4Lj34k7hNAlAaIW/nCgOOjiI4k489PxR/OLcrzilvqxbfRurAw6XFJl6GGbrcCS5JDw8uhC98rPCPma8UBb+4gcrEZEVDDjiVEDylcXrmzrgkCQpIvUwzFYalSTg7ovC71GJJyz8RUQUOg6pJAir17ckl/KllUzuJcVt/CfRkIdU4i2HgxMCiCiRMOCIU+rZIFbHRrV6M8z0cASb3WL2Iuf4IQeHnx4RUaQx4IhTU95dG9bj1QGKZDKH44LOzfH45d3COjbgzCmb8VZplIgokTDgSBDh5gxIkOAy8WpLkoQbB7YN72AOxbVUiIhCx4AjQUTi+ibv9Uh2a+9RAqe26ZF0fiYiouAYcCSIsIMASXmRvO3CjjrHcc7qrnZiUEZEZA0DjgQR7uVN/Xh5Psf79xQpt+XFVBuHVIgM8W1BRhhwJIiwOzgkZTExeS+GS1XQKpQejr/cdF7IbUsUyqTR2HNeWi4RORkDjgZKXhhMUUFTsj4Do1vLDFzTr3XE2havmDRKRBQ6BhxxSr34mplhDqNNJCi7OOT7kx9KCrIf8uGTRKTGXjcywoAjTnVXrexq5vIWbBt5bQy3zrogZut1NER8VoiIQseAI07VeKwvA2/UCyJJyiqhPfMzFfept6VATKYlMsZ3CBnh4m1xyuNVlzYP/pigPRyyXY4ozMX/3dgXhS0zVPuQ2MOhI94Wb2uSyrcvESUO9nDEqYCAw8R3B/VF8Pr+bWSPV46vul0SruvfBj3zs5T7DqEOB7/5x8ZtF3aIdROIiExjwBGn1KuymuvhUG70/I19ZY+XFEveq2em1O+DAYSeeHta0lPYw0FEiYMBR5yqrrWew2FE3cOhuE8x3TPOrqpxhAu2ERnjxwcZYcARp85U1yp+N/U+DrKRmaXlfbuYPFS79HmDJmn+SETnmPmMoYaLAUcc+mjdQZRVqgIOM3U4gt6p/Wmg9c39iSt6BD2eEX7wEBGRHAOOOPTAf9ZHfJ9GwQgraJoTj8/NHcPYE0VEiYEBR4IwlTQa4pCKHdM94/HiHK54PKVfjeuBvm2bxroZRERBMeBIEOEv3iaZSxoN4bIajxdiO8RrQm2qm29jig9x+hahOMFPqgRhqg6HwTaSBDRvnKL7SPl2FFxcBR9x1BQiIj2cyB8jQghLF62wezgADOqYjftHdkGXFk0ium8tTkwajdfr+theLbFy9ym0btoo1k0hItLFgCMGPtt4GI9/uAkv3dIfQ7vkmHqMBKBzXhPsOFYe8nElScIvLu2que9wxNOXfTvF63n+pKgD2jdPx3ltm8W6KUREujikEgNT3l2Lkooa3PrGCkuPe+v2QbiuX+uQjmm8sJtkaruGLl4Lf7ldEi7p1gLZukNmRESxZyngePnll9GnTx9kZmYiMzMTRUVFmDt3rl1tczz1sIMwGIeQJAmtmzbCHcML9PenmxZqTNL5ORxOj1usrjdDRNTQWQo42rRpg+nTp2PNmjVYvXo1LrnkElx99dXYsmWLXe1rMN5bvR+Dnl2IzQdLNO/3Xd9CvZDbWYfD6cGFT7K7/kTbN28cw5YQESUeSzkcV155peL3Z555Bi+//DK+/fZb9OzZM6INa2ge+e9GAMDPZ67T3kDy/Wft6p6S5EJ1rdd0rYZIDRs4MWk0ye3CZz8fCiGArEbJsW4OEVFCCTlp1OPx4L333sOZM2dQVFSku11VVRWqqqr8v5eWloZ6yAahxmu8aJtRb4LWRX7Oz4dh9ur9uMtgKCZecxPiUc/8rFg3gYgoIVkOODZt2oSioiJUVlaiSZMm+PDDD9Gjh/66G9OmTcPTTz8dViMbEr2eAde5SMNqaNA5rwkev7y74TYsbU5EkcAvL2TE8iyVwsJCrF+/HitWrMA999yDiRMn4rvvvtPdfurUqSgpKfH/279/f1gNdrpg5ccNezgicPxQPi74IUNERMFY7uFISUlB586dAQADBgzAqlWr8MILL+DVV1/V3D41NRWpqanhtdLhSipq/D97gyQ/2DFtVa+H46q++fhkw6Gw90lEDQPf92Qk7DocXq9XkaPRUOw7eRbT5mzFsdLKsPd16ky1/+dar0DP/EwAwC/HFPpvl/xJo5GnDGLqf372ut74791F6NEq0/I+nZg0SkREobMUcEydOhVLlizBnj17sGnTJkydOhWLFi3ChAkT7Gpf3HrgP+vw6pJduP2tVRHdr9crkJbsBgB0yq0vQe4btjDs4QjxIq+3WmyT1CQM7JCNNyadb/x4fqshIqIgLAUcx44dw09+8hMUFhZi5MiRWLVqFebPn49LL73UrvbFrbX7igEAmw+GP+tGfr2u8Xjh8dZFDi5FBdDAbeUu6pob+vGDBAwts9Lw4yHtQt4/ERGRpRyO119/3a52NGjyC35pZa0/j8PtAlLcLlR7vBjcMTtgW583bzsfF3XJRbcn5oV2fPlqsTrbuCx2Y3BEhajh4dIIZISLt0WI1dVfjfgCDpck4fNfDMfn3x3Bj4e0B6AMDkb3bAG3S8JFXXLhckmhlzZXJI1qn4NRwKF1D0t/ExGRHBdvi4BNB0pw/jML8d7q0Kb8qhMsPedqf7kkCR1yGuPO4Z2QnlIXG8qv+4+O6Ya/TxgAV5hXdzNrqViNpThVloiI5BhwRMD9s9bhRHmVvzy5VR5VxCH8QyrGF227p8jKWR1SYbxB1PDwbU9GGHCEaUD7ZgEBg1Ver/LxvqRRrWu8YvgjrKPKdxqRTYiogWMKBxlhwBEmlxT+xVgdsPh+d2u8eyWNmSvhUiSN6vVwGPW2aLUz3EYREZGjMOAIk0uSrA83qHh0eji0hlSU+RbK+0PtaFH2mmifi/oU/3u3/oJ9WtsTkfPxfU9GGHCEQD4E4pIk7DpxJsz9KX+v9fiGVLR6OLR/DoeZ3aiDqoEdsi1tT0TOx2RxMsKAIwTyIZA9J60HG/KOC49XBAyp1JybpqLdw2HDkIqJHRmOqJi8jYiIGi4GHCGQD4EcLrG+lkqSu/5pr671Bgyp1PorjQY+Vh4bqAOSUFNXzQQH/OZCREThYOEvC4QQqKjxwBtmGc0klwTfcm1VtR7/NFgfXw+H1rCE/JZIDVuY2Y3RTByte1hxkIiI5NjDYcFtM1ahx2/mY0+YORvyS3FVrRfrzq3L4uPL4dCsw2FLDkfwHamn7gbfJxE1NPyeQUYYcOg4VlqJ2av3o7LG479t0bbjAID/rAqtoigALNx6FGeq6/dZVePFM3O2Krap9er3cMipp82qe0pMM/EhUePR37fmw/nBQ9Tg8G1PRjikouPavy/DweIK/HCkDL++oofivtowxlQmv7Va8XtVrSdgm1r/tFiNHcgOHalhCzO7OXWmyto+Q2wLERE5E3s4dBwsrgAALNh6NOA+q8MLRqpqvQG3+ToqgvVw2HFR11sArsbqkAr7VokanNyM1Fg3geIYA44gajWGEqo9gUFCqORDNmpaAYe8Neq77ZylYjRcI2/HNeflAwDuu7hziK0hokTz5qTzcVHXXPz+mt6xbgrFMQ6paCitrPH/7NW40FZr9EqEymh4RitpVN6eSE1VNdMbYTY95Pkbz8O9F3dGl7wmYbaKiBLFxd3ycHG3vFg3g+IcezhUNh0oQZ+nPvf/rq6RAZjv4ais8WDKO2sx2yDJ1Gh4RisOkF/4pQi9embCFq3AS4vbJaFriwwOqRARkQIDDpVXluxU/B5OD8d7q/fjs02H8cv365at/+FoWcA2RvUttHo4ktz1tyWpC39FYC0VPUb7ZmhBRETBcEhFRX0R1+qA0JpZoqW0slbx+yfrDwVsYzikohEJ5GWk4a7hBUhJciE9JTIvn6k6HJHLkyUiogaIAYdKQLnwMHo45MHL9qNlmj0JxkMq2oHA1Mu7mzq+WfLD6PVkhFzjg4iICBxSCaDuVdCKB8zmcMjXTLn0z0uw63hghVKtHBF/W4xWTIsyszkcREREWhhwqARMNQ2xh2PzwRL87n/fKW5bsftkwHZeITCsS47mPqIVb5hZ8p7hBhERhYMBRxBaF9pgAYfXK3DFi0sD96WxM48XyExL1tyPK0oRR7g5HEa9NERERAADjqC0goQ9J88aPmbtvtOat2sliHqE0K3uqZU0agdzs1T0g4pqg3VWiIiIAAYcAdTX1UjmLtRo5H7U1Hrh1ekwidTy81bona7R86B1XkRERHIMOIIIJd5I0lx1DThbHTid9qH3Nuj2cLii9OqYK20eeFvrpo0AAJewwiAREQXBabFBhNLDYfUxeikQ0RtSMZPDEdjID++9AIt/OI4r++bb0SwiInIQBhwqvlVifULp4bBas0Jv82gNqZg5SrJGr01eZhpuGNg28g0iIiLH4ZCKzMYDxVi2Uzl1VW+4w4jVlAa9ACVqs1RMFP56+qqeUWkLERE5EwMOmY81So+HMuOz1mLEoXWIaNb8kg+p6BUbK8htgr5tsqLVJCIichgGHDJa19pQSnrXWIxStPIjYlVlNFq9KkRE1LAw4JDRypkIZVKs5R4OjYPEYkosECRRlUvOExFRiBp8wPHJhkNYfi5vQ2u2hhDmF2sDgJW7T2He5iOW2qDVwxGzgMOgh8NooTkiIiIjDXqWyr6TZ/HzmesAANufGaubN3HDK8tw+9COpvY58Y2VqKgxt3y9kVgNqRgdlwW+iIgoVA2uh2PbkTK8v+YAhBA4cabKf/vZKo9ur8KGAyW4f9Z6U/uv9nhxfodmltqktRZJrEYvGHAQEZEdGlwPx+i/LAEANElLQk6TFP/t9/9nHZo3Tg17/0II9MzPwqo92uupaNEKOKLdw9EzPxOHiivQMz9Td5sueRnYefxMFFtFRERO0eACDp91+4oxXLYs/KJtxyOyX68IpdKoRsARQhfHHcM64rWvd+PO4QWWH/vpfUPhEUKzwJfP76/thewmKbhlUDvL+yciooatwQYc1bVeVFlIBrXC6kxa7SEV6wHHY2O74+rzWqN7K/1eCj0ulwRXkJqjOU1S8ey1vS3vm4iIqMEGHDUeL6pqw0/u1GK1h0N7SMX6cd0uCb1aszgXERHFnwaXNOrjFcK2Hg6rs0c9ERpSISIiilcNNuAQAKpq7BpSsdrDEXhbKEMqRERE8apBBBwerwgo3iWEQJVN0zwtJ43GwSwVIiIiOzWIgOOKF5diyLSFqJQV5PJ6rVUQtSLYkMojowvx2c+H+n+v9Qa2g/EGERE5ieOTRoUQ2Hq4FADww9Gy+tsR2OsRKcF6OG4tao9GyW7Z9oHbxKq0ORERkR0c38MhTwyVZNM+Z68+gM0HS+w5aJAeDpckKQIKrVkqXLWViIicxPEBh3wYRd1p8Nmmw7YcM1gPh0tSDplo1uGIdKOIiIhiqAEEHPU9HFoXdjsEO4xLkhSzUDR7ODikQkREDuL4gEO+cqtWvQs7nKmqNbzfF0v4/j9bHbg94w0iInISxwccS3ec8P9stT5GqI6XVxne7+u98P1fWqkVcDDiICIi53B8wJGaVH+K0VpdPVgOhy+UMMoLZbhBRERO4viAIzOtfuav1YJcoQoW2Ph6Nox6MVyOf2WIiKghcfxlrcZTH2RoVfS0g0ejkJecP4fDaBv2cRARkYM4PuCQzwCpiXDAccewjkGPqUVS5XBoYRkOIiJyEscHHLWyi39thJM4MtKSNW83O/3WMKhg0igRETmI4wMO+fCGfHglEpLc2kFBremAgz0cRETUMDhqLZXXluzCa1/vwnltm+KWwe0wojBP2cMRJLfCqmSdzE6zPRxGnRiMN4iIyEkc08Ph9Qo8M2crjpVV4fPvjmLSm6sAKC/+tXHWw6GepdKrdab/Z1YaJSIiJ3FMwPHGN7s1b5cPo1RHOIcj2R1eD4d82KRbywy4ZT0mjDeIiMhJHBNwfL7lqObt8hyOSPdwpOgEHGaTU+W9GC5JORGWlUaJiMhJHBNwpCTpXPxlvQ2Pf7gposfUG1Ixn8MhCzhcyl4NhhtEROQkjgk4kvXyKSLcqyGXpDekYrKiqUsRYEgBPR5ERERO4ZiA46ttxzVv/78FP9h2zJQweziUAQZUQyrhtIyIiCi+WAo4pk2bhvPPPx8ZGRnIy8vDNddcg23bttnVtrjXKEV7VrH5WSrynyXF7+zhICIiJ7EUcCxevBhTpkzBt99+iwULFqCmpgaXXXYZzpw5Y1f7ouK3V/cM6XFJsjERt+xns2vEyYMKSVLmdDDeICIiJ7FU+GvevHmK32fMmIG8vDysWbMGw4cPj2jDouniwjwAWyw/LtxeCHWPBmepEBGRU4VVabSkpAQAkJ2drbtNVVUVqqqq/L+XlpaGc0hb6M02MbJ72uVYvfe0/3eXBHgs7iMgh4OzVIiIyKFCThr1er144IEHcOGFF6JXr166202bNg1ZWVn+f23btg31kLZxh7BwiSRJqiER6/swnqVieXdERERxK+SAY8qUKdi8eTNmzZpluN3UqVNRUlLi/7d///5QD2no7cmDQ35sks6aKMHIC3+FEh+oczbUSaREREROEdKV9r777sP//vc/fPXVV2jTpo3htqmpqcjMzFT8s8PQLjm4rn/rkB4bSg8HACQnmUvy7NEqEwPbNwu4PTCHgz0cRETkTJYCDiEE7rvvPnz44Yf48ssv0bFjR7vaFRKzs0PUklySbplyI/K1VIwSSAW0A5LAWSryexlxEBGRc1i6yk6ZMgVvv/023n33XWRkZODIkSM4cuQIKioq7GqfJV5ZxJGe4obXZD0Mt0vC3AeGWT6e2SEVoRMJudQ9HMzhICIih7IUcLz88ssoKSnBiBEj0KpVK/+///znP3a1zxL5db1RslsRgBhJcknolNvE8vHks1uC5VxIGiFJQA+HfHsGHERE5CCWpsXqfVOPFz8a0AafbDgEoK63w+yaJqHmcLhNFuoSAkFHSCRJMj1EQ0RElGgcs5YKAAzvmot/3DoAAOAV5nM6Qp0RophlYrCdgN6QinIIJS1ZNkTDeIOIiBzEUQEHAHRpkQEAKKmowU3/+NbWY5ntGLm4W55mQCKfjeuSJDRKdvt/57RYIiJykrAqjcYjeRCwYX+xrcdSrJ+icf/LE/qjrKoWV/XNx6Q3Vwbc71L1kCSFWdeDiIgoXjkw4IjepVrRC6ERcbTISsPYdnX1N7SSRtVrp3C1WCIicirHDalEk8s43ggaNKgrjSpKnTPeICIiB3FcwOGKYgELeUChNYMnWACh7NFQJ5Ey4iAiIudwXsARxnXaarVReQ7HtRpl1dV1Nozud6kXg7PUEiIiovjmwIAj9Ev10scuxrt3mF8ETn6onvlZhvdrUfeAcPE2IiJyKscljYZznc7LSENeRprp7ZVDKsb3ayeNKiuVBusRISIiSlTs4QhBq6w0U8cKFkCoZ6W4VDkdRERETsGAIwT/ubPo3LHqb9OqJhqsKYohFKgrlzLiICIi53BgwGHv/ru2aIJ2zdMBKAME7SGV+p+1cjLUpc0VPR6Oe2WIiKghc9xlLVbJllp1OORtyU5PDrjfaJYK56kQEZGTODDgCLwt2W3t4v3X8f3Qo1Wm5n26C8Jp1uGoP+5jY7ujTbNGeGR0oXZbVYW/mMNBRERO4rhZKlo5HEkuF2o8HtP7uKpvPsb2aokuv5obcJ/JBWjPtaX+55ZZaVj66CWK+yWDHg7OUiEiIidxYMAReFuSWwJqAm8/v0Mz3DGsQHM/etd7rYqiQIilzRXbBgYgRERETuHAgCPwQq1XQXT2XUW6OR96t5dX1eoeu6igOZbvOinbh1FLVUmlkFS/ExEROUeDyOFI0snhMEow1bvnTJX20IxLkvDPiQNN79/3GP/PLnXdDoYcRETkHI4LOPRyOKzSu97rDam4JAmNU5OQnuKW3Wb+GJIkqZarN9lQIiKiBNAgAg6rs1QA/R4Gj07AkXQuunDLHhcs0FEW+lKudMscDiIichIHBhyBtyVZXAXWiFdnmoovWCiT5Xi4g3RxKKfBSgGVR4mIiJzCcQGHVs9EUiSLWugEHH3bBK4WGyzgkA+iuCR1TgdDDiIicg7HzVLREsnhCfWQyjePXYJjpZXo0iIjYNtggY58xEWSOEuFiIicq0EEHNUeb8T25VUFHK2bNkLrpo00t7XSwyFJnKVCRETO5bghFS07jpVHbF+6pc01BA04AnI4WGmUiIicqUEEHLHitlCHQ0Jg5VEiIiKnYMBho2CJn8rl6NU5HIw4iIjIORwZcORmpMa6CaaoeziUdThi0CAiIiKbODLgyExLjFxYZWVRZQ4HkziIiMhJHBlwaOV1/nRox6i3Ixjl6rDqQmAxaBAREZFNHBlweDXKgf76ih5RbUPXFk2CbuNSdWgoh1gYcRARkXM4M+CwMHXVLmYChoBpsZD/Hvk2ERERxYpDA47YRxxmUjDUhb6YwkFERE7lyIBjUIfsoNsMbN8sCi0xpl6sTWKlUSIicihHBhxPXtUTgzvqBx2DO2Zj1p1DotgibcqkUUk1ayX67SEiIrKLIwOOrEbJ+NklXXTvb5WVFtEl60OlnpWiDkCIiIicIvZXXZsIvXXkEZ2LuZkhkcDF2+T3EREROYdjAw6jmSpm8yOGdckJ+fjBlqYH1NNipYBZK0RERE7h2IBDGMxUMXstf+u2QXjo0q4hHT/YSrF17VDncHC1WCIicibnBhyyn9XX/n7tmprah8sloUmIZdLNBRzKnyVVjwcREZFTJMaiI6GQRRxJrrq46osHh2Pl7tO46fy2pncTahGxQQazZHxcih4O1bTY0A5LREQUlxwbcNR4vP6fz8Ub6JyXgc55GZb2YzQ0o+XLhy7Cl98fw4+HtA+6rbKyKCuNEhGRczk24PBEqL651aKlBblNUJAbfB0VQLkcPcAhFSIici7H5nDUyAKOcBZCM5peGy51D4d6iIWIiMgpHBtwZIaY7Klm50Jw6uXpJeWd9h2YiIgoyhwbcAzvkuv/uXFq6MGHnQvBGdfhsO2wREREUefYgMPlkvDstb2RmZaE134yIOT92LnwrDrAUM5SYcRBRETO4dikUQC4ZXA7jB/UNqwETKuzVKwIWJ5ecZ9thyUiIoo6x/Zw+IQ72yPZxkXelMvRB/5ORETkFI4POMI1YUh79MzPtGXfgbNUZPcx4iAiIgdhwBFEk9QkfPbzYWjdtBEAoHfrrIjt26WqLKqowxGxoxAREcWeo3M4ImnWnUPw72/34rYLO0Rsn+rVYeVZHFwtloiInIQBh0lts9Px+OXdI7pP5RAKFN0ajDeIiMhJOKQSQ5Jqloqy0igjDiIicg4GHDEUUIdD5z4iIqJEx4AjhtQ9Gly8jYiInIoBRwypezTUs1aIiIicggFHDBlXGmXIQUREzsGAI4YCFmvjLBUiInIoBhwxpF6sjUMqRETkVAw4YshlOEuFIQcRETkHA44YUgcYkmLWSvTbQ0REZBfLAceSJUtw5ZVXIj8/H5Ik4aOPPrKhWQ2DyyVPGlX3eDDiICIi57AccJw5cwZ9+/bFSy+9ZEd7GhTJsA5HDBpERERkE8trqYwdOxZjx461oy0NjnIarPIW9nAQEZGT2L54W1VVFaqqqvy/l5aW2n3IhKGswwH2cBARkWPZnjQ6bdo0ZGVl+f+1bdvW7kMmDHUpc3UhMCIiIqewPeCYOnUqSkpK/P/2799v9yEThjpJNHCIhYiIyBlsH1JJTU1Famqq3YdJSJKq0Jeix4Olv4iIyEFYhyOG1GunuFiHg4iIHMpyD0d5eTl27Njh/3337t1Yv349srOz0a5du4g2zumMAgzmcBARkZNYDjhWr16Niy++2P/7gw8+CACYOHEiZsyYEbGGNQSKmIKzVIiIyMEsBxwjRoyAEMKOtjQ4LlXhL/XvRERETsEcjhhSJolqLFdPRETkEAw4YkhR2twlKWamsIODiIichAFHDLlUPRyK3xlxEBGRgzDgiCF5j4bLJamGVBhwEBGRczDgiCF5j4ZbkiCvzMFwg4iInIQBRwypl6dXlzonIiJyCgYcMaSuu6Eodc54g4iIHIQBRwzJezHcLuXqKQw4iIjISRhwxFDAarEcUiEiIodiwBFD6kJfrDRKREROxYAjhtSFv5T3Rbs1RERE9mHAEUPq5elZ2pyIiJyKAUcMuV3yIRTlMAorjRIRkZMw4IghZcAhBSzmRkRE5BQMOGIoyVX/9KsXb2PSKBEROQkDjhgKHFKB7HcGHERE5BwMOGIoSTWkIh9HYbxBREROkhTrBjRk6hwOOQYcRETkJAw4Ykg9pCLHWSpEROQkDDhiSF1ZtHFq/cuRzEIcRETkIAw4YkiRJOqSkJLkwhuTBuJEeTXyMtNi1zAiIqIIY8ARQ1qVRS/p1iI2jSEiIrIRZ6nEkMTF2oiIqIFgwBFD8iDDzZwNIiJyMAYcMcQYg4iIGgoGHDHUtFGK/2f2cBARkZMxaTSGstKT8cakgUh2u5DsZuxHRETOxYAjxjgrhYiIGgJ+rSYiIiLbMeAgIiIi2zHgICIiItsx4CAiIiLbMeAgIiIi2zHgICIiItsx4CAiIiLbMeAgIiIi2zHgICIiItsx4CAiIiLbMeAgIiIi2zHgICIiItsx4CAiIiLbRX21WCEEAKC0tDTahyYiIqIQ+a7bvuu4VVEPOMrKygAAbdu2jfahiYiIKExlZWXIysqy/DhJhBqqhMjr9eLQoUPIyMiAJEkh7aO0tBRt27bF/v37kZmZGeEWxpeGcq48T+dpKOfK83SWhnKegPVzFUKgrKwM+fn5cLmsZ2REvYfD5XKhTZs2EdlXZmam4/8gfBrKufI8naehnCvP01kaynkC1s41lJ4NHyaNEhERke0YcBAREZHtEjLgSE1NxZNPPonU1NRYN8V2DeVceZ7O01DOlefpLA3lPIHon2vUk0aJiIio4UnIHg4iIiJKLAw4iIiIyHYMOIiIiMh2DDiIiIjIdjELOF5++WX06dPHX3CkqKgIc+fO9d8/YsQISJKk+Hf33Xcr9rFv3z6MGzcO6enpyMvLwyOPPILa2lrFNosWLUL//v2RmpqKzp07Y8aMGdE4PV3Tp0+HJEl44IEH/LdVVlZiypQpaN68OZo0aYLrr78eR48eVTwu0c5V6zyd8po+9dRTAefRrVs3//1OeT2DnadTXk8AOHjwIH784x+jefPmaNSoEXr37o3Vq1f77xdC4De/+Q1atWqFRo0aYdSoUdi+fbtiH6dOncKECROQmZmJpk2bYvLkySgvL1dss3HjRgwbNgxpaWlo27Yt/vjHP0bl/HyCneekSZMCXtMxY8Yo9pEI59mhQ4eA85AkCVOmTAHgnPdosPOMu/eoiJFPPvlEfPbZZ+KHH34Q27ZtE48//rhITk4WmzdvFkIIcdFFF4k77rhDHD582P+vpKTE//ja2lrRq1cvMWrUKLFu3ToxZ84ckZOTI6ZOnerfZteuXSI9PV08+OCD4rvvvhMvvviicLvdYt68eVE/XyGEWLlypejQoYPo06ePuP/++/2333333aJt27Zi4cKFYvXq1WLIkCHiggsu8N+faOeqd55OeU2ffPJJ0bNnT8V5HD9+3H+/U17PYOfplNfz1KlTon379mLSpElixYoVYteuXWL+/Plix44d/m2mT58usrKyxEcffSQ2bNggrrrqKtGxY0dRUVHh32bMmDGib9++4ttvvxVff/216Ny5sxg/frz//pKSEtGiRQsxYcIEsXnzZjFz5kzRqFEj8eqrr8bNeU6cOFGMGTNG8ZqeOnVKsZ94P08hhDh27JjiHBYsWCAAiK+++koI4Zz3aLDzjLf3aMwCDi3NmjUT//znP4UQdU+U/GKlNmfOHOFyucSRI0f8t7388ssiMzNTVFVVCSGE+OUvfyl69uypeNxNN90kRo8eHfnGB1FWVia6dOkiFixYoDi34uJikZycLN577z3/tlu3bhUAxPLly4UQiXWueucphHNe0yeffFL07dtX8z4nvZ5G5ymEc17PRx99VAwdOlT3fq/XK1q2bCmee+45/23FxcUiNTVVzJw5UwghxHfffScAiFWrVvm3mTt3rpAkSRw8eFAIIcTf//530axZM/+5+45dWFgY6VPSFOw8hagLOK6++mrd+xPhPLXcf//9olOnTsLr9TrqPaomP08h4u89Ghc5HB6PB7NmzcKZM2dQVFTkv/2dd95BTk4OevXqhalTp+Ls2bP++5YvX47evXujRYsW/ttGjx6N0tJSbNmyxb/NqFGjFMcaPXo0li9fbvMZBZoyZQrGjRsX0J41a9agpqZGcXu3bt3Qrl07fzsT6Vz1ztPHKa/p9u3bkZ+fj4KCAkyYMAH79u0D4LzXU+88fZzwen7yyScYOHAgbrjhBuTl5aFfv3547bXX/Pfv3r0bR44cUbQzKysLgwcPVrymTZs2xcCBA/3bjBo1Ci6XCytWrPBvM3z4cKSkpPi3GT16NLZt24bTp0/bfZpBz9Nn0aJFyMvLQ2FhIe655x6cPHnSf18inKdadXU13n77bdx+++2QJMlx71Ef9Xn6xNN7NOqLt8lt2rQJRUVFqKysRJMmTfDhhx+iR48eAIBbbrkF7du3R35+PjZu3IhHH30U27ZtwwcffAAAOHLkiOJJAuD//ciRI4bblJaWoqKiAo0aNbL7FAEAs2bNwtq1a7Fq1aqA+44cOYKUlBQ0bdo0oJ3BzsN3n9E20TxXo/MEnPOaDh48GDNmzEBhYSEOHz6Mp59+GsOGDcPmzZsd9XoanWdGRoZjXs9du3bh5ZdfxoMPPojHH38cq1atws9//nOkpKRg4sSJ/rZqtVN+Hnl5eYr7k5KSkJ2drdimY8eOAfvw3desWTNbzs8n2HkCwJgxY3DdddehY8eO2LlzJx5//HGMHTsWy5cvh9vtTojzVPvoo49QXFyMSZMm+dvglPeonPo8gfj7zI1pwFFYWIj169ejpKQE//3vfzFx4kQsXrwYPXr0wJ133unfrnfv3mjVqhVGjhyJnTt3olOnTjFstTX79+/H/fffjwULFiAtLS3WzbGNmfN0yms6duxY/899+vTB4MGD0b59e8yePTvqHzJ2MjrPyZMnO+b19Hq9GDhwIJ599lkAQL9+/bB582a88sor/guxE5g5z5tvvtm/fe/evdGnTx906tQJixYtwsiRI2PS7nC9/vrrGDt2LPLz82PdFFtpnWe8vUdjOqSSkpKCzp07Y8CAAZg2bRr69u2LF154QXPbwYMHAwB27NgBAGjZsmVAVrHv95YtWxpuk5mZGbULw5o1a3Ds2DH0798fSUlJSEpKwuLFi/HXv/4VSUlJaNGiBaqrq1FcXBzQzmDn4bvPaJtonWuw8/R4PAGPSdTXVK1p06bo2rUrduzYgZYtWzri9dQiP08tifp6tmrVyt+z6tO9e3f/8JGvrVrtlJ/HsWPHFPfX1tbi1KlTll53OwU7Ty0FBQXIyclRvKbxfp5ye/fuxRdffIGf/vSn/tuc+B7VOk8tsX6PxkUOh4/X60VVVZXmfevXrwdQ96YBgKKiImzatEnxx79gwQJkZmb631RFRUVYuHChYj8LFixQ5InYbeTIkdi0aRPWr1/v/zdw4EBMmDDB/3NycrKindu2bcO+ffv87UyEcw12nm63O+AxifqaqpWXl2Pnzp1o1aoVBgwY4IjXU4v8PLUk6ut54YUXYtu2bYrbfvjhB7Rv3x4A0LFjR7Rs2VLRztLSUqxYsULxmhYXF2PNmjX+bb788kt4vV7/h3xRURGWLFmCmpoa/zYLFixAYWFhVIYZgp2nlgMHDuDkyZOK1zTez1PuzTffRF5eHsaNG+e/zYnvUa3z1BLz96jlNNMIeeyxx8TixYvF7t27xcaNG8Vjjz0mJEkSn3/+udixY4f47W9/K1avXi12794tPv74Y1FQUCCGDx/uf7xvOs9ll10m1q9fL+bNmydyc3M1p/M88sgjYuvWreKll16K6bRYH3Xm8N133y3atWsnvvzyS7F69WpRVFQkioqK/Pcn6rnKz9NJr+lDDz0kFi1aJHbv3i2++eYbMWrUKJGTkyOOHTsmhHDO62l0nk56PVeuXCmSkpLEM888I7Zv3y7eeecdkZ6eLt5++23/NtOnTxdNmzYVH3/8sdi4caO4+uqrNafF9uvXT6xYsUIsXbpUdOnSRTFdtLi4WLRo0ULceuutYvPmzWLWrFkiPT09atNFg51nWVmZePjhh8Xy5cvF7t27xRdffCH69+8vunTpIiorKxPmPH08Ho9o166dePTRRwPuc8p7VAj984zH92jMAo7bb79dtG/fXqSkpIjc3FwxcuRI8fnnnwshhNi3b58YPny4yM7OFqmpqaJz587ikUceUcwfFkKIPXv2iLFjx4pGjRqJnJwc8dBDD4mamhrFNl999ZU477zzREpKiigoKBBvvvlmtE5RlzrgqKioEPfee69o1qyZSE9PF9dee604fPiw4jGJeK7y83TSa3rTTTeJVq1aiZSUFNG6dWtx0003KWoZOOX1NDpPJ72eQgjx6aefil69eonU1FTRrVs38Y9//ENxv9frFU888YRo0aKFSE1NFSNHjhTbtm1TbHPy5Ekxfvx40aRJE5GZmSluu+02UVZWpthmw4YNYujQoSI1NVW0bt1aTJ8+3fZzkzM6z7Nnz4rLLrtM5ObmiuTkZNG+fXtxxx13KKZMCpEY5ymEEPPnzxcAAl4nIZzzHhVC/zzj8T3K5emJiIjIdnGVw0FERETOxICDiIiIbMeAg4iIiGzHgIOIiIhsx4CDiIiIbMeAg4iIiGzHgIOIiIhsx4CDiIiIbMeAg4iIiGzHgIOIiIhsx4CDiIiIbMeAg4iIiGz3/98GpKLsaeQWAAAAAElFTkSuQmCC", "text/plain": [ "
" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "spec = miles.interpolate(\n", " age=5.7*u.Gyr, met=-0.45*u.dex, imf_slope=1.3,\n", ")\n", "plt.plot(spec.spectral_axis, spec.flux)\n", "plt.show()" ] }, { "cell_type": "markdown", "id": "2b06a6ca-1792-4c52-b2ab-47a89db26d2e", "metadata": {}, "source": [ "2. Given a range of values" ] }, { "cell_type": "code", "execution_count": 24, "id": "83f6ebc2-ec0a-4045-ad64-a3203d538778", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "112" ] }, "execution_count": 24, "metadata": {}, "output_type": "execute_result" } ], "source": [ "spec2 = miles.in_range(\n", " age_lims=[7.0, 10.0]<" ] }, "execution_count": 25, "metadata": {}, "output_type": "execute_result" } ], "source": [ "np.unique(spec2.age)" ] }, { "cell_type": "code", "execution_count": 26, "id": "e39abcba-17ee-420e-ad85-c414e1c88c07", "metadata": {}, "outputs": [ { "data": { "text/latex": [ "$[0,~0.22] \\; \\mathrm{dex}$" ], "text/plain": [ "" ] }, "execution_count": 26, "metadata": {}, "output_type": "execute_result" } ], "source": [ "np.unique(spec2.met)" ] }, { "cell_type": "markdown", "id": "66a9518f-6634-4266-a007-d89c68c68015", "metadata": {}, "source": [ "3. From a list of parameters" ] }, { "cell_type": "code", "execution_count": 27, "id": "6b85a6f9-e48c-43d2-b072-243f64f54e71", "metadata": { "editable": true, "slideshow": { "slide_type": "" }, "tags": [] }, "outputs": [ { "data": { "text/latex": [ "$[0.2512,~0.0708,~1.4125] \\; \\mathrm{Gyr}$" ], "text/plain": [ "" ] }, "execution_count": 27, "metadata": {}, "output_type": "execute_result" } ], "source": [ "spec3 = miles.in_list(\n", " age=[0.2512, 0.0708, 1.4125] << u.Gyr,\n", " met=[0.22, 0.0, -1.71] << u.dex,\n", " imf_slope=np.array([1.3, 1.3, 1.3]),\n", ")\n", "spec3.age" ] }, { "cell_type": "markdown", "id": "bc4da328-0835-4d79-8f8c-42114e82ef17", "metadata": { "editable": true, "slideshow": { "slide_type": "" }, "tags": [] }, "source": [ "### Analyze the spectra\n", "\n", "You can load filters from the database and compute derived quantities from those and the SSP spectra" ] }, { "cell_type": "code", "execution_count": 28, "id": "0cdd5823-514c-4525-a750-3fa3058c48f8", "metadata": { "editable": true, "slideshow": { "slide_type": "" }, "tags": [] }, "outputs": [ { "name": "stderr", "output_type": "stream", "text": [ "milespy.magnitudes: Filter SLOAN_SDSS.i [6430.0,8630.0] is outside ofthe spectral range [3540.5, 7409.6]\n", "WARNING:milespy.magnitudes:Filter SLOAN_SDSS.i [6430.0,8630.0] is outside ofthe spectral range [3540.5, 7409.6]\n", "milespy.magnitudes: Filter SLOAN_SDSS.u [2980.0,4130.0] is outside ofthe spectral range [3540.5, 7409.6]\n", "WARNING:milespy.magnitudes:Filter SLOAN_SDSS.u [2980.0,4130.0] is outside ofthe spectral range [3540.5, 7409.6]\n", "milespy.magnitudes: Filter SLOAN_SDSS.z [7730.0,11230.0] is outside ofthe spectral range [3540.5, 7409.6]\n", "WARNING:milespy.magnitudes:Filter SLOAN_SDSS.z [7730.0,11230.0] is outside ofthe spectral range [3540.5, 7409.6]\n" ] }, { "name": "stdout", "output_type": "stream", "text": [ "5.7 Gyr -0.45 dex 0.6223571750923741 nan 1.9871261886373415 1.6757822396390838 nan nan\n" ] } ], "source": [ "import milespy.filter as flib\n", "filts = flib.get( flib.search(\"sloan\") )\n", "outmls = spec.mass_to_light(filters=filts, mass_in=\"star+remn\")\n", "\n", "print(\n", " spec.age,\n", " spec.met,\n", " spec.Mass_star_remn,\n", " outmls[\"SLOAN_SDSS.u\"],\n", " outmls[\"SLOAN_SDSS.g\"],\n", " outmls[\"SLOAN_SDSS.r\"],\n", " outmls[\"SLOAN_SDSS.i\"],\n", " outmls[\"SLOAN_SDSS.z\"],\n", ")" ] }, { "cell_type": "markdown", "id": "25ac7c14-6db4-4de7-a731-1fef385c1f2c", "metadata": { "editable": true, "slideshow": { "slide_type": "" }, "tags": [] }, "source": [ "## Compute mags" ] }, { "cell_type": "markdown", "id": "e9ea90dc-f654-40fe-8fd8-14bc9d0cfd94", "metadata": {}, "source": [ "Note that we can use fancy regex for selecting the filters!" ] }, { "cell_type": "code", "execution_count": 29, "id": "3356c56f-5bcc-4554-bad7-70f46d74eabd", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "{'SLOAN_SDSS.g': 6.401271533249293, 'SLOAN_SDSS.r': 5.723690182070875}\n" ] } ], "source": [ "import milespy.filter as flib\n", "filts = flib.get( flib.search(\"sdss.(r|g)\") )\n", "outmags = spec.magnitudes(\n", " filters=filts, zeropoint=\"AB\"\n", ")\n", "print(outmags)" ] }, { "cell_type": "markdown", "id": "516eef1f-438d-4cfc-87ac-b9bc6bea7e4c", "metadata": {}, "source": [ "They can be saved as astropy tables using `write`. By default, the `basic` format is used and rather than saving to a file it prints the output to `stdout`.\n", "This can be changed with the `format` and `output` parameters, respectively." ] }, { "cell_type": "code", "execution_count": 30, "id": "1b2774ab-ca31-447f-90b6-f34ca5563a2a", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "| SLOAN_SDSS.g | SLOAN_SDSS.r |\n", "| 6.401271533249293 | 5.723690182070875 |\n" ] } ], "source": [ "outmags.write(format='fixed_width')" ] }, { "cell_type": "markdown", "id": "3d054b3f-e9c0-498c-8cc6-638ddbb74087", "metadata": {}, "source": [ "## Compute LS indices" ] }, { "cell_type": "code", "execution_count": 31, "id": "cb850f49-fa86-4e40-8067-e84505b4b960", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "['H10Fe', 'H_K', 'HdA', 'HdF', 'HgA', 'HgF', 'Hg_sigma_275', 'Hg_sigma_200', 'Hg_sigma_125', 'Hg_sigma_130', 'Hbeta_o', 'Hbeta', 'Halpha', 'Ha_Gregg94']\n", "[2.41603984]\n" ] } ], "source": [ "import milespy.ls_indices as lslib \n", "names = lslib.search(\"^H.*\")\n", "print(names)\n", "indeces = lslib.get(names)\n", "\n", "outls = spec.line_strength(indeces)\n", "print(outls['Halpha'])" ] }, { "cell_type": "markdown", "id": "14ed1952-b2b9-4892-938a-62ab57c0b2ce", "metadata": {}, "source": [ "## Compute mass-to-light ratios" ] }, { "cell_type": "code", "execution_count": 34, "id": "4583fa7f-66ee-4ed5-819d-a4f6b292fa52", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "{'SLOAN_SDSS.r': 1.6757822396390838}\n" ] } ], "source": [ "import milespy.filter as flib\n", "filts = flib.get( flib.search(\"sdss.r\") )\n", "outmls = spec.mass_to_light(filters=filts, mass_in=\"star+remn\")\n", "print(outmls)" ] }, { "cell_type": "code", "execution_count": null, "id": "937956ac-c361-4fda-beb4-d106a688ffea", "metadata": {}, "outputs": [], "source": [] } ], "metadata": { "kernelspec": { "display_name": "Python 3 (ipykernel)", "language": "python", "name": "python3" }, "language_info": { "codemirror_mode": { "name": "ipython", "version": 3 }, "file_extension": ".py", "mimetype": "text/x-python", "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", "version": "3.10.12" } }, "nbformat": 4, "nbformat_minor": 5 }