{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "84841feb",
   "metadata": {},
   "source": [
    "# Source model interface\n",
    "\n",
    "`icecube_tools` has a simple source modelling interface built in that we demonstrate here."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "id": "7a3b57df",
   "metadata": {
    "execution": {
     "iopub.execute_input": "2023-10-27T13:42:27.853515Z",
     "iopub.status.busy": "2023-10-27T13:42:27.852855Z",
     "iopub.status.idle": "2023-10-27T13:42:28.812218Z",
     "shell.execute_reply": "2023-10-27T13:42:28.811193Z"
    }
   },
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "from matplotlib import pyplot as plt\n",
    "\n",
    "from icecube_tools.source.flux_model import PowerLawFlux, BrokenPowerLawFlux\n",
    "from icecube_tools.source.source_model import PointSource, DiffuseSource"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "48159555",
   "metadata": {},
   "source": [
    "## Spectral shape"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "5ce7493a",
   "metadata": {},
   "source": [
    "We start by defining a spectral shape, such as a power law or broken power law. Let's start with the definition of a simple power law flux."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "id": "b2885c16",
   "metadata": {
    "execution": {
     "iopub.execute_input": "2023-10-27T13:42:28.817042Z",
     "iopub.status.busy": "2023-10-27T13:42:28.816518Z",
     "iopub.status.idle": "2023-10-27T13:42:28.821262Z",
     "shell.execute_reply": "2023-10-27T13:42:28.820664Z"
    }
   },
   "outputs": [],
   "source": [
    "# Parameters of power law flux\n",
    "flux_norm = 1e-18  # Flux normalisation in units of GeV^-1 cm^-2 s^-1 (sr^-1)\n",
    "norm_energy = 1e5  # Energy of normalisation in units of GeV\n",
    "spectral_index = 2.0  # Assumed negative slope\n",
    "min_energy = 1e4  # GeV\n",
    "max_energy = 1e8  # GeV\n",
    "\n",
    "# Instantiate\n",
    "power_law = PowerLawFlux(flux_norm, norm_energy, spectral_index, min_energy, max_energy)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "id": "ebb8bb44",
   "metadata": {
    "execution": {
     "iopub.execute_input": "2023-10-27T13:42:28.824792Z",
     "iopub.status.busy": "2023-10-27T13:42:28.824208Z",
     "iopub.status.idle": "2023-10-27T13:42:29.710966Z",
     "shell.execute_reply": "2023-10-27T13:42:29.710256Z"
    }
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "Text(0, 0.5, 'F $[GeV^-1 cm^-2 s^-1 (sr^-1)]$')"
      ]
     },
     "execution_count": 3,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAksAAAG1CAYAAADpzbD2AAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjguMCwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy81sbWrAAAACXBIWXMAAA9hAAAPYQGoP6dpAABwF0lEQVR4nO3deVxU9f4/8NcsDDu4oCAK4oImsgwioia54b6XZWqFVt4WXIqsq91b6r2Vld1yYcqb5ZJb7mZmXa9cSXAJZXNBxQUFRUBEtgEGmDm/P/w13whBZmA4M/B6Ph7zSM45c+ZFB4b3fD6f8/lIBEEQQEREREQPJRU7ABEREZE5Y7FEREREVAcWS0RERER1YLFEREREVAcWS0RERER1YLFEREREVAcWS0RERER1YLFEREREVAe52AEsnU6nQ1ZWFhwdHSGRSMSOQ0RERPUgCAKKi4vh7u4OqbTutiMWSw2UlZUFDw8PsWMQERGRETIzM9GpU6c6j2Gx1ECOjo4AHvzPdnJyEjkNERER1UdRURE8PDz0f8frwmKpgX7venNycmKxRNVUVFRg1apVAIAFCxZAoVCInIiIiP6sPkNoJFxIt2GKiorg7OyMwsJCFktUjVqthoODAwCgpKQE9vb2IiciIqLfGfL3my1LRlKpVFCpVNBqtWJHITMll8sRHh6u/zcREVkmtiw1EFuWiIiILI8hf785zxIRERFRHVgsEREREdWBxRKRiajVarRq1QqtWrWCWq0WOw4RERmJo06JTKiwsFDsCERE1EAslohMxNbWFmlpafp/ExGRZWpR3XBTpkxB69atMXXq1Br70tPTMXToUPj4+MDPz4/dJtRgUqkU3t7e8Pb2fuS6Q0REZL5a1Dv4ggUL8N133z1036xZs/CPf/wDqamp+PXXX2Ftbd3E6YiIiMgctahiaciQIQ9dA+bChQuwsrJCaGgoAKBNmzacRJAarLKyUj95aWVlpdhxiIjISGZTLB07dgwTJkyAu7s7JBIJ9u/fX+MYlUoFLy8v2NjYICQkBPHx8Y3y2leuXIGDgwMmTJiAPn364KOPPmqU8zZUpVYndgRqgIqKCsydOxdz585FRUWF2HGIiMhIZtN8olarERAQgBdffBFPPvlkjf07duxAZGQk1q5di5CQEKxcuRKjRo3C5cuX0b59ewCAUqlEVVVVjecePnwY7u7utb52VVUVYmNjkZycjPbt22P06NEIDg7GiBEjahyr0Wig0Wj0XxcVFRnz7T5ScXklpnx5As8Ge+ClQV3qtdAfmReZTKYfHyeTyUROQ0RExjKbYmnMmDEYM2ZMrfs///xzzJkzB7NnzwYArF27Fj/99BPWr1+PRYsWAQCSk5ONeu2OHTuib9++8PDwAACMHTsWycnJDy2Wli9fjmXLlhn1OobYnXALV3NL8MFPF3Hq+j189nQAWtlx1XpLYmNjg127dokdg4iIGshsuuHqUlFRgYSEBISFhem3SaVShIWF4eTJkw0+f3BwMHJzc3H//n3odDocO3YMvXr1euixixcvRmFhof6RmZnZ4Nd/mFkDvfDPSb2hkElx5GIuxq6KRcLNfJO8FhEREdXOIoqlvLw8aLVauLq6Vtvu6uqK7Ozsep8nLCwMTz/9NA4dOoROnTrpCy25XI6PPvoITzzxBPz9/eHt7Y3x48c/9BzW1tZwcnLC5s2b0b9/fwwfPtz4b6wOEokEzw/wwt7XB8KrrR2yCsvxzL9PYe2v16DTce1jIiKipmI23XBN4ciRI7Xue1Q34J9FREQgIiJCv2qxqfh2dMaP8wbh3X3n8WNKFj7++RJOXb+Hfz0dgLYOnN7AnJWWlsLb2xvAg5sI7OzsRE5ERETGsIiWJRcXF8hkMuTk5FTbnpOTAzc3N1EyqVQq+Pj4IDg42OSv5WhjhdXPKrH8ST9Yy6WIuXwXY1fH4rfr90z+2mQ8QRCQlZWFrKwsCAJbA4mILJVFFEsKhQJBQUGIjo7Wb9PpdIiOjsaAAQNEyRQREYHU1FScPn26SV5PIpFgej9P7I94HF3b2SOnSIPp604h6n9X2C1npmxsbJCUlISkpCTY2NiIHYeIiIxkNt1wJSUluHr1qv7r9PR0JCcno02bNvD09ERkZCTCw8PRt29f9OvXDytXroRardbfHddS9OrghB/nDsJ7+89jb9JtfHY4Daeu5+OLaUq0c2S3nDmRyWRQKpVixyAiogaSCGbSPxATE4OhQ4fW2B4eHo6NGzcCAKKiorBixQpkZ2dDqVRi9erVCAkJaeKkD/w+M7NWq0VaWhoKCwvh5OTUpBl2ncnEez+cR3mlDu0crbFqmhIDu7s0aQYiIiJL9PuY4/r8/TabYslSGfI/2xSu5BQjYlsi0nJKIJEA84d5Y/5wb8iknMRSbJWVldi6dSsAYObMmbCyshI5ERER/Y7FUhMSu1gCgLIKLZYcOI+dZ24BAPp3bYPVzwaivRPHyYhJrVbDwcEBwINuZnt7e5ETERHR7wz5+20RA7zNUVPeDfcotgoZPp0agC+mBcBOIcOp6/kYsyoWx9Luih2tRZPJZBg7dizGjh3L5U6IiCwYW5YayBxalv7o2t0SRGxNxKXsYkgkwOtDuuHNsB6Qy1gXExER/Y4tSy1Yt3YO2B/xOGaGeEIQANXRa5ix7jfcKSwTOxoREZFFYrFkJHPqhvszGysZPpzihzXTA+FgLUf8jXyMXRWLo5dyxY5GRERkcdgN10Dm1g33Zzfy1Ji7PRHnbxcBAF55oisWjuoJK3bLmVxpaSkCAgIAACkpKVzuhIjIjLAbjvS8XOyx57WBmDXQCwDw72PXMe3fJ3G7gN1ypiYIAq5evYqrV69yuRMiIgvGYqkFsJbLsHRib6x9rg8cbeRIzCjA2FWx+G9qzqOfTEazsbFBXFwc4uLiuNwJEZEFYzeckcxhBm9jZOaXYu62RKTcKgQAvPh4Fywa8xgUctbNRETUcnBSyiZk7mOWHqaiSodPf7mEb+LSAQABnZwRNaMPPNpwTA0REbUMHLNEdVLIpfj7eB9880JfONtaIeVWIcaujsXP5+6IHa1Zqaqqwq5du7Br1y5UVVWJHYeIiIzElqUGssSWpT+6XVCGedsSkZhRAAAIH9AZi8f2go0VZ5xuKC53QkRkvtiyRPXWsZUtdrwyAK8M7goA2HTyJp766gRu5KlFTmb5pFIpBg8ejMGDB0Mq5a8aEZGlYsuSkSx1gHddjl7OxVs7U5CvroCDtRzLn/TDhAB3sWMRERE1Og7wbkKW3g33Z9mF5Zi/PQnxN/IBANP7eWLJBB92yxERUbPCbjgympuzDbbNCcG8Yd0hkQDb4zMwWXUc1+6WiB2NiIhIFCyWqAa5TIq3RvbEdy/2g4uDApeyizFhTRz2Jd0SO5pFKSsrg1KphFKpRFkZZ0wnIrJULJaoVqHe7XBofigGdG2L0got3tyRgrd3paC0grfB14dOp0NKSgpSUlKg0+nEjkNEREZisUR1au9kgy0vh+CNMG9IJMCuhFuYFHUcaTnFYkczezY2Njh8+DAOHz7M5U6IiCwYB3g3UHMb4F2XE9fy8Mb3ycgt1sDGSop/TPLF00GdIJFIxI5GRERkEA7wJpMY2M0FhxaEItTbBeWVOryz+ywid6ZArWG3HBERNV8sloykUqng4+OD4OBgsaM0KRcHa2ya3Q9vj+oJqQTYl3QbE9bE4eKdIrGjmZ2qqir89NNP+Omnn7jcCRGRBWM3XAO1pG64P4tPz8f87UnILiqHQi7Fkgk+mNHPk91y/x+XOyEiMl/shqMm0a9LGxxaEIqhPduhokqHv+07j3nbk1BcXil2NLMglUrRt29f9O3bl8udEBFZMLYsNVBLbln6nU4nYF3sdaz4z2VU6QR0bmsH1Yw+8O3oLHY0IiKih2LLEjUpqVSCVwZ3w45XBqBjK1vcvFeKJ788gU0nboC1OBERWToWS9Rogjq3xk/zB2GEjysqtDosOXABr29NRGEZu+WIiMhysViiRtXKToGvnw/C++N9YCWT4Ofz2Ri/JhYpmQViR2tyZWVlePzxx/H4449zuRMiIgvWooqlKVOmoHXr1pg6dWq17ZcvX9av4aVUKmFra4v9+/eLE7IZkEgkeHFQF+x+dSA82tgiM78MU9eewLdx6S2qW06n0+HEiRM4ceIElzshIrJgLWqAd0xMDIqLi7Fp0ybs3r37oceUlJTAy8sLN2/erNet3hzgXbfCskos2nMWP5/PBgCE9XLFZ0/7o5WdQuRkpldVVYWDBw8CAMaPHw+5XC5yIiIi+h0HeNdiyJAhcHR0rPOYAwcOYPjw4ZwTp5E421rhy5l98I9JvaGQSXHkYg7GrY5Dws37YkczOblcjsmTJ2Py5MkslIiILJjZFEvHjh3DhAkT4O7uDolE8tBuMJVKBS8vL9jY2CAkJATx8fGNnmPnzp2YNm1ao5+3JZNIJHhhgBf2vj4QXm3tcLugDNP+fRL//vUadLoW07BJREQWymyKJbVajYCAAKhUqofu37FjByIjI7FkyRIkJiYiICAAo0aNQm5urv4YpVIJX1/fGo+srKx6ZSgqKsKJEycwduzYWo/RaDQoKiqq9qD68e3ojB/nDcKEAHdU6QQs//kSXv7uDPLVFWJHMwmtVouYmBjExMRAq9WKHYeIiIxklmOWJBIJ9u3bh8mTJ+u3hYSEIDg4GFFRUQAeDJ718PDAvHnzsGjRonqfOyYmBlFRUQ8ds7R582b85z//wZYtW2p9/tKlS7Fs2bIa2zlmqf4EQcD2+Ews+/ECNFU6uDnZYM2MQAR7tRE7WqPicidEROar2Y1ZqqioQEJCAsLCwvTbpFIpwsLCcPLkyUZ7nfp0wS1evBiFhYX6R2ZmZqO9fkshkUgwI8QT+yMeR9d29sguKsezX5+C6ujVZtUtJ5FI4OPjAx8fH66XR0RkwSyiWMrLy4NWq4Wrq2u17a6ursjOzq73ecLCwvD000/j0KFD6NSpU7VCq7CwEPHx8Rg1alSd57C2toaTkxM2b96M/v37Y/jw4YZ9M6TXq4MTfpw7CE8GdoRWJ2DFfy4jfEM88ko0YkdrFHZ2drhw4QIuXLgAOzs7seMQEZGRLKJYaixHjhzB3bt3UVpailu3bmHAgAH6fc7OzsjJyYFCUb9b2iMiIpCamorTp0+bKm6LYG8tx7+eCcCnU/1hYyVF7JU8jF0Vi5PX7okdjYiICICFFEsuLi6QyWTIycmptj0nJwdubm6iZFKpVPDx8UFwcLAor9+cSCQSPNPXAz/OHQTv9g7ILdZg5jensPJIGrTNqFuOiIgsk0UUSwqFAkFBQYiOjtZv0+l0iI6OrtY61JTYstT4vF0dcWDuIDzTtxN0ArDyyBU8/+1vyC0uFzuaUcrKyjBixAiMGDGCy50QEVkws5kpr6SkBFevXtV/nZ6ejuTkZLRp0waenp6IjIxEeHg4+vbti379+mHlypVQq9WYPXu2iKmpsdkqZPh0agAGdGuLv+07jxPX7mHsqlisnBaIQd4uYscziE6nw5EjR/T/JiIiy2Q2xdKZM2cwdOhQ/deRkZEAgPDwcGzcuBHTpk3D3bt38f777yM7OxtKpRK//PJLjUHfTUWlUkGlUnH+HBOZEtgJ/p1aIWJrIi5lF+P59b8hYkh3vBHmDbnMIhpEYW1trZ+GwtraWuQ0RERkLLOcZ8mScG040yqv1OIfB1Ox7bcMAEA/rzZYPT0Qbs42IicjIiJL1uzmWaKWy8ZKho+m+GHN9EA4WMsRfyMfY1fH4ujl3Ec/mYiIqBGwWDIS74ZrWhMC3HFw3iD0dndCvroCszecxvKfL6JSa75jgbRaLU6fPo3Tp0+zu5aIyIKxG66B2A3XtMortfjo0EV8d/ImAKCPZyusmdEHHVvZipysJi53QkRkvtgNR82WjZUM/5jki69m9oGjjRyJGQUYuyoW/03NefSTm5hEIkHnzp3RuXNnLndCRGTB2LJkpD/eDZeWlsaWJRFk3CvFvO2JSLlVCAB4aVAX/HX0Y1DI+RmAiIjqZkjLEoulBmI3nLgqqnT4+OdLWH88HQAQ4NEKUdMD4dGGa7EREVHt2A1HLYZCLsX7E3zw9fNBcLKRIyWzAGNXx+KX83fEjkZERM0EiyVqFkb2dsOhBaEI9GyF4vIqvLolEUt+OA9NlXh3oZWXl2Py5MmYPHkyysstc8kWIiJiN5zROGbJPFVqdfjs8GX8+9frAADfjk6Imt4HXi5Nfyca74YjIjJfHLPUhDhmyTwdvZSLyJ3JuF9aCQdrOT5+yg/j/d2bNENlZSU2btwIAJg1axasrKya9PWJiKh2LJaaEIsl83WnsAzztyfh9I37AIAZIZ54f7wPbKxkIicjIiKxcYA3EYAOzrbYPqc/IoZ2g0QCbPstA5NVx3HtbonY0YiIyIKwWKJmTS6T4u1Rj+G7F/uhrb0Cl7KLMWFNHPYl3TL5a+t0Oly4cAEXLlyATme+y7IQEVHdWCwZiWvDWZZQ73b4eUEoBnRti9IKLd7ckYJ3dqegrMJ0d8uVlZXB19cXvr6+KCsrM9nrEBGRaXHMUgNxzJJl0eoErI6+gtX/uwJBAHq4OkA1ow+8XR0b/bXUajW8vLwAADdu3ODdcEREZoQDvJsQiyXLdOJqHhbsSMbdYg1srKT4xyRfPB3UiWu4ERG1EBzgTfQIA7u74ND8UIR6u6C8Uod3dp/FWztToNZUiR2NiIjMDIslarHaOVpj0+x+WDiyB6QSYG/SbUyIisPFO0ViRyMiIjPCYolaNKlUgrnDvPH9XwbA1cka1++qMVl1HNt+y0BDe6jLy8sxc+ZMzJw5k8udEBFZMI5ZaiCOWWo+7pVo8NauFMRcvgsAmBDgjo+m+MLRxriZt7ncCRGR+TLk77e8iTIRmb22DtZYHx6Mr2OvY8V/LuPHlCycu1WAqBl94NvR2eDzKRQKfPHFF/p/ExGRZWLLkpG4kG7zlnAzH/O2JSGrsBwKmRTvje+F5/p35t1yRETNBKcOaELshmu+CkorsHDXWRy5mAMAGOvnhuVP+sPZlgviEhFZOpMUSwcOHDA4yIgRI2Bra2vw8ywJi6XmTRAEfBuXjk9+uYRKrQCPNraImt4HAR6tHvlcnU6HjIwMAICnpyekUt5PQURkLkxSLBn6Ri+RSHDlyhV07drVoOdZGhZLLUNyZgHmbkvErftlsJJJsHhML8x+3KvObjkO8CYiMl8mm5QyOzsbOp2uXg87O7sGfRNE5kTp0Qo/zQ/F6N5uqNQK+MfBVPxlcwIKSivqfJ6dnR1/F4iILFy9i6Xw8HCDutSee+45trRQs+Jsa4WvnuuDf0zqDYVMiv+m5mDc6jgkZtx/6PH29vZQq9VQq9VsVSIismAc4N1A7IZrmc7fLkTEtkTcvFcKuVSCd0b3xMuDukIq5d1yRESWgGvD1WLKlClo3bo1pk6dWmPfF198gd69e8PHxwfz589v8OzN1Lz5dnTGwXmDMN6/A6p0Aj46dAkvf3cG+eq6u+WIiMjyGF0sVVZWIjMzE5cvX0Z+fn5jZjKZBQsW4Lvvvqux/e7du4iKikJCQgLOnTuHhIQEnDp1SoSEZEkcbaywZnogPpriB4Vciv9dysW41bE4fePB74NGo8GcOXMwZ84caDQakdMSEZGxDCqWiouL8dVXX2Hw4MFwcnKCl5cXevXqhXbt2qFz586YM2cOTp8+baqsDTZkyBA4Ojo+dF9VVRXKy8tRWVmJyspKtG/fvonTkSWSSCSYEeKJHyIeR9d29rhTWI5nvz4F1dGrqKioxDfffINvvvkGVVVVYkclIiIj1btY+vzzz+Hl5YUNGzYgLCwM+/fvR3JyMtLS0nDy5EksWbIEVVVVGDlyJEaPHo0rV64YFOTYsWOYMGEC3N3dIZFIsH///hrHqFQqeHl5wcbGBiEhIYiPjzfoNWrTrl07LFy4EJ6ennB3d0dYWBi6detm0DnUanW1rruKigqo1eoaLQq/D/jV6XT6bZWVlVCr1TUWWzXk2NLSUqjVami1Wv22qqoqqNVqlJWVGX1sWVkZ1Gp1tT/2Wq3W4GNLS0urHVteXg61Wo3KykqjjtXpdPr/P3+k0WigVqtRUVFh1LGCIOiPfdj1rO3Yx9wc8ePcQXgysCOqKivwyY8peHlTPBa/vxQffPABrKysHno9G+Pn5GHXszF+Tn6/ng39Ofnz9Wzoz0lt17OhPyd/vJ6GHGvI7z3fI1rue4Sx157vEXUf29D3iHoT6unZZ58Vzp8//8jjysrKhK+++kr49ttv63tqQRAE4dChQ8Lf/vY3Ye/evQIAYd++fdX2f//994JCoRDWr18vXLhwQZgzZ47QqlUrIScnR39MQECA0Lt37xqP27dv6485evSo8NRTT1U7d35+vjBy5Ejh3r17QmlpqTB48GDh119/fWjO8vJyobCwUP/IzMwUAAgAhNzcXP1xH3zwgQBAePnll6s9387OTgAgpKen67d98cUXAgBhxowZ1Y51cXERAFT7//71118LAIRJkyZVO7Zz584CACE+Pl6/bcuWLQIAISwsrNqxPj4+AgDh6NGj+m379u0TAAgDBw6sdmzfvn0FAMLBgwf12w4fPiwAEAICAqodO3jwYAGAsHPnTv22uLg4AYDQvXv3aseOHTtWACBs2LBBvy0pKUkAILi7u1c7durUqQIAISoqSr8tLS1NACA4OztXOzY8PFwAIHz66af6bbdu3RIACHK5vNqxr7/+ugBAWLJkiX7b/fv39dezoqJCv33hwoUCAGHhwoX6bRUVFfpj79+/LwiCIOh0OmHqnDcEAIJD4Dgh+IP/Cieu5gmCIAhyuVwAINy6dUt/jk8//VQAIISHh1fL5uzsLAAQ0tLS9NuioqIEAMLUqVOrHevu7i4AEJKSkvTbNmzYIAAQxo4dW+3Y7t27CwCEuLg4/badO3cKAITBgwdXOzYgIEAAIBw+fFi/7eDBgwIAoW/fvtWOHThwYI3f26NHjwoABB8fn2rHhoWFCQCELVu26LfFx8cLAITOnTtXO3bSpEkCAOHrr7/Wbzt//rwAQHBxcal27IwZMwQAwhdffKHflp6eLgAQ7Ozsqh378ssvCwCEDz74QL8tNzdXfz3/aMGCBQIA4d1339VvKykp0R9bUlKi3/7uu+8KAIQFCxZUOwffIx7ge8QDS5YsEQAIr7/+erXX43vEA031HlFYWCgAEAoLC4VHqfdCutu3b6/XcTY2Nnj11Vfre1q9MWPGYMyYMbXu//zzzzFnzhzMnj0bALB27Vr89NNPWL9+PRYtWgQASE5ONvh1AeDIkSPo3r072rRpAwAYN24cTp06hSeeeKLGscuXL8eyZcuMeh1q/iQSCXq7O2M3gFa2Vsgt1mDmN6ewYHgPsaMREZGRGnXqgMzMTCxZsgTr169v0HkkEgn27duHyZMnA3jQDGlnZ4fdu3frtwEP5n4qKCjADz/8UO9zx8TEICoqCrt379ZvO3XqFF577TWcPHkSVlZWmDhxIv7yl79g0qRJNZ6v0WiqNYcWFRXBw8MDWVlZcHNz08/oXFFRgcrKSsjlclhbW+uP/71Z0NbWVj8remVlJSoqKiCTyWBjY2PUsaWlpRAEATY2NpDJZAAeNIdqNBpIpdJqc2QZcmxZWRl0Oh2sra0hlz+orbVaLcrLyw06ViKRVJucsby8HFqtFgqFAlZWVgYfq9Pp9M28f5zDSKPRoKqqClZWVlAoFAYfKwiCvpnXzs6uxvU05NgKLfDernj8kHIHUlsnhHjY45Op/ujcvrX+ejbGz8nDrmdj/Jz8fj0b+nPy5+vZ0J+T2q5nQ39O/ng9G/pzUtv1NORYvkc0//eI+lx7vkeY7j1CtIV0U1JS0KdPn2r9l8b4c7GUlZWFjh074sSJExgwYID+uHfeeQe//vorfvvtt3qdNywsDCkpKVCr1WjTpg127dqlP9/f/vY37N27F1KpFMOHD8eqVavqXMpCpVJBpVJBq9UiLS2N8yxRDX9c7qTnX/ehHFZwcVBg5bRADPJ2ETkdEVHLZkixVO9uOODRi+lev37dkNM1uSNHjtS678MPP8SHH35Y73NFREQgIiJC/z+bqC47XxmAd364jEvZxXh+/W+IGNIdb4R5Qy5rUVOdERFZJIOKpcmTJ0MikdQ5YWNdrTHGcnFxgUwmQ05OTrXtOTk5cHNza/TXq48/tiwRPYy9vX2135X9ES5Y9mMqtsdnIOroVcSn52P19EC4OdvUcRYiIhKbQR9rO3TogL1799a6eG5iYqJJQioUCgQFBSE6Olq/TafTITo6ulq3XFOKiIhAamqqWc8rRebFxkqG5U/6YfX0QNgrZIi/kY+xq2Nx9HKu2NGIiKgOBhVLQUFBSEhIqHX/o1qd6lJSUoLk5GT9HW3p6elITk5GRkYGACAyMhLr1q3Dpk2bcPHiRbz22mtQq9X6u+OILMXEAHccnB+K3u5OyFdXYPaG01j+80VUanWPfjIRETU5gwZ4x8bGQq1WY/To0Q/dr1arcebMGQwePNjgIDExMRg6dGiN7eHh4di4cSMAICoqCitWrEB2djaUSiVWr16NkJAQg1+rMXCANz2KRqPBX//6VwDAJ598Uu1OFgAor9Tio0MX8d3JmwCAoM6tsXp6IDq2sq1xLiIialyi3Q3XEhnyP5talj/eDVdSUlLtFtY/+vncHbyz5yyKy6vgbGuFz54OwAgf16aMSkTU4pjsbjgiqj8rKyu8++67+n/XZoxfB/R2d8a87YlIuVWIOd+dwUuDuuCvox+DQs675YiIxMaWJSOxG44aW0WVDh//fAnrj6cDAAI6OSNqRh94tLF7xDOJiMhQ7IZrQuyGo8Z2+EI2Fu5KQVF5FRxt5Fgx1R+jfTuIHYuIqFkx5O+3QW382dnZDQpG1JIItaw4/igje7vh0IJQBHq2QnF5FV7dkoglP5yHpopzehERicGgYmnkyJGmymFxVCoVfHx8EBwcLHYUMlOlpaVwcHCAg4ODfn2o+urU2g47XxmAV57oCgDYdPImnvrqBG7kqU0RlYiI6mBQN5yfnx/OnTtnyjwWh91wVJv63g33KEcv5SJyZzLul1bCwVqOj5/yw3h/98aMSkTU4pisG84US5kQNVd2dnYoKSlBSUlJtVWyDTX0sfY4tCAUwV6tUaKpwtxtSXh33zmUV7JbjoioKfC+ZCITkUgksLe3h729fYM/aHRwtsX2Of0RMbQbJBJg228ZmKw6jmt3SxopLRER1YbFkpE4ZomamlwmxdujHsOm2f3Q1l6BS9nFmLAmDvuTbosdjYioWTNozFJgYCCSkpJMmcficMwS1aaiogLLli0DACxZsgQKhaLRzp1bVI753yfh1PV8AMC0vh5YOrE3bBWyRnsNIqLmjPMsNSEWS1SbxhrgXRutTsDq6CtY/b8rEASgh6sDVDP6wNvVsVFfh4ioOeJyJ0RmQC6XY8GCBfp/NzaZVII3R/RASJc2WLAjGWk5JZgYdRz/mNQbT/f1aPTXIyJqqQxuWdJqtfjxxx8xfPhwODryEyxblsgc3C3WIHJnMmKv5AEAnuzTEf+c5At7a34eIiJ6GJNNHQAAMpkM06dPx927d40O2BxwgDeZk3aO1tg0ux8WjuwBqQTYm3gbE6PicCm7SOxoREQWz6i74YKDg5Gent7YWSxKREQEUlNTcfr0abGjEAEApFIJ5g7zxvY5/eHqZI1rd9WYFHUc2+MzDFpuhYiIqjOqWJo3bx7effddZGZmNnYeomZDrVZDIpFAIpFArW66ZUpCurbFofmhGNyjHTRVOizeew4Lvk9GcXllk2UgImpOjLobTip9UGM5ODhg4sSJGDJkCAIDA+Hn59eot0dbAo5ZotqY+m64R9HpBHwdex0r/nMZWp0Ar7Z2iJrRB74dnZs0BxGROTL51AE3btzA2bNnkZycjJSUFCQnJ+PGjRuQy+Xo2bMnzp49a3R4S8NiiWojCALy8h4MuHZxcRFtuaCEm/mYty0JWYXlUMileG9cLzzXvzOXLyKiFk2UeZaKi4uRnJyMs2fPIiIiojFOaRFYLJEluK+uwNu7U3DkYi4AYKyfGz5+yh9ONlYiJyMiEodJ74YDgBdffBEbN27Uf33z5k3ExcXB39+/RRVKRJaitb0C617oi7+P6wW5VIJD57IxfnUczt4qEDsaEZHZM6pYOnToEB577DEAQEFBAYKCgjB58mT4+PggLS2tUQMSWaqKigp8+OGH+PDDD1FRUSF2HEgkErwc2hW7XxuITq1tkZFfiqe+OoH1cem8W46IqA5GFUuFhYXo2LEjAGDPnj1wc3NDUVERpk2bhkWLFjVqQHPFeZboUSorK/H3v/8df//731FZaT53oik9WuGn+aEY1dsVlVoB/ziYilc2J6Cw1HwyEhGZE6OKJQ8PD/08S7t27cKsWbNgbW2NV199FcePH2/UgOaK8yzRo8jlcrz88st4+eWXTbLcSUM421ph7XNBWDaxNxQyKQ6n5mDs6lgkZtwXOxoRkdkxaoD3Rx99hJ07d2LChAn4+OOPcfHiRXTv3h2XLl1CUFBQk84pIzYO8CZLd+5WIeZuT8TNe6WQSyX46+jH8NKgLpBKebccETVfJh/gvXjxYjz99NM4duwYPv74Y3Tv3h0AcPr0aXh6ehpzSiISiV8nZxycNwjj/DugSifgw0MX8fJ3Z3BfLf44KyIic9BoUwcAwIoVK1BeXo733nuvsU5p9tiyRM2FIAjYFp+BZT+moqJKhw7ONlg9PRDBXm3EjkZE1OhEmWeppWKxRLVRq9Vo3749ACA3N7fJZ/A2VmpWEeZuS8T1PDVkUgkiR/TAa4O7sVuOiJoVk3TDZWRkGBTi9u3bBh1P1ByVlpaitLRU7BgG8XF3woF5gzBZ6Q6tTsCK/1zGrI2nkVeiETsaEZEo6l0sBQcH45VXXqnz7q/CwkKsW7cOvr6+2LNnT6MEbExTpkxB69atMXXq1Br7PvvsM/Tu3Ru+vr7YsmWLCOmoubG1tUV6ejrS09Nha2srdhyDOFjL8cU0JT59yh82VlIcS7uLsaticfLaPbGjERE1uXp3w927dw8ffvgh1q9fDxsbGwQFBcHd3R02Nja4f/8+UlNTceHCBfTp0wfvvfcexo4da+rsBouJiUFxcTE2bdqE3bt367efO3cO4eHhOHHiBARBwNChQ/HLL7+gVatWjzwnu+GoubucXYyIbYm4mlsCqQRYMLwH5g7rDhm75YjIgpmkG65t27b4/PPPcefOHURFRcHb2xt5eXm4cuUKAGDmzJlISEjAyZMnzbJQAoAhQ4bA0dGxxvaLFy9iwIABsLGxga2tLQICAvDLL7+IkJDI/PR0c8SBuY9jalAn6ATgiyNpeGH9b8gtLhc7GhFRkzB46gBbW1tMnToVK1euxL59+/DLL79gy5YteOutt+Dr62t0kGPHjmHChAlwd3eHRCLB/v37axyjUqng5eUFGxsbhISEID4+3ujX+yNfX1/ExMSgoKAA9+/fR0xMDMdcUYNVVlZi5cqVWLlypVnN4G0MO4Ucnz0dgH89HQBbKxmOX72HsavicPxqntjRiIhMzmymFVar1QgICMCLL76IJ598ssb+HTt2IDIyEmvXrkVISAhWrlyJUaNG4fLly/o7jpRKJaqqqmo89/Dhw3B3d6/1tX18fDB//nwMGzYMzs7O6N+/P2Qy2UOP1Wg00Gj+b6BrUVGRod8qtRAVFRV48803AQBz5syBlZWVyIka7qmgTgjwcEbE1iRczinGc9/+hnlDu2P+cG/IZUZN20ZEZPYa9d0tMzMTL774olHPHTNmDD744ANMmTLlofs///xzzJkzB7Nnz4aPjw/Wrl0LOzs7rF+/Xn9McnIyzp8/X+NRV6H0u1deeQWJiYk4evQorKys4O3t/dDjli9fDmdnZ/3Dw8PDqO+Xmj+ZTIYZM2ZgxowZtRbflqh7e0f8MPdxTO/nAUEAVv/vKmZ88xtyitgtR0TNU6MWS/n5+di0aVNjnhLAg0/oCQkJCAsL02+TSqUICwvDyZMnG+U1cnNzAQCXL19GfHw8Ro0a9dDjFi9ejMLCQv0jMzOzUV6fmh8bGxts3boVW7duhY2NjdhxGpWNlQzLn/THqmeVsFfIEJ+ejzGrYhFzOVfsaEREjc6gbrgDBw7Uuf/69esNClObvLw8aLVauLq6Vtvu6uqKS5cu1fs8YWFhSElJgVqtRqdOnbBr1y4MGDAAADBp0iQUFhbC3t4eGzZsqHXhU2tra1hbW0OlUkGlUkGr1Rr/jRFZuEnKjvDr6Iy525KQeqcIszacxmtDuiFyRA9YsVuOiJoJg4qlyZMnQyKRoK7ZBiQS872d+MiRI7XuM7SFKiIiAhEREfpbD4laqq7tHLD39YH48KeL2HzqJr6KuYbT6flYPT0Q7q0sa34pIqKHMeijX4cOHbB3717odLqHPhITE00S0sXFBTKZDDk5OdW25+TkwM3NzSSv+SgqlQo+Pj4IDg4W5fXJ/KnVarRr1w7t2rWDWq0WO45J2VjJ8M/JvlDN6ANHaznO3LyPsatjEX0x59FPJiIycwYVS0FBQUhISKh1/6NanYylUCgQFBSE6Oho/TadTofo6Gh9N1pTi4iIQGpqap0zmhPl5eUhL6/l3F4/zr8DfpofCv9OzigorcRLm87gg4MPFuYlIrJUBnXDvf3223V+Qu7evTuOHj1qVJCSkhJcvXpV/3V6ejqSk5PRpk0beHp6IjIyEuHh4ejbty/69euHlStXQq1WY/bs2Ua9HpGp2dra4vz58/p/txSebe2w69UB+PjnS9hw/Aa+iUvH6Zv3ETU9EB5t7MSOR0RksHovd2JqMTExGDp0aI3t4eHh2LhxIwAgKioKK1asQHZ2NpRKJVavXo2QkJAmTvrAHwd4p6WlcbkToof4z4VsvL0rBUXlVXCykePTqQEY7StO1zkR0R8ZstyJ2RRLloprwxHV7db9UszdloTkzAIAwKyBXlg89jFYy5vP3FNEZHlMsjYcERmmsrIS69atw7p16yx+uZOG6NTaDjtfGYA5oV0AABtP3MBTX53AjbzmPeidiJoPg1qWsrOzRbv7zNywG44eRa1Ww8HBAcCDMXn29vYiJxJf9MUcvLUrBQWllXCwluPjp/ww3v/RM+wTETU2k3XD+fv74+zZsw0O2JywG45qU15ejmeffRYA8P333ze7WbyNlVVQhvnbk3Dm5n0AwMwQT7w33gc2VuyWI6KmY7Jiyc/PD+fOnWtwwOaExRKR4Sq1Onzx3zR8GXMNANCrgxNUMwLRtZ2DyMmIqKUw2Zglc56du6lxUkoi41nJpHhn9GPY9GI/tLFX4OKdIkxYE4cfkm+LHY2IqAZ2wzUQW5aIGianqBzztyfht/R8AMCzwR5YMqE3bBXsliMi0+HdcERmoLS0FF5eXvDy8kJpaanYccyWq5MNtr4cgvnDvSGRAN+fzsRk1XFczS0WOxoREQADiyWZjJ/0iOpLEATcvHkTN2/eNMkyQM2JXCZF5Ige2PJSCFwcrHE5pxgT1hzH7oRbYkcjIjKsWEpKSjJVDovDMUv0KDY2NoiPj0d8fDzvhKunx7u74NCCQXi8e1uUVWqxcFcK3tqZgtKKKrGjEVELxhm8G4hjlogan1Yn4MujV/HFkTToBKB7eweoZvRBTzdHsaMRUTPBMUtEZNFkUgnmDffGtjn94epkjau5JZgYFYfv4zPYpUlETY7FEpGJVFVVYevWrdi6dSuqqtiNZIz+Xdvi0PxQDO7RDpoqHRbtPYc3diSjRMP/n0TUdIzqhnvxxRfxxBNPYNasWQCAmzdvIjU1FQMHDoSzs3NjZzRr7Iaj2nC5k8aj0wn497Hr+OzwZWh1Arq62GPNjED0dm9Z7zdE1HhM3g136NAhPPbYYwCAgoICBAUFYfLkyfDx8cHly5eNOaXF4QBvehSpVIqwsDCEhYVBKmUjbkNIpRK8NqQbdvylPzo42+B6nhpTvjyBzad4pyERmZ5RLUu2trZIS0uDh4cHvv32W3zxxRdISEjA4sWLcePGDezdu9cUWc0SW5aImtZ9dQUW7kpB9KVcAMA4vw5Y/pQfnGysRE5GRJbE5C1LHh4eSE9PBwDs2rULs2bNgrW1NV599VUcP37cmFMSEdVLa3sFvgnvi7+P6wW5VIKfzt3B+NVxOHurQOxoRNRMGVUszZo1C/Pnz8d7772H6OhoTJ48GQCg0+lQUlLSmPmIiGqQSCR4ObQrdr06AB1b2SIjvxRPfXUCG46ns1uOiBqdUcXS4sWL8fTTT+PYsWP4+OOP0b17dwDA6dOn4enp2agBiSxVaWkpevfujd69e3O5ExMJ9GyNQ/NDMdLHFZVaAct+TMWrWxJQWFopdjQiakYadVLKFStWoLy8HO+9915jndLsccwS1YZ3wzUdQRCw8cQNfHToIiq1Ajq2skXUjEAEerYWOxoRmSlD/n5zBu8GYrFEtdFqtYiNjQUAhIaGcm3FJnD2VgHmbktCRn4p5FIJ/jr6Mbwc2gUSiUTsaERkZkw6wLusrAxxcXFITU2tsa+8vBzfffedoackapZkMhmGDBmCIUOGsFBqIv6dWuHg/EEY59cBVToBHx66iJc3ncF9dYXY0YjIghlULKWlpaFXr1544okn4Ofnh8GDB+POnTv6/YWFhZg9e3ajhzRHnGeJyDw52VghakYgPpjsC4VciuhLuRi3OhZnbuSLHY2ILJRBxdJf//pX+Pr6Ijc3F5cvX4ajoyMef/xxZGRkmCqf2YqIiEBqaipOnz4tdhQyU1VVVdi/fz/279/P5U6amEQiwXP9O2Pf6wPRxcUeWYXlmPb1KXwZcxU6HUceEJFhDBqz5OrqiiNHjsDPzw/Ag0GVr7/+Og4dOoSjR4/C3t4e7u7u0Gq1JgtsbjhmiWrDAd7moURThb/tO4cfkrMAAIN7tMPnzwSgrYO1yMmISEwmG7NUVlYGuVyu/1oikeCrr77ChAkTMHjwYKSlpRmXmKgZkkqlGDhwIAYOHMjlTkTkYC3HymlKfPKUH6zlUvyadhdjV8fit+v3xI5GRBZC/uhD/s9jjz2GM2fOoFevXtW2R0VFAQAmTpzYeMmILJytrS1ntDcTEokE04I9ofRojYhtibiaW4Lp607hzbAeeH1od8ikvFuOiGpn0MfdKVOmYPv27Q/dFxUVhenTp3P2XCIyWz3dHHFg7uN4qk8n6ATgX/9Nwwvrf8PdYo3Y0YjIjHGepQbimCUiy7Q74Rbe238eZZVauDhYY9WzSjze3UXsWETUREy+kK4lyszMxJAhQ+Dj4wN/f3/s2rWr2v6DBw+iZ8+e8Pb2xjfffCNSSmpOysrKEBwcjODgYJSVlYkdh/5kalAn/DjvcfR0dUReiQbPffsbPv9vGrS8W46I/qTFtCzduXMHOTk5UCqVyM7ORlBQENLS0mBvb4+qqir4+Pjg6NGjcHZ2RlBQEE6cOIG2bds+8rxsWaLa8G44y1BWocWyHy/g+9OZAICQLm2wenogXJ1sRE5GRKZkspal7OzsBgUTU4cOHaBUKgEAbm5ucHFxQX7+g0nq4uPj0bt3b3Ts2BEODg4YM2YMDh8+LGJaag6sra1x8OBBHDx4ENbWvE3dXNkqZPj4KX+selYJe4UMv6XnY+yqWPyadlfsaERkJgwqlkaOHGmqHDh27BgmTJgAd3d3SCQS7N+/v8YxKpUKXl5esLGxQUhICOLj4416rYSEBGi1Wnh4eAAAsrKy0LFjR/3+jh074vbt20adm+h3crkc48aNw7hx46pNuUHmaZKyI36cNwi9OjjhnroC4evj8ckvl1Cl1YkdjYhEZlCxZMoeO7VajYCAAKhUqofu37FjByIjI7FkyRIkJiYiICAAo0aNQm5urv4YpVIJX1/fGo+srCz9Mfn5+XjhhRfw9ddfG5VTo9GgqKio2oOImoeu7Ryw7/WBeK6/JwDgq5hrePbrU8gq4JgzopbMoI+7ply5e8yYMRgzZkyt+z///HPMmTNHv/bc2rVr8dNPP2H9+vVYtGgRACA5ObnO19BoNJg8eTIWLVqEgQMH6re7u7tXa0m6ffs2+vXr99BzLF++HMuWLavvt0UtmFarxf/+9z8AwLBhw7iYroWwsZLhg8l+6N+1LRbvOYczN+9j7OpYfP5MAIY95ip2PCISgUXcDVdRUYGEhASEhYXpt0mlUoSFheHkyZP1OocgCJg1axaGDRuG559/vtq+fv364fz587h9+zZKSkrw888/Y9SoUQ89z+LFi1FYWKh/ZGZmGv+NUbNWXl6OkSNHYuTIkSgvLxc7DhlovL87Ds4fBL+OzigorcSLG8/gw59SUcluOaIWxyKKpby8PGi1Wri6Vv9U5+rqWu9B58ePH8eOHTuwf/9+KJVKKJVKnDt3DsCDsSX/+te/MHToUCiVSrz11lu13glnbW0NJycnbN68Gf3798fw4cMb9s1RsyWVShEQEICAgAAud2KhOre1x+7XBmDWQC8AwLrYdDy99iQy80vFDUZETcqgbjhL7kYYNGgQdLraPxFOnDjRoOVaIiIiEBERob/1kOjPbG1tH9k1TObPWi7D0om9MaBbW7y9KwXJmQUYtzoWK54OwKjebmLHI6ImYNDH3aSkJFPlqJOLiwtkMhlycnKqbc/JyYGbmzhvViqVCj4+PggODhbl9YmoaY3q7Yaf5odC6dEKReVVeGVzApYeuABNlVbsaERkYhbRN6BQKBAUFITo6Gj9Np1Oh+joaAwYMECUTBEREUhNTcXp06dFeX0ianoebeyw85UBmBPaBQCw8cQNTP3qJG7eU4ucjIhMyWyKpZKSEiQnJ+u7LdLT05GcnIyMjAwAQGRkJNatW4dNmzbh4sWLeO2116BWq/V3xxGZm7KyMgwZMgRDhgzhcifNiEIuxd/G+eDb8L5oZWeFc7cLMX51HH46e0fsaERkImaz3ElMTAyGDh1aY3t4eDg2btwIAIiKisKKFSuQnZ0NpVKJ1atXIyQkpImTPqBSqaBSqaDVapGWlsblTqgGLnfS/GUVlGH+9iScuXkfAPBcf0/8fZwPbKwsd3wnUUthyHInRhVLRUVF2LBhA7Kzs9GlSxcEBATAz88PdnZ2Roe2VFwbjmpTVVWFffv2AQCmTJnCWbybqUqtDl/8Nw1fxlwDAPTq4ATVjEB0becgcjIiqovJi6WwsDCkpKQgODgYGRkZuHz5MgCgW7duCAgIwI4dO4xLboFYLBERAPyadhdv7khGvroC9goZPnrSD5OUHR/9RCIShSF/v436qHvy5EnExMTo7wTTaDQ4d+4ckpOTkZKSYswpLc4fu+GIiAb3aIdD80Mx//skxKfnY8H3yTh57R6WTuzNbjkiC2dUy9KAAQPw5ZdfIjAw0BSZLApblqg2Wq0Wp06dAgD079/foucpo/qr0uqwOvoK1hy9CkEAero6QjUzEN3bO4odjYj+wOTdcLGxsfj000+xe/duWFtbGx20OWCxRLXhAO+WLe5KHt7YkYy8Eg1srWT4YLIvngrqJHYsIvr/DPn7bdTUAV5eXigqKoKPjw/effddHDhwoMWtkcZJKelRJBIJunfvju7du5t0EWoyT4O8XXBowSAM7NYWZZVavLUrBQt3paC0okrsaERkIKNalvr164ecnBwMHjwYGRkZSElJQVFREdq0aYPAwEAcPnzYFFnNEluWiKguWp0A1dGrWHkkDToB6N7eAaoZfdDTjd1yRGIy+QDv8+fP4+TJkwgICNBvu3HjBpKSknD27FljTklE1CzJpBLMH+6Nfl3aYP72JFzNLcEkVRyWTeyNZ/p6sNWRyAIY1bI0ePBgLF++HAMHDjRFJovCliUiqq+8Eg0id6bgWNpdAMBkpTs+mOIHB2vOwUXU1Ew+ZmnBggVYunQpCgoKjHl6s8AxS/Qo5eXlGDduHMaNG4fy8nKx45AZcHGwxsZZwXhndE/IpBLsT87CxDVxSM0qEjsaEdXBqJYlqfRBjdW2bVtMmTIFISEhCAwMhK+vLxQKRaOHNGdsWaLa8G44qsuZG/mYtz0JdwrLoZBL8f54H8wM8WS3HFETMfnUATdv3kRKSop+Esrk5GTcuHEDcrkcPXv2bFHjllgsUW0qKyuxdetWAMDMmTNhZWUlciIyN/fVFVi4KwXRl3IBAOP8O+DjJ/3gaMOfFSJTM3mx9DDFxcVITk7G2bNnERER0RintAgsloioIQRBwDex6fjkl0uo0gno3NYOUdP7wK+Ts9jRiJo1kxdLJ06cgJOTE3x9fY0O2VywWCKixpCYcR/ztiXhdkEZFDIp/jauF14Y0JndckQmYvIB3hEREfjtt99qbL927RqKi4uNOaXF4QBvehStVovk5GQkJydzDUF6pD6erXFofihG+riiQqvDkgMX8NqWRBSWVYodjajFM6plyc7ODufOnUO3bt2qbf/3v/+NH3/8EQcPHmy0gOaOLUtUGw7wJmMIgoCNJ27go0MXUakV0Km1LaJm9IHSo5XY0YiaFZO3LDk5OeH+/fs1toeGhuoXDiVq6SQSCdzd3eHu7s6uFKo3iUSC2Y93wZ7XBsKzjR1u3S/D1K9O4JvY62ikIaZEZCCjiqXRo0fjs88+q3kyqRQVFRUNDkXUHNjZ2eH27du4ffs27OzsxI5DFsa/UyscnD8I4/w6oEon4IOfLmLOd2dQUMr3WKKmZlSx9M9//hO//vornnrqKZw7dw7Agwn4PvnkE/j7+zdqQCKilsrJxgpRMwLxz8m+UMilOHIxF2NXxSLhZr7Y0YhaFKOKJQ8PD5w6dQplZWUICAiAra0tHB0d8eOPP2LFihWNnZGIqMWSSCR4vn9n7Ht9ILq42COrsBzP/PsU1v56DTodu+WImkKD51nKyMhAcnIyrKysEBISgjZt2jRWNovAAd5Um/Lycjz//PMAgM2bN8PGxkbkRGTpSjRV+Nu+c/ghOQsAMKRnO/zr6QC0dbAWORmR5RFlUsqWisUS1YZ3w5EpCIKAHaczseTABWiqdHB1ssbqZwMR0rWt2NGILIrJ74arTWZmJl588cXGPKXZ4jxL9CgKhQJRUVGIiopqcWsmkulIJBI8288TP8x9HN3a2SOnSIPp605hTfQVaNktR2QSjdqylJKSgj59+rSoCfjYskREYimtqMJ7+y9gT+ItAMCg7i74YpoS7RzZLUf0KIb8/ZYbcuIDBw7Uuf/69euGnI6IiBrATiHHv54JwIBubfHe/vOIu5qHMatisfpZJQZ2dxE7HlGzYVDLklQqhUQiqXNiNIlEwpYlIgA6nQ7Xrl0DAHTr1g1SaaP2ehNVcyWnGBHbEpGWUwKJBJg3zBsLhntDJuWEqEQPY7IxSx06dMDevXuh0+ke+khMTGxQcKLmpKysDD169ECPHj1QVlYmdhxq5rxdHfFDxCA8G+wBQQBWR1/BzG9OIaeoXOxoRBbPoGIpKCgICQkJte5/VKsTUUvj7OwMZ2dnsWNQC2GrkOHjp/yx6lkl7BUynLqej7GrYnEs7a7Y0YgsmkHdcLGxsVCr1Rg9evRD96vVapw5cwaDBw9utIDmjt1wRGSOrt8tQcS2JFy8UwQAeH1IN0SO6AG5jN3BRIAJu+FCQ0NrLZQAwN7e3mwLpczMTAwZMgQ+Pj7w9/fHrl27qu2fMmUKWrdujalTp4qUkIio8XRt54B9rw/EzBBPAMCXMdcwfd0pZBWwS5jIUC1mUso7d+4gJycHSqUS2dnZCAoKQlpamn6iwJiYGBQXF2PTpk3YvXt3vc/LliUiMncHz2Zh0Z5zKNFUobWdFf71TACGPeYqdiwiUYk2KaU569ChA5RKJQDAzc0NLi4uyM//v8UohwwZAkdHR5HSUXOk0Wgwa9YszJo1CxqNRuw41IKN93fHT/MHwbejE+6XVuLFjWfw0aGLqNTqxI5GZBEMKpays7NNlQPHjh3DhAkT4O7uDolEgv3799c4RqVSwcvLCzY2NggJCUF8fLxRr5WQkACtVgsPD48GpiaqXVVVFTZt2oRNmzahqqpK7DjUwnVua489rw3ErIFeAICvj13HM/8+iVv3S8UNRmQBDCqWRo4caaocUKvVCAgIgEqleuj+HTt2IDIyEkuWLEFiYiICAgIwatQo5Obm6o9RKpXw9fWt8cjKytIfk5+fjxdeeAFff/21UTk1Gg2KioqqPYgexsrKCp9++ik+/fRTWFlZiR2HCNZyGZZO7I21z/WBo40cSRkFGLsqFocvmO6DMFFzYNCYJT8/P5w7d86UeQA8mIJg3759mDx5sn5bSEgIgoODERUVBeDBhH8eHh6YN28eFi1aVK/zajQajBgxAnPmzNGvBv9HMTExiIqKqnPM0tKlS7Fs2bIa2zlmiYgsSWZ+KeZuT0JKZgEAYPbjXlg8phcU8hYzOoNaOJONWZJIxJkJtqKiAgkJCQgLC9Nvk0qlCAsLw8mTJ+t1DkEQMGvWLAwbNuyhhVJ9LV68GIWFhfpHZmam0eciIhKLRxs77HplAF4e1AUAsOH4DUxdewIZ99gtR/RnFvERIi8vD1qtFq6u1e/ecHV1rfc4quPHj2PHjh3Yv38/lEollEpltVaysLAwPP300zh06BA6depUaxFmbW0NJycnbN68Gf3798fw4cON/8aoWdPpdLh9+zZu374NnY4Dacn8KORS/H28D755oS+cba1w9lYhxq2OxaFzd8SORmRWDFpI15INGjSozj9YR44cMeh8ERERiIiI0DfjEf1ZWVkZOnXqBAAoKSnRT1NBZG7CfFxxaEEo5m9PQsLN+3h9ayKe798ZfxvXCzZWMrHjEYnOoJYlmUycXxoXFxfIZDLk5ORU256TkwM3NzdRMqlUKvj4+CA4OFiU1yfLIJfLIZe3mM8kZME6trLF93/pj1cHdwMAbD51E099dQLpeWqRkxGJz6BiKSkpyVQ56qRQKBAUFITo6Gj9Np1Oh+joaAwYMECUTBEREUhNTcXp06dFeX0yf/b29qisrERlZSVblcgiWMmkWDTmMWycHYw29gpcyCrC+NWxOJCS9egnEzVjZjNmqaSkBMnJyUhOTgYApKenIzk5GRkZGQCAyMhIrFu3Dps2bcLFixfx2muvQa1WY/bs2SKmJiJqfob0bI9D80PRr0sbqCu0mL89CYv3nkN5pVbsaESiMGjqgPv376N169YmCRITE4OhQ4fW2B4eHo6NGzcCAKKiorBixQpkZ2dDqVRi9erVCAkJMUmeR1GpVFCpVNBqtUhLS+PUAUTU7FRpdVgVfQVRR69CEIDH3BwRNaMPurd3EDsaUYMZMnWAQcVSu3bt8M9//hOvvPKKaNMImBuuDUe10Wg0iIyMBAB8/vnnsLa2FjkRkXHiruThjR3JyCvRwE4hwweTffFkn05ixyJqEJMVS8uXL8fy5cvRtWtXrFmzBqGhoQ0Oa+lYLFFt1Go1HBwefALn3XBk6XKLy/HG98k4ce0eAODpoE5YNqk37BS8gYEsk8kmpVy8eDEuX76MwMBADB06FNOnT8ft27cbFNZS8W44ehQrKyssWbIES5Ys4XInZPHaO9pg80sheDOsB6QSYFfCLUyKOo60nGKxoxGZnEEtS3+UkJCAN954A0lJSVi0aBHefvvtFtnNwJYlImppTl67hwXfJyG3WAMbKyn+MdEXT/ftxOEZZFFM1rL0R0FBQYiNjcW3336Lb7/9Fr169cK+ffuMPR0REVmIAd3a4tCCUIR6u6C8Uod39pxF5M4UqDVVYkcjMokGTx0wbdo0XLp0CS+99BLCw8MxYsSIxshl9tgNR48iCAIKCgpQUFAAIxtwicyWi4M1Ns3uh7dH9YRMKsG+pNuYsCYOF+8UiR2NqNEZ3Q1XUVGBS5cu4fz58/rHb7/9pl/HraVgNxzVhgO8qaU4fSMf87YlIbuoHAq5FEsm+GBGP092y5FZM+Tvt0G3MSxbtkxfGF27dg1VVVVwdnaGr68v/P39MXbsWPj7+zcoPBERWZZgrzY4tCAUC3el4H+XcvG3fedx8to9LH/SD442vLmBLJ9BLUu+vr7w8/ODv7+//r+enp6mzGf22LJEtREEAVVVD8ZwyOVyfsqmZk+nE/BN3HV8+stlVOkEdG5rB9WMPvDtyMXGyfyYbJ4l+j+cwZuI6OESM+5j3rYk3C4og0Imxd/G9cILAzrzAwOZFRZLTYgtS0RENRWUVuDt3Wfx39QcAMDo3m74ZKo/nG3ZLUfmweRTBxQWFuIvf/kLunfvjl69euHOnTtGBSVqzioqKvD222/j7bffRkVFhdhxiJpUKzsFvn4+CO+P94GVTIJfLmRj3OpYJGcWiB2NyGBGFUsRERE4d+4cPv30U9y8eRNlZWUAgDfffBNRUVGNGpDIUlVWVuKzzz7DZ599hsrKSrHjEDU5iUSCFwd1we5XB8KjjS1u3S/D02tP4JvY65xOgyyKUcXSzz//jC+//BJPPvkkZDKZfvuoUaOwadOmRgtHZMmsrKywcOFCLFy4kMudUIsW4NEKB+eFYoyvGyq1Aj746SLmfHcGBaVscSXLYFSxJAgCHB0da2z39vbGlStXGhzKEnBSSnoUhUKBFStWYMWKFVAoFGLHIRKVs60VvpzZB/+c1BsKmRRHLuZi7KpYJNzMFzsa0SMZVSyNGTMGW7durbFdrVa3mLsdIiIikJqaitOnT4sdhYjIIkgkEjw/wAt7Xx8Ir7Z2yCosxzP/PoW1v16DTsduOTJfBk1K+bvly5ejb9++AB60MkkkEpSXl+Of//wn+vTp06gBiSwV51kiejjfjs74cd4gvLvvPH5MycLHP1/Cqev38K+nA9DWoeUtyE7mz6iWJU9PT5w4cQInTpxAaWkp+vXrh1atWuHXX3/FJ5980tgZiSxSaWkpFAoFFAoFSktLxY5DZFYcbayw+lkllj/pB2u5FDGX72Ls6ljEp7NbjsxPg+dZysjIQEpKCqysrBASEoLWrVs3VjaLwHmWqDZcG46ofi7eKULEtkRcv6uGVAJEjuiB14d0h1TK1lgyHU5K2YRYLFFtBEFAYWEhAMDZ2ZndcER1UGuq8N7+89ibdBsAEOrtgs+fUaKdI7vlyDRMNinlnj17oFQq9V8vWrQI69evR0JCAjQajVFhiZoriUSCVq1aoVWrViyUiB7B3lqOz6cpsWKqP2yspIi9koexq2Nx4mqe2NGIDGtZGj9+PMLCwvDGG28AABwdHaHValFeXg6ZTIZevXrh2LFjaNWqlYnimh+2LBERNa4rOcWI2JaItJwSSCTA/GHemD/cGzJ2y1EjMlnL0oULFzBy5Mhq286dO4fr169j7969sLKywtq1aw1PbIE4zxI9SkVFBZYuXYqlS5dyuRMiA3i7OuKHiEGY1tcDggCsir6C5775DblF5WJHoxbKoJYlGxsbXLt2DR07dgQAtGrVComJiejatSsA4Pvvv8eaNWtw/Phx06Q1Q2xZotpwgDdRw+1Puo13951DaYUWbe0V+GKaEk/0aCd2LGoGDPn7bdA8Sy4uLrhx44a+WMrOzq62jINSqURqaqoRkYmaH7lcjtdff13/byIy3OTAjvDr5IyIrYm4lF2M8A3xeH1IN7wZ1gNymVGz3xAZzKCWpRdeeAFyuRzr169/6P60tDT06dMHJSUljRbQ3LFliYjI9MortfjnwVRs/S0DANDPqw1WTVeig7OtyMnIUplszNLbb7+Nbdu2YdWqVQ/df/z4cX2XHBERUWOxsZLhwyl+WDM9EA7WcsTfyMfYVbE4eilX7GjUAhhULPn5+WHLli14++23ERYWhj179iAjIwNZWVnYuXMnFi9ejJkzZ5oqKxERtXATAtxxcN4g+HZ0wv3SSszeeBrLD11EpVYndjRqxoyalDIpKQlvvvkmjh07pp8/RhAETJgwAbt37642jqm5Yzcc1UatVuun0SgoKOAAb6JGpKnSYvmhS9h44gYAINCzFdZMD0Sn1nbiBiOLYbJuuN8FBgYiJiYGN27cwIEDB7B161akpKTghx9+MNtCKTMzE0OGDIGPjw/8/f2xa9eueu0jaoiqqir9YrpE1His5TIsndgba5/rA0cbOZIyCjBudRwOX8gWOxo1Q/VuWTp79ix8fX0hldavvrpw4QJ69uxpNncB3blzBzk5OVAqlcjOzkZQUBDS0tJgb29f575HYcsS1Uan0+HOnTsAgA4dOtT7d4eIDJOZX4q52xKRcuvB8kIvPt4Fi8Y8BoWcv3NUO5O0LAUGBuLevXv1DjFgwABkZGTU+3hT69Chg36pFjc3N7i4uCA/P/+R+4iMJZVK0bFjR3Ts2JGFEpEJebSxw65XB+KlQV0AAOuPp+PptSeQmV8qcjJqLurd7CMIAt577z3Y2dWvP9jQGYuPHTuGFStWICEhAXfu3MG+ffswefLkaseoVCqsWLEC2dnZCAgIwJo1a9CvXz+DXgcAEhISoNVq4eHhYdA+IiIyTwq5FO+N90H/rm2xcFcKUm4VYuzqWKyY6o/Rvh3EjkcWrt7F0hNPPIHLly/X+8QDBgyArW39579Qq9UICAjAiy++iCeffLLG/h07diAyMhJr165FSEgIVq5ciVGjRuHy5cto3749gAeTYj5sfMjhw4fh7u4OAMjPz8cLL7yAdevW1Tiurn2/02g01RYNLioqqvf3SC1LRUWFfpqNBQsWQKFQiJyIqPkb4eOKQwtCMW9bIhIzCvDqlkSED+iMd8f1grVcJnY8slBG3Q1nahKJpEbLUkhICIKDgxEVFQXgwXgQDw8PzJs3D4sWLarXeTUaDUaMGIE5c+bg+eefr/e+P1q6dCmWLVtWYzvHLNGfcbkTIvFUanX47PBl/PvX6wAA345OiJreB14u/D2kB0x+N1xTq6ioQEJCAsLCwvTbpFIpwsLCcPLkyXqdQxAEzJo1C8OGDatRDNW1788WL16MwsJC/SMzM9Pwb4haBLlcjvDwcISHh5vNjQ5ELYWVTIrFY3phw+xgtLFX4PztIoxfE4cfU7LEjkYWyCKKpby8PGi1Wri6ulbb7urqiuzs+t0mevz4cezYsQP79++HUqmEUqnEuXPnHrnvz6ytreHk5ITNmzejf//+GD58eMO+OWq2rK2tsXHjRmzcuBHW1tZixyFqkYb2bI9D80PRr0sblGiqMG97Et7ddw7llVqxo5EFaTEfdwcNGgSd7uEzvNa1rzYRERGIiIjQN+MREZF5cnO2wbaXQ7Aq+gqijl7Ftt8ykHjzPlQz+6BbOwex45EFsIiWJRcXF8hkMuTk5FTbnpOTAzc3N1EyqVQq+Pj4IDg4WJTXJyKi+pPLpHhrZE9892I/uDgocCm7GBPWxGFf0i2xo5EFMKhYun79OsQYD65QKBAUFITo6Gj9Np1Oh+joaAwYMKDJ8wAPWpZSU1Nx+vRpUV6fzN/vy520atUKarVa7DhEBCDUux0OzQ/FgK5tUVqhxZs7UvDO7hSUVbBbjmpnULHk7e2Nu3fv6r+eNm1ajdYeY5WUlCA5ORnJyckAgPT0dCQnJ+sntoyMjMS6deuwadMmXLx4Ea+99hrUajVmz57dKK9PZAq/3whAROajvZMNtrwcgjfDekAqAXaeuYWJUXFIyykWOxqZKYOmDpBKpcjOztbPa+To6IiUlBR07dq1wUFiYmIwdOjQGtvDw8OxceNGAEBUVJR+UkqlUonVq1cjJCSkwa9tDJVKBZVKBa1Wi7S0NE4dQDXodDpcu3YNANCtWzfO4k1khk5eu4f53yfhbrEGNlZS/GOSL54O6qRfJJ6aL0OmDjCbYslScW04IiLLlleiwZs7khF7JQ8AMCWwIz6Y7At76xZzD1SLZLJ5liQSSY1qm9U3ERFZMhcHa2ya3Q9vj+oJqQTYl3QbE6LicPEOV2igBwxuWRozZox+zpgff/wRw4YNqzEz8d69exs3pRliNxw9SmVlJb7++msAwF/+8hdYWVmJnIiIHiU+PR/ztychu6gcCrkUSyf0xvR+HmwYaIZM1g1X38HUGzZsqO8pLR674ag2XO6EyDLlqyvw1s5kHL384IamCQHu+GiKLxxt+IGnOTFZsUQ1sVii2pSXl+uXz9m8eTNsbGxETkRE9aXTCVgXex0r/nMZVToBXm3tEDWjD3w7chLi5oLFUhNgNxwRUfOXcPM+5m9Pwu2CMihkUvxtXC+8MKAzu+WaARZLTYgtS0REzVtBaQUW7jqLIxcfzCs4xtcNHz/lD2dbdstZMpPdDUdERNTStLJTYN0LQXhvvA+sZBL8fD4b49fEIiWzQOxo1ERYLBGZSGlpKTp27IiOHTuitLRU7DhE1AASiQQvDeqC3a8OhEcbW2Tml2Hq2hP4Ni5dlGXAqGmxWDISF9KlRxEEAVlZWcjKyuKbKVEzEeDRCgfnhWKMrxsqtQL+eTAVc75LQEFphdjRyIQ4ZqmBOGaJaqPVanHu3DkAgJ+fH2QymciJiKixCIKAzadu4oODF1Gh1aFjK1usnh6IoM6txY5G9cQB3k2IxRIRUct1/nYh5m5LxI17pZBLJXh7VE/MCe0KqZR3y5k7DvAmIiJqAr4dnfHjvEGYEOCOKp2A5T9fwkubTiNfzW655oTFEpGJVFZWYuPGjdi4cSMqKyvFjkNEJuJoY4XVzyrx0RQ/WMulOHr5LsauikV8er7Y0aiRsBvOSJyUkh6Fy50QtTwX7xQhYlsirt9VQyaVIHJED7w2uBu75cwQxyw1IY5ZotqUl5fjqaeeAgDs2bOHy50QtRBqTRXe238ee5NuAwBCvV3w+TNKtHO0FjkZ/RGLpSbEYomIiP5MEATsSriF9384j/JKHdo5WmPVs0oM7OYidjT6/zjAm4iISEQSiQTP9PXAj3MHwbu9A+4Wa/DcN79h5ZE0aHVso7A0LJaIiIhMxNvVEQfmDsIzfTtBJwArj1zBc9/8htyicrGjkQFYLBGZSGlpKby9veHt7c3lTohaMFuFDJ9ODcAX0wJgp5Dh5PV7GLs6FrFX7oodjeqJxRKRiQiCgKtXr+Lq1atc7oSIMCWwEw7MHYTH3ByRV1KBF9bHY8V/LqFKqxM7Gj0CiyUiE7GxsUFcXBzi4uJ4JxwRAQC6t3fA/ojHMSPEE4IAqI5ew4x1v+FOYZnY0agOvBvOSJxniYiIGuLHlCws3nsOJZoqtLazwufPKDH0sfZix2oxOHVAE+LUAUREZKwbeWpEbEvEhawiAMArT3TFwlE9YSVjx4+pceoAIjNQVVWFXbt2YdeuXaiqqhI7DhGZIS8Xe+x5bSDCB3QGAPz72HVM+/dJ3C5gt5w5YctSA7FliWrD5U6IyBA/n7uDd/acRXF5FZxtrfDZ0wEY4eMqdqxmiy1LRGZAKpVi8ODBGDx4MKRS/qoRUd3G+HXAofmhCOjkjMKySsz57gz+eTAVFVW8W05sbFlqILYsERFRY6qo0uGTXy7h27h0AECARytETQ+ERxs7kZM1L2xZIiIislAKuRTvjffBuhf6wtnWCimZBRi7Oha/nL8jdrQWq8UUS5mZmRgyZAh8fHzg7++PXbt26fcVFBSgb9++UCqV8PX1xbp160RMSkREBIzwccWhBaHo49kKxeVVeHVLIpb8cB6aKq3Y0VqcFtMNd+fOHeTk5ECpVCI7OxtBQUFIS0uDvb09tFotNBoN7OzsoFar4evrizNnzqBt27aPPC+74ag2ZWVlGDBgAADg5MmTsLW1FTkREVmiSq0Onx2+jH//eh0A4NvRCVHT+8DLhTeNNAS74R6iQ4cOUCqVAAA3Nze4uLggPz8fACCTyWBn96AvWKPRQBAELk9BDabT6ZCSkoKUlBTodBygSUTGsZJJsXhML2yYFYzWdlY4f7sI49fE4eDZLLGjtRhmUywdO3YMEyZMgLu7OyQSCfbv31/jGJVKBS8vL9jY2CAkJATx8fFGvVZCQgK0Wi08PDz02woKChAQEIBOnTrh7bffhouLi7HfChGAB8udHD58GIcPH+ZyJ0TUYEMfa49DC0IR7NUaJZoqzN2WhHf3nUN5JbvlTM1siiW1Wo2AgACoVKqH7t+xYwciIyOxZMkSJCYmIiAgAKNGjUJubq7+mN/HHP35kZX1f9V3fn4+XnjhBXz99dfVzt+qVSukpKQgPT0d27ZtQ05OzkNzaDQaFBUVVXsQPYxMJsOIESMwYsQIyGQyseMQUTPQwdkW2+f0R8TQbpBIgG2/ZWCy6jiu3S0RO1qzZpZjliQSCfbt24fJkyfrt4WEhCA4OBhRUVEAHnRxeHh4YN68eVi0aFG9zqvRaDBixAjMmTMHzz//fK3Hvf766xg2bBimTp1aY9/SpUuxbNmyGts5ZomIiJrSsbS7eHNHMu6pK2CnkOHDKb6YEthJ7FgWo9mNWaqoqEBCQgLCwsL026RSKcLCwnDy5Ml6nUMQBMyaNQvDhg2rUSjl5OSguLgYwIOi59ixY+jZs+dDz7N48WIUFhbqH5mZmUZ+V9TcVVVV4aeffsJPP/3E5U6IqNE90aMdfl4Qiv5d26C0Qos3d6Tgnd0pKKtgt1xjs4hiKS8vD1qtFq6u1ad9d3V1RXZ2dr3Ocfz4cezYsQP79++HUqmEUqnEuXPnAAA3b95EaGgoAgICEBoainnz5sHPz++h57G2toaTkxM2b96M/v37Y/jw4Q375qjZ0mg0GD9+PMaPHw+NRiN2HCJqhto72WDry/2xYLg3JBJg55lbmKSKw5WcYrGjNStysQM0lUGDBtV6R1K/fv2QnJxs0PkiIiIQERGhb8Yj+jOpVIq+ffvq/01EZAoyqQRvjuiBkC5tsGBHMtJySjAhKg7/nOSLp/t6PPoE9EgW8Q7u4uICmUxWY9B1Tk4O3NzcRMmkUqng4+OD4OBgUV6fzJ+trS1Onz6N06dPc44lIjK5gd1dcGh+KAZ1d0F5pQ5v7z6LyJ3JUGs4DKChLKJYUigUCAoKQnR0tH6bTqdDdHS0ftK/phYREYHU1FScPn1alNcnIiL6s3aO1vjuxX5YOLIHpBJgb+JtTIyKw6Vs3rndEGZTLJWUlCA5OVnfHZaeno7k5GRkZGQAACIjI7Fu3Tps2rQJFy9exGuvvQa1Wo3Zs2eLmJqIiMi8SKUSzB3mje1z+sPVyRrX7qoxKeo4tsdncMJlI5nN1AExMTEYOnRoje3h4eHYuHEjACAqKgorVqxAdnY2lEolVq9ejZCQkCZO+oBKpYJKpYJWq0VaWhqnDqAaysrK9HdwHjlyhF1xRNTk7pVo8NauFMRcvgsAmBDgjo+m+MLRxkrkZOIzZOoAsymWLBXXhqPaqNVqODg4AHjQcmpvz3WciKjp6XQC1sVex6f/uQytToBXWztEzegD344t++akZjfPEpElsra2xr59+7Bv3z5YW1uLHYeIWiipVIJXBnfDzlcGwN3ZBjfuleLJL09g88kb7JarJ7YsGYndcEREZGkKSiuwcNdZHLn44O7ysX5u+Pgpfzi1wG45dsM1IXbDERGRJREEAd/GpeOTXy6hUivAo40toqb3QYBHK7GjNSl2wxGZAa1Wi5iYGMTExECr5fIDRGQeJBIJXg7til2vDkSn1rbIzC/D1LUn8G1cOrvlasGWJSOxG44ehQO8icjcFZZV4q+7z+KXCw+WDhvh44oVU/3Ryk4hcjLTYzdcE2I3HNWmtLRUP8P76dOnYWdnJ3IiIqKaBEHAdydv4sOfLqJCq0PHVrZYPT0QQZ1bix3NpFgsNSEWS0RE1Bycv12IiG2JuHmvFHKpBG+P6ok5oV0hlUrEjmYSHLNEREREBvHt6IyD8wZhvH8HVOkELP/5El7+7gzy1RViRxMdiyUjcSFdIiJqbhxtrLBmeiA+nOILhVyK/13KxdhVsTh9I1/saKJiN1wDsRuOalNWVoaJEycCAA4cOMDlTojIoqRmFWHutkRcz1NDJpUgckQPvDa4W7PpluOYpSbEYolqw7vhiMjSqTVV+Pv+89iXdBsAEOrtgi+mKeHiYPmrEnDMEpEZsLa2xpYtW7BlyxYud0JEFsneWo7PnwnAp0/5w8ZKitgreRizKhYnruWJHa1JsWWpgdiyRERELUFaTjEitibiSm4JpBJgwfAemDusO2QW2i3HlqUmwAHeRETUkvRwdcQPcx/H00GdoBOAL46k4flvf0NucbnY0UyOLUsNxJYlqo1Wq0ViYiIAoE+fPpDJZCInIiJqHHsTb+Hv+8+jtEILFwcFVk4LxCBvF7FjGYQDvJsQiyWqDQd4E1FzdjW3BHO3JeJSdjEkEmDu0O5YMNwbcplldFqxG47IDEgkEnTu3BmdO3eGRGKZffpERLXp3t4B+yMex/R+nhAEYM3/rmLGN78hu7D5dcuxZamB2LJEREQt3YGULCzecxbqCi3a2Cvw+TMBGNKzvdix6sSWJSIiImoyEwPccXB+KHq7OyFfXYFZG07j458voVKrEztao2CxRERERA3WxcUee14biBcGdAYArP31Gp79+hSyCspETtZwLJaITKS8vByTJ0/G5MmTUV7e/PrwiYj+zMZKhn9M8sVXM/vA0VqOhJv3MXZ1LI6k5ogdrUE4ZslIKpUKKpUKWq0WaWlpHLNENfBuOCJqyTLulWLu9kScvVUIAHh5UBe8M/oxKOTm0U7DqQOaEAd4U20qKyuxceNGAMCsWbNgZWUlbiAioiZWUaXDxz9fwvrj6QCAAI9WiJoeCI82diInY7HUpFgsERER1e3whWws3JWCovIqONrIsWKqP0b7dhA1E++GIyIiIrMxsrcbDi0IRaBnKxSXV+HVLYlY8sN5aKq0YkerFxZLRCai0+lw4cIFXLhwATpd87h9lojIWJ1a22HnKwPwyhNdAQCbTt7EU1+dwI08tcjJHo3dcA3EbjiqDQd4ExE93P8u5eCtnSm4X1oJB2s5Pn7KD+P93Zs0A7vhiMyEi4sLXFwsa3FJIiJTG/aYKw4tCEWwV2uUaKowd1sS3t13DuWV5tkt12KKpczMTAwZMgQ+Pj7w9/fHrl27ahxTWlqKzp07Y+HChSIkpObG3t4ed+/exd27d9mqRET0Jx2cbbF9Tn+8PqQbAGDbbxmYrDqOa3dLRE5WU4spluRyOVauXInU1FQcPnwYb7zxBtTq6v2kH374Ifr37y9SQiIiopZFLpPindGPYdOL/dDWXoFL2cWYsCYO+5Nuix2tmhZTLHXo0AFKpRIA4ObmBhcXF+Tn5+v3X7lyBZcuXcKYMWNESkhERNQyDe7RDocWhKJ/1zYordDijR3J+OvusyirMI9uObMplo4dO4YJEybA3d0dEokE+/fvr3GMSqWCl5cXbGxsEBISgvj4eKNeKyEhAVqtFh4eHvptCxcuxPLly42NT1RDeXk5Zs6ciZkzZ3K5EyKiR3B1ssHWl/tj/nBvSCTAjjOZmKSKw5WcYrGjmU+xpFarERAQAJVK9dD9O3bsQGRkJJYsWYLExEQEBARg1KhRyM3N1R+jVCrh6+tb45GVlaU/Jj8/Hy+88AK+/vpr/bYffvgBPXr0QI8ePR6ZU6PRoKioqNqD6GG0Wi22bduGbdu2Qas1j09HRETmTCaVIHJED2x9KQTtHK2RllOCiVHHsetMpqi5zHLqAIlEgn379mHy5Mn6bSEhIQgODkZUVBSAB3PYeHh4YN68eVi0aFG9zqvRaDBixAjMmTMHzz//vH774sWLsWXLFshkMpSUlKCyshJvvfUW3n///RrnWLp0KZYtW1ZjO6cOoD+rrKzUF/8RERFc7oSIyAB3izV4c0cy4q7moaerI36cN6hR15Wz+OVO/lwsVVRUwM7ODrt3765WQIWHh6OgoAA//PDDI88pCAJmzJiBnj17YunSpbUet3HjRpw/fx6fffbZQ/drNBpoNBr910VFRfDw8GCxRERE1Mi0OgFrf72GUb1d0b29Y6Oeu9nNs5SXlwetVgtXV9dq211dXZGdnV2vcxw/fhw7duzA/v37oVQqoVQqce7cOYOzWFtbw8nJCZs3b0b//v0xfPhwg89BREREjyaTShAxtHujF0qGkov66k1o0KBB9VpyYtasWfU6X0REBCIiIvSVKdGf6XQ6ZGRkAAA8PT0hlVrEZxMiIvoTiyiWXFxcIJPJkJOTU217Tk4O3NzcRMmkUqmgUqk4cJdqVVZWhi5dugDgcidERJbMIj7qKhQKBAUFITo6Wr9Np9MhOjoaAwYMECVTREQEUlNTcfr0aVFenyyDnZ0d7OzsxI5BREQNYDYtSyUlJbh69ar+6/T0dCQnJ6NNmzbw9PREZGQkwsPD0bdvX/Tr1w8rV66EWq3G7NmzRUxNVDt7e/sas8QTEZHlMZti6cyZMxg6dKj+68jISAAP7njbuHEjpk2bhrt37+L9999HdnY2lEolfvnllxqDvpsKu+GIiIhaBrOcOsCSGHLrIREREZmHZjd1AJEl0mg0mDNnDubMmVNtbi4iIrIsbFky0h+74dLS0tiyRDWo1Wo4ODgA4N1wRETmxuJn8LYk7Iaj2lRUVGDFihUAgLfffhsKhULkRERE9DsWS02IxRIREZHl4ZilJqBSqeDj44Pg4GCxoxAREZEJsWWpgdiyRLURBAF5eXkAHsxCL5FIRE5ERES/M+Tvt9nMs0TU3JSWlqJ9+/YAOMCbiMiSsVhqoN8b5oqKikROQubmj7N3FxUVcQJTIiIz8vvf7fp0sLFYMtLvUwdUVFQAADw8PERORObM3d1d7AhERPQQxcXFcHZ2rvMYjllqIJ1Oh6ysLDg6OlYbkxIcHFzrIru17fvz9qKiInh4eCAzM1P08VB1fT9NeT5DnveoYxuy/2H7HrbNXK5hc7x+9Tmmvr9rtW03l+sHtMxryPdR05yP1/ABQRBQXFwMd3d3SKV13+/GlqUGkkql6NSpU43tMpms1ota277atjs5OYn+S17X99OU5zPkeY86tiH7H7avruPFvobN8frV5xhDf9f4O9i4z2voNeT7qGnOx2v4fx7VovQ7Th1gIhEREQbvq+s5YmvsbMaez5DnPerYhux/2D5ev8Z9Xn2ONfYa8nfQMq4h30dNcz5eQ8OxG86McVoCy8draNl4/Swfr6HlM4dryJYlM2ZtbY0lS5bA2tpa7ChkJF5Dy8brZ/l4DS2fOVxDtiwRERER1YEtS0RERER1YLFEREREVAcWS0RERER1YLFEREREVAcWS0RERER1YLFkwUpLS9G5c2csXLhQ7ChkBC8vL/j7+0OpVGLo0KFixyEjpKenY+jQofDx8YGfn1+1xZPJvF2+fBlKpVL/sLW1xf79+8WORQb64osv0Lt3b/j4+GD+/Pn1WhTXGFzuxIJ9+OGH6N+/v9gxqAFOnDgBBwcHsWOQkWbNmoUPPvgAoaGhyM/P51w+FqRnz55ITk4GAJSUlMDLywsjRowQNxQZ5O7du4iKisKFCxdgZWWFJ554AqdOncKAAQMa/bXYsmShrly5gkuXLmHMmDFiRyFqkX5/gw4NDQUAtGnTBnI5P39aogMHDmD48OGwt7cXOwoZqKqqCuXl5aisrERlZSXat29vktdhsSSCY8eOYcKECXB3d4dEInlo069KpYKXlxdsbGwQEhKC+Pj4avsXLlyI5cuXN1Fi+rPGuIYSiQSDBw9GcHAwtm7d2kTJ6XcNvYZXrlyBg4MDJkyYgD59+uCjjz5qwvTUGL+Dv9u5cyemTZtm4sT0Zw29hu3atcPChQvh6ekJd3d3hIWFoVu3bibJymJJBGq1GgEBAVCpVA/dv2PHDkRGRmLJkiVITExEQEAARo0ahdzcXADADz/8gB49eqBHjx5NGZv+oKHXEADi4uKQkJCAAwcO4KOPPsLZs2ebKj6h4dewqqoKsbGx+PLLL3Hy5En897//xX//+9+m/BZatMb4HQQerDt24sQJjB07tili0x809Brev38fBw8exI0bN3D79m2cOHECx44dM01YgUQFQNi3b1+1bf369RMiIiL0X2u1WsHd3V1Yvny5IAiCsGjRIqFTp05C586dhbZt2wpOTk7CsmXLmjI2/YEx1/DPFi5cKGzYsMGEKakuxlzDEydOCCNHjtTv//TTT4VPP/20SfJSdQ35Hfzuu++EmTNnNkVMqoMx13Dnzp3C66+/rt//6aefCp988olJ8rFlycxUVFQgISEBYWFh+m1SqRRhYWE4efIkAGD58uXIzMzEjRs38Nlnn2HOnDl4//33xYpMf1Kfa6hWq1FcXAzgweDS//3vf+jdu7coeamm+lzD4OBg5Obm4v79+9DpdDh27Bh69eolVmT6g/pcv9+xC8481ecaenh44MSJEygvL4dWq0VMTAx69uxpkjwcjWhm8vLyoNVq4erqWm27q6srLl26JFIqMkR9rmFOTg6mTJkCANBqtZgzZw6Cg4ObPCs9XH2uoVwux0cffYQnnngCgiBg5MiRGD9+vBhx6U/q+z5aWFiI+Ph47Nmzp6kj0iPU5xr2798fY8eORWBgIKRSKYYPH46JEyeaJA+LJQs3a9YssSOQEbp27YqUlBSxY1ADjRkzhnekWjBnZ2fk5OSIHYMa4MMPP8SHH35o8tdhN5yZcXFxgUwmq/ELnJOTAzc3N5FSkSF4DS0fr6Fl4/WzfOZ2DVksmRmFQoGgoCBER0frt+l0OkRHR5tkoi1qfLyGlo/X0LLx+lk+c7uG7IYTQUlJCa5evar/Oj09HcnJyWjTpg08PT0RGRmJ8PBw9O3bF/369cPKlSuhVqsxe/ZsEVPTH/EaWj5eQ8vG62f5LOoamuQeO6rT0aNHBQA1HuHh4fpj1qxZI3h6egoKhULo16+fcOrUKfECUw28hpaP19Cy8fpZPku6hhJBMNGqc0RERETNAMcsEREREdWBxRIRERFRHVgsEREREdWBxRIRERFRHVgsEREREdWBxRIRERFRHVgsEREREdWBxRIRERFRHVgsEREREdWBxRIRtSizZs2CRCKBRCLB/v37RcsRExOjzzF58mTRchDRo7FYIiKL9cfC54+P0aNH1/m80aNH486dOxgzZky17UePHsX48ePRrl072NjYoFu3bpg2bRqOHTtW70x+fn549dVXH7pv8+bNsLa2Rl5eHgYOHIg7d+7gmWeeqfe5iUgcLJaIyKL9Xvj88bF9+/Y6n2NtbQ03NzdYW1vrt3355ZcYPnw42rZtix07duDy5cvYt28fBg4ciDfffLPeeV566SV8//33KCsrq7Fvw4YNmDhxIlxcXKBQKODm5gZbW9v6f7NEJAoWS0Rk0X4vfP74aN26tUHnyMjIwBtvvIE33ngDmzZtwrBhw9C5c2f4+/tjwYIFOHPmTLXj4+LiEBoaCltbW3h4eGD+/PlQq9UAgOeeew5lZWXYs2dPteekp6cjJiYGL730UsO+YSJqciyWiKjF27NnDyorK/HOO+88dL9EItH/+9q1axg9ejSeeuopnD17Fjt27EBcXBzmzp0LAHBxccGkSZOwfv36aufYuHEjOnXqhJEjR5ruGyEik2CxREQW7eDBg3BwcKj2+Oijjww6R1paGpycnODm5qbftmfPnmrnPHfuHABg+fLlmDlzJt544w14e3tj4MCBWL16Nb777juUl5cDeNAVFxMTg/T0dACAIAjYtGkTwsPDIZXybZfI0sjFDkBE1BBDhw7FV199VW1bmzZtDD7PH1uPAGDUqFFITk7G7du3MWTIEGi1WgBASkoKzp49i61bt+qPFQQBOp0O6enp6NWrF0aMGIFOnTphw4YN+Mc//oHo6GhkZGRg9uzZRnyHRCQ2FktEZNHs7e3RvXv3Bp3D29sbhYWFyM7O1rcuOTg4oHv37pDLq79NlpSU4JVXXsH8+fNrnMfT0xMAIJVKMWvWLGzatAlLly7Fhg0bMHToUHTt2rVBOYlIHGwPJqIWb+rUqbCyssInn3zyyGP79OmD1NRUdO/evcZDoVDoj5s9ezYyMzOxd+9e7Nu3jwO7iSwYW5aIyKJpNBpkZ2dX2yaXy+Hi4lLvc3h6euJf//oXFixYgPz8fMyaNQtdunRBfn4+tmzZAgCQyWQAgL/+9a/o378/5s6di5dffhn29vZITU3Ff//7X0RFRenP2aVLFwwbNgx/+ctfYG1tjSeffLIRvlsiEgNblojIov3yyy/o0KFDtcegQYMMPs+8efNw+PBh3L17F1OnToW3tzfGjh2L9PR0/PLLL/Dz8wMA+Pv749dff0VaWhpCQ0MRGBiI999/H+7u7jXO+dJLL+H+/fuYMWMGbGxsGvy9EpE4JIIgCGKHICJqKrNmzUJBQYGoS538kbnlIaKa2LJERC3O79MNHDx4ULQMsbGxcHBwqHZXHRGZJ7YsEVGLkpubi6KiIgBAhw4dYG9vL0qOsrIy3L59G8CDO+/+OMcTEZkXFktEREREdWA3HBEREVEdWCwRERER1YHFEhEREVEdWCwRERER1YHFEhEREVEdWCwRERER1YHFEhEREVEdWCwRERER1eH/AU86zNs6eI1IAAAAAElFTkSuQmCC",
      "text/plain": [
       "<Figure size 640x480 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "energies = np.geomspace(min_energy, max_energy)\n",
    "\n",
    "fig, ax = plt.subplots()\n",
    "ax.plot(energies, [power_law.spectrum(e) for e in energies])\n",
    "ax.axhline(flux_norm, color=\"k\", linestyle=\":\")\n",
    "ax.axvline(norm_energy, color=\"k\", linestyle=\":\")\n",
    "ax.set_xscale(\"log\")\n",
    "ax.set_yscale(\"log\")\n",
    "ax.set_xlabel(\"E [GeV]\")\n",
    "ax.set_ylabel(\"F $[GeV^-1 cm^-2 s^-1 (sr^-1)]$\")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "9a05928e",
   "metadata": {},
   "source": [
    "We can also use the `PowerLawFlux` class to perform some simple calculations, such as integration of the flux."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "id": "cc227264",
   "metadata": {
    "execution": {
     "iopub.execute_input": "2023-10-27T13:42:29.714639Z",
     "iopub.status.busy": "2023-10-27T13:42:29.714174Z",
     "iopub.status.idle": "2023-10-27T13:42:29.722985Z",
     "shell.execute_reply": "2023-10-27T13:42:29.722202Z"
    }
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([9.999e-13])"
      ]
     },
     "execution_count": 4,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "total_flux = power_law.integrated_spectrum(min_energy, max_energy)  # cm^-2 s^-1 (sr^-1)\n",
    "total_flux"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "id": "654e5335",
   "metadata": {
    "execution": {
     "iopub.execute_input": "2023-10-27T13:42:29.726514Z",
     "iopub.status.busy": "2023-10-27T13:42:29.726001Z",
     "iopub.status.idle": "2023-10-27T13:42:29.733279Z",
     "shell.execute_reply": "2023-10-27T13:42:29.732654Z"
    }
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "9.210340371976183e-08"
      ]
     },
     "execution_count": 5,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "total_energy_flux = power_law.total_flux_density()  # GeV cm^-2 s^-1 (sr^-1)\n",
    "total_energy_flux"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "66081f5f",
   "metadata": {},
   "source": [
    "Sampling from the power law shape is also possible:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "id": "904c45fa",
   "metadata": {
    "execution": {
     "iopub.execute_input": "2023-10-27T13:42:29.737033Z",
     "iopub.status.busy": "2023-10-27T13:42:29.736437Z",
     "iopub.status.idle": "2023-10-27T13:42:30.743748Z",
     "shell.execute_reply": "2023-10-27T13:42:30.742801Z"
    }
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<matplotlib.legend.Legend at 0x7f4f82ad4e50>"
      ]
     },
     "execution_count": 6,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAjUAAAG1CAYAAAABTQXdAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjguMCwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy81sbWrAAAACXBIWXMAAA9hAAAPYQGoP6dpAABXYUlEQVR4nO3deVxU9f7H8dcMMKAiqJkYLlGaJqaQhmilSWKIiqlZtrtlaph1uVbaolkprWblFF1TIdvMTCs1syjD9boUbS6ZYZImaioIKjPMzO8PfnEjN8CBwwzv5+Mxj3vPMue8jyfh4/kux+RyuVyIiIiIeDiz0QFERERE3EFFjYiIiHgFFTUiIiLiFVTUiIiIiFdQUSMiIiJeQUWNiIiIeAUVNSIiIuIVVNSIiIiIV/A1OkBVcTqd7N27l7p162IymYyOIyIiImXgcrk4evQooaGhmM1nfhZTY4qavXv30qxZM6NjiIiISAVkZ2fTtGnTM+5TY4qaunXrAsV/KEFBQQanERERkbLIy8ujWbNmJb/Hz6TGFDV/NTkFBQWpqBEREfEwZek6oo7CIiIi4hW8vqixWq2Eh4cTFRVldBQRERGpRCaXy+UyOkRVyMvLIzg4mNzcXDU/iYjUEA6HA7vdbnQMOQuLxXLakU3l+f1dY/rUiIhIzeFyudi3bx9HjhwxOoqUgdls5qKLLsJisZzTcVTUiIiI1/mroGnUqBG1a9fW/GTV2F/zyP3xxx80b978nO6VihoREfEqDoejpKA577zzjI4jZXD++eezd+9eioqK8PPzq/BxvL6jsIiI1Cx/9aGpXbu2wUmkrP5qdnI4HOd0HBU1IiLildTk5Dncda88rqg5duwYF154IePHjzc6ioiIiFQjHlfUTJ06lc6dOxsdQ0RERKoZj+oovGPHDrZt20ZCQgI//vij0XFERMTDhE1YWmXn2vV0nyo7V1UwmUwsWrSI/v37Gx3ltNz2pCYjI4OEhARCQ0MxmUwsXrz4pH2sVithYWEEBAQQHR3Nhg0bynWO8ePHk5yc7KbEbuTQxE4iInLuDhw4wJgxY2jevDn+/v40btyYuLg41qxZY3Q0j+C2JzUFBQVEREQwfPhwBg4ceNL2+fPnk5SUREpKCtHR0cyYMYO4uDi2b99Oo0aNAIiMjKSoqOik765YsYKNGzfSqlUrWrVqxdq1a8+ap7CwkMLCwpLlvLy8c7i6Myiywdx4uLg7dJ8IPh718EtERKqRG264AZvNRlpaGhdffDE5OTmkp6fz559/Gh3NI7jtSU18fDxPPfUUAwYMOOX26dOnM3LkSIYNG0Z4eDgpKSnUrl2bOXPmlOyTmZnJjz/+eNInNDSU9evX89577xEWFsb48eOZNWsWTzzxxGnzJCcnExwcXPJp1qyZuy61tO3LYM8mWPU8pCVA7p7KOY+IiHi1I0eOsGrVKp555hliYmK48MIL6dSpExMnTqRfv35A8e/Sdu3aUadOHZo1a8Y999xDfn5+yTFSU1OpV68eS5YsoXXr1tSuXZtBgwZx7Ngx0tLSCAsLo379+owbN67U8OmwsDCefPJJbrnlFurUqUOTJk2wWq1nzJudnc1NN91EvXr1aNCgAddffz27du0q2b5y5Uo6depEnTp1qFevHldddRW//fabe//Q/qFKOgrbbDY2b95MbGzs/05sNhMbG8u6devKdIzk5GSys7PZtWsXzz//PCNHjmTSpEmn3X/ixInk5uaWfLKzs8/5Ok6pbX+4YTZY6sLutZByNfy8onLOJSIiXiswMJDAwEAWL15cqqXh78xmMy+//DI//fQTaWlpfPnllzz44IOl9jl27Bgvv/wy7733HsuXL2flypUMGDCAZcuWsWzZMubNm8frr7/OBx98UOp7zz33HBEREXz77bdMmDCB++67j88///yUOex2O3FxcdStW5dVq1axZs0aAgMD6dWrFzabjaKiIvr3788111zD999/z7p167j77rsrfZh9lbSVHDx4EIfDQUhISKn1ISEhbNu2rVLO6e/vj7+/P1arFavVes4T+pxRu0EQejksGAr7vod3boQrx0GPSeBT8ZkRRUSk5vD19SU1NZWRI0eSkpJChw4duOaaa7j55ptp3749APfff3/J/mFhYTz11FOMHj2aV199tWS93W7ntddeo0WLFgAMGjSIefPmkZOTQ2BgIOHh4cTExPDVV18xePDgku9dddVVTJgwAYBWrVqxZs0aXnzxRXr27HlS1vnz5+N0OnnjjTdKCpW5c+dSr149Vq5cyRVXXEFubi59+/YtydGmTRv3/oGdgscN6QYYOnQozz//fJn2TUxMZMuWLWzcuLFSM4U9t43Wu5JILbqueMXal9k8pQtXTkgjbMLSKu1xLyIinumGG25g7969fPzxx/Tq1YuVK1fSoUMHUlNTAfjiiy/o0aMHTZo0oW7dutxxxx38+eefHDt2rOQYtWvXLikkoPgBQlhYGIGBgaXW7d+/v9S5u3TpctLy1q1bT5nzu+++45dffqFu3bolT5gaNGjAiRMn2LlzJw0aNGDo0KHExcWRkJDASy+9xB9//HGufzxnVSVFTcOGDfHx8SEnJ6fU+pycHBo3blyp57ZarYSHhxMVFVWp5wEoxMLjRUMZbbufPFdtOpp3sMx/IrHmzZV+bhER8Q4BAQH07NmTxx57jLVr1zJ06FAmT57Mrl276Nu3L+3bt2fhwoVs3ry5pN+LzWYr+f4/351kMplOuc7pdFY4Y35+Ph07diQzM7PU5+eff+bWW28Fip/crFu3jiuvvJL58+fTqlUr1q9fX+FzlkWVFDUWi4WOHTuSnp5ess7pdJKenn5SZehuVfWk5u+WOzvRxzaVTOfF1DMV8IblBR7znVc8UkpERKQcwsPDKSgoYPPmzTidTl544QU6d+5Mq1at2Lt3r9vO88+CY/369adtMurQoQM7duygUaNGtGzZstQnODi4ZL/LL7+ciRMnsnbtWi677DLeeecdt+U9FbcVNfn5+SWVGkBWVhaZmZns3r0bgKSkJGbNmkVaWhpbt25lzJgxFBQUMGzYMHdFqFayXSHcaHucN4riARjh+ynMiYPDu4wNJiIi1dKff/7Jtddey1tvvcX3339PVlYWCxYs4Nlnn+X666+nZcuW2O12XnnlFX799VfmzZtHSkqK286/Zs0ann32WX7++WesVisLFizgvvvuO+W+t912Gw0bNuT6669n1apVZGVlsXLlSsaNG8fvv/9OVlYWEydOZN26dfz222+sWLGCHTt2VHq/Grd1FN60aRMxMTEly0lJSQAMGTKE1NRUBg8ezIEDB5g0aRL79u0jMjKS5cuXn9R52N2qpKPwadjx5amiO1jvDOd5vxTq7f0GUrrB9TMhvF+V5xERqemq8yy/gYGBREdH8+KLL7Jz507sdjvNmjVj5MiRPPzww9SqVYvp06fzzDPPMHHiRLp160ZycjJ33nmnW87/73//m02bNjFlyhSCgoKYPn06cXFxp9y3du3aZGRk8NBDDzFw4ECOHj1KkyZN6NGjB0FBQRw/fpxt27aRlpbGn3/+yQUXXEBiYiKjRo1yS9bTMblcLlelnqGayMvLIzg4mNzcXIKCgtx+/LN1BA7lIGtbvgW///8syp3uhuueAl9/t2cREanJTpw4QVZWFhdddBEBAQFGx/EIYWFh3H///aVGV1WlM92z8vz+9sjRT55oLw1h2DK46v8f5W34D8zuCX/uNDaYiIiIl/D6oqYqRz+dlY8f9HwCbl0AtRrAH9/B69fAjwuNTiYiIuLx1PzkJuWdh6Yxf/KSxUq0uXjywbeLevBE0R1sf/rUr5kQEZGyUfOT51Hzk4fbx3ncanuEV4r643SZuM03ncWWSXBwh9HRREREPJLXFzXVqvnpHxz48ELRTdxpn8ABVxBtzLuLm6O+f9/oaCIiIh7H64saIybfK6/Vznb0LkxmrSMc7AXw4Uj4aCzYjp39yyIiIgLUgKLGUxygPrfbH4buEwETfDsPZl0L+yvnhZ8iIiLeRkVNNeLEDN0nwJCPITAEDmyF/3SHb982OpqIiEi15/VFTXXuU3NaF3WD0avh4hgoOg4f3QOLRkNhvtHJRETEg61cuRKTycSRI0fK/J2wsDBmzJhRaZncyeuLGk/oU3NKgY3g9g/h2sfAZIbv3oVZMbDvR6OTiYhIJRk6dCgmk4nRo0eftC0xMRGTycTQoUOrPpiHcNu7n8Q9Tp7vpg2dTI/wsmUmjQ/+zInXuhOQ8Bx0HAomkxERRUSkEjVr1oz33nuPF198kVq1agHF87i88847NG/e3OB01ZvXP6nxBhtcbehdmMxXjggCTHZYcj8sHAEn8oyOJiIibtahQweaNWvGhx9+WLLuww8/pHnz5lx++eUl6woLCxk3bhyNGjUiICCAq6+++qRWiWXLltGqVStq1apFTEwMu3btOul8q1evpmvXrtSqVYtmzZoxbtw4CgoKKu36KpOKGg9xiCCG2x8g2X4LmHyKX63wn2uKX7UgIiJn5nKBraDqPxWctH/48OHMnTu3ZHnOnDkMGzas1D4PPvggCxcuJC0tjW+++YaWLVsSFxfHoUOHAMjOzmbgwIEkJCSQmZnJXXfdxYQJE0odY+fOnfTq1YsbbriB77//nvnz57N69WrGjh1bodxG8/rmJ6vVitVqxeFwGB3lnLkw87ojgY0nWvOK5RWaHPqVwpRrearoduY5egKlm6N2Pd3HmKAiItWN/RhMC6368z68Fyx1yv2122+/nYkTJ/Lbb78BsGbNGt577z1WrlwJQEFBAa+99hqpqanEx8cDMGvWLD7//HNmz57NAw88wGuvvUaLFi144YUXAGjdujU//PADzzzzTMl5kpOTue2220rezn3JJZfw8ssvc8011/Daa6953GsmvP5Jjcd2FD6Db1yt6F2YzOeODvibinjSLxWr30sE4ZmPC0VEpLTzzz+fPn36kJqayty5c+nTpw8NGzYs2b5z507sdjtXXXVVyTo/Pz86derE1q1bAdi6dSvR0dGljtulS5dSy9999x2pqakEBgaWfOLi4nA6nWRlZVXiFVYOr39S461yCWSk/d+McH7KQ77v0sdnA+1MWYy1j+N7Vwuj44mIVC9+tYufmhhx3goaPnx4STOQ1Wp1V6JS8vPzGTVqFOPGjTtpmyd2SlZR49FMzHb0ZpOzFTP9XqG5+QAfWB7n6aJbmePoZXQ4EZHqw2SqUDOQkXr16oXNZsNkMhEXF1dqW4sWLbBYLKxZs4YLL7wQALvdzsaNG0uaktq0acPHH39c6nvr168vtdyhQwe2bNlCy5YtK+9CqpDXNz/VBN+5WtLHNo1PHVFYTA4m+c1jlt90OHbI6GgiIlJBPj4+bN26lS1btuDj41NqW506dRgzZgwPPPAAy5cvZ8uWLYwcOZJjx44xYsQIAEaPHs2OHTt44IEH2L59O++88w6pqamljvPQQw+xdu1axo4dS2ZmJjt27OCjjz7y2I7CKmq8RB51GGO/n0n2IRS6fOnpsxle7wbZG4yOJiIiFRQUFERQUNAptz399NPccMMN3HHHHXTo0IFffvmFzz77jPr16wPFzUcLFy5k8eLFREREkJKSwrRp00odo3379nz99df8/PPPdO3alcsvv5xJkyYRGmpAp2o3MLlcFRxv5mHy8vIIDg4mNzf3tP+BnIuTJ80zTltTFla/lwkz54DZF3pMgi73glk1rIh4vxMnTpCVlcVFF13kcaN3aqoz3bPy/P7Wbzkv9JPrIvrapkLbgeAsgs8nwbuDoeBPo6OJiIhUGq8vajzyhZZukE9tGDQH+s4AH3/YsQJSrobf1hodTUREpFJ4fVHjjfPUlJnJBFcMg5FfwnmXwNG9kNoXMp4Hp9PodCIiIm6lId1e7O/9fGozkaf85jDQZzV8+SQZny/mX/Z72Pz0rQYmFBERcR+vf1IjxY4RQJJ9DA/Y7+a4y0I3nx/41H8iZGUYHU1ERMQtVNTUKCYWOLrTz/YUPzub0Mh0BN68HlY+DU7PfzeWiMjf1ZDBvV7BXfdKRU0NtMPVlH62p3i/6BpwOWFlMszrD0dzjI4mInLO/Pz8ADh27JjBSaSsbDYbwEmTDJaX+tTUUCfw58GiUdx0462wJKm4GSrlKhg4C1rEGB1PRKTCfHx8qFevHvv37wegdu3amEwmg1PJ6TidTg4cOEDt2rXx9T23ssSjipqwsDCCgoIwm83Ur1+fr776yuhIni/iZgjtAAuGwv6fYN4A6Ppv6D4RfDzqPw8RkRKNGzcGKClspHozm800b978nItPj/uttXbtWgIDA42O4V3ObwUj02H5RNg8F1Y9D7vXwQ1vQJBnTpUtIjWbyWTiggsuoFGjRtjtdqPjyFlYLBbMbpj13uOKGnGv0q936EmCuQ7T/GZT97c1/PlCJ5Ls95A27WHD8omInAsfH59z7qchnsNtHYUzMjJISEggNDQUk8nE4sWLT9rHarUSFhZGQEAA0dHRbNhQvpctmkwmrrnmGqKionj77bfdlFz+7hPnlSTYnuJHZxjnmY6SZnkGPp8MDv1LR0REqje3PakpKCggIiKC4cOHM3DgwJO2z58/n6SkJFJSUoiOjmbGjBnExcWxfft2GjVqBEBkZCRFRUUnfXfFihWEhoayevVqmjRpwh9//EFsbCzt2rWjffv2p8xTWFhIYWFhyXJeXp6brtT77XJdwA22x3nY922G+H4Oa2b8f3PUbKjXzOh4IiIip1Qpb+k2mUwsWrSI/v37l6yLjo4mKiqKmTNnAsW9nZs1a8a9997LhAkTyn2OBx54gLZt2zJ06NBTbn/88ceZMmXKSetrwlu63amXeQMpdedCYS7Uqg/9X4PW8UbHEhGRGqLavaXbZrOxefNmYmNj/3dis5nY2FjWrVtXpmMUFBRw9OhRAPLz8/nyyy9p27btafefOHEiubm5JZ/s7Oxzu4gaarmzE4z6uniE1PHD8O7N8NkjUGQzOpqIiEgpVdJR+ODBgzgcDkJCQkqtDwkJYdu2bWU6Rk5ODgMGDADA4XAwcuTIM75529/fH39/f6xWK1arFYdDM+ZWVNizW/Djfib4vssI309h3Uwy13zKWPu9/O4qbjrc9XQfg1OKiEhN5zGjny6++GK+++67cn8vMTGRxMTEksdXUjF2fHmy6A7WO9vwvF8KkeadLLM8zAP2UXzmPH1xKSIiUlWqpPmpYcOG+Pj4kJNTehr+nJyckgmSKovVaiU8PPyMT3Wk7D53XkHvwmS+cbYkyHSM1y0vMtk3DYoKz/5lERGRSlQlRY3FYqFjx46kp6eXrHM6naSnp9OlS5dKPXdiYiJbtmxh48aNlXqemmQP53OTbRIpRX0BGOb7GczuCYd+NTiZiIjUZG4ravLz88nMzCQzMxOArKwsMjMz2b17NwBJSUnMmjWLtLQ0tm7dypgxYygoKGDYsGHuiiBVqAhfni66lWG2BzjkCoQ/voOUbvDjh0ZHExGRGsptfWo2bdpETMz/XoSYlJQEwJAhQ0hNTWXw4MEcOHCASZMmsW/fPiIjI1m+fPlJnYfdTR2FK9dXzsvpXZjM+lbvFM9l88Ew2LUK4qaBXy2j44mISA1SKfPUVEflGedeEd46T01Z7ZoaByunwarpgAtCLoMbU6HhJUZHExERD1bt5qmRGsDHF3pMgtsXQu2GkPMjvH4NfP++0clERKSG8PqiRqOfqljLHjB6NYR1BXsBfDgSPhoLtmNGJxMRES+n5ic3qenNT/9kxsk43w8Z57MIs8nFdmdTWo9dCI0uNTqaiIh4EDU/ieGcmJlRNIjb7A+z31WP1ubfYVYMfKu3q4uISOXw+qJGzU/GWudsS+/CZFY5LgP7MfjoHlg0GgrzjY4mIiJeRs1PbqLmpzMz4SSr9zb4ahq4nNCwVfHoqJDTv5RUREREzU9S7bgwQ7cHYMgSqHsBHPwZZl0Lm1OhZtTVIiJSyVTUSNUKu6p4dFTLWCg6AZ/cBwvvgsKjRicTEREP5/VFjfrUVEN1GsKtCyB2Cph84McP4PVuxa9aEBERqSCvL2r0QstqymyGq++HYZ9CUNPil2G+0RM2zFJzlIiIVIjXFzVSzTWPhtGroFU8OAph2XhYMARO5BqdTEREPIzbXmgpcjZnHiF2OyN8zuMh33exbPmouClq0Fxo0qHK8omIiGfTkxqpJkzMdvTmRttkqNccDu+C2dfB+tfUHCUiImXi9UWNOgp7lu9cLWHUKmiTAE47LJ8A82+H44eNjiYiItWc1xc16ijsgWrVg5vmQfxz4GOBbUsgpRv8vsnoZCIiUo15fVEjHspkgui7YcQKqH8R5O6GOXGw9hVwOo1OJyIi1ZCKGqneQi+HURnQdgA4i2DFo/DuzXDskNHJRESkmtHoJ6l2Tj1KahC3+jRgsu+b+O/4DFKuhkFzoHnnKs8nIiLVk57UiIcw8Y6jB/1tT7DTeQHk7YG5vWHVC2qOEhERQEWNeJitrgvpZ3sK2t0ELgekPwFvD4L8A0ZHExERg6moEY9TQC0Y+B/oNxN8a8HO9OLmqF2rjY4mIiIG8vqiRvPUeCmTCTrcASO/hIatIX8fpCXAymfA6TA6nYiIGMDrixrNU+PlQsLh7q8g8nZwOWHlNJjXH47mGJ1MRESqmNcXNVIDWOpAfysMeB38akNWBqRcBTu/MjqZiIhUIRU14j0iboa7v4ZGbaHgAMwbAF8+BY4io5OJiEgV0Dw14pHO9MZvfx5gsu+b3Or7JWQ8B7+thRvegKDQKkwoIiJVTU9qxOsUYuHhorvghtlgCYTf1hSPjtrxhdHRRESkEqmoEe/VblDxKxYat4Njf8LbN8AXj4PDbnQyERGpBB5V1GRlZRETE0N4eDjt2rWjoKDA6EhS3Z3XAkZ8AVEji5dXvwipfSD3d2NziYiI23lUUTN06FCeeOIJtmzZwtdff42/v7/RkcQT+AVAn+fhxjTwD4Ls/xY3R21fbnQyERFxI48pan766Sf8/Pzo2rUrAA0aNMDXV/2cpRza9i9ujgq9HI4fhncHw2ePQJHN6GQiIuIGbitqMjIySEhIIDQ0FJPJxOLFi0/ax2q1EhYWRkBAANHR0WzYsKHMx9+xYweBgYEkJCTQoUMHpk2b5q7oUpM0uAiGfwad7yleXjeTb5/owtUTUwmbsPS0HxERqf7c9qijoKCAiIgIhg8fzsCBA0/aPn/+fJKSkkhJSSE6OpoZM2YQFxfH9u3badSoEQCRkZEUFZ08p8iKFSsoKipi1apVZGZm0qhRI3r16kVUVBQ9e/Y8ZZ7CwkIKCwtLlvPy8tx0peIpzlyMXE1PcwDP+6VwufkXllom8qB9FJ859ToNERFP5baiJj4+nvj4+NNunz59OiNHjmTYsGEApKSksHTpUubMmcOECRMAyMzMPO33mzRpwhVXXEGzZs0A6N27N5mZmactapKTk5kyZUoFr0Zqgs+dV9C7MJlXLK/QwfwLr1teZG5RHMlFt2LDz+h4IiJSTlXSp8Zms7F582ZiY2P/d2KzmdjYWNatW1emY0RFRbF//34OHz6M0+kkIyODNm3anHb/iRMnkpubW/LJzs4+5+sQ77OH87nJNonXi/oAMMz3Mz6wPE5zk94dJSLiaaqkqDl48CAOh4OQkJBS60NCQti3b1+ZjuHr68u0adPo1q0b7du355JLLqFv376n3d/f35+goCDmzZtH586d6dGjxzldg3ivInxJLrqNYbYHOOQKpL05i6WWh+ljXm90NBERKQePGf0ExU1cP/zwAz/++CPTp08v03f0lm4pq6+cl9O7MJmNzlbUNR3HanmZp3xn449GR4mIeIIqKWoaNmyIj48POTmlH+nn5OTQuHHjSj231WolPDycqCh1AJWz28d53Gx7DGtRPwBu901nkWUyHPzF4GQiInI2VVLUWCwWOnbsSHp6esk6p9NJeno6Xbp0qdRz60mNlJcDH54rupk7bQ9x0BVEuPk3eL0bfP++0dFEROQM3FbU5Ofnk5mZWTKCKSsri8zMTHbv3g1AUlISs2bNIi0tja1btzJmzBgKCgpKRkOJVDcZzgh6FyazzhEO9gL4cCR8NBZsx4yOJiIip+C2Id2bNm0iJiamZDkpKQmAIUOGkJqayuDBgzlw4ACTJk1i3759REZGsnz58pM6D7ub1WrFarXicDgq9TzinfZTn9vsD/PrtT/A18/Ct/Ngz2a4MRXOb210PBER+RuTy+VyGR2iKuTl5REcHExubi5BQUFuP75mnfVuu57uA7+uhIUjoWA/+NWGPi9A5K1GRxMR8Wrl+f3tUaOfRAx1cXcYs6b4f+3HYPEYWDQabHpbvIhIdeD1RY1GP4lbBTaC2z+EmEfBZIbv3oX/dIecn4xOJiJS46n5yU3U/FTzRJu28pJlJo1Nhznh8iOg3wvQ4U4wmYyOJiLiNdT8JFIF/utqQ+/CZFY6Iggw2eGTccUjpAqPGh1NRKRG8vqiRs1PUpkOEcQw+wM8bb8ZTD7wwwJ4/Rr443ujo4mI1DhqfnITNT/JrnsawgfDIe938PGHXtPgihGlmqPK8t/Jrqf7VGZMERGPouYnESM0j4bRq6BVPDgKYem/YcFQOJFrdDIRkRpBRY2IO9VuALe8C9dNBbMvbFlc/IqFPd8YnUxExOt5fVGjPjVS5UwmuHIsDP8MgpvD4V0w+zpYnwLUiNZeERFDqE+Nm6hPjZxKEPk86zeLXj7FL1T9zHEFD9jvJo/A035HfWpERP5HfWpEqok8Ahltv5/J9iEUunyJ89nEMv+HiTT9YnQ0ERGvo6JGpNKZSHPEcYPtcXY5Q2hqOsgCyxTu8lmKmqNERNxHRY1IFfnRdTEJtqkscXTGz+TgUb+3ecPveeqhyfpERNzB64sadRSW6uQotRlrv5dH7MMpdPkR6/MtS/0fpqNpu9HRREQ8ntcXNYmJiWzZsoWNGzcaHUXk/5l42xHLANsUfnU2ponpT+ZbnmSMz8eYcBodTkTEY3l9USNSXW1xhZFgm8pix5X4mpw85Pcec/2eg4KDRkcTEfFIvkYHEKnJCqjF/fZE1jrb8oRvKt19vmPfs1cwzjaWDa42p/2ehn2LiJxMT2pEDGfifUcM/WxPscPZhMamw7xreYqxPoswqzlKRKTMVNSIVBM/u5rRz/YkHzi64WNyMd5vAW/6JXM+R4yOJiLiEVTUiFQjxwlgvH00/7aN5pjLn6t9fmKZ/0SuNP9odDQRkWpPRY1INbTQ2Y0E21NsczbjfFMub/kl8y/fD9QcJSJyBl5f1GieGvFUO11N6G97gneLYjCbXNzn+yHvWKbSiMNGRxMRqZb0Qks30QstpTL1M69hmt9sAk0n+NNVl/PuSIWWsUbHEhGpdHqhpYiX+dh5FQm2qWxxXsh5pqPw1g3wxePgKDI6mohItaGiRsRDZLkuYIBtCvOK/v8JzeoXIbUP5P5ubDARkWpCRY2IBynEwmNFw2HQXPAPguz1kHI1bF9udDQREcOpqBHxRJcNhFFfwwWRcPwwvDsYPnsEimxGJxMRMYyKGhFP1eBiGLECokcXL6+bCXPj4fBvxuYSETGIxxQ127dvJzIysuRTq1YtFi9ebHQsEWP5+kP8MzD4bQgIhj2b4PWusHWJ0clERKqcxxQ1rVu3JjMzk8zMTFavXk2dOnXo2bOn0bFEqoc2fWHUKmhyBZzIhfm3wacPQVGh0clERKqMxxQ1f/fxxx/To0cP6tSpY3QUkeqj/oUw7FPoMrZ4+b8pMPs6OPSrsblERKqI24qajIwMEhISCA0NxWQynbJpyGq1EhYWRkBAANHR0WzYsKFC53r//fcZPHjwOSYW8UK+FoibCrfMh1r14Y9M8l66knsefoywCUtP+RER8Ra+7jpQQUEBERERDB8+nIEDB560ff78+SQlJZGSkkJ0dDQzZswgLi6O7du306hRIwAiIyMpKjp5MrEVK1YQGhoKFM8suHbtWt57770z5iksLKSw8H+P3vPy8s7l8kQ8S+teMHo1G18YQJT5Z161vMy8oi08VXQ7hVhK7VqWwmbX030qK6mIiNu4raiJj48nPj7+tNunT5/OyJEjGTZsGAApKSksXbqUOXPmMGHCBAAyMzPPep6PPvqI6667joCAgDPul5yczJQpU8p+ASLeJrgpt9ge5V++H5Do+zF3+H5BB/MOEu3j2OW6wOh0IiJuVyV9amw2G5s3byY29n/vqjGbzcTGxrJu3bpyHausTU8TJ04kNze35JOdnV3u3CKerghfniu6mTttD3HQFURb828ssTxCP/Nao6OJiLhdlRQ1Bw8exOFwEBISUmp9SEgI+/btK/NxcnNz2bBhA3FxcWfd19/fn6CgIObNm0fnzp3p0aNHuXOLeIsMZwS9C5NZ72xDoOkEL1tmkuw7C380WZ+IeA+PGv0UHBxMTk4OFovl7Dv/v8TERLZs2cLGjRsrMZlI9bef+txme5iXigbidJm4xfcrPrI8RgvTHqOjiYi4RZUUNQ0bNsTHx4ecnJxS63NycmjcuHGlnttqtRIeHk5UVFSlnkfEEzjw4cWiQdxun8gBVzCXmrP5xPIoA80ZRkcTETlnbusofCYWi4WOHTuSnp5O//79AXA6naSnpzN27NhKPXdiYiKJiYnk5eURHBxcqecSqSrnOhR7rfMy4gufZobfTK72+YnplhS6FG1hUtFQjnPmTvgiItWV257U5Ofnl8z4C5CVlUVmZia7d+8GICkpiVmzZpGWlsbWrVsZM2YMBQUFJaOhRKRqHSSYO+0Ted5+Iw6XiRt9M/jY8hitTOpULyKeyW1PajZt2kRMTEzJclJSEgBDhgwhNTWVwYMHc+DAASZNmsS+ffuIjIxk+fLlJ3Uedjer1YrVasXhcFTqeUQ8kRMzMx0D2Oi8lJcsM7nEvIePLI8xuWgI7zu6AyajI4qIlJnJ5XK5jA5RFf5qfsrNzSUoKMjtx9fMrOLpGpDHi36vco3P9wAsclzFo/bhFFBLk++JiGHK8/vbo0Y/iUjlOUQQQ+0P8oz9ZopcZgb4rOFjy6O0Mf1mdDQRkTLx+qJGo59Eys6Fmdcc/Rhse4y9rga0MP/BYssk2DgbasZDXRHxYGp+chM1P4m3qcdRXvBLoYfPt8Ur2g6AhJcgQKMIRaTqqPlJRM7ZEepyl/3fPGW/Dcy+8NMieP0a2Put0dFERE7J64saNT+JVJwLM284+sCw5RDcHA5nwezr4L+vqzlKRKodry9q9JoEETdoFgWjM+DSvuCwwacPwvzb4fhho5OJiJTw+qJGRNykVn0Y/Bb0egbMfrBtCbzeDX7fbHQyERFARY2IlIfJBJ1Hw4gVUD8MjuyGOdfB2plqjhIRw3l9UaM+NSKVoEkHGJUB4deDswhWPALv3gLHDhmdTERqMA3pdhMN6RZvdtoZhV0u2DQblj8MjkL2uM7jXtu9fONqVb7jiIichoZ0i0jVMJkg6i646wt+dTamielP3rc8wWifjzHhNDqdiNQwbnuhpYh4r7I8iazDVKb6zaa/z1om+L1HZ/NWkuxjOIT7n4yKiJyKntSIiFsUUIv77Yk8ZB/JCZcf3X2+Y5n/RDqZthodTURqCK8vatRRWKQqmZjviOF625P84gylsekw71qeYqzPIsxqjhKRSqaOwm6ijsIipdXiBE/6pTLIJwOA1Y62XP3QYghsZGwwEfEo6igsIoY7TgDj7aP5t200x1z+XO3zE7x2Ffy60uhoIuKlVNSISKVa6OxGP9uTbHc2hYL98GZ/+GoaOB1GRxMRL6OiRkQq3S+uplxvexI63Am44Otn4M3rIe8Po6OJiBdRUSMiVeIE/tDvFRj4BlgCYdcqSLkafvnC6Ggi4iVU1IhI1Wp/I9z9NYS0g2MH4a0b4Isp4CgyOpmIeDgVNSJS9Rq2hLu+gCuGFy+vng5pfSF3j7G5RMSjeX1Ro3lqRKopvwDo+yIMmguWurB7XXFz1M+fGZ1MRDyU1xc1iYmJbNmyhY0bNxodRURO5bKBMDoDLoiA44fgnZtgxaPgsBudTEQ8jNcXNSLiARpcDCM+h06jipfXvgJz4+HIbmNziYhHUVEjItWDrz/0fhYGvwUBwfD7xuLmqG2arVtEykZFjYhUL20SYNQqaNIRTuTCe7fCpxOgyGZ0MhGp5lTUiEj1U/9CGLYcuowtXv7vazDnOjiUZWwuEanWVNSISPXka4G4qXDLe1CrPuz9Fl7vBj8tNjqZiFRTvkYHKI8XX3yRN954A5fLRWxsLC+99BImk8noWCJSmVrHFzdHLRwB2f+FBUNg111w3dTiYeFnEDbh7P1xdj3dx11JRcRgHlPUHDhwgJkzZ/LTTz/h5+dHt27dWL9+PV26dDE6moi42amKEV8SSfIN4R7fj2HjG8UFzo1pcF4LAxKKSHXkUc1PRUVFnDhxArvdjt1up1GjRkZHEpEqUoQvzxbdzBDbQ/zpqgv7fihujvrhA6OjiUg14baiJiMjg4SEBEJDQzGZTCxevPikfaxWK2FhYQQEBBAdHc2GDRvKfPzzzz+f8ePH07x5c0JDQ4mNjaVFC/0LTaSm+doZQe/CZLjwKrDlFzdLfTwO7MeNjiYiBnNb81NBQQEREREMHz6cgQMHnrR9/vz5JCUlkZKSQnR0NDNmzCAuLo7t27eXPHGJjIykqOjkl9qtWLGCWrVqsWTJEnbt2kWtWrWIj48nIyODbt26nTJPYWEhhYWFJct5eXluulIRMVoODWixfTT3+Z7PWJ+PMH+TxrZNX5JoH8dOVxOj44mIQdxW1MTHxxMfH3/a7dOnT2fkyJEMGzYMgJSUFJYuXcqcOXOYMGECAJmZmaf9/oIFC2jZsiUNGjQAoE+fPqxfv/60RU1ycjJTpkyp4NWISHXnwIfpRTfxX2cbZvi9yqXmbD6xPMqj9mF86Dz1zwUR8W5V0qfGZrOxefNmYmNj/3dis5nY2FjWrVtXpmM0a9aMtWvXcuLECRwOBytXrqR169an3X/ixInk5uaWfLKzs8/5OkSk+lnjbEfvwmTWONpS21TIdEsKz/mmUIsTRkcTkSpWJUXNwYMHcTgchISElFofEhLCvn37ynSMzp0707t3by6//HLat29PixYt6Nev32n39/f3JygoiHnz5tG5c2d69OhxTtcgItXXAepxh30i0+2DcLhM3OibwceWx2hl0j9mRGoSjxnSDTB16lSmTp1aru8kJiaSmJhIXl4ewcHBlZRMRIzmxMzLjoFscF3KS34zucS8h48sjzG5aAjvO7oDnjOnlebXEamYKnlS07BhQ3x8fMjJySm1Picnh8aNG1fqua1WK+Hh4URFRVXqeUSkeljvDKd3YTIZjnbUMtl41m8WL/q9Sh00OkrE21XJkxqLxULHjh1JT0+nf//+ADidTtLT0xk7dmylnltPakSqj7I8gXCHPwlmiP0hxjg/Icl3AQN81tDe9Ctj7ePY6rqwSjKISNVz25Oa/Px8MjMzS0YwZWVlkZmZye7duwFISkpi1qxZpKWlsXXrVsaMGUNBQUHJaCgREXdyYeZVx/XcbHuUva4GtDD/wWLLJG7z+QJwGR1PRCqB257UbNq0iZiYmJLlpKQkAIYMGUJqaiqDBw/mwIEDTJo0iX379hEZGcny5ctP6jzsblarFavVisPhqNTziEj1tMl1KX0Kp/G83+v08PmWqX5z6GLewgT7XeRT2+h4IuJGJpfLVSP+yfJX81Nubi5BQUFuP35VPVYXkYox4eQun2U86DsfP5ODXc4QEu3jWJpcuU3gFaGOwiL/U57f3x717icRkYpyYWaWoy832Sbxu6shYeYcPrRMhv++DjXj33YiXs/rixqNfhKRv/vWdQm9C6fxmeMK/E1F8OmD8P4dcPyI0dFE5Bx5fVGTmJjIli1b2Lhxo9FRRKSayCOQUfZ/8bj9TjD7wdZP4PWu8Ptmo6OJyDnw+qJGROTUTKQ6esGIFVA/DI7shjnXwdqZao4S8VBeX9So+UlEzqhJBxiVAeHXg7MIVjwC794Cxw4ZnUxEysnrixo1P4nIWQUEw41p0OcF8PGHnz+FlK6w+79GJxORcvD6okZEpExMJoi6C+76Ahq0gLzfYW48rH4RnE6j04lIGaioERH5uwvaw6iv4bJB4HLAF4/DOzdBwUGjk4nIWXh9UaM+NSJSbv514YY3IOEl8A2AXz6HlKth1xqjk4nIGXh9UaM+NSJSISYTdBwKI7+Ehq3g6B+Q1he+fg6ceu2KSHVUJW/pFhGprsr0SoIpX8Gy8fDdu/DVU/Dbahg4CwIbVUFCESkrr39SIyJyzvwDYUAKXP8q+NWGX1cWN0f9+rXRyUTkb1TUiIiU1eW3wciv4Pw2kJ8Db14PXyWrOUqkmvD6okYdhUXErRpdWtzP5vI7ABd8/XRxcXN0n9HJRGo8ry9q1FFYRNzOUhuun1ncr8avDuxaBa9dBb+kG51MpEbz+qJGRKTStL+p+BULIe3g2EF46wZIfwIcRUYnE6mRVNSIiJyLhi2LZyG+YgTgglUvFA/9zt1jdDKRGkdFjYjIufILgL7TYdBcsNSF3euKR0f9vMLoZCI1iooaERF3uWwgjM6ACyLg+CF450ZY8Rg47EYnE6kRVNSIiLhTg4thxOfQaVTx8tqXi1+MeWS3sblEagDNKCwi4m6+/tD7WQi7Gj4aC79vhJSu0P9VwlKNDifivbz+SY3mqRERw4T3K26OatIRThyB925lku+b+KHRUSKVweuLGs1TIyKGqh8Gw5ZDl7EADPddzgeWx2lmyjE2l4gX8vqiRkTEcL4WiJsKt7zHEVcdIsy/stTyMPHm/xqdTMSrqKgREakqrePpXZjMJmcrgkzHec3yEk/4zsUfm9HJRLyCihoRkSq0l4bcbHuU14oSALjT93MWWh4nzPSHwclEPJ+KGhGRKlaEL88U3cJQ24P86arLZeZdLLE8Qj/zWqOjiXg0DekWETHISmckvQuTedkyk2jzNl62zKRz0U9MKRpCIZYzfjdswtKzHn/X033cFVXEI3jUk5rnn3+etm3bctlll/HWW28ZHUdE5Jzl0IBbbY/wclF/nC4Tt/p+xWLLY7Qw6d1RIuXlMUXNDz/8wDvvvMPmzZvZuHEjM2fO5MiRI0bHEhE5Zw58mF50E3faJ3DAFUQbczYfWx5lgHmV0dFEPIrHND9t3bqVLl26EBAQAEBERATLly/n5ptvNjiZiIh7rHa2o3dhMi/5WbnSZwsvWl6jS9EWJhcN4TgB5T5eWZqoQM1U4j3c9qQmIyODhIQEQkNDMZlMLF68+KR9rFYrYWFhBAQEEB0dzYYNG8p8/Msuu4yVK1dy5MgRDh8+zMqVK9mzR49nRcS7HKA+t9sfZrp9EA6XiZt8v+Yjy2NcYvrd6Ggi1Z7bntQUFBQQERHB8OHDGThw4Enb58+fT1JSEikpKURHRzNjxgzi4uLYvn07jRo1AiAyMpKiopOnD1+xYgXh4eGMGzeOa6+9luDgYDp37oyPj89p8xQWFlJYWFiynJeX54arFBGpfE7MvOwYyAbXpbzsN5NW5j18bHmUSUVDWeC4BjAZHVGkWjK5XC6X2w9qMrFo0SL69+9fsi46OpqoqChmzpwJgNPppFmzZtx7771MmDCh3Oe46667GDBgAH36nPqx6eOPP86UKVNOWp+bm0tQUFC5z3c2ZX3MKyKex53NM+X9WXEeubzo9yrdfH4A4EPH1TxqH86xCjRHnY6an6Q6y8vLIzg4uEy/v6uko7DNZmPz5s3Exsb+78RmM7Gxsaxbt67Mx9m/fz8A27dvZ8OGDcTFxZ1234kTJ5Kbm1vyyc7OrvgFiIgY5E+CGWJ/iGftN+FwmRjos5pPLI9wqWm30dFEqp0q6Sh88OBBHA4HISEhpdaHhISwbdu2Mh/n+uuvJzc3lzp16jB37lx8fU8f39/fH39/f6xWK1arFYfDUeH8IiJGcmHmVUd/Njov5WXLTFqY/+Ajy2NMKbqTdxzXouYokWIeM6QbYN26dSVv3O7YsWOZvqO3dIuIt9joupTehdP40hGJv8nONL/ZvOL3CoEcMzqaSLVQJUVNw4YN8fHxIScnp9T6nJwcGjduXKnntlqthIeHExUVVannERGpCocJYoR9PNPst2B3+ZDgs55PLI/Q1pRldDQRw1VJUWOxWOjYsSPp6ekl65xOJ+np6XTp0qVSz60nNSLibVyY+Y8jgcG2x/jd1ZCLzDl8aJnMnT6fAW4f+yHiMdxW1OTn55OZmUlmZiYAWVlZZGZmsnt3cWe2pKQkZs2aRVpaGlu3bmXMmDEUFBQwbNgwd0UQEalRvnG1ok/hND53dMTfVMQTfmm86vcSQRQYHU3EEG7rKLxp0yZiYmJKlpOSkgAYMmQIqampDB48mAMHDjBp0iT27dtHZGQky5cvP6nzsLupo7CIeLNcAhlpT2KYczkTfd+ht88G2pmyGGu/l+9cLY2OJ1KlKmWemuqoPOPcK0Lz1Ih4LyPnqSmP9qadzPR7mebmA9hcPjxTdAuzHfGcbXSU5qmR6qzazVMjIiKV73tXC/raprHM0QmLycFjfm8xy286weQbHU2kSnh9UaPRTyJSk+RRh3vs9/GYfSiFLl96+mxmqf/DdDD9bHQ0kUrn9UWNRj+JSM1jYp7jOgbaniDLGUJT00HmW57kbp9PMOE0OpxIpfH6okZEpKb6yRVGgm0qHzu64Gdy8LDfu8z2e5766AW/4p28vqhR85OI1GT51GacfSwT7HdxwuXHtT6ZLPN/mChT2V9RI+IpvL6oUfOTiIiJ9xzX0t/2JDudF3CB6RDvWp7iHp/Fao4Sr+L1RY2IiBTb5mpOgm0qCx1X42ty8qDf+6T5PQP5B4yOJuIWKmpERGqQYwTwb/sYHrDfzXGXhW4+P0DKVZCVYXQ0kXPm9UWN+tSIiPyTiQWO7iTYnuJnZxPIz4E3r4eVT4NTs6+L5/L6okZ9akRETu0XV1P62Z6Cy28HlxNWJhcXN0f3GR1NpEK8vqgREZHTO4E/XG+FAf8BvzqwaxWkXA07vzQ6mki5qagRERGIGAx3r4SQy6DgAMwbCOlPgqPI6GQiZaaiRkREip3fCu76AjoOA1yw6nlIS4DcPUYnEykTry9q1FFYRKQc/GpBwgy4YTZY6sLutcXNUT+vMDqZyFl5fVGjjsIiIhXQbhCM+hoat4fjh+CdG2HFY+CwG51M5LS8vqgREZEKOq9FcXNUp7uLl9e+DHN7w5FsY3OJnIaKGhEROT1ff+j9HNz0JvgHw+8bipujti0zOpnISVTUiIjI2YVfX9wcFdoBThyB926B5Q9Dkc3oZCIlVNSIiEjZNLgIhn8GnROLl9dbYU4cHN5laCyRv6ioERGRsvO1QK9pcPO7EFAP9n4DKd1gy8dGJxNRUSMiIhVwaW8YvQqadoLCXHj/Dlj2ABQVGp1MajCvL2o0T42ISCWp1xyGLYOr7ite3vAfmN0T/txpbC6psby+qNE8NSIilcjHD3o+Abd9ALUawB/fwevXwI8LjU4mNZDXFzUiIlIFLukJo1dD8yvBdhQ+GA6f3A/240YnkxpERY2IiLhHcBMY8gl0HQ+YYPNceCMWDu4wOpnUECpqRETEfXx8ocdjcMeHUOd8yPmxuDnqu/lGJ5MaQEWNiIi4X4tri5ujwrqCvQAW3Q0fJYLtmNHJxIupqBERkcpRtzHc+RF0nwiY4Nu3YNa1sH+b0cnES1XLombAgAHUr1+fQYMGnbRtyZIltG7dmksuuYQ33njDgHQiIlJmZh/oPgGGfAyBIXBgK/ynO3z7ttHJxAtVy6Lmvvvu48033zxpfVFREUlJSXz55Zd8++23PPfcc/z5558GJBQRkXK5qBuMXgMXx0DRcfjoHlg0GgrzjU4mXqRaFjXdu3enbt26J63fsGEDbdu2pUmTJgQGBhIfH8+KFSsMSCgiIuUWeD7c/iFc+xiYzPDduzArBnJ+MjqZeIlyFzUZGRkkJCQQGhqKyWRi8eLFJ+1jtVoJCwsjICCA6OhoNmzY4I6s7N27lyZNmpQsN2nShD179rjl2CIiUgXMZug2HoYuhbqhcPDn4n42m1PB5TI6nXi4chc1BQUFREREYLVaT7l9/vz5JCUlMXnyZL755hsiIiKIi4tj//79JftERkZy2WWXnfTZu3dvxa/kHwoLC8nLyyv1ERGRauLCK4tHR7XsCUUn4JP7YOEIOKGf1VJxvuX9Qnx8PPHx8afdPn36dEaOHMmwYcMASElJYenSpcyZM4cJEyYAkJmZWaGwoaGhpZ7M7Nmzh06dOp1y3+TkZKZMmVKh84iISBWocx7c+j6sewW+mFL8aoW938KNqXBBhNHpxAO5tU+NzWZj8+bNxMbG/u8EZjOxsbGsW7funI/fqVMnfvzxR/bs2UN+fj6ffvopcXFxp9x34sSJ5Obmlnyys7PP+fwiIuJmZnPxCzGHL4fgZnDo1+JZiDfMUnOUlJtbi5qDBw/icDgICQkptT4kJIR9+/aV+TixsbHceOONLFu2jKZNm5YURL6+vrzwwgvExMQQGRnJv//9b84777xTHsPf35+goCDmzZtH586d6dGjR8UvTEREKlezTjAqA1r3AYcNlo2HBUPgRK7RycSDlLv5qSp88cUXp93Wr18/+vXrV+ZjJSYmkpiYSF5eHsHBwe6IJyIilaF2A7j5bVj/Gnw+CbZ8BHsz4ca50KSj0ek8UtiEpWfdZ9fTfaogSdVw65Oahg0b4uPjQ05OTqn1OTk5NG7c2J2nKjOr1Up4eDhRUVGGnF9ERMrBZIIu98CIz6BeczjyG8yOKy501BwlZ+HWosZisdCxY0fS09NL1jmdTtLT0+nSpYs7T1VmiYmJbNmyhY0bNxpyfhERqYAmHWHUKmiTAE47LJ8A790Gxw4ZnUyqsXIXNfn5+WRmZpaMYMrKyiIzM5Pdu3cDkJSUxKxZs0hLS2Pr1q2MGTOGgoKCktFQIiIiZVKrHtw0D3o/Dz4W2L4UXu8G2fpHqpxaufvUbNq0iZiYmJLlpKQkAIYMGUJqaiqDBw/mwIEDTJo0iX379hEZGcny5ctP6jxcVaxWK1arFYfDYcj5RUTkHJhM0GkkNI2CBUPhcBbM7QU9JkGXe4tHT4n8P5PLVTMaKf/qKJybm0tQUJDbj1+Wzlgi4pnc2ZGyOv6s8JiOoifyiifp++nD4uVLroP+KcXz3cgpeUNH4fL8/laJKyIiniEgCAbNgb4vgo8/7FgBKVfDb2uNTibVhNcXNRr9JCLiRUwmuGI4jEyH81rC0b2Q2hcyngen0+h0YjCvL2o0+klExAs1bgd3fw3tB4PLAV8+CW/fAPkHjE4mBvL6okZERLyUfyAMeB36zQTfWrDzS0i5CrIyjE4mBvH6okbNTyIiXsxkgg53wN1fwfmXQn4OvHk9rHwanBr1WtN4fVGj5icRkRqgURsY+SVE3g4uJ6xMLi5ujpb9vYPi+by+qBERkRrCUgf6W4ubpPzqwK5VxaOjdn5pdDKpIipqRETEu0TcDHevhEZtoeAAzBsI6U+Co8joZFLJvL6oUZ8aEZEa6PxWxcO+Ow4FXLDqeUhLgLy9RieTSuT1RY361IiI1FB+tSDhJbhhNlgCYffa4uaoHZ8bnUwqidcXNSIiUsO1GwSjMqBxezj2J7w9CD6fBA670cnEzVTUiIiI9zuvBYz4HKJGFi+veQlS+8CRbGNziVupqBERkZrBLwD6PA83vQn+wZD93+LmqG3LjE4mbuL1RY06CouISCnh18OoryG0A5w4Au/dAp89AkU2o5PJOfL6okYdhUVE5CQNLoLhn0Hne4qX182Eub3g8C5DY8m58fqiRkRE5JR8LdArGW5+FwLqwZ7NkNINtnxsdDKpIBU1IiJSs13aG0avgqZRUJgL798Byx6AokKjk0k5qagRERGp1xyGfQpX3Ve8vOE/MLsn/LnT2FxSLipqREREAHz8oOcTcOsCqNUA/vgOXr8GfvzQ6GRSRipqRERE/q7VdTB6NTTvAraj8MEwWPIvsB83OpmchYoaERGRfwpuAkOWQNd/AybYNAfeiIWDO4xOJmfg9UWN5qkREZEK8fGFHpPg9oVQuyHk/FjcHPX9+0Ynk9Pw+qJG89SIiMg5adkDxqyBsK5gL4APR8JHY8F2zOhk8g9eX9SIiIics7qN4c6P4JoJgAm+nQezroX924xOJn+jokZERKQszD4QM7G4uAkMgQNbYVYMfPu20cnk/6moERERKY+LrykeHXVxd7Afg4/ugUWjoTDf6GQ1nooaERGR8gpsBLcvgmsfBZMZvnu3+KlNzk9GJ6vRVNSIiIhUhNkM3R4oHvpd9wI4+HNxP5vNaeByGZ2uRqqWRc2AAQOoX78+gwYNKtc2ERGRKhd2VXFzVMueUHQCPhkHC++CwqNGJ6txqmVRc9999/Hmm2+We5uIiIgh6jSEW9+H2Clg8oEfP4DXuxW/akGqTLUsarp3707dunXLvU1ERMQwZjNcfX/xizGDmsKhX+GNnrBhlpqjqki5i5qMjAwSEhIIDQ3FZDKxePHik/axWq2EhYUREBBAdHQ0GzZscEdWERGR6q95NIxeBa3iwVEIy8bDgiFwItfoZF6v3EVNQUEBERERWK3WU26fP38+SUlJTJ48mW+++YaIiAji4uLYv39/yT6RkZFcdtllJ3327t1b8Sv5h8LCQvLy8kp9REREqkTtBnDLuxA3Dcy+sOUjSOkKe74xOplX8y3vF+Lj44mPjz/t9unTpzNy5EiGDRsGQEpKCkuXLmXOnDlMmDABgMzMzIqlLYfk5GSmTJlS6ecRERE5JZMJuiRCs87wwVA48hvMvg6uexKiRxdvF7dya58am83G5s2biY2N/d8JzGZiY2NZt26dO091VhMnTiQ3N7fkk52dXaXnFxERAaBpRxi1Ci7tC047LJ8A82+H44eNTuZ13FrUHDx4EIfDQUhISKn1ISEh7Nu3r8zHiY2N5cYbb2TZsmU0bdq0VEF0pm1/5+/vT1BQEPPmzaNz58706NGjYhclIiJyrmrVg8FvQfxz4GOBbUsgpRtk62XL7lTu5qeq8MUXX1Ro26kkJiaSmJhIXl4ewcHB5xpNRESkYkwmiL4bmkXBgmFwOAvm9oIek6HL2OLRU3JO3Pon2LBhQ3x8fMjJySm1Picnh8aNG7vzVGVmtVoJDw8nKirKkPOLiIiUEno5jPoa2g4AZxF8/hi8ezMcO2R0Mo/n1qLGYrHQsWNH0tPTS9Y5nU7S09Pp0qWLO09VZomJiWzZsoWNG/WIT0REqomAYBg0F/q+CD7+sOMzSLkafqva/qfeptxFTX5+PpmZmSUjmLKyssjMzGT37t0AJCUlMWvWLNLS0ti6dStjxoyhoKCgZDSUiIiIUNwcdcVwGJkO57WEvD2Q2gdWvQBOp9HpPFK5+9Rs2rSJmJiYkuWkpCQAhgwZQmpqKoMHD+bAgQNMmjSJffv2ERkZyfLly0/qPFxVrFYrVqsVh8NhyPlFRETOqHE7uHslLEmCH96H9Cdg1xoY8DoEnm90Oo9icrlqxtzNf3UUzs3NJSgoyO3HD5uw1O3HFJHqYdfTfdx2rOr4s8Kd1yfnwOWCb9+CZQ9A0XEIbAyDZkPY1RU+ZFn+e6vu9788v7/V1VpERKQ6MJmgwx1w91fQsDXk74O0BFj5DDjV2lAWXl/UaPSTiIh4lEZtigubyNvB5YSV02DeADiac/bv1nBeX9Ro9JOIiHgcSx3oby3uV+NXG7K+Lh4d9etKo5NVa15f1IiIiHisiJvh7q+hUVso2A9v9ocvp4KjyOhk1ZLXFzVqfhIREY92fqviYd8dhgAuyHgW3uwHeX8Ynaza8fqiRs1PIiLi8fxqQb+X4YbZYAmE39ZAylWwo3yvDvJ2Xl/UiIiIeI12g2BURvHcNsf+hLdvgC8eB4fd6GTVgooaERERT3JeCxjxBUTdVby8+sXimYhzfzc2VzXg9UWN+tSIiIjX8QuAPi/AjWngHwTZ/y0eHbV9udHJDOX1RY361IiIiNdq27+4OSr0cjh+GN4dDJ89AkU2o5MZwuuLGhEREa/W4CIY/hlEjyleXjcT5vaCw78Zm8sAKmpEREQ8na8/xD8Ng9+GgGDYsxle70qcuWa1UqioERER8RZt+sLo1dDkCjiRy+uWF5nsm4aFmjE6yuuLGnUUFhGRGqVecxi+HK68F4Bhvp/xgeVxmpu8/91RXl/UqKOwiIjUOD5+cN1TDLM9wGFXIO3NWSyxPExv83qjk1Uqry9qREREaqqvnJfTuzCZjc5WBJmO86rlZZ70nYM/3jk6SkWNiIiIF/uD87jZ9hjWon4A3OH7BYssk7nI5H3vjlJRIyIi4uUc+PBc0c3caXuIg64gws2/8YnlEfqZ1xgdza1U1IiIiNQQGc4Iehcms84RTqDpBC9brPDRWLAdMzqaW6ioERERqUH2U5/b7A/zUtFAnC4TfDsP3ugBB7YbHe2cqagRERGpYZyYebFoELfZH4Y6jWD/FvhPd8h8x+ho58TrixrNUyMiInJq65xtYcwauLg72I/B4jGwaDTYCoyOViFeX9RonhoREZEzCGwEt38IMY+CyQzfvVv81CbnJ6OTlZvXFzUiIiJyFmYfuOYBGPIJ1L0ADv4Ms66FzWngchmdrsxU1IiIiEixsKuL3x3VMhaKTsAn42DhXVB41OhkZaKiRkRERP6nTkO4dQHEPg4mH/jxA3j9Gvjje6OTnZWKGhERESnNbIar/wXDlkFQEzi0E96IhY1vVOvmKBU1IiIicmrNOxc3R7WKB0chLP03LBgKJ3KNTnZK1bKoGTBgAPXr12fQoEGl1mdnZ9O9e3fCw8Np3749CxYsMCihiIhIDVG7AdzyLlw3Fcy+sGUxvN4N9nxjdLKTVMui5r777uPNN988ab2vry8zZsxgy5YtrFixgvvvv5+CAs8cSy8iIuIxTCa4ciwM/wyCm8PhXTD7OlifUq2ao6plUdO9e3fq1q170voLLriAyMhIABo3bkzDhg05dOhQFacTERGpoZpeAaMz4NK+4LTD8odg/u1w/LDRyYAKFDUZGRkkJCQQGhqKyWRi8eLFJ+1jtVoJCwsjICCA6OhoNmzY4I6spWzevBmHw0GzZs3cfmwRERE5jVr1YfBbEP8s+Fhg2xJI6Qa/bzI6WfmLmoKCAiIiIrBarafcPn/+fJKSkpg8eTLffPMNERERxMXFsX///pJ9IiMjueyyy0767N27t0wZDh06xJ133sl//vOf0+5TWFhIXl5eqY+IiIi4gckE0aNgxAqoHwa5u2FOHKx9xdDmKN/yfiE+Pp74+PjTbp8+fTojR45k2LBhAKSkpLB06VLmzJnDhAkTAMjMzKxYWoqLlf79+zNhwgSuvPLK0+6XnJzMlClTKnweEREROYvQy2FUBnxyH/y0qLiPzeV3QK16hsRxa58am83G5s2biY2N/d8JzGZiY2NZt27dOR/f5XIxdOhQrr32Wu64444z7jtx4kRyc3NLPtnZ2ed8fhEREfmHgGAYNBf6TIdBsw0raMDNRc3BgwdxOByEhISUWh8SEsK+ffvKfJzY2FhuvPFGli1bRtOmTUsKojVr1jB//nwWL15MZGQkkZGR/PDDD6c8hr+/P0FBQcybN4/OnTvTo0ePil+YiIiInJ7JBFEjiue1MVC5m5+qwhdffHHK9VdffTVOp7Ncx0pMTCQxMZG8vDyCg4PdEU9ERESqIbc+qWnYsCE+Pj7k5OSUWp+Tk0Pjxo3deaoys1qthIeHExUVZcj5RUREpGq4taixWCx07NiR9PT0knVOp5P09HS6dOnizlOVWWJiIlu2bGHjxo2GnF9ERESqRrmbn/Lz8/nll19KlrOyssjMzKRBgwY0b96cpKQkhgwZwhVXXEGnTp2YMWMGBQUFJaOhRERERCpDuYuaTZs2ERMTU7KclJQEwJAhQ0hNTWXw4MEcOHCASZMmsW/fPiIjI1m+fPlJnYeritVqxWq14nA4DDm/iIiIVI1yFzXdu3fHdZaJdcaOHcvYsWMrHMqd1FFYRESkZqiW734SERERKS+vL2o0+klERKRm8PqiRqOfREREagavL2pERESkZvD6okbNTyIiIjWD1xc1an4SERGpGby+qBEREZGaoVq+0LIy/DW3Tl5eXqUc31l4rFKOKyLGc+fPjer4s6Kyfi6K8cry31t1v/9/5TvbHHkAJldZ9vJgf80obLPZ2Llzp9FxREREpAKys7Np2rTpGffx+qLmL06nk71791K3bl1MJlPJ+qioqNP2tzndtn+uz8vLo1mzZmRnZxMUFOT+8OVwpuupyuOV53tl2bci9+l02061zlvvoSfcvzNt199B3UMj1MR7WJ1/F7pcLo4ePUpoaChm85l7zdSY5iez2XzKCs/Hx+e0f/in23a69UFBQYb/ZTzT9VTl8crzvbLsW5H7dLptZ9rf2+6hJ9y/M23X30HdQyPUxHtY3X8XlvU1RzW+o3BiYmK5t53pO0Zzd7aKHq883yvLvhW5T6fbVp3vH7g3nyfcvzNt199B3UMj1MR76C2/C2tM81Nl+utlmbm5uYb/C0MqRvfQs+n+eT7dQ89XHe5hjX9S4w7+/v5MnjwZf39/o6NIBekeejbdP8+ne+j5qsM91JMaERER8Qp6UiMiIiJeQUWNiIiIeAUVNSIiIuIVVNSIiIiIV1BRIyIiIl5BRU0VOHbsGBdeeCHjx483OopUQFhYGO3btycyMpKYmBij40gFZGVlERMTQ3h4OO3ataOgoMDoSFJG27dvJzIysuRTq1YtFi9ebHQsKacXX3yRtm3bEh4ezrhx48r0csqKqDGvSTDS1KlT6dy5s9Ex5BysXbuWwMBAo2NIBQ0dOpSnnnqKrl27cujQIc2F4kFat25NZmYmAPn5+YSFhdGzZ09jQ0m5HDhwgJkzZ/LTTz/h5+dHt27dWL9+PV26dHH7ufSkppLt2LGDbdu2ER8fb3QUkRrprx+kXbt2BaBBgwb4+urfc57o448/pkePHtSpU8foKFJORUVFnDhxArvdjt1up1GjRpVyHhU1Z5CRkUFCQgKhoaGYTKZTPvK0Wq2EhYUREBBAdHQ0GzZsKLV9/PjxJCcnV1Fi+Sd33EOTycQ111xDVFQUb7/9dhUll7+c6z3csWMHgYGBJCQk0KFDB6ZNm1aF6cUdfwf/8v777zN48OBKTiz/dK738Pzzz2f8+PE0b96c0NBQYmNjadGiRaVkVVFzBgUFBURERGC1Wk+5ff78+SQlJTF58mS++eYbIiIiiIuLY//+/QB89NFHtGrVilatWlVlbPmbc72HAKtXr2bz5s18/PHHTJs2je+//76q4gvnfg+LiopYtWoVr776KuvWrePzzz/n888/r8pLqNHc8XcQit8rtHbtWnr37l0VseVvzvUeHj58mCVLlrBr1y727NnD2rVrycjIqJywLikTwLVo0aJS6zp16uRKTEwsWXY4HK7Q0FBXcnKyy+VyuSZMmOBq2rSp68ILL3Sdd955rqCgINeUKVOqMrb8TUXu4T+NHz/eNXfu3EpMKWdSkXu4du1a13XXXVey/dlnn3U9++yzVZJXSjuXv4Nvvvmm67bbbquKmHIGFbmH77//vuuee+4p2f7ss8+6nnnmmUrJpyc1FWSz2di8eTOxsbEl68xmM7Gxsaxbtw6A5ORksrOz2bVrF88//zwjR45k0qRJRkWWfyjLPSwoKODo0aNAcSfFL7/8krZt2xqSV05WlnsYFRXF/v37OXz4ME6nk4yMDNq0aWNUZPmbsty/v6jpqXoqyz1s1qwZa9eu5cSJEzgcDlauXEnr1q0rJY96y1XQwYMHcTgchISElFofEhLCtm3bDEol5VGWe5iTk8OAAQMAcDgcjBw5kqioqCrPKqdWlnvo6+vLtGnT6NatGy6Xi+uuu46+ffsaEVf+oaw/R3Nzc9mwYQMLFy6s6ohyFmW5h507d6Z3795cfvnlmM1mevToQb9+/Solj4qaKjJ06FCjI0gFXHzxxXz33XdGx5BzFB8frxGIHiw4OJicnByjY8g5mDp1KlOnTq3086j5qYIaNmyIj4/PSX/RcnJyaNy4sUGppDx0Dz2f7qFn0/3zfNXtHqqoqSCLxULHjh1JT08vWed0OklPT6+UCYXE/XQPPZ/uoWfT/fN81e0eqvnpDPLz8/nll19KlrOyssjMzKRBgwY0b96cpKQkhgwZwhVXXEGnTp2YMWMGBQUFDBs2zMDU8ne6h55P99Cz6f55Po+6h5UypspLfPXVVy7gpM+QIUNK9nnllVdczZs3d1ksFlenTp1c69evNy6wnET30PPpHno23T/P50n30ORyVdJbpURERESqkPrUiIiIiFdQUSMiIiJeQUWNiIiIeAUVNSIiIuIVVNSIiIiIV1BRIyIiIl5BRY2IiIh4BRU1IiIi4hVU1IiIiIhXUFEjItXS0KFDMZlMmEwmFi9ebFiOlStXluTo37+/YTlE5OxU1IhIpft7gfL3T69evc74vV69evHHH38QHx9fav1XX31F3759Of/88wkICKBFixYMHjyYjIyMMmdq164do0ePPuW2efPm4e/vz8GDB7nyyiv5448/uOmmm8p8bBExhooaEakSfxUof/+8++67Z/yOv78/jRs3xt/fv2Tdq6++So8ePTjvvPOYP38+27dvZ9GiRVx55ZX861//KnOeESNG8N5773H8+PGTts2dO5d+/frRsGFDLBYLjRs3platWmW/WBExhIoaEakSfxUof//Ur1+/XMfYvXs3999/P/fffz9paWlce+21XHjhhbRv35777ruPTZs2ldp/9erVdO3alVq1atGsWTPGjRtHQUEBALfffjvHjx9n4cKFpb6TlZXFypUrGTFixLldsIhUORU1IuIxFi5ciN1u58EHHzzldpPJVPL/d+7cSa9evbjhhhv4/vvvmT9/PqtXr2bs2LEANGzYkOuvv545c+aUOkZqaipNmzbluuuuq7wLEZFKoaJGRKrEkiVLCAwMLPWZNm1auY7x888/ExQUROPGjUvWLVy4sNQxf/jhBwCSk5O57bbbuP/++7nkkku48sorefnll3nzzTc5ceIEUNwEtXLlSrKysgBwuVykpaUxZMgQzGb9eBTxNL5GBxCRmiEmJobXXnut1LoGDRqU+zh/fxoDEBcXR2ZmJnv27KF79+44HA4AvvvuO77//nvefvvtkn1dLhdOp5OsrCzatGlDz549adq0KXPnzuWJJ54gPT2d3bt3M2zYsApcoYgYTUWNiFSJOnXq0LJly3M6xiWXXEJubi779u0reVoTGBhIy5Yt8fUt/eMsPz+fUaNGMW7cuJOO07x5cwDMZjNDhw4lLS2Nxx9/nLlz5xITE8PFF198TjlFxBh6vioiHmPQoEH4+fnxzDPPnHXfDh06sGXLFlq2bHnSx2KxlOw3bNgwsrOz+fDDD1m0aJE6CIt4MD2pEZEqUVhYyL59+0qt8/X1pWHDhmU+RvPmzXnhhRe47777OHToEEOHDuWiiy7i0KFDvPXWWwD4+PgA8NBDD9G5c2fGjh3LXXfdRZ06ddiyZQuff/45M2fOLDnmRRddxLXXXsvdd9+Nv78/AwcOdMPViogR9KRGRKrE8uXLueCCC0p9rr766nIf595772XFihUcOHCAQYMGcckll9C7d2+ysrJYvnw57dq1A6B9+/Z8/fXX/Pzzz3Tt2pXLL7+cSZMmERoaetIxR4wYweHDh7n11lsJCAg452sVEWOYXC6Xy+gQIiL/NHToUI4cOWLoKxL+rrrlEZGT6UmNiFRbfw0DX7JkiWEZVq1aRWBgYKlRVCJSPelJjYhUS/v37ycvLw+ACy64gDp16hiS4/jx4+zZswcoHmn19zlyRKR6UVEjIiIiXkHNTyIiIuIVVNSIiIiIV1BRIyIiIl5BRY2IiIh4BRU1IiIi4hVU1IiIiIhXUFEjIiIiXkFFjYiIiHiF/wM4VMB4fJLMZAAAAABJRU5ErkJggg==",
      "text/plain": [
       "<Figure size 640x480 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "samples = power_law.sample(1000)\n",
    "\n",
    "fig, ax = plt.subplots()\n",
    "ax.hist(samples, bins=energies, density=True, label=\"Samples\")\n",
    "ax.plot(energies, [power_law.spectrum(e) / total_flux for e in energies], label=\"Model\")\n",
    "ax.set_xscale(\"log\")\n",
    "ax.set_yscale(\"log\")\n",
    "ax.set_xlabel(\"E [GeV]\")\n",
    "ax.legend()"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "fe1dcf83",
   "metadata": {},
   "source": [
    "The `BrokenPowerLaw` class is also available and behaves in a very similar way."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "c3a9ec40",
   "metadata": {},
   "source": [
    "## Diffuse and point sources\n",
    "\n",
    "Once the spectral shape is defined, we can specify either a `DiffuseSource` or a `PointSource`. It is assumed that diffuse sources are isotropic and the flux model describes the per-steradian flux over the entire $4\\pi$ sky. We also specify a redshift of the source such that adiabatic neutrino energy losses can be accounted for. Naturally, `PointSource` objects also have a direction specified in (ra, dec) coordinates."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "id": "4c3848b8",
   "metadata": {
    "execution": {
     "iopub.execute_input": "2023-10-27T13:42:30.748238Z",
     "iopub.status.busy": "2023-10-27T13:42:30.747674Z",
     "iopub.status.idle": "2023-10-27T13:42:30.753626Z",
     "shell.execute_reply": "2023-10-27T13:42:30.752976Z"
    }
   },
   "outputs": [],
   "source": [
    "diffuse_source = DiffuseSource(power_law, z=0.0)\n",
    "\n",
    "ra = np.deg2rad(50)\n",
    "dec = np.deg2rad(-10)\n",
    "point_source = PointSource(power_law, z=0.5, coord=(ra, dec))"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "75ab2a27",
   "metadata": {},
   "source": [
    "The original flux model can now be accessed from within the source along with its other properties:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "id": "e92b1383",
   "metadata": {
    "execution": {
     "iopub.execute_input": "2023-10-27T13:42:30.757243Z",
     "iopub.status.busy": "2023-10-27T13:42:30.756759Z",
     "iopub.status.idle": "2023-10-27T13:42:30.764325Z",
     "shell.execute_reply": "2023-10-27T13:42:30.763365Z"
    }
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<icecube_tools.source.flux_model.PowerLawFlux at 0x7f4f8500da90>"
      ]
     },
     "execution_count": 8,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "diffuse_source.flux_model"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "e3116482",
   "metadata": {},
   "source": [
    "Sources and lists of sources can be used as input to simulations, as demonstrated in the simulation notebook."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "ebd41fee",
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "icecube_dev",
   "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.9.18"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}