{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "\n", "\n", "
\n", "**This is a fixed-text formatted version of a Jupyter notebook.**\n", "\n", "- Try online [![Binder](https://mybinder.org/badge.svg)](https://mybinder.org/v2/gh/gammapy/gammapy-webpage/v0.8?urlpath=lab/tree/spectrum_simulation_cta.ipynb)\n", "- You can contribute with your own notebooks in this\n", "[GitHub repository](https://github.com/gammapy/gammapy/tree/master/tutorials).\n", "- **Source files:**\n", "[spectrum_simulation_cta.ipynb](../_static/notebooks/spectrum_simulation_cta.ipynb) |\n", "[spectrum_simulation_cta.py](../_static/notebooks/spectrum_simulation_cta.py)\n", "
\n" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "# Spectrum simulation for CTA\n", "\n", "A quick example how to simulate and fit a spectrum for the [Cherenkov Telescope Array (CTA)](https://www.cta-observatory.org).\n", "\n", "We will use the following classes:\n", "\n", "* [gammapy.spectrum.SpectrumObservation](..\/api/gammapy.spectrum.SpectrumObservation.rst)\n", "* [gammapy.spectrum.SpectrumSimulation](..\/api/gammapy.spectrum.SpectrumSimulation.rst)\n", "* [gammapy.spectrum.SpectrumFit](..\/api/gammapy.spectrum.SpectrumFit.rst)\n", "* [gammapy.irf.CTAIrf](..\/api/gammapy.irf.CTAIrf.rst)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Setup" ] }, { "cell_type": "code", "execution_count": 1, "metadata": {}, "outputs": [], "source": [ "%matplotlib inline\n", "import matplotlib.pyplot as plt" ] }, { "cell_type": "code", "execution_count": 2, "metadata": {}, "outputs": [], "source": [ "import numpy as np\n", "import astropy.units as u\n", "from gammapy.irf import EnergyDispersion, EffectiveAreaTable\n", "from gammapy.spectrum import SpectrumSimulation, SpectrumFit\n", "from gammapy.spectrum.models import PowerLaw\n", "from gammapy.irf import CTAIrf" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Simulation" ] }, { "cell_type": "code", "execution_count": 3, "metadata": {}, "outputs": [], "source": [ "# Define simulation parameters parameters\n", "livetime = 1 * u.h\n", "offset = 0.5 * u.deg\n", "# Energy from 0.1 to 100 TeV with 10 bins/decade\n", "energy = np.logspace(-1, 2, 31) * u.TeV" ] }, { "cell_type": "code", "execution_count": 4, "metadata": {}, "outputs": [], "source": [ "# Define spectral model\n", "model = PowerLaw(\n", " index=2.1,\n", " amplitude=2.5e-12 * u.Unit(\"cm-2 s-1 TeV-1\"),\n", " reference=1 * u.TeV,\n", ")" ] }, { "cell_type": "code", "execution_count": 5, "metadata": {}, "outputs": [], "source": [ "# Load IRFs\n", "filename = (\n", " \"$GAMMAPY_DATA/cta-1dc/caldb/data/cta/1dc/bcf/South_z20_50h/irf_file.fits\"\n", ")\n", "cta_irf = CTAIrf.read(filename)" ] }, { "cell_type": "code", "execution_count": 6, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "NDDataArray summary info\n", "energy : size = 42, min = 0.014 TeV, max = 177.828 TeV\n", "offset : size = 6, min = 0.500 deg, max = 5.500 deg\n", "Data : size = 252, min = 0.000 m2, max = 5371581.000 m2\n", "\n" ] }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYgAAAEOCAYAAACTqoDjAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMi4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvIxREBQAAIABJREFUeJzt3Xd4VGX6xvHvk4TQCSIgvUmTIigRe28oAnaxrGWtu2v5re5aVteu4OpacVdxsa/orthQ7Loi6qr0XiKghCIGIRAgpD2/PyZoNk4yk2QmZya5P9eVK5kzZ865k3NNnnnPOe/7mrsjIiJSXkrQAUREJDGpQIiISFgqECIiEpYKhIiIhKUCISIiYalAiIhIWCoQIiISlgqEiIiEpQIhIiJhqUCIiEhYaUEHqInWrVt7t27dgo4hIpJUZsyYkePubSKtl9QFolu3bkyfPj3oGCIiScXMvo1mPZ1iEhGRsFQgREQkLBUIEREJSwVCRETCUoEQEZGwVCBERCSspL7NVUQkGXS7/q2o1ls5dnick1RNUhYIMxsBjOjZs2fQUUSkDkrWf+ixlpQFwt0nA5MzMzMvDjqLiNSe7I3b6JDRmJQUi9k2N20rIGt9HsvW5/30PZZ+3FrANUf35tn/fssPW3YA0LpZQ3K3F1BY7Ez6zf4M6doqpvuMlaQsECJSP8T6k3y02wNIT02he+um9NytGT3bNCMnbwdT5q1l47bCn9Y58dHPGDmoAycMak/b5o3+5/VZ67cwYdpKXpmZzY6iEgD6tmvOhQd1Z+TgDjzyYRbjPs7ixlfnM/mKg2iQmniXhFUgRETK+eiaQ+nSqglp5f5p3zqyP9OW5fD67NW8t/B7Zq/axOxVm7jzrYUcsHtrRg7qQJsWDXnm85X8Z8kPP73uiL5tufCg7hyw+66YhVo/lx/RkzfmrGHxui089dkKLjlk91r9HaNh7h50hmrLzMx0jcUkUndtzi/k0mdn8MXyDTRNT+W4ge15eUY27TMa8f7Vh9KsYfU+4y5cs5kR46ZR4s6k3xzA3l12qfI2thcUs8fN70S1bkUtnP8sWc/5T31N4wapfHDNoXRs2bjKOarDzGa4e2ak9RKvTSMiAqzLzef0x77gi+UbaNO8IS9duj/3nLInAztmsDY3nwfeX1qt7RaXODe8MpfiEufc/bpWqzgANE5PrdbryjqsT1uOH9iO7YXF3PbGghpvL9Z0iklEEs7S77dw/pNfsSY3nx5tmvLMBUPp3KoJAGNOHsjIcdN46rMVnLRXRwZ0zKjStp/5fCVzsnNpn9GIPw7rW6OcsbiL6eYT+vPJkh94b+H3vL/we47ut1uNtxkrakGISEL5cvkGTv3756zJzWdI112YdNkBPxUHgAEdM7jgwO6UONzwyjyKS6I/TZ69cRv3vbcEgDtGDaj2KapYapfRiGuO6QPArW8sYFtBUcCJfqYCISIJ4625a/nVhK/YnF/EMf12458X7csuTdN/sd7VR/emQ0Yj5q3O5dkvVka1bXfn5tcXsK2gmOED23NUAn1SP3f/rvTv0ILVm7bz0IfLgo7zExUIEUkIT05bweUTZ1JQXMKv9uvK388ZQqMG4c/zN22Yxm2jBgBw37tLWJu7PeL235y7lo8Wr6d5ozRuGdEvptlrKi01hbtOGogZTPh0BUvWbQk6EqACISIBKylx7p6yiNvfXIg7XDusD7eP6k9qhM5wR/fbjWP778bWgmJujXCBd9O2Am6bHFrnT8fvQdsWjSpdPwiDO7fknH27UlTi3PTaPEqqcOosXoI/ASci9Uqkzmp/eWcJf3kndJ0g0kXgnf0S3l3wPe8tWMcx/duFXW/MlMXk5BUwtHsrzsjsXL3gteAPx/bh7fnr+HrlRl6ekc3p+wSbNSlbEGY2wszG5+bmBh1FRALUPqMxfzg2dIH3ljcWkLfjlxd4v/hmAy9NX0V6agp3nzQwpsN0xFpG4wb8+YQ9ALj77UX8uLUg0DzqKCciVbJ+cz7rt+ygb7vmv+hpHMmPWwt4+MNlPP/fbykqcRo3SOXSQ3twySE9aJJevRMaxSXOiY9+xrzVuVx4UHf+fMLP1xfyC4s57qFPWZGzlauP7s2VR/aq1j5qk7tzzoQv+SxrA6cN6cS9pw2K+T6i7SinU0wiUqloxy9aMeb4n4aRKC+/sJhnPl/JuI+z2JJfhBmckdmZa47pXePrAakpVmHfiHEfZbEiZyu92jbjskMTbyiLcMyMO0YNYNiDn/LvGdmcltmZod2DGcwvKU8xiUh8rd+cz7+mr+K3/5wR9Wv2vftDrn5pNpNmZLMuNx8IfRp+Y84ajrr/E8a8vZgt+UUc3Ks1U648mHtO3TNmF4sHdMzg1+X6Rixet5nHPvkGMxh7ykDS05Ln312PNs34zWGhgnbTa/MoKB3sr7bpFJOI4O7MWrWJDxd9z8eLf2Dh2s3/83yP1k05rE9bDuvThqHdW9GoQSprNm3ns6wcPsvKYVrWBnLydkS1r3jNobB1RxFH3/8Ja3LzuWn4Hrw5dy2zV23iV/t15Y4TB8Rln7FWldFma/J31CkmEalUtP+Mpv7xcLrs2uQXyzu0bMxpmZ05LbMz7k73G6bEOmKVNG2Yxu2jBnDRs9O5a8oi3KFdi0ZcO6xPoLmSmQqESD2TtX4Lb8xZG/X64YpDeWb2P59oC4pKmJO9iTWbtnPkHrvFfUiL8sVu54mRdZvzGXjrez8tT/QZ4MLly8nbQetmDQNIowIhUi98t2Ebk+euYXLp/AM77do0neMGtmPEnh3Yp1urmN0Cmp6Wwj7dEnOWtGQTVHEAFQiROqOgqITsjds44q+fRLX+cxcOZf8eu1b5VtVElOgtg2SlAiGShPILi3l5RjYL1uTy7YZtfLthG2tzt1OV0RkO7tUmfgGlTlCBEEkS0V5UPmD3Xem6axO6tGpa+r0JPds2q3DgO5GKqECI1DEvXLxf0BGkjlCBEElwm7YV8MD7S0lNMYpLnIzGDbjmmN6cNbRLnbh+IIlLBUIkQRUVlzDxq+/46/tL2bStkBQLTSzz+6N6h51ERyTWVCBEEtAX32zgtskLfroldf8eu3LLyH70bdci4GRSn6hAiAQsmovPXyzfoOIgtS4pT2BqPggRkfhLyhaEu08GJmdmZl4cdBaR6pqXncsLX31Lk/RUthUUA7BLkwacltmZCw7sRvuMxgEnlPouKQuESLLauqOIyXPW8MJX3zE3++cW8NDurTh73y4MG9COhmnqryCJQQVCJA6qMmzzB1cfQs+2zeOYRqR6kvIahEhdouIgiUotCJE4WDl2ODuKivnLO0uYMG0FAAf1bM1NJ+yhu5EkaahAiMTB8h/yuPLFWcxfvZm0FOOaY/pw6SE9YjactkhtUIEQiSF3Z9LM1dz8+ny2FRTTuVVjHh69F3t12SXoaCJVpgIhEiNb8gv582vzeW32GgBGDurAnScNoEWjBgEnE6keFQiRGJizahNXTJzFdz9uo3GDVG4f1Z9Th3TCTKeUJHmpQIjUQEmJ88Sny7n33SUUlTj92rfgkbP2Yvc2zYKOJlJjKhAiVRCpf8PCtZs58q+faApMqRPUD0JERMJSC0KkClaOHY67M/btxTw+dTkN01J45tdD2a/HrkFHE4k5tSBEqujRj7N4fOpy0lKMx84ZouIgdZYKhEgVPP3ZCu57bylm8MAZgzm8b9ugI4nEjQqESJRenpHNrZMXAjDmpIGMGNQh4EQi8aUCIRKFd+av5dqX5wBw0/A9GD20S8CJROJPBUIkgqlLf+CKibMocbjyyF5cdHCPoCOJ1AoVCJFKTF/5I5c8N53CYueCA7vx+6N6BR1JpNaoQIhUYP7qXC546mvyC0s4bUgn/jy8n4bOkHpFBUIkjKz1eZz75Fds2VHE8QPbMfaUPTVUt9Q76ignQuVDaEyZt44p86YAaAgNqVeSsgVhZiPMbHxubm7klUVEpFrM3St+0uzqKLax1d0fj12k6GVmZvr06dOD2LXUURO/+o4bXplHk/RU3rziIHpoVFapg8xshrtnRlovUgvij0AzoHklX9fULKpIYliybgu3vrEAgLtOGqDiIPVepGsQz7n77ZWtYGZNY5hHJBDbC4q5/IWZ7Cgq4dQhnThpr05BRxIJXKUtCHe/NtIGollHJNHd+sYClq3PY/c2Tbl9VP+g44gkhIgXqc2sr5kdaWbNyi0fFr9YIrXn9dmreWn6KtLTUhh31t40SdfNfSIQoUCY2ZXA68AVwHwzG1Xm6bvjGUykNqzM2cqfXpkHwM0n9GOP9i0CTiSSOCJ9VLoYGOLueWbWDXjZzLq5+0OAeg1JUttRVMzlE2eytaCY4we24+x9NQCfSFmRCkSqu+cBuPtKMzuMUJHoigqEJLl73l7C/NWb6bRLY8acvKeG0RApJ9I1iHVmNnjng9JicQLQGhgYz2Ai8fTBwu958rMVpKUY487am4zGDYKOJJJwIhWIc4F1ZRe4e5G7nwscErdUInG0ZtN2/lA6t8O1w/owuHPLgBOJJKZKTzG5e/bOn81sF6Bzmddsj2MukbgoKi7hqhdnsWlbIYf1acNFB2luB5GKRHU/n5ndAZwPfAPsHJvDgSPiE0skNiobhO8/S36gx580CJ9IRaK94ft0YHd3L4hnGBERSRzRFoj5QEtgfRyziMTczpZB7vZCjntwKmty87niiJ5cc0yfgJOJJL5oC8QYYJaZzQd27Fzo7iPjkkokhtydG1+dx5rcfAZ1bsmVR2raUJFoRFsgngHuAeYBJfGLIxJ7r85azZtz19IkPZWHzhhMg9SknAZFpNZFWyBy3P3huCYRiYNVP27j5tdDQ3jfOqI/3Vpr8GGRaEVbIGaY2RjgDf73FNPMuKQSiYGi4hL+76XZ5O0o4rgB7TgtU0N4i1RFtAVir9Lv+5VZpttcJaE9+vE3zPh2I+1aNGLMyQM1lIZIFUVVINz98HgHEYmlGd9u5OGPlmEG958+iJZN0oOOJJJ0Ig33fUKkDUSzjkht2pJfyP+9NIviEueSg3twQM/WQUcSSUqRWhD3mtlqKh+59W7gzdhFEqmZW99YyKoft9O/QwuuPqZ30HFEklakAvE9cH+EdZbFKItIjU2es4ZJM7Np1CCFh0YPpmFaatCRRJJWpMH6DqulHCI1tmbTdm58NTQ73I3D+9GzbfOAE4kkN/UYkjqhuMT5/Uuz2ZxfxJF923KOZocTqTEVCKkTxk9dzpcrfqR1s4bcc6pmhxOJhWj7QYgklIqG8c7J20HmnR/89FjDeItUX9QFwswGAP2ARjuXufuz8QglIiLBi3bCoFuAwwgViCnAccA0QAVCArFy7HCKS5zzn/qKT5flMKTrLrx4yX4aiE8khqJ9N50KHAmsc/cLgEFAw7ilEonCIx8t49NlObRqms64s/ZScRCJsWjfUdvdvQQoMrMWhCYO0mS+EpipS3/goQ9DQ2k8NHow7TMaBx1JpM6J9hrEdDNrCTwBzADygK/ilkqkEms2beeqF2fhDlcf3ZuDe7UJOpJInRTtYH2/Lf3xMTN7B2jh7nPjF0skvIKiEn73wkw2bivkkN5tuPzwnkFHEqmzojrFZCHnmNnN7r4S2GRmQ+MbTeSXxry9iFnfbaJDRiMePGMwKSnq7yASL9Feg/gbsD9wZunjLcCjsQxiZilmdpeZPWJm58Vy21I3vDV3LU99tpIGqca4s/emVVMN4S0ST9EWiH3d/XdAPoC7bwQivjvN7EkzW29m88stH2ZmS8wsy8yuL108CugIFALZUf8GUi9880Me1748B4Abj9+DvbvsEnAikbov2gJRaGaphGaRw8zaACVRvO5pYFjZBaXbeZRQX4p+wJlm1g/oA3zh7lcDv4kyl9QD2wuK+e3zM9laUMzwPdtz3gHdgo4kUi9EWyAeBl4F2prZXYQ6yd0d6UXuPhX4sdzioUCWuy939wLgRUKth2xgY+k6xVHmkjrO3bnxtXks+X4LPdo05Z5TNM6SSG2J9i6mf5rZDEKd5Qw40d0XVXOfHYFVZR5nA/sCDwGPmNnBwNSKXmxmlwCXAHTpohE767qXvl7FKzNX06hBCn8/ewjNGmr4MJHaEvHdZmYpwFx3HwAsjsE+w338c3ffBlwY6cXuPh4YD5CZmekxyCMJpKJB+PILSzj2wZ8/N2gQPpH4i3iKqbQH9Rwzi9XH9Wygc5nHnYA1Mdq2iIjESLTt9fbAAjP7Cti6c6G7j6zGPr8GeplZd2A1MBo4qxrbkTpo5djhPPXZCm6bvBCASw/pwfXH9dV1B5EARFsgbqvOxs1sIqFRYFubWTZwi7tPMLPLgXeBVOBJd19Qne1L3eLu3P/+Uh75KAuAG47ry6WH7h5wKpH6K9qL1J+UfWxmBxL61P9J+Ff89LozK1g+hdCw4dViZiOAET17apiFuqK4xPnz6/N54cvvSE0xxpw8kNMzO0d+oYjETdTjI5vZYDP7i5mtBO4EqnsXU425+2R3vyQjIyOoCBJDO4qKuWLiTF748jvS01J47JwhKg4iCaDSFoSZ9SZ0jeBMYAPwEmDufngtZJN6IG9HEZc+N53PsjbQvGEa/zgvk3177Bp0LBEh8immxcCnwAh3zwIws9/HPZXUCxvydnDB018zNzuX1s0a8syv96F/B7UKRRJFpFNMpwDrgI/N7Akz29lRTqRGVm/azmmPf8Hc7Fw6t2rMpN/sr+IgkmDMPXJfMzNrCpxI6FTTEcAzwKvu/l5841UuMzPTp0+fHmQEiVJFHeDKUwc4kfgzsxnunhlpvaguUrv7Vnf/p7ufQKhj22zg+ggvExGRJBZVCyLRlLnN9eJly5YFHUeilLV+C6PHf0lO3g4O2H1XJpy3D43TU4OOJVLvxLQFkWh0m2vyWfb9FkaP/y85eTs4qGdrFQeRJKChMSXulqzbwllP/JcNWws4uFdrnjg3k0YNVBxEEl1VOsp1NbOjSn9ubGbN4xdL6orF6zb/VBwO6d1GxUEkiURVIMzsYuBl4PHSRZ2A1+IVSuqGRWs3c9YTX7JhawGH9m7D+F8NUXEQSSLRtiB+BxwIbAZw92VA23iFkuS3cE2o5fDj1gIO79OGx1UcRJJOtNcgdrh7wc4hl80sjdL5qUXKW7Aml7P/8SWbthVyRN+2/P2cvWmYpuIgkmyibUF8YmZ/Ahqb2dHAv4HJ8YtVOTMbYWbjc3Nzg4ogFZi/OpeznggVh6P2UHEQSWbR9qROITQd6DGEhtp4F/iHB9yJQj2pgxVt72hQD2mRRBJtP4hoTzGNAp519ydqFktERJJFtAViJPCgmU0FXgTedfei+MWSZLBy7HA+y8rhuklzyd64nbQU44ojevHbw3enQWpS9sEUkTKinVHuAjNrABxHaCa5v5nZ++5+UVzTScLakl/ImLcX88KX3wHQv0ML7j11EP06tAg4mYjEStQ9qd290MzeJnT3UmNCp51UIOqhT5f9wPWT5rF603YapBpXHtGLyw5Tq0GkromqQJjZMEIzyx0O/Af4B3B6/GJJItqSX8jdUxYx8atVAAzo2IL7ThtE33ZqNYjURdG2IM4ndO3hUnffEb84kqg+WfoDN0yay5rcfNJTU7jqqF5cckgPtRpE6rBor0GMjncQSRyRbl8tKC7h3neX8LvDe9ZSIhEJQqUf/8xsWun3LWa2uczXFjPbXDsRw+ZSRzkRkThLygmDdlJHufhxd8a+vZjHpy4nLcV49Oy9ObZ/u6BjiUgMxHTCIDN7LpplUjeoOIgIRD8WU/+yD0oH6xsS+zgSNBUHEdkp0jWIG8xsC7Bn2esPwPfA67WSUGqNioOIlFVpgXD3Me7eHLjX3VuUfjV3913d/YZayii1QMVBRMqL9hTTV2aWsfOBmbU0sxPjlElqmYqDiIQTbYG4xd1/uqfU3TcBt8QnktQmFQcRqUi0BSLcelGP4ySJScVBRCoT7T/56WZ2P/AoocH6rgBmxC2VxEVlPaSLSpxLnwsdUk3uIyIQfQviCqAAeAn4F7Ad+F28QkWintQiIvFXpZ7UZtbM3fPimKdK1JO66tyd+95bwqMff0NaijHurL0ZNkCnlUTqk1j3pD7AzBYCC0sfDzKzv9UwowTggQ+W8ejH35CaYjxy5l4qDiJSoWhPMT0AHAtsAHD3OcAh8Qol8fHQB8t4+MNlpKYYD40ezHED2wcdSUQSWNSD+bv7qnKLimOcReJo3EfLeOCDpaQYPHDGYE7Ys0PQkUQkwUV7F9MqMzsAcDNLB64EFsUvlsTS3/6TxX3vhYrD/acPZuQgFQcRiSzaFsRlhO5a6ghkA4MJ8C4mid7jn3zDX95Zghnce+ogTtyrY9CRRCRJVNqCMLN73P064HB3P7uWMkmM/OPT5Yx5ezFmcM8pe3LKkE5BRxKRJBKpBXG8mTUANDBfkpkwbQV3vhU6Czj25IGcntk54EQikmwiXYN4B8gBmpZOMWqEelIb4O7eIs75JAqR5pC+btI8rps0Tz2kRaRKIrUgbnL3DOCtMkN9//S9NgKKiEgwIrUgvgD2BjbXQhappnm3HsP97y/lmc9XUuLQpnlDbh3Rn+F7qp+DiFRfpAKRbmbnAQeY2cnln3T3V+ITS6Lh7kyeu5Y731zI+i07SDG44MBu/P7o3rRo1CDoeCKS5CIViMuAs4GWwIhyzzkQSIEwsxHAiJ49ewax+4TwzQ953Pz6fD7L2gDAXl1acueJA+jfISPCK0VEohPVYH1mdqG7T6iFPFVSHwfr215QzLiPlzF+6nIKi52WTRpw/bC+nJ7ZmZQUCzqeiCSBmAzWZ2bXArj7BDM7rdxzd9csolTVx0vWc9T9n/Dox99QWOyM3qczH11zGKOHdlFxEJGYq7QFYWYz3X3v8j+HexyEut6CiHT76k66fVVEqiJWw31bBT+HeywiInVIpIvUXsHP4R5LjK0cOxx355LnZvD+wu/J7LoLL16yH2mpUQ/CKyJSbZEKxKAyPagbl/5M6eNGcU0mADz/5Xe8v/B7mjdK48HRg1UcRKTWVFog3D21toLILy1Zt4U731wIwJiTB9JplyYBJxKR+kQfRxNUfmExV0ycyY6iEs7I7KwJfkSk1qlAJKi73lrE0u/z6NGmKbeM7Bd0HBGph1QgEtC7C9bx3H+/JT01hYdH70WT9Ggn/hMRiR0ViASzNnc7102aC8C1w/owoKOGzhCRYKhAJJDiEuf3L81m07ZCDuvThl8f2D3oSCJSj6lAJJDHPvmG/y7/kdbNGnLfaYM0fIaIBEoFIkHM+HYj97+/FID7Tx9E62YNA04kIvWdCkQC2JxfyFUvzqK4xLn44O4c0rtN0JFERJKzQJjZCDMbn5ubG3SUGnN3bnx1PtkbtzOgYwv+eGzfoCOJiABRzgeRqJJ1NFeN0ioiQYrVaK4iIlJPqQdWAFaOHc6aTdsZOW4aOXkFXHxwd24crt7SIpJY1IIIQH5hMZc+N4OcvAIO6tma64bpuoOIJB4ViFrm7tzwyjzmrc6lc6vGPHLmXhrCW0QSkv4z1bIJ01bw6qzVNG6QyhPnZrJL0/SgI4mIhKUCUYumLcvh7imLAPjr6YPo265FwIlERCqmAlFLvtuwjcsnzqTE4XeH787xA9sHHUlEpFIqELVgW0ERlzw3nU3bCjmib1uuPrpP0JFERCJSgYgzd+eP/57L4nVb6NG6KQ+cMZhUDcInIklABSLO/vafb3hr3lqaNUxj/LlDyGjcIOhIIiJRUYGIo48Xr+e+95YA8OAZg+nZtnnAiUREoqee1DFU2RhLFz3785hRGmNJRJKBWhAiIhKWWhAxtHLscDZuLeC0x78ga30eQ7ruwvMX7kvj9NSgo4mIVJlaEDG0raCIC57+mqz1efTerRkTzstUcRCRpKUCESMFRSVc9vxMZq/aRMeWjXn21/vSsomG0RCR5KUCEQMlJc4f/j2HqUt/oFXTdJ67cCjtMhoFHUtEpEZUIGrI3bn9zYW8MWcNTdNTefqCfejRplnQsUREakwFooYe/TiLpz9fSXpqCuPPzWTPTi2DjiQiEhMqEDXwwpffcd97SzGDB0cP5sCerYOOJCISM0lZIMxshJmNz83NDSzD2/PWctNr8wC4Y9QAjc4qInVOUhYId5/s7pdkZGQEsv/Ps3K46sXZlDhcfXRvztmvayA5RETiydw96AzVlpmZ6dOnT4+8Yg1UNnxGeRpCQ0SSgZnNcPfMSOslZQtCRETiT0NtRFC2VTDuo2Xc995Suu7ahHeuOkS9pEWkTlMLIkpZ6/N4+MMsAMacNFDFQUTqPBWIKJSUODe8MpeC4hJOz+zEAbqdVUTqARWIKPzzq+/4euVGWjdryI3H9ws6johIrVCBiGBt7nbueXsxALeP6k9GE00ZKiL1gwpEJdydP782n7wdRRzTbzeOG9Au6EgiIrVGBaISb81byweL1tO8YRq3jxqAmQUdSUSk1qhAVGDj1gJufWMBANcf31fDd4tIvaMCUYG7piwiJ6+Aod1bceY+XYKOIyJS61Qgwpi2LIeXZ2STnpbC2JMHkpKiU0siUv+oQJSzraCIG16dC8BVR/bS5D8iUm+pQJTzwPtLWfXjdvZo34JLDukRdBwRkcCoQJQxZ9UmJkxbQYrBPacMpEGq/jwiUn/pP2CpwuISrps0lxKHCw/qrqlDRaTeq7ejuVY2z8MTn67giU9XAJrjQUTqL7UgREQkrHrbgijfMigucT7/JoeDe7UJKJGISGJRC6JUaoqpOIiIlKECISIiYalAiIhIWCoQIiISlgqEiIiEpQIhIiJhqUCIiEhYKhAiIhKWuXvQGarNzHKBZeUWZwC5YVYvv7w1kBOnaJFUlDHe24l2/UjrVfZ8tH//ipYFdVyCOiZVeU11j0tNl+u9Uv31EvW90tXdI3f8cvek/QLGR7Ms3HJgeiLlro3tRLt+pPUqez7av38lywI5LkEdk9o4LjVdrvdK7I9JVY9LUO+VZD/FNDnKZZUtD0KsslR1O9GuH2m9yp6vyt9fx6Rqr6nucYnV8iDovRLdfuIiqU8x1YSZTXf3zKBzyP/ScUk8OiaJqTaOS7K3IGpifNABJCwdl8SjY5KY4n5c6m0LQkREKlefWxCZN/jYAAAE/0lEQVQiIlIJFQgREQlLBUJERMJSgQjDzHqY2QQzeznoLPWZmTU1s2fM7AkzOzvoPBKi90fiMbMTS98nr5vZMbHabp0rEGb2pJmtN7P55ZYPM7MlZpZlZtdXtg13X+7uF8Y3af1UxeNzMvCyu18MjKz1sPVIVY6L3h+1o4rH5LXS98n5wBmxylDnCgTwNDCs7AIzSwUeBY4D+gFnmlk/MxtoZm+W+2pb+5HrlaeJ8vgAnYBVpasV12LG+uhpoj8uUjuepurH5KbS52MiLVYbShTuPtXMupVbPBTIcvflAGb2IjDK3ccAJ9RuwvqtKscHyCZUJGZTNz/MJIwqHpeFtZuufqrKMTGzRcBY4G13nxmrDPXlTdeRnz+JQugfT8eKVjazXc3sMWAvM7sh3uGkwuPzCnCKmf2dxBr+ob4Ie1z0/ghURe+VK4CjgFPN7LJY7azOtSAqYGGWVdhD0N03ADH7I0tEYY+Pu28FLqjtMPKTio6L3h/BqeiYPAw8HOud1ZcWRDbQuczjTsCagLLIL+n4JCYdl8RTq8ekvhSIr4FeZtbdzNKB0cAbAWeSn+n4JCYdl8RTq8ekzhUIM5sIfAH0MbNsM7vQ3YuAy4F3gUXAv9x9QZA56ysdn8Sk45J4EuGYaLA+EREJq861IEREJDZUIEREJCwVCBERCUsFQkREwlKBEBGRsFQgREQkLBUIqbPMrNjMZpf5qnSY99piZivNbJ6ZZZrZq6XZsswst0zWAyp47UVm9ly5ZbuVDgvdwMxeMrMfzezE2vltpC5TPwips8wsz92bxXibaaWdlWqyjZVAprvnlFl2GPAHd690dGEz2wVYBnRy9/zSZZcDA9390tLHzxOaR+O1muQUUQtC6p3ST/C3mdnM0k/yfUuXNy2dpOVrM5tlZqNKl59vZv82s8nAe2aWYmZ/M7MFpXOITDGzU83sSDN7tcx+jjazV2qQcx8z+8TMZpjZ22a2m7tvBD4HhpdZdTQwsbr7EamICoTUZY3LnWIqO9NWjrvvDfwd+EPpshuBj9x9H+Bw4F4za1r63P7Aee5+BKGZ7roBA4GLSp8D+AjYw8zalD6+AHiqOsHNrCHwEHCKuw8BngfuKH16IqGigJl1Ls0ytTr7EalMfRnuW+qn7e4+uILndn6yn0HoHz7AMcBIM9tZMBoBXUp/ft/dfyz9+SDg3+5eAqwzs48hNOZy6fWBc8zsKUKF49xqZt8D6A98YGYAqYRG8oTQ4GwPm1kzQtNL/qs0i0hMqUBIfbWj9HsxP78PjNAn9iVlVzSzfYGtZRdVst2nCE1ulE+oiFT3eoUBc9394PJPuPtWM/uA0Oxuo4HfVHMfIpXSKSaRn70LXGGlH9nNbK8K1ptGaKa7FDPbDThs5xPuvobQ+Pw3EZpTuLoWEpq9bWhplnQz61/m+YnAH4GW7v51DfYjUiEVCKnLyl+DGBth/TuABsBcM5vPz+f8y5tE6HTPfOBx4Esgt8zz/wRWuXu152529x3AqcD9ZjYHmAXsW2aVdwid/nqxuvsQiUS3uYpUg5k1c/c8M9sV+Ao40N3XlT43Dpjl7hMqeO1Kyt3mGuNsus1VYkItCJHqedPMZgOfAneUKQ4zgD0J3XVUkR+AD80sM9ahzOwl4EBC10BEakQtCBERCUstCBERCUsFQkREwlKBEBGRsFQgREQkLBUIEREJSwVCRETC+n8mfp7qOb8Y3QAAAABJRU5ErkJggg==\n", "text/plain": [ "
" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "aeff = cta_irf.aeff.to_effective_area_table(offset=offset, energy=energy)\n", "aeff.plot()\n", "plt.loglog()\n", "print(cta_irf.aeff.data)" ] }, { "cell_type": "code", "execution_count": 7, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "NDDataArray summary info\n", "e_true : size = 30, min = 0.112 TeV, max = 89.125 TeV\n", "e_reco : size = 30, min = 0.112 TeV, max = 89.125 TeV\n", "Data : size = 900, min = 0.000, max = 0.926\n", "\n" ] }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAZcAAAEUCAYAAADnQnt7AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMi4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvIxREBQAAFNFJREFUeJzt3XuwXWV9xvHnIRCIARKBUCUBSQIGQ4KEBPDGSLloqIabWgmXFkSiIzjtVK2oVNEZr6XWUrA0AsZamoCASJgwtFq5KZUg3hJDbKAgh0tjCCcFBELCr3/sjWwO52S/6+Tde629zvczk8nZa//2Oj99Z/PkXZd3OSIEAEBO25TdAACgfggXAEB2hAsAIDvCBQCQHeECAMiOcAEAZEe4AACyI1wAANn1XLjYPt72N2x/z/bbyu4HAPBylQgX25fbXmt7xYDtc22vtr3G9rmSFBHXRcRZkk6X9N4S2gUAtFGJcJG0SNLc1g22R0m6WNIxkqZLmm97ekvJec33AQAVU4lwiYhbJa0fsPkQSWsi4r6I2ChpiaTj3PBlSTdGxN3d7hUA0N62ZTewBRMlPdjyuk/SoZI+LOkoSeNs7xMRlwz2YdsLJC2QpLFjx86ett++HW4XALpr0/ObC9Wv6Pt9cu3zv7t3XURMKNrTC6ocLh5kW0TEhZIubPfhiFgoaaEkzZ4zK370k5vzdgcAJVv/bH+h+td+PP1gz1MXnfhA0X5aVeKw2BD6JO3Z8nqSpIdL6gUAUECVw2W5pH1tT7Y9WtJJkq4vsgPb82wv7O/f0JEGAQCDq0S42F4s6Q5J02z32T4zIjZJOkfSTZJWSboqIlYW2W9ELI2IBePHj8vfNABgSJU45xIR84fYvkzSsi63AwDYSpWYuQAA6qUSM5dOsT1P0rwpUyeX3QoAJHngyd8m1x541o8L7fvTn9ovufaTFxXa9cvUeubCORcAKEetwwUAUA7CBQCQXa3DhftcAKActQ4XzrkAQDlqHS4AgHIQLgCA7Gp9nwsAVMHjG9PP+x78N2uSa7/4mWKPEtlnp42F6rdGrWcunNAHgHLUOlw4oQ8A5ah1uAAAykG4AACyI1wAANkRLgCA7Gp9KTJL7gPolPufeCC5dtYJ1yXXfuDzb0+uPXXqnsm1krTjdmML1W+NWs9cuFoMAMpR63ABAJSDcAEAZEe4AACyI1wAANkRLgCA7LgUGQCGYdZ7bkiu/eGVRyTX7rJ9+uXC3by0uKhaz1y4FBkAylHrcAEAlINwAQBkR7gAALIjXAAA2REuAIDsan0pMgAUWr34jFuTa//zisOSaw/cdWZybV0wcwEAZFfrcLE9z/bC/v4NZbcCACNKrcOFmygBoBy1DhcAQDkIFwBAdoQLACA7wgUAkB33uQDoOY/+fm1y7axTv59ce8FXXp++310PSK4diZi5AACyI1wAANkRLgCA7AgXAEB2hAsAILtahwtriwFAOWp9KXJELJW0dPacWWeV3QuAoT301MOF6meceG1y7VWXH5xce/TEOYX6wNBqPXMBAJSDcAEAZEe4AACyI1wAANkRLgCA7AgXAEB2tb4UGUB5HnrqkeTaGe9dWmjf37/isOTaaeP2KbRv5MHMBQCQHeECAMiOcAEAZEe4AACyI1wAANkRLgCA7LgUGUCyIqsXzzj+6uTaZQUuLZakfXaenFy743ZjC+0beTBzAQBkR7gAALLruXCxPcX2ZbbT59wAgK6qRLjYvtz2WtsrBmyfa3u17TW2z5WkiLgvIs4sp1MAQIpKhIukRZLmtm6wPUrSxZKOkTRd0nzb07vfGgCgqEqES0TcKmn9gM2HSFrTnKlslLRE0nFdbw4AUFglwmUIEyU92PK6T9JE27vavkTSLNufGOrDthfYvsv2Xb/73WOd7hUA0KLK97l4kG0REY9J+mC7D0fEQkkLJWn2nFmRuTegNgrdu3Lcd5Jrb7nqyOTa3cfsnlwrSeNG71yoHt1X5ZlLn6Q9W15PkpT+LQAAlKbK4bJc0r62J9seLekkSdeX3BMAIEElwsX2Ykl3SJpmu8/2mRGxSdI5km6StErSVRGxsuB+59le2N+/IX/TAIAhVeKcS0TMH2L7MknLtmK/SyUtnT1n1lnD3QcAoLhKzFwAAPVS63DhsBgAlKMSh8U6hcNiGKn+54n7k2sPevcNybW3f+fo5Nr9X8mCGiNZrWcuAIBytJ252N4lYT/PR0R/hn4AADWQcljs4eafwe6Yf8EoSXtl6Sgj2/MkzZsyNf2pdQCArZcSLqsiYtaWCmz/LFM/WXHOBQDKkXLO5Y2ZagAAI0RKuFxg+81bKoiIZzL1AwCogZTDYv+tRsC8WtKVkhZHxM872xaAgX6x/lfJtYfPvyW59pbFhyfXTtl5SnItRra2M5eI+IeIeKOkt6rxQK9v2l5l+9O2X9vxDrcCN1ECQDmS73OJiAci4svNk/snSzpBjQUlKysilkbEgvHjx5XdCgCMKMnhYnu75kzgCkk3SvqNpHd1rDMAQM9KuYnyaEnzJb1D0p1qPMt+QUQ81eHeAAA9KuWE/icl/Zukj0bE+g73AwCogbbhEhF/LEluOFXSlIj4nO29JL0qIu7sdJMAgN5SZFXkr0t6XtIRkj4n6QlJ10g6uAN9ZcHyL6iy1f2rC9UffsptybW3LzkiufY1O74muXbMqB2SazGyFVkV+dCIOFvSM5IUEY9LGt2RrjLhajEAKEeRcHnO9ihJIUm2J6gxkwEA4CWKhMuFkr4raXfbn5d0u6QvdKQrAEBPS7kUeduI2BQRV9j+qaQj1Vh+//iIqPRNlACAcqSc0L9T0kGSFBH3SLqnox0BAHpeymGxLT0kDACAl0mZuUyw/VdDvRkRX83YT1ZcigwA5UgJl1GSdlQPzmB4EiW6bfWG3yTXvmH+Dwrte/mVb0+unTR2UnLtDqO2L9QHkCIlXB6JiM91vBMAQG1wzgUAkF1KuBzZ8S4AALWSEi7fb1dg++4MvQAAaiLlnMvrbP9yC+9bEot3AQD+ICVc9kuo2by1jQAA6iPleS4PdKMRoKpWPr4yufYt829Orr198eGF+pj4ij2Sa7m8GGUrsnBlz7E9z/bC/v4NZbcCACNKrcOF57kAQDlqHS4AgHIMK1xsH2v7KttLbB+XuykAQG8b7szlnRHxpxFxkqS5ORsCAPS+lEuRBzPG9l7Nn8fmagYAUA/DDZfzJX24+TOLWqLn3L1uS/cFv9SRp9+RXHvHkvTVkqbsPDW5VpJGb7NdoXqgTMMNlz+KiI9Jku03SFqTryUAQK8b7jmXE1p+PjZHIwCA+hj2zMX2VEkhKf22YQDAiDDccDlP0tnNnz+TqRcAQE1szWGxXSLi45L+LGM/AIAaGG64TJX0YPPnnTL1AgCoieEeFgs17nWZoQqfc7E9T9K8KVMnl90KuuAX63+VXHvkGf+VXPvDb70puXa/8SlPqADqb7gzl79T4yFhp0n6RL528mLhSgAox7BmLhHxW0nnSpLtd+jFQ2QAABQPF9t/LelASTeosa7YT3I3BQDobcOZubwuIk62fbukoyLimdxNAQB623DOuexm+08krZN0RPNnAAD+oG242N5/wKarJU2Q9N3m3xM60BcAoIelHBb7tqSDJMn2+yPi0hfesP2KiPh9p5oDAPSmlHBxy88fknRpy+vbJM3O2hHQYuXjv06uPfzP05fGv63AvSszdpmRXLvp+U3JtdtuM9zbzIDqSznnEi0/e8B7w71PBgBQYyn/dHqV7dMl/UIvD5d4eTkAYKRLCZfzJc2RdIakSbZXSrqn+We3zrUGAOhVbcMlIha2vrY9SdIBkmZKurVDfQEAeljhM4oR0SepT9Ky/O0AAOqAE/IAgOy4FhJddU//PYXq33JK+pHXm//1zcm1nVoan8uLgQZmLgCA7AgXAEB2PTeHtz1W0tclbZR0c0RcUXJLAIABKjFzsX257bW2VwzYPtf2attrbJ/b3HyipKsj4ixJx3a9WQBAW5UIF0mL1Hjw2B/YHiXpYknHSJouab7t6ZIm6cUnX27uYo8AgESVCJeIuFXS+gGbD5G0JiLui4iNkpZIOk6Ne2wmNWsq0T8A4KWqfM5lol6coUiNUDlU0oWSLrL9DklLh/qw7QWSFkjSnnvt2cE262nj5o3Jtfc+cW9y7ZtOu61QH7dfcVhy7b7jXptcyyXDQGdV+Rs2cJFMSYqIeEqNdc62qLlszUJJmj1nFgtsAkAXVfmwUp+k1inHJEkPl9QLAKCAKofLckn72p5se7SkkyRdX2QHtufZXtjfv6EjDQIABleJcLG9WNIdkqbZ7rN9ZkRsknSOpJskrZJ0VUSsLLLfiFgaEQvGjx+Xv2kAwJAqcc4lIuYPsX2ZWH0ZAHpOJWYuAIB6qcTMpVNsz5M0b8rUyWW3UgnPPf9ccu2vHl+VXHvUgp8l1/7H5Ycm10rStHHTkmu5vBiojlrPXDjnAgDlqHW4AADKQbgAALKrdbhwnwsAlKPW4cI5FwAoR63DBQBQDsIFAJAdNwZUUJH7UZ7Z/Gxy7V3r0u9dOfEvVifXXnPR/sm1B+w6PblW4t4VoFfVeubCCX0AKEetw4UT+gBQjlqHCwCgHIQLACA7wgUAkB3hAgDIrtbXefbqkvubYnNy7U/X3ZNce+JH70uuveZr6UvdH7HHwcm1AEaGWs9cuFoMAMpR63ABAJSDcAEAZEe4AACyI1wAANkRLgCA7LgUuQue3vxMofqb+n6aXHvGZ9cl1157wZTk2jftfkByLQAMVOuZC5ciA0A5ah0uAIByEC4AgOwIFwBAdoQLACA7wgUAkF2tL0XupCeeezK59kf/u6LQvt//988m11726QnJtW+YMCO5dvtR2yfXAsBAzFwAANnVOlxsz7O9sL9/Q9mtAMCIUutw4SZKAChHrcMFAFAOwgUAkB3hAgDIjnABAGTHfS4t+jf+X3LtbY+m37tyxleLLbn/lQ/tlFz7tonTk2vHbDumUB8AMFzMXAAA2REuAIDsCBcAQHaECwAgO8IFAJBdrcOFtcUAoBy1vhQ5IpZKWnrg7Nef9eRzT7WtL3J58dmXbE6uPX/BK5NrJemEvack1+6w7Q6F9g0A3VDrmQsAoByECwAgO8IFAJAd4QIAyI5wAQBkR7gAALKr9aXIL3j82ad1zf0r29Z96tJNyfs85+T0y4vnT52UXCtJO2+3Y3LtKI8qtG8A6AZmLgCA7AgXAEB2hAsAIDvCBQCQHeECAMiOcAEAZDciLkV+dL30hSXtVzH+yGk7J+/zlH1enVw7brudkmslLi8G0PuYuQAAsiNcAADZ9Vy42J5i+zLbV5fdCwBgcF0NF9uX215re8WA7XNtr7a9xva5W9pHRNwXEWd2tlMAwNbo9gn9RZIukvQvL2ywPUrSxZKOltQnabnt6yWNkvTFAZ9/X0Ss7U6rAIDh6mq4RMSttvcesPkQSWsi4j5Jsr1E0nER8UVJ7+xmfwCAPBwR3f2FjXC5ISJmNF+/W9LciHh/8/Vpkg6NiHOG+Pyukj6vxkzn0mYIDVa3QNKC5ssZklYMVtdF4yRtqMD+inyuXe1w3y+yfTdJ67bwO7ol5/hVYeza1QznvaqOXx2/e+1qcnz3pkVEsfsoWkVEV/9I2lvSipbX71EjJF54fZqkf8z8O+/q9v/OQXpYWIX9Fflcu9rhvl9kexXGLvf4VWHs2tUM572qjl8dv3s5x6hTY1eFq8X6JO3Z8nqSpIdL6qWTllZkf0U+1652uO8X3V4FOXurwti1qxnOe1Udvzp+99rVlP7dq8JhsW0l/UbSkZIekrRc0skR0f7pXum/866ImJNrf+gexq63MX69a2vHrtuXIi+WdIekabb7bJ8ZEZsknSPpJkmrJF2VM1iaFmbeH7qHsettjF/v2qqx6/rMBQBQf1U45wIAqBnCBQCQHeECAMhuRIcLi2D2HttjbX/L9jdsn1J2P0jH96232T6++b37nu23tavv2XBhEcz6KDiWJ0q6OiLOknRs15vFSxQZO75v1VNw/K5rfu9Ol/Tedvvu2XBRYxHMua0bWhbBPEbSdEnzbU+3PdP2DQP+7N79ljGERUocSzVusn2wWdb+8aLotEVKHztUzyIVH7/zmu9vUc8+5jhYBLM2ioylGis6TJL0c/X2P45qoeDY/bq73aGdIuNne5WkL0m6MSLubrfvun05J+rFf9VKjf8QTRyq2Pauti+RNMv2JzrdHAoZaiyvlfQu2/+k6i43MtINOnZ833rGUN+9D0s6StK7bX+w3U56duYyBA+ybci7RCPiMUlt/09CKQYdy4h4StIZ3W4GhQw1dnzfesNQ43ehpAtTd1K3mctIWQRzJGAsexdj19uyjF/dwmW5pH1tT7Y9WtJJkq4vuScMD2PZuxi73pZl/Ho2XEpcBBOZMZa9i7HrbZ0cPxauBABk17MzFwBAdREuAIDsCBcAQHaECwAgO8IFAJAd4QIAyI5wAQBkR7gAALIjXIAh2P6A7Uds/7zlz8yM+9/b9tPN/e7a8jsetf1Qy+vRQ3z+ZttvH7DtL21/3faY5mc32t4tV89AqrqtigzkdICk8yLisg7+jnsj4sDmzwdKku3zJT0ZERe0+exiNdZ9uqll20mSPhYRT0s60Pb9edsF0jBzAYY2U42HkpXO9qm272zORv65+bTAqyW90/b2zZq9Je0h6fbyOgUaCBdgaPtL+mbL4akFZTRh+3VqPLP8zc1ZzmZJpzSfj3KnXnxM7UmSrgwWDEQFcFgMGITtPSWtjYgDWraNaT5JcQ9Jr5S0UtLfRsS9treJiOc71M6RkmZLWm5bksZIWtt874VDY99r/v2+DvUAFEK4AIM7QNI9rRua5zE+aPtwSTMi4iLbp9v+rKS7bPdLWhcRNzSfO/5xSR9R48l+90bE14bZiyV9KyIGezTwdZK+avsgSWNSnm0OdAOHxYDBzdSAcNmCG4cIjg9JelrSY839DdcP1Hhu+e6SZHsX26+RpIh4UtLNki5XYxYDVAIzF2BwMyW91fYxzdch6bDmf8wH2tD8+1m9+J0aq8Y/3r4dEb/cmkYi4te2z5P077a3kfScpLMlPdAsWSzpWjUOiwGVQLgAg4iIU4bxsVskfcX2ZEnjJV0k6Qu2H5H0RER8NvF3nz/ItislXTlE/XfVOHQGVAZPogRK0rxo4MeSHmu51yXXvseo8fjaCZJmRsT6nPsH2iFcAADZcUIfAJAd4QIAyI5wAQBkR7gAALIjXAAA2REuAIDsCBcAQHaECwAgu/8H1+4m8eDJV8kAAAAASUVORK5CYII=\n", "text/plain": [ "
" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "edisp = cta_irf.edisp.to_energy_dispersion(\n", " offset=offset, e_true=energy, e_reco=energy\n", ")\n", "edisp.plot_matrix()\n", "print(edisp.data)" ] }, { "cell_type": "code", "execution_count": 8, "metadata": {}, "outputs": [], "source": [ "# Simulate data\n", "sim = SpectrumSimulation(\n", " aeff=aeff, edisp=edisp, source_model=model, livetime=livetime\n", ")\n", "sim.simulate_obs(seed=42, obs_id=0)" ] }, { "cell_type": "code", "execution_count": 9, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "*** Observation summary report ***\n", "Observation Id: 0\n", "Livetime: 1.000 h\n", "On events: 411\n", "Off events: 0\n", "Alpha: 1.000\n", "Bkg events in On region: 0.00\n", "Excess: 411.00\n", "Excess / Background: inf\n", "Gamma rate: 411.00 1 / h\n", "Bkg rate: 0.00 1 / min\n", "Sigma: nan\n", "energy range: 0.10 TeV - 100.00 TeV\n" ] }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAmcAAAJkCAYAAABUJoa/AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMi4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvIxREBQAAIABJREFUeJzs3Xl8lPW5///XRdj3JYAIgYAidUeluLRWqtaKiHpq677+bLHVntZTbV2+9njaaqWtbTme1oVWC9pWpVrrilatoBY3UNxFtigRZN8DZLt+f8w9IYRJMklm7nvmnvfz8ZhHMp+5lytDPsOVz/25P5e5OyIiIiKSG9pFHYCIiIiI7KTkTERERCSHKDkTERERySFKzkRERERyiJIzERERkRyi5ExEREQkhyg5ExGRgmBmN5rZGjP7LHj+H2a2zMy2mNkhGTzP0Wa2IFPHk8Kj5EzaxMzOMbO5wYfbCjObaWZfzPI53cz2zuY5RCT/mFmZmW0LPo+Sj98Fr5UAVwL7ufsewS63AN919+7u/mYbzrvLZ5K7v+juo9ryszRzvouCc56RrXNItJScSauZ2Q+AKcDPgYHAUOA24NQo4xKRgjYxSLaSj+8G7cOAte6+qt62w4D3wg+xzS4E1gVfG2Vm7cMJRzJNyZm0ipn1An4KXO7uf3f3re5e5e6PufsPzayTmU0xs+XBY4qZdQr2vcjMXmpwvLq/PM1smpn93syeMLPNZvaqme0VvPZCsMtbwV/FZ5pZsZk9bmYbzGydmb1oZvrdFhEAzOx44Blgz+Bz4z4z2wIUkfgsWRxst6eZPWRmq81sqZl9r94xiszsOjNbHHwuzTOzkkY+k8aZWXmw3zVm9mCDeP7XzG4Nvu9lZncFVx4+DS69FjXxswwDjgEmAV81s4H1XhtnZuVmdnVw6fZPQfvJZjY/+IycY2YH1dvnmno/0/tm9h9teKslQ/QfmLTWkUBn4OFGXv9/wBHAaOBgYCxwfQuOfzbwE6APsAi4CcDdvxS8fnDwV/EDJC5VlAP9SYzgXQeoLpmIAODuzwLjgeXB58bZ7t49ePlgd98r+IPuMeAtYDBwHHCFmX012O4HJD6XTgJ6Av8fUNHIZ1J99wEnmVlPSCR5wBnAX4PXpwPVwN7AIcAJwDeb+HEuAOa6+0PAB8C5DV7fA+hLYlRwkpkdCtwNXAr0A+4EHk3+sQwsBo4GepH4zP2zmQ1q4vwSAiVn0lr9gDXuXt3I6+cCP3X3Ve6+mkSnP78Fx/+7u78WHP8vJJK8xlQBg4Bhwejdi66isSKF6h/BCFHy8a009/s80N/df+rule6+BPgDcFbw+jeB6919gSe85e5rmzuou38MvAGcFjQdSyKpeyUY9RoPXBFcfVgF/LbeOVO5gJ2J3V/Z/dJmLXCDu+9w923At4A73f1Vd69x9+nADhJ/POPuf3P35e5eGySWC0n8MS0RUnImrbUWKG5iTsOewMf1nn8ctKXrs3rfVwDdG9sQ+BWJ0bV/mtkSM7umBecRkXg5zd1713v8Ic39hpG47FmX2JEYhU9eNiwhMcrUGn8lMeoGcA47k6thQAdgRb1z3gkMSHUQM/sCMBy4v95xDzSz+n+8rnb37Q1+risb/FwlBJ/HZnZBvUueG4ADgOJW/pySIZosKK31MrCdxF+DD6Z4fTm7TrYdGrQBbAW6Jjc0sz1oA3ffTOLS5pVmtj/wvJm97u7PteW4IlJQlgFL3X1kE6/vBbzbimP/Dfi1mQ0B/oPEtJDkMXcAxU1chajvQsCA+WZWv/0CYH7wfcOrBsuAm9z9poYHC+av/YHEJdyX3b3GzOYH55AIaeRMWsXdNwL/DfzezE4zs65m1sHMxpvZL0nMs7jezPqbWXGw7Z+D3d8C9jez0WbWGfifFp5+JTAi+SSY7Lq3JT6tNgE1wUNEJF2vAZuCyfRdghsADjCzzwev/xH4mZmNtISDzKxf8Noun0kNBVM7ZpGYoL/U3T8I2lcA/ySRuPU0s3ZmtpeZHdPwGMFn5RkkbgQYXe/xn8C5TVzF+APwbTM7PIi7m5lNMLMeQDcSydzq4BwXkxg5k4gpOZNWc/ffkJgkez2Jzr0M+C7wD+BGYC7wNvAOiTkXNwb7fUTiTs9nScxveKnhsZvxP8D0YBj+DGBkcKwtJEb0bnP3WW340UQkfz1mu65z1thNS7tw9xpgIomEZymwhkRC1ivY5DfADBLJ1CbgLqBL8Nr/sOtnUip/BY5n5yXNpAuAjsD7wHoSVyJSTcg/DdgG3OPunyUfQRxFwImN/FxzScw7+11w/EXARcFr7wO/JvG5uRI4EPh3I/FLiEzzpkVERERyh0bORERERHKIkjMRERGRHKLkTERERCSHKDkTERERySFKzkRERERySF4sQltcXOylpaVRhyEibTBv3rw17t4/6jhymT7rROKhrZ93eZGclZaWMnfu3KjDEJE2MLOPm9+qsOmzTiQe2vp5p8uaIiIiIjlEyZmIiIhIDlFyJiIiIpJD8mLOmUgcVFVVUV5ezvbt26MOJas6d+7MkCFD6NChQ9ShiIjkJSVnIiEpLy+nR48elJaWYmZRh5MV7s7atWspLy9n+PDhUYcjIpKXdFlTJCTbt2+nX79+sU3MAMyMfv36xX50UEQkm5SciYQozolZUiH8jCIi2aTLmiIRuCVLCcxV7lk5roiIhEcjZyIiIiI5JFYjZ+PGjQNg1qxZkcYhkq5MjXRlayQuXQsWLABg1KhRkcZRSPR5JxK+sPpdrJIzEWlaWVkZ48eP54tf/CJz5sxh8ODBPPLII3Tp0mW3befPn8+3v/1tKioq2Guvvbj77rvp06cP48aN4/DDD+f5559nw4YN3HXXXQwYMCCCn6awXXnllVGHIFJwwup3uqwpUmAWLlzI5ZdfznvvvUfv3r156KGHUm53wQUX8Itf/IK3336bAw88kJ/85Cd1r1VXV/Paa68xZcqUXdolPBMnTmTixIlRhyFSUMLqd0rORArM8OHDGT16NACHHXYYZWVlu22zceNGNmzYwDHHHAPAhRdeyAsvvFD3+te+9rUm95fsW7BgQd3lZBEJR1j9Tpc1RQpMp06d6r4vKipi27ZtrT5GUVER1dXVGYtN0nfppZcCmnMmEqaw+l2skrPvfOc7UYcg0iJRT+RvTK9evejTpw8vvvgiRx99NPfee2/dKFoq/fv3DzE6EZF4i1VyduaZZ0YdgkhsTJ8+ve6GgBEjRvCnP/2p0W379u0bYmQiIvEWq+Rs2bJlAJSUlEQciUjTolostrS0lHfffXdnHFdd1ei2o0eP5pVXXtmtvf5wfnFxMWVlZVRWVgLQsWPHzAUrIlKgYpWcnX/++YDmYIiEbenSpYDWORMRyYRYJWci0nKXX345//73v3dp+/73v8/FF18cUUSSjuuvvz7qEEQKTlj9TsmZSIjcPecKg//+97/P6PFc9T1Dcfzxx0cdgkjBCavfaZ0zkZB07tyZtWvXxjp5cXfWrl1L586dow4l9ubPn8/8+fOjDkOkoITV7zRyJhKSIUOGUF5ezurVq6MOJeM+++wzAGpra+ncuTNDhgyJOKL4u+KKKwDNsRUJU1j9LuvJmZkVAXOBT939ZDMbDtwP9AXeAM5398pMnEu15iSXdejQgeHDh0cdRlYsWrQIgH333TfiSERE8l8YI2ffBz4AegbPfwH81t3vN7M7gEuA2zNxItWZE4mG+p6ISOZkdc6ZmQ0BJgB/DJ4bcCzwYLDJdOC0TJ1PteZEopELfc/MyszsHTObb2Zzg7a+ZvaMmS0MvvYJ2s3MbjWzRWb2tpkdWu84FwbbLzSzC+u1HxYcf1Gwr7X2HCIiTcn2DQFTgB8BtcHzfsAGd08W4ysHBqfa0cwmmdlcM5ub7hydSy+9tK7ulYiEJ4f63pfdfbS7jwmeXwM85+4jgeeC5wDjgZHBYxLB6L2Z9QVuAA4HxgI3JJOtYJtJ9fY7sTXnEBFpTtYua5rZycAqd59nZuOSzSk2TXnrmrtPBaYCjBkzJr63t4lINp0KjAu+nw7MAq4O2u/xxK2zr5hZbzMbFGz7jLuvAzCzZ4ATzWwW0NPdXw7a7yEx6j+zpedw9xWZ+MF+/vOfZ+IwItICYfW7bM45+wJwipmdBHQmMedsCtDbzNoHo2dDgOVZjEFECocD/zQzB+4M/sAbmEyG3H2FmQ0Ith0MLKu3b3IUv6n28hTttOIcGUnOjjrqqEwcRkRaIKx+l7XkzN2vBa4FCEbOrnL3c83sb8DXSdyxeSHwSLZiEJGC8gV3Xx4kR8+Y2YdNbNvYKH5L25uS1j5mNonEZU+GDh3azCF3mjNnDqAkTaS1Sq95otHXyiZPSNkeVr+LYp2zq4H7zexG4E3grghiEJGYcfflwddVZvYwiTljK5OXEoPLlquCzcuBknq7J0fxy9l5iTLZPitoH5Jie1pxjoZxt2oKx3XXXQdonTORMIXV70JJztx9FokPONx9CYkPzYxrrubVLWmUzbkqxqu3i2RL1HUezawb0M7dNwffnwD8FHiUxAj9ZHYdqX8U+K6Z3U9i8v/GILl6Gvh5vZsATgCudfd1ZrbZzI4AXgUuAP6v3rHSPkeW3gIRaaGyyRO4Y/ZiJs9MDLL/z8T9uOgLubEWZawqBKjWnEg0cqDvDQQeDla3aA/81d2fMrPXgRlmdgnwCfCNYPsngZOARUAFcDFAkIT9DHg92O6nyZsDgO8A04AuJG4EmBm0T27JOUQk87ZX1bBo1RYOGNyrRfvc9dLSuudlayuyEVqrxCo5S9a7Gj16dJPbpRodS2dUTURSS7fvZUswIn9wiva1wHEp2h24vJFj3Q3cnaJ9LnBAJs4hIpnRmnljSX9/41NWb95BhyKjqsb5eO3WTIfXarFKzlRrTiQa6nsikk9qap07X1gMwKVf2ovfPb9II2ciItI2U6ZMiToEkcjU1Dq3/HNnVZIJBw2iYkc1zy9YzWmj92TKWYc0uf/Md1fw8doKhvbtynfGJZKzZesqqK6ppX1R4+vzh9XvlJyJiOShqC4hi0Rt3dZKvn//m7y4cA1F7Yxrx3+OS744nPL12zj+N7P5x/zlnDV2KEeM6Jdyf3fnjtmJUbNJXxpBt07tGdSrMys2bmf5hu0M7de10XOH1e+yXb5JRESy4Nlnn+XZZ5+NOgyRFilbs5V5H69v9f7vlG9k4v+9xIsL19CvW0f+fMnhfPPoEZgZJcEoGMANj7xHVU1tymO8tGgN7366ieLuHfn6YYkVcoYFCdnSZuadhdXvNHImIpKHbrzxRiAn7pQVaVJTk/ah+Yn7Sf9481N+9NDbVFbXcnBJb+4471AG9eqyyzbfPmYvHnqjnAUrN3PPyx9zyRd3Xxrj9lmJUbOLvzCczh2KABhe3I1XlqwLbgro32gMYfW7WCVnqjUnEg31PRHJhsYSu7eWbeDIm/+1W2LXuUMRN5y8P9+8Zy5TnvmIiQcPYkCPznWvz1+2gTmL19K9U3vOO2JYXfuwft0AKFuTGzcFxCo5UxkTkWio74lIY5IJ1H2vfcK1f38HyO6Cr8fvN5BjPzeAf324islPfshvztw5T+yOYNTs3COG0qtLh7r20uCyZlmOLKcRq+RMteZEoqG+JyLNeXnx2rrvl65JLwlKJnaPzP+U798/n6/uP5A7zx/T7H43TNyPlxat4e9vfspZY4cydnhfFq/ewtPvf0bHonZc0iAxLC0ORs5yJDmL1Q0B1113XV3dKxEJj/qeiDTF3Xl5yc7kbEmayVnd9qsT24/o3z2t7Yf168a3j0ncHPDfj7xLdU0tU2cvwR1OP2wIA3p23mX7oX0TI2fL1lVQUxt9GcdYjZyJiBSKO++8M+oQRNK2ePXWXVbjTyZb6UomcyOCEa50XDZuL/7+RjkffraZW/75EX9/s5x2Bpd+acRu23bt2J6BPTuxctMOlm/YRknf1MtphNXvYjVyJiJSKEaNGsWoUaOiDkMkLclRs2M/N4CidsbyjdvYXlWT9v5LVm8B0h85g8TNAf998n4A3DF7MVU1zvgDBtVdwmyotF/zlzbD6ndKzkRE8tBjjz3GY489FnUYIml5JUjOjh7Zn6F9u+Ke/vwud6+bo7ZX//RHzkqveYJJ987bpe2Jd1ZQes0TKe8CrUvOmrjkGla/02VNEZE89Otf/xqAiRMnRhyJSNPcnVeD5OzIvfrx/IerWLpmK0tXb+Vze/Rsdv/PNm2norKGvt060rtrx6zFOaw4ecdm48tphNXvYpWcqdacSDTU90SkMQtXbWHNlkoG9OjEiOJuDA8uK6Z7U8DiVS0fNYNdF7etqqmlfTvDzBrdfngwcvZxDtyxGavkTLXmRKKhvicijUkuoXHkXv0wM4YHSVa6NwUsWRPMNytOf75ZQx2aKGaelFyINt1lPrIpVnPOVGtOJBrqeyLSmGRylixEnkyylgZJV3N2LqPRspGzlkrW11y2blvky2nEauRMteZEoqG+JyKp1NY6rywNRs6SyVn/lo1QLW7FnZqt0a1Tewb06MSqzTtYsXEbQ/qkXk4jDLFKzkRECsW9994bdQgizVqwcjMbKqoY1Ktz3cjUgB6d6NaxiPUVVazfWkmfbk1P8g9r5AwSd2yu2ryDsjUVKZOzsPpdrC5riogUipKSEkpKSqIOQ6RJdfPNRvSrm4y/y7yzZkbPtlfVsHzjNtq3s7pV/LOptLjpGpth9TuNnGXALU3c/ZF0lUdfDkJE4uOBBx4A4Mwzz4w4EpHGJRefPWKvfru0Dy/uzrufbmLpmq0cNqxPo/svXbMVdxjar2tak/rbalgzd2yG1e+UnImI5KHbb78dUHImuaumtt76ZiMaJmfJOzabvikgzEuasHMh2qVrUq91Fla/i1VyFnWtuVSjY+mMqonku6j7nojkng9WbGLT9mqG9OmyW63KvdK8KaA1ZZvaInlZM+q1zmKVnKnOnEg01PdEpKGGS2jUlxw5azY5a0XZpraou6y5roLaWqddu2gGWGJ1Q4BqzYlEQ31PRBp6pZFLmrBrclbbxJpiYY+cde/UnuLunaisrmXFpu2hnDOVWI2cqdacSDTU90SkvuqaWl5bug5IVAZoqEfnDvTv0YnVm3ewvJE1xdx955yz4nBGzgCGF3dlzZYdfLxmK4N7dwntvPXFKjkTESkUDz74YNQhiDTqveWb2LyjmmH9urJnIwnO8OJurN68g6VrtqZMzlZv3sHmHdX06tKBvs2shZZJw/p14/Wy9ZStreCovXd9Lax+F6vLmiIihaK4uJji4uKowxBJ6eUmLmkmjShuusbm4np3ajZVsDzTSvs1vtZZWP1OyZmISB6aNm0a06ZNizoMkZTqFztvTHNlnDJR8Lw1SoOksSxFXGH1OyVnIiJ5SMmZ5KqqmlpeL0vMN0t1p2bS8CDpWtzIWmdhr3GWVFq3EO3ua52F1e9iNedMteZEoqG+JyJJb5dvpKKyhhHF3RjYs3Oj2zW3nEbyTs2wltFIGlbvsmZUy2nEKjnLZr0rLSYr0jjVeBSRpFcaKdnU0NC+XSlqZ3y6YRvbq2ro3KFol9eTa5yFtYxGUo/OHSju3pE1WypZuXk7g3qFf8dmrC5rPvDAA3V1r0QkPOp7IpLU1Ppm9XVs346SPl1w3/0S4o7qGpatq6Cd7RzJCtOwfuktkpstsUrObr/99rq6V9lylXujD5FCFUbfS4eZFZnZm2b2ePB8uJm9amYLzewBM+sYtHcKni8KXi+td4xrg/YFZvbVeu0nBm2LzOyaeu0tPodIXFVW1zK3bD3Q9HyzpJ2XNnedd/bJ2gpqHUr6dqVT+6JUu2ZVMiFMNe8sDLG6rCkiBe/7wAdAz+D5L4Dfuvv9ZnYHcAlwe/B1vbvvbWZnBdudaWb7AWcB+wN7As+a2T7BsX4PfAUoB143s0fd/f2WniNTP+iTTz6ZqUOJtFnpNU/s1vb5m54FoGzyhEb3G9G/O88vWF13CTNpcQSLz9Y3PBg5a7icRlj9LlYjZyJSuMxsCDAB+GPw3IBjgeSqkdOB04LvTw2eE7x+XLD9qcD97r7D3ZcCi4CxwWORuy9x90rgfuDUVp4jI7p27UrXruFf7hHJpOGNrHVWt4xGyPPNkoY1spxGWP1OI2ciEhdTgB8BPYLn/YAN7l4dPC8HBgffDwaWAbh7tZltDLYfDLxS75j191nWoP3wVp5jTdt+zITbbrsNgMsuuywThxNpk7LJE1i0agtf+e1sOhS146Wrv8yAHo3fqZnU2FpnUS2jkVTayGXNsPqdRs5EJO+Z2cnAKnefV785xabezGuZam/u/HXMbJKZzTWzuatXr06xS2ozZsxgxowZaW8vkm1/eGEJ7nD6oUPSSsxg5wKzuydn0SxAmzSs3mVNrzenPKx+F6uRM9WaE4lGDvS9LwCnmNlJQGcSc86mAL3NrH0wsjUEWB5sXw6UAOVm1h7oBayr155Uf59U7WtacY5duPtUYCrAmDFjdGeR5KWVm7bz8JufYgbfOnp42vsN7NmJrh2LWLe1kg0VlfTu2hF3r5tzFvYaZ0nJep7rtlayctMO9uiVXrKZKbEaOVOtOZFoRN333P1adx/i7qUkJvT/y93PBZ4Hvh5sdiHwSPD9o8Fzgtf/5Yk/jx8FzgrutBwOjAReA14HRgZ3ZnYMzvFosE9LzyESO3/6dxmVNbWcuP8eLZonZmY7550Fo2frtlaycVsV3Tu1p3+PTlmJNx3DmqixmW2xSs5UzkQkGjnc964GfmBmi0jM97oraL8L6Be0/wC4BsDd3wNmAO8DTwGXu3tNMCr2XeBpEneDzgi2bfE5ROJm8/Yq/vLKxwBM+tKIFu9ft5xGMFq2c/HZcAue7xZXXRmn8JOzWF3WTP7ncNFFF0Uah0ihyaW+5+6zgFnB90tI3GnZcJvtwDca2f8m4KYU7U8Cu91H35pziMTJfa99wuYd1Rw+vC+HDO3T4v1H1I2cJeaZ7ZxvFs0lzaSdC9GGv9ZZrJIzEclfZvaDNDbb6u53Zj2YPDBr1qyoQ5CYSbVWWX2p1iurrK7lrpeWAvDtY/Zq1XmTl0GTNwUsqZtvFs3NAEmlxck7NneOnIXV72J1WVNE8toPge4klsJo7HFlZNGJyG7+Mf9TVm7awaiBPRg3qn+rjtFwrbO6BWijTs7q7tjUyJmIFK573f2nTW1gZtFe58ght9xyCwBXXXVVxJFIXJRNnsBdLy3lZ4+/X9fWoci447zDOG7fgbttX1vrTH1hCQCXHjOi1fPDhvffuWxFba3XW4A22u5eWm/OmbtjZqH1u6yNnJlZZzN7zczeMrP3zOwnQXvKOnQiUtjc/UeZ2KZQPP744zz++ONRhyEx8tC88rrE7JdfP4iLv1BKVY3znT+/wfMfrtpt+399uIpFq7awZ6/OTDx4z1aft2fnDhR378T2qlqWra/gk7UVmO0cUYtKr64d6N21AxWVNazevAMIr99lc+RsB3Csu28xsw7AS2Y2k8RdS6nq0LWZas2JRCNTfc/MPkdiZf1X3X1LvfYT3f2pjJxERHbz7Psr+dFDbwNw/YR9OWNMCe6OO0ybU8al987jzgsO48ujBtTtc+cLiwH4/744nA5FbRvrGVHcjTVbdvDCR6uprnUG9+5C5w7hFzxvqLRfN+ZXbKBsbQUDeoa31lnWRs48Ifnh2iF4OI3XoWsz1ZoTiUYm+p6ZfY/EGmH/CbxrZqfWe/nnbTq4iDTq1SVrufyvb1BT61z+5b345tGJ5TDMjBsm7sdFR5VSWVPLpffOY9aCxAjavI/X8XrZenp2bs9ZY4e2OYbkKNmzHySOH/UlzaRkGaeGNTazLas3BJhZkZnNB1YBzwCLabwOXcN9W1zS5LbbbqureyUi4clQ3/sWcJi7nwaMA35sZt8PXotusSORGHv30418c/pcdlTXcs7hQ7nqhFG7vJ5M0C44chiV1bVMuncesz9azZ2zE3PNzj9yGN07tf0iXDIZe3nxWiD6OzWTSot3zocLU1ZvCHD3GmC0mfUGHgb2TbVZI/u2uKRJst6VCgGLhCtDfa8oOdru7mVmNg540MyGoeRsN126dIk6BMlzS1Zv4cK7X2PzjmomHDSIn516QMpJ/WbGT07ZH3e495WP+dY9c6mqqaVj+3ZcdFT6pZqakhw5q6ypBaIfOWu4rMhtsxZz26zFrPxkS6vvSm2JUO7WdPcNZjYLOILG69CJSGH7zMxGu/t8gGC+6snA3cCB0YaWe2bOnBl1CJLHPtu4nfPveo21Wys5emQxvz1jNEXtGv8byMy4N6gCUFldW/f18zc9C6ReA60lGiZjURU8b87AM37CzDb+rOnIWnJmZv2BqiAx6wIcD/yCnXXo7mfXOnQiUtguAKrrNwR/xF1gZlp4ViQDUi00++LCNexz/cw2J1htMbRvN9oZ1AbXyaIeOYvyvYDsjpwNAqabWRGJuW0z3P1xM3sfuN/MbgTeZGcdOhEpYO5envzezPoAJez8jNoWSVA57Gc/+xkAP/7xjyOORApFNhOWju3bUdK3Kx+vraBrxyL2CPHOyJYIq99lLTlz97eBQ1K0p6xDJyICYGY/Ay4icQNRcr5p8k5vCTz33HOAkjNpmR+dOIpfPrWAbh2LePjyL7DPwB5Rh1RneHE3Pl5bwfDibrRr4hJrlMLqd7GqEKBacyLRyHDfOwPYy90rM3lQkUL37Psr+dXTCwD47ZmjcyYxa3ip9b3lm+raor68GBXV1hSRXPMu0DvqIETiZOHKzVzxwHzc4cqv7MMJ++8RdUjShFiNnKnWnEg0Mtz3bgbeNLN3SVQaAcDdT8nEwUUKzcaKKr51z1y27KhmwoGD+O6xe0cd0i7qj47V1DrtjFbX6YyLWCVnyXpXSs5EwpXhvjedxJ3d7wC1mThgHPXr1y/qECQPVNfU8t373qBsbQX7DerJr75xUE4nPk0t55ELwup3sUrORCQW1rj7rVEHkeseeuihqEOQPHDzzA95ceEa+nXryNQLDqNrR/233xZh9Tv9K4lIrplnZjcDj7LrZc03ogtJJL+4O/e9toy7XlpK+3bG7ecdxpA+qj2dL5SciUiuSS7Bc0S9Ni2l0cC1114LwM2xXg5hAAAgAElEQVQ33xxxJJJLUi0yW13rnHHny0Dh3v2YKWH1u1glZ6o1JxKNTPY9d/9yxg4WYy+//HLUIUgO2VFdwyNvqhpitoXV72KVnKnWnEg0MtH3zOxkd3+8rduIFJItO6q579VP+ONLS1i5KTELYHDvLnzr6OGc8fkSzTHLU/pXE5Fc8Ssz+xRo6natnwNKzqSgpbp0mfTbMw/m5IP2pEORljHNZ7FKzlRrTiQaGep7K4HfNLPNwracQCTu/uOQIVGHIBkQq+RMteZEopGJvufu4zIUTkEYMkT/CReqsskTeOb9lXzrnrn07daR2T8cR4/OHaIOqyCE1e9ilZyJiBSKP//5z1GHIBGpqXV++dSHAPznsXsrMQtRWP1OF6VFRETyyENvlLNw1RZK+nbhnMOHRh2OZIGSMxGRPHTFFVdwxRVXRB2GhGx7VQ2/feYjAK78yig6tS+KOKLCEla/i9VlTdWaE4lGpvuemR0A7Ad0Tra5+z0ZPUmemz9/ftQhSASmzyljxcbt7DuoJ6ccvGfU4RScsPpdrJIz1ZoTiUYm+56Z3QCMI5GcPQmMB14ClJxJQdtYUcVtsxYDcM34z9Eux4uES+vpsqaI5JqvA8cBn7n7xcDBQKdoQxKJ3u2zF7NxWxVHjujHl0YWRx2OZFGskrNrr722ru6ViIQnw31vm7vXAtVm1hNYBYzI1MFF8tGKjdv407+XAolRMzONmsVZrC5rqtacSDQy3Pfmmllv4A/APGAL8FomTxAH++yzT9QhSIj+99mF7KiuZcKBgzi4pHfU4RSssPpdi5MzM+sDlLj721mIR0QKnLtfFnx7h5k9BfTU583upk6dGnUIEpJFqzYzY+4yitoZV56gpDxKYfW7tJIzM5sFnBJsPx9YbWaz3f0HWYytxZbNng3ALRruFclblrhecy4wwt1/amZDzWysuzc6emZmnYEXSMxNaw886O43mNlw4H6gL/AGcL67V5pZJxI3GBwGrAXOdPey4FjXApcANcD33P3poP1E4H+BIuCP7j45aG/xOURa4pdPLaDW4dyxJYzo3z3qcCQE6c456+Xum4CvAX9y98OA47MXlogUsNuAI4Gzg+ebgd83s88O4Fh3PxgYDZxoZkcAvwB+6+4jgfUkki6Cr+vdfW/gt8F2mNl+wFnA/sCJwG1mVmRmRUEM40ncRXp2sC0tPUemTJo0iUmTJmXykJKD5n28nn++v5IuHYr4/nEjow6n4IXV79K9rNnezAYBZwD/L4vxtEmv4OtV7pHGIVJoMlxv7nB3P9TM3gRw9/Vm1rGpHdzdScxNA+gQPBw4FjgnaJ8O/A9wO3Bq8D3Ag8DvghG7U4H73X0HsNTMFgFjg+0WufsSADO7HzjVzD5o6TmCWNvso48+ysRhJAeVXvPEbm3bqmoY+/NEDduyyRPCDkkCYfW7dJOznwBPAy+5++tmNgJYmL2wWuec5jcRkSzIcL25qmCkygHMrD9Q29xOwT7zgL1JjHItBja4e3WwSTkwOPh+MLAMwN2rzWwj0C9of6XeYevvs6xB++HBPi09x5oGcU8CJgEMHapSPCKSfnK2wt0PSj5x9yVm9pssxSQihe1W4GFggJndRGLds+ub28nda4DRwZ2eDwP7ptos+JpqYqo30Z5qCkhT2zd1jl0b3KcCUwHGjBmjYX+hbPIEfvPPBdz6r0V079Seh75zFKP26BF1WBKidJOz/wMOTaMtUo8EX69qwzF0M4FIyyVrzU2ZMqXNx3L3v5jZPBIL0Rpwmrt/0IL9NwQ3MR0B9Daz9sHI1hBgebBZOVAClJtZexKzItbVa0+qv0+q9jWtOIdIkx5+s5xb/7WIdga/O+cQJWYFqMnkzMyOBI4C+ptZ/Tsze5K4YymnLG9+ExHJgkzVmzOzdsDb7n4A8GEL9usPVAWJWRcSNyz9AniexMjb/cCF7Pwb7tHg+cvB6/9ydzezR4G/BlcG9gRGklhjzYCRwZ2Zn5K4aeCcYJ8WnaMVb0tKo0ePztShJIe8XraOqx98B4AbJu7PuFEDIo5I6gur3zU3ctYR6B5sVz9130TiwyY2dBOBSPTcvdbM3jKzoe7+SQt2HQRMD+adtQNmuPvjZvY+cL+Z3Qi8CdwVbH8XcG8w4X8diWQLd3/PzGYA7wPVwOXB5VLM7Lsk5t4WAXe7+3vBsa5uyTkyJROjlJJbPllbwaX3zqOyppYLjxzGhUeVRh2SNBBWv2syOXP32cBsM5vm7h+HEpGIFLpBwHtm9hqwNdno7qc0tkOwSO0hKdqXsPNuy/rt24FvNHKsm4CbUrQ/SaIQe5vPIdLQxm1VXDztNdZtrWTcqP78+OT9mt9JYivdOWedzGwqUFp/H3c/NhtBiUhB+0nUAeSD8847D8j4nbISgaqaWi7/yxssXr2VUQN78H9nH0L7oliVvo6NsPpdusnZ34A7gD+SWDU7J/WPOgCRApXJenPBiH0dM/sCiZVyZqfeozCVl5dHHYJkgLtzw6Pv8dKiNRR378hdF42hR+cOUYcljQir36WbnFW7++1ZjSQDYjUJjvTvHNV8OYlapuvNmdloEgnZGcBS4KGMnkAkYqkWml2zpZIv/uJ5LTIraSdnj5nZZSTWDtqRbHR33RYuIhlhZvuQmDR/NolalA8A5u5fjjQwEZGQpZucXRh8/WG9NgdGZDactnkw+NqWdc5yUWMjY1qTTXJFstZcG0fQPgReBCa6+yIAM/uvNgcnkoPKJk/g0beW87373mRQr848f9U4OnfIuRWqJCJpJWfuPjzbgWTC6qgDEClQGao3dzqJkbPnzewpEuuG6S+QRhx55JFRhyBtUF1Ty5RnEv3me8eNVGKWJ8Lqd2klZ2Z2Qap2d78ns+GISKFy94eBh82sG3Aa8F/AQDO7HXjY3f8ZaYA55uabb446BGmDv7/xKUvWbGVYv658/bAhUYcjaQqr36V7WfPz9b7vTKKsyhuAkrM0ZfMSZFPH1s0Ckm/cfSvwF+AvZtaXxFph1wBKziQWdlTX8L/PLQTgB1/Zhw5aNkMaSPey5n/Wf25mvYB7sxKRiEgguOnozuAh9Zx++ukAPPSQbmTNN/e9+gmfbtjGqIE9mHjQnlGHIy0QVr9Ld+SsoQoSNedySi7+imdz5KqpY+tmAQmT6jyGb+3atVGHIK1QUVnN755fDMAPTtiHdu30WZ1Pwup36c45e4zE3ZmQqCu3LzAjW0G11qlRByBSoFTnUSQ90+d8zJotOzh4SC9O2G9g1OFIjkp35OyWet9XAx+7u5anFpGsMLNhwEh3f9bMugDt3X1z1HGJtMWm7VXcMTsxanblCaMwXeGQRqQ1CzEop/Ih0APoA1RmM6jW+mvwEJFwnXfeeXU159rKzL5FYtnC5DyzIcA/MnJwkQj98cWlbNxWxdjhfTl6ZHHU4UgOS/ey5hnAr4BZJNYd+j8z+6G7P9jkjiHbGHUAIgUqw/XmLgfGAq8CuPtCMxuQyRPEwXHHHRd1CNIC67ZWcteLSwD44Vc1apavwup36V7W/H/A5919FYCZ9QeeZeei/CIimbLD3SuT/3mZWXt2znmVwI9//OOoQ5AWuGP2YrZW1nDMPv35fGnfqMORVgqr36W7uEq7ZGIWWNuCfUVEWmK2mV0HdDGzrwB/Ax6LOCaRVlu5aTvT55QBcNUJo6INRvJCuiNnT5nZ08B9wfMzgSeb2sHMSkgsUrsHUAtMdff/DRaVfAAoBcqAM9x9fctDF5GYuga4BHgHuJTEZ80fI40oB40fPx6AmTNnRhyJNKb0mid2a5v4u5eARG1NyT9h9bsmkzMz2xsY6O4/NLOvAV8kMefsZRIreDelGrjS3d8wsx7APDN7BrgIeM7dJ5vZNSQ+iK9u488BwLBMHEREWizD9eZOBe5x9z9k8qBxs23btqhDECk4YfW75kbOpgDXAbj734G/A5jZmOC1iY3t6O4rgBXB95vN7ANgMIkP3nHBZtNJ3GSQkeTspEwcRERaLMP15k4BppjZCySKnz/t7tWZPIFIGMomT2D6nDJuePQ9Dh/elwcuVbF6SU9z88ZK3f3tho3uPpfEZcm0mFkpcAiJu68GBolbMoFLeReWmU0ys7lmNnf16tXpnkpE8py7XwzsTWKu2TnAYjPTZU3JO+7OPS+XAXDhUaVRhiJ5prnkrHMTr3VJ5wRm1h14CLjC3TelG5i7T3X3Me4+pn///mntMz14iEi4Tj/99Lqac5ng7lXATBIjZ/NQARDJQ3MWr2Xx6q0M7NmJr6gagLRAc5c1XzezbzWc+2Fml5D4wGySmXUgkZj9JbgsCrDSzAa5+wozGwSsavwILVORqQOJSItkst6cmZ0InAV8mcS0hz8CZ2TsBDFx8sknRx2CNCM5anbO2GF0KNICB3EQVr9rLjm7AnjYzM5lZzI2BugI/EdTO1pikaK7gA/c/Tf1XnoUuBCYHHx9pBVxi0h8XURixOxSd98RcSw566qrroo6BGnCpxu28cz7K+lQZJx9eEnU4UiGhNXvmkzO3H0lcJSZfRk4IGh+wt3/lcaxvwCcD7xjZvODtutIJGUzgtG3T4BvtCpyEYkldz8r6hhE2uqvr35MrcOEAwYxoEdTM4REdpfWOmfu/jzwfEsO7O4vkVh2IxXVHRGRXZjZS+7+RTPbzK4VAQxwd+8ZUWg5ady4cQDMmjUr0jhkdzuqa7j/tWUAXHCkFnmKk7D6XbqL0OaFvaMOQKRAZaLenLt/Mfjao80HE4nQzHc+Y+3WSvYd1JMxw/pEHY7koVglZ1+JOgCRApXJenNmdq+7n99cm0iumv5yGZAYNVOBc2kN3T4iIrlm//pPgsLnh0UUi0iLvFO+kTc/2UCPzu05dfSeUYcjeSpWydkfUQE+kSiMHz++ruZca5nZtcF8s4PMbFPw2AysRHd1S55ILp9xxpgSunaM1cUpCVGsfnOqog6gANySgSH6q9yb30jySibqzbn7zcDNZnazu1/b9qji7YwztPRbrlm/tZJH31oOwHlH6EaAOAqr38UqORORWHjNzHq5+0YAM+sNjHP3f0QcV0657LLLog5BGvjbvGXsqK7lmH36M7y4W9ThSBaE1e+UnEmrtGb0KxOjblIQbnD3h5NP3H2Dmd0AKDmrp6IiUROla9euEUciADW1zr2vfAxo+Yw4C6vfxWrOmYjEQqrPpSb/kDSzEjN73sw+MLP3zOz7QXtfM3vGzBYGX/sE7WZmt5rZIjN728wOrXesC4PtF5rZhfXaDzOzd4J9bg2qoLTqHJlw0kkncdJJJ2XykNIGsz9axbJ12xjSpwvjRg2IOhzJkrD6XaySs32Dh4iE6+STT85kzbm5ZvYbM9vLzEaY2W9pvpZvNXClu+8LHAFcbmb7AdcAz7n7SOC54DnAeGBk8JgE3A6JRAu4ATgcGAvckEy2gm0m1dvvxKC9ReeQeLrn5cSo2flHDKOona4SSNvE6rLmuKgDyCPpXGLUxH1JV4brzf0n8GPggeD5P4Hrm9rB3VcAK4LvN5vZB8Bg4FR2fjRMJ1FI/eqg/R53d+AVM+ttZoOCbZ9x93UAZvYMcKKZzQJ6uvvLQfs9wGnAzJaeI4hVYqD0mid2a7t55ofcPPNDAMomTwg7JImJWCVnIpL/3H0rcI2ZdXf3LS3d38xKgUOAV4GByWTI3VeYWfJ602BgWb3dyoO2ptrLU7TTinPskpyZ2SQSI2sMHTq0ZT+siMRSrJKz5DWDcGrGx0Oq0TFN3JeWymS9OTM7isSShd2BoWZ2MHCpuzd7m5SZdQceAq5w901NrM6e6gVvRXuT4aSzj7tPBaYCjBkzRsPVeaRs8gS2V9Vw5M3Psb6iiocvO4pDhqpck7RdrJIzEYmF3wJfBR4FcPe3zOxLze1kZh1IJGZ/cfe/B80rk5cSg8uWq4L2cqCk3u5DgOVB+7gG7bOC9iEptm/NOTLioosuytShpA1mvruC9RVV7L9nT0aX9I46HMmysPpdrG4IEJF4cPdlDZpqmto+uHPyLuADd/9NvZceBZJ3XF7IzkoDjwIXBHdUHgFsDC5NPg2cYGZ9ghsBTgCeDl7bbGZHBOe6oMGxWnKOjLjooouUoOWAP7/yCZBYdFZ1NOMvrH6nkTMRyTXLgkubbmYdge8BHzSzzxeA84F3zGx+0HYdMBmYYWaXAJ8A3wheexI4CVgEVAAXA7j7OjP7GfB6sN1PkzcHAN8BpgFdSNwIMDNob9E5MmXNmjUAFBcXZ/Kw0gIfrNjEvI/X071Te045WHU0C0FY/U7JWczl0/wx3UEqgW8D/8vOSfj/BC5vagd3f4nUc7wAjkuxvTd2THe/G7g7Rftc4IAU7Wtbeo5M+PrXvw5kZp6ftM5fXk0sn/G1QwfTrZP+Oy0EYfW7WP02HRx1ACIFKhP15szsF+5+NfBldz+37VGJZM+WHdU8/MangOpoSubFKjk7KuoAckg+jzDpDtL8k6F6cyeZ2fXAtcDfMnFAkWz5x5ufsrWyhrGlfdlnYI+ow5GYiVVyVhl1ACIFKkP15p4C1gDdzGwTicuUyWUs3N17tjFMkYxwd/4c1NE89witTSeZF6u7Ne8KHiISrgzVm7ve3XsBT7h7T3fvUf9rBsIUyYg3PlnPh59tpl+3jpx4wB5RhyMxFKuRMxHJay8DhwKbog4kH3znO9+JOoSClVw+44zPl9CpfVHE0UiYwup3Ss5EJFd0NLMLgaPM7GsNX6y3sKwAZ555ZtQhFKR1Wyt54u0VmME5Y3VJs9CE1e+UnIlIrvg2cC7QG5jY4DUHlJzVs2xZYp3ekpKSZraUTHpw3jIqa2oZN6o/JX3bNMdS8lBY/U7JmYjkhGCtspfMbK67a/poM84//3xA65yFqbbW+curQUWAw7V8RiEKq9/FKjkbE3UAIgUqE+VMzOxH7v5Ld7/LzL7h7n+r99rP3f26Np9EpA1eWrSGj9dWMLh3F778uQFRhyMxFqvk7PNRByBSoDJUa+4s4JfB9w3XOjuRRDkmkcgkl884e2wJRe209qJkT6ySs61RByBp0YKy8ZOhenPWyPepnouEasXGbTz7wUratzPO+Lzm+Ul2xSo5uyf4ekOkUYgUngzVm/NGvk/1XCQUpdc8scvzWnfG3vQcAGWTJ0QRkhSAWCVnktvyuaSUhOLgepUBugTfEzzvHF1YuenKK6+MOgSRghNWv1NyJiI5wd21mmcLTJzYcLURyYayyROYPPND7pi9mHGj+jPt4rFRhyQRCqvfxap8k4hIoViwYAELFiyIOozYW7NlB9PnlAHwX8fvE20wErmw+p1GzqRRmrgvkrsuvfRSQOucZdudsxezraqG4/cdwMElvaMORyIWVr+LVXJ2ZNQBiBQo1XmUOFq1eTv3BstnXKFRMwlRrJKz0VEHEBOauC8tpTqPEkd3zFrC9qpaTthvIAcM7hV1OFJAYpWcbYg6AJECpTqPEjcrN23nz69q1EyiEavk7L7g642RRiFSeFTnUeLmtucXUVldy/gD9mC/PXtGHY4UmFglZyIiheL666+POoTYWr5hG/e9tgwzjZrJrsLqd0rORETy0PHHHx91CLF126xFVNbUcvJBgxi1R4+ow5EcEla/0zpnIiJ5aP78+cyfPz/qMGKnfH0FD7yeGDX7/nEjow5HckxY/U4jZyIieeiKK64ANM8v037//CKqapxTR+/JyIEaNZNdhdXvYpWcHRN1ACIFSnUeJQ6Wravgb3PLaWfwPY2aSYRilZztF3UAIgVKdR4lX5Ve80TK9uN+PRtI1NYUCVuskrNVUQcgOSmdMlRaeLdtkrXmRo0aFXEkIiL5L1bJ2UPB119GGoVI4VGdR8lXZZMnsHFbFeOnvMDyjdv5/nEj+a+vaPkMiVaskjORpqQaHVNxd8lXP//5z6MOITb++5F3Wb5xOweX9Oa7x+4ddTiSw8Lqd1lLzszsbuBkYJW7HxC09QUeAEqBMuAMd1+frRhEROLqqKOOijqEWHhk/qc8Mn85XToUMeXM0XQo0gpT0riw+l02fwunASc2aLsGeM7dRwLPBc9FRKSF5syZw5w5c6IOI68t37CN6//xLgA/Pnk/hhd3izgiyXVh9busjZy5+wtmVtqg+VRgXPD9dGAWcHW2YhARiavrrrsO0Dy/1qqtda6c8Rabt1dz3OcGcPbYkqhDkjwQVr8Le/x2oLuvAAi+DsjkwY8LHiISruuvvz7SWo9mdreZrTKzd+u19TWzZ8xsYfC1T9BuZnarmS0ys7fN7NB6+1wYbL/QzC6s136Ymb0T7HOrWWKyYmvOIbnh7n8v5eUla+nXrSOTTz8I0/xTySE5e3HdzCaZ2Vwzm7t69eq09tkneIhIuI4//vioaz1OI/1pFOOBkcFjEnA71M2JvQE4HBgL3JBMtoJtJtXb78TWnENyw4efbeKXTyWWf/nF6QfRv0eniCMS2VXYydlKMxsEEHxtdGkyd5/q7mPcfUz//v3TOvinwUNEwhV1nUd3fwFY16D5VBLTJwi+nlav/R5PeAXoHXwefRV4xt3XBTcqPQOcGLzW091fdncH7mlwrJacQyK2o7qGK+6fT2VNLWePHcrx+w2MOiSR3YSdnD0KJC8VXAg8kumDP5rJA4pIWq644oq6mnM5pLFpFIOBZfW2Kw/ammovT9HemnPspjVXCaT1fv3Pj/jws82U9uvK9RP2jTockZSyuZTGfSQm/xebWTmJywWTgRlmdgnwCfCNbJ1fRKQRqSYXeSvaW3OO3RvdpwJTAcaMGZN2qYopU6aku2nBS1WiqWxtBfvf8LTKM0mLhNXvsnm35tmNvKQ5+yIShpVmNsjdVzSYRlEO1L81bwiwPGgf16B9VtA+JMX2rTlHxowePTqThxORNITV73L2hgARkTZqbBrFo8AFwR2VRwAbg0uSTwMnmFmf4EaAE4Cng9c2m9kRwV2aFzQ4VkvOkTHPPvsszz77bCYPGUsbK6o4emQxAB2KjJu/diBlkyfUPURaIqx+p/JNIpL3WjiN4kngJGARUAFcDODu68zsZ8DrwXY/dffkTQbfIXFHaBdgZvCgpefIpBtvvBEg6rtkc9qiVZv55vS5lK2toF+3jtx+3mGMHd436rAkj4XV72KVnI2POgCRAhV1nceWTKMI7ri8vJHj3A3cnaJ9LnBAiva1LT2HhONfH67ke/fNZ8uOavYb1JOpFxzGkD5dow5LJC2xSs5Kow5ApECpzqPkCnfnjtlL+OXTH+IOEw4cxK++cRBdO8bqvzuJuVj9tpZFHYBIgUrWmlOSJlFJdUcmwBPvrOD356pAg+SXWCVnM5vfRGLoFpVdiZzqPIqIZE6skjMRkUJx5513Rh1CTnB37nn5YzoWtaOyppZ9B/Xk9+ccwoj+3aMOTWIorH6n5Exi4ypPe/1Okbw3atSoqEOI3MZtVVz94Ns89d5nAJx3xFCun7AfnTsURRyZxFVY/U7JmYhIHnrssccAmDhxYsSRRGP+sg18969vUL5+Gz06tefm0w/k5IP2jDosibmw+p2SMxGRPPTrX/8aKLzkbHtVDdPnlPGrpxdQXescOLgXvzvnEIb16xZ1aFIAwup3sUrOTok6AMk6Tf7PTarzKNnW2N2Y73y6kWN+NUur/UusxCo5Gxx1ACIFSnUeJRvcndfL1jNtztKoQxEJVaySs4+iDkCyRpP9c1uy1pxKCUkm7Kiu4ZE3l/OnOWV8sGITAO3bGScesAcXf6GUQ4f2wTSKLjEWq+TsuagDEClQqvMomdDYpUuAl64+lj16dQ4xGpHoxCo5ExEpFPfee2/UIWTUq0vWNvm6EjPJBWH1OyVnIiJ5qKSkJOoQMmLJ6i1Mnvkh/3x/JQD9e3TiqhP24euHlVDUTpcuJbeE1e+UnImI5KEHHngAgDPPPDPiSFpn7ZYd3PrcQv7y6idU1zpdOhRx6TEj+NbRI+jWSf81SW4Kq9+pB4iI5KHbb78dyK/krKk5Zduqarji+H1CjEak5cLqd7FKzk6POgCRAqU6j9KcRas2Rx2CSN6IVXI2IOoARAqU6jxKY974ZD13zFpcN6esncFJBw7i28fsxQGDe0UcnUhuilVy9n7UAYgUqEKv8yg7NXXp8tzDhzLpSyNUakmkGbFKzmZHHYDkrebKQmkR3KYVap1HSSwY+075Rl5duo7Xlq5rctub/uPAkKISyW+xSs5ERArFgw8+GMl5mxoZSxrWrytjS/sydnhfTjxgD3p07hBCZCLZF1a/U3ImBa25ETEVWpdcVVxcHPo5a2ubH0F+5drjtGCsxFZY/U7JmYhIHpo2bRoAF110UdbPtb2qhn+8+SlTX1xS19azc3suOLKUC44cxoCeSsakMITV75SciYjkoTD+k1i2roJH31rOtDllrN68A4DBvbtwyReHc+bnS7RYrBQcJWetcHbUAYgUqLjVeSxkS1Zv4dhfN3571awfjqNDUbsQIxIpPLFKznpHHYDEVlNzzxqbt5bOfLVM3wWa7hy5TJ83LnUeC5G789HKLTz5zgqeevczFqxserFYJWYi2Rer5Gx+1AGIFKh8r/NYSKpqavlgxSZO+d2/G93ma4cOZvwBgzh6ZDGdOxSFGJ2IQMySs5ejDkBip6kRpraMUmX7LtC2jOa1Rj7WeSwUFZXVzC1bz+tl65hbtp75yzawraqmyX1+c8bokKITkVRilZyJiBSKJ598MmV7ZXUt85dt4Iw7m/5z9RuHDWFMaR8OG9aXvfp3w7RsjEizGut3mabkTEQkD3Xt2rXu+w9WbGL2R6v596I1zC1b3+zIGMCvvnFwNsMTiaX6/S6blJyJtFGmLxVGNalf8sttt93Gj//xLj0OnZDy9YuOKuWovfpx+Ih+9OqiFfpFMuG2224D4LLLLsvqeZSciYjkoRkzZrB1ydpGk7P/OWX/kCMSib8ZM2YASs5a5IKoA5CCku2Rq7An9bdFVKZlK2YAACAASURBVHUeC90RI/oxa3Lq5ExE8lesFqzpFjxEJFzFxcWR1HrMB2Z2opktMLNFZnZN1PGISO6L1cjZ61EHIFKgwqzzmE/MrAj4PfAVoBx43cwedff3o41MRHJZrJKzuVEHIFKglJw1aiywyN2XAJjZ/cCpgJIzEWlUrJIzEUlPa8pRSasMBpbVe14OHF5/AzObBEwCGDp0aNoHnjVrVtujE5EWCavfxWrOmYhIjkmVBe+S/br7VHcf4+5j+vfvH1JYIpLLNHImUkAyUY5KWqQcqF8VfgiwPKJYRCRPaORMRCR7XgdGmtlwM+sInAU8GnFMIpLjYjVydknUAYgUqLDqzeUbd682s+8CTwNFwN3u/l7EYYlIjotVctYx6gBEWqgtlxKzdRmyueOmujQaVr25fOTuTwLKXkUkbbFKzuZEHYBIgQqr3pyISCGIVXL2VtQBiKSpLctVZGupi+aO29SIWlj15kRECkEkNwSonImIiIhIaqEnZ/XKmYwH9gPONrP9wo5DREREJBdFMXJWV87E3SuBZDkTERERkYIXxZyzZsuZiEjuSjX3bBlQcswx4QcjIhJDUSRnzZYzgV3rzQHbzaz+2kC9gI2NPC82szWZCDSFhufN1D5NbdPYa6nam3pfUj0vBrLxXuXb+9SwLc7vU1OvN/f70uQ2S2bPbq7vDWviNQHmzZu3ycwWNmgu9N/N1rbt8v8C+fU+Nbddtv5fyKX3Kd39cvV9atvnnbuH+gCOBJ6u9/xa4Npm9pma7nNgbhZjn5qNfZraprHXUrW35H3K5nuVb+9Tit+h2L5PLXmvWvo7ls2+VygP/W62vL828d5k/XczW+9Tc9tl6/+FXHqf0t0vTu9T/UcUc85aU87ksRY+z5bWnCedfZraprHXUrXrfUqvvbm2OL9PTb2ezu9LVL9ThUK/m5lrC+O9ytb71Nx2+fb/QmvPUYj/fwJgQRYYKjM7CZjCznImN2Xw2HPdfUymjhdneq/So/cpPXqfwqf3PD16n9Kj9yk9YbxPkSxC69ktZzI1S8eNI71X6dH7lB69T+HTe54evU/p0fuUnqy/T5GMnImIiIhIapFUCBARERGR1JSciYiIiOQQJWciIiIiOaSgkjMzG2Fmd5nZg1HHkmvMrJuZTTezP5jZuVHHk6v0O5Q+Mzst+H16xMxOiDqeQqLf06bp8y49+j1KTzY+6/ImOTOzu81slZm926D9RDNbYGaLzOyapo7hiXqel2Q30tzRwvfsa8CD7v4t4JTQg41QS96nQvsdaqiF79U/gt+ni4AzIwg3L+mzrnX0eZcefd6lJ+rPurxJzoBpwIn1G8ysCPg9MB7YDzjbzPYzswPN7PEGjwHhhxy5aaT5ngFD2FnztCbEGHPBNNJ/nwrdNFr+Xl0fvC7pmYY+61pjGvq8S8c09HmXjmlE+FkXyTpnreHuL5hZaYPmscAid18CYGb3A6e6+83AyeFGmHta8p6RKEA/BJhPfiXtbdbC9+n9cKPLLS15r8zsA2AyMNPd3wg10Dymz7rW0eddevR5l56oP+vy/ZdyMDv/+oFEhxvc2MZm1s/M7gAOMbNrsx1cjmrsPfs7cLqZ3Y5K8kAj75N+h1Jq7HfqP4Hjga+b2bejCCxG9FnXOvq8S48+79IT2mdd3oycNcJStDW6qq67rwUK/T+JlO+Zu28FLg47mBzW2Puk36HdNfZe3QrcGnYwMaXPutbR51169HmXntA+6/J95KwcKKn3fAiwPKJY8oXes/Toffr/27v3aLvK8t7j358BJAZJlFuVoFxEFIESEkWxDqmXCmrEWysoWi2SWounHce2amsVtdbWWi8ctDYKYtFyEW8E8WBtRUU5CngFgYoUSgSKXIKAICQ854+1Iott9txrX+fcK9/PGHtkr/nMudazX/ZY/PY753rn8Byr2ecYT43jNhzHaThzNk7zPZxdAOyZZLckWwGHA2e23FPXOWbDcZyG51jNPsd4ahy34ThOw5mzcZo34SzJKcD5wF5J1iY5qqrWA8cA5wCXAqdX1SVt9tkljtlwHKfhOVazzzGeGsdtOI7TcNoeJ298LkmS1CHzZuZMkiRpc2A4kyRJ6hDDmSRJUocYziRJkjrEcCZJktQhhjNJkqQOMZzpV5JsSPK9ga83tt0TQJKrkvwwyYokn+33dkWSWwd6PWicY1+d5OQx23ZKckOSLZOcluTmJM+fm59GUtt8r1PXuc6ZfiXJ7VW1zQw/5xb9hfum8xxXASuq6saBbQcDf1ZVz53g2IcAPwaWVtVd/W3HAPtW1R/2H38COKOqPjedPiXND77X+V7Xdc6caUL9v+beluQ7/b/qHtPfvijJiUkuSPLdJIf1t78yyaeSrAG+lOQBST6U5JIkZyU5O8mLkzw9yWcHXueZST4zjT4fn+SrSS5K8sUkO1XVLcA3gecM7Ho4cMpUX0fSaPK9Tl1hONOghWOm+l8yULuxqg4A/gn4s/62vwL+o6oeD/w28A9JFvVrTwJ+v6qeBrwQ2BXYF3h1vwbwH8Bjk+zQf/wq4GNTaTzJA4EPAC+qquXAJ4B39Mun0HuTIsku/V6+NpXXkTQSfK9Tp23RdgPqlDurav9xahv/yruI3hsQwO8Az0uy8Q1sa+AR/e//rapu7n//W8Cnqupe4PokXwGoqupfI3Fkko/ReyN7xRR7fyzwOODLSQAWAGv7tTOB45JsA7yE3v3Q7p3i60ia/3yvU6cZzjSsX/b/3cB9vzeh99fb5YM7JjkQuGNwU8PzfgxYA9xF701tqtdsBPhBVT1lbKGq7kjyZeAwen9V/tEUX0PS6PO9Tq3ztKam4xzgden/+ZZk2Tj7nQe8qH89xk7AwRsLVXUtcC3wZuCkafTyI2DnJE/o97JVkscN1E8B/hxYUlUXTON1JG1+fK/TnDKcadDY6zD+boL93wFsCfwgycXcd93DWJ+mN+1+MfDPwLeAWwfqnwSuqaofTbXxqvol8GLgvUm+D3wXOHBgl/9L7zTEqVN9DUkjw/c6dZpLaWhOJNmmqm5Psh3wbeDJVXV9v3Y88N2qOmGcY69izMfLZ7g3P14uaUb4XqeZ4MyZ5spZSb4HfB14x8Cb1UXAfvQ+cTSenwH/nmTFTDeV5DTgyfSuA5Gk6RrqvS7JXyb56JhjW3uvS/KyJF+a6dfV1DhzJknqpP5M0k70Ls7f6KSqOqadjoaT5FzgicA9QNFbHPZTwPv6pyWlRs6cSZK6bGVVbTPwNePBLMlsrFxwTFU9GHgY8Hp6n548e+OHCuZaevx//jzhfyhJ0rzTX53/vCTvSXJLkv9KcuhAfXGSE5Jcl+SnSf4myYKBY7+R5H1JbgaOTbIgyT8mubH/XMckqSRbJPnd/mnJwdd/fZIJr92qqjuq6lzgefTWN3tO//hj+9eAkWTrJJ9IclOSdendiWCnfu3cJO9K8u307rH5+SQPHejjiUm+2T/u++nd7omBY9+Z5BvAL4Dd+z/7lUlu6/+cLxscz4FjD+r3cWv/34PGPO87+mN4W5IvJdl++P96mojhTJI0Xx0IXA5sD7wbOGFgZurjwHrgUcAyegvJvnrMsVcCOwLvBI4GDgX2Bw4ABm8QfiawW5LHDmw7ErjfjcabVNV/AxcCv7Y+GfD7wGJgF2A74DXAnQP1VwB/ADy8/zMdB5BkZ+ALwN8AD6V3R4NP5747EQC8HFgFPJjeNW3HAYf2Z/UOAr43tpl++PtCf9/tgPcCX+h/yGGjl9K708GOwFbcdzcFzQDDmSSpyz7XnxXa+HX0QO3qqvpIVW2gF8YeBuzUn3U6FPjT/szVDcD76N/aqO/aqvo/VbW+qu4Efg/4QFWt7d+n8lfLa/SvEzuNXiCjv67YrsBZk/xZrqUXosa6h14IelRVbaiqi6rq5wP1k6vq4qq6A/hr4Pf6s4BHAmdX1dlVdW9V/Ru9APjsgWNPqqpL+overgfuBfZJsrCqrquqSzbRz3OAH1fVyf3xOQW4DFg5sM/Hquo/+2N3Or1QqxliOJMkddnzq2rJwNdHBmrXb/ymqn7R/3Yb4JH01iW7bmOoo7fu2I4Dx14z5nUePmbb2PrHgZf2Z+ZeTu/WSJO9uH9n4OZNbD+Z3kK3pya5Nsm7k2w5Ti9X0/vZtqf3c/7uYHildwuph23q2H64ewm9mbnrknwh/Zu7j/Hw/usMurrf/0bXD3z/C3rjrhliOJMkjZpr6N2GafuBULdtVQ2upD92qYLrgKUDj3cZLFbV/wPupnda8qVM4pQm/OpG5MvpLbFxP1V1T1W9rar2pneq8bnc/96bg708gt5M2430fs6Tx4TXRVU1uKju/X7Oqjqnqp5JL8BdBgyG3Y2upRf8Bj0C+OkQP6pmgOFMkjRSquo64EvAPybZNr3bKe2R5KkNh50O/EmSnZMsAd6wiX3+BTgeWF9V522i/muSPKj/up+ntyjt2ZvY57eT7Ns/VflzeuFrcPmQI5PsneRBwNvpLSS7gd6aaSuTPKv/gYatkxycZOnY1+i/zk5JnpdkEb3wevuY19nobODRSV7a/0DES4C9mfxpXE2R4UyS1GVrktw+8PXZIY97Bb0L1X8E3AKcwf1P9431EXqB7gf0bol0Nr1rtAbDy8nAPgw3a3Z8ktuA/wHeT+/WTodU1b2b2Pc3+v39HLgU+Cr3X5j7ZHr347we2Br4XwBVdQ29m5z/Jb2L/a+hd1/N8f7f/gB6y3pcS+/06lOB147dqapuojd793rgJuAvgOfO1p0L9OtchFaSpDH6y3J8uKoeObBtIXADcEBV/XiO+jgX+ERVjb2bgEaYM2eSpM1ekoVJnt0/jbcz8FZg7CzdHwEXzFUw0+ZrNlZFliQBSZ5Pb1mCHYEPVpX3LuyuAG+jt2TGnfTW+XrLr4q9W0mF+69/Js0KT2tK0iQkOZHe9Tg3VNU+A9sPAT4ALAA+OviJuSQPAd5TVUfNdb+S5h9Pa0rS5JwEHDK4of8puw/SW/h0b+CIJHsP7PLmfl2SJuRpTUmahKr6WpJdx2x+AnBFVV0JkORU4LAkl9Jbaf6LVfWdTT1fklX0bq/DokWLlu/1mD1nq3VJQ1h/76ZWF7nPxWt/0VgHuPdnP7mxqnaYcMdxGM4kafp25v6ruK+ld+/G1wHPABYneVRVfXjsgVW1GlgNsHzFsvrGt86d/W4ljevmX65rrD/6DZv8O+t+7jj+hWPvsDAphjNJmr5sYltV1XH0b1ItScPymjNJmr613P8WO0vpLfQpSZNmOJOk6bsA2DPJbkm2Ag4Hzhz24CQrk6xet+7WWWtQ0vxhOJOkSUhyCnA+sFeStUmOqqr1wDHAOfRuv3N6VV0y7HNW1ZqqWrVkyeLZaVrSvOI1Z5I0CVV1xDjbz2YTN7WWpMly5kySJKlDDGeS1DKvOZM0yNOaktSyqloDrFm+YtnRbfcijbqrb//vxvr+R3+zsf6Wv3rMhK/xl8dPqqVf48yZJElShxjOJEmSOsRwJkmS1CGGM0lqmR8IkDTIcCZJLXMRWkmDDGeSJEkdYjiTJEnqENc5kyRJI+OWu5uv3Xz8X1/RWH/XW/dsrD/qwXdPuqfJcuZMklrmBwIkDTKcSVLL/ECApEGGM0mSpA4xnEmSJHWI4UySJKlDDGeSJEkdYjiTJEnqENc5k6SWJVkJrNx9j93abkXqvKtuu7qxvuwFn2us/+E7n9VYP3KPXRrr22y5qLE+E5w5k6SWuZSGpEGGM0mSpA4xnEmSJHWI4UySJKlDDGeSJEkdYjiTJEnqEMOZJElShxjOJKllSVYmWb1u3a1ttyKpA1yEVpJaVlVrgDXLVyw7uu1epK5b9rtnNda/ctrTGusPfWDzIrJzscjsRJw5kyRJ6hDDmSRJUocYziRJkjrEcCZJktQhhjNJkqQOMZxJkiR1iOFMkiSpQ1znTJIkzZmrbru6sb7sVV9rrP/HJ5/SWN9/u30n3VPXOHMmSS3zDgGSBhnOJKllVbWmqlYtWbK47VYkdYDhTJIkqUMMZ5IkSR1iOJMkSeoQw5kkSVKHGM4kSZI6xHXOJEnSjLn+Fzc01pcd+eXG+nve/ZvNx2+336R7mm+cOZMkSeoQw5kkSVKHGM4kSZI6xHAmSZLUIYYzSZKkDjGcSZIkdYjhTJIkqUNc50ySWpZkJbBy9z12a7sVaUI/vePaxvo+L/xMY/30Ex/fWH/mzism3dOoceZMklpWVWuqatWSJYvbbkVSBxjOJEmSOsRwJkmS1CGGM0mSpA4xnEmSJHWI4UySJKlDDGeSJEkd4jpnkiQJgJ/ecd2E++zzkjWN9S9/8imN9b0WP2pSPW2OnDmTJEnqEMOZJElShxjOJEmSOsRwJkmS1CGGM0mSpA4xnEmSJHWI4UySJKlDDGeSJEkd4iK0kiRtJn56x7WN9X2ef8aEz3H2BIvMPmrb3Rrr22y5aMLX2Nw5cyZJktQhhjNJkqQOMZxJkiR1iOFMkmZJkt2TnJBk4gt5JKnPcCZJk5DkxCQ3JLl4zPZDklye5IokbwSoqiur6qh2OpU0XxnOJGlyTgIOGdyQZAHwQeBQYG/giCR7z31rkkaB4UySJqGqvgbcPGbzE4Ar+jNldwOnAocN83xJViW5MMmFP/vZTTPcraT5yHXOJGn6dgauGXi8FjgwyXbAO4FlSd5UVe8ae2BVrQZWAyxfsazmolmNrgnXMTvsU431r57+9AlfY8eFOzbWF2+17YTPoWaGM0mavmxiW1XVTcBr5roZSfObpzUlafrWArsMPF4KNE9hSNI4DGeSNH0XAHsm2S3JVsDhwJnDHpxkZZLV69bdOmsNSpo/DGeSNAlJTgHOB/ZKsjbJUVW1HjgGOAe4FDi9qi4Z9jmrak1VrVqyZPHsNC1pXvGaM0mahKo6YpztZwNnz3E7kkaQM2eSJEkdYjiTpJZ5zZmkQZ7WlKSWVdUaYM3yFcuObrsXddt/3XZVY/2AF5/VWD/vU89srD/uId7YogucOZMkSeoQw5kkSVKHGM4kSZI6xHAmSS3zAwGSBhnOJKllLkIraZDhTJIkqUNcSkPSyEjy0CF2u7eq1s16M5I0RYYzSaPk2v5XGvZZADxibtqRJuf7N/+wsX7wEV9trH/1lIMb67tvu/tkW1ILDGeSRsmlVbWsaYck352rZiRpKrzmTNIoedIM7TOn/LSmpEGGM0mj5D1Jnty0Q1XdNVfNDMtPa0oaZDiTNEp+TC+gXZXk75Ps33ZDkjRZhjNJI6OqPlBVTwKeCtwMfCzJpUnekuTRLbcnSUMxnEkaOVV1dVX9ff/DAS8FXgBc2nJbkjQUw5mkkZNky/5F9p8Evgj8J/CiltuSpKG4lIakkZHkmcARwHOAbwOnAquq6o5WG5tAkpXAyt332K3tViR1QKqq7R4kaUYk+Qrwr8Cnq+rmtvuZrOUrltU3vnVu221oFl2+7vLG+hOP+PfG+nn/enBj/ZHbPLKxvs2WixrrmhkLt1hyUVWtmOrxzpxJGhlV9dsA6TkS2L2q3p7kEcBvVNW32+1QkibmNWeSRtGH6C02e0T/8W3AB9trR5KG58yZpFF0YFUdsPFWTVV1S5Kt2m5KkobhzJmkUXRPkgVAASTZAbi33ZYkaTiGM0mj6Djgs8COSd4JnAf8bbstSdJwPK0paWQk2aKq1lfVJ5NcBDwdCPD8qnIRWknzguFM0ij5NnAAQFVdBlzWbjvDcZ0zSYMMZ5JGSdpuYCqqag2wZvmKZUe33Yum5/Jb/7OxPtE6Zhec9qzG+tJFSxvrWy94YGNd84PhTNIo2SHJ/x6vWFXvnctmJGkqDGeSRskCYBvm6QyaJIHhTNJoua6q3t52E5I0HS6lIWmUOGMmad4znEkaJU9vuwFJmi7DmaRR8uWJdkjynbloRJKmymvOJI2Sxyb5QUM9wOK5akaSpsJwJmmUPGaIfTbMehcaWZfccklj/beOOLexft4pBzfWd37QwxvrrmO2eTCcSRoZVXV12z1I0nR5zZkktSzJyiSr1627te1WJHWA4UySWlZVa6pq1ZIlXg4nyXAmSZLUKYYzSSMtyfOSnJ7k1CSHtd2PJE3EcCZp1D23qn6vqg4HDmm7GUmaiJ/WlDTqFiZ5RP/7Ra12IklDMJxJGnXHAq/rf+9N0dXoOzc2rWEMT3/l+Y31809tvoPY7tvu0Vjf6gFbNta1eTCcSRp1O1XVnwMkeSJwRcv9SFIjrzmTNOpeMPD981rrQpKG5MyZpFG3U5I9gAKa740jSR1gOJM06t4M/HH/+7e22YgkDcPTmpJG3QuAh1bVG4BXtN2MJE3EcCZp1O0BXNP//sFtNiJJwzCcSRp1RW+ts33wmjNJ84DhTNKo+0cgwMuBN7XcyyYlWZlk9bp1t7bdiqQOSFW13YMkzYkkz6mqL7Tdx3iWr1hW3/jWuW23MdK+f/MPG+sH/37zIrNfOelJjfX9t9t30j1p9CzcYslFVbViqsf7aU1JIyvJXwD7A2fRu6/mt9rtSJImZjiTNMoeW1UvTXIe8IyquqvthiRpIl5zJmmUbZ/k2cCNwNP630tSpxnOJI2MJI8bs+kMYAfgs/1/d5jzpiRpkjytKWmUnAwcAJDk1VX10Y2FJA+qql+01pkkDcmZM0mjJAPfv3ZM7etz2YgkTZXhTNIoGVwbKGNqvt9Jmhc8rSlplPxGklcC3+fXw5mLOm4GLrnlR431idYx+/rHD2qs7/PQfRrr6+9d31jf4gH+b1cT87dE0ig5FlgBvApYmuQS4LL+1/Yt9iVJQzOcSRoZVbV68HGSpcB+wL7A11ppSpImyXAmaWRV1VpgLXB2271I0rC8QFaSJKlDDGeSJEkdYjiTJEnqEMOZJElSh/iBAEnSvHHZussa67/1suYP5Z77iSc31h+z5DGT7mmQ65hpJjhzJkmS1CGGM0mSpA4xnEmSJHWIJ8claRYkWQR8CLgbOLeqPtlyS5LmCWfOJGlISU5MckOSi8dsPyTJ5UmuSPLG/uYXAmdU1dHA8+a8WUnzluFMkoZ3EnDI4IYkC4APAocCewNHJNkbWApc099twxz2KGmeM5xJ0pCq6mvAzWM2PwG4oqqurKq7gVOBw+jd03Npfx/fayUNzWvOJGl6dua+GTLohbIDgeOA45M8B1gz3sFJVgGrAHZ5xC6z2Gb77t5wd2P9J7f9ZMLnOOjlX2+sn/fJpzTW91z86Ma665SpC/wtlKTpySa2VVXdAbxqooOrajWwGmD5imU1w71Jmoecapek6VkLDE55LQWubakXSSPAcCZJ03MBsGeS3ZJsBRwOnNlyT5LmMcOZJA0pySnA+cBeSdYmOaqq1gPHAOcAlwKnV9Ulk3zelUlWr1t368w3LWne8ZozSRpSVR0xzvazgbOn8bxrgDXLVyw7eqrPIWl0OHMmSZLUIYYzSZKkDjGcSVLLvOZM0iCvOZOklo3KNWf33HtPY/2Ht1zaWH/Gqu9O+Br/duKBjfW9Fu/VWHeRWc0HzpxJkiR1iOFMkiSpQwxnkiRJHWI4k6SW+YEASYMMZ5LUsqpaU1WrlixZ3HYrkjrAcCZJktQhhjNJkqQOccEXSRIw8Tpld234ZWP9whub1zF74Z9c3lj/9PGPa6wD7Lfd3o111zHTKHDmTJJa5gcCJA0ynElSy/xAgKRBhjNJkqQOMZxJkiR1iOFMkiSpQwxnkiRJHWI4kyRJ6hAXhJGkliVZCazcfY/dWu1jfW1orF9042WN9Rf+2ZWN9U+/f6/G+tMe/vjGurS5cOZMklrmUhqSBhnOJEmSOsRwJkmS1CGGM0mSpA4xnEmSJHWI4UySJKlDDGeSJEkd4jpnktSyuVrn7M4NdzXWz1l7UWP9VW+7sbH+mffs3lg/aMf9GuuSepw5k6SWuc6ZpEGGM0mSpA4xnEmSJHWI4UySJKlDDGeSJEkdYjiTJEnqEMOZJElShxjOJEmSOsRFaCVpRNx2z+2N9W/8z8WN9Ve/75eN9RPeskNj/Yk77NNYf+CCBzbWJfU4cyZJLUuyMsnqdetubbsVSR1gOJOklnmHAEmDDGeSJEkdYjiTJEnqEMOZJElShxjOJEmSOsRwJkmS1CGucyZJ88S6u3/eWP/69c3rmL3qvXc11t/92gc31n9n570b6wu3WNhYlzQcZ84kSZI6xHAmSZLUIYYzSZKkDjGcSZIkdYjhTJIkqUMMZ5IkSR1iOJOkliVZmWT1unW3tt2KpA5wnTNJallVrQHW7L/8N4++/Z47xt1vonXM/vjDGxrrx656SGP9Bbvu3ljfeoutG+uSZoYzZ5IkSR1iOJMkSeoQw5kkSVKHGM4kSZI6xHAmSZLUIYYzSZKkDjGcSZIkdYjrnElSR9zyyzv59FWXjFv/q4+ubzz+mJc2r2N2xB5LG+vbbrlNY31BFjTWJc0MZ84kSZI6xHAmSZLUIYYzSZKkDjGcSZIkdYjhTJIkqUMMZ5IkSR1iOJMkSeoQ1zmTpI64/mb421M3jFt//cu3bTz+ZY96WGN98ZYPbqy7jpnUDc6cSZIkdYjhTJIkqUMMZ5IkSR1iOJMkSeoQw5kkzZIkuyc5IckZbfciaf4wnEnSJiQ5MckNSS4es/2QJJcnuSLJG5ueo6qurKqjZrdTSaPGpTQkadNOAo4H/mXjhiQLgA8CzwTWAhckORNYALxrzPF/UFU3zE2rkkZJqqrtHiSpk5LsCpxVVfv0Hz8JOLaqntV//CaAqhobzMY+zxlV9eJxaquAVf2H+wAXb2q/GbQYuHUWj5tov/Hqw26f6PH2wI1D9DlVUx2/yRw7k2M4zLbNbQwnW5vKGO5VVc0LCzapKr/88ssvvzbxBewKXDzw+MXARwcevxw4vuH47YAPAz8B3jTE6104Bz/T6tk8bqL9xqsPu32Ix7M6hlMdv7bGcJhtm9sYTrbWxhh6WlOShpdNbBv39ENV3QS8ZvbamZI1cOAnZAAABMxJREFUs3zcRPuNVx92+0SPZ9t0Xq+NMRxm2+Y2hpOtzfkYelpTksYxU6c1J/F6F1bVipl4rs2VYzh9juH0TXcM/bSmJA3vAmDPJLsl2Qo4HDhzBp9/9Qw+1+bKMZw+x3D6pjWGzpxJ0iYkOQU4mN7F0f8DvLWqTkjybOD99D6heWJVvbO9LiWNIsOZJElSh3haU5IkqUMMZ5IkSR1iOJMkSeoQw5kkzQPeRH1qkixK8vEkH0nysrb7mY/83Zu+JM/v/w5+PsnvTLS/4UySZpk3UZ9ZkxzPFwJnVNXRwPPmvNmOmswY+ru3aZMcw8/1fwdfCbxkouc2nEnS7DsJOGRww8BN1A8F9gaOSLJ3kn2TnDXma8e5b7nTTmLI8QSWAtf0d9swhz123UkMP4batJOY/Bi+uV9v5O2bJGmWVdXX+ncbGPQE4IqquhIgyanAYf27DTx3bjucXyYznsBaegHtezgh8SuTHMMfzW1388NkxjDJpcDfAV+squ9M9Nz+okpSO3bmvhkd6IWIncfbOcl2ST4MLNt42yjdz3jj+RngRUn+ibm/h+R8s8kx9HdvUsb7PXwd8AzgxUkmvN+uM2eS1I5RuIl6l2xyPKvqDuBVc93MPDXeGPq7N7zxxvA44Lhhn8SZM0lqx1pgl4HHS4FrW+plFDie0+cYTt+MjKHhTJLaMds3Ud/cOJ7T5xhO34yMoeFMkmZZ/ybq5wN7JVmb5KiqWg8cA5wDXAqcXlWXtNnnfOF4Tp9jOH2zOYbe+FySJKlDnDmTJEnqEMOZJElShxjOJEmSOsRwJkmS1CGGM0mSpA4xnEmSJHWI4UySJKlDDGeSJPUl+cMk1yX53sDXvjP4/LsmubP/vNsNvMb1SX468HircY4/N8mzxmz70yQfSrKwf+zdSbafqZ4197zxuSRJ99kPeHNVnTCLr/GTqtq///3+AEmOBW6vqvdMcOwp9G4JdM7AtsOBP6+qO4H9k1w1s+1qrjlzJknSffYFvtd2EwBJjkzy7f5s2D8nWQCcATw3yQP7++wKPBw4r71ONdMMZ5Ik3edxwMcGTi+uaqOJJI8FXgI8uT/LtgF4WVXdBHwbOKS/6+HAaeW9GEeKpzUlSQKS7ALcUFX7DWxbmOTD9GanHgJcAvxDVf0kyQOq6t5ZaufpwHLggiQAC4Eb+rWNpzY/3//3D2apB7XEcCZJUs9+wGWDG/rXcb0mycHAPlV1fJJXJnkbcGGSdcCNVXVWklOBNwCvB0Lv2rL3T7GXAB+vqjdtovY54L1JDgAWVtV3pvga6ihPa0qS1LMvY8JZgy+OE7xeC9wJ3NR/vqn6d+DFSXYESPLQJI8EqKrbgXOBE+nNomnEOHMmSVLPvsBTkxzaf1zAU/phaKxb+//+kvv+X7qI3qTHyVX1g+k0UlU/SvJm4EtJHgDcA/wxcHV/l1OAz9A7rakRYziTJAmoqpdN4bCvAu9OshuwBDge+Nsk1wG3VdXbhnztYzex7TTgtHH2/yy9U58aQfEDHpIkzY3+hw6+Cdw0sNbZTD33QuB8YAdg36q6eSafX3PHcCZJktQhfiBAkiSpQwxnkiRJHWI4kyRJ6hDDmSRJUocYziRJkjrEcCZJktQhhjNJkqQOMZxJkiR1yP8HbV/jiVzYskoAAAAASUVORK5CYII=\n", "text/plain": [ "
" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "sim.obs.peek()\n", "print(sim.obs)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Spectral analysis\n", "\n", "Now that we have some simulated CTA counts spectrum, let's analyse it." ] }, { "cell_type": "code", "execution_count": 10, "metadata": {}, "outputs": [], "source": [ "# Fit data\n", "fit = SpectrumFit(obs_list=sim.obs, model=model, stat=\"cash\")\n", "fit.run()\n", "result = fit.result[0]" ] }, { "cell_type": "code", "execution_count": 11, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "\n", "Fit result info \n", "--------------- \n", "Model: PowerLaw\n", "\n", "Parameters: \n", "\n", "\t name value error unit min max\n", "\t--------- --------- --------- --------------- --------- ---\n", "\t index 2.116e+00 3.446e-02 nan nan\n", "\tamplitude 2.382e-12 1.205e-13 1 / (cm2 s TeV) nan nan\n", "\treference 1.000e+00 0.000e+00 TeV 0.000e+00 nan\n", "\n", "Covariance: \n", "\n", "\t name index amplitude reference\n", "\t--------- ---------- ---------- ---------\n", "\t index 1.188e-03 -9.160e-16 0.000e+00\n", "\tamplitude -9.160e-16 1.451e-26 0.000e+00\n", "\treference 0.000e+00 0.000e+00 0.000e+00 \n", "\n", "Statistic: -1591.397 (cash)\n", "Fit Range: [ 0.1 100. ] TeV\n", "\n" ] } ], "source": [ "print(result)" ] }, { "cell_type": "code", "execution_count": 12, "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAaMAAAEOCAYAAAAkF3jEAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMi4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvIxREBQAAIABJREFUeJzt3Xd4W/eV5//3lwALwN5BkRS7rOJuR+69yY5jO47j2JKcMklmf9lJduf5TdvsZHZ2N5tJ9jc7eTKzM5vZzKRLtuOS6okludux427HLbZYxA6wgJ0oRDm/PwDKNEWK7YIAyPN6Hj8RgIuLQyPUx/fec8/XiAhKKaVUMmUkuwCllFJKw0gppVTSaRgppZRKOg0jpZRSSadhpJRSKuk0jJRSSiWdhpFSSqmk0zBSSimVdBpGSimlkk7DSCmlVNLZk11AuigrK5P6+vpkl6GUUmnllVdeGRaR8qW20zBapvr6el5++eVkl6GUUmnFGNO1nO30NJ1SSqmk0zBSSimVdBpGSimlkk7DSCmlVNJpGCmllEo6DSOllFJJp2GklFIq6TSM1sHU1BQikuwylFIqZelNr+vA4/EwMzNDbW0t+fn5yS5HKaVSjh4ZrRO/38/Ro0fp6OhgZmYm2eUopVRK0SOjdTY6Osr4+DgulwuXy4UxJtklKaVU0umRURJEo1H6+/t5++23GRsbS3Y5SimVdBpGSRQMBmlvb6e1tZVAIJDscpRSKmk0jFLAxMQE77zzDr29vUQikWSXo5RS607DKEWICAMDA7z99tt4vd5kl6OUUutKwyjFhEIhOjs7effdd5menk52OUoptS40jFLU9PQ07777Lp2dnYRCoWSXo5RSCaWt3SnO6/UyNjZGVVUVFRUV2gqulNqQ9MgoDUQiEXp7e3n77bcZHx9PdjlKKWU5DaMEe6Z1iK8/0c+7wzNrnk8XDAZpa2ujra1NW8GVUhuKnqZLsP4xP7/tnuax9igNRXb2NDm5ZKuDbPvqT7eNj48zMTFBRUUFVVVV2Gw2CytWSqn1p0dGCfaJej+Pnvsif3RmFhGBb78ywecfGuT7r0/QPxle9X5nW8HfeusthoeHLaxYKaXWn9GlDZbn3HPPlZdffnnlb3zqb+GJ/0Eoq4iBxtv4Td71/KrL8HxvgIjAGZVZXN/s5OyqbGxraE5wOp3U1taSl5e36n0opZTVjDGviMi5S26nYbQ8qw4joPe5+8l/7TsUDr1IODOPwYZbObrlZg71ZHKkw8eIP0q5M4NrG51c1eikMHv1B6zFxcXU1NSQlZW16n0opZRVNIwstpYwamtrY3x8HOfYe7haD1DseZaIzcFQ/U3019/GcyO5HGrz8dbQDPYMuLAmh+ubnbSUZK6qlTsjI4PKykpcLhcZGXomVimVPBpGCzDGNAJ/CRSKyG3GmFuADwMVwD+JyJHF3mtFGM3KmTiGq+1uSvqeQDLsDG+9AU/zHXSESjjc7uPJTj/+sNBYZGdPs5OLa1fX8JCZmUl1dTWlpaWrqlsppdYq6WFkjMkBngayiXXtPSAif73KfX0PuBEYFJFT5722B/h7wAb8q4h8Yxn7e0BEbpvzuBj4XyLy2cXeY2UYzcqe6sXVdg+lvUcQDN7aa/E038l4VhVPdQc41OajZyJMXqbhigYH1zU5qcpbeQNkbm4uNTU1ej1JKbXuUiGMDJArIlPGmEzgN8B/FJHn52xTAfhFZHLOc80i0jZvX5cCU8CP5oaRMcYGHAWuAXqBl4A7iQXT1+eV9AciMhh/3/ww+jvgoIi8utjPk4gwmpXl81DZ/hPKun+NiUYYqb4ST8te/Hl1vDMc4lCbjxf6Yg0PZ8YbHs5aRcNDSUkJ1dXVej1JKbVulhtGCbvPSGIpNxV/mBn/Z37yXQZ8wRhzg4gEjDGfBz4K3DBvX08bY+oX+JjdQJuIdAAYY+4FbhaRrxM7kjqpeGB+A3j4ZEGUaDNOFz2n/UfcLfupbL+P8q5fUdL3GGNVl+Bs2ceuC1oY8Ud4tMPPkQ4fX392jAqnjWubHFzV4KRgmQ0PIyMjjI2N6fUkpVTKSehNr/Ejl1eAZmLXZF6Y+7qI3G+MaQDuNcbcD/wBsaOc5aoGeuY87gXOO0k9pcDXgLOMMV8GpoGrgcL4Edk/L/CejwAfaW5uXkFZqxPOKaVv1xfwNO+l8tiDVBz7GcXupxmvOA93y35u37WLW3fk8mJfkEPt0xx4c4qfvD3FRbU5XNe0vIaHaDSK2+1meHiYLVu2UFZWlvCfSymllrIuDQzGmCLgZ8CXROStBV6/l9jRUJOIDC2yj3rgoXmn6T4OXCcin4s/vgvYLSJfsvpnSORpusVkhKaoOPZzKjsewB6aYKLsLNwt+5kqPROMoXs8xOF2P092+QmEhabi2ISHi7Y6yLYt7xSe0+mkpqaG/Pz8FdenlFJLWe5punU5TyMiY8CTwJ75rxljLgFOJRZWK21w6AVq5zyuAfpXV2XqiWbm4dm2nzevvoeenV/AMdnFKb/9E0559ksUDDzP1gI7nz+7gH+5sZzPnZVPMCz808sT/OFDg/zwdxN4ppae8ODz+Th69Cjt7e0Eg8F1+KmUUupEiWxgKAdCIjJmjHEAR4D/KSIPzdnmLOAeYu3Vx4ADQIeIfGWB/dVz4pGRnVgDw1VAH7EGhr0i8rbVP08yjozmM5EZSnsextV2D9n+QXwFzbhb9jNWdTGYDESEt4ZmONTm48X+ICJwpive8ODKJmOJU3jGGMrLy6mqqsJu17GFSqm1S4VuutOBHxLrbMsA7hOR/z5vm4uACRF5M/44E/i0iPzLvO3uAS4HyoAB4K9F5Lvx124AvhX/nO+JyNcS8fOkQhgdFw1T2vsIrrZ7yJnuxZ9Xh6dlHyNbroCM2NBUrz/CIx0+HunwMxaIUpFr47omJ1fVO8hfouHBZrPp+klKKUskPYw2mpQKo1kSobj/KapaD+KYPEbQWYWn+U68Ndcitlj7djgqvNAXu2fpneEQWRlw0VYHe5qcNJdknnT32dnZVFdXU1xcbH3tSqlNQcPIYikZRrMkSuHAb6k6eoDc8feYySnH0/wJhrd+GLFlH9+sezzEoXYfT3UGCESE5uJM9jQ7uag2h6yTNDzk5eVRU1NDbm5u4n4GpdSGpGFksZQOo1kiFAy9jKv1APkjbxLKKmag6eMM1d9E1O48vpkvFOXJLj+H23z0TkbIzzJc1eDk2iYHlbmLXysqLi6murqa7OzsRbdRSqm5NIwslhZhNEee93dUHT1AwfArhDPzGWy4lcGGW4lkvd/CvVDDw9lV2expcnKmK2vBhgdtclBKrYSGkcXWEkZDQ0N4PB5mZmYsrmppztHfU9V6kKKB54jYnQzW38xg422Esz94Hcjri3Ckw8ejHX7GglFc8YaHKxoc5Ged2PBgs9lwuVxUVFToJAel1KI0jCy2ljCaNTk5ycjICKOjo0QiEYsqWx7HRDuu1rsp7n8SychiqO7DDDR9gpCj/APbhaLCC70BDrX7+H284eHirQ6ub3bSWHxiw0NWVhZVVVWUlpZq551S6gQaRhazIoxmiQjj4+N4vV7Gx8dZz+8ge6obV+s9lPY9ghgb3prr8DTfwUzulhO27RwLcbjdx9NdsYaHlpJYw8OFNSc2PDgcDqqrqyksLFyvH0UplQY0jCxmZRjNFYlEGB0dxev1MjU1tfQbLJLl88SWr+g5hJEII9VX427eSzB/6wnbToeiPNnp51C7j/7JCAVZhisbnFzX5KBiXsODdt4ppebSMLJYosJormAwyMjICCMjIwQCgYR+1qzMwHBsUnjnrzDRGUarLsXTsh9/YdMJ24oIbwzGGh5e7g8ixBoerm92ckblBxseioqKqK6uJicnZ11+DqVUatIwsth6hNFc09PTeL1eRkdHCYeXnjG3VvbgGBUdD1DR+XNsYR9jlefjbrkLX/GOBbcfntPwMB6M4sqLNTxcWe8gL97wYIyhpKSELVu26BpKSm1SGkYWW+8wmrXe15dsM5OUd/6cyo4H45PCz4lPCj8dFmhQCEWE5/sCPNzm4z1viCwbXBKf8DDb8KDt4EptXhpGFktWGM0VDoePX1+anp5O6GdlhP2Ud/2Syvb7yAyOMlV8Ku5t+5ko/9CCoQRwbCy2Ku0z3QGCEeGU0kz2NDm5oCaHTJvBZrNRUVFBZWUlNpstofUrpVKDJWEUH3a6lJCI/H4lxaWjVAijuQKBAF6vF6/XSygUStjnmEiQsu5f42r7CVmBQaYLt+Fu2c+460IwC99fND0T5Yl4w4N7KkJBdgZXNzi4tslJudOG3W6nsrJS71FSahOwKowmgdeAk91AUisi9SuuMM2kWhjNNTExgdfrZWxsjGg0mpDPMNEQJb2P4Gq9mxxfP/78etwt+xjdcjmYhY9yoiK8MTDDoXYfr/TH1ko6Z0tswsPplVlkZ2XhcrkoLy/Xe5SU2qCsCqOnReTSJT5oyW02glQOo1mRSISRkZHEnsaLRijuf5Kq1gM4proI5Fbjad7LSM3VSMbiU8CHfBGOtPt49JifiWCUqjwbe5qcXFHvoDgvR2+cVWqD0mtGFkuHMJor4afxJEqR+zdUtR7AOdFG0FHBQNMdDG+94fjyFQsJRYTfxic8vOcNkW0zXLI1hz3NTrZX5lFVVUVJSYmGklIbhKVhZGJ/M5wKbAH8wNsi4l1zlWkk3cJologwMTHB8PBwYrrxRCgYfJGq1gPkjb5NKLsET9PtDNd9hKjdcdK3dozGJzx0+5mJwCmlmVzf7OSyxkK21mzRUFJqA7DqNF098OfAHmLLgg8BOUALMAb8M3BANsHhVbqG0VzhcBiv18vw8LD1N9WKxCaFtx6gYPhVwpkFDDR+jKGGjxLJzDvpW6fmNDx4piIUZmdwVYODm3YWc3pzLcXFxRpKSqUpq8LoPuDbwFMiEp33WhWwDxgWkR+srdzUtxHCaK7p6WmGh4cZGRmxvOkhd/QdXK0HKBp4nog9l8H6Wxhs/Bjh7KKTvi8qwu8GYhMeXnEHMcC5W7K5aWcR15/VoEdKSqUhvWZksY0WRrOi0SgjIyMMDw9b3vTgGG+jqvUARe5niNqyGa77CANNHyeUU7bkewenwxzp8PNYh4+JGWFLvo0bTylg74VN1FVVaCgplSasvmZ0K/CIiEwaY/4TcDbwNyLy+tpLTQ8bNYzm8vv9DA8P4/V6LV3iImeyC1fb3ZT0PYYYG8O11zPQfAczTteS752JCM/1BjjU5qN1JNbwcEVDLp+6qIHzt9dqKCmV4qwOozdE5HRjzIXA3wLfBP5MRM5fe6npYTOE0axoNMro6CjDw8OWThLPmu7D1XYvpT2HMQjemqvxNO8lmFe7rPe3j8YmPPym289MFHaWZ7P3Q9Xcdn4LOVk6ZkipVGR1GL0mImcZY/6GWCfdwdnnrCg2HWymMJorEUdLmf7B2KTwrocw0TCjWy7D07IXf8GJk8IXMjkT5Yljfg63+/BMRyjKyeCjp1fw2cu3U1OiS1colUqsDqNfE+um2wOcC0wDL4nIGWstNF1s1jCaNXu0NDQ0ZNm1JXtwhMqOByg/9gtsET9jlRfh3rYPX9H25dUUb3h4uM3Hq+4gxsCljYV89tIWLt6m15WUSgVWh1EecAPwhoi8a4zZApwhIg+vvdT0sNnDaC6fz8fQ0JBlnXi2mQkqjv2UimM/xR6aYrz8XDyzk8KXaWA6zJF2P48d8zE5I9QVZ3PXBfXcvruOgpzFJ0MopRJLu+kspmF0okgkgtfrZWhoyJL7ljJC0/FJ4feTOTPGZMlpuFv2M1l+7qKTwuebiQjP9sQaHtpGQ+TYM7j5jCo+c0kj210Fa65RKbUyGkYW0zA6ucnJSQYHBy2Z8mDCgdik8PZ7yQoMM110SmxSeOWFyw4lgLaREIfaYw0PoSicU1vIpy9u5LpdLrLsOi1cqfWgYWQxDaPlmZmZYWhoiOHh4TWvUGsiM5T2HsHVdg/ZPje+gkY8zfsY3XLpopPCFzIZjPJYZ6zhYXA6QlluJnvPr2fv7q24CnVZdKUSScPIYhpGKyMijIyMMDg4iM/nW9vOohFK+h/H1XoQx1Q3gdwaPC178VZfDRnLb+mOiPC6Jzbh4TVPkAwD1+5y8ckL6jm/Uac7KJUIVo0Dqgb+J1ANPAx8U0TC8dceFJGPWVRvytMwWr2pqSkGBwcZGxtb2yk8iVLkfiY+KbydoKMST/MdeGuvP+mk8IV4pmITHh4/5mdyJkpzeS6fvLCeW8+uIS9b71lSyipWhdFh4FfA88BniU3uvklERvU+I7VSs6fwhoaG1nbPkgiFg8/jOnqAvLHfM5NdykDT7QzX3bjkpPD5ghHh2R4/h9r8tI+GcGbZ+NjZNdx1QR3bKvNXX6NSCrAujD4QOMaYTwF/CtwEPCgiZ1tRbDrQMLJONBrF6/UyODi4ti48EfKHX6Oq9QD53tcJZRUy2Hgbg/U3E11iUvhCWkdip/Ce7Q0SigjnNZTwyQvquXZXJZk2bXhQajWsCqN3gLNEJDjnuT3APwJOEdliRbHpQMMoMcbHxxkYGGBycnJN+8kdeYuq1gMUDr5I2J7LYMOtDDbeSiSrcMX7mghGefyYjyMdfgamI1TmZ3PneVu5c/dWKgu04UGplbAqjP6M2KSFJ+c9fy7wtyJyxVoLTRcaRonl8/kYGBhgdHR0TdeVHGNHqWo9SLHnGSK2HIbqb2Kg8XbCOSUr3ldEhNfcQQ53+HnVHcSeYbhul4u7LqjjvAZteFBqObSbzmIaRutjZmaGwcFBhoaG1jTdIWfyGK7WuynpewLJsDG89QY8TXcQclauan+eqTCH23080RVgMhjllMp89l9Qx0fPqtaGB6VOwupxQFuBLwL1wPHfPBG5dQ01phUNo/UViUQYGhpicHCQUCi06v1kT/fharuH0p7DAHhrr8PTfCfB3OpV7S8YFn7T4+dwR4D2kRnysm3cenYNd51fR4s2PCh1AqvD6HXgR8CbwPH/XBWRx9ZSZDrRMEqO2WaHgYEBgsHg0m9YRKZvAFf7Tyjr/jdMNMJI9RV4WvYSyG9Y1f5EhNaREIfb/TzbEyAUFS5oLOWuC+q4Zqc2PCg1y+owelFEdltSWZrSMEouEWFsbAy3243f71/1fuwBL5Ud91Pe+UtskQCjrktwt+zDX7Rt1fscD0Z57JiPR48FGJgKU1mQzd7dddy5u5YKbXhQm5zVYXQXUAccBo7/56mIvLGWItOJhlHqGB8fx+PxrGnhP9vMOBUd8Unh4WnGK3bjbtnPdMmpq95nRIRX3UEOdwR4zR3AnmHYc6qLu86vY7c2PKhNyuow+irwOaCD90/TiYhcuqYq04iGUeqZnJzE4/EwMTGx6n1khKao6PwlFR33kzkzzmTpmbFJ4WVnrWgo63z9k7EJD090+pmaibLdlc/+82MND7na8KA2EavD6D3g9Ln3G202Gkapa3p6Grfbzfj4+Kr3kRH2U9b1EJXt95EV9DJVtAPPtv2MV5y/plAKhoVnevwcbvfTMRoiL9vObefUsP/8rTRXaMOD2visDqP7gS+IyLAVxaUjDaPU5/P5cLvdjI2NrXofJjJDac8hXG33ku334Ctowt2yj7GqS1Y0KXw+EeE9b4jD7T6e6w0QjsKFTaV88oI6rt5RiV0bHtQGZXUYPQ6cDrzAB68ZaWu3Sjl+vx+3283o6OjqdxINU9L3GFWtd5Mz3YM/byue5r2MVF8FGasPJYDxQIRHj/l5pMPPkC+CqyCHfedt5Y7dWynPz17TvpVKNVaH0VULPa+t3SqVWRJKEqG4/2lcbQdxTnQQdFbhab4Tb821K54UPl8kKrziDnKo3cfvBmawZxhuOK2Kuy6o49y6Ym14UBtCIm56HRSRQPyxAygTkZ41V5omNIzSl9/vp7+/f02n7xChcOC3VLX+mNyx95jJKWOg6RMMbf0wYl97+3b/ZJhD7T6e7PQzHRK2u/L45AUN3HLWFpxZ2vCg0pfVYfQycKGIzMQfZwPPbKZ7jzSM0p/P56O/v39NjQ6xSeGvUHX0APkjbxDKKmKg6eMM1d1ENDN3zTUGwlGe6Q5wqN1H51iY/Gwbt51by13n19FYvvJJ5Eolm+UTGETkzHnP/U5EzlhDjWlFw2jjmJ6epr+/f00t4QB53jdwtR6kcOglwpn5DDZ8lMGGW4lkFay5xtmGh4fbfDzfF2t4uKiplE9eWM9V2yu04UGlDavD6DHg70Tk1/HHNwJ/olO7VTqbmpqir69vTTfPAjjH3qXq6EGKBp4lYnMwVH8zA023Ec5e+aTwhYzFGx6OtPvw+qO4CrLYf349d+zeSlmeNjyo1GZ1GLUA9wClgADDwD4RaV1roelCw2jjmpiYoK+vD5/Pt6b95Ex0xJav6H8SychkqO5GBpo+QchRbkmdkajwsjvIoTYfbwzOkJlh2HNqJZ++qIGzt2rDg0pNCVlCwhhTBCAia7gSnJ40jDa+0dFR+vv717b6LJA91RObFN77CII5Pil8Jte6tSh7J2JLWjzZ6ccXFk6pcPLpixq55awaHFlraz1XykpWLa53B/ATWWQjY0w9sEVEnltlnWlDw2hzEBG8Xi9ut5uZmZk17SvL56Gy7V7Keh7GSISRLVfiadlHIL/OomrBH47yTFeAh9t9dI+HycvO4KNnbuEPLmmmoWztDRVKrZVVYfQnwF3Ai8ArwBCQAzQDlwMTwF+IyHsW1JzSNIw2l2g0ytDQEB6Ph3A4vKZ9ZQaGqWy/n7KuX5ERCTJWdQnulv34C5stqjYWou/ONjz0BogInFdXwKcvbuLaXVXYMvQUnkoOy07TGWPswDXARUAV4Ad+D/xaRI5ZUGta0DDanCKRCB6Ph8HBwTWtPAtgD45RcexBKo79HFt4mrHK8/G07Ge6eKdF1caMBiI81uHncIePEX+Uijw7d5xTzScvbqYsX5e0UOtLlx23mIbR5hYKhejv78fr9bLW3xlbaIryYz+jsuNB7KEJJsrOxt2yj6nSM9c0lHW+SFR4qT824eHNwRnsGXBVSzGfubiJ85ortOFBrQsNI4tpGCmAQCBAb2/v2m6cjYtNCv8Vrvb7yAyOMFW8C3fLPiYqzrM0lCDW8DA74cEfFppKsrjzQ9XccX4zeY61jTVS6mQ0jCymYaTmmpycpLe3d83t4DA7KfxhXG33kO0fxFfQjHvbfsZcF4Ox9uZWfzjK010BDrX56J4Ik5dluH57MZ+6sJFT6yst/SylwLoGhg+JyEuWVpamNIzUQkZGRujr61tz5x2AiYYo6X0UV9vd5Ez34c+rw92yn9Etl695Uvh8IsI7wyEOtfl4oS/W8HB2VQ53nLOFG89pwOnQa0vKGlaF0b8C5wBvA4eAwyIyZFmVaUTDSC1GRBgcHMTtdhOJRCzYYYTi/iepaj2IY7KTgHMLnpa9jNRcg2Rkrn3/84z4Izza4edIh4/RQJQKp40bdxSx97w6GmsqycjQ0UNq9ayewHAqcD1wHbHW7seJhdPzIrK2FqM0oWGklhIOh3G73QwNDa25yQEAiVLkeRZX60Fyx48yk1OBp/kTDG+9AbFZPwYoHBVe7AtyqH2at4dCZGbARbUObj/bxcU7aigoWPvMPbX5JOyakTEmF7iKWDidJyJnr67E9KJhpJYrGAzS29u7tiUr5hKhYOglqo4eIG/0LULZxQw0fpyh+puI2p3WfMY83eMhDrf7ebLLTyAsNBXb+fApBXz07Fq2VJaTna0z8dTyaAODxTSM1EpZ2eQAgAh53jeoaj1AwfArhDMLGGj8GEMNHyWSmZjlJfyhKE/FJzz0ToTJyzRc2eDgo6eVc1qDi+LiYmw2HT+kFqdhZDENI7Vaw8PD9Pf3EwqFLNunc/T3VLUeoGjgt0TsTgbrb2Gw8TbC2UWWfcZcIsLbQzMcavfxQl8QETjTlcUNLXlcsb2SivIyPY2nFqRhZDENI7UW0WgUj8fDwMDAmic5zOUYb8fVeoBi99NEbdkM193IQNPthHLKLPuM+bz+CI92+Hikwx9reMi1cV2Tkz3bCqivKqe0tJScHO3GUzFWddN9C7hbRF60srh0pGGkrDAzM0Nvby+jo6OW7jd7sjs2KbzvEcTY8NbuiU0Kd7os/Zy5wlHhhb7YPUvvDIfIyoCLtjrY0+TkjNoiSktLKSkp0dN4m5yVg1LvAEqAe4F7ROQty6pMIxpGykpTU1P09PRYdz0pLsvnxtV2L6U9hzASwVt9DZ6WOwnmbbX0c+brHg9xqN3HU50BAhGhuTiTPc1OLt7qoKK0mNLSUgoKCnQE0SZkdWt3E7FQugMwwN3AvSLSsdZC15MxphH4S6BQRG4zxtwCfBioAP5JRI4s9l4NI5UIibieBJDpH6Ky/SeUd/0bJjrD6JbL8LTsw1/QZOnnzDcdivJUp59D7T76JiPkZxmubHByXZOD6sIcSkpKKC0txelMTBegSj2JbO0+B/hX4HQRWfT42xhTC/wIcAFR4Dsi8vcr+rD39/U94EZgUEROnffaHuDvARvwryLyjWXs7wERuW3O42Lgf4nIZxd7j4aRSpTZyeADAwPW3J80hz04SkXHA1R0/gJb2MdY5YW4W/bhK95h6efMJyK8NTTDw20+XuqPNTycXZXNniYnZ7qyyHU6j5/Gy8y0/kZelTqsPjKyAdcSOzK6DniW2Cm7B07yniqgSkReNcbkE1sP6RYReWfONhWAX0Qm5zzXLCJt8/Z1KTAF/GhuGMXrOkpsiYte4CXgTmLB9PV5Jf2BiAzG3zc/jP4OOCgiry7282gYqUQLBoP09PRYMoR1PtvMJBWdP6Oi40HsoUkmys7BvW0/U6VnWP5Z83l9EY50+Hi0w89YMEplvOHhygYHBdk2CgoKKC0tpbCwUKc9bEBWXTO6gthf7jcBrxG7bvTTueGxgoJ+AfyjiDwy57mPA18AbhCRgDHm88BHReSGBd5fDzw0L4wuAP6riFwXf/xlABGZH0Tz9/VA/DSdAb4BPCIij57sPRpGar1MTEzQ09Oz5uXPF5IR9lHe+Usq2+8nc2aUyZLT8LTsZ6L8XMsnhc8Xigov9gV4uM3H7+MNDxdvdbCn2UlTcSboUUhwAAAezklEQVQ2m43i4tj1pby8xNw3pdafVWH0DLHrQw+sZSZdPEieBk4VkYl5r/05cCFwP/BF4BoRmVpkH/PD6DZgj4h8Lv74LmJTIb64SB2lwNeIHUn9KzANfIrYEdXrIvLPC7znI8BHmpubP9/a2rqin1up1bJ83t08JhKkrPvXuNruJSswxHThKbhb9jHuutDySeEL6RwLcbjdx9NdsYaHlpJM9jQ5ubA2hyybITs7+/j1JZ32kN4sv2ZkjDkf2CYiP4r/pZ4rIt3LeF8e8BTwNRH56SLb3AvcADQtFnqLhNHHgevmhdFuEfnSsn6oFdAjI5UM4XCYvr4+hoeHE7J/Ew1R0nOEqra7yfa58ec34G7Zx+iWy8AkviV7oYaHq+INDxW5dgDy8vIoLS3VaQ9pyuprRl8htux4k4hsM8ZUAz8RkYuXeF8m8BCxad/fXGSbS4BvE7umNHmSo5p6LDpNtxoaRiqZpqen6enpYXp6OjEfEI1Q0v84rta7cUx1EcitwdN8J96aayDDnpjPnENEeHMwNuFhtuHhnKps9jQ7OaMyiwxjyMjIoLCwUNvE04zVYfQ6cBbwqoicFX/uDRE5/STvMcAPgRER+eNFtjkLuIdYe/Ux4ADQISJfWWDbek4MIzuxBoargD5ip9v2isjbS/5QK6RhpFKB1+ulr6/P8lbw4yRKkfs3VLUewDnRRtBRyUDzHQzXXo/Y1mdF2GFfhEfmNDy48uIND/UO8rJipxAzMzOPX1/SNvHUZnUYvSAi5xljXhWRs40xTmLLR5wsjC4GngHeJNbaDfCfReTXc7a5CJgQkTfjjzOBT4vIv8zb1z3A5UAZMAD8tYh8N/7aDcC3iHXQfU9EvrbkD7QKGkYqVUQiEdxuN4ODg5a3gh8nQsHgC1S1HiBv9B1msksZaPo4w3UfIWp3JOYz5wlFhed7YxMe3vWGyLLBxbUOrm920lj8fju4w+HQNvEUZnUY/QWwFdgD/A/gs8SaGr611kLThYaRSjWBQIDu7m4mJ1fc3Lp8IuR7X8PVepCC4dfik8JvY7DhFqIJmhS+kGNjsVVpn+kOEIwI20oyua7ZyYU1sYYHAGMM+fn5lJaWUlRUpG3iKSIRDQzXE7vXyBC7BvTw2kpMLxpGKlWNjo7S29trydLnJ5M78jZVrQcoHHyBsD2XoYZbGGi4jUh2YUI/d67pmShPdPk51ObDPRWhIMtwVaOTaxudVOS+39yQkZFx/DRefn7+utWnTmRVa/cREbnW0srSlIaRSmXRaBS3252QKQ7zOcaOUtV2N0XuZ+KTwj+Cp+l2wjmlCf3cuaLxhoeH23y80h8E4JwtsQkPp8cbHmZlZWUdbxPXaeLrz6owem22YWGz0zBS6SAQCNDT08PExMTSG69RzmQnrta7Kel7HMmwMbz1BgaaPpHQSeELGfJFONLu49FjfiaCUbbEGx6uqHeQm/XBU3XOOWOI7PbEdwkq68KoA/jTxV5f7L6hjUjDSKWTsbExenp6En7qDiBrui+2fEXPEQyCt+ZaPM13EsyrSfhnzxWKCL/tDXCo3cd73hDZNsMlW3PY0+ykoeiDjQ3GmONjiIqKirRNPIGsCiMv8Ati14nmExH5g9WXmF40jFS6mV3Qz+PxJPzUHUCmfxBX+08o6/o3TDTMaPXluJv3EShoSPhnz9cxGlvS4pluPzMROKU0k+ubnZxfk0Nmxgf/OtMxRIllVRi9KiJnW1pZmtIwUukqGAzS3d29LqfuAOzBESrbH6C88xfYIn5GXRfhadmPr+iUdfn8uaZmojwRn/DgmYpQmJ3B1Y0Orm10UuY8cZqDjiGynl4zspiGkUp363nqDsA2M0HFsZ9Sceyn2ENTjJfvxr1tP9Mlpy79ZotFRXhjIDbhYbbh4dwtsQkPp1dkLXiaLi8vj5KSEl2tdo2sCqNTN+vKrvNpGKmNYD277mZlhKYp75qdFD7GZOkZuFv2M1l2dsInhS9kcDq2pMVjHT4mZiTW8NDs5Iq6ExseIHZ9qaioiJKSEgoLC/X60gpZFUYPiciNS3zQkttsBBpGaiNZlxtm5zHhQGxSePu9ZAWGmS7ajrtlP+OVFyQllGYiwnPxCQ+tI+83PFzf7KS+aOFJDna7/fj1pdzc3HWuOD1ZFUZjxJZ+WHQTYJeINK68xPSiYaQ2otHRUXp6ehI3624BJjJDac9hXO33ku1z4ytoxNOyn9GqS9ZlUvhC2kdjEx5+0+1nJgrbSzPZs0jDw6ycnPeXUc/KWp+5fenIqjC6bBmfNSMiv11JcelIw0htVOsy624h0QglfY/harsbx1Q3gdxa3C37GKm+cl0mhS9kMt7wcLjNh2c6QtGchofSBRoeZuXn51NSUqLLXCzA8nFAm52Gkdro/H4/3d3dTE2dsLZlYklkzqTwdoIOV2z5itrr1m1S+HxREX43EJvw8Ko7iDHwofiEh9MWaXiA2Bii2etLusxFjIaRxTSM1Gbh9Xrp7e0lHA6v7weLUDjwW6paD5A79i4zOaUMNH2Coa03IvbkjfEZnA5zuN3PY8d8TM4I1fk29jQ5uazeQW7m4sNY7Xb78dN4m3mZCw0ji2kYqc0kEonQ19fH0NCCCy8nlgj5w69SdfTH5I+8QSiriMHG2xisv5loZvKaBmYiwnM9sQkPrSMhcmyGS+tiEx7qCk++dIXD4TgeTJttmQurl5CoEJHBec+dIiLvraHGtKJhpDaj6elpuru78fl8Sfn8XO+bVLUepHDoRcKZeQw23Mpgw0eJZK3fpPCFtI3EJjw8G2942FkWa3jYXb14w8Os2WUuiouLN8UyF1aH0XvAX4nIffHHfwJ8VkR2rrnSNKFhpDYrEWFoaIj+/n4ikUhSanCOvYer9SDFnt8QseUwVH8zA40fJ5xTkpR6Zk0Gozze6edwu4+B6QhFORlc3eDg2iYnpY6TNzLMXl+aXUZ9o7I6jKqA7wABoBL4PfAnIrLOVzqTR8NIbXahUIje3l5GRkaSVkPOxDGq2g5S3PckkmFneOsNeJrvIOSoSFpNEGt4eM0zw6E2H695Yg0Pu+MTHk4tX7zhYVZmZubx03gOx/qspLteErG43h8BXya2hPidIvLs2kpMLxpGSsVMTEzQ3d1NMBhMWg3ZU72xSeG9RxAM3trYpPCZ3Oqk1TTLMxXmcLuPxzv9TM0INfmxCQ+X1zlwnqThYdZGW0bd6iOjRwA38B+AGuB7wNMisujyEhuNhpFS7xMRPB4Pbrd7fe9NmifTN4Cr/V7Kun+NiUYYqb4ST8teAvn1SatpVjAiPNsTW5W2fTRMjt1wWV0Oe5qcbF2i4QE2zjLqVofRLSLy8zmP7cCXReSrayszfWgYKXWi9Z4Ivhh7wEtlx/2Ud/6SjEiQMdfFuFv24S/altS6ZrWOxE7hPdsTIBSFXeWZ7GnKZXd1NvYlGh4gtszF7PWldFtGXVu7LaZhpNTikjFWaCG24DiVxx6k/NjPsIenGa84D3fLvqRMCl/IRDDK48d8HG73M+iLUJKTwTWNTq5udFCyRMPDrHRbRt3qI6NJYHbDLCATmBKR5PZXriMNI6VOLhKJ0N/fz+Dg4NIbJ1hGaIqKYz+nsuMB7KEJJkrPxLNtP5OlZyVlKOt8ERFecwc51O7jNc8MNgO7q2NDWneWZS57ckM6LKOe0CMjY8wtwG4R+c+rKS4daRgptTw+n4+urq6k3Zs0V0bYT1nXQ1S230dW0MtU8U7cLfuZqDgvJUIJwD0V5ki7j8eP+ZkKCbUF9tiEh7ocHMtoeIAPLqNeWFiYUteXEn6azhjzvIicv6o3pyENI6WWLxXuTZorNin8YVxt95DtH8RX0Iy7ZT9jVReDSY2/uINh4TfxhoeOsTAOu+GyOgd7mp3UFiz/qCfVllG3+jTdrXMeZgDnApeJyAWrLzG9aBgptXKhUIienh5GR0eTXUpMNExp7yO42u4hZ7oXf14dnpa9jGy5EjJSY9q2iNA6O+GhJ0A4CqeWZ7Gn2cmHtiyv4WFWVlYWpaWlSV1G3eow+v6ch2GgE/iX+SOCNjINI6VWb3x8nJ6enqTem/QBEqG4/ylcrXfjnOwg6KzC3byXkdprkYzUubdn/HjDg48hXzTW8NDk5JoGB8XLbHiYlZube3wM0XpeX9JuOotpGCm1NslY8nxJEo1NCj96gNzx95jJKcfT/AmGt34YsSXnSGIhERFedQc51Obj9YFYw8P5NTlc17SyhgeIXV8qLCw8fn0p0ctcWLW43v/m/S66E4jIf1hdeelHw0gpayRt3aSTEaFg6GVcrQfIH3mTUFYxA00fZ6j+JqL21Fr+oX8yzJGOWMPDdEjYWhhreLi0LgeHfWXXv9ZjGXWrwuhTJ3uziPxwFbWlJQ0jpaw1PDxMX1/f+q+btIQ87++oOnqAguFXCGfmxyeF30okK7VuNp1teHi4zcexsTBOu+HyegfXNTmpWUHDw6zs7OzjbeJWXl+yKozsIpJa/09JEg0jpawXDofp7e3F6/Umu5QTOEd/T1XrQYoGniNidzJYfzODjbcRzi5OdmkfICIcHQlxqM3Hcz0BwgKnVWSxpynW8GBbQcPDrLy8vOPXl9a6jLpVYfSqiJwd//P/FpEvramqNKZhpFTiTE5O0tXVlToNDnM4Jtpxtd5Ncf+TSEYWQ3UfZqDpE4Qc5cku7QTjgQiPHvNzpN3HsD9KqSODaxudXNXooDhn5aGSkZFBVVUVLpdr1TVZFUavichZ8T8fD6bNSMNIqcQSEdxuNx6PJ3UaHObInurG1XoPpX2PIMaGt3YPnuY7mHFWJbu0E0SiwivxCQ+/G5jBbuC8mtiEh+2lK2t4KC4uprGxcdW1JOLISMNIw0iphAsEAnR3dzM5OZnsUhaU5fNQ2XYvZT0PYyTCSPXVuJv3EszfmuzSFtQ3GVvS4olOP76QUBdveLhkmQ0PqRJGPqANMEBT/M/EH4uInL7qCtOMhpFS68vr9dLb25tyDQ6zMgPDVLbfR3nnrzDRGUarLsXTsh9/YVOyS1tQIBzlme4Ah9p8dI6/3/Cwp9lJdf7iDQ+pEkZ1J3uziHStora0pGGk1PpL5QaHWfbgGBUdD1DR+XNsYR9jlefjbrkLX/GOZJe2IBHhPW9swsNv5zQ8XN/s5NyqExseUiKM1Ps0jJRKnlRucJhlm5mkvPPnVHY8GJsUXnYO7pb9TJWenjJDWecbm9Pw4PVHKXPEJjxc3eCgKN7woGGUYjSMlEquVG9wmJUR9lPe9Usq2+8jMzjKVPGpuLftZ6L8QykbSpGo8HJ8wsMbg7GGhwtqY6vSnt9cQVPT6k89ahhZTMNIqdQQCATo6upKrQkOCzCRIGXdv8bV9hOyAoNMF27D07KPMddFKTMpfCG9E7GGhyc7/fjCQnNpNv/l5jO4dNvqWtktCyNjzOki8oYx5jQReXNV1WwAGkZKpZbh4WF6e3tTYomKkzHRECW9j+BqvZscXz/+/HrcLfsY3XI5mNSYFL4QfzjKM10BHukM8lc3nc4V2ytWtR8rw+hbwD8BfyQif7yqajYADSOlUk84HKanp4eRkZFkl7K0aITi/iepaj2AY6qLQG41nua9jNRcnVKTwucrKiqisbFx1QNVlxtGJz1WNMb8dXyb54EMY8x/WVU1SimVAHa7nYaGBlpaWpK2Xs+yZdgYrbmKdy7/Lu3n/jci9lzqf/e37Hr8LsqP/RwTmUl2hQsyxiR8sjcs78joJmAPcEhEfpnwilKUHhkpldpScomKkxGhYOglqo7+mLzRtwlll+Bpup3huo8QtTuSXd1x69VNt5yraOeJyL8HPrTqapRSKsEyMjKorq5mx44dCVsOwVLGMFGxm/cu+gfeu+Cb+PPrqX3nnznt0TtxHf0xtlBqN2hYTbvplkmPjJRKL4ODg/T396d8g8NcuaPv4Go9QNHA80TsuQzW38Jg48cIZxclrab1OjJav7VnlVJqHVVUVFBUVERPTw9jY2PJLmdZpot30r77b3CMt1LVehBX291UHHuQ4bobGWi6nVBOWbJLTBgNI6XUhpWVlUVTUxNjY2N0d3cTCoWSXdKy+Atb6Dj3v5Iz2RUPpJ9S3vkLhmuvZ6D5Dmacq1/SIVVpGCmlNryioiLy8/Pp6+tjaGgo2eUsWyC/js6zvkz/tk/haruXsu5fU979b3hrrsbTvJdgXm2yS7TMUq3dNmPMvzPGfNUYc9G8176S2NKUUso6NpuNrVu3sn37dhyO1OlWW46Z3C10n/H/8tZVBxmsv5mSvifZ9cRnaHjlq+RMdCS7PEss1U33f4HLAC/wD8aYb8557daEVaWUUgmSm5vLjh072LJlCxkZqTuWZyEhRzm9p36RN68+yEDz7RQOPM+upz5H04t/hXPs3WSXtyZLLSHxxuyaRcYYO/B/gDLgTuD52VVgNwPtplNq4wkGg3R1daXsQn5Lsc1MUHHsZ1Qc+yn20CTj5R/C07IvNincIqlyn1HW7B9EJCwifwi8DjwO5K26OqWUSgHZ2dls27aNuro67Pb0u4QeySrAfcqnePPqe+jd8Yc4x9s45bk/Ztuzf0z+4EuQRrfuLBVGLxtj9sx9QkT+O/B9oD5RRSml1HoqKytj165dlJSUJLuUVYnanQw038GbVx2ke9cXyfb1s+2Fv2D7b/49hZ5nQaLJLnFJetPrMulpOqU2h/Hxcbq7u5mZSc1ZccthIjOU9h7B1XYP2T43vvxGPC37GN1y6YonhafEaTpjzJ/P+fPH5732N6uuTimlUlRhYSG7du2isrJyXQaEJoLYshiuu5G3rvgRx876MkbCNL76VXY98RlKew5BNJzsEk+w1Gm6O+b8+cvzXtuDUkptQBkZGdTU1LB9+3acTmeyy1m9DBsjNdfwzuXfo/2cvyZqy6b+9f+PUx//JGWdv0ypSeFLhZFZ5M8LPVZKqQ3F6XSyfft2ampq0q4N/ANMBmNbLuP3l36Htt1fI5RTQt2b3+LUx/dR0X4/GWF/sitccgKDLPLnhR4rpdSGY4yhsrKSoqIiuru7mZiYSHZJq2cM45UXMF5xPvne13AdPUDtO9/G1XY3g423MVh/M9HM5DRKLxVGZxhjJogdBTnifyb+OCehlSmlVArJzs6mpaUFr9dLb28v4XDqXXdZNmOYLDubybKzyR15i6rWg1S/+10q2+5lsOFWBhtvJZJVuK4lnTSMRCR1F2hXSqkkKC0tpbCwMH2WO1/CdMmptJ33dRxjR6lqPciW1h9T2XE/Q/U3MdB4O1C8LnVoa/cyaWu3Umq+iYkJurq60roNfL6cyWO4Wu+mpO8JJMPGZMutFF7/V1C0uqGsVq70qpRSagEFBQVp3wY+XyC/gc6z/5K3r/whIzXXUND6ILQeSfjn6pHRMumRkVLqZKanp+nq6sLvT35nmpXKMgPU7TgHMlfXJqBHRkoptY5mp4FXV1endxv4PJG8qlUH0UpsnH9jSimVZMYYXC4XO3fuJD8/P9nlpBUNI6WUstjcaeA2mzYlL4eGkVJKJcjsNPDi4vVpj05nGkZKKZVAmZmZNDY20tTURGZmZrLLSVnpt5qUUkqloaKiIvLz8+nr62NoaCjZ5aQcPTJSSql1YrPZ2Lp1K6eccgo5OTpRbS4NI6WUWmd5eXns3LmTqqqqDXOz7FppGCmlVBIYY9iyZQs7duwgNzc32eUknYaRUkolkcPhYPv27dTW1m6om2VXavP+5EoplUIqKirYtWsXBQUFyS4lKTSMlFIqRWRlZdHS0kJ9fT12++ZqdtYwUkqpFFNaWrrpbpbVMFJKqRRkt9tpbGykubmZrKysZJeTcBpGSimVwgoLC9m5cyfl5eXJLiWhNIyUUirFbYabZTWMlFIqTeTl5bFjxw5cLteGu1lWw0gppdJIRkYG1dXVbN++HafTmexyLKNhpJRSacjpdLJ9+3Zqamo2xM2y6f8TKKXUJmWMobKyckOsLKthpJRSaW4jrCyrYaSUUhvE7MqyRUVFyS5lxTSMlFJqA8nMzKSpqYnGxsa0Wll2cw0/UkqpTaK4uJj8/Hx6e3vxer3JLmdJemSklFIblN1up76+npaWlpQfKaRhpJRSG1xBQQG7du2ioqIi2aUsSsNIKaU2gYyMDGpra1N2pJCGkVJKbSJ5eXns3Lkz5UYKaRgppdQmY4yhurqaHTt2pMxIIQ0jpZTapBwOB9u3b6e6ujrpI4U0jJRSahMzxuByudixYwd5eXlJq0PDSCmlFDk5OZxyyils3bo1KSOFNIyUUkodV15ezs6dOyksLFzXz9UJDEoppT4gKyuL5uZmvF4v09PT6/KZGkZKKaUWVFpaSmlp6bp8lp6mU0oplXQaRkoppZJOw0gppVTSaRgppZRKOg0jpZRSSadhpJRSKuk0jJRSSiWdhpFSSqmk0zBSSimVdEZEkl1DWjDGjAOt854uBMYX2Hz+82XAcIJKW8piNSZ6P8vdfqntTvb6cv/9L/Zcsr6XZH0nK3nPar+XtT6vvyur3y5Vf1fqRKR8ya1ERP9Zxj/Ad5bz3ELPAy+nUt3rsZ/lbr/Udid7fbn//k/yXFK+l2R9J+vxvaz1ef1dsf47Wen3kqzfFT1Nt3y/WuZzJ3s+GayqZaX7We72S213stdX8u9fv5OVvWe134tVzyeD/q4s73MSQk/TrQNjzMsicm6y61AfpN9L6tHvJDWtx/eiR0br4zvJLkAtSL+X1KPfSWpK+PeiR0ZKKaWSTo+MlFJKJZ2GkVJKqaTTMFJKKZV0GkZJZoxpNMZ81xjzQLJr2cyMMbnGmB8aY/7FGLMv2fWoGP39SD3GmFvivye/MMZca9V+NYzWwBjzPWPMoDHmrXnP7zHGvGeMaTPG/KeT7UNEOkTks4mtdHNa4fdzK/CAiHweuGndi91EVvK96O/H+ljhd/Lz+O/Jp4FPWFWDhtHa/ADYM/cJY4wN+CfgemAncKcxZqcx5jRjzEPz/qlY/5I3lR+wzO8HqAF64ptF1rHGzegHLP97UevjB6z8O/lK/HVL2K3a0WYkIk8bY+rnPb0baBORDgBjzL3AzSLydeDG9a1wc1vJ9wP0Eguk19H/SEuoFX4v76xvdZvTSr4TY8zvgW8AD4vIq1bVoL901qvm/f/ChthfctWLbWyMKTXG/DNwljHmy4kuTi36/fwU+Jgx5tuk1oiazWLB70V/P5Jqsd+VLwFXA7cZY/4fqz5Mj4ysZxZ4btE7i0XEC1j2haolLfj9iMg08Jn1LkYdt9j3or8fybPYd/IPwD9Y/WF6ZGS9XqB2zuMaoD9JtagT6feTmvR7ST3r+p1oGFnvJaDFGNNgjMkC7gB+meSa1Pv0+0lN+r2knnX9TjSM1sAYcw/wW+AUY0yvMeazIhIGvggcBn4P3Ccibyezzs1Kv5/UpN9L6kmF70QHpSqllEo6PTJSSimVdBpGSimlkk7DSCmlVNJpGCmllEo6DSOllFJJp2GklFIq6TSMlLKAMSZijHl9zj8nXTpkvRhjOo0xbxpjzjXG/CxeW5sxZnxOrRcu8t7PGWN+PO+5yvhSA5nGmJ8YY0aMMbesz0+jNjK9z0gpCxhjpkQkz+J92uM3Hq5lH53AuSIyPOe5y4E/FZGTTpE3xhQDrUCNiATiz30ROE1E/l388QFi60D9fC11KqVHRkolUPzI5L8ZY16NH6Fsjz+fG1/Q7CVjzGvGmJvjz3/aGHO/MeZXwBFjTIYx5v8YY96Or4H1a2PMbcaYq4wxP5vzOdcYY366hjo/ZIx5yhjzijHmYWNMpYiMAs8BH56z6R3APav9HKUWo2GklDUc807TzV0Bc1hEzga+Dfxp/Lm/BB4XkQ8BVwB/a4zJjb92AfApEbmS2Aq09cBpwOfirwE8DuwwxpTHH38G+P5qCjfGZAN/D3xMRM4BDgBfjb98D7EAwhhTG6/l6dV8jlIno0tIKGUNv4icuchrs0csrxALF4BrgZuMMbPhlANsjf/5EREZif/5YuB+EYkCHmPMExCb4x+/nrPfGPN9YiH1yVXWvgPYBTxqjAGwEZvYDLHBmP9gjMkjtsT0ffFalLKUhpFSiReM/2+E93/nDLEjkffmbmiMOQ+YnvvUSfb7fWILAQaIBdZqry8Z4A0RuWT+CyIybYx5lNiqq3cAX1jlZyh1UnqaTqnkOAx8ycQPRYwxZy2y3W+IrUCbYYypBC6ffUFE+omtL/MV4AdrqOUdYquq7o7XkmWM2TXn9XuAPwOKROSlNXyOUovSMFLKGvOvGX1jie2/CmQCbxhj3uL9azTzPUjslNlbwP8FXgDG57x+EOgRkXdWW7iIBIHbgG8aY34HvAacN2eTQ8ROId672s9Qaina2q1UijPG5InIlDGmFHgRuEhEPPHX/hF4TUS+u8h7O5nX2m1xbdrarSyhR0ZKpb6HjDGvA88AX50TRK8ApxPrflvMEPCYMeZcq4syxvwEuIjYNSul1kSPjJRSSiWdHhkppZRKOg0jpZRSSadhpJRSKuk0jJRSSiWdhpFSSqmk0zBSSimVdP8/8T/3RKQT0zAAAAAASUVORK5CYII=\n", "text/plain": [ "
" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "energy_range = [0.1, 100] * u.TeV\n", "model.plot(energy_range=energy_range, energy_power=2)\n", "result.model.plot(energy_range=energy_range, energy_power=2)\n", "result.model.plot_error(energy_range=energy_range, energy_power=2);" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Exercises\n", "\n", "* Change the observation time to something longer or shorter. Does the observation and spectrum results change as you expected?\n", "* Change the spectral model, e.g. add a cutoff at 5 TeV, or put a steep-spectrum source with spectral index of 4.0" ] }, { "cell_type": "code", "execution_count": 13, "metadata": {}, "outputs": [], "source": [ "# Start the exercises here!" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## What next?\n", "\n", "In this tutorial we simulated and analysed the spectrum of source using CTA prod 2 IRFs.\n", "\n", "If you'd like to go further, please see the other tutorial notebooks." ] } ], "metadata": { "kernelspec": { "display_name": "Python 3", "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.6.0" }, "nbsphinx": { "orphan": true } }, "nbformat": 4, "nbformat_minor": 2 }