{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\n",
    "<script type=\"text/javascript\" src=\"../_static/linksdl.js\"></script>\n",
    "<div class='alert alert-info'>\n",
    "**This is a fixed-text formatted version of a Jupyter notebook.**\n",
    "\n",
    " You can contribute with your own notebooks in this\n",
    " [GitHub repository](https://github.com/gammapy/gammapy-extra/tree/master/notebooks).\n",
    "\n",
    "**Source files:**\n",
    "[spectrum_models.ipynb](../_static/notebooks/spectrum_models.ipynb) |\n",
    "[spectrum_models.py](../_static/notebooks/spectrum_models.py)\n",
    "</div>\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Spectral models in Gammapy"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "## Introduction\n",
    "\n",
    "This notebook explains how to use the functions and classes in [gammapy.spectrum.models](http://docs.gammapy.org/dev/spectrum/#module-gammapy.spectrum.models) in order to work with spectral models.\n",
    "\n",
    "The following clases will be used:\n",
    "\n",
    "* [gammapy.spectrum.models.PowerLaw](http://docs.gammapy.org/dev/api/gammapy.spectrum.models.PowerLaw.html)\n",
    "* [gammapy.utils.modelling.Parameter](http://docs.gammapy.org/dev/api/gammapy.utils.modeling.Parameter.html)\n",
    "* [gammapy.utils.modelling.ParameterList](http://docs.gammapy.org/dev/api/gammapy.utils.modeling.ParameterList.html)\n",
    "* [gammapy.spectrum.models.SpectralModel](http://docs.gammapy.org/dev/api/gammapy.spectrum.models.SpectralModel.html)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Setup\n",
    "\n",
    "Same procedure as in every script ..."
   ]
  },
  {
   "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.spectrum import models\n",
    "from gammapy.utils.modeling import Parameter, ParameterList"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Create a model\n",
    "\n",
    "To create a spectral model, instantiate an object of the spectral model class you're interested in."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "PowerLaw\n",
      "\n",
      "Parameters: \n",
      "\n",
      "\t   name     value   error       unit      min max frozen\n",
      "\t--------- --------- ----- --------------- --- --- ------\n",
      "\t    index 2.000e+00   nan                 nan nan  False\n",
      "\tamplitude 1.000e-12   nan 1 / (cm2 s TeV) nan nan  False\n",
      "\treference 1.000e+00   nan             TeV nan nan   True\n"
     ]
    }
   ],
   "source": [
    "pwl = models.PowerLaw()\n",
    "print(pwl)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "This will use default values for the model parameters, which is rarely what you want.\n",
    "\n",
    "Usually you will want to specify the parameters on object creation.\n",
    "One way to do this is to pass `astropy.utils.Quantity` objects like this:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "PowerLaw\n",
      "\n",
      "Parameters: \n",
      "\n",
      "\t   name     value   error       unit      min max frozen\n",
      "\t--------- --------- ----- --------------- --- --- ------\n",
      "\t    index 2.300e+00   nan                 nan nan  False\n",
      "\tamplitude 1.000e-12   nan 1 / (cm2 s TeV) nan nan  False\n",
      "\treference 1.000e+00   nan             TeV nan nan   True\n"
     ]
    }
   ],
   "source": [
    "pwl = models.PowerLaw(\n",
    "    index=2.3 * u.Unit(''),\n",
    "    amplitude=1e-12 * u.Unit('cm-2 s-1 TeV-1'),\n",
    "    reference=1 * u.TeV,\n",
    ")\n",
    "print(pwl)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "As you see, some of the parameters have default ``min`` and ``values`` as well as a ``frozen`` flag. This is only relevant in the context of spectral fitting and thus covered in [spectrum_analysis.ipynb](https://github.com/gammapy/gammapy-extra/blob/master/notebooks/spectrum_analysis.ipynb). Also, the parameter errors are not set. This will be covered later in this tutorial."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Get and set model parameters\n",
    "\n",
    "The model parameters are stored as ``ParameterList`` on the spectal model. Each model parameter is ``Parameter`` instance. It has a ``value`` and a ``unit`` attribute, as well as a ``quantity`` property for convenience."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "ParameterList\n",
      "Parameter(name='index', value=2.3, unit='', min=nan, max=nan, frozen=False)\n",
      "Parameter(name='amplitude', value=1e-12, unit='1 / (cm2 s TeV)', min=nan, max=nan, frozen=False)\n",
      "Parameter(name='reference', value=1.0, unit='TeV', min=nan, max=nan, frozen=True)\n",
      "\n",
      "Covariance: \n",
      "None\n",
      "Parameter(name='index', value=2.3, unit='', min=nan, max=nan, frozen=False)\n",
      "Parameter(name='index', value=2.6, unit='', min=nan, max=nan, frozen=False)\n",
      "Parameter(name='amplitude', value=1e-12, unit='1 / (cm2 s TeV)', min=nan, max=nan, frozen=False)\n",
      "Parameter(name='amplitude', value=2e-12, unit='1 / (m2 s TeV)', min=nan, max=nan, frozen=False)\n"
     ]
    }
   ],
   "source": [
    "print(pwl.parameters)\n",
    "\n",
    "print(pwl.parameters['index'])\n",
    "pwl.parameters['index'].value=2.6\n",
    "print(pwl.parameters['index'])\n",
    "\n",
    "print(pwl.parameters['amplitude'])\n",
    "pwl.parameters['amplitude'].quantity = 2e-12 * u.Unit('m-2 TeV-1 s-1')\n",
    "print(pwl.parameters['amplitude'])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## List available models\n",
    "\n",
    "All spectral models in gammapy are subclasses of ``SpectralModel``. The list of available models is shown below."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "[gammapy.spectrum.models.ConstantModel,\n",
       " gammapy.spectrum.models.CompoundSpectralModel,\n",
       " gammapy.spectrum.models.PowerLaw,\n",
       " gammapy.spectrum.models.PowerLaw2,\n",
       " gammapy.spectrum.models.ExponentialCutoffPowerLaw,\n",
       " gammapy.spectrum.models.ExponentialCutoffPowerLaw3FGL,\n",
       " gammapy.spectrum.models.PLSuperExpCutoff3FGL,\n",
       " gammapy.spectrum.models.LogParabola,\n",
       " gammapy.spectrum.models.TableModel,\n",
       " gammapy.spectrum.models.AbsorbedSpectralModel,\n",
       " gammapy.spectrum.crab.MeyerCrabModel]"
      ]
     },
     "execution_count": 6,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "models.SpectralModel.__subclasses__()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Plotting\n",
    "\n",
    "In order to plot a model you can use the ``plot`` function. It expects an energy range as argument. You can also chose flux and energy units as well as an energy power for the plot"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<matplotlib.axes._subplots.AxesSubplot at 0x10f326240>"
      ]
     },
     "execution_count": 7,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAZIAAAEOCAYAAACjJpHCAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMS4wLCBo\ndHRwOi8vbWF0cGxvdGxpYi5vcmcvpW3flQAAIABJREFUeJzt3Xd8FWXa//HPlYQk1EiX3kKHgEgv\nwUIXBEFX1/azK4og7K79Wd1ldXVLwIIFC1hWRAQUEKSomwDSQTrSkV6kCoR6//7IYZ8sj5BDTpmc\nk+/79crLc+bMmflmd/TKPffMNeacQ0REJLdivA4gIiKRTYVEREQCokIiIiIBUSEREZGAqJCIiEhA\nVEhERCQgKiQiIhIQFRIREQmIComIiAREhURERAIS53WAcChVqpSrWrWq1zFERCLKokWL9jnnSue0\nXlQXEjPrAfRITk5m4cKFXscREYkoZrbFn/Wi+tSWc26ic+6BpKQkr6OIiEStqC4kIiISeiokIiIS\nEBUSEREJiAqJiIgERIVEREQCokJyEc45xi7axsnTZ72OIiKSZ0V1ITGzHmY2/NChQ7n6/rxN+/nd\nmKX0eG0WP2w9GOR0IiLRIaoLSaD3kbSsXpJ372zKoeOn6P3GbF74ahXHT54JckoRkcgW1YUkGDrU\nK8u0Qanc0rwy78zcROehGXy/YZ/XsURE8gwVEj8USyzAizc0ZNT9LTGDW9+Zx1PjlnE485TX0URE\nPKdCcgla1SjJ1wNSeSC1OqMXbKVjWjozVu32OpaIiKdUSC5RwfhYnu5Wl/EPt6F4oXju+3Ah/Uct\n4edfTngdTUTEEyokudSo0mVM6NeWQR1rMWXFTjqkpfPlD9txznkdTUQkrFRIAhAfF0P/a2vyVf92\nVClZmAGf/sC9Hyxk56HjXkcTEQkbFZIgqFW2KGP7tuZ/utdjzoaf6ZiWwb/mbeHsWY1ORCT6qZAE\nSWyMcW/bakx9LJVGlZJ4ZvwKbn13Lpv3HfU6mohISEV1IQn0zvbcqFyyEB/f24KX+zRk5Y7DdB6a\nwfCMDZw+ozYrIhKdLD9MDjdt2tR58ajd3YczefaLFUxftZuUikm83CeFuuWKhT2HiEhumNki51zT\nnNaL6hGJ18oWS2T4HVfy+q1XsP3AcXq8Nou06Ws5cVptVkQkeqiQhJiZ0T2lPDMGtef6RuV59Zt1\ndH91Fot/OuB1NBGRoFAhCZPiheNJu7kxI+5uxtETp+nz5vf8eeIqjp087XU0EZGAqJCE2dW1yzB1\nYCq3t6jC+7OzmkDOXq8mkCISuVRIPFA0sQCDezVg9AMtiYuJ4bZ35/HE58s4dFxNIEUk8qiQeKhF\n9ZJMGdCOh9rX4PPF2+iYls7Ulbu8jiUicklUSDyWWCCWJ7vW4YuH21CySAIPfrSIRz5ZzN4jagIp\nIpFBhSSPaFgxiQn92vD7TrWYvnI3HYekM27xNjWBFJE8T4UkDykQG0O/a2oyeUBbqpcqzKDPlnL3\nyAVsP6gmkCKSd0V1IfGiRUowJJcpypiHWvN8j3rM37SfTmnpfDRns5pAikiepBYpedzW/cd4evxy\nZq7bR/OqJXipT0Oqly7idSwRyQfUIiVKVCpRiA/vac7fbkxhza7DdHllJm/+W00gRSTvUCGJAGbG\nb5pWYsag9lxduzQvf72GXm/MZuWOyDplJyLRSYUkgpQplsjbdzTlzduasOvQCa5/fTZ/n7qGzFNq\nAiki3lEhiUBdG5ZjxqBUejWuwLDvNnDdqzNZtGW/17FEJJ9SIYlQlxWK55+/acQH9zQn89RZbnxr\nDs9PWMnRE2oCKSLhpUIS4drXKs3Uganc2bIKH8zZTKchGWSs3et1LBHJR1RIokCRhDj+1LMBYx5s\nRUKBGO58fz6/H7OUQ8fUBFJEQk+FJIo0rVqCyf3b8cjVNRi/ZDsdhqTz9YqdXscSkSinQhJlEgvE\n8ofOdZjQrw1liibw0MeL6fvxIvYcyfQ6mohEKRWSKFW/fBJfPNKGx7vU5ps1e+iYlsGYhVvVBFJE\ngk6FJIoViI3h4auSmTKgHbXKFuEPny/jzvfns3X/Ma+jiUgUiepCEqlNG4OtRukijH6gFX/uWZ/F\nWw7QeWgGI2dvUhNIEQkKNW3MZ7YdOMbT41eQsXYvTasU56U+KSSXURNIEfm/1LRRflXF4oX44O5m\n/POmRqzb8wvdXpnJsO/Wc0pNIEUkl1RI8iEzo8+VFZkxqD0d65Xl71N/pOfrs1mxPX+fAhSR3FEh\nycdKF01g2G1NeOv2K9n7ywl6DpvNy1+rCaSIXJq4i31oZil+bOOUc251kPKIB7o0uJxW1UvywuRV\nvPnvDUxdsYuX+qTQvFoJr6OJSAS46GS7mR0BlgB2kW1Ucs5VDXKuoNJku/9mrtvLU+OWs+3Ace5s\nVYXHu9ShSMJF/94QkSjl72R7Tv+FWOKcS81hRxmXlEzytHY1SzP1sVT+OW0tI77fxIxVu3mhd0Ou\nrl3G62gikkdddI4kpyLi7zoSWQonxPHHHvX4/KHWFEqI4+4RCxg0+gcOHD3pdTQRyYP8mmy3LA3N\nrLOZpZpZyVAHE+9dWaU4X/Vvy6PXJDNh6Q46Dkln0rIdarMiIv/looXEzKqa2RvABmAocDcwCMgw\ns9lmdoeZXWz+RCJcQlwsv+tUmwn92nJ5UiL9PlnCgx8tYs9hNYEUkSw5TbZ/BrwJpDvnzp73WTng\nNmCfc25kKEMGSpPtwXH6zFnenbWJIdPXEh8Xw/9cV4+bmlZEf0uIRCd/J9vVIkUu2ca9v/DkuOXM\n37SftsmlePGGhlQuWcjrWCISZEFtkWJmvc2sqO/1k2b2mZk1DjSkRKbqpYvw6f0t+UuvBvyw9SCd\nh2bw3qxNnFETSJF8yd872593zh0xs9ZAD2A08FboYkleFxNj3N6yCtMGptKiegkGT1rFjW99z7rd\nR7yOJiJh5m8hOdczozvwhnNuLJAQmkgSScpfVpARdzVj6M2N2bzvKNe9OotXv1nHydNqAimSX/hb\nSHaa2TDgZmCymcVfwnc9o+eRhIeZ0euKCkwf1J5O9cuSNn0t178+i2XbDnodTUTCwK/JdjMrAnQD\nljnn1phZeaCRc25KqAMGgybbw2v6qt08+8Vy9h45wf3tqvNYh1oUjI/1OpaIXCJdtZWNCkn4HTp+\nir9OXs2nC7ZStWQhXuqTQsvquo9VJJLowVbiqaSCBXipTwqf3NeCsw5uGT6XZ8Yv50jmKa+jiUiQ\nqZBISLVOLsXXj7Xj3rbVGDX/JzoNyeDbNbu9jiUiQaRCIiFXKD6O/+lej7F9W1MkIY57Ri7ksU+X\nsF9NIEWiQk69tiqY2cdm9p2ZPW5mcdk+Gxv6eBJNrqhcnEn92zLg2ppMWraTDmnpTFiqJpAikS6n\nEcn7wFzgD0A14DszK+77rHoog0l0SoiLZWDHWkzq35ZKxQvSf9QS7v9wIbsOqQmkSKTKqZCUcc69\n7pxb6JzrC7xLVuffaoD+jJRcq3N5McY93IZnutVl1vp9dExLZ9T8nzQ6EYlAORWSBDP7zx3szrkP\nyBqdTAcuD2UwiX6xMcb9qdX5ekAq9SsU46lxy7n1nXls+fmo19FE5BLkVEhGAK2yL3DOfQ3cAvwY\nqlCSv1QtVZhP7mvJX3s3ZMX2Q3QemsG7MzeqCaRIhNANiZKn7DqUybNfLGfG6j00qpjE325sRO3L\ni3odSyRfCuqd7WZWGegHVAX+c+WWc653ABnDRoUksjjnmLhsJ3+asJLDmad45OpkHr4qmfg4Xa0u\nEk7+FpK4nFbwmQB8SNbciNq6SkiZGdc3Kk/b5FL8eeJKhs5Yx+TlO/nbjY1oXOkyr+OJyHn8HZHM\nd841D0OekNCIJLJ9u2Y3z4xfwe7DmdzbthqDOtZWE0iRMAj2qa07gCrAVODEueXOuWWBhAwXFZLI\ndyTzFC9NWcO/5v1E5RKFeKlPQ1rXKOV1LJGoFuxCMhi4D9jI/57acs651IBShokKSfSYu/Fnnhy7\njM0/H+O3zSvzVLc6FEss4HUskagU7ELyI5DinDuR48p5kApJdDl+8gxDZ6zlnZkbKV00gRd6NaRD\nvbJexxKJOsFuI78M0DWYkicUjI/lqW51+eKRNhQvFM99Hy7k0VFL+PmXiPw7RyTi+XvVVklgjZnN\n47/nSCLi8l+JTikVL2NCv7a8lb6B175dx6x1e3muR316Ni6PmXkdTyTf8PfU1rW/ttw5903QEwWR\nmfUAeiQnJ9+/bt06r+NICK3dfYTHP1/GD1sPck2dMrxwQwPKJRX0OpZIRAvFDYl7nHOZvvcFgVLO\nua0BJw0DzZHkD2fOOkZ+v5l/TP2R2BjjqW51+G2zysTEaHQikhvBniMZx3/fiHgW0PNIJE+JjTHu\nbVuNqY+l0qhSEs+MX8Fv35nLpn1qAikSSv4Wkjjn3H8eZ+e7eivhIuuLeKZyyUJ8fG8L/tYnhVU7\nD9NlaAZvp2/g9Bk1ZRAJBX8Lyc9m1u3cGzPrDuwPTSSRwJkZv2lWiRmD2tO+Vmn+OmUNvd/8ntU7\nD3sdTSTq+DtHUhMYRdbVWw7YB9zmnIuIGWzNkeRvzjkmL9/FcxNWcPDYKR6+qgaPXJNMQpzarIhc\nTFAn27Nt9DIA59zBALKFnQqJABw4epLBX61i3OLtJJcpwst9UriySvGcvyiSTwVlst3MbrFsF+Q7\n5w5mLyJmVtXMWgcWVSQ8iheOJ+03jRlxdzOOnTjNjW99z58nruLYydNeRxOJaDndkFgBWGJm84FF\nwF4gEUgGrgIOA0+EMqBIsF1duwzTBrXnb1+v4f3Zm5i+ehcv9U6hTbKaQIrkRo6ntswsDugItAHK\nAceB1cBk59ymkCcMAp3akguZv2k/T4xdxqZ9R7m5aSWevq4uSQXVBFIEQjRHEqlUSORiMk+d4ZVv\n1jE8YyMlC8czuFcDOte/3OtYIp4L9g2JIlErsUAsT3SpwxcPt6FkkQQe/GgRj/xrMXuPqAmkiD9U\nSER8GlZMYkK/Nvyhc22mr9pNxyHpjFu8jfwwahcJRE5XbTULVxCRvKBAbAyPXJ3M5AFtqV6qMIM+\nW8pdIxaw/eBxr6OJ5Fk5jUgeNLMlZvaxmd1uZqXDkkrEY8llijLmodY836MeCzbvp1NaOh/O2czZ\nsxqdiJzP3zvbGwBdgc5kXf77LfA1MNc5l+cbGGmyXQKxdf8xnh6/nJnr9tGsanFe6pNCjdJFvI4l\nEnIhu2rLzAoD15JVWFo455rkLmL4qJBIoJxzjF28ncGTVnH81Bke61CT+9tVp0Csphkleuny32xU\nSCRY9hzJ5LkvVzJlxS7qly/Gy31SaFAhyetYIiGhy39FQqBM0UTevP1K3rytCbsPn6DnsNn8feoa\nMk+d8TqaiGdUSERyoWvDcswYlMoNV1Rg2Hcb6PbqTBZu1pMVJH/K6fLfoWbWPFxhRCLJZYXi+cdN\njfjwnuacOHWWm96ew3NfruCXE2oCKflLTiOSrcAwM9tgZi/4rt4SkWxSa5Vm2sBU/l+rqnw4dwud\nh2SQsXav17FEwuaihcQ590/nXDOgE3AMGGVmK8zsaTOrHpaEIhGgcEIcz19fnzEPtiKhQAx3vj+f\n349ZysFjJ3P+skiEy83lv1cC7wIpzrmIeMScrtqScMo8dYbXvl3HW+kbKV4onr/0qk+XBuW8jiVy\nyYJ61ZaZxZpZVzP7APgK2AjcHGBGkaiUWCCWP3Suw4R+bShbLIGHPl5M348XsedIptfRREIip8n2\nq81sOLAd6E/WHe01nXN9nHOfhyOgSKSqXz6JLx5pw+NdavPNmj10TMtgzMKtagIpUSenEcmfgSVA\nQ+dcV+fcB865I2HIJRIVCsTG8PBVyUwZ0I5aZYvwh8+Xcef789m6/5jX0USCJqfJ9nbOuTedc3vN\nrKWZ3QlgZiXNrHJ4IopEvhqlizD6gVYM7lmfxVsO0HloBiNnb1ITSIkK/s6RPAs8BzzrW5QIfBKq\nUCLRKCbGuKNVVaYNak+zqiV4fuIqbnp7Duv3aJAvkc3fO9tvBLoBRwGcc9uBYqEKFSxm1sPMhh86\ndMjrKCL/UeGygoy8uxlpv2nEhr2/0O2VWbz+7TpOncnzjbRFfpW/heSEy5ohdABmVih0kYLHOTfR\nOfdAUpKa6kneYmb0blKR6QPb07FeWf4xbS3Xvz6bFdv1R49EHn8LyTgzGwYkmdndwDTg/dDFEskf\nShdNYNhtTXjr9ivZ90tWE8iXpqgJpEQWv29INLOuZN3hbsBU59yUUAYLJt2QKJHg0LFTvDh5NaMX\nbqVaqcK83CeF5tVKeB1L8rGgPI/EzKY55zoFNZkHVEgkksxat4+nxi9j6/7j3NGyCk90rUORhDiv\nY0k+FKw72/WMdpEwa1uzFFMfS+WeNtX4eN4WOqWl892Pe7yOJXJBOY1INgK/v9DnzrlxoQgVbBqR\nSKRatOUAT4xdxvo9v9D7igr8T/d6FC8c73UsySf8HZHkNF5OArqTNS9yPgdERCERiVRXVinOV/3b\nMuzb9bzx7w1krNvL89fX57qG5TD7tX8tRcIvpxHJYudckzDmCQmNSCQarN55mCfGLmPZtkN0qleW\nwb0aULZYotexJIoFa45Ef/KI5BF1yxVjXN/WPNW1Dulr99IhLZ3RC35SE0jxXE6F5I6wpBARv8TF\nxvBg+xp8/VgqdcsV44mxy7n9vXn89LOaQIp3ciokL+W0ATObFKQsIuKnaqUK8+n9LflLrwYs3XqI\nzkMzeG/WJs6oCaR4IKc5koNAxsW+D9R3zuXpx+5qjkSi2Y6Dx3lm/HK++3EvV1S+jL/1SaFm2aJe\nx5IoEKwbEtv7sa+Tzrk5lxIu3FRIJNo55/jyhx38aeJKjp44Q79rknmofQ3i4/ztgiTyfwWlkEQL\nFRLJL/b9coLnJ6xk0rKd1Lm8KH+7MYWUipd5HUsiVFCf2S4ikaFUkQRev7UJ79zZlAPHTtJr2Gxe\nnLya4yfVBFJCR4VEJAp1rFeW6YPac3OzSgzP2EjXVzKYu/Fnr2NJlPL3CYllfmVZ7eDHEZFgKZZY\ngL/2TuGT+1pw1sEtw+fyzPjlHMk85XU0iTL+jkhmmtlvzr0xs98B40MTSUSCqXVyVhPI+9pWY9T8\nn+g0JINv1+z2OpZEEX8LyVXAHWY2xswygFpA85ClEpGgKhgfy7Pd6zG2b2uKJsZxz8iFDPh0CfuP\nnvQ6mkQBvwqJc24n8DXQCqgKfOic+yWEuUQkBK6oXJxJj7ZjwLU1mbx8Jx3S0pmwdIfarEhA/J0j\nmQ60ABoA3YAhZvaPUAYTkdCIj4thYMdaTHq0HZVKFKL/qCXc/+FCdh3K9DqaRCh/T20Nc87d6Zw7\n6JxbAbQGDoUwl4iEWO3LizKub2ueva4us9bvo2NaOqPmqwmkXDrdkCgibPn5KE+OXc6cjT/TqnpJ\nXurTkColC3sdSzwW1BsSzeyImR32/WSa2Rkz04hEJEpUKVmYT+5vwYs3NGTF9qwmkO9kbFQTSPGL\nv5PtRZ1zxXw/iUAfYFhoo4lIOJkZt7aozPRB7WmbXIoXJq+m9xuzWbPrsNfRJI/L1Z3tzrkvgGuC\nnEVE8oDLkxJ5586mvHJLY7YeOE6P12YxZPpaTp4+63U0yaNyemY7AGbWO9vbGKApWc9sF5EoZGb0\nbFyBdjVL8+eJK3nlm3VMWbGTl/ukcEXl4l7HkzzG3xFJj2w/nYEjQM9QhRKRvKFE4XiG3nIF79/V\nlCOZp+n95vcMnrSKYydPex1N8hBdtSUifjmSeYqXv17Dx3N/onKJQrzUuyGtk0t5HUtCKFgPtnqN\ni5zCcs71z1288FIhEQmeuRt/5smxy9j88zFuaVaJp6+rS7HEAl7HkhDwt5DkNEei//qKyH9pWb0k\nXz+WypAZa3knYyPf/biHv/RqSMd6Zb2OJh7JaUQS55yL+JOhGpGIhMaybQd5/PNlrNl1hO4p5Xj+\n+vqUKpLgdSwJkmDdkDg/2wZfCziViESVlIqXMaFfWwZ1rMXUlbvomJbO+CXb1GYln8mpkFi2121C\nGUREIlN8XAz9r63JV/3bUbVUYQaOXso9Ixew4+Bxr6NJmORUSPRnhYj4pVbZonz+UGv+2L0eczfu\np9OQDD6au4WzarMS9XKaIzkGrCdrZFLD9xrfe+ecSwl5wiDQHIlIeG3df4ynxi1n1vp9NK9Wgpf7\npFCtlJpARppgXf5b5WJfds5tyUW2sFMhEQk/5xxjFm5j8FerOHn6LAM71uK+ttWIi81VZybxQFAK\nSbRQIRHxzu7DmTz7xQqmr9pNwwpJvNwnhXrli3kdS/wQ1DbyIiK5VbZYIsPvuJJhtzZh56HjXP/6\nLP457UdOnD7jdTQJEhUSEQk5M+O6lHJMH9ie6xuV57Vv13Pdq7NYtOWA19EkCHIsJGaW4vtnw9DH\nEZFoVrxwPGk3N2bE3c04duI0N771PX+auFJNICOcPyOSe8ysJnBvqMOISP5wde0yTBvUnttbVGHE\n7M10GpLBrHX7vI4luXTRQmJmz/nWmQvEmNkfw5JKRKJekYQ4BvdqwGcPtiI+Nobb35vHH8Ys5dCx\nU15Hk0t00ULinPsTMAMYDcxwzv05LKlEJN9oXq0Ekwe0o+9VNRi3ZDsdhqQzdeUur2PJJfDn1FYL\n59zDQLNQhxGR/CmxQCxPdKnDl4+0oXSRBB78aBGP/Gsxe4+c8Dqa+CHP30diZtWBZ4Ak59yNvmUx\nwGCgGLDQOffBxbah+0hEIsepM2cZnrGRV2aso2B8LH/sXo/eTSpgZjl/WYIqT9xHYmbvm9keM1tx\n3vIuZvajma03sycvtg3n3Ebn3PkT/T2BCsApYFtwU4uIlwrExvDI1clMHtCO5DJF+N2Ypdw1YgHb\n1QQyzwr1fSQjgS7ZF5hZLDAM6ArUA35rZvXMrKGZTTrvp8wFtlsbmOOcGwT0DWF+EfFIcpkijHmw\nFc/3qMeCzfvplJbOh3M2qwlkHhTSQuKcywD2n7e4ObDeN9I4CXwK9HTOLXfOdT/vZ88FNr0NOHcn\nk26PFYlSMTHGXW2qMfWxVJpUKc4fv1zJzcPnsGHvL15Hk2xyuvw31sweNLPBZtbmvM+ezeU+KwBb\ns73f5lt2oQwlzewt4Aoze8q3eBzQ2fewrYwLfO8BM1toZgv37t2by6gikhdUKlGID+9pzt9vTOHH\nXUfo+spM3vj3ek6dOet1NCHn7r/vAoXIelLiHUC673QSZrbYOdckxx2YVQUmOeca+N7fBHR2zt3n\ne38H0Nw592hgv8qFabJdJHrsOZLJc1+uZMqKXdQvX4yX+6TQoEKS17GiUrAm25s75251zg0FWgBF\nzGycmSXw309PvBTbgErZ3lcEduRyWyKSz5Qpmsibt1/Jm7c1YffhE/QcNpu/T11D5imd5fZKToUk\n/twL59xp59wDwA/At0CRXO5zAVDTzKqZWTxwCzAhl9sSkXyqa8NyzBiUyg1XVGDYdxvo9upMFm4+\nf0pWwiGnQrLQzP7rqivf3e0jgKo5bdzMRgFzgNpmts3M7nXOnQb6AVOB1cBnzrmVuQkvIvnbZYXi\n+cdNjfjwnuacOHWWm96ew3NfruCXE2oCGU55/obEYNAciUj0O3riNH+f+iMfzNlM+aSC/LV3Q1Jr\nlfY6VkQLyhyJmT2e7fVN5332Yu7jiYgEV+GEOJ6/vj6fP9SKxAIx3Pn+fH4/ZikHj530OlrUy+nU\n1i3ZXj913mddyOPMrIeZDT906JDXUUQkTK6sUoKv+rej39XJjF+ynQ5pGUxZvtPrWFEtp0JiF3j9\na+/zHOfcROfcA0lJujRQJD9JLBDL7zvXZkK/NlyelEDffy3moY8WsedwptfRolJOhcRd4PWvvRcR\nyVPql0/ii4fb8ESXOnz74x46pKXz2cKt5Ie54XDKqZA0MrPDZnYESPG9Pvdej94VkTwvLjaGvlfV\nYMqAdtS+vCiPf76MO9+fz9b9x7yOFjVyerBVrHOumHOuqHMuzvf63PsC4QopIhKoGqWLMPqBVgzu\nWZ/FWw7QeWgGI2dv4oyaQAYs1N1/RUTyjJgY445WVZk2qD3Nqpbg+YmruOmt71m/54jX0SKaComI\n5DsVLivIyLubkfabRmzcd5Rur8zi9W/XqQlkLqmQiEi+ZGb0blKR6QPb07F+Wf4xbS09XpvF8m26\nXeBSqZCISL5WumgCw25twtt3XMn+oyfp9cZsXpqiJpCXIqoLiW5IFBF/da5/OdMHtefGJhV5K30D\nXV+ZybyNP3sdKyKo15aIyHlmr9/Hk+OWsXX/ce5oWYXHu9SmaGL+u1A1WM8jERHJd9okl2LqY6nc\n06YaH8/bQuchGXz344We/C0qJCIiv6JQfBx/7FGPsX1bUzghjrtHLGDQ6B84cFRNIM+nQiIichFN\nKhdnUv+29L+2JhOW7qBDWjqTlu1Qm5VsVEhERHKQEBfLoI61mPhoWyoUL0i/T5bw4EeL2K0mkIAK\niYiI3+qWK8a4vq15ulsd0tfupUNaOqMX/JTvRycqJCIilyAuNoYHUmsw9bFU6pUrxhNjl3P7e/P4\n6ef82wRShUREJBeqlirMqPtb8sINDVi69RCdh2bw7syN+bIJZFQXEt2QKCKhFBNj3NaiCtMHpdKq\nRkn+8tVq+rz5PWt3568mkLohUUQkCJxzTFi6gz9NXMWRzFM8ek1NHmpfg/i4yP17XTckioiEkZnR\ns3EFpg9MpUuDcqRNX8v1r89i6daDXkcLORUSEZEgKlkkgdd+ewXv3NmUA8dOcsMbs3lx8mqOn4ze\nJpAqJCIiIdCxXlmmD2rPLc0rMzxjI11fyWDOhuhsAqlCIiISIsUSC/DiDQ355P4WOOC378zl6fHL\nOZJ5yutoQaVCIiISYq1rlOLrAak8kFqdT+f/RKchGXy7ZrfXsYJGhUREJAwKxsfydLe6jHu4DcUS\nC3DPyIUM+HQJ+6OgCaQKiYhIGDWudBkTH23LwA61mLx8Jx3S0vnyh+0R3WZFhUREJMzi42IY0KEm\nkx5tR6UShRjw6Q/c98FCdh467nW0XFEhERHxSO3LizKub2ueva4uszfso1NaBp/M+4mzEdZmJaoL\niVqkiEheFxtj3NeuOlMfS6UdpruoAAAH5klEQVRBhSSeHr+cW9+dy+Z9R72O5je1SBERySOcc4xe\nsJUXJq/m1Jmz/K5jbe5pW43YGPMkj1qkiIhEGDPjluaVmT6wPW2TS/PC5NX0fmM2a3Yd9jraRamQ\niIjkMZcnJfLOnVfy+q1XsO3Acbq/Oou06Ws5cTpvtllRIRERyYPMjO4p5Zk+qD09GpXn1W/W0eO1\nWSz56YDX0f4PFRIRkTysROF4htzcmBF3NeNI5ml6v/k9gyet4tjJ015H+w8VEhGRCHB1nTJMG5jK\nbS0q896sTXQZOpPv1+/zOhagQiIiEjGKJhbgL70a8ukDLYkxuPXdeTw5dhmHjnvbBFKFREQkwrSs\nXpKvH0vlwfbV+WzhVjqmpTNt5S7P8qiQiIhEoMQCsTzVtS5fPNKGEoXjeeCjRfT7ZDH7fjkR9iwq\nJCIiESylYlYTyN91rMW0lbvpkJbO+CXbwtoEUoVERCTCFYiN4dFra/JV/7ZUK1WYgaOXcs/IBew4\nGJ4mkCokIiJRombZonz+UGue61GPuRv30zEtnY/nbgn5fqO6kKhpo4jkN7Exxt1tqjFtYCpXVC7O\nut1HQr5PNW0UEYlSzjlOnXHEx+VuzOBv08a4XG1dRETyPDMjPi70nYOj+tSWiIiEngqJiIgERIVE\nREQCokIiIiIBUSEREZGAqJCIiEhAVEhERCQg+eKGRDPbCxwELnaLe9JFPi8F5I0nyPjvYr9PXt5X\nINu61O/6u74/6+W0TrQdXxC+Y0zHl3fHVxXnXOkc13LO5YsfYHhuPwcWep0/2L9vXt1XINu61O/6\nu74/6+W34yvY/7+Haz86vkLzk59ObU0M8PNIE87fJ5j7CmRbl/pdf9f3Z738dnxB+H4nHV95/PjK\nF6e2AmVmC50f/WZEckPHl4RSOI6v/DQiCcRwrwNIVNPxJaEU8uNLIxIREQmIRiQiIhIQFRIREQmI\nComIiAREheQSmVkvM3vHzL40s05e55HoYmZ1zewtM/vczPp6nUeik5kVNrNFZtY9GNtTIQHM7H0z\n22NmK85b3sXMfjSz9Wb2JIBz7gvn3P3AXcDNHsSVCHOJx9dq59xDwG8AXRIsfrmUY8znCeCzYO1f\nhSTLSKBL9gVmFgsMA7oC9YDfmlm9bKs86/tcJCcjuYTjy8yuB2YB34Q3pkSwkfh5jJlZB2AVsDtY\nO1chAZxzGcD+8xY3B9Y75zY6504CnwI9LcvLwBTn3OJwZ5XIcynHl2/9Cc651sBt4U0qkeoSj7Gr\ngZbArcD9ZhZwHYgLdANRrAKwNdv7bUAL4FGgA5BkZsnOube8CCcR71ePLzO7CugNJACTPcgl0eNX\njzHnXD8AM7sL2OecOxvojlRILsx+ZZlzzr0KvBruMBJ1LnR8/Rv4d3ijSJT61WPsPy+cGxmsHenU\n1oVtAyple18R2OFRFok+Or4k1MJ2jKmQXNgCoKaZVTOzeOAWYILHmSR66PiSUAvbMaZCApjZKGAO\nUNvMtpnZvc6500A/YCqwGvjMObfSy5wSmXR8Sah5fYypaaOIiAREIxIREQmIComIiAREhURERAKi\nQiIiIgFRIRERkYCokIiISEBUSCRfM7MzZvZDtp8nc/5W6JnZZjNbbmZNfe/jzOxFM1uXLeszOWxj\npJk9eN6yXmY22cwK+rZx0sxKhfJ3keinXluS3x13zjUO5gbNLM53M1igrnbO7fO9/gtwOdDQOZdp\nZkWB3+Xw/VHAk8Db2ZbdAoxyzh0HGpvZ5iDklHxOIxKRX+EbEfzJzBb7RgZ1fMsL+x4itMDMlphZ\nT9/yu8xsjJlNBKaZWYyZvWFmK81skm8UcKOZXWtm47Ptp6OZjcshSyHgfuBR51wmgHPuiHPu+Wzr\n3G5m832jjLd9z6KYAdQxs3LZttMB+CKY/1uJqJBIfnfuFM+5n+xPvdznnGsCvAn83rfsGeBb51wz\nsp7r8HczK+z7rBXw/5xz15DVCr4q0BC4z/cZwLdAXTMr7Xt/NzAih4zJwE/OuSO/9qGZ1SXraZ1t\nfKOrM8BtzrkzwDiynrYIcD3w3YW2I5JbKiSS3x13zjXO9jM622fnRgqLyCoKAJ2AJ83sB7LavScC\nlX2fTXfOnXu4UFtgjHPurHNuF/AdZPWJBz4Cbjezy8gqMFMuJbCZ3e0relvNrBJwLXAlsMCX61qg\num/1UWSdzsL3z1GXsi8Rf2iOROTCTvj+eYb//XfFgD7OuR+zr2hmLYCj2RddZLsjgIlAJlnFJqf5\nlPVAZTMr6julNQIY4Xs+d6xvXx845576le/OBsqZWSOgNf9bVESCRiMSkUszFXjUzAzAzK64wHqz\ngD6+uZKywFXnPnDO7SDruRDPkvWs7Ytyzh0D3gNeN7NE335jgXjfKt8AN5pZGd9nJcysiu+7DvgM\n+ACYfG6ORSSYVEgkvzt/juSlHNYfDBQAlvlGBIMvsN5Ysh4stIKsq6bmAYeyff4vYKtzbpWfOZ8B\ndgIrzGwJMJOs4rDDt41nyZrkXwZMB8pl++4ooBFZz+wWCTq1kRcJETMr4pz7xcxKAvPJmgzf5fvs\ndWCJc+69C3x3M9A02+W/ocoYlv1IdNOIRCR0Jvkmv2cCg7MVkUVACvDxRb67F/jm3A2JwXbuhkSy\nRldnQ7EPyT80IhERkYBoRCIiIgFRIRERkYCokIiISEBUSEREJCAqJCIiEhAVEhERCcj/B164BYhS\nJtPWAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x10f30fcc0>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "energy_range = [0.1, 10] * u.TeV\n",
    "pwl.plot(energy_range, energy_power=2, energy_unit='GeV')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Parameter errors\n",
    "\n",
    "Parameters are stored internally as covariance matrix. There are, however, convenience methods to set individual parameter errors on ``ParameterList``."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "PowerLaw\n",
      "\n",
      "Parameters: \n",
      "\n",
      "\t   name     value     error        unit      min max frozen\n",
      "\t--------- --------- --------- -------------- --- --- ------\n",
      "\t    index 2.600e+00 2.000e-01                nan nan  False\n",
      "\tamplitude 2.000e-12 2.000e-13 1 / (m2 s TeV) nan nan  False\n",
      "\treference 1.000e+00 0.000e+00            TeV nan nan   True\n",
      "\n",
      "Covariance: \n",
      "\n",
      "\tname/name index amplitude\n",
      "\t--------- ----- ---------\n",
      "\t    index  0.04       0.0\n",
      "\tamplitude   0.0     4e-26\n"
     ]
    }
   ],
   "source": [
    "errors = dict(\n",
    "    index = 0.2 * u.Unit(''),\n",
    "    amplitude = 0.1 * pwl.parameters['amplitude'].quantity\n",
    ")\n",
    "pwl.parameters.set_parameter_errors(errors)\n",
    "print(pwl)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "You can access the parameter errors like this"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[[  4.00000000e-02   0.00000000e+00   0.00000000e+00]\n",
      " [  0.00000000e+00   4.00000000e-26   0.00000000e+00]\n",
      " [  0.00000000e+00   0.00000000e+00   0.00000000e+00]]\n",
      "0.2\n"
     ]
    }
   ],
   "source": [
    "print(pwl.parameters.covariance)\n",
    "print(pwl.parameters.error('index'))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "You can plot the butterfly using the ``plot_error`` method."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<matplotlib.axes._subplots.AxesSubplot at 0x10fd1eb00>"
      ]
     },
     "execution_count": 10,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAZIAAAEOCAYAAACjJpHCAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMS4wLCBo\ndHRwOi8vbWF0cGxvdGxpYi5vcmcvpW3flQAAIABJREFUeJzt3Xec1PW1//HXWWABAZEOUgQERcQk\nKCom1sQYuIoYsWAsUVGwYOyiKTcxxmsssaNI7A1ENLFEJbkmxnjj9YIdUASxsIJSRJoUgfP74+z+\ndtg6u1N39v18PPYh853vzHx4ZLKHz+d8PueYuyMiIlJfRbkegIiINGwKJCIikhIFEhERSYkCiYiI\npESBREREUqJAIiIiKVEgERGRlCiQiIhIShRIREQkJQokIiKSkqa5HkA2dOzY0Xv37p3rYYiINCiv\nv/76cnfvVNt9BR1IzGwEMKJfv37MmjUr18MREWlQzOyTZO4r6KUtd3/G3ce2bds210MRESlYBR1I\nREQk8xRIREQkJQokIiKSEgUSERFJiQKJiIikJO8DiZn1NbN7zGx6wrXdzGySmU03s7NzOT4RkcYu\no4HEzO41s6VmNrvC9WFmNs/MFpjZ5TW9h7svdPcxFa695+5nAccBQ9I/chERSVamZyT3A8MSL5hZ\nE2AiMBwYCJxgZgPNbA8ze7bCT+fq3tjMjgReAV7M3PDDggWwfn2mP0VEpGHK6Ml2d3/ZzHpXuLwP\nsMDdFwKY2VRgpLtfAxxRh/d+GnjazP4CPJqeEVdt3Tp47z3o0gW6dYOivF8QFBHJnlz8SuwOLEp4\nXFJ6rUpm1sHMJgGDzeyK0msHm9mtZnYX8Fw1rxtrZrPMbNayZctSHrQ7fP45zJ0Lq1en/HYiIgUj\nF7W2rIprXt3N7r4COKvCtZeAl2r6EHefbGZLgBHFxcV71X2YVdu4EebPh/btoWdPaFrQ1cpERGqX\nixlJCdAz4XEPYHEmPijVWltbtsBVV8GXX1Z+7ssvYfZsWL48xUGKiDRwuQgkM4H+ZtbHzIqB0cDT\nmfggMxthZpNXrVpVr9fPnAm//S0cdRQ8/ngElkRbtsAnn8C8ebBhQxoGLCLSAGV6++8U4FVgVzMr\nMbMx7r4ZGA/MAN4Dprn7nEx8fqozkqFD4e23YcAAuPZaOOUUeOedyvetXRu5k8WLI5ciItKYmBfw\nb76EfiRnzp8/v97v89Zb8PzzcPPNsHQpjBgB550XeZKKWrSAXr2gTZv6j1tEJB+Y2evuXutZvYLe\nyJqufiRmcNhhMH16zEqeew5GjYJp0yovd23YAB98EEtemzen9LEiIg1CQQeSdNtuO/jZz+Cxx2C3\n3eC66+Dkk2P5q6Lly2HOnKoT9SIihaSgA0mqyfbq9O4NEyfC738Pq1bBmDHwm9/AihXb3rd5M3z0\nUWwX3rgxrUMQEckbBR1IMtlq1wwOPTR2c516KrzwQix3TZ1aeUlr9epIxn/+uZLxIlJ4CjqQZGpG\nkmi77WD8+AgggwbBDTfEctebb25739at8NlnUWpl3bqMDUdEJOsKOpBkckZSUe/ecNttkTdZswbO\nPBN+9avKBxbXr4f334dFiyon6kVEGqKCDiTZZgbf/34sd512Gvz3f8dy16OPVl7uWro0lru++io3\nYxURSRcFkiS0bl23+1u2hHPPjd1d3/oW3HgjnHgivP76tvdt2gQffhg/33yTvvGKiGRTQQeSdOVI\ndt4Zdt217ocMe/WCW2+NvMnXX8O4cfDLX0LFYsRffRVbhZcuTWmYIiI5UdAn28sMGTLEZ82alZb3\nWr06SqHUNWG+YQPcfz888AA0awZjx8Lo0ZWrB7dqBTvtFLMaEZFc0sn2DNl++6i91a9f7NhKVosW\ncNZZcRp+8OAot/KTn0DF+FbWROuzz2Knl4hIvlMgqae2beN0e9++ESSS1bNnBJE//CFmKWedBT//\n+bbLWmqiJSINSUEHkmycI2nXDnbfHfr0gebNkx0XHHRQzE7OPBNeegmOOQYeemjb3V1lTbQ++kh1\nu0QkfylHkkbuUSZlyZLYkZWskpKYofzrXxGQLr0U9tln23uaNoXu3aFjx/SOWUSkOsqR5IBZ/KIf\nNKhubXh79ICbboqfTZvgnHPgiivgiy/K79m8OSoKf/CBmmiJSH5RIMkAM+jcGfbYI2YRTZok97oD\nDojlrnHj4OWXY7nrgQe2PWOyZk0k45csUd0uEckPCiQZVFQEXbtGQOnaNR7XpnnzyJtMmwZ77x1l\nV044AV57rfyerVtjC/LcudGdUUQklxRIsqBJk5iZ7LFHzFTMan9N9+5xIv7mm2NZ69xzYcKE2M1V\nZsOG6Bf/6aeq2yUiuVPQyfZ0tdpNt02bYmlqxYrklqc2bowdXffdF0FozJgouVJcXH5Ps2aRl2nX\nLnPjFpHGJdlke0EHkjLZ2rVVVxs2xBLVypXJ3b9kScxS/vGPKL9y2WUwdOi297RtG88lBhkRkfrQ\nrq0GoEWLONA4cGAEgNp06wbXXx95E4g+KJdeuu1y16pV5XW7GsG/EUQkDyiQ5IGWLaPkSrKFIffb\nLxppnXsu/PvfUar+nnvKz65s3Rr9Tt5/P4pFiohkkgJJHmndGnbZBfr3r72OV3Fx9Dx54gn43vfg\nzjvh+OPhf/6n/J6vv45gUlKiul0ikjkKJHlo++2Tr+PVtWt0Zbz99thefP75cPHFkXuBWN764otY\n7spgpRgRacQUSPJYu3aRP9lpp9qT50OHxnLX+PFx5uTYY+Huu2PHF8Sy14IFUbdLTbREJJ3yPpCY\nWV8zu8fMple43srMXjezI3I1tmyoS9mVZs3g1FNh+nQ48ECYNCmWu155pfyeL7+M2UnFXvIiIvWV\n0UBiZvea2VIzm13h+jAzm2dmC8zs8prew90XuvuYKp6aAExL53jzWVnZlUGDYMcday670rUrXHMN\nTJwYgeeCC+DCC6PHCcThxU8+icOMqtslIqnK9IzkfmBY4gUzawJMBIYDA4ETzGygme1hZs9W+Olc\n1Zua2aHAXOCLqp4vZE2axDbgQYOgS5eaT8nvuy9MmQLnnRcNtI47DiZPLg8ea9dGmRXV7RKRVGQ0\nkLj7y8CXFS7vAywonWlsAqYCI939XXc/osJPdV3MDwGGAj8BzjSzvF+iS7emTaNq8KBBsfRVXUBp\n1gx++tPy5a7Jk2O56+WX43l31e0SkdTk4hdwd2BRwuOS0mtVMrMOZjYJGGxmVwC4+y/c/QLgUeCP\n7l5pc6uZjTWzWWY2a9myZen9G+SR4uJIxg8cWHN5lC5dYrnrzjvjNRddFEteJSXxfFndrk8+Ud0u\nEambJDtmpFVV/3audmHF3VcAZ1Xz3P01vG6ymS0BRhQXF+9V10E2NGWn5Neti1zImjVV37f33rHc\nNWUK/PGPsdz105/GT4sWkYRftUp1u0QkebmYkZQAPRMe9wAWZ+KD3P0Zdx/bNpn6IwWiVavaDzU2\nbQonnxzLXYccUh5QXnoplrq++QYWLoztwnXp9CgijVMuAslMoL+Z9TGzYmA08HQmPigbPdvzVTKH\nGjt3hquvjm3CLVrAJZfEgcZFpQuPZXW7vvhCyXgRqV6mt/9OAV4FdjWzEjMb4+6bgfHADOA9YJq7\nz8nE5zfGGUlFiYcamzWr+p4hQ+DRR2OL8Ntvx+zkzjsjb7J1a+RRVLdLRKpT0GXk87UfSa5s3RpV\ngT//vPqE+vLlcMst8Pzzsc344ovhoINiV1jZWZYdd0yu26OINGzqR5IgX/uR5MqWLRFMli6tvpjj\nG2/AtdfChx/Cd78by169esVzZTvFtt8+e2MWkexTPxKpVlnr35rOoOy5JzzySGwTfvvtOHsycSKs\nXx8J+Pnzo27X5s3ZH7+I5JeCDiSNOdmejGbNys+g7LBD5eebNoWf/CRK1f/wh9Hq95hj4MUXI/mu\nul0iArUsbZnZG0m8xzJ3/1H6hpR+WtpKTm1nUN56K5a75s+P8iuXXgq9e8dzbdrE0ldtZe9FpOFI\nS47EzOYAR9b0euBJd/9W3YeYeUq218+qVRFQ1q+v/NzmzXH+5M47o0T9iSfCmDFxZsUsEvRdu9Zc\nA0xEGoZ0BZKD3f2lWj6o1ntyTTOS+vnyywgoVR1KXLEimmk980yUX7nwQvjBDyKAtGgRS2atW2d/\nzCKSPtq1lUCBpP7cYdmyqBBcVWL97bdjueuDD2CffWK5q0+feK5Tp0jq11TyXkTyV9p2bZnZ3mZ2\ni5m9YWZLzGyhmT1tZuPMrE16hiv5KrEPSrdulc+PfPvb8NBDMGECvPcejB4d51DWrYsANGcOrFyZ\nm7GLSHbUtrT1LLACeAqYBSwFWgC7EKXcDweuc/dnMz/UulOOJP2++SZmJ8uXVy6b8uWXcNttsdzV\nqVNUFz7ssAhGbdtGMr62lsEikj/SlSPp4u41No8ys8419A3JC1raSr+NGyN/UtVs4513Yrlr3rwo\nv3LppbDzzjGb6d49Zjgikv/StbR1hZntU9MN+R5EJDOaN4+CkAMGxNbfRN/6Fjz4IFx+eeROfvIT\nuPlmWL06CkK+/37VO8JEpGGqLZAsAiaa2YdmdrWZDcrGoKThKCtb368ftGxZfr1Jkzi8+MQTMGIE\nPPxwPH7hhejE+N57UQyyuhItItJwJLVry8x2Jsq9jybOjjwKTHX3hZkdXmqUI8m+FSuidW/FLcOz\nZ8N110VL3z33hMsui+DTvHnkTlS3SyT/ZGz7r5ntBdwNfMvdG8TGTuVIsmvr1vItw4lVhrdsgT//\nGe64I2Ylxx8PY8fGeZP27aMrY9Nc9OwUkSqltWijmTUxs+Fm9gDwF2AhcHyKY5QCVVQUhxT32CP+\nW7ZluEkTGDUqlruOPDLa/Y4aBc89FzOZOXPivyLSsNS2a+sQ4ASiTMqbwFSiJEo11Zjyk2YkubVp\nUyx3VQwSc+bE7q65c2Hw4Fju6t8/kvc77RTLXiKSO+na/vsvIh8y3d2XpXF8WaVAkh/Wr48E++rV\n5de2boWnn47zJ2vXRnfGceMiZ9K1q+p2ieRS2nMkZjYU2MXdHzSzDkArd/80xXFmhQJJflmzJgJK\nYuveVasid/Lkk5EvOe88OPzwKAbZq5fqdonkQloDiZn9EvgesLO772Jm3YHH3H3/1IeaOdq1ld+q\nKgr53nux3DV7dpRfmTAhtherbpdI9qU7kLwFDAbecPfBpdfeydfy8RVpRpK/3Mv7yJcVhdy6Ncqs\n3HZbLIMdcwycfXb5zq527XI7ZpHGIt2tdjd6RBwvffPtUhmcSBmz2Nk1aFD5Dq+iIhg5MnZ3jRoV\n/U+OPjqWvRYsiD7yVZW2F5HcSDaQPGlmE4G2ZnYa8Ffg3swNSxqbJk2gRw/YfXfo0CGutW0bS1sP\nPhgzkSuvhDPOgP/939jxtVTFeUTyQl2S7cOBw4iT7TPc/flMDiydtLTV8Hz9deRPynZ4bd0Kf/kL\n3HprJOZHjYrlrm7dIhm/nebIImmXru2/f3X3w9I6shxQIGm4Vq+OHV5lRR7XrIFJk+Dxx2OL8Hnn\nxeHGrl1hxx0r90sRkfpLV46kU5rGI1Iv228PAwdC797Ry6RNmyhL/9BDcWjxqqvgtNPgpZfiYGPi\nGRURyY7aZiQLgUuqe97dn8zEoCqMoS/wC6Ctux9Teu1g4CpgDlE88qWa3kMzksKwdWv5Dq8tW2LH\nV9ly18qVkZA/55xo9au6XSKpS9eMpC1wBDCiip8jkhjEvWa21MxmV7g+zMzmmdkCM7u8pvdw94Xu\nPqbiZWAt0a2xpLZxSGEoKoolrEGDojlWUREccUTs7jr+ePjTnyKY3HNPNNdavjzXIxZpHGqbkbzh\n7nvW+83NDiR+4T/o7oNKrzUBPgB+SASBmUQ9rybANRXe4vSyxllmNj1hRlLk7lvNrAtwo7ufWNM4\nNCMpTBs3Rv7kq6/i8fz5cZjxrbdiOWzCBNh331gCa9Eit2MVaYiSnZHUNvlPqcqRu79sZr0rXN4H\nWFDWy8TMpgIj3f0akpjllL5vWTuklYBK+zVSzZtHC9+1ayOg9O8Pf/wjPP883HILnHoqHHUUjB8f\nnRy7dVPdLpFMqG1p6+QMfGZ3ovNimZLSa1Uysw5mNgkYbGZXlF472szuAh4Cbq/mdWPNbJaZzVq2\nrMHWm5QktG4dgaJv35h5/Md/xHLXCSdEQcijj4bbb4d3341dXyKSXjXOSNx9dk3P11NV/yasdn3N\n3VcAZ1W49iRQY6Lf3Seb2RJgRHFx8V71Gag0LO3awQ47REJ+yRK46KLYGnzddXDNNfDUU1Gq/qCD\n4vCjkvEi6ZGLXfclQM+Exz2AxZn4IHd/xt3Htm3bNhNvL3morORKWVOt/v3hrrvgd7+LAHPaaXDh\nhfDKK2qiJZIudQ4kZtbOzFIp1jgT6G9mfcysmOgD/3QK71ctMxthZpNXrVqVibeXPJZYcqV9exg2\nLGp2nXhiFIQcORJuuCGqDW/cmOvRijRsybbafcnMtjez9sDbwH1mdmMSr5sCvArsamYlZjbG3TcD\n44EZwHvANHefU/+/QvU0I5HmzSN3suuuMUO54IJo8bvLLvD730f+ZOrUWApLslqQiFSQbBn5N919\nsJmdAfR09183hDLy6kciFa1cGTW8NmyAv/4Vbr4Zli2LXMoll0TLXzXREgnpLiPf1My6AccBz6Y0\nsizSjEQqatculrt69ozdXdOnw8knxwn5ww+PkisLF8bJeRFJTrKB5LfEUtQCd59ZWrZE/8SXBimx\nB0qfPrHcNXUq7LZb7PAaPhwefjg6OIpI7ZIuI98QaWlLkrFhQ/kJ+RdfhJtugi++iPIrV1wBe+0V\nuRaRxiatrXYbOpVIkWSsWQOLFsW24HvvjVlJixbR9+S882IXmE7GS2OS7hxJg6Ttv1IXbdpEja7d\ndouzJlOnRj7l+uvhBz+IwKKT8SKVFXQgUbJd6qNjx8ifDB0Kd9wRhSBXr4ZTTomyKzNnwubNuR6l\nSP6oNZCY2QAz+4GZta5wfVjmhiWSW0VF0L17BJRjj43dXaedBjNmwCGHwM9/Hn1RRKSWQGJmPwOe\nAs4DZpvZyISn/yuTAxPJB8XFsbNr8OCo0/XYY1F+5frr4cAD4ZFHytsAizRWtc1IzgT2cvejgIOB\nX5nZ+aXP5X3aUTkSSZdWraLC8MEHR+2u666LfMlJJ8Exx8Abb0QHR5HGqLZA0sTd1wK4+8dEMBle\nWh4l7wOJciSSbu3bx3LXiSdGqfrTT4e//S1mJ5dfHqfkRRqb2gLJ52b2nbIHpUHlCKAjsEcmByaS\nr4qKYMcdYcgQ+OUvY7nr29+O5a6hQ2O5S4UgpTGpLZCcAmyTUnT3ze5+CnBgxkYl0gCU5U9++MPo\nzHjDDXG48aSTojPjm2+qEKQ0DrU1tiop+7OZtSP6iJS9Ju9TjAkn23M9FClgrVvH+ZMuXeCAA+Du\nu+GBB2D//eGcc+J0fPv2uR6lSOYkW/33KuBU4EPKuxm6u38/c0NLH51sl2zZsiW2Bc+aBX/4A/zr\nX1HG/re/heOOg2bNcj1CkeQle7I92WajxwE7u/um1IYlUtiaNInzJx07xqn4Z56JgHLSSfDgg3G4\n8dvfVqkVKSzJnmyfDeyQyYGIFJLmzaFfv9jV9ec/w9ix8M9/wne/G31PVFlYCkmygeQa4E0zm2Fm\nT5f9ZHJgIoWgbduoHnzllfD447HT68Yb49rDD8M33+R6hCKpS3Zp6wHgWuBdQMeuROrADLp2jb7x\ngwbB00/HVuGTTy5f7vrOd7TcJQ1Xssn2f7r7QVkYT1qpH4nko7VrYf58mDwZ7rsvAshZZ8WZlA4d\ncj06kXJp7UdSepJ9I/B06X8BcPc3UhlktmjXluQb9zgFP2tWzE5eegl69YpWvyecoN1dkh/SHUj+\nUcVlbf8VSdHmzfDZZ/DUU3GgcdEi+P73Y7lrr7203CW5pQ6JCRRIJN+tW1e+3HXvvXFt3Dj4xS+g\nc+fcjk0ar7R2SDSz/zKzHRIetzOz36UyQBEp16pVJNyvvBKefDK2Cd96K+y9d+RRVLtL8lmy23+H\nu/tXZQ/cfSXwH5kZkkjj1akTHHZYzEpuuy2unX46DB8Or72mUvWSn5INJE3MrHnZAzNrCTSv4X4R\nqaemTaF3bzj11MidnHsuvPoqHHQQ/OxnsHhxrkcosq1kA8nDwItmNsbMTgf+RpwtyTgz62tm95jZ\n9IRrRWZ2tZndZmY/zcY4RLKtdesop/KrX0Xvk/33h4kTYd99I5fy9de5HqFISCqQuPt1wO+A3YDd\ngatKr9XIzO41s6VmNrvC9WFmNs/MFpjZ5bV89kJ3H1Ph8kigO/ANUFL5VSKFwSyqCh96aHRmvP32\nqOc1blwsd73ySuz8EsmlGndtmZl5Ldu6arrHzA4E1gIPuvug0mtNgA+AHxJBYCZwAtCEKMWS6HR3\nX1r6uunufkzpny8HVrr7XYnXq6NdW1IoVq2ChQsjAX/33ZEzOe206M7Yq5e2C0t6pWvX1j/M7Dwz\n61XhzYvN7Ptm9gBQ7dKSu78MVCxPtw+woHSmsQmYCox093fd/YgKP0ureesSYGXpn7fU8ncQKRht\n28Zy14QJsdx14IEwaRJ873twxx0RaESyrbZAMoz4RT3FzBab2Vwz+wiYT8wibnL3++v4md2BRQmP\nS0qvVcnMOpjZJGCwmV1RevlJ4EdmdhvwcjWvG2tms8xs1jI10pYCUlQUpeoPOSR2dt1xR3RrHD8e\nRoyAv/89OjWKZEttHRI3AHcAd5hZM6JX+/rErcD1UNXku9rlM3dfAZxV4drXQMW8ScXXTQYmQyxt\n1X2YIvmtRQvYddfofbL33vDQQ9Hyd/jw2PF1ySXRVKtJk1yPVApdsru2cPdv3H1JikEEYgbSM+Fx\nDyAjGxrNbISZTV6l+b4UsA4d4jDjhRfGctfBB8eurkMOieT80qXqHS+ZlXQgSaOZQH8z62NmxcBo\nohikiNRT2dmT/fePfieTJkVzrQsugKOPhhkzYM2aXI9SClVGA4mZTQFeBXY1sxIzG+Pum4HxwAzg\nPWCau8/JxOe7+zPuPrZt27aZeHuRvNOmDQwcCEceCY8+GoHkzTdh5Ei46CKYPVvlViT9atv+OwN4\nAXje3d/P2qjSRP1IpDHbsAE+/TS2C99yC7zwAnTrBhdfDMcdBzvuqPyJ1Cwt1X/NrCuxc2sYsAvw\nGhFYXnT3tWkaa8bpHIk0ZsuXR6n6116D666DDz+MopATJsB++0WyXudPpCppLyNvZkXAvsBw4AfA\neuCvyZxwzzUFEmnsNm+OXidLl8K0aXFKftMmOOkkOOec2P3Vpk2uRyn5JuP9SMysI/Ajd3+kXm+Q\nBVraEtnW6tWx3PXZZ1Gm/rnnop/8RRfBj38MPXtGkl4E1NhqG5qRiJTbujUqCC9dGon4a6+NplpD\nh8Jll8WZlG7dlD+RNDe2aqh0jkSksqIi6NEDBgyIXMlDD8XhxXffjST8r34FM2dGT/lG8O9MSQPN\nSEQaMfeYmSxeHIHjttvg2Wej4vCFF8Lhh0cxyO23z/VIJRcyPiMxs9Pq+1oRyQ9lZep33x369IHf\n/CaqCrdtGxWFx4yBv/0NFixQ/S6pXipLW1embRQiklPFxdC/fwSTIUPgwQfh0kth7lwYPRquugpm\nzYpEvfqfSEW1nSN5p7qngF3cPa/3d2jXlkjdbd4MJSWwYgV8+WUsdz3zDHTuHMtdP/pRHGbs3Fnn\nTwpdug4kfgH8iPLeH///KeDf7r5jSqPMEuVIROpu9Wr45JM4b/LOO7G7a9682NV12WWRrO/RA3bY\nIdcjlUxJV47kWaC1u39S4edj4KU0jFNE8tT220fupEuXaKb14INxGv7992O567rrIsB88IH6xzd2\n2rUlIrX6+mv4+GNYvx5Wrozy9E89BZ06RWHIww6LP++4IzRrluvRSrroQCLKkYikkzt88QUsWRKH\nGt99N5a73n8f9torlrv694+T8l26xHkVadjSlSN5w933rOWDar0n1zQjEUmfDRsid7J2LWzZAn/6\nU7T7XbculrzOPBPat492wO3b53q0kop0BZL1RH/2am8B2rp7r7oPMXsUSETSb9myqNm1ZQt89RVM\nnAh//nN0bDz/fBg2DFq3jvpdrVrlerRSH+kKJDsl8Vlb3L2kLoPLNgUSkcz45ps4W/JVaQPu2bMj\nCT93Luy5Zyx39esH7drFDq/i4tyOV+pGOZIECiQimbVyZflhxS1bIhE/cWIsfx13HIwbF7vAOneO\nHIoKQjYMKtooIlnTrl1sFe7QIYLE0UfDE09Ei9+pU2HUqKjhtWQJzJkTzbakcBT0jES7tkSyL/Eg\nI8Qy17XXRgAZPLh8d1fLlrHcpYKQ+SutS1tmNtDd51a4drC7v1T/IWaPlrZEsmvr1kjEL11a/vjp\np+P8yZo1cOyxsdzVpk0UiOzRA1q0yO2YpbJ0L21NM7MJFlqa2W3ANakNUUQKVVFR7NYaMCACRFER\nHHVULHf9+Mfw2GPly10rV8asZdEiFYRsqJINJPsCPYF/AzOBxcD3MjUoESkMrVrBwIHRcdGsvDz9\nQw/FOZPf/AbOOCMONS5dGru+vvhCDbUammQDyTfAeqAl0AL4yN23ZmxUIlIwzKJ0ym67wXbbxbUB\nA+Cee+A//zN2e518cuRRVq6MysNz5pRvKZb8l2wgmUkEkr2B/YETzGx6xkYlIgWnZcvyisFFRfFz\n5JHw5JOxzPXEE/Hfp56Kml4ffqiCkA1Fssn2Ie4+q8K1k939oYyNLI2UbBfJLxs3xs6uNWvKr82b\nF7OSd96BPfaISsMDBsRzHTrEUpgKQmZXundtVVkCxd0/rcfY6sTM+gK/IEqxHFN67QDgRKApMNDd\nv1vTeyiQiOSnxDIrELu7nnsObr01lrlGjYKzz47cSlGRCkJmW7oDybuAE7W1WgB9gHnuvnstr7sX\nOAJY6u6DEq4PA24BmgB3u/vvkxjD9LJAknDtKKCLu99V02sVSETy16ZNkSdZtar82po1cNddMG1a\nnDMZPz6WwYqKosyKCkJmR1q3/7r7Hu7+rdL/9gf2AV5J4qX3A8MqDKwJMBEYDgwk8i0DzWwPM3u2\nwk/nWt7/J8CUZP4OIpKfiou1w8bqAAAQ3ElEQVSjHlefPtC0aVxr0wYuuQQefhh22gl+9zs4/XR4\n770IPB99FDu91q7N7dglNK3Pi9z9DTPbO4n7Xjaz3hUu7wMscPeFAGY2FRjp7tcQs5eklC63rXL3\n1UkPXETyVvv2Mfv49NNY1gLYZRe4++7y5a5TTolzKOecE8/PmxflWbp3h+bNczf2xi6pQGJmFyU8\nLAL2BJbV8zO7A4sSHpcQ51Sq++wOwNXAYDO7ojTgAIwB7qvhdWOBsQC9euV1lXsRKdW0KfTtG1t/\nP/00qgubweGHw0EHlS93vfhiLHeNHBlB56uvIneigpC5kWzKqk3CT3PgL8DIen6mVXGt2kSNu69w\n97PcfeeEIIK7/9rd/13D6ya7+xB3H9KpU6d6DlVEcmGHHcqLQJZp3RouvhgeeQR23hmuvhpOPTXO\nnLjD55/HgcZly3SgMduSmpG4+5Vp/MwS4pR8mR7ESfm0SyjamIm3F5EMatIEeveOJa/EIpD9+sXM\n5IUX4OabI5gcdRSce24EoE8/jVPyPXuqIGS21BhIzOwZap4tHFmPz5wJ9DezPsBnwGgiaS4iUsn2\n28fspKQkZhsQy13Dh8MBB8Af/xil6v/+99gq/OMfRzvg+fPjtT17qiBkptXWIfGgml7s7v+s8c3N\npgAHAx2BL4Bfu/s9ZvYfwM3E9t973f3qOo67TrT9V6QwrF0LH38cBxoTffhhdGZ8/fUoxTJhAgwq\nPXBgBh07RpmWpvXaXtR4pavVbq9sHDrMFPUjESk8W7fC4sWxfJX468sd/vrXWO5atiwS8ePHx64u\niKWybt2iS6NVlamVStIVSN5w9z1L//yEu49K4xizRjMSkcKzbl3kTtavr3z97rvh0Uej+vDZZ0fH\nxrLdXM2bx3bhsgAj1UvXgcTEuN03tSGJiKRPq1axjFVWoj7x+vnnw5QpsOuuUb/rlFOihhfEstjC\nhXEGRQUh06O2QOLV/LlBMLMRZjZ5VWLtBREpGGUl6gcMKC9RX6ZvX7jjDrjmmjhrcvrpcOWV8OWX\n8fzatXFS/uOP47yK1F9tS1tbgHXEzKQlUBa/DXB3bxCb67S0JVL4ys6SLFlS+RzJ119H/5NHHokd\nXGefHQUhy5LvRUXlBxpVELJcWos2NlRKtos0Phs2xCxj3brKz338cezu+r//i/Irl10G3/lO+fPN\nmkX+JPEgZGOmQJJAMxKRxsU9dnUtXhy7vCo+9+KLcNNN0db38MPhZz/bNnhst12cP2ndOrvjzjdp\nrf4rItKQmMVS1cCBlYOBGRx6KEyfDqedBjNmxK6uqVNh8+a45+uvIxm/cGHlMytSmQKJiBSs5s1j\n51bPnpVzHy1bRlmVxx6Ljow33BC94998s/yelSujlldJSXnzLamsoAOJdm2JCMQhxIEDo89JRTvt\nBLfdFrmT1avhzDPhP/8Tli+P591jCUwFIaunHImINCoV2/smWr8e7r03GmoVF8O4cXDccduWVmnR\nAnr0iPa/hU45EhGRKnTqFLOTqioDJy53ffvbcOONcOKJ8MYb5fds2AALFkRRyA0bsjfufKZAIiKN\nTnEx9O8fy1pVNcLq1QtuuSXyJl9/DWPHwi9/Wb7cBbEMNndulK0vS9I3VgUdSJQjEZGadOwYs5Oq\nlqnM4OCD4fHH4Ywzokz9qFGx7FUWONxjqWz27MijNIJMQZWUIxERAVasgEWLqt+dtWhRzFD+53+i\n/Mpll8GQCtmD5s0jf7LDDpkfbzYoRyIiUgcdOkQDreqS6D17Ron6P/whciNnnQU//3kcfCyzcWP0\nRvngg8ZVEFKBRESkVLNm0cq3T5+qm2CZwUEHwbRpsU34pZfgmGPgoYe2zZOsWdO4CkIqkIiIVNC+\nfeROqluiatEitgZPmwZ77RWJ+RNOgJkzt71vxYrInyxZUrlUSyFRIBERqUKzZrDzztXPTiDyITfd\nFD+bNkVV4SuuiMR7mbKOjnPmlJewLzQFHUi0a0tEUtW+feROakqgH3BAzE7GjYOXX47lrgce2HZZ\na9Mm+OgjeP/96IVSSLRrS0QkSStX1n5u5LPP4iDjP/8JvXvDpZfCvvtWvq9du5jRFBdnbLgp064t\nEZE0a9cuZic19Xvv3j12dt18cwScc8+FCROi6VailSsjf1JduZaGRIFERKQOmjaNcyR9+1afOwHY\nf/8otXLWWfDKK7Hcdd99scRVpqyr45w5256ab2gUSERE6iGZ2Unz5nEqfvp02G8/mDgRRo+GV1/d\n9r5vvoFPPomSK2vWZHbcmaBAIiJST8nOTrp1g+uvh1tvjVnIeedF7qTictf69XGYccGChlUQUoFE\nRCRFycxOAL773VjuOucc+Pe/o3bXPfdsu9wFsGpVzE5qKtmST/I+kJhZXzO7x8ymJ1zrZWZPm9m9\nZnZ5LscnIgLJz06Ki+H00+GJJyKw3HknHH981PBKVNZ3fvbs+G8+b7DNaCAp/UW/1MxmV7g+zMzm\nmdmC2gKBuy909zEVLu8C/MXdTwcGpnnYIiL1VjY7qa1wY9eusdx1++1ReuX88+Hii+PwYqLNm2Nm\nMnduzFTyUaZnJPcDwxIvmFkTYCIwnAgCJ5jZQDPbw8yerfDTuZr3fRMYbWZ/B/6RwfGLiNRZ06a1\nn4ovM3QoTJ0K48fDa6/BscfC3XdHAchEiQ211q/P3NjrI6OBxN1fBioWBdgHWFA609gETAVGuvu7\n7n5EhZ+lld40nAb82t2/Dxyeub+BiEj91Vazq0xxMZx6auzu2n9/mDQpWvy+8krle1evjoKQ+dRQ\nKxc5ku7AooTHJaXXqmRmHcxsEjDYzK4ovfwC8LPS6x9X87qxZjbLzGYtW7YsPSMXEamjZGp2lena\nFa69NrYJN20KF1wAF14IJSXb3pdvDbUyXiLFzHoDz7r7oNLHxwI/cvczSh+fDOzj7udlagwqkSIi\n+aDsvEgyuY5vvoFHH41lri1bYsZyyilRebiiTDXUyucSKSVAz4THPYDF1dybEhVtFJF8UtbvpHfv\nqnvFV7z3pz+N5a6DDoLJk2N317/+VfneXDfUykUgmQn0N7M+ZlYMjAaezsE4RERyorZujIm6dIFr\nroltwsXFsdRV1XIX5K6hVqa3/04BXgV2NbMSMxvj7puB8cAM4D1gmrvPycTnu/sz7j62bTL/a4mI\nZFHZ7GSnnWqfnQDsvTdMmRJ5k9dfj2T8XXdVfQI+2w21CrqMvJmNAEb069fvzPnz5+d6OCIiVdq0\nKXInq1cnd/+yZVFdeMYM2HFHuOiiWP4yq3xvjx4xq6mPfM6RZI1mJCLSEBQXQ//+0KsXFCXxW7lT\nJ7j66tgm3KIFXHJJHGhctKj212ZCQQcSEZGGpFOnOHfSpk1y9w8ZEju7LrwQ3n47lrvuvDP7BR8L\nOpBo15aINDTNm8Muu0DPnsnNTpo2hRNPjNpdhx4aRSCPPRb+8Y/snS8p6BxJGZ0jEZGGaOPG2IFV\nlx7vb7wB110X5VT22w9uuCGKQ9aHciRoRiIiDVvz5rDrrpEwryqRXpU994SHH44CkO+8A5ddltkx\ngmYkIiINwoYN8NFHdTtwuHx5nFUZOrR+n6kZiYhIAWnRAgYMiO2+yc5OOnaMGl+ZpkAiItJAmEXb\n3gEDoGXLXI+mXEEHEuVIRKQQbbcd7LZbVAtOdnaSSQUdSHQgUUQKlRl07x7J+KoqAmdTQQcSEZFC\n16pVzE46V9dPNgsUSEREGriiojjAuMsuUW4l65+f/Y/MHuVIRKQxadMmSqx07Jjdzy3oQKIciYg0\nNk2aRGn6fv2iVH021NJBWEREGqK2bWN2smlT5j9LgUREpEA1bRo/mVbQS1siIpJ5CiQiIpKSgg4k\n2rUlIpJ5BR1ItGtLRCTzCjqQiIhI5imQiIhIShRIREQkJQokIiKSkkbRatfMlgFfATVt32pbw/Md\ngeXpHleG1fT3yefPSuW96vraZO9P5r7a7im07xdk7zum71fuvl87uXunWu9y90bxA0yu7/PArFyP\nP91/33z9rFTeq66vTfb+ZO5rbN+vdP/vnq3P0fcrMz+NaWnrmRSfb2iy+fdJ52el8l51fW2y9ydz\nX2P7fkH2/k76fuX596tRLG2lysxmufuQXI9DCpO+X5JJ2fh+NaYZSSom53oAUtD0/ZJMyvj3SzMS\nERFJiWYkIiKSEgUSERFJiQKJiIikRIEkBWbW18zuMbPpuR6LFAYza2VmD5jZH83sxFyPRwpPJn5v\nNdpAYmb3mtlSM5td4fowM5tnZgvM7PKa3sPdF7r7mMyOVBq6On7Xjgamu/uZwJFZH6w0SHX5jmXi\n91ajDSTA/cCwxAtm1gSYCAwHBgInmNlAM9vDzJ6t8NM5+0OWBup+kvyuAT2ARaW3bcniGKVhu5/k\nv2Npl4W28PnJ3V82s94VLu8DLHD3hQBmNhUY6e7XAEdkd4RSKOryXQNKiGDyFo37H3pSB3X8js1N\n9+fri7qt7pT/axDi/9Tdq7vZzDqY2SRgsJldkenBSUGp7rv2JDDKzO6kMMuqSPZU+R3LxO+tRjsj\nqYZVca3aE5vuvgI4K3PDkQJW5XfN3dcBp2V7MFKQqvuOpf33lmYk2yoBeiY87gEsztFYpLDpuyaZ\nlrXvmALJtmYC/c2sj5kVA6OBp3M8JilM+q5JpmXtO9ZoA4mZTQFeBXY1sxIzG+Pum4HxwAzgPWCa\nu8/J5Til4dN3TTIt198xFW0UEZGUNNoZiYiIpIcCiYiIpESBREREUqJAIiIiKVEgERGRlCiQiIhI\nShRIpNEzsy1m9lbCT43tA7LFzD42s3fNbIiZ/al0bAvMbFXCWL9bzWvPMLOHKlzrUlpqvJmZPWZm\nX5rZUdn520gh0zkSafTMbK27t07zezYtPRCWynt8DAxx9+UJ1w4GLnH3GqtRm1k7YD7Qw903lF4b\nD+zh7uNKHz9M9D75cyrjFNGMRKQapTOCK83sjdKZwYDS661KGwnNNLM3zWxk6fVTzexxM3sG+KuZ\nFZnZHWY2p7SHzXNmdoyZ/cDM/pTwOT80sydTGOfeZvZPM3vdzJ43sy7uvhL4N3B4wq2jgSn1/RyR\n6iiQiEDLCktbxyc8t9zd9wTuBC4pvfYL4O/uvjdwCHC9mbUqfW4/4Kfu/n2i22FvYA/gjNLnAP4O\n7GZmnUofnwbcV5+Bm1lz4BZglLvvBTwMXFX69BQieGBmPUvH8nJ9PkekJiojLwLr3f071TxXNlN4\nnQgMAIcBR5pZWWBpAfQq/fPf3P3L0j/vDzzu7luBz83sHxB1vEvzFyeZ2X1EgDmlnmPfDdgd+G8z\nA2hCVH2FKNB3q5m1Bo4nai1trefniFRLgUSkZhtL/7uF8v+/GDEDmJd4o5ntC6xLvFTD+95HNK7a\nQASb+uZTDHjH3Q+o+IS7rzOz/ya64o0Gzq7nZ4jUSEtbInU3AzjPSqcAZja4mvteIbodFplZF+Dg\nsifcfTHRG+KXRL/t+ppLdL3bp3QsxWa2e8LzU4BLgR3cfWYKnyNSLQUSkco5kt/Xcv9VQDPgHTOb\nTXlOoqIniGWm2cBdwGvAqoTnHwEWuXu9e2i7+0bgGOBGM3sbeBPYN+GWF4hlt6n1/QyR2mj7r0gG\nmVlrd19rZh2A/wO+5+6flz53O/Cmu99TzWs/psL23zSPTdt/JS00IxHJrGfN7C3gX8BVCUHkdeBb\nxC6r6iwDXjSzIekelJk9BnyPyNGIpEQzEhERSYlmJCIikhIFEhERSYkCiYiIpESBREREUqJAIiIi\nKVEgERGRlPw/k66QW0ULSJMAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x10fce6b00>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "ax = pwl.plot_error(energy_range, color='blue', alpha=0.2)\n",
    "pwl.plot(energy_range, ax=ax, color='blue')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Integral fluxes\n",
    "\n",
    "You've probably asked yourself already, if it's possible to integrated models. Yes, it is. Where analytical solutions are available, these are used by default. Otherwise, a numerical integration is performed."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "1.2186014196061302e-12 1 / (m2 s)\n"
     ]
    }
   ],
   "source": [
    "integral = pwl.integral(emin=1 * u.TeV, emax= 10*u.TeV)\n",
    "print(integral)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "## User-defined model\n",
    "\n",
    "Now we'll see how you can define a custom model. To do that you need to subclass ``SpectralModel``.\n",
    "\n",
    "As an example we will use a PowerLaw plus a Gaussian (with fixed width)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [],
   "source": [
    "class UserModel(models.SpectralModel):\n",
    "    def __init__(self, index, amplitude, reference, mean, width):\n",
    "        self.parameters = ParameterList([\n",
    "                Parameter('index', index, parmin=0),\n",
    "                Parameter('amplitude', amplitude, parmin=0),\n",
    "                Parameter('reference', reference, frozen=True),\n",
    "                Parameter('mean', mean, parmin=0),\n",
    "                Parameter('width', width, parmin=0, frozen=True)\n",
    "            ])\n",
    "    def evaluate(self, energy, index, amplitude, reference, mean, width):\n",
    "        pwl = models.PowerLaw.evaluate(energy=energy, index=index, amplitude=amplitude, reference=reference)\n",
    "        gauss = amplitude * np.exp(-1 *(energy - mean) ** 2/( 2 * width ** 2))\n",
    "        return pwl + gauss"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "UserModel\n",
      "\n",
      "Parameters: \n",
      "\n",
      "\t   name     value   error       unit      min max frozen\n",
      "\t--------- --------- ----- --------------- --- --- ------\n",
      "\t    index 2.000e+00   nan                 nan nan  False\n",
      "\tamplitude 1.000e-12   nan 1 / (cm2 s TeV) nan nan  False\n",
      "\treference 1.000e+00   nan             TeV nan nan   True\n",
      "\t     mean 5.000e+00   nan             TeV nan nan  False\n",
      "\t    width 2.000e-01   nan             TeV nan nan   True\n"
     ]
    }
   ],
   "source": [
    "model = UserModel(\n",
    "    index=2 * u.Unit(''),\n",
    "    amplitude=1e-12 * u.Unit('cm-2 s-1 TeV-1'),\n",
    "    reference=1 * u.TeV,\n",
    "    mean=5 * u.TeV,\n",
    "    width=0.2 * u.TeV,\n",
    ")\n",
    "print(model)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<matplotlib.axes._subplots.AxesSubplot at 0x10fe31470>"
      ]
     },
     "execution_count": 14,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAZIAAAEOCAYAAACjJpHCAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMS4wLCBo\ndHRwOi8vbWF0cGxvdGxpYi5vcmcvpW3flQAAIABJREFUeJzt3Xd4VGXax/HvnUlPSAIBQu9dEIEI\nSIsrXSkusoK6ioooNsC8W/RVd33ddV230BRR7LoKKqKgIM2SUJSu0nsVAgklQEgCIc/7x0w0iwkz\nSWbmzJncn+vKJTkzc+aO6PzydDHGoJRSSpVXiNUFKKWUsjcNEqWUUhWiQaKUUqpCNEiUUkpViAaJ\nUkqpCtEgUUopVSEaJEoppSpEg0QppVSFaJAopZSqEA0SpZRSFRJqdQH+UL16ddOoUSOry1BKKVtZ\nt25dljGmhrvnVYogadSoEWvXrrW6DKWUshUR2e/J87RrSymlVIVokCillKoQDRKllFIVokGilFKq\nQjRIlFJKVYgtgkREmojIayIyu9i1G0XkFRGZKyL9rKxPKaUqM58HiYi8LiLHRGTTJdcHiMh2Edkl\nIo9e7h7GmD3GmNGXXPvEGDMGuBMY4fXClVI+senHbC5cLLS6DOVF/miRvAkMKH5BRBzANGAg0Aa4\nRUTaiEg7Efnskq+abu7/hOtePqFn2ivlPct2ZjLo+eW8uWKf1aUoL/J5kBhj0oETl1zuDOxytTTO\nA7OAocaYjcaYQZd8HSvpvuL0HPC5MWa9L2rPPneBm6avJH1Hpi9ur1SlcuFiIU/N2wzAR+sPWVyN\n8iarxkjqAgeLfX/Ida1EIpIoIi8BHUTkMdflh4E+wHARGVvCa+4VkbUisjYzs3xBcOxMHqdyL3DH\n66v5nw++59S58+W6j1IK3lq5j92ZOfRuVZNtGWfYeuS01SUpL7EqSKSEa6X2IRljjhtjxhpjmhpj\nnnVdm2qM6eS6/lIJr5lhjEk2xiTXqOF2q5gSNU+qwoJxPXnwV02Z+92P9JmYxmc/HNbuLqXK6NiZ\nPCYv3cmvWtbgH8OvxBEifPLdj1aXpbzEqiA5BNQv9n094LBFtVxWZJiD3/dvxbyHelAnIYqH3tvA\nmLfXkZGdZ3VpStnGc59vJ7/gIk8OakNibAQpLWowd8NhCgv1l7JgYFWQrAGai0hjEQkHRgLzLKrF\nI23qxDHn/m48fn1rlu/KpO/ENN5dtV//R1DKjYMnzvHR+kPc3aMxTWrEAvDrDnXJOJ3Ht3uPW1yd\n8gZ/TP+dCXwDtBSRQyIy2hhTADwELAK2Ah8YYzb7upaKCnWEMKZXExZN6EW7evE8/vEmRr7yLXsy\nz1pdmlIBa9OP2QAMalfnp2t9WicRGxHKx+u1eysY+GPW1i3GmNrGmDBjTD1jzGuu6wuMMS1c4x7P\n+OK9RWSwiMzIzs726n0bJsbw7j1deO6mdmw9cpoBU5Yx7atdOjdeqRJsP3oGEWhWM/ana1HhDga0\nrcXnmzLIu3DRwuqUN9hiZXt5GWM+NcbcGx8f7/V7iwgjrm7AF6kpXNeyJv9ctJ0hL6xg4yHvhpZS\ndrc94wyNEmOICnf81/VhHepyNr+AL7eVOMNf2UhQB4k/1IyL5KXbO/HSbztx/Gw+Q6ct59kFW8k9\nr79lKQXOFkmLpNhfXO/cuBqOEGHLYZ0GbHcaJF4yoG0tlqSmcHNyfV5O38OAKems3J1ldVlKWSrv\nwkX2ZeXQMqnKLx4LdYRQNyGKAyfOWVCZ8iYNEi+Kjwrj7zddyXtjugBw6yurePSjH8jOvWBxZUpZ\nY9exsxQaaFkrrsTHGyZGs1+DxPaCOkh8NdjuTrem1Vk4vhf39WrCB2sP0ndiGgs3Zfi1BqUCwY6j\nZwBoWeuXXVsA9atFc1CDxPaCOkh8OdjuTlS4g8eub83cB3uQGBvB2P+sY+w76zh2Whcyqspje8YZ\nwh0hNEyMKfHxBtWiOZFznjN52mq3s6AOkkDQrl488x7qzh8GtOTL7cfoMzGN99cc0G1WVKWw/egZ\nmtaMJcxR8kdNw2rRADpOYnMaJH4Q5gjhgWubsXB8T1rVjuOPH23ktldXsf94jtWlKeVT2zPO0LKE\nGVtF6ruCRLu37E2DxI+a1Ihl1piuPPPrtmw8lE3/yenMSN9NgS5kVEEoO/cCR7LzSh1oB2iQqC2S\nYBDUQWLVYPvlhIQIt3VpyJLUFHo0q8HfFmzj1y+u1Ln0KujsdDPQDhAXGUbV6DD2H9cgsbOgDhIr\nB9vdqRUfySt3dOKFWztwJDuXIS8s55+Ltul2ESpobMtwBkmLEtaQFNegWrS2SGwuqIMk0IkIg66s\nw5JHUhhyVR2mfbWb66cuY/XeSw+UVMp+dhw9Q2xEKHUToi77vPoaJLanQRIAqsaEM/Hmq3j77s6c\nLyjk5pe/4YlPNuqUSGVr2zKcW6OIlHSO3c8aJkbz48lcHSu0MQ2SANKrRQ0WTejF3d0b8+6qA/Sd\nmM7SLUetLkupMjPGsOPoGVrWuny3Fji7tgoKDUf0sDjb0iAJMDERofxpcBvm3N+N+Kgw7nl7LQ/P\n3EDW2XyrS1PKY8fO5HPq3IUS99i6VH1dS2J7GiQBqkODqnz6cA9S+7Zg0aYM+kxMY876Q7qQUdnC\n3iznGqmiExEvp4EGie0FdZAE4vTfsggPDWFc7+bMH9eDpjViSf3ge0a9sUYXb6mAl+HqpqrjZqAd\noHZ8FGEO0SCxsaAOkkCe/lsWzZOq8OF91/B/Q65g3b4T9J+czuvL93JRz4tXAapovKNWfKTb5zpC\nhHpVozmga0lsK6iDJJiEhAijujVicWoKnRtX4+nPtjD8pZU/7a6qVCDJyM6lSmQosRGhHj1fpwDb\nmwaJzdRNiOKNO69m8oir2JeVww1TlzFpyQ7yC3QhowocR7LzqO1Ba6RIQw0SW9MgsSER4cYOdVma\nmsL17Woz5YudDJq6nPUHTlpdmlKAM0hqxbsfHynSoFo02bkXyD6na6fsSIPExhJjI5gysgOv35lM\nTn4BN01fyVPzNpOTX2B1aaqSO5KdR50ytEh0CrC9aZAEgetaJbE4NYXbuzbkrW/20W9SOl9vP2Z1\nWaqSOl9QSNbZfI8G2os01F2AbU2DJEjERoTy9NC2fHjfNUSGhXDnG2tIff87TuSct7o0VckcdZ0C\nWpYxEm2R2FtQB4nd15GUR3KjaiwY35Nx1zVj3veH6Tsxjbnf/agLGZXfZJwumvrr+RhJbEQo0eEO\njusODrYU1EESLOtIyioi1EFqv5Z8Nq4H9apFM37Wd9zz1loOn8q1ujRVCRT9d1aWFglAtZhwjmsL\n2paCOkgqu1a14phzfzeeuKE1K3Zn0W9SOu98s49CXciofCijDIsRi0uMjdA95WxKgyTIOUKEe3o2\nYfGEFK6qn8CTczczYsY37Dp21urSVJA6kp1HbEQocZFhZXpd9ZhwHdOzKQ2SSqJBYjTvjO7MP4df\nyY6jZ7l+yjJe+HInF/QMCOVlGdl5ZW6NgKtr66wGiR1pkFQiIsJvkuuzJLUXfdsk8a/FOxj8/HJ+\nOHTK6tJUEDlyumyr2oskxkZwPCdfJ4bYkAZJJVSzSiTTbuvIjNs7cfLceW6ctoJn5m8h97xus6Iq\nLiM7l1pxZQ+S6rHhXLhoOJ2nC2rtRoOkEut3RS2WpKYwsnMDXlm2l/6T01mxK8vqspSNXbhYyLEz\n+eVqkVSLCQfQcRIb0iCp5OIiw/jbr9sx696uOEKE215dxR9mf697HqlyOXYmH2PKtoakSGJsBICu\nJbGhoA6Syrggsby6Nknk8/E9uf/apny0/kf6TErj841HtL9alUlGtmsNSUI5xkhcLZIsHXC3naAO\nksq6ILG8IsMc/HFAK+Y+2J2kuAjuf3c9972z7qctL5Ryp+hAq/INtmvXll0FdZCo8mlbN55PHujO\nowNbkbYjkz4T05i5+oC2TpRbRYsRa8eVvWuraIxEu7bsR4NElSjUEcLYlKYsnNCLK+rE8dicjdzy\nyrfsy8qxujQVwI5k5xEV5iAuyrOTEYuLCHVQJTJUt0mxIQ0SdVmNq8fw3j1deXZYOzYfPk3/yem8\nlLabAl3IqEpwJDuX2vGRiEi5Xl89NkKDxIY0SJRbISHCLZ0bsDQ1hZQWNfj759u48cUVbD6skxjU\nfztSzlXtRZyr27Vry24u2/4UkfUe3CPTGNPfS/WoAJYUF8nLt3di4aYMnpy7mSEvrODeXk0Y37s5\nkWEOq8tTASAjO49uTauX+/WJMeHsP65nktiNu47MCGDIZR4XYI73ylGBTkQY2K423ZpW55kFW5j+\n9W4Wbsrg2WHt6Nok0erylIUKKrAYsUhibATrD+iWPXbjrmvrQWPM7st87QLG+aNQFVjio8P4x/D2\nvHtPFwoKCxk541v+9+ONnM7ThYyVVdbZ81wsNBXq2kqMCedETr4edWAzlw0SY8zX7m7gyXNU8Ore\nrDqLJvRiTM/GzFp9gL4T01i8OcPqspQFjrgWI5Znn60iibHhFBo4lau/kNiJ28F2EblaRKaIyHoR\nOSIie0RknojcJyJV/FGkCmzR4aE8fkMbPn6gO1Wjw7n3nXU8+O56Ms/ooGllcsz1912hFolrm5QT\nOfrfjp1cNkhE5DPgISANuBFoDHQE/gokAPNFZJCvi1T20L5+Ap8+3IPf9WvBki1H6TMxjQ/XHtSF\njJVE0S8ONapElPseuk2KPblrkYw2xowyxswxxhwwxuQZY04ZY1YbY54zxvQCVvuj0PLQvbb8L8wR\nwkPXNWfB+J40rxnL72f/wO2vreaAzsQJepln8hH5eYV6eRRtk6IHXNmLuyB5TEQ6X+4JxphjXqzH\nq3SvLes0qxnLB/ddw19ubMt3B0/Rf3I6ry7bw0UdRA1amWfzqRYdTpij/MvTEmNcOwBr15atuPsb\nPwhME5HdIvKMiLT1R1EqOISECLd3bcjiR3rRrWkif52/lWEvrmBbxmmrS1M+kHkmn+qx5e/WAqga\n7TznXVsk9uJu1ta/jTFXA/2Ac8BMEdkkIv8rIk38UqGyvToJUbw6Kpmpt3Tg0MlcBk1dzr8Xbye/\nQE9kDCaZZ/IrND4Czj3eqkaHaYvEZjxqg7rWjDxjjGkHjAJ+A+z0aWUqqIgIQ9rXYUlqCoPb1+H5\nL3dx/ZRlrN13wurSlJd4I0jAdXa7tkhsxaMgERGHiAwUkbeA+cAeYIRPK1NBqVpMOJNGXMWbd11N\n3oVCfvPyN/xp7ibO5us53XZmjCHzrHeCpFpMuG7caDPupv/+SkRmAD/iXMH+JdDcGHOTMWa2PwpU\nwenaljVZ/EgvRl3TiHe+3U+/iWl8tT1g520oN07nFXC+oJAaFRwjAageqxs32o27FsnTwAagnTFm\noDHmLWPMGT/UpSqBmIhQnhpyBbPHdiMmIpS73ljDhFkb9EPEhryxhqRIYoxuJW837gbbexpjphtj\nMkWkq4jcASAiiSLSwD8lqmDXqWFVPhvXg/G9mzN/4xH6Tkrnkw0/6kJGG/FmkFSLCefUuQt65o2N\neDpG8gTwZ+AJ16VI4D1fFaUqn4hQB4/0bcH8cT1pUC2aCe9/x11vruHQSV3IaAeZZ70XJNWLzm4/\np60Su/B05dBw4HogB8AY8yMQ56uiVOXVIqkKH93fjT8NasPqvSfoNymdN1fs1d1gA1xWUYvEC2Mk\nRftt6cwt+/A0SPKNs5/BAIhItO9KUpWdI0S4u0djFk3oRXKjajz16RaGv7SSnUd1eC5QZZ7NJ8wh\nxEeFVfheRfttndBxEtvwNEjmiMg0IF5E7gIWA6/7riyloH61aN6662om3tyePVk53DB1OVOW7uR8\ngfadB5qiVe0hIeU7q724ov22snTShW14uiDxOeAzYB7QHnjGGDPZl4UpBc6FjMM61mNpagr929Zi\n0tIdDH5+ORsOnLS6NFWMtxYjQrH9trRryzbcrSNZXPRnY8znxphHjDETjDGf+740pX5WPTaC52/p\nwGujkjmdd4Fh01fy9KdbOHdeFzIGgswz+V4ZHwGIjwrDESLatWUj7lokNfxShVIe6t06icWP9OK2\nLg14fcVe+k1KJ31HptVlVXreWtUOzs0+46PCOKmztmwj1M3j8SIyrLQHjTFzvFyPUm5ViQzjrze2\nY+hVdfnjRz9wx+urualjPZ64oTVVK3AWhiqfi4WG414MEoCE6DA9btdG3AYJMAgoaQTNABokyjJX\nN6rGgnE9ef7Lnbyctoe0Hcd4asgV3NCuNiIVH/RVnjmRc55C4501JEUSosI4pS0S23AXJPuNMXf7\npRIfEJHBwOBmzZpZXYrykcgwB7/v34ob2tXh0Tk/8NB7G/ik9WH+emPbCp0drjyX6cU1JEWqRoeT\ncTrPa/dTvuVujMTWv9bpCYmVR5s6ccy5vxuPX9+a5bsy6TsxjXdX7deFjH7gzVXtReKjwzh1Tru2\n7MJdkNzulyqU8oJQRwhjejVh0YRetKsXz+Mfb2LkK9+yJ/Os1aUFNW/us1UkISpcu7ZsxN2mjZv8\nVYhS3tIwMYZ37+nCcze1Y+uR0wyYsowXv97FBd0E0CeKgqSix+wWVzU6jJzzF3XxqU14urJdKVsR\nEUZc3YAvUlPo3aom/1i4naEvrGDjoWyrSws6WWfziQ53EBPhbsjVcwmus9tP5WqrxA7KHCQiUlVE\nrvRFMUp5W824SKb/thMv/bYTmWfzGTptOc8u2ErueT0v3lu8uaq9SEK0cxp3to6T2IKn28h/LSJx\nIlIN+B54Q0Qm+rY0pbxnQNtaLE1N4ebk+rycvocBU9JZuTvL6rKCgjdXtRcpapGc1CCxBU9bJPHG\nmNPAMOANY0wnoI/vylLK++Kjwvj7TVfy3pguANz6yioe/egHsnXhW4V4c1V7kaquFokOuNuDp0ES\nKiK1gZtxbt6olG11a1qdheN7cV+vJnyw9iB9J6axcFOG1WXZli+6toq2o9cpwPbgaZA8DSwCdhlj\n1ohIE2Cn78pSyreiwh08dn1r5j7Yg8TYCMb+Zx33/2cdx87oIriyyC+4SHbuBZ91belguz14uo38\nh8aYK40xD7i+32OMucm3pSnle+3qxTPvoe78YUBLvth2jD7/TuP9NQf0vHgPZbm2evd2iyQ2IpTQ\nENExEpvQ6b+q0gtzhPDAtc34fHxPWtWK448fbeS2V1ex/3iO1aUFPF8sRgTn9O0EXd1uGxokSrk0\nrRHLrHu78syv27LxUDb9J6czI303BbqQsVTHXPtheTtIwDkFOFu7tmxBg0SpYkJChNu6NGRJago9\nmtXgbwu2MWz6SrYcPm11aQHpSLYzSGrHR3n93glRYZzM0RaJHbgNEhFpJSK9RST2kusDfFeWUtaq\nFR/JK3d04oVbO3D4VC5DXljOPxdtI++CLmQs7nB2LuGOEBJ9cA5MQnS4nkliE+6O2h0HzAUeBjaJ\nyNBiD//Nl4UpZTURYdCVdVjySApDr6rLtK92c/3UZazee8Lq0gLGkVN51IqPJCTE+xuFO8dItGvL\nDty1SMYAnYwxNwLXAk+KyHjXY7beYl4pT1WNCeffN7fn7bs7c76gkJtf/oYnPtnImTz9bflIdi61\nfXTui/NwK/13bAfugsRhjDkLYIzZhzNMBrq2R9EgUZVKrxY1WDShF6N7NObdVQfoNymdL7Yetbos\nSx0+lUedBO+Pj4AzwHMvXNTuRBtwFyQZInJV0TeuUBkEVAfa+bIwpQJRTEQoTw5qw5z7uxEXGcbo\nt9by8MwNZLkOd6pMLhYajp7O81mLpGh1u25hE/jcBckdwH/tHWGMKTDG3AH08llVSgW4Dg2q8unD\nPUjt24JFmzLoMzGNOesPVaqFjFln8ykoNNT2VYvEtd/WSR0nCXjuDrY6ZIzJgJ+3jxeRjiLSEcj1\nS4VKBajw0BDG9W7O/HE9aFojltQPvmfUG2s4eOKc1aX5xeFTzo+AOr4aI4nW/bbswqOTaETkL8Cd\nwG6g6FcuA1znm7KUso/mSVX48L5reOfb/fxj4Tb6T07nd/1aMqpbIxw+mM0UKHy5hgQ0SOzE0yPN\nbgaaGmO0jalUCUJChFHdGtGnTRKPf7yRpz/bwqc/HOa5m66kRVIVq8vziZ9aJAm+apHoVvJ24enK\n9k1Agi8LUSoY1E2I4o07r2byiKvYl5XDDVOXMWnJDvILgm/m0ZHsPKLCHD8Nintb1Z92ANYWSaDz\ntEXyLLBBRDYBP01PMcYM8UlVStmYiHBjh7r0bF6dpz/bwpQvdrJg4xGeG34lHRtUtbo8rzmSnUvt\nhEhEfNN9FxXmINwRooPtNuBpkLwFPAdsBPy6g53r7JPHcZ7SONx1rTUwHuc05C+MMdP9WZNSnkiM\njWDKyA4MvaoOT3y8iZumr2TUNY34ff+WxER4+r9e4Dp8Ko86PhofAWcgx0eH6bntNuBp11aWMWaq\nMeYrY0xa0Ze7F4nI6yJyzNWSKX59gIhsF5FdIvLo5e7hOvtk9CXXthpjxuIcu0n28GdQyhLXtUpi\ncWoKt3dtyFvf7KPfpHTSdmRaXVaF+XJVe5Gq0WHaIrEBT4NknYg8KyLXFE3/dU0BdudN4L82dxQR\nBzANGAi0AW4RkTYi0k5EPrvkq2ZpNxaRIcBy4AsPfwalLBMbEcrTQ9vy4X3XEBkWwqjXV5P6/nec\nzLHnh+SFi4UcO5PvszUkRRKiwnXWlg142r7u4Ppn12LX3E7/Ncaki0ijSy53xnlk7x4AEZkFDDXG\nPItz1bxHjDHzgHkiMh94z9PXKWWl5EbVmD+uJy9+tYsXv95N2o5M/jS4DUPa1/HZWIMvHD2dhzG+\nW0NSJCE6jAOVZF2OnXl61O6vSvgq7xqSusDBYt8fcl0rkYgkishLQAcRecx17VoRmSoiLwMLSnnd\nvSKyVkTWZmbavxtBBY/IMAep/Vry2bge1KsWzfhZ33HPW2t/mk5rBz+tIfF1i0S7tmzBoyARkb+J\nSEKx76uKyF/L+Z4l/dpV6r4SxpjjxpixxpimrlYLxpivjTHjjDH3GWOmlfK6GcaYZGNMco0aNcpZ\nqlK+06pWHHPu78aTg9qwcvdx+k1K551v9lFYGPjbrPh6VXuRqtHatWUHno6RDDTGnCr6xhhzEri+\nnO95CKhf7Pt6wOFy3kspW3OECKN7NGbxI73o0CCBJ+duZsSMb9h17KzVpV2Wv1ok8dFh5BcUkns+\n+NbhBBNPg8QhIj8dyiwiUUB5D2leAzQXkcYiEg6MBOaV815KBYX61aJ5++7O/Os37dlx9CzXT1nG\nC1/u5EKAnhd/5FQuVSJDifXxNOaEKNfqdj27PaB5GiT/Ab4QkdEicjewBOfakssSkZnAN0BLETkk\nIqONMQXAQ8AiYCvwgTFmc/nKVyp4iAjDO9VjaWoKfdsk8a/FOxj8/HJ+OHTK/Yv97HC2b9eQFCla\n3a5ntwc2j36dMMb8Q0R+APrgHOP4izFmkQevu6WU6wsoZZDcm0RkMDC4WbNmvn4rpbymRpUIpt3W\nkaGbM3hy7iZunLaC0T0ak9q3JVHhDqvLA35e1e5r8T9tk6ItkkDm7sz2nwbGjTELjTG/M8b8T/EQ\nkQCes2iM+dQYc298fLzVpShVZv2uqMWS1BRGXN2AV5btpf/kdFbsyrK6LMB5Vruvdv0trupPGzdq\niySQueva+kpEHhaRBsUviki4iFwnIm8Bo3xXnlKVW1xkGM8Oa8fMMV0JEbjt1VX8Yfb3lm4bknfh\nIsdzzvt8xhboVvJ24S5IBgAXgZkiclhEtojIXmAncAswyRjzpo9rVKrSu6ZpIgsn9GJsSlM+Wv8j\nfSal8fnGI5bUkuGnGVtQrEWiXVsBzd0JiXnGmBeNMd2BhkBvoIMxpqExZowx5ju/VKmUIjLMwaMD\nWzH3we7UrBLB/e+u57531nL0dJ5f6zic7Z81JOD8mSNCQ2y7lUxl4emsLYwxF4wxR4qvJwl0IjJY\nRGZkZ2dbXYpSXtO2bjxzH+zOowNb8fX2TPpMTGPm6gN+Oy9+51HnGpcGidF+eb8aVSLIOqtBEsg8\nDhI70sF2FaxCHSGMTWnKwgm9uKJOHI/N2cgtr3zLvqwcn7/3qr3HqRMfSV0/dG0BJMVF+r3Vpcom\nqINEqWDXuHoMM8d05e/D2rH58Gn6T07npbTdFPhoIaMxhtV7T9C1SaLfNplMiovQIAlw7qb/LhKR\nR0Sklb8KUkqVjYgwsnMDlqamcG3LGvz9823c+OIKNh/2fpfu7syzZJ09T5cm1bx+79LUrBLJsdP5\n7p+oLOOuRTIKOAk8JSLrRWS6iAwVkVg/1KaUKoOkuEhevj2Z6bd1JCM7nyEvrOC5hdvIu+C9faq+\n3XMCgC6NE712T3eS4iI5k19ATn6B395TlY27WVsZxpg3jTEjcZ5E+DbQCVgkIktF5A/+KLK8dLBd\nVUYD29Xmi9QUbupYl+lf72bglGV8u+e4V+797Z7j1IqLpKGfBtrB2bUFcOyMtkoCVVlmbRUaY74x\nxvzJNR14JPCj70qrOB1sV5VVfHQY/xjennfv6cLFQsPIGd/y2JyNnM4r/8I+Ywyr9p6gS5Nqfj2E\nKynOOc1Yx0kCV7kH240xWcaYd71ZjFLKu7o3q86iCb0Y07Mx7685QN+JaSzenFGue+3NyiHzTL5f\nu7Xg5xaJBkng0llbSgW5qHAHj9/Qho8f6E7V6HDufWcdD767nswydhUVjY909eNAO0BNV4tEB9wD\nlwaJUpVE+/oJfPpwD37XrwVLthylz8Q0Plh70OOFjKv2HqdGlQgaV4/xcaX/rUpEKFFhDm2RBLBy\nB4mI3OXNQpRSvhfmCOGh65qzYHxPWiZV4Q+zf+C3r61i//HLL2Q0xrBqzwm6NPbv+Ag4pzcnxUVw\nVAfbA1ZFWiT/57UqlFJ+1axmLLPu7cozv27LDwez6TcpnSlLd5Y6VXj/8XNknM6jSxP/jo8Uqamr\n2wPaZQ+2ch1mVeJDQJL3y/EuPdhKqdKFhAi3dWlI71ZJ/GX+FiYt3cHHGw7x5KA2/KplTUJCnC2P\nH0/l8sePnB8F3ZpaEyRJcZFsDMCTIpWTuxMSk4D+OBclFifASp9U5EXGmE+BT5OTk8dYXYtSgapW\nfCTTbu3IiORM/jxvM6PfWkuATijgAAAQgklEQVTDxGhuTq5P1ehwnl2wlUJj+MfwK2law5q1yElV\nIlh6Oh9jjN+71pR77oLkMyC2pO3iReRrn1SklLJErxY1WDihJ/N/OML7aw7yz0XbAejcqBr/vrk9\n9av5bxHipZLiIsm9cJEz+QXERYZZVocq2WWDxBgz+jKP3er9cpRSVooIdTCsYz2GdazH3qwc9mXl\n0KtFDRwh1rYCahatbj+dp0ESgHT6r1KqRI2rx/CrVjUtDxEovrpdZ24FIne7/653dwNPnqOUUhWh\n26QENndjJK0vM3MLnIPuupGVUsqnalYp2iZFWySByF2QeHIOiff2qFZKqRLERIRSJSJUWyQByt1g\n+35/FeILuo5EqeBRMy6CY2c0SAJRUA+26zbySgUP59nt2rUViII6SJRSwSNJt0kJWB4FiYi0KeHa\ntV6vRimlSlEzLoJjrtXtKrB42iL5QET+KE5RIvI88KwvC1NKqeKSqkRy/mIhp86V/5RH5RueBkkX\noD7O/bXWAIeB7r4qSimlLvXTWhIdcA84ngbJBSAXiAIigb3GmEKfVaWUUpf4+chdHXAPNJ4GyRqc\nQXI10AO4RURm+6wqpZS6hK5uD1zuFiQWGW2MWev6cwYwVERu91FNSin1CzWq/LxxowosngbJMRFp\ncMm1NG8X4226IFGp4BEZ5iAhOowfT+VaXYq6hKdBMh8wOPfWigQaA9uBK3xUl1fowVZKBZcO9RNY\nseu4HnAVYDwaIzHGtDPGXOn6Z3OgM7Dct6UppdR/6906iQMnzrE786zVpahiyrWy3RizHufAu1JK\n+U3v1jUBWLr1mMWVqOI86toSkdRi34YAHYFMn1SklFKlqB0fRZvacXy59RhjU5paXY5y8bRFUqXY\nVwTOMZOhvipKKaVK06d1TdbuP8HJnPNWl6JcPGqRGGP+z9eFKKWUJ65rncTUL3fx9Y5j/LpDPavL\nUbgJEhH5FOdsrRIZY4Z4vSKllLqMK+vGUz02gi+2apAECnctkn/5pQqllPJQSIhwXasafL4pgwsX\nCwlz6GkYVnMXJHuNMQf8UolSSnmod+skPlh7iDX7TtCtaXWry6n03EX5J0V/EJGPfFyLUkp5pEez\n6oQ7Qpi99pCeTxIA3AVJ8aWjTXxZiFJKeSomIpQ7uzdizoYfefbzbRomFnPXtWVK+bMt6F5bSgWv\nxwa2Iv/CRWak7+FioeGJG1rrtikWcRck7UXkNM6WSZTrz7i+N8aYOJ9WV0G615ZSwUtEeGrIFYgI\nry3fy+ncC/x5yBXERni6haDylsv+GzfGOPxViFJKlZWI8OfBbagSGcoLX+1i5e7j/HP4lXRrpgPw\n/qTz5pRStiYi/E+/lsweew3hoSHc+uoqHv94I2fy9Gx3f9EgUUoFhU4Nq7FgXE/u6dGYmasP0G9S\nOl9sPWp1WZWCBolSKmhEhTt4YlAb5jzQnbjIMEa/tZaHZ24g66ye8+5LGiRKqaBzVf0EPn24B4/0\nacHCTUfoMzGNOet1zYmvaJAopYJSeGgI4/s0Z8G4njSpHkPqB98z6o01HDxxzurSgo4GiVIqqDVP\nqsKHY7vx1OA2rNt3gv6T03l9+V4uFmrrxFs0SJRSQc8RItzZvTGLU1Po3LgaT3+2heEvrWTH0TNW\nlxYUNEiUUpVG3YQo3rjzaiaNaM++rBxumLqMyUt3cL6g0OrSbE2DRClVqYgIv+5Qj6WpKQxsW5vJ\nS3cy6PllrD9w0urSbEuDRClVKSXGRjD1lg68fmcyZ/MKuGn6Sp6at5mc/AKrS7MdDRKlVKV2Xask\nFqemcHvXhry5ch/9JqWTtiPT6rJsRYNEKVXpxUaE8vTQtsweew2RYSGMen01qe9/x8mc81aXZgsa\nJEop5ZLcqBrzx/Xk4euaMe/7w/SZmMa87w/rQkY3NEiUUqqYyDAH/9OvJZ8+3IN6VaMYN3MD97y1\nlsOncq0uLWBpkCilVAla145jzgPdeeKG1qzYnUW/Sem88+1+CnUh4y8EdZCIyGARmZGdnW11KUop\nG3KECPf0bMLiCSlcVT+BJz/ZxIgZ37Dr2FmrSwsoUhn6/pKTk83atWutLkMpZWPGGGavO8Rf528l\n9/xFxvVuxn0pTQlzBO/v4yKyzhiT7O55wftvQCmlvEhE+E1yfZak9qJvmyT+tXgHg59fzg+HTlld\nmuU0SJRSqgxqVolk2m0dmXF7J06eO8+N01bwzPwt5J6/aHVpltEgUUqpcuh3RS0WP5LCiKsb8Mqy\nvfSfnM6KXVlWl2UJDRKllCqn+Kgwnh3Wjln3dsURItz26ip+/+H3ZJ+rXOfFa5AopVQFdW2SyOfj\ne3L/tU2Zs+FHek9MY8HGI5VmIaMGiVJKeUFkmIM/DmjF3Ae7kxQXwQPvrue+d9Zx9HSe1aX5nAaJ\nUkp5Udu68cx9sDuPDmxF2o5M+kxMY9bqA0HdOtEgUUopLwt1hDA2pSkLJ/TiijpxPDpnI7e+sop9\nWTlWl+YTGiRKKeUjjavHMHNMV54d1o5Nh7PpPzmdl9J2U3AxuE5k1CBRSikfEhFu6dyApakpXNuy\nBn//fBs3vriCzYeDZ+smDRKllPKDpLhIXr49mem3dSQjO58hL6zguYXbyLtg/4WMGiRKKeVHA9vV\n5ovUFG7qWJfpX+9m4JRlrNpz3OqyKkSDRCml/Cw+Oox/DG/Pu/d04WKhYcSMb/nfjzdyOs+eCxk1\nSJRSyiLdm1Vn4YSejOnZmFmrD9B3YhpLthy1uqwy0yBRSikLRYeH8vgNbfj4ge5UjQ5nzNtrefC9\n9WSeybe6NI9pkCilVABoXz+BeQ/14Hf9WrBk81H6TExj9rpDtljIqEGilFIBIjw0hIeua86C8T1p\nkRTL7z78njteX83BE+esLu2yNEiUUirANKsZy/v3XsNfhl7B+v0n6TcpnVeX7eFigJ4Xr0GilFIB\nKCREuP2aRixJTeGapon8df5Whk1fybaM01aX9gsaJEopFcDqJETx2qhkpoy8ioMnzjFo6nImLt5O\nfkHgLGTUIFFKqQAnIgy9qi5LU1MY3L4OU7/cxQ1Tl7Nu/wmrSwM0SJRSyjaqxYQzacRVvHnX1eSe\nv8jwl77hT3M3cTa/wNK6NEiUUspmrm1Zk8WP9GLUNY1459v99JuYxlfbjllWT8AHiYg0EZHXRGT2\nJddjRGSdiAyyqjallLJKTEQoTw25gtljuxETEcpdb65h/KwNHD/r/4WMPg0SEXldRI6JyKZLrg8Q\nke0isktEHr3cPYwxe4wxo0t46I/AB96sVyml7KZTw6p8Nq4H43s3Z8HGI/SZmMYnG37060JGX7dI\n3gQGFL8gIg5gGjAQaAPcIiJtRKSdiHx2yVfNkm4qIn2ALYD9NqVRSikviwh18EjfFswf15OGiTFM\neP877npzDYdO+mcho0+DxBiTDlw6raAzsMvV0jgPzAKGGmM2GmMGXfJVWqffr4CuwK3AGBEJ+C46\npZTytRZJVfjo/m78aVAbVu05Qb9J6bz9zT6fv68VH8B1gYPFvj/kulYiEUkUkZeADiLyGIAx5nFj\nzATgPeAVY8wvzq0UkXtFZK2IrM3MzPTuT6CUUgHKESLc3aMxix/pRaeGVdmeccbn7xnq83f4JSnh\nWqmdecaY48DYUh578zKvmwHMAEhOTg7MfQWUUspH6leL5u27O3PeD+fDW9EiOQTUL/Z9PeCwBXUo\npVRQExEiQh0+fx8rgmQN0FxEGotIODASmGdBHUoppbzA19N/ZwLfAC1F5JCIjDbGFAAPAYuArcAH\nxpjNvqxDKaWU7/h0jMQYc0sp1xcAC3z53gAiMhgY3KxZM1+/lVJKVVpBPW3WGPOpMebe+Ph4q0tR\nSqmgFdRBopRSyvc0SJRSSlWIBolSSqkKsWJBot8UDbYDp0VkJxAPZJfjVtWBLG/WpkpV3r+jQBeo\nP5cVdfn6PX1xf2/cs6L3sOLzq6EnTxJ/7hBpNRGZYYy5txyvW2uMSfZFTeq/lffvKNAF6s9lRV2+\nfk9f3N8b96zoPQL586uydW19anUByq1g/TsK1J/Lirp8/Z6+uL837lnRewTqf0OVq0VSXtoiUUrZ\nlbZIAscMqwtQSqly8vnnl7ZIlFJKVYi2SJRSSlWIBolSSqkK0SBRSilVIRokZSQiMSLyloi8IiK3\nWV2PUkqVhYg0EZHXRGS2t+6pQQKIyOsickxENl1yfYCIbBeRXSLyqOvyMGC2MWYMMMTvxSql1CXK\n8hlmjNljjBntzffXIHF6ExhQ/IKIOIBpwECgDXCLiLTBeTTwQdfTLvqxRqWUKs2beP4Z5nUaJIAx\nJh04ccnlzsAuV3qfB2YBQ3GeOV/P9Rz996eUslwZP8O8Tj8IS1eXn1se4AyQusAc4CYRmU4Ab1mg\nlKr0SvwME5FEEXkJ6CAij3njjYJ6998KkhKuGWNMDnCXv4tRSqkyKu0z7Dgw1ptvpC2S0h0C6hf7\nvh5w2KJalFKqrPz2GaZBUro1QHMRaSwi4cBIYJ7FNSmllKf89hmmQQKIyEzgG6CliBwSkdHGmALg\nIWARsBX4wBiz2co6lVKqJFZ/hummjUoppSpEWyRKKaUqRINEKaVUhWiQKKWUqhANEqWUUhWiQaKU\nUqpCNEiUUkpViAaJqvRE5KKIfFfs61H3r/I9EdknIhtFJFlEPnbVtktEsovV2q2U194jIu9cci3J\ntdV4mIi8LyInRORG//w0KpjpOhJV6YnIWWNMrJfvGepaEFaRe+wDko0xWcWuXQv8zhgzyM1rqwI7\ngXrGmDzXtYeAdsaY+1zf/wfn2TqfVKROpbRFolQpXC2C/xOR9a6WQSvX9RjXQUJrRGSDiAx1Xb9T\nRD4UkU+BxSISIiIvishmEflMRBaIyHAR6S0iHxd7n74iMqcCdV4tImkisk5EPheRJGPMSWAlcEOx\np44EZpb3fZQqjQaJUhB1SdfWiGKPZRljOgLTgd+5rj0OfGmMuRr4FfBPEYlxPXYNMMoYcx3O0zQb\nAe2Ae1yPAXwJtBaRGq7v7wLeKE/hIhIBTAFuMsZ0Av4D/MX18Eyc4YGI1HfVkl6e91HqcnQbeaUg\n1xhzVSmPFbUU1uEMBoB+wBARKQqWSKCB689LjDFFBwz1AD40xhQCGSLyFTj38XaNX/xWRN7AGTB3\nlLP21sAVwFIRAXDg3PUVnBv0TRWRWGAEzr2WCsv5PkqVSoNEqcvLd/3zIj///yI4WwDbiz9RRLoA\nOcUvXea+b+A8GC0PZ9iUdzxFgB+MMT0vfcAYkyMiS3GeijcSuL+c76HUZWnXllJltwh4WFxNABHp\nUMrzluM8TTNERJKAa4seMMYcxnk2xBM4z9sury04T73r7KolXESuKPb4TOD3QIIxZk0F3kepUmmQ\nKPXLMZK/u3n+X4Aw4AcR2cTPYxKX+ghnN9Mm4GVgFZBd7PF3gYPGmC3lLdwYkw8MByaKyPfABqBL\nsacsxNntNqu876GUOzr9VykfEpFYY8xZEUkEVgPdjTEZrsdeADYYY14r5bX7uGT6r5dr0+m/yiu0\nRaKUb30mIt8By4C/FAuRdcCVOGdZlSYT+EJEkr1dlIi8D3THOUajVIVoi0QppVSFaItEKaVUhWiQ\nKKWUqhANEqWUUhWiQaKUUqpCNEiUUkpViAaJUkqpCvl/oPoPfo3B0P4AAAAASUVORK5CYII=\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x10fd19438>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "energy_range = [1, 10] * u.TeV\n",
    "\n",
    "fig, ax = plt.subplots(1)\n",
    "model.plot(ax=ax, energy_range=energy_range)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## What's next?\n",
    "\n",
    "In this tutorial we learnd how to work with spectral models. Go to [gammapy.spectrum](http://docs.gammapy.org/dev/spectrum/index.html) to see what else you can do with gammapy.spectrum."
   ]
  }
 ],
 "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.3"
  },
  "nbsphinx": {
   "orphan": true
  }
 },
 "nbformat": 4,
 "nbformat_minor": 1
}