toy_spectral_lines.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "import yt\n",
    "from yt.units import mp, clight\n",
    "import numpy as np\n",
    "import matplotlib.pyplot as plt"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "ckms = clight.to_value(\"km/s\")\n",
    "sigma_inst = 4.5e-3/2.355 # in keV"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "E0 = 6.7\n",
    "dE = 0.2\n",
    "Emin = E0-0.5*dE\n",
    "Emax = E0+0.5*dE\n",
    "emid = np.linspace(Emin, Emax, 1000)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "T, N = np.loadtxt(\"iron_XXV.dat\", unpack=True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "yt : [INFO     ] 2020-01-30 14:32:46,941 Omega Lambda is 0.0, so we are turning off Cosmology.\n",
      "yt : [INFO     ] 2020-01-30 14:32:46,974 Parameters: current_time              = 1.5380632172269546\n",
      "yt : [INFO     ] 2020-01-30 14:32:46,974 Parameters: domain_dimensions         = [1 1 1]\n",
      "yt : [INFO     ] 2020-01-30 14:32:46,975 Parameters: domain_left_edge          = [0. 0. 0.]\n",
      "yt : [INFO     ] 2020-01-30 14:32:46,976 Parameters: domain_right_edge         = [40000. 40000. 40000.]\n",
      "yt : [INFO     ] 2020-01-30 14:32:46,977 Parameters: cosmological_simulation   = 0\n"
     ]
    }
   ],
   "source": [
    "ds = yt.load(\"ArepoBullet/snapshot_150.hdf5\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [],
   "source": [
    "cyl = ds.disk(\"c\", [1.0, 0.0, 0.0], (50.0,\"kpc\"), (2.0,\"Mpc\"))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "yt : [INFO     ] 2020-01-30 14:32:47,033 Allocating for 2.653e+07 particles\n",
      "Loading particle index: 100%|██████████| 51/51 [00:00<00:00, 1036.02it/s]\n"
     ]
    }
   ],
   "source": [
    "v_x = cyl[\"gas\",\"velocity_x\"].to_value(\"km/s\")\n",
    "kT = cyl[\"gas\",\"kT\"]\n",
    "EM = cyl[\"gas\",\"emission_measure\"]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [],
   "source": [
    "norm = np.interp(kT, T, N)*EM.d"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [],
   "source": [
    "sigma_v = np.sqrt(kT/(56.0*mp)).to_value(\"km/s\")\n",
    "e0 = E0*(1.0+v_x/ckms)\n",
    "sigma_e = e0*sigma_v/ckms\n",
    "sigma_E = E0*sigma_v/ckms"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [],
   "source": [
    "f = np.zeros(emid.size)\n",
    "f0 = np.zeros(emid.size)\n",
    "for i in range(sigma_v.size):\n",
    "    f += norm[i]*np.exp(-0.5*((emid-e0[i])/sigma_e[i])**2)\n",
    "    f0 += norm[i]*np.exp(-0.5*((emid-E0)/sigma_E[i])**2)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [],
   "source": [
    "instrument = np.exp(-0.5*((emid-E0)/sigma_inst)**2)\n",
    "instrument /= instrument.sum()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [],
   "source": [
    "f0_conv = np.convolve(f0, instrument, mode='same')\n",
    "f_conv = np.convolve(f, instrument, mode='same')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<matplotlib.lines.Line2D at 0x128d51400>"
      ]
     },
     "execution_count": 13,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAmIAAAJwCAYAAAAqSA12AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjIsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8li6FKAAAgAElEQVR4nOzde7hcdX3v8c93z77lnkDCJQGSGEACApoJiBesCl5OvVZsqx61WpvZ9ZRaqD09pzcv1do+p9J6qbaZrW2V1kvpQStqPYIKXtBCJqIgBAJJIJAEkpDrvu89v/PHmjUzmcxlrZm1Zvae9X49T541mXX7xT7P9Mv39/19f+acEwAAANqvp9MDAAAASCoCMQAAgA4hEAMAAOgQAjEAAIAOIRADAADoEAIxAACADiEQAwAA6JA5HYiZ2R+Z2U1mtsPMnJntivj5rzKz28zskJmNmtlDZvZ3FddsMLOPmtnWwnWHzOxuM/sfZtYX5XgAAEB3sbnc0NXMnKSnJW2VlJZ01Dm3JqJnv1/SByT9P0nflDQq6RxJlzjnXl923ZckXS3pq5JyklKSXi3pFZK+LemVbi7/jwwAAGIz1wOxZzjndhQ+3ydpYRSBmJldLelWSe9zzn2owbUvkJRzzo1XfP8vkv67pNc4577e6pgAAED3mdNTk34QFpSZXW1m3zazw2Y2bmY/N7PfrnLpH0t6StJfFu5baGZV/7dyzv2oMggr+HLh+KwwYwQAAMkxpwOxMMwsI2+qcKGkv5D0+5IekfT3ZvbXZdctkPQiSf8l6V1m9oSkY5KOm9mXzOz0gK88q3B8MqJ/AgAA6DJzemqyXL2pSTM7U9JOSTc7595Sce7jkq6VdL5z7hEzu1TSPZL2S1os6a8k/UzSlZJ+T9KDkjY650brjGWhpJ9LOkXSM5xzT7f+LwQAAN2mt9MDaJM3ShqQ9FkzW15x7hZJ75F0lbwM2aLC9yskbXLOfabw96+Y2VFJ75f0G5L+vtqLzCwl6V8krZX0FoIwAABQS1KmJtcXjrfJy3SV/7m1cM6fchwrHPOSbqx4zucKxxdXe0mhjuwfJb1O0p84577Y6sABAED3SkpGzArHt0vaW+Mav/D/8cLxkHNuouIa/95lJ73AzCR9pvCODzrnPtL8cAEAQBIkJRDbXjgecM7dVu9C59yTZvaYpLPNbH5FLZhfgP9U+T1lQdg7JX3YOfeBaIYNAAC6WVKmJv9N0oSkD5rZvMqTZrbEzAbKvrpRXhZtqOLSdxeO3yy71yQNS/pNSR9xzv1ZlAMHAADda06vmjSzt0laXfjr70rql3RD4e+POuduLLv2nfKyVrvlBVqPyivIv1jS6yVd6JzbVbh2sbz2FedLyspbNflCeQ1avyvp5c65mcK1H5X03sI1H60yzEeccz+O5l8MAAC6yVwPxG6X9Es1Tt/hnHtxxfUvkPQHkl4gaamkA/LaUXxd0qfKG7MWVld+SF7h/XJ5tWNflPShiuvqjUGSPuece0eIfxYAAEiIOR2IAQAAzGVJqREDAACYdebcqkkzI4UHAADmFOecVfuejBgAAECHzLmMmI/aNgBR8jrR8NsCIFr+b0stczYQA4Aobd68udNDAJBAc27VpF8jNtfGDQAAkqcs206NGAAAwGxCIAYAkrLZrLLZbKeHASBhmJoEAFGsDyAeTE0CAADMUgRiAAAAHUIgBgAA0CEEYgAAAB1CIAYAANAhBGIAAAAdwhZHACDaVgDoDDJiAAAAHUIgBgAA0CEEYgAgKZ1OK51Od3oYABKGLY4AQGxxBCAebHEEAAAwSxGIAQAAdAiBGAAAQIcQiAEAAHQIgRgAAECH0FkfACS94zffpW37jumJw2NatXRep4cDICHIiAGApNP+23u09+K367dvzHV6KAAShEAMQOI553RTbrck6d4njujpkckOjwhAUhCIAUi8Q6NTGtv7sCb2PSxJ2rH/eIdHBCApCMQAJN7OAyPa97nrtO9z10mSdhwY6fCIACQFgRiAxNtZEXhV/h0A4kIgBiDx9h4eq/t3AIgLgRiAxDtYUZxf+XcAiAuBGIDE23984sS/H5uocSUARItADEDiHagIvMiIAWgXAjEAiXegIiP29Mik8nnXodEASBICMQCJd3BkUmf8xsd02/fv1JJ5fZrJOx0aJSsGIH4EYgASzTmno2NTGjjjXL3o+c/Vsvl9kqQjY1MdHhmAJCAQA5Boo5MzyjtpsK9HfakeLZ5HIAagfQjEACTa8YlpSdLBb/2dMpmMlhQCsaPj050cFoCEIBADkGjHxr3M14Et39Tw8LAWD5IRA9A+BGIAEq0y8+VPTR4lEAPQBgRiABLt+EmBWK8kMmIA2oNADECiHasIxEo1YgRiAOJHIAYg0Y5PnBhw+TViTE0CaAcCMQCJVpkRWzToTU2OTMx0YjgAEqa30wMAgE7yA7Gzzr1Ipy0e0Ly+lCRpdJL2FQDiR0YMQKL5gdj7//FryuVyWjBARgxA+xCIAUg0v0ZsYSEAm99PRgxA+xCIAUg0PyO2qFCkX8yITZIRAxA/AjEAieZvcfTqS1fKzEoZsQkyYgDiRyAGINEqO+sv6PcyYqNTZMQAxI9ADECiHa9o3Dp/wM+IEYgBiB+BGIBEq+wj1p/qUW+PaXImr8npfIdGBSApCMQAJNrxilowM9O8Qp3YGAX7AGJGIAYgsZxzGq0SbPl1YiO0sAAQMwIxAIk1UZh6HOg98aewWCdGIAYgZmxxBCCx/GzYvP6UNm/eXPy+mBGjYB9AzAjEACTWWKFFxby+lDKZTPF7v5cYU5MA4sbUJIDEGivLiJXzu+vTwgJA3AjEACTWeFlGLJvNKpvNSiIjBqB9mJoEkFjFGrG+lIaGhiRJmUym1F2f9hUAYkZGDEBiFWvEKqYm/VWTI+w3CSBmgQIxM/sjM7vJzHaYmTOzXWFeYmarCs+4w8z2mtmImf3CzP7azE5tauQA0KKxsoxYufk0dAXQJkEzYh+R9FJJj0g61MR7XiPpA5IOSvprSddJurNw/KmZndHEMwGgJeO1MmLFhq4EYgDiFbRGbJ1zbockmdl9khaGfM8PJK12zu0r+27YzP5L0rCkPyj8AYC2Ga2REVvQT0NXAO0RKCPmB2HNcs79oiII8325cHxWK88HgGb4NWKDlVOTAzR0BdAenS7WP6twfLKjowCQSLWmJv3AbHyaQAxAvDrdvuKDhePnqp00M9fGsQBIGL8Yf35fSs6Vfm4GC3tPjlMjBiBmHcuImdl7Jf2qpKxz7rudGgeA5Bqt0VmfjBiAdulIRszMfkve6slvSLq21nXOOatyL1kyAJGoVSNWDMSm8m0fE4BkaXtGzMx+U1JW0rclXeOcm2r3GABAOnGLo3Q6rXQ6LUka7Os54TwAxKWtGTEze6e8dhW3SXq9c26ine8HgHLFGrH+lLZu3Vr8vpQRIxADEK/IM2Jmdo6ZXWBmfRXfv0PSZyR9T9LrnHPjUb8bAMIoTk1W1oj1MjUJoD0CZcTM7G2SVhf+ukJSv5n9aeHvjzrnbiy7/POSfknSWkm7Cve/VtJnJR2V1zvsGrMTyr+OO+e+2uS/AQCaUmuLI39qcoJifQAxCzo1+S55wVW5DxWOd0i6UfVtkJd9WyqvPqzSo5IIxAC01dhU9UBsgGJ9AG0SKBBzzr046AOrXeuc+4C8vSYBYNYYq9nQlWJ9AO3R6c76ANAxtaYm+1M9MpOm807TM2TFAMSn0531AaBjyjNimzZtKn5vZhrsTWlsakbj03ktTPHfrADiQSAGILHKM2LZ7Inlq4N9PV4gNjWjhQP8VAKIB/+ZByCRnHPFLYwGek/+KaSXGIB2IBADkEhTM07OSb09pt5Uj3K5nHK5XPH8PFZOAmgD8u0AEmmiIhu2ceNGSV6mTCpvYUFGDEB8yIgBSKSJaS/TNVCxYtJHU1cA7UAgBiCRioFYlfowiW2OALQHgRiARJqYql2oL9HUFUB7EIgBSKRSRqzW1CQZMQDxIxADkEilGrFaGTGK9QHEj0AMQCJNNqoR86cmKdYHECPaVwBIpFL7Ci/ztWXLlhPOD1CsD6ANCMQAJNLE1IkZsXQ6fcJ5piYBtANTkwASqXGNGKsmAcSPQAxAIlVOTWYyGWUymeJ5MmIA2oFADEAi+Rmx/pT3Mzg8PKzh4eHi+cFePyNGjRiA+BCIAUikYkNX2lcA6CACMQCJ1HCLo0Ig5l8HAHEgEAOQSI066/cXArRJAjEAMSIQA5BIpWL96j+D/vcTNHQFECMCMQCJVOwjVqNGrJgRmyEjBiA+NHQFkEiVU5MbNmw44bz//QSrJgHEiEAMQCJVTk3mcrkTzpMRA9AOTE0CSKTJBp31izViZMQAxIhADEAiBV41SUYMQIwIxAAkUmUfMTOTmRXP+x33aV8BIE4EYgASqXKvyUr+lCXtKwDEiUAMQCI1al8xkKKzPoD4EYgBSKRGWxyVMmIEYgDiQyAGIJH8Kcf+GoFYeY2Yc65t4wKQLARiABKp0arJnh5Tb49XvM/KSQBxIRADkEjFGrEaGbHyc6ycBBAXOusDSKTKzvqbN28+6Zr+3h6NTM4QiAGIDYEYgEQqTk32eVOTmUzmpGu8acspCvYBxIapSQCJ1GjVpFTWXZ9ADEBMCMQAJM70TF4zeaceU7EgP5vNKpvNnnBdcb9JAjEAMWFqEkDi+KsgB3pTxW2NhoaGJJ04RUlGDEDcyIgBSJxGXfV9pY2/2eYIQDwIxAAkTpD6sPLzfuAGAFEjEAOQOI02/Pb1F85P0NAVQEwIxAAkDhkxALMFgRiAxAlfI0YgBiAeBGIAEqe44XeqQUYs5WfEKNYHEA/aVwBInGobfjvnTrrOz5iREQMQFzJiABKnWKzfaGoyRR8xAPEiEAOQOMUasUbF+oV9KOmsDyAuBGIAEqfa1GQ6nVY6nT7hOjJiAOJGjRiAxCn1ESv9t+jWrVtPuq601yTF+gDiQUYMQOIUM2JB21eQEQMQEwIxAIkzWWVqshoCMQBxIxADkDjBO+tTrA8gXgRiABLHb9Da3yAQIyMGIG4EYgASx9/Eu1EgVirWJxADEA9WTQJIHD/DVb7F0aZNm066rp9ADEDMCMQAJM5UlYxYNps96bp+2lcAiBlTkwASp1pGrJoBasQAxIxADEDiFAOxsoxYLpdTLpc74bpiIMam3wBiwtQkgMTxA6u+sozYxo0bJUnOueJ3xfYVUwRiAOJBRgxA4kxOe8FW4PYVZMQAxIRADEDiTIZuX0GxPoB4BArEzOyPzOwmM9thZs7MdjXzMjN7u5n91MzGzOxJM/uMma1o5lkA0KxJf9PvBsX6NHQFELegGbGPSHqppEckHWrmRWZ2vaTPSToi6fckbZb0Jkm3m9mCZp4JAM2YmvGmJvsaTU2mCMQAxCtosf4659wOSTKz+yQtDPMSM1su6cOS7pZ0lXNupvD93ZK+Ji8w+0iYZwJAswK3r+hjr0kA8QqUEfODsBa8XtJ8SZ/0g7DCc2+RtEPSW1t8PgAEVq19RTVkxADErV3tKy4rHH9c5dxPJL3ZzBY65463aTwAEmyqSvuKLVu2nHRdX8pkJk3nnWbyTqkea9sYASRDuwKxlYXjE1XOPSHJCtc8VH7CzFyV6wGgJf5U40BZRiydTp90nZmpP9Wjiem8JqfzmtefatsYASRDu9pXzC8cJ6qcG6+4BgBiFbR9hUQLCwDxalcgNlo4DlQ5N1hxTZFzzir/xDZCAIlRbWoyk8kok8mcdG1/obs+dWIA4tCuQGxP4biqyrlVklzZNQAQq2rF+sPDwxoeHj7p2lJGjEAMQPTaFYjdXTg+r8q550p6kEJ9AO0StH2FRCAGIF6RB2Jmdo6ZXWBmfWVf/4ekMUnXmlmq7NrXSFon6V+jHgcAVJPPO03nCw1dU42rHeiuDyBOgVZNmtnbJK0u/HWFpH4z+9PC3x91zt1YdvnnJf2SpLWSdkmSc26/mf2ZpI9Kus3MvihvSvK9krZJ+liL/w4ACKRYqJ/qkVmIQIyNvwHEIGj7infJC67KfahwvEPSjWrAOXeDmR2UdL2kT0g6KunfJP1vpiUBtEuYFZNSqaCfjBiAOAQKxJxzLw76wHrXOuf+WdI/B30WAEQtaFd9n19HNkVGDEAM2tXQFQBmhVLrihOnJTds2FD1emrEAMSJQAxAotTKiOVyuarX97NqEkCM2tW+AgBmhTCtKySK9QHEi0AMQKJMVumqX88AxfoAYkQgBiBRJqts+C15G3xXa2fBqkkAcSIQA5AooVdN9rJqEkB8CMQAJMrUjN9VP2SNGBkxADEgEAOQKJMzM5LCZ8Qo1gcQBwIxAIkSetVkivYVAOJDIAYgUfyAqi9sRoxADEAMCMQAJIpfIzYQMiNGIAYgDnTWB5AotVZNbt68uer1pRqxmXgHBiCRCMQAJMrkdPVi/UwmU/X6YvuKaRfvwAAkElOTABIldPuKFKsmAcSHQAxAovgBVWVGLJvNKpvNnnQ9xfoA4sTUJIBEmajRvmJoaEjSyVOUfiBG+woAcSAjBiBRpmpkxGqhoSuAOBGIAUiUZhu6+kX+ABAlAjEAidL8pt+smgQQPQIxAIniB2KhV01SIwYgBgRiABKl6RoxAjEAMSAQA5AoExTrA5hFaF8BIFFqFes7V70GjKlJAHEiIwYgUUpTkxboevqIAYgTgRiARCllxFKBrqd9BYA4EYgBSJRa7SvS6bTS6fRJ19O+AkCcqBEDkCj+1GRf6sSpya1bt1a9nmJ9AHEiIwYgUSZCNnTt7TGZSTN5p5k8WTEA0SIQA5AofmZrIGAgZmasnAQQGwIxAIlSmpoM/vNHU1cAcSEQA5AoYfealEorJydmWDkJIFoEYgASpVZD13pYOQkgLqyaBJAoxU2/KzJimzZtqnkPU5MA4kIgBiBR/KxWZUYsm83WvIdifQBxYWoSQGI454qrJpuZmiQQAxA1AjEAiTFZ1sy1p+fEhq65XE65XK7qfaWmrhTrA4gWU5MAEsOflqzWumLjxo2SvKxZJf96Nv4GEDUyYgASo5nWFVKp+StTkwCiRiAGIDGaaV1Rfj3tKwBEjUAMQGI001VfolgfQHwIxAAkhl/jFXSfSR/F+gDiQiAGIDGarRGjjxiAuBCIAUiMZqcm+5iaBBAT2lcASIxiM9cqGbEtW7bUvK+f9hUAYkIgBiAx6q2aTKfTNe8bYNNvADFhahJAYtTa8LsRVk0CiAuBGIDEqLfPZCaTUSaTqXpfsVifVZMAIkYgBiAxJuu0rxgeHtbw8HDV+8iIAYgLgRiAxGi6fQWBGICYEIgBSIxS+woLdV9fcWqSQAxAtAjEACRGvfYV9ZQyYqyaBBAtAjEAiVFqX5EKdd9ALxkxAPEgEAOQGH4g1dcbbmqytMURqyYBRIuGrgASo7hqskr7ig0bNtS8j2J9AHEhEAOQGPVWTeZyuZr39TM1CSAmTE0CSIymN/1OkREDEA8CMQCJQR8xALMNgRiAxKjXvsLMZFa9iL+0xRHtKwBEi0AMQGJMTNfea7KeYvsKVk0CiBiBGIDEmCpktJqemqRYH0DECMQAJIaf0QqbEaNGDEBcCMQAJEbTxfqsmgQQk8C/RmbWY2bXm9k2Mxs3s91mdoOZLQh4/0Iz+2Mzu9fMjpnZATO708zeYbUqZAEgQv7UZOj2FWTEAMQkzK/R30r6G0n3S/pdSTdJeo+kW8ys7nMK5/9T0ock3S3pvZI+LCkl6Z8k/VXokQNASK1mxKZYNQkgYoE665vZRfKCr5udc9eUfb9T0ickvUnSF+o84rmSXijpY86568vu/7SkbZKGJP2v0KMHgBAm6rSv2Lx5c837Su0r8nLO1WxzAQBhBd3i6M2STNLHKr4flpfNeqvqB2KLC8c95V865ybN7ICkgYDjAICmTdVpX5HJZGre19Nj6kuZpmacJmfyGuhNxTZGAMkSNBC7TFJe0l3lXzrnxs3snsL5eu6SdFjSH5rZLkn/JWmepHdISkv67Wo3mRnzAAAiU6+hayP9qR5NzcxocppADEB0ggZiKyUdcM5NVDn3hKTnm1m/c26y2s3OuUNm9lpJn5H0b2Wnjkm6xjn31TCDBoBmTNbJiGWzWUm1M2P9vT0amZyhYB9ApIL+Z+F8SdWCMEkaL7umnuOS7pP0UUlvkPRbkh6W9AUze1m1G5xzVvkn4HgB4CTFTb+rZMSGhoY0NDRU896+FE1dAUQvaEZsVNJpNc4Nll1TlZldLOlOSdc75/6h7PsvygvOhs1snXOO/UMAxKZeRqwRmroCiEPQX6M9kpabWbWi+lXypi2rTksWXC8vYLup/Evn3Kikb0haLWlNwLEAQFOabV9Rfs8UGTEAEQr6a3R34drLy780s0FJz5a0pcH9qwrHahWuvRVHAIiF375ioMlifam0cTgARCHor9GXJTlJ11V8v0lebdi/+l+Y2Tozu6DiuvsLx3eUf2lmSyW9TtIhSY8EHAsAhOacK9WINTE1OcDUJIAYBMpCOefuNbNPSbrWzG6W9E1J6+V11r9DJ/YQ+468qcbywvqPSXq7pL8q1Iv9SNIp8gK5MyX9jnNuusV/CwDUNJ13ck5K9ZhSPeHX/VAjBiAOYaYDr5O0S1JG0qskHZD0SUnvc87V/WVyzj1qZpdLep+kq+R14h+TdI+k9zrnbg4/dAAIzg+g+lLNLb4uBmLUiAGIUOBArLCi8YbCn3rXranx/SOSfiPM4AAgKv60ZK0Vk87V7x9dbF9BRgxAhMIXSgDAHFRaMdlcV/zSxt8EYgCiQyAGIBH81Y7NrJiUSlOTrJoEECUCMQCJ0GifyXQ6rXQ6XfN+ivUBxIHeXQASoVGx/tatW+veP0CxPoAYkBEDkAhTDTJijfRTrA8gBgRiABKhlX0mJVZNAogHgRiARGhln8ny+1g1CSBKBGIAEmFipsX2FRTrA4gBgRiARChNTbbWWX+CjBiACLFqEkAiNCrW37RpU937KdYHEAcCMQCJ0KhYP5vN1r1/gKlJADFgahJAIrRarM+qSQBxIBADkAiNOuvncjnlcrma9/fT0BVADJiaBJAIpc761QOxjRs3SpKcc1XP074CQBzIiAFIhEYZsUYo1gcQBwIxAIngB1ADTXbWL7avIBADECECMQCJEFVnfTJiAKJEIAYgERrViDUyQLE+gBgQiAFIhEYNXRuhfQWAOBCIAUiElov1WTUJIAa0rwCQCBMNOutv2bKl7v2smgQQBwIxAInQqFg/nU7XvZ9ifQBxYGoSQCIUa8RabF9BsT6AKBGIAUiERhmxTCajTCZT8/6BVEoSfcQARItADEAiNCrWHx4e1vDwcM37+3rNew6BGIAIEYgBSITJBsX6jRSL9WfyNfejBICwCMQAJEKxoWuT7St6Uz3qMck5aSZPIAYgGgRiABJhcsYLnprNiEkU7AOIHoEYgESYnJ6RVNqqqBn0EgMQNQIxAInQamd9715v5SSBGICo0NAVQCI0al+xYcOGhs/ws2m0sAAQFQIxAIlQLNavUSOWy+UaPqMvVWhhQY0YgIgwNQkgEab8Yv2WpibZ+BtAtAjEACRCq33EJPabBBA9AjEAXc85VyrWrxGImZnMrO5zWDUJIGoEYgC6nh+E9faYenrqB1v1kBEDEDUCMQBdL4r6MO/+wsbf1IgBiAiBGICu16h1RVD9KTb+BhAtAjEAXS+KQn2JqUkA0SMQA9D1osuI0b4CQLQIxAB0vckZb59JMmIAZhs66wPoepPTjYv1N2/e3PA5xUCMjBiAiBCIAeh6QTb8zmQyDZ/Tn2LTbwDRYmoSQNeLqli/r9dbNcmm3wCiQiAGoOsFKdbPZrPKZrN1nzNAZ30AEWNqEkDX81c59tXJiA0NDUmqP0XJpt8AokZGDEDXm4iqfQWrJgFEjEAMQNcLUqwfRHHTbzJiACJCIAag6/kZrIGW+4ixahJAtAjEAHQ9P3CqVyMWBFOTAKJGIAag601FNDXZV9j0e4KpSQARIRAD0PWi2mtywF81SUYMQERoXwGg6wUp1nfONXwOWxwBiBoZMQBdbyKizvpscQQgagRiALpeVDViFOsDiBqBGICuF2SvyXQ6rXQ6Xfc5TE0CiBo1YgC6XpBi/a1btzZ8jr9qkowYgKiQEQPQ9aJeNUkgBiAqBGIAut5kgE2/gygW6zM1CSAiBGIAul5ke02SEQMQMQIxAF0vSLF+EBTrA4ha4F8lM+sxs+vNbJuZjZvZbjO7wcwWhHjGKWb2UTN7uPCM/Wb2PTO7srnhA0BjxU2/yYgBmGXCrJr8W0nvkfQVSTdIWl/4+3PM7GrnXN1fJjNbLel2SQslfVbSQ5KWSLpE0qrQIweAgIIU62/atKnhc/yMGoEYgKgECsTM7CJJvyvpZufcNWXf75T0CUlvkvSFBo/5l8L7LnHO7W1uuAAQXpBi/Ww22/A5fvuK6bxTPu/U02PRDBBAYgXN079Zkkn6WMX3w5JGJb213s1m9iJJL5T0f5xze82sz8zmhx0sADQjqs76ZlbKilEnBiACQX+VLpOUl3RX+ZfOuXFJ9xTO1/PLheNjZnaLpDFJI2b2kJnVDeIAoFVBivVzuZxyuVzDZ1GwDyBKQWvEVko64JybqHLuCUnPN7N+59xkjfufWTgOS9ou6TckDUj6fUk3mlmfc+6fKm8yMxdwfABQU5AasY0bN0qSnKv/s9Pf2yNNUCcGIBpBA7H5kqoFYZI0XnZNrUBsUeF4TNJL/IDNzL4iaYekj5jZ5xoV/ANAMyYial9R/gwCMQBRCPqrNCovg1XNYNk1tYwVjl8sz5o55w5J+pqkM1TKmqnsvFX+CTheACiKqkas/BkEYgCiEPRXaY+k5WZWLRhbJW/aslY2TJIeLxz3VTnnr6BcFnAsABBKVJ31pbKNv6kRAxCBoL9Kdxeuvbz8SzMblPRsSVsa3O8X+Z9V5Zz/3VMBxwIAoUS16bf3jNQJzwSAVgT9VfqyJCfpuorvN8mrDftX/wszW2dmF1Rc91V59WFvNbOFZdeeKen1krY75x4OOXYACCSqLY4kVk0CiFagYn3n3L1m9ilJ15rZzZK+qVJn/Tt0YjPX70haLa/vmH//IQkOBXkAACAASURBVDP7A0mbJf3EzP5RUr+kdxeO10bwbwGAk+TzTtN5byWkP63YigGK9QFEKMwWR9dJ2iUpI+lVkg5I+qSk9wVZ7eicy5rZAUl/KOlD8vqS/VjSW5xzPwo5bgAIpFgfluqRWe1AbMuWRhUWHor1AUQpcCDmnJuRt8fkDQ2uW1Pn3M2Sbg76TgBoVdBC/XQ6Heh5BGIAotR6wQQAzGJRFupLrJoEEC0CMQBdzW/mOtAgEMtkMspkMg2fx6pJAFEiEAPQ1SYDBmLDw8MaHh5u+Dw2/QYQJQIxAF1tYnpGUnRTk9SIAYgSgRiArjYx5WfEUpE8b4BADECECMQAdDV/CrHR1GRQNHQFECUCMQBdzc+IRTY1SUNXABEiEAPQ1fwasagyYn0EYgAiFKazPgDMOaVVk/VrxDZs2BDoeUxNAogSgRiArjYRsKFrLpcL9DxWTQKIElOTALpa1FOTZMQARIlADEBXK05N9kXzczdAjRiACBGIAehqEwFrxMxMZtbweUxNAogSgRiArha0RiwoVk0CiBKBGICuFnTT76CoEQMQJQIxAF2tVKwfzRZHfiA2RSAGIAIEYgC6Wlyd9SeYmgQQAQIxAF0ttr0mCcQARIBADEBX8zNiUQViAwRiACJEZ30AXc2vEWs0Nbl58+ZAz6NYH0CUCMQAdLWge01mMplAz6N9BYAoMTUJoKtNRNxZn1WTAKJEIAagqxXbV6Tq/9xls1lls9mGz+snIwYgQkxNAuhqQfeaHBoaktR4ipJVkwCiREYMQFcLutdkUP6qyQmmJgFEgEAMQFeLc69J51wkzwSQXARiALraZMR7TaZ6TKkekyRNzRCIAWgNgRiArhb1XpNSWcE+05MAWkQgBqCrRT01Wf6sKQr2AbSIQAxAV4t6alKiuz6A6NC+AkBXmwgYiIUpvKeXGICokBED0LWmZ/KayTuleky9DRq6hlFsYUEgBqBFBGIAupY/ddgfYRAmsd8kgOgQiAHoWhNTwfeZTKfTSqfTgZ5LjRiAqFAjBqBrBa0Pk6StW7cGfi7bHAGIChkxAF1rMobWFVJ5jdhMpM8FkDwEYgC6VhzNXCVpsM97nj/1CQDNIhAD0LXCTE2GMVioORsnIwagRQRiALpWHF31pVKGjYwYgFYRiAHoWqWpSTJiAGYnVk0C6FqlqcnGNWKbNm0K/FwyYgCiQiAGoGuFWTWZzWYDP3eAjBiAiDA1CaBrxVWs72fExsmIAWgRgRiArjUxFbx9RS6XUy6XC/Rcv0aMPmIAWsXUJICuVdxrMkBGbOPGjZIk51zDawepEQMQETJiALpWca/JqKcmyYgBiAiBGICuVawRC7DpdxiD1IgBiAiBGICu5a+aHEiREQMwOxGIAehaxYaufRHvNUlGDEBECMQAdK3Y2lf4fcSmyIgBaA2BGICuFaahaxiDhQybH+gBQLNoXwGga/kZq8EAfcS2bNkS+LmlqUkyYgBaQyAGoGuNh1g1mU6nAz+3VKxPRgxAa5iaBNC1xiYLGbHYivXJiAFoDYEYgK7lr5oMEohlMhllMplAzyUjBiAqBGLAbLD359J3/0La9o1Oj6SrlGrEGv/UDQ8Pa3h4ONBzyYgBiAo1YkCnPfwd6Qu/LuWnvL9f9X7pyt/v7Ji6hN/na15/tFOTZMQARIWMGNBJkyPSV9/tBWFrrvS+++6Hpad3dHZcXaKYEYu4RszvSzY5nVc+33iTcACohUAM6KS7Pysdf1Ja+Rzp7V+TLn2z5GakH3+60yPrCuPTwdtXhGFmxd5kZMUAtIJADOiU6Unpzk94n1/yp1JPj/S83/H+ft+/e+fREn9qcjDiTb+lUt0Z+00CaAWBGNAp278tjeyXVqyXzr3K++6Mi6UVF0hjh6TdP+ns+LrA+GQ8e01KpelO9psE0IpAgZiZ9ZjZ9Wa2zczGzWy3md1gZgvCvtDM5pvZTjNzZvZ34YcMdIl7vuAdn/PfJbPS9+sKQdmO29s+pG5TnJoMkBHbsGGDNmzYEPjZpYJ9MmIAmhc0I/a3kv5G0v2SflfSTZLeI+kWMwubVftzSctD3gN0l4ljXkbMeqSLf+3Ec+te4h0f+V77x9VFZvJOUzNOZlJ/qvHPVC6XUy6XC/z8UgsLMmIAmtewfYWZXSQv+LrZOXdN2fc7JX1C0pskfSHIy8xsg6TrJP2hpBuaGTDQFXbc4a2UPPsKadHpJ55b/Xypp0/a81NvinLess6McY4r32fSyjOOESEjBiAKQbJZb5Zkkj5W8f2wpFFJbw3yIjNLFe75lqSbQ4wR6D7bv+0dz7v65HP9C6RVGyQ56YmtbR1WN/EDsah7iPnIiAGIQpBA7DJJeUl3lX/pnBuXdE/hfBDXS7pA0rVhBgh0Heekh2/zPp/38urXnPls77j3nvaMqQv5G34H6aoveS0pwmTO/IwY3fUBtCLIL9RKSQeccxNVzj0habmZ9dd7gJmtlfRBSX/unNsVdHCFgv4T/gS9F5i1nrpfOvqEtPB06YxLql+zshCI7SEQa1ZczVx9fkaMPmIAWhEkEJsvqVoQJknjZdfU8/eSdsor+AeSbfut3vHcq09cLVmOjFjL/EAsjtYVUnn7CjJiAJoXZK/JUUmn1Tg3WHZNVWb2Vkkvl/Qi59xUmME55076/1JkxTDn+YHYeS+rfc3y86XeedLhx6TRp6X5p7RnbF2klBGLp13iAJ31AUQgyC/UHnnTjwNVzq2SN21ZtQV44Z6/kfRNSfvM7FwzO1fS6sIlSwrfLW1i7MDcM37Ea9RqKekZL6l9XapXOv1C7/NT97dnbF2m2FU/4u2NfANkxABEIEggdnfhusvLvzSzQUnPlrSlzr3zJK2Q9CpJ28v+3F44/9bC338rzKCBOWvH7VJ+Wjr7udK8Bv/9sfyZ3nH/g7EPqxuREQMwFwSZmvyypD+W1//rB2Xfb5JXG/av/hdmtk5Sn3NuW+GrEUm/WuWZKyR9Wl4ri89K+nnokQNzUXFaskrbikorzveOBx6KbzxdrLTPJDViAGavhoGYc+5eM/uUpGvN7GZ504zr5XXWv0MnNnP9jrxpRyvcOyXp3yufaWZrCh8fcc6ddB7oSuVtK86tUx/mIyPWkmIfsYCB2ObNm0M9v5gRIxAD0IIgGTHJy4btkpSRN814QNInJb3POUdeHgjiyfukY3ulhWd4m3s3sqIQiJERa4q/z2TQVZOZTCbU8/2MGFOTAFoRKBBzzs3I25Ko7rZEzrk1AZ+3S4WsGZAY5dOSQRqHLl0tpfq9nmMTx6SBRfGOr8uUpibjqREbpKErgAjE8wsF4GTFQKxGN/1KqV7p1HO9zwcfjmdMXSxsQ9dsNqtsNhv4+QM0dAUQAQIxoB3GDku7/0vq6ZWe8eLg9y1b6x0P7YphUN2tfNPvIIaGhjQ0NBT4+WTEAESBQAxohx3fk9yMdPYV0uCS4PctW+MdCcRCi799BRkxAK0jEAPaYbu/yXeAthXlTilkxJ7eGe14EiD+9hVkxAC0jkAMiFs+Lz3s7y8ZoG1FOTJiTWtXRswP+ACgGQRiQNz23iMdf1JatFI6/aJw9xKINW18Ot6M2Lz+QkZsmowYgOYRiAFxe/A/veMz/1uwthXllpwtyaQjj0szU5EPrZuFXTUZlv/csUkCMQDNIxAD4vaQH4j9cvh7+walxSu9Qv8jj0c7ri4XdyDmd+wfo0YMQAsIxIA4Hd4t7btX6l8orb2yuWfQwqIppfYVwX7mnHNyzgV+/vx+rx82GTEArSAQA+L00Le847qXSr0DzT2jWCfGyskw4l41Oa+fqUkArSMQA+L04De9YzPTkj4K9pvC1CSAuYBADIjL+BFp5w8k6wm+rVE1BGJN8VczBm1fkU6nlU6nAz+/v7dHvT2m6bzTJE1dATQp0KbfAJpw/9ek/JS05kppwanNP2fpOd6RYv1Qwk5Nbt26NfQ75vWldGxiWmNTM+oPWIsGAOX45QDicu9N3vGSX2vtOUtWeUcCsVDinpqUqBMD0DoCMSAOR/dKO78vpfql9a9t7VkLz/CmN48/KU1PRDO+BJgoZsTi+5krBmLUiQFoEoEYEIf7/q8k59WGzVva2rNSvV5Xfkk6uqfloSXBTN5pciYvM6k/FWMgVsi2jU5Ox/YOAN2NQAyIQ1TTkj5/evLoE9E8r8uVeoilZGF3MwiBqUkArSIQA6J2YLu3v+TAYum8V0TzzCVneUfqxAKJe8Nv33ymJgG0iFWTQNR+/m/ecf1rvS2KorCYgv0wmtnwe9OmTaHfU5qaJBAD0BwCMSBKzpVNS/5qdM/1M2JMTQYyVqjZ8qcOg8hms6HfM6+wzdE4GTEATWJqEojSEzlvK6KFZ3j9w6LC1GQofoZqfohArBnzClOfZMQANItADIiSnw171jVST4RBQHFqkoxYEMVArC940j+XyymXy4V6Dxt/A2gVU5NAVGamC20rJF38xmifveRs70hGLBA/MAozNblx40ZJknMu8D2D7DcJoEVkxICo7PqBNLJfOmWdtPI50T57/ilS76A0cUQaPxrts7tQu6Ym59O+AkCLCMSAqPziZu/4rGukqHtXmVGwH8JoE8X6zWDVJIBWEYgBUZiZkh64xfv8rDfE8w7qxALzVzHGXqxPHzEALSIQA6Kw4w5p7JC0Yr102vp43lFcObk7nud3kdLUZLxlsH5GbIwtjgA0iUAMiII/LXnRr8T3DqYmA/MDsTANXZtBZ30ArSIQA1o1PSk98HXvc1zTkhJTkyGMtWlqcrCfGjEAraF9BdCqR3/orWY87SJp+Xnxvcff+JupyYb8Yv0wgdiWLVtCv2d+IeNGZ30AzSIQA1q1/TbveH5EG3zXstifmtwT73u6gJ+hmhdiajKdTod+j1+DRkYMQLOYmgRatf3b3vG8l8X7Hj8jdvQJb09L1DTWrmL9fu8nlBoxAM0iEANacWiXdHC7NLBEOuvyeN81sMh7z/S4NPp0vO+a45pp6JrJZJTJZEK9Zx5bHAFoEYEY0Irtt3rHdS+WUm2Y6V+80jseZaujeprZ4mh4eFjDw8Oh3kNDVwCtIhADWvHId73juVe3533F6UnqxOoZnQpfrN8M2lcAaBWBGNCsfF569E7v89oXteedxRYWZMTqaddekwO9PTKTJqfzmslTtwcgPAIxoFn7t0njh6VFK6Wlq9vzzsVlBfuoqTQ1Ge90sZmVuuuTFQPQBAIxoFmPFbJhq58X/SbftTA1GUgxIxZzZ32pvE6MbY4AhEcgBjTr0R97x3Oe17530l0/kGaK9Zvlv2N8Mh/7uwB0Hxq6As1wTnqsEIitfn773lvcb5IasVqmZ/KanMmrx7warqA2bNjQ1PuKGbEpMmIAwiMQA5pxZLdXpzW4VFqxvn3vLbav2OMtFughqV1pdKrUzNVCTBnncrmm3ldcOUkLCwBN4FccaMbjhX0Jz7qsvcFQ/wIv+JuZlEYPtu+9c0g7pyUlabCPQAxA8wjEgGbsvcc7rnxO+9/N9GRd7Wpd4fPfQ1NXAM0gEAOasccPxJ7d/nf705MU7Fc1MuHVai0I2brCzEJNZfoWDHjvGWHVJIAmEIgBYTnX2YzYYlpY1HO8EIgtHGhPCaz/Hv+9ABAGgRgQ1qGd0vgRacFp0qIz2//+Yi+xWTY1+bMvSZ96rvRXq6V/e7t0eHdHhlHMiA20Z2qymBEjEAPQBAIxIKzyacl2NXItNxt7if3wY9JXhkq7Ddz/H9JnX96RMR4vBmLtyYgtKGbEqBEDEB6BGBDWnp96xzM7UB8mzb5tjh6+Tbrt/ZJM+uWPSr/3M+nsK6Rje6QvvUXKtzdAGSkERO2bmkwV3ktGDEB4BGJAWHt/5h07Uagvla2anAWB2OSIdMv13ueX/ol0+SZp2RrpzV+Ulpzt1dJt+ce2DmmkQxkxAjEAzSAQA8J66n7vePqzOvP+YlPXvV5T1066+7PSkcekMy6WXnB96fv5p0iv/Evv8/f+Qpo43rYhFacm29S+wl+dSbE+gGbQWR8I4/h+aWS/1L/Qy/h0Qt88ad4p0tjT0shT0qIzOjOOyVHpzk94n696v5Sq+Dm54NVew9vH75Z+eqN0xbvbMqxmM2KbN29u6n1kxAC0gowYEMb+B7zjaes7u73QkllQJ3bf//WC0pXPkc69+uTzZtILrvM+//jTbasV8/t5hQ3EMpmMMplM6PctKNaIUawPIDwCMSCMJwvTkqe1cX/JahYX6sQ6uXIy98/e8fJM7dWjz/xlaelqb/pyx+1tGdbxthfrMzUJoHkEYkAYfn3YaRd2dhzFOrEOBWL77pOe2CINLJEufH3t63p6pOe81fv8039py9BGm5yazGazymazod/nv2eUzvoAmkAgBoTxlD812eFAzJ+aPNKhpq5bP+cdL/k1qX9+/Wuf/RZJJm37hjRxLPahHW+yoevQ0JCGhoZCv28hfcQAtIBADAjKudkTiPlTk53Y5mhmWrrvZu/zhrc1vn7JWdI5V0gzE9JD/y/esalUI9auqUmK9QG0gkAMCOrIbmnymLRghbRwRWfH0smpyUd/KI0ekE5ZJ51xSbB71r/WOz7wtfjGVeAXzberj9j8Pi/zNjY1o5m8a8s7AXQPAjEgqKe2eccVF3R2HFLZ1GQHArFffNU7XvQrwbd4Wv8a77j9Vq/tRYzavel3T48Ve5aNUCcGICQCMSCogw97x+Xnd3YcUmmbo2N727uF0Mx0Kat1UZ0i/UpLz5ZWpaWpUW9LpBi1u7N++buYngQQFoEYEJQfiJ16bmfHIUm9A94UqZuRjj/Zvvc+dqc0etCblgy7s8AFr/aO2+OrE8vnnUYnvcDUnzJsh4UEYgCaRCAGBHVwu3ecDYGYVKoTa+f0pF9sf8Grgk9L+vymrw9/11v4EINiM9f+lHp6Qo6vBQtYOQmgSQRiQFAHH/GOy2dJIOZvsXTksfa9c/ut3vH8V4S/94yLpYWnS8f2lFafRqyVQn3nnFyTAWKpuz4ZMQDhBArEzKzHzK43s21mNm5mu83sBjNbEODe883sz83sJ2a238yOmdk9ZvYnQe4HZoXJEW+FYk+ftOScTo/Gs3S1dzzcpkDs0C7pwINeE9eznxv+fjNp3Uu9z498J9Kh+dpdqO/z33dsnEAMQDhBM2J/K+lvJN0v6Xcl3STpPZJuMbNGz/hNSddLekTSn0v6n5IelPRhSXea2bwmxg20l58NO2XtyZtbd8qyQiB26NH2vM/Phq17iZTqa+4Z667yjg/HE4gdHZ+SJC0abO//jRYPev97HCu8HwCCavhrZWYXyQu+bnbOXVP2/U5Jn5D0JklfqPOIf5f0l865I2Xf/YOZbZf0J5LeJenvmhg70D7FQv3zOjuOcsWMWJsCMb8+7LyXN/+MdS+RZNKjd3ptLBp15Q/Jz0gtnhc+UEyn05KkXC4X+l4/8CMjBiCsIBmxN0sySR+r+H5Y0qikt9a72Tm3pSII8325cAy59AroAD8jduq6zo6jXDszYpOj0q4feJ/9ovtmLFgurXy212X/0R9FM7YyR8e8jJSfoQpj69at2rp1a1Pv9QO/o2TEAIQUJBC7TFJe0l3lXzrnxiXdUzjfjMIeLaq59t7MXOWfJt8FtGY2ta7wLS3Uqh3ZLeXz8b5r1w+l6XFp5XOkRae39ix/evKR77U+rgp+RqrdU5NkxAA0K0ggtlLSAefcRJVzT0habmb9YV5qZilJ75M0rfrTmsDs4LeuWD6Lpib7F3i9xGYmvcaucdoewbSkb+2LvKOfYYuQn5FqZmqyFX4Gzs/IAUBQQQKx+ZKqBWGSNF52TRgfk3SFpPc55x6sdZFzzir/hHwP0DrnZmdGTGpPnZhz0vZve5+jCMTOvlxK9Uv77pVGn279eWX8QGhRm1dN+oEfGTEAYQUJxEYlDdQ4N1h2TSBm9iFJ10rKOuf+Muh9QMeMHpTGj0gDi70M1GzSjjqx/Q96LTLmL5dWbmj9eX3zpLMuk+Skx37c+vPKtFKs3wp/apIaMQBhBQnE9sibfqwWjK2SN205GeRlZvYBSX8q6Z8k/XbQQQIddcDvqL8ufDf5uLUjI+Znw869WuqJqAf0mhd6x53RTk92vn0FGTEA4QT5Vb27cN3l5V+a2aCkZ0vaEuRFZvZ+Se+X9HlJv+WabWENtNtsbF3ha0dGzA/Ezo9gWtK35krvuOuH0T1TZRmxJlZNbtq0SZs2bWrqvWTEADQryH82flnSH0u6TlL5f75uklcb9q/+F2a2TlKfc25b+QPM7H2SPiDpRknvdM7FvMQLiNBsrQ+T4s+IjR/xpg+tp9QVPwpnXSalBqQnC3Vi80+J5LHFGrEmMmLZbLbp91IjBqBZDX+tnHP3mtmnJF1rZjdL+qak9fI669+hE1c9fkfSanl9xyRJZvY7kj4o6TFJt0l6i504vfOkc+7WFv8dQHyKgdgs6iHmO2Wtd/T7nEXtke9J+WnpnOdJ85ZF99y+Qa9of9cPvH5i618TyWM7XiM2NiXnnGy2TWEDmLWC/mfjdZJ2ScpIepWkA5I+KW/VY6Pslt9n7BxJn6ty/g5JBGKYvfxAbDa1rvAtOdvLLB3fJ40flQYXR/t8f1ujKFZLVlrzQi8Q2/mDyAKxVtpX+B31/Q77YQz0pjTQ26OJ6bzGp/Ka158K/QwAyRSo8tY5N+Ocu8E590zn3IBzbpVz7vedc8crrltT2WLCOfeOam0oyv68OMJ/DxCt/Iz09A7v8ymzMCPWkypl6vyAMSr5fLRtKyrFUCfWytTkxo0btXHjxqbfTXd9AM2IaAkU0KUOP+Y1TF10pjSwsNOjqc6vXYs6ENt7jzTylLR4lXT6RdE+W5LO2ij1DkpP/UIaOdDy46Zn8hqZnJGZtLC//Ruzl7rrE4gBCI5ADKinuMfkLCzU9/lTpn6bjaiUT0vGUfPUO+DViUmR7Dt5fMKrD1s40KuenvbXaPkrNY+MUbAPIDgCMaCe2bxi0ue31TgYdSBW2Nbo/FdE+9xyawrbHUXQT6yV1hVRKK2cJCMGIDgCMaCeuRCIFTNiEU5NHt8vPbHVWwjg7w0ZB7+xawT7Th4e7cw+k77FhanJI+w3CSAEAjGgntm42Xel8hqxfEQt+h6+VZLzAqX+BdE8s5pVaalvvrR/mxf8teDQqLfBx7L5nQnEls3v98YxEmijEQCQRCAG1DcXasTmLfX2wJwek44+Ec0zi930Y5yWlKTefuns53qfW8yKlQKx/lZH1RQ/ADw0SkYMQHAEYkAtU2PSkd1ST6+09JxOj6a+KOvEpiekh7/jfT7vZa0/r5GIpif9qcmlTWbEtmzZoi1bAu3YVtXSQgB4eJSMGIDgCMSAWvxs2LK1Uqoz012BrXimd3zqgdafteMOaeKodPqzpFOe0frzGlkbTcF+qxmxdDrdVDNX37IFZMQAhEcgBtQyFwr1fX6fr333tf6sB77mHSPqdt/QyudIfQu8bN6xfU0/ptWMWKv8jNghMmIAQiAQA2qZzXtMVjrjYu/4ZIuB2My0tO0b3uf1r23tWUGl+qRzrvA+t9Blv9WMWCaTUSaTafr9y4pTk2TEAARHIAbUMpv3mKzkZ8T2b5NmWggEHrtTGnvaywKetj6asQWxtrDd0c7vN/0If0rQnyIMa3h4WMPDw02/v1SsT0YMQHAEYkAtc2lqcmCRtGyNtx1TKx32H7jFO65/TTzd9GvxG7u2ULDvF8kv7dCqyaVkxAA0gUAMqMa5UkAzFwIxySuul6R9P2/u/vyMdH+b68N8Z14q9S/yNlg/0lwLjk63r1g82KtUj+n4xLQmpyPq5wag6xGIAdWMHpTGD0sDi6WFp3d6NMGsKqz4e7zJFgw775CO7/NWSq7cEN24gkj1Squf531usk7s8EhharJDxfpmpqWFrv6Hx5ieBBAMgRhQTXk2rJ1TdK046zLv+Phdzd3/sy95x0t+vTP/5jWFOrFd4evEpmbyOjYxrR7r3F6TUmnFJtOTAIIiEAOqmQtbG1VatUGyHq+FxeRouHsnjpfqwy75tejHFkSxYD98nZgf+CyZ16eens4FzmxzBCAsAjGgmgMPecdT51Ag1r/AWz3pZqQ9Pw137wNfk6ZGve2G2tHEtZozLpEGlkiHH5UO7Qp16+EI6sM2bNigDRtam5It9RIjIwYgGAIxoJoDfuuKOVKo7zu70I/r0R+Fu+/uz3jH57w12vGE0ZOSnlFYPelvsRTQgeNeIHbqwuYDsVwup1wu1/T9krS88P6DIxMtPQdAchCIAdX4U5NzKSMmSete4h0f+W7we57IeX8Gl0rPemM84wrq3Ku9Y5jxS9p/3At8li8ciHpEofjv33+MQAxAMARiQKWZqcLUmM2Nrvrl1rxQspS0+y5p/Giwe+4qy4b1z49vbEGsu8o77rhDmg5eZ3WgEPisWNTZQMx//4HjBGIAgiEQAyod2iXlp6UlZ0t98zo9mnAGl3irJ92M146ikSNPSPf9uySTLntX7MNraOnZ0ooLpMljoVZ/+hmxFS1kxMxM1uJqUTJiAMIiEAMqHZiDKybLnf8K7/iLrzS+9gc3eN34L/qVzhXpV/KzYiHqxPyM2PJZkxFj1SSAYAjEgEpzsXVFuYsLdV7bvum1pajl8G5p6+clmfRL/6stQwvkXD8Quy3wLVFkxKLgF+uTEQMQFIEYUGmubW1Uaek50jnPk6bH6mfFbvuAlJ+SnnWNdNoFbRteQ6tfIPXO87ZqOvZkoFv8mqzZkhEjEAMQFIEYUMnf7HuuZsQkacNveMcffdzbQ7LSw9/xasN6B6Wr/qy9Y2ukb1Ba8wLvc8DVkweOeVOBnS7WXzjQq4HeHo1NzWhkYrqjYwEwNxCIAeWck/Y/6H1efn5nx9KKi9/oZcYObpd+9sUTzx3bJ33lt73PL/qf0rI1bR9eQ34bi+3fbnhpPu+KGbFTF3Rmw2+fmZEVAxAKQWyfigAAIABJREFUgRhQ7vhT0tjTXof3RWd2ejTNS/VJL/kT7/N//m/pyV94n488Ln3+ddLIU97eji+8vnNjrOf8V3rH7bdK0/UDmiNjU5rOOy0e7NVgX6oNg6vPXzlJCwsAQfR2egDArPLU/d7xtPVzZ7PvWi75denBb0r3/4c0fJXXY+yxH0uTx6UV66U3/pPXzX42OmWtdPrF0pP3ej3Fzn95zUv3R1Qftnnz5pbu9/kZsafIiAEIgEAMKLd/m3c8bX1nxxEFM+n1/+AVvv/8S9LDt3rfX/Bq6TUflxYs7+z4Gln/Gi8Qe+BrdQOxvUfGJUlnLB5s6XWZTKal+32nL/YCsX2FcQFAPQRiQLnyjFg36J8vvWGz9PxrpYOPeJuCz5VFCOtfLd3+EenB//QWHNTI3u09PCZJWrl0djTf9cex98hYh0cCYC6gRgwo91QXZcTKnXGxdNHr504QJkmnXSgtWyuNHpAe+0nNy/b4gdiS1jJi2WxW2Wy2pWd445hXGBcZMQCNEYgBPuekpx7wPq/oskBsLjLzpiclb3qyhj2FKcBWM2JDQ0MaGhpq6Rnl49hDRgxAAARigO/I494eh/OXSwtXdHo0kKQLX+8d77tZmqnel8vPiJ05S6Ymzyxk5vxxAUA9BGKAr5sK9bvFqg3eDgcjT0k7vlf1Er9Yf9XS1qYmo3LGkkGZeasmp2bynR4OgFmOQAzw7bvXO552YWfHgRIz6dI3eZ9/9qWTTjvnShmxJbMjI9aX6tFpiwbknPTkUerEANRHIAb49v7MO555aWfHgRNd/GvecdvXpfGjJ5x6emRSE9N5LZnXpwUDs2cReLFOjIJ9AA0QiAE+ArHZadlqafULpelxrzltGT/QObPFFZNRK62cpE4MQH0EYoAkjR2WDu2UUgPSimd2ejSodOmve8ef/ssJX+86OCJJWn3q/HaPqK5zCuPxxwcAtRCIAVKpPuz0i7x9GjG7XPQGaWCxtPsnpf9bSdp5wAt01ixf0PIrnHNyzrX8HElae6o3nl0HCMQA1EcgBkhMS852AwulS9/sfb5ruPi1H+j4gc9s4QeGOwnEADRAIAZI0r6fe8czL+nsOFDbZb/lHe+9SRo7JEnaWZj6WxtBRixKa8sCsaiybAC6E4EYIEl77vGOZMRmrxXnS+teKk2NSnd9RlIp4xRFIJZOp5VOp1t+jiQtX9ivhQO9Ojo+rUOjU5E8E0B3IhADxg5JBx6UUv3SaRd1ejSo54XXe8effFqHDx/S4dEpLehPacWigZYfvXXrVm3durXl50iSmWnNcq9gn+lJAPUQiAG77/aOKzdIfbOrDQIqrLlSOusyaexpHbvTy4qtPnWBzKzDAzvZ2uULJUk79h/v8EgAzGYEYsBjP/aO51zR2XGgMTPpyvdKkpbf82kt0JguOGNRhwdV3TNP9wKxbfuOdXgkAGYzAjHgsZ94x3Oe19lxIJjzXymddZnmTT6tTO83dOHKxZ0eUVX+uO7fc7TBlQCSjEAMyTY9IT2R8z6ffXlnx4JgzKSXfUiStCn1DV26ZHZ2r19/pheIPbDvKCsnAdREIIZk23OPNDMhrVgvzT+l06NBQDNnX6Fb3WWabxO65Bd/1enhVHXG4kEtm9+nw6NT2nuEPScBVEcghmR79Ife8ZzndnYcCOXRgyN6/8TbNKpBDTz4NenBb7X8zE2bNmnTpk0RjM5jZsXpyV8wPQmgBgIxJNvD3/WOz3hJZ8eBUH762GHt0XLdcso7vS9u+T1p5GBLz8xms8pmsxGMruRZq5ZIku7ZfSjS5wLoHgRiSK7xo97ehZaSnvHiTo8GIdy962lJ0tFL3imd83zp+D7pP/6HNMtqsS5b7U13372TQAxAdQRiSK6d35fy015fqnlLOz0ahHBXIRDb+IzTpDdkpcGl0kPfku74P00/M5fLKZfLRTVESdJla06RmXTP7sMan5qJ9NkAugOBGJJr2ze843lXd3YcCGX/sQnt2D+ieX0pb+pv6dnSr2yWZNLtH5Hu/femnrtx40Zt3Lgx0rEumd+nZ56+SJMzef388SORPhtAdyAQQzJNT0oPFgKxC1/f2bEglB89fECSlF69TH2pwk/YM18pveIj3uevDEkP3NKh0Z3suWu96ckfbN/f4ZEAmI0IxJBMO78vjR+RTrtQWn5ep0eDEG69/0lJ0lXrTzvxxBXvll7we950803vkO75QvsHV8VV60+XVBo3AJQjEEMy/azw/6TJhs0p41Mzuv3BpyRJL7vw9BNPmklXf1B6wXVeMPbVd0u3vk/Kd7Y264pnnKpFA73atu+YHjs42tGxAJh9CMSQPCMHvakr65Ge/ZZOjwYh3P7gfo1Mzmj9mYt11rL5J19gJr3sg9Kr/1bq6ZV+9HHpsy+XntrW/sEW9Pf26MUXeNm7r9+7p2PjADA7EYgheX76eWlmUjr3aq/QG3PGF+56TJJ0zYZV9S/c+JvS274qLVopPbFF2nyl9K0/lo53pk7rDYXxfvGux5TPz64WGwA6i0AMyTI1Jv34U97ny4c6OxaEsuvAiL7/0H719/bojemzGt+w9krpd34ibXi7F3j/5FPSxy+Vvn69t7VVG3uOvei8FTpr2TztfnpMd1C0D6AMgRiS5e7PSCP7pTOfLZ17VadHgxA+/p3tkqTXXbpSS+f3B7tpcIn02k9KQz+Qzvv/7Z15mF1VlbffVUNSlXkgQAgxwYAgOADKlBbhs8EWlEb6A7uhEYIgojaDkbYVhxaHbh8R2xZa2sgkQcAAoqAig4CMAomAEYWQYIQkIElISEhS8+o/9r7Jyalzh0rduqdu5fc+z3nOvXuvvfdad59z9rp7On8HnRtg/pUw5zD47n5w+2dh0Z3w+krmz5/P/PnzB0T3xgbj5IOnAfCdu5/TS8CFEJuxensgmJkDepCJvrPuJbj0AOhYDyfdCG96b94aiQr5/bK1HPs/D9HUYNzz6cOZOiFjflglvPInWHA1/H4ebHp167ixU2Gnt8CE3WD8dBi/G4yfBqMnQ8uY/prAxo4u3v3N+1j1ejv//U/7cuy+ZYZXhRBDAjMDwN0tM75Sh8bMGoBzgY8B04GVwDzgS+6+YaDTJ/KRIyb6Tk83/Oh4WHIP7Hk0nHh93hqJCtnY0cUxlzzIkpUbOONdu/GFD+zd/0x7umHZ47DoDnjx0TBU2VniMTRsNIyZDGN2CfPOxuwSvm/+vAuM2AEaSg8yXP/YC3zuJwsZ29rMr847lMljW/tvixBiUFNNR+y/gXOAW4DbgTcDZwMPAEe4e89Apk/kI0dM9A13+NXn4NHLYMTEMEw1Vr0R9UBbZzdnzl3A/YtWsvuOo/j52e+ipbmx+gX1dHPmqSdC21rmnHs0rPkzrFkKa18IPaldm8rn0dAces8KTtqYKfH7ls8+emdOm/sU9z27kj13Gs31Zx7MhJEVDrMKIeqSqjhiZrYPsBC4xd3/fyL8bOC7wD+7e9HdE/ubPpWXHDFROR0b4Bfnh33DGprhn2+EGf8vb61EBSx+ZT2fnvcUTy17jYkjhzHvrEOYMWnUgJWXeFhuHeEObWuDQ7ZuBaxfEc6FY30MTw91FqFnxCQWt43hL51jeX3YJN7x1rfwhmkzggM3chIMH73laGyutplCiBpTLUfsa8DngXe7+wOJ8BZgNfAbdz96oNKn8pIjJsrz2nL4w03w28tCQ9nUCidcBXselbdmogTr2zp59PlXueWJ5dzx9Mt09ThTxrVyxax3stfO/Z+nVYqijlildG7a4pStewnWLU98jw7b+pfB+7DBbFPrFqds2EhoboWmFmgeAc0tIb65NREe45paEuGtvcPTeZQZUhVCbDvVcsTuAI4ARrh7eyruIeBN7j5poNKn5AfcEVvx52d48bFbEyHpsjL+Mfcj3lLxno5PK9jP8irSKStNifTlbeirjmn5Xgps9a2xp5OWzjWM6FjNxA2LGdW+5XUyr4x6M/fs9e+8OmrLq4z6e/mUu/7K5V9J8eXzGGAdytnYv+R0u7O+rZO1GztZs7GDpas2snztliHAxgbjQ+/clX97316Vr5LsB/12xCqhpxtefwXWraBjzTIefnIhzy1exMSe1UzmVcbbesY2bGK0tdHqG2mkohkb/aa7YRhdDcPpbmyhu2E4XY0tic/D6bFmehqacGumxxrj5yZ6GprosSZ6rDHExe9ujfRYdO7MCE+AcKbwO4fIEBdlQrjFj5ltlhBVxxqaOPD42QOXf5UcsYXAju6+U0bcPOAEYLi7d1QzfcHpymIgH5ZP3Hkt+z38yQHLXww8672Vh3v2YV73YdzTsx+unVrqgqYG4+1Tx3HYmybxoXdOZeexLTUruyaOWAavbujgtqdWcMfTL7Nw2Wusb++KMU4LHYxmE6NsEyPZRAsdtFhnONNBq7UznPC9lXZarIOW+L3F2qNc51Zxw2O6QtwIay+pnxBDnU0+jNYLB25/v3KOWFOF+YwAit2tbQmZTEesCulrypid38ijE8u8g9C2bth7P7pTv7elf/+tv/dKb6UdB+/1bzGdn/Uqslya3jr3MX0vG1Ppe8WXcY6st03F4ntoZGPzODY0j2dNy1RWt0zDrJE9gGKv9C6jbtn/4+XTlxYol74SHcplkrcNpaLNYExrM2Pj8YYJI5g6YQTNjduX0zxh5DBOnTmdU2dOp6fHeXHNRl56rY2/rmtj7cZO2jq7aevsoa2rm86uHpzQ2+j45l7HTnc6gLWp8Ipwp8k7aOppp9k7aO5po7mnnWZv3/y5wbtp9C4avIvGeDTQHb93bx0ev4c+MIfC2ZM94775s7lnhgtRMxqaOCjH4it1xDYCOxaJa0nIVDV9lvdYqpesWsx420xmvG3mQBcjhBBb0dBgTJs4kmkTR+atihCiRlTqiK0A9jaz4ek5XsAUYFWxYckqpRdCiAFl//33z1sFIcR2SKVjAI9H2QOTgXHV475AufeC9De9EEIMKAsWLGDBggV5qyGE2M6o1BH7MWEa03mp8I8S5nb9qBBgZjPMbK9tTS+EEEIIsb3Ql531LwH+hbAz/i8JO+OfAzwEvKewM76ZLQWmped3VZq+Aj20j5gQQggh6oJqvuKokdCjdSbhXZGrCD1dX3L31xNyS8l2xCpKX4EecsSEEFUnr+0rhBBDm6o5YoMFOWJCiIFAjpgQYiAo54htXxv2CCGEEEIMIuSICSGEEELkhBwxIYQQQoickCMmhBBCCJETcsSEEEIIIXKi0lccCSHEkOb73/9+3ioIIbZDtH2FEEIIIcQAoe0rhBBCCCEGKXLEhBACmDNnDnPmzMlbDSHEdoaGJoUQAu2sL4QYGDQ0KYQQQggxSJEjJoQQQgiRE3LEhBBCCCFyQo6YEEIIIUROyBETQgghhMiJut1Zv7AKQQghqomeLUKIWqIeMSGEEEKInKi7fcRqRWK/su3u7/H2bDts3/bL9u3Tdti+7d+ebYft2/7BYLt6xIQQQgghckKOmBBCCCFETsgRE0IIIYTICTliQgghhBA5ocn6QgghhBA5oR4xIYQQQoickCMmhBBCCJETcsSEEEIIIXJCjpgQQgghRE7UvSNmZhPM7FtmttjM2sxspZnda2aHVph+hJl9ycyeNrNNZvaqmT1iZsdlyO5pZj81szVmtsHMHjCz9xTJd6yZXWJmy6NeT5vZx63KL7Krhf0WONnMbojlbDSzF8zsVjM7qEi+XuR4vZ5sj3JfLmHP+Rn5NpjZp8zsmajXi2Z2sZmNrEPbi9ldOD5foXzV6r0/9pvZrApsmpJKs4uZXRPL2GRm883shCL5Dzezr5jZn82s3cyWmNkXzKy53mw3s783s6vidbzBzFaY2d1m9r4i+S8tke8OdWh/KflLi5Rxipk9Ea+Tv5rZ5WY2qQ5tL1WXbmY/SOU/qOs+pm02s0+Y2QIzWxuP35nZuWY2LEO+Zu193b70G8DMpgH3AaOAK4BFwFjgbcCU4ik3px8P/BrYA7gK+DYwEngzMC0lOwN4GOgCvgm8BnwUuMPMjnL3uxOyw4C7gP2AS4A/AUcB3wN2Ar68bRb30r9W9g8H5gJPAjcAfwYmA2cBj5jZKe5+bUYRDwBzUmGdlVlXVvea1X2CTwGrUmELMuT+CzgHuAW4OOZ5DrCfmR3h7j3l9Cujey1t/3CRbL4MzABuy4gbsHqHftt/P9k2TSbc10+6+/JEWROAB4EdCb/TMuAkYJ6ZfcTdr0rl82PgWOBK4BHgEOCrwO7ArEptLEYtbSfU4TrgZ8CzwATgNOB2M/uCu389I69ngKzw9WV0q4ga21/gPwjP8CTPZuj2KcI18hvgXGBXYDZwiJkd6O4byuhXkhrbfl4sJ80ngYPJvu8Hc90DXE24d28GLgcagWOA7wAzgX9MlFXb9t7d6/YgPPBfBCZvY/q5hAfN3hXIzgO6gX0TYaOAvxBuSkuEfwJw4OxUHjcDHcC0erKf4LAflhG+E8Ex+SvQkIpz4OohUvdfjvZMr0B2H6AHuDkVfnbM46R6sr1I+l3jvfB4RtyA1ns17C+S5+ei7p9MhX8zhh+TCGsEHgNWA6MS4UdH2YtTeVwcw2fWme3vyZAdEZ93HcD4VNxS4L4hVPezYvjhFeSxA7AhXheNifBjYh4X1JPtRWRbgbXACqCpnuoe2CXaeUsq3GK+PcnrmRq39wP2ow30Abw7aTzQDIzoQ/rp8Ye+OH5vTD5UU7IjgTbg1xlxX4x6HJgIezDelC0p2UOj7Gfqyf4y+dwc9dg5Fe6EfyDDtiXfwWQ7CUcMGJN+CKVkvxZlD02Ft8Rr4pf1ZHuRPArX/Mcy4gas3qthf5E8DXgO2AiMS8UtAxZnpPlw1ONDibBrY9jUlOzUGP69erK9RJqCY3lwKnwpoceiCRgzBOp+VizvcGA0MKxEPmdE2Q9nxC0B/lhPtheRL1zz/5kRN6jrHtgzpv+fjLh5hJ6v1vi95u19Pc8ROzqeXzCz24BNwAYzW2RmJ1eQ/n2EOXJ/NLO5hItxvZkti13MSd5GGJ57JCOf38bzARDmBwH7A0+4e1tK9jGC531ABfqVo5b2l2JXgte/NiPu+ES+r8Qx9LF9yLsYedn+e0IXdZuZPWxmR2XIHECo48eSgfFaeJL+132u9R7nPJwW011fRGyg6h36b38WhxGGDm92983XsZlNJgx5/DYjzVb3feLzcnd/MSkYv68g/7rPItP2Muwaz69kxB1EqPvX4hycH5rZLtuoW5q87L+V0IPcZmZPFSmrULfF2oi9zCxrqK9SBkPdn05wLK4oEj+Y635JPD5iZmeY2XQzm2Fms4F/IDiXm6Jszdv7ep4jtmc8/4Dg1Z9K+PFmA3PNrNl7z9/ISv+fhOG1swgOxVnAt81snLv/e5QpXExZ8wcKYYUx6vGELtxesu7ebmarqWw8uxy1tD8TMzsaOBCYW+QivBFYTOhFOhr4F+AwM5vp7v2ZvF1r29cS5ss8DKyJ6c8DfhHnCV2dkN0FWOXu7RnlLgdmmtkwd++o2Nps3fOq9/cAuxGGH9dlxA9kvSf131b7szg9ni9Phfflvi/I/7FIGcvZ4sBsK7W0PRMzezuh4XrA3Z9PRT8d83mG0LYcTugp+ts4R2pFH3VLU2v7NwLXAfcQnM7dCHOk5prZDHe/MCFb7lqxKLOoj/oVyLXuzWx3Qq/Ub9x9cYbIoK57d+8ys78HfhjzKNBJ6GW7LBFW+/a+2l2ItTqAuwne+RISXcbxh1kDvERq3lIq/eUx/VpgYiK8mVDRm4hjxmzpkv1IRj5vjHHfid8LwxDXFCn3BcLEyLqxv0j6PQhzZJYBkyrU+YJY5ufr2fYoOzGWs4at5wktAV4okuaaWG5FQ0CD0XZCL5gD7+qDzlWp92rYn5HfOEKD+xyJeR8xrjC08JWMdA0x7qeJsG7g/iLl3A+srRfbi8hPimWvA95UYRknRZ1/UE91XyLNcGAhoQGfngj/ddStV/nAV2LcvpXqNthsJ/xxc+DkPpQxqOqe4EhfR5hEfzxwSqy3buCUhFzN2/t6HposdCNe74neBXdfQ+hK3pktXnSp9D9399WJ9J2EymohrA6BcMFCuAnTtKRkSskW5DcWiesLtbR/K8xsN7Y8eI5y95UV6nwRoffl/RXKFyM32xOyq4H/JTzQZiaiNlK67gsy20qe9T4eOA54xt0f7IPO1ap36L/9aU4i/KO9wuOTM0Ff7vvC54G872tp+1bE1aN3EXoLPujuFfXsuPt1hPlD9Vb3mXjo6f4WodfnvYmovl4rfSXPum8k9ECtBW6qtIDBVPdmtjPwOOFP8ifc/SZ3vwY4EngUuDQ+3yCH9r6eHbFl8fxyRtxL8Tw+I25b0he6VbO6GAthha7JNYSLppesmQ0n9KRkdXn2lVravxkzmw7cS1hBcqS7LyynaIHY2K8grDDqD7nYnsHSeE7aswLYIdZ1mimEYcttHZaEfG0/mfDAKTZHJJMq1jtUr+4KnE6YqHt1Rlxf7vuCfLFhiCn0/76vpe2biU7Y3cBewHHufk8fyoBwn9Rb3ZdiaTyn73sofq14QmZbyNP2ownbXPzIe09BKcdSBkfdn0loe29MBnrYSuhmwmKM/WNwzdv7enbECpOhs+ZdlJpMui3pFwLthD2B0hR6D+bD5or9HWHPqHRjfCDhN59fQq9KqaX9wOZ9XO4l7N1ypLs/UZmqm9O3xLz/2pd0GdTc9iLsEc9Jex4n1PGBScFo+770v+7ztP10wpDMNaUUTFPFeofq1R1mti/h4fsLd+/1gHf3lwgP0awewq3u+8jjwBQzm5oqZyqhJynvuk/qVNL2hNx4Qk/YPgQn7FeVq7uZ3amzui9DsfsestuIg4BnvX/zI/O0/Yx4rmgeYYrBUvcFR6kxI64pda59e9/fsdu8DoL3u47gKSfn6EwGXgcWJcLeQPg315wIayR46xuAKYnwkYRx3TXAyET4jYSx5Lcnwgr7iixi631FPknxfUU6gd3q0P5phI1c1wIHlNFtYpHwi6jC9h21tJ1wc47N0GEqYY7cKuKy5xj+VkrvI1bxHIu8bU+V+86o/80ldBvQeq+G/am8LiW1R1gJ3bP2EVsDjE6Ev5/S+4hVPK9ukNg+nrBhcTvwgTK6TSgSXngW9mvrjpzs73U9E/6EPhd/k10T4ZMIQ1CPkr2P2BfqyfaE7M6ENmtBPdc9YYPdXroQ5sY+FW3cMRFe0/a+Xz9O3gehu9GBPxBWT3w2/lAdwHsTcveRsSEnYffbLkID9NmYx0IyJuoRPPtXCd79ZwmbuD0R0/9dSnYYwQvuJDyEzwB+EvP9ar3ZT+i2fT6Gf5cwRJU+dkrI/xdh6e9/EFbjnU9YeeSE5b+tdWT7OEKDexXwGcLuyt8iOKRdwAkZul0S8/lJrPuL47VwH32YUJu37ak0l7FlTmAxvQa83qthf4xrIdzPy0k0nBlyEwmO63rgwlj2vTHf0zPkb4txlxN6EAuLI+bWoe3zYx7XkX3PvzEhe168hi4iNEznEt4s4YQVtBUt6Blk9q8gLE65gHAffy2GOTA7Q/7TMe7eqOeFBCfhT1RhT71a2p6Q/7eY18dLyAz6uies4C60YbdGPf+V4IQ58M1UWTVt7/v94+R9EJZS/5bwD389cCfwNymZUhfmYYTGYj3hH81DFPmnQHgFzM8IjfBGwkZuRxSRHUf457GC8O/pj4Rl/BWtzhlM9hM2AfUyx+EJ+WOBO+LN3hZ1e5LwQGupM9uHExrThQSHrJMwJ+EmEpv6pdI0Eh7Kz8a6X0549UnVNjit8XVf2FH7RUqvyKxJvVfJ/sKKrq9XUNYUwtsIVkW7fgf8YxHZFkKDvTTW/fOETSAzeycGs+2Uv+dnJWT/htDAvUCYM9NGcEC+QT9WCeds/8WEHsHVhPt+FfBLUg1xKs0sQuPeRhgqu5JET0u92J6Qf5bwfOg1KlBvdU94TdmlhJGdDra8CeGjZLTL1LC9t5iJEEIIIYSoMfU8WV8IIYQQoq6RIyaEEEIIkRNyxIQQQgghckKOmBBCCCFETsgRE0IIIYTICTliQgghhBA5IUdMCCGEECIn5IgJIYYsZvYWM+sysyMTYbPMzM3s8BxVw8zOM7PV8X2OQojtFDliQohBhZkdHh2lYkdXH7L7NvCQu981gPruE/X6aRm506LcBTHofwm7kH9xoHQTQgx+msqLCCFELlxPeKVMmp5KEpvZIcCRwAerqVQad3/azB4F3m9mO7r7K0VETyO8SPiHMV2bmX0fuMDMvu7uqwdSTyHE4ESOmBBisPI7d7+2H+k/QXhPYJYzV22uAA4ivAz72+lIM9sdOBT4pbsvT0RdS3g59CzCuw2FENsZGpoUQgw5zKyJ0BN2l7t3Vpjm83Ho8BIza0iEH2Fmd5rZWjNrM7Pfm9lZqeQ3EF4ifFqR7AvhVyYD3f15wouVT6hERyHE0EOOmBBisDLCzHbIOMZUkPYdwCjgsXKCZtZgZt8DvgZ8zt3PdveeGHcmcGfM6+vAbGAJcJmZXVTIw93XAzcBbzGzA9L5A6cAK4FbM1R4BHiHmY2qwC4hxBBDjpgQYrByIcF5SR/XVZB273heUkrIzFqBm4EzgFPd/RuJuMnAd4Eb3H2mu1/k7t9z9+Ni+Gwzm5HI7op4TveKHQnsClxbpHduCWGayJ4V2CWEGGJojpgQYrAyB7gxI3xlBWknxfOrJWQmAHcBbweOcfc7UvHHA8OBK8xsh1TcbcA5wN8SnT13f8DMFgEnmtlsd2+LspnDkgkKk/R3LKGrEGKIIkdMCDFYec7d797GtB7PVkLmasKQ47vd/cGM+DfHcykddkp9vxL4BnAccH3cI+yDwGPu/ocieRR09CLxQoghjBwxIcRQpNBrNqGEzI8JvVVfMrNj3X1TKr7gIJ0CvFQkj+dT339ImGt2GmH7jZOIvWol9CjoWElPnxBiiCFOkE/0AAABxklEQVRHTAgxFCn0Pu1RQuZHwK+BucDPzewYd9+YiH8unldV2jPn7i+b2e2EPcWmEhyyjYRVlcXYHegirJ4UQmxnaLK+EGIo8gSwDji4lJC73wCcSNjj6/bUysV5QDtwYZzUvxVmNtbMhmdkewXh2XoxYfXmTe6+roQaBwML3P31UroKIYYm6hETQgxW9jezk4vE/bSU4+Lu3Wb2E+BYMxvu7u0lZG8ys06C43WHmR3l7uvcfZmZfRy4HPiTmc0F/kJYCPBWwtyvvYGlqSx/AbzMlr3Bik3SJ6663BM4v5iMEGJoI0dMCDFYOTEeWewBLC6T/jLCjvUfIGxRURR3/5mZ/UOUu9PM3ufua939qrgS8nzgY8A4YBVhGPGLBIcrnVeXmV0DfIawovL+EkWfTOh1u7qMLUKIIYq5a6GOEGJoYma/Aka6+6F565LGzFoIk/1vcPfZeesjhMgHzRETQgxlPg0cYmbvzVuRDM4CWoCv5q2IECI/1CMmhBBCCJET6hETQgghhMgJOWJCCCGEEDkhR0wIIYQQIifkiAkhhBBC5IQcMSGEEEKInJAjJoQQQgiRE3LEhBBCCCFyQo6YEEIIIURO/B8VxFuhvU5EtAAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 720x720 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "plt.rc(\"font\", size=18)\n",
    "plt.rc(\"axes\", linewidth=2)\n",
    "fig, ax = plt.subplots(figsize=(10,10))\n",
    "ax.plot(emid, f0_conv*1.0e12, lw=2)\n",
    "ax.plot(emid, f_conv*1.0e12, lw=2)\n",
    "ax.tick_params(width=2, length=6)\n",
    "ax.set_xlabel(\"E (keV)\")\n",
    "ax.axvline(6.7, ls='--', color='k', lw=2)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(1000,)"
      ]
     },
     "execution_count": 14,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "instrument.shape"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(1000,)"
      ]
     },
     "execution_count": 15,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "f.shape"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "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.8.1"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}