{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Ultrarelativistic right-handed neutrino production"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### In this notebook we reproduce the $2\\leftrightarrow 2$ part of [1202.1288](https://arxiv.org/pdf/1202.1288)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Importing the modules"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "#start by importing the controller and manipulate modules\n",
    "from analytical.controller import *\n",
    "from numerical.manipulate import *\n",
    "#reimport numpy (though it is pulled by numerical) for smarter syntax highlighting in vscode\n",
    "import numpy\n",
    "#import matplotlib too\n",
    "import matplotlib.pyplot as plt"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Field-theoretical component"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The relevant part of the config file is as follows\n",
    "```ini\n",
    "[Model]\n",
    "modelpath = /path/to/autotherm/analytical/models/rhn.fr\n",
    "# Symbol for the Lagrangian in the model file\n",
    "lagrangian = Ltot\n",
    "# \"Name\" of the particle whose production rate must be computed\n",
    "produced = F[6]\n",
    "# List of the particles in the thermal bath (or leave empty for SM assumption)\n",
    "inbath = \n",
    "assumptions = Element[hhh,Reals]\n",
    "replacements =\n",
    "includeSM = yes\n",
    "noneq = hhh\n",
    "#comma-separated list of particles to be treated with a flavor expansion\n",
    "flavorexpand = \n",
    "```"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "rhn_analytical=analytical_pipeline(\"../../MyModels/RHN/RHN.cfg\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Inspect the output. The `hhh` coupling in the model file stands for $\\sqrt{\\mathrm{Tr}h^\\dagger h}$, with $h$ the Yukawa matrix between the right-handed neutrinos and the left-handed lepton doublet. In this models there are three generations of massless right-handed neutrinos."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/markdown": [
       "-  $$ \\left\\vert \\mathcal{M}(-1;-1,-1)\\right\\vert^2=36 \\mathrm{Tr}[h^\\dagger h] \\vert h_t\\vert^2$$\n",
       "-  $$ \\left\\vert \\mathcal{M}(1;-1,1)\\right\\vert^2=- \\frac{\\mathrm{Tr}[h^\\dagger h] \\left(g_{1}^{2} + 3 g_{2}^{2}\\right) \\left(s^{2} + t^{2}\\right)}{s t}$$\n",
       "-  $$ \\left\\vert \\mathcal{M}(1;1,-1)\\right\\vert^2=- \\frac{\\mathrm{Tr}[h^\\dagger h] \\left(g_{1}^{2} + 3 g_{2}^{2}\\right) \\left(s^{2} + u^{2}\\right)}{s u}$$\n",
       "-  $$ \\left\\vert \\mathcal{M}(-1;1,1)\\right\\vert^2=\\frac{\\mathrm{Tr}[h^\\dagger h] \\left(g_{1}^{2} + 3 g_{2}^{2}\\right) \\left(t^{2} + u^{2}\\right)}{t u}$$\n"
      ],
      "text/plain": [
       "<IPython.core.display.Markdown object>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "from IPython.display import display, Markdown\n",
    "out=\"\"\n",
    "s = sympy.Symbol(\"s\")\n",
    "t = sympy.Symbol(\"t\")\n",
    "u = sympy.Symbol(\"u\")\n",
    "g1 = sympy.Symbol(\"g1\")\n",
    "g2 = sympy.Symbol(\"g2\")\n",
    "g3 = sympy.Symbol(\"g3\")\n",
    "ht = sympy.Symbol(\"ht\")\n",
    "hhh=sympy.Symbol(\"hhh\")\n",
    "for key, item, in rhn_analytical[0].items():\n",
    "    sitem=sympy.simplify(sympy.sympify(item).subs(2*t*u,(s*s-t*t-u*u)))\n",
    "    sitemsubst= sitem\n",
    "    if sitemsubst.as_coeff_Add()[0] ==0:\n",
    "        sitemsubstcollect = sitemsubst\n",
    "    else:\n",
    "        sitemsubstcollect=0\n",
    "        for exprtemp in sympy.factor(sitemsubst.expand().as_independent(ht)).subs(t+u,-s)\\\n",
    "            .subs(2*t*t+2*t*u,2*t*t+(s*s-t*t-u*u)).subs(t*t+2*t*u,t*t+(s*s-t*t-u*u)).subs(t*t+t*u,t*t+(s*s-t*t-u*u)/2):\n",
    "            sitemsubstcollect += exprtemp.simplify()\n",
    "        # beautify the output\n",
    "    out+=f\"-  $$ \\\\left\\\\vert \\\\mathcal{{M}}({key[0]};{key[1]},{key[2]})\\\\right\\\\vert^2={sympy.latex(sitemsubstcollect).replace('ht^{2}',r'\\vert h_t\\vert^2').\\\n",
    "                                    replace('hhh^{2}',r'\\mathrm{Tr}[h^\\dagger h]')}$$\\n\"\n",
    "    # display(sympy.pprint(sitem))\n",
    "display(Markdown(out))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Numerical part"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### First, we must call NumRate"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We set the degeneracy to unity, since we are treating the right-handed neutrinos as massless and thus distinct from their antiparticles. The rate for the latter is then equal to that for the former that is determined here"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "rate_rhn=NumRate(*rhn_analytical,1) "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The leading-log term reads"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/latex": [
       "$\\displaystyle \\frac{T^{2} hhh^{2} \\left(g_{1}^{2} + 3 g_{2}^{2}\\right) \\log{\\left(\\frac{64 k^{2}}{T^{2} \\left(g_{1}^{2} + 3 g_{2}^{2}\\right)} \\right)}}{256 \\pi k \\tanh{\\left(\\frac{k}{2 T} \\right)}}$"
      ],
      "text/plain": [
       "T**2*hhh**2*(g1**2 + 3*g2**2)*log(64*k**2/(T**2*(g1**2 + 3*g2**2)))/(256*pi*k*tanh(k/(2*T)))"
      ]
     },
     "execution_count": 5,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "leadlog=rate_rhn.get_leadlog().simplify()\n",
    "leadlog"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "reexpress this in terms of the lepton asymptotic mass"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/latex": [
       "$\\displaystyle \\frac{\\mathrm{Tr}[h^\\dagger h] m_{\\infty\\,L}^2 \\left(2 n_{B}{\\left(k \\right)} + 1\\right) \\log{\\left(\\frac{4 k^{2}}{m_{\\infty\\,L}^2} \\right)}}{16 \\pi k}$"
      ],
      "text/plain": [
       "\\mathrm{Tr}[h^\\dagger h]*m_{\\infty\\,L}^2*(2*n_B(k) + 1)*log(4*k**2/m_{\\infty\\,L}^2)/(16*pi*k)"
      ]
     },
     "execution_count": 6,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "minf = sympy.Symbol(r\"m_{\\infty\\,L}^2\")\n",
    "htrace=sympy.Symbol(r'\\mathrm{Tr}[h^\\dagger h]')\n",
    "T = sympy.S(\"T\")\n",
    "k= sympy.S(\"k\")\n",
    "nB= sympy.Function(\"n_B\")\n",
    "leadlog.subs(g2*g2,16*(minf/T/T-g1*g1/16)/3).subs(sympy.tanh(k/(2*T)),1/(1+2* nB(k))).subs(hhh*hhh,htrace).simplify()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "This agrees with Eqs. (24) and (28) of of [1202.1288](https://arxiv.org/pdf/1202.1288) and with (4.23) of [1605.07720](https://arxiv.org/abs/1605.07720)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now prepare to compare with the integration results in (29) of of [1202.1288](https://arxiv.org/pdf/1202.1288)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [],
   "source": [
    "#prepare to compute the total rate in the strict LO scheme\n",
    "def ratevec(x,g1,g2,ht):\n",
    "    #add factor of 2 for the \"antiparticle\"\n",
    "    #factor out factor of 2/(2\\pi^2)\n",
    "    statfac=x*x/(numpy.exp(x)+1.)\n",
    "    evalf=rate_rhn.rate(x,1,tuple([g1,g2,ht]),1)[1][0]\n",
    "    return statfac*evalf"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "we first test $c_Q$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/markdown": [
       "### We find $c_Q$=2.52136"
      ],
      "text/plain": [
       "<IPython.core.display.Markdown object>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "from scipy.integrate import quad\n",
    "#re-instate factor of 2/(2\\pi^2)\n",
    "yukawaint=1/numpy.pi**2*quad(ratevec,0.,20.,args=(1e-7,1e-7,1),epsrel=1e-5)[0]\n",
    "display(Markdown(r'### We find $c_Q$=%.5f'%(yukawaint*1536*numpy.pi)))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "It checks out! Let us then test $c_V$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/markdown": [
       "### We find $c_V$=3.16693"
      ],
      "text/plain": [
       "<IPython.core.display.Markdown object>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "g1test=1\n",
    "g2test=1\n",
    "gaugeint=1/numpy.pi**2*quad(ratevec,0.,20.,args=(g1test,g2test,0.),epsrel=1e-5)[0]\n",
    "display(Markdown(r'### We find $c_V$=%.5f'%(gaugeint*1536*numpy.pi/(g1test**2+3*g2test**2)+numpy.log(g1test**2+3*g2test**2))))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Checks out as well!"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Plot the rate at $T=200$ GeV. Couplings chosen along the procedure of [2004.11392](https://arxiv.org/abs/2004.11392)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "k=numpy.geomspace(0.01,15,100)\n",
    "coups=(numpy.sqrt(0.13031),numpy.sqrt(0.406827),numpy.sqrt(0.718285))\n",
    "rateplt=rate_rhn.rate(k,1,coups,0)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAjAAAAGdCAYAAAAMm0nCAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjkuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8hTgPZAAAACXBIWXMAAA9hAAAPYQGoP6dpAABVKklEQVR4nO3deXhU5f3//+csmeyZQPZAICxhkxj2CFJRQVFRi1ZFPioUrf7aogVxVxBtVUTFYhWl2K9LqyjiglYQxVREAdnDIsi+BMhCgOz7zPz+GBKNBMjAhMkcXo/rmivJmfucec/Rel697/uc2+RyuVyIiIiI+BGzrwsQERER8ZQCjIiIiPgdBRgRERHxOwowIiIi4ncUYERERMTvKMCIiIiI31GAEREREb+jACMiIiJ+RwFGRERE/I4CjIiIiPid0wowM2bMIDk5maCgINLT01m5cuUJ2/7444/87ne/Izk5GZPJxPTp08/4mCIiInJu8zjAzJkzhwkTJjB58mTWrl1LWloaQ4cOJS8vr8H2ZWVltG/fnmeffZb4+HivHFNERETObSZPF3NMT0+nb9++vPLKKwA4nU6SkpK45557ePjhh0+6b3JyMuPHj2f8+PFeO6aIiIice6yeNK6qqmLNmjU88sgjddvMZjNDhgxh+fLlp1XA6RzT6XSyZ88eAgICMJlMddsDAwMJDAw8rTpERETEuyorK6msrKz72+VyUV1dTXJyMmbzmU3D9SjA5Ofn43A4iIuLq7c9Li6On3766bQKOJ1j7tmzhw4dOpzW54mIiIhv7dy5k/bt25/RMTwKMM1FQEAAACtXriQhIaFue20PzJw5cNddacAeFi1aRL9+/XxUafNWVFREUlISWVlZRERE+Locv6Zz6T06l96h8+g9Open79c9MNnZ2fTr16/uOn4mPAow0dHRWCwWcnNz623Pzc094QTdpjhm7bBRQkICrVu3Pu797t0hIiKYoiKwWq36F+4UIiIidI68ROfSe3QuvUPn0Xt0Lr3nl9M/TpdHA1A2m43evXuTkZFRt83pdJKRkUH//v1Pq4CmOGb//tC+vXsuzC+Tn4iIiBiDx0NIEyZMYPTo0fTp04d+/foxffp0SktLGTNmDACjRo2iVatWTJkyBXBP0t28eXPd7wcOHCAzM5OwsDA6duzYqGOejtrJvAowIiIixuNxgBkxYgSHDh3i8ccfJycnhx49erBw4cK6Sbj79u2rN7P44MGD9OzZs+7vF154gRdeeIFBgwaxePHiRh3z12rDycnuOFKAObXAwEAmT56sO7e8QOfSe3QuvaO5nUen00lVVZWvyzhtzz//PAAVFRU+rqR5CwgIwGKxnPD9xly/G8vj58A0B0VFRdjtdgoLCxscj1y/Hvr2vZzq6kX85z//4dZbb/VBlSIiAu7e9927d+N0On1dipwFkZGRxMfHNzjP5VTXb0/45V1Ip2KxQHW1emBERHzN5XKRnZ2NxWIhKSnpjJ/9Ic2Xy+WirKys7in6v7xLuCkYMsB07AiXXRbIokUKMCIivlRTU0NZWRmJiYmEhIT4uhxpYsHBwQDk5eURGxt70uGkM2XIKBwUBDEx6oEREfE1h8MBuO84lXNDbVCtrq5u0s8xZIABTeIVEWlOvPHcD/EPZ+uftSEDjMMBu3a5A0xpqQKMiIiI0RgywJjN8O237gBTVKQAIyIiYjSGDDAmE3Tq5A4wVVW6Z19ERHzjrbfeIjIy0tdlGJIhAwzAzTcHAWA2qwdGRES85/e//z3Dhw9vVNsRI0awbdu2Rh87OTmZ6dOnn3G7rKwsbr/9dhITE7HZbLRt25Zx48Zx+PDhRtfS3Bk2wAQFuQOMnpooIiK+UF1dTXBwMLGxsWf1c3ft2kWfPn3Yvn077733Hjt27GDmzJl1awweOXLkrNbTVAwbYGrvRVeAERFphkpLT/z69X+3T9a2vLxxbT304YcfkpqaSnBwMFFRUQwZMoTS0lKeeOIJ3n77bT799FNMJhMmk4nFixezZ88eTCYTc+bMYdCgQQQFBfHuu+82OIT03//+l759+xIUFER0dDTXXXcdABdffDF79+7l3nvvrTv26Rg7diw2m42vvvqKQYMG0aZNG6688kq+/vprDhw4wGOPPXZax21uDBtg/vlPdw/MgQMKMCIizU5Y2Ilfv/td/baxsSdue+WV9dsmJzfczgPZ2dmMHDmS22+/nS1btrB48WKuv/56XC4X999/PzfddBNXXHEF2dnZZGdnM2DAgLp9H374YcaNG8eWLVsYOnTocceeP38+1113HVdddRXr1q0jIyODfv36AfDxxx/TunVr/vrXv9Yd21NHjhzhyy+/5M9//nPd/5GvFR8fzy233MKcOXPww1WEjmPIJ/EClJe7A0xZmQKMiIg0XnZ2NjU1NVx//fW0bdsWgNTU1Lr3g4ODqaysJD4+/rh9x48fz/XXX3/CYz/99NPcfPPNPPnkk3Xb0tLSAGjZsiUWi4Xw8PAGj90Y27dvx+Vy0bVr1wbf79q1K0ePHuXQoUNnfWjL2wwbYMaODeLBByE4uPzUjUVE5OwqKTnxe79+/PyxtXUa9Ou1lfbsOe2SaqWlpTF48GBSU1MZOnQol19+OTfccAMtWrQ45b59+vQ56fuZmZnceeedZ1zjqRihh+VUDDuElJLi7oGprlYPjIhIsxMaeuLXsZswGtX2V8MkJ2znAYvFwqJFi/jiiy/o1q0bL7/8Mp07d2b37t2N+Fon/6xfD+t4W8eOHTGZTGzZsqXB97ds2UKLFi2IiYlp0jrOBsMGGN2FJCIip8tkMnHhhRfy5JNPsm7dOmw2G5988gngXtepdo0nT51//vlkZGSc8P0zOTZAVFQUl112Ga+++irlv5rgnJOTw7vvvsuIESMMsbSDYQPM9u3uAFNYqAAjIiKNt2LFCp555hlWr17Nvn37+Pjjjzl06FDdvJLk5GQ2bNjA1q1byc/P92jRwsmTJ/Pee+8xefJktmzZwsaNG5k6dWrd+8nJySxZsoQDBw6Qn59/0mMdOHCAzMzMeq+jR4/yyiuvUFlZydChQ1myZAlZWVksXLiQyy67jFatWvH000+f3olpZgwbYBYscAeYggIFGBERabyIiAiWLFnCVVddRadOnZg4cSLTpk3jymN3PN1555107tyZPn36EBMTw9KlSxt97Isvvpi5c+fy2Wef0aNHDy699FJWrlxZ9/5f//pX9uzZQ4cOHU45zPPCCy/Qs2fPeq/58+eTkpLC6tWrad++PTfddBMdOnTgrrvu4pJLLmH58uW0bNny9E5MM2Ny+eFMn6KiIux2O4WFhURERDTY5v771zNtWg9atEjgyJGDZ7lCEREB9zD+7t27adeuXd3Qvhjbyf6ZN+b63ViG7YG5887ak6YeGBEREaMxbICpTX2/nsQkIiIi/s/wAaaiouKcuB9eRETkXGLYAPPhhz+Pu1VVVfmwEhEREfE2wwaYqqqfA4yeBSMiImIshg0wN95oq/tdAUZERMRYDBtgWrc26Wm8IiIiBmXYAAM/rzmhACMiImIshg0whw4BqAdGRETEiAwbYH76CY4e1bNgRETk7EtOTmb69Om+LqNJPPHEE/To0cPXZRg3wERHQ0iIemBERKR5MplMzJs376x8VnMJHd5k2ADTtSt07qwAIyIi/kvPMTsxwwYY0HICIiLNVWmp+/XLB6VXVbm3VVY23Nbp/HlbdbV726///+mJ2p6ODz/8kNTUVIKDg4mKimLIkCGUlpZy8cUXM378+Hpthw8fzu9///t624qLixk5ciShoaG0atWKGTNm1L2XnJwMwHXXXYfJZKr7u7an5F//+le9xRAXLlzIwIEDiYyMJCoqiquvvpqdO3fW+7z9+/czcuRIWrZsSWhoKH369GHFihW89dZbPPnkk6xfvx6TyYTJZOKtt94CoKCggD/84Q/ExMQQERHBpZdeyvr16+sd99lnnyUuLo7w8HDuuOOOZtMpYOgAo7uQRESap7Aw9ys//+dtzz/v3nb33fXbxsa6t+/b9/O2GTPc2+64o37b5GT39i1bft527FrtkezsbEaOHMntt9/Oli1bWLx4Mddff71HS9M8//zzpKWlsW7dOh5++GHGjRvHokWLAFi1ahUAb775JtnZ2XV/A+zYsYOPPvqIjz/+mMzMTABKS0uZMGECq1evJiMjA7PZzHXXXYfzWFIrKSlh0KBBHDhwgM8++4z169fz4IMP4nQ6GTFiBPfddx/nnXce2dnZZGdnM2LECABuvPFG8vLy+OKLL1izZg29evVi8ODBHDlyBIAPPviAJ554gmeeeYbVq1eTkJDAq6++6vkJbQJWXxfQVMrKYPVqd4ApKCjzcTUiIuJPsrOzqamp4frrr6dt27YApKamenSMCy+8kIcffhiATp06sXTpUv7+979z2WWXERMTA0BkZCTx8fH19quqquLf//53XRuA3/3ud/XavPHGG8TExLB582a6d+/O7NmzOXToEKtWraJly5YAdOzYsa59WFgYVqu13md9//33rFy5kry8PAIDAwF44YUXmDdvHh9++CF33XUX06dP54477uCOY0nxqaee4uuvv24WHQOG7YGx2aCoKASAoiINIYmINCclJe5XdPTP2x54wL3tlVfqt83Lc29v0+bnbWPHurf9v/9Xv+2ePe7tXbv+vO1XIzuNkpaWxuDBg0lNTeXGG2/k9ddf5+jRox4do3///sf9veWXXUMn0LZt23rhBWD79u2MHDmS9u3bExERUTfktO9Yt1RmZiY9e/asCy+NsX79ekpKSoiKiiIsLKzutXv37rrhqS1btpCenn7S7+Urhu2BsVphyJBgvv4anE4FGBGR5iQ09PhtNpv71Zi2AQHuV2PbespisbBo0SKWLVvGV199xcsvv8xjjz3GihUrMJvNxw0lVZ/uRJsGhDbwJa655hratm3L66+/TmJiIk6nk+7du9dN8q2dMuGJkpISEhISWLx48XHvRUZGeny8s82wPTAAHTq4/4FWVWkISUREPGMymbjwwgt58sknWbduHTabjU8++YSYmBiys7Pr2jkcDjZt2nTc/j/88MNxf3f9RddQQEAADofjlHUcPnyYrVu3MnHiRAYPHkzXrl2P6w06//zzyczMrJu78ms2m+24z+rVqxc5OTlYrVY6duxY7xV9rGusa9eurFix4qTfy1cMHWBqE6nuQhIREU+sWLGibuLqvn37+Pjjjzl06BBdu3bl0ksvZf78+cyfP5+ffvqJP/3pTxQUFBx3jKVLl/Lcc8+xbds2ZsyYwdy5cxk3blzd+8nJyWRkZJCTk3PS4akWLVoQFRXFrFmz2LFjB//73/+YMGFCvTYjR44kPj6e4cOHs3TpUnbt2sVHH33E8uXL6z5r9+7dZGZmkp+fT2VlJUOGDKF///4MHz6cr776ij179rBs2TIee+wxVq9eDcC4ceN44403ePPNN9m2bRuTJ0/mxx9/9MIZPnOGDjD5+e45MEeOKMCIiEjjRUREsGTJEq666io6derExIkTmTZtGldeeSW33347o0ePZtSoUQwaNIj27dtzySWXHHeM++67j9WrV9OzZ0+eeuopXnzxRYYOHVr3/rRp01i0aBFJSUn07NnzhLWYzWbef/991qxZQ/fu3bn33nt5/vnn67Wx2Wx89dVXxMbGctVVV5Gamsqzzz6LxWIB3JOAr7jiCi655BJiYmJ47733MJlMLFiwgIsuuogxY8bQqVMnbr75Zvbu3UtcXBwAI0aMYNKkSTz44IP07t2bvXv38qc//ckbp/iMmVye3BPWTBQVFWG32yksLCQiIuKE7eLiniIvbxLDht3J55/POosViogIuB9jsXv37nrPNBFjO9k/88ZevxvD0D0wrVu7h5BqajQHRkRExEgMHWD+8Af3EFJoqIaQREREjMTQAUaTeEVERIzpnAgwZWUaQhIRETESQweYzz93B5h9+9QDIyIiYiSGDjAFBe45MKWlCjAiIiJGYugA87vfuXtggoMVYERERIzE0AHm/PPdAaa6WnNgREREjMTQASYkxD2EpLuQREREjMXQAaa4uPYuJAUYERHxT4sXL8ZkMjW43tK5zNABZulSd4CprKzA6XT6uBoREfEXF198MePHj/d1GXIShg4wtUsJgHttBhERETEGQweY4cN/DjCaByMi4nsul4vS0lKfvBq7dvHvf/97vv32W1566SVMJhMmk4m33nqLyMjIeu3mzZuHyWSq+/uJJ56gR48e/Oc//yE5ORm73c7NN99McXFxXRun08mUKVNo164dwcHBpKWl8eGHH9Y77oIFC+jUqRPBwcFccskl7Nmz57TPt5FZfV1AU7JarQQEBFBdXa0AIyLSDJSVlREWFuaTzy4pKSE0NPSU7V566SW2bdtG9+7d+etf/wrA/PnzG/UZO3fuZN68eXz++eccPXqUm266iWeffZann34agClTpvDOO+8wc+ZMUlJSWLJkCbfeeisxMTEMGjSIrKwsrr/+esaOHctdd93F6tWrue+++07/SxuYoQMMuO9EKiwspLS01NeliIiIH7Db7dhsNkJCQoiPjwfAYrE0al+n08lbb71FeHg4ALfddhsZGRk8/fTTVFZW8swzz/D111/Tv39/ANq3b8/333/PP//5TwYNGsRrr71Ghw4dmDZtGgCdO3dm48aNTJ06tQm+qX8zdIDZswfKy0OAQq2HJCLSDISEhFBSUuKzz25qycnJdeEFICEhgby8PAB27NhBWVkZl112Wb19qqqq6NmzJwBbtmwhPT293vu1YUfqM3SAcTigqsrdXagAIyLieyaTqVHDOM2N2Ww+bg5NdXX1ce0CAgLq/W0ymerugq0NbvPnz6dVq1b12gUGBnqz3HPCaU3inTFjBsnJyQQFBZGens7KlStP2n7u3Ll06dKFoKAgUlNTWbBgQb33S0pKuPvuu2ndujXBwcF069aNmTNnnk5p9SQkQNu27v+haAhJREQay2az4XA46v6OiYmhuLi43rUkMzPTo2N269aNwMBA9u3bR8eOHeu9kpKSAOjatetx19Qffvjh9L+IgXkcYObMmcOECROYPHkya9euJS0tjaFDh9Z1kf3asmXLGDlyJHfccQfr1q1j+PDhDB8+nE2bNtW1mTBhAgsXLuSdd95hy5YtjB8/nrvvvpvPPvvs9L8ZEBICiYm1CzoqwIiISOMkJyezYsUK9uzZQ35+Punp6YSEhPDoo4+yc+dOZs+ezVtvveXRMcPDw7n//vu59957efvtt9m5cydr167l5Zdf5u233wbgj3/8I9u3b+eBBx5g69atp/U55wqPA8yLL77InXfeyZgxY+p6SkJCQnjjjTcabP/SSy9xxRVX8MADD9C1a1f+9re/0atXL1555ZW6NsuWLWP06NFcfPHFJCcnc9ddd5GWlnbKnp3GqO2q1BCSiIg01v3334/FYqFbt27ExMRQVFTEO++8w4IFC0hNTeW9997jiSee8Pi4f/vb35g0aRJTpkyha9euXHHFFcyfP5927doB0KZNGz766CPmzZtHWloaM2fO5JlnnvHytzMGk6uxN8bjnmgUEhLChx9+yPDhw+u2jx49moKCAj799NPj9mnTpg0TJkyo90TDyZMnM2/ePNavXw/AXXfdxbp165g3bx6JiYksXryYa6+9lvnz53PRRRcdd8yioiLsdjtZWVlERETUbQ8MDDxuHHHQoOEsWfIpL7/8T+6++67GflUREfGCiooKdu/eTbt27QgKCvJ1OXIW/PKfuclkorKysu69oqIikpKSKCwsrHf9Ph0e9cDk5+fjcDiIi4urtz0uLo6cnJwG98nJyTll+5dffplu3brRunVrbDYbV1xxBTNmzGgwvPxSUlISdru97jVlypTj2ixf7h5C2r9fQ0giIiJn05QpU+pdp2vn+nhDs7gL6eWXX+aHH37gs88+o23btixZsoSxY8eSmJjIkCFDTrhfQz0wv9aiRSh5eVBeriEkERGRs+mRRx5hwoQJdX/X9sB4g0cBJjo6GovFQm5ubr3tubm5dQ/7+bX4+PiTti8vL+fRRx/lk08+YdiwYQCcf/75ZGZm8sILL5w0wERERJyyC2rkyFBeeglCQ9UDIyIicjY1NLXDWzwaQrLZbPTu3ZuMjIy6bU6nk4yMjBM+aKd///712gMsWrSorn11dTXV1dWYzfVLsVgsXllBunYSr+5CEhERMQ6Ph5AmTJjA6NGj6dOnD/369WP69OmUlpYyZswYAEaNGkWrVq3q5qOMGzeOQYMGMW3aNIYNG8b777/P6tWrmTVrFuDuRRk0aBAPPPAAwcHBtG3blm+//ZZ///vfvPjii2f8BWufvKi7kEREfMeD+0XEz3mj86ExPA4wI0aM4NChQzz++OPk5OTQo0cPFi5cWDdRd9++ffV6UwYMGMDs2bOZOHEijz76KCkpKcybN4/u3bvXtXn//fd55JFHuOWWWzhy5Aht27bl6aef5o9//OMZf8GlS909MNu2qQdGRORsCwgIwGQycejQIWJiYuqt3izG4nK5qKqq4tChQ5jNZmw2W5N+nke3UTcXtbdRN+Y2rL59X2f16rvo3v1aNm48/jZvERFpWiUlJezfv1+9MOeIkJAQEhISGgwwnly/T6VZ3IXUlK64IoTVqyEiQkNIIiK+EBYWRkpKSoNrB4mxWCwWrFbrWelpM3yA6d3bPYTkcmkISUTEVywWCxaLxddliIGc1mKO/kR3IYmIiBiP4QNMZaX7LqSCAg0hiYiIGIXhA8x337l7YPLy1AMjIiJiFIYPMHFx7gDjdKoHRkRExCgMH2BGjHAPITmd6oERERExCsMHmLCwMABqamrqLektIiIi/svwAab2LiTQnUgiIiJGYfgAs3+/FbM5CIDi4mIfVyMiIiLeYPgAU1kJTqd7GKmkpMTH1YiIiIg3GDfA7NgBBQUkJkJUlAKMiIiIkRgzwGzfDl27wsUXEx7mIjExHFCAERERMQpjBpiMDKipgfXrISOj7k4kzYERERExBmMGmO3bf/79008xm90B5sgR9cCIiIgYgTEDzC9XPN2xg9Wr3QFm1y4FGBERESMwZoB57jn4+mv37zt3EhqqSbwiIiJGYswAA5CeDhs3wrp13HyzexKv3a45MCIiIkZg9XUBTSYsDLp3P/aremBERESMxHg9MEVFcN550L8/VFcDCjAiIiJGY7wemMJC2LwZAgJgzhz4/ns2/NQCgC1bFGBERESMwJg9MAAREfDll/DPf5K3zwZAXp7mwIiIiBiBcQOM3Q4JCQBcFJsFQHi4emBERESMwLgBJiIC4uMBSA/cCYDJpAAjIiJiBMYOMMd6YMIKCwEtJSAiImIUxgswx8IKdntdD4zlaAUABQXqgRERETEC4wUYqxVat3aHl2M9MKtzuwGQn68AIyIiYgTGu4369793vwAKCgCIqTwCgMNRjMvlwmQy+aQ0ERER8Q7jBZhfstthyxaGh4dD69a4XA7Ky8sJCQnxdWUiIiJyBowdYEwm6NKF0GO9Li6Xi6KiIgUYERERP2e8OTBTp7qXEXjzzbpNJpOJiIgIAIpq71ISERERv2W8ALN9O/zwA2Rnu//+z38ou+NuaipDAThyRAFGRETE3xkvwFS4b5kmMND9c9EirG/MorTCvR5STo4CjIiIiL8zXoCprHT/rA0wsbHYqKZNuAOAigoFGBEREX9n/AATEwNA13B3cKmqUoARERHxd8YNMEFB7p/HAkyE0wlAYe2TekVERMRvGS/A/HoOTGys+88q98Pr8vLUAyMiIuLvjBdgQkPdD7CrfdbLsR6YJYUXAbBunQKMiIiIvzPeg+wWLKj/97EAE+YeQaKsTAFGRETE3xmvB+bX2rSB7du546nzAUhIUIARERHxd8YPMFYrdOyIPS4O0JN4RUREjMB4AeaGG+Cyy2DnznqbtZSAiIiIcRhvDsy330J+/s93IwG8+ip73l4HwM6dCjAiIiL+zngB5tcPsgOYN4+jK20AHD2qACMiIuLvjDeE1FCAiY6mL1sBsFoVYERERPydsQKMywVVVe7ffxVg+rADgPLyQlwulw+KExEREW8xVoCp7X2B+gEmJobIY79WV1dT8cv5MSIiIuJ3jBVgantf4LgemBDAjHs5gYKCgrNaloiIiHiXsQJMdbV7CQGbDQICft4eE0M+sThpAcDhw0d9VKCIiIh4g7ECTFQUlJa6h5Islp+3R0cTSQEcCzDZ2QW+qE5ERES8xFgB5kT69SNw+4/0SgsHoLq6wLf1iIiIyBk5NwJMSAh07EiL6ChAc2BERET8nbECzIEDMGwYjBzZ4NuRkZGAAoyIiIi/M9aTeIuKYMECaNHi+PeeeYbs7/IB2LCh4OzWJSIiIl5lrB6Ymhr3T2sDuezf/yYnzx1sdu/WXUgiIiL+zJgB5pe3UNeKjqYP7hWqQ0MLzl5NIiIi4nWnFWBmzJhBcnIyQUFBpKens3LlypO2nzt3Ll26dCEoKIjU1FQWLFhwXJstW7Zw7bXXYrfbCQ0NpW/fvuzbt8+zwk7WAxMTw0VsBMBiKfDsuCIiItKseBxg5syZw4QJE5g8eTJr164lLS2NoUOHkpeX12D7ZcuWMXLkSO644w7WrVvH8OHDGT58OJs2baprs3PnTgYOHEiXLl1YvHgxGzZsYNKkSQQFBXlW3CkCTOSxXzWJV0RExL+ZXB6ubJienk7fvn155ZVXAHA6nSQlJXHPPffw8MMPH9d+xIgRlJaW8vnnn9dtu+CCC+jRowczZ84E4OabbyYgIID//Oc/jaqhqKgIu91OYWEhERERP7/x3Xdw0UXQqRNs3Vp/p0cf5bMpU/ktTtLS+pCZucqTry0iIiJn6ITX79PgUQ9MVVUVa9asYciQIT8fwGxmyJAhLF++vMF9li9fXq89wNChQ+vaO51O5s+fT6dOnRg6dCixsbGkp6czb968U9ZTVFRU71VVVuZ+4wQ9MCu5HIBt2zSJV0REpKlVVlYed632Fo8CTH5+Pg6Hg7i4uHrb4+LiyMnJaXCfnJyck7bPy8ujpKSEZ599liuuuIKvvvqK6667juuvv55vv/32pPUkJSVht9vrXs8sXw5OJ6xff3zjmBjica9CXVNT0MhvLCIiIqdrypQp9a7TSUlJXju2z58D43Q6Afjtb3/LvffeC0CPHj1YtmwZM2fOZNCgQSfcNysrq14XVGBgIJhMDffAXHMN1y7pwD0XDcDpLMDlcmEymbz7ZURERKTOI488woQJE+r+Lioq8lqI8SjAREdHY7FYyM3Nrbc9NzeX+Pj4BveJj48/afvo6GisVivdunWr16Zr1658//33J60nIiKi8WNodjvRvdMAcDgclJSUEB4e3rh9RURExGOBgYHuzoUm4NEQks1mo3fv3mRkZNRtczqdZGRk0L9//wb36d+/f732AIsWLaprb7PZ6Nu3L1t/Nel227ZttG3b1pPy4Icf4Kab4MknG3w7ODgYm80GwNGjmgcjIiLirzweQpowYQKjR4+mT58+9OvXj+nTp1NaWsqYMWMAGDVqFK1atWLKlCkAjBs3jkGDBjFt2jSGDRvG+++/z+rVq5k1a1bdMR944AFGjBjBRRddxCWXXMLChQv573//y+LFiz0rbt8+mDsX8vJg8uT67zmdOB6bTIAzlCqq2L//MG3atPH064uIiEgz4HGAGTFiBIcOHeLxxx8nJyeHHj16sHDhwrqJuvv27cNs/rljZ8CAAcyePZuJEyfy6KOPkpKSwrx58+jevXtdm+uuu46ZM2cyZcoU/vKXv9C5c2c++ugjBg4c6Flx1dXunw09iddsxvLydMpq2gBH2bPnCAMGePrtRUREpDnw+DkwzcEJ7yN/+234/e/hiivgiy+O37FDB5J31bCXffzrX3O4446bzlrNIiIi5zqfPQem2TvZk3gBYmPpiXt5gqqqw2epKBEREfG2cyvAxMQQdezXI0eOnJWSRERExPuMGWAamgMDEBtLC9zPfjl4UD0wIiIi/sqYAeYkPTBLGA7A0qUKMCIiIv7KWAFm7FgoLoZ//avh92NjicI9Z7mkRENIIiIi/spYAcZqhbAwCAlp+P1Roxj1kntBx9hY9cCIiIj4K2MFmFOJiiIhzb1kweHDCjAiIiL+ylgB5uOPYcwYePfdEzaJinLfh6S7kERERPyXsQLM2rXw1luwYkXD75eVUTnNHW4OHz5StxK2iIiI+BdjBZhT3YVktWJ760MAXC4nhYWFZ6kwERER8SZjBZjatZBOFGBsNtq2qCYA94rUhw9rGElERMQfGSvAnKoHBoiIDyGeKgCOHNFEXhEREX90zgUY4uKIPvZrfn5+k5ckIiIi3mesAONwuH+eLMDEx9cFmP37DzV5SSIiIuJ9xgowtXcVmU/yteLiyMT9MLvvv887C0WJiIiItxkrwEybBrm5MH78idvExxNBIACHDqkHRkRExB8ZK8CEh0NsrHs5gRO54w5uG9cVgPh4BRgRERF/ZKwA0xgxMbRO7QRAXp6GkERERPyRsQLMG2+4V6RevPikzWJjYwENIYmIiPgrYwWYr76CV1+FDRtO3KasjNJ/fQ3A9u3qgREREfFHxgowtbdRWywnbmO1UvPZSgAKC9UDIyIi4o+MFWAacxu1zUbvCPcaSC5XGaWlpWehMBEREfEmYwWYxvTAAF1aQdCx3zUPRkRExP8YM8CcrAcGMCUmEnPsd92JJCIi4n+MFWBqh5BO0QPDLwLMgQPqgREREfE3xgowjRxCIiGBbVwIwLJl6oERERHxNydZ9dAPvf02lJVBdPTJ2yUmEkooJUBOjnpgRERE/I2xemDi4qBdO/eSAidz663cNLotALGxuWehMBEREfEmYwWYxoqKom1qZwBycnJ8XIyIiIh4ylgB5tVX4cEHYdOmUzZNSEgAIDs7u6mrEhERES8zVoB59114/nnYvv3k7VwuKt5dBsCmjQfPQmEiIiLiTcYKMI15Ei+AyUTNN+6Qc/iIhpBERET8jbECTGNvowYuiK8CwOkspLy8vCmrEhERES87ZwNMajtz3XICmgcjIiLiX4wZYE41hASYWrUi4djvCjAiIiL+xVgBprFLCQAkJhJ/7Nc9exRgRERE/ImxAowHPTC0bs1PXALAt98qwIiIiPgTYy0l8OmnUFEBycmnbtu6NXZCOIqGkERERPyNsXpgOnaE7t0hLOzUbQcPZtRfugEQHa0AIyIi4k+MFWA8ER5Oco+ugJYTEBER8TfGGkL6+9+hqAj+8Ado1eqUzbWcgIiIiH8yXoDJyoJhwxoVYGo+XAHA9q37m7oyERER8SJjDSE1dimBY4JW7gSgrOIwFRUVTVWViIiIeJmxAowHT+IF6JEShJUAAA4cONBUVYmIiIiXndMBJjqlJe2pBiArK6upqhIREREvM2aAaeQQEq1bk3TsVwUYERER/2GsAOPJUgIASUm0Pvbr1q0KMCIiIv7CWAHmNHpg1vBbAP73P92JJCIi4i+MdRv14sVQXQ1JSadsCkDbtrTCzCYg/9C+pqxMREREvMhYPTC9ekF6OgQFNa59VBR/fO1qAEJC1QMjIiLiL4wVYDxlMtFxYD9Ak3hFRET8ibECzPPPwwsvQElJo3dJOjbcdOTIEcrKypqqMhEREfEiYwWYhx6CBx6A0tJG7xL4v6UEmIMB2LVLw0giIiL+wDgBxuVyv6DxdyEBgTs2U+NsC8C6dRpGEhER8QfGCTC1t1BD458DA5jaJdOBSgCOHtWdSCIiIv7AmAHGgx4YkpO5lN0AHD68x7s1iYiISJMwToCpfQoveNQDQ3Iy7Y/9umv7dq+WJCIiIk3jtALMjBkzSE5OJigoiPT0dFauXHnS9nPnzqVLly4EBQWRmprKggULTtj2j3/8IyaTienTp3tW1GkOIREVRVtbIACbN2zz7DNFRETEJzwOMHPmzGHChAlMnjyZtWvXkpaWxtChQ8nLy2uw/bJlyxg5ciR33HEH69atY/jw4QwfPpxNmzYd1/aTTz7hhx9+IDEx0fNvcrpDSCYTOyNvBGDTVt2FJCIi4g88DjAvvvgid955J2PGjKFbt27MnDmTkJAQ3njjjQbbv/TSS1xxxRU88MADdO3alb/97W/06tWLV155pV67AwcOcM899/Duu+8SEBDQqFqKiop+fjkcVH31Ffzvf2CzefSd+rQNB6CqJpfy8nKP9hUREZGGVVZW1r9WFxV57dgeBZiqqirWrFnDkCFDfj6A2cyQIUNYvnx5g/ssX768XnuAoUOH1mvvdDq57bbbeOCBBzjvvPMaXU9SUhJ2u939iorimaVL4ZJLPOuBAQZ/9iQREREA7Nmzx6N9RUREpGFTpkz5+Tptt9c9PNYbPLrS5+fn43A4iIuLq7c9Li6OnJycBvfJyck5ZfupU6ditVr5y1/+4kk5ZGVlUVhYWPd65JFHPNq/ljU+hvbt3VN5d+3adVrHEBERkfoeeeSRetdpby7b4/PVqNesWcNLL73E2rVrMZlMHu0bERFR13NCWRn861/u3pc//cnjOtq1a0dmZia7d+/2eF8RERE5XmBgIIGBgU1ybI96YKKjo7FYLOTm5tbbnpubS3x8fIP7xMfHn7T9d999R15eHm3atMFqtWK1Wtm7dy/33XcfycnJjS+uqAjuvtv98lRhIeUr3Osnfblwp+f7i4iIyFnlUYCx2Wz07t2bjIyMum1Op5OMjAz69+/f4D79+/ev1x5g0aJFde1vu+02NmzYQGZmZt0rMTGRBx54gC+//LLxxdU+B8bD+S8AhIZScdAOwE8/7vB8fxERETmrPB5CmjBhAqNHj6ZPnz7069eP6dOnU1paypgxYwAYNWoUrVq1YsqUKQCMGzeOQYMGMW3aNIYNG8b777/P6tWrmTVrFgBRUVFERUXV+4yAgADi4+Pp3Llz4wurDTAeDkMBYLVyZVQ+iw8D1Vs9319ERETOKo8DzIgRIzh06BCPP/44OTk59OjRg4ULF9ZN1N23bx/mX/SCDBgwgNmzZzNx4kQeffRRUlJSmDdvHt27d/fet4Az64EBru3m4KHvICd/Hy6Xy+P5OCIiInL2mFyu2iWc/UdRURF2u53CwsKfJ/Hu2QPt2kFwsHtCr4cqx44l5NVXcQIHDx4kISHBqzWLiIic6xq8fp8m46yFVJvDTrPnJLBzZ5KOnY7MTA0jiYiINGfGCTBnOIRESgq5XATA//6nACMiItKc+fw5MF4THw/z559+gOnYkTgi2Qvs3PmTV0sTERER7zJOD0xoKFx1FVxxxent36EDD06/BICqqu1eLExERES8zTgB5kyZzXQ9PxWArVs1hCQiItKcGSfAHD0Kb74J779/2oeofe7M7t27qaqq8lZlIiIi4mXGCTAHDsDtt4OHC0L+UtiSdVhNwTgcDrZu1ZICIiIizZVxAsyZ3oUEhLlKcLi6APDddxpGEhERaa4UYH7BfF5X0ikEIP+Q7kQSERFprowTYM7wQXYAdOrElabdAOzZknnmNYmIiEiTME6A8UIPDEFBdD22ptOP69d7oSgRERFpCgowv9Kps/tW6vU7duOsPaaIiIg0K8YJMLVDSGcYYJzthgJBVNZUsnPnrjOvS0RERLzOOAGmQwf44AN4+eUzOsx5F8UTQAoAK1Zs9EZlIiIi4mXGCTBRUXDjjXD11Wd0GNvokdzy+94A7Ny5wRuViYiIiJcZJ8B4i9lMaqp7HszGjeqBERERaY6Msxr1oUPwzTcQEXH6Czoec/755wOwcYMCjIiISHNknB6YzZthxAi4994zPpT1ky0AbNu+nbKysjM+noiIiHiXcQKMl26jBugYFwrEAC4yMzef8fFERETEu4wTYLx0GzVAq0EdST32+5YtGkYSERFpbowTYGp7YM5kKYFjTD3SGMIhANavWH7GxxMRERHvMl6A8UIPDHY7vWNiAFi9dOmZH09ERES8SgHmBDp3uwCA1T/toLq62ivHFBEREe8wToDx4hwYgMDuVwB2qp1VbNjwo1eOKSIiIt5hnACTmgpvvgmTJ3vlcN2GdyXS0hmApUtXeeWYIiIi4h3GCTCtW8Pvfw/XXOOVw1mGXMKfHhwMwMaNK71yTBEREfEO4wSYJtC3b18AVq1SD4yIiEhzYpwAk5MD8+fDcu/d9lwbYDZt2qQn8oqIiDQjxgkwy5a5V6J+8EGvHTLg85WYiMXhcLByZabXjisiIiJnxjgBxosPsqsVm94eKz0BWPiF5sGIiIg0F8YJMF6+jRrAlNqdu63rAdiZ+YXXjisiIiJnxjgBxssPsgPAauX68+IAWPLDD7hqQ5KIiIj4lPECjBeHkAD6Dh5MIJBXVMS2bdu8emwRERE5PcYJME0whAQQOGQIrWgNwMcfL/HqsUVEROT0GCfANMUQEsCFF1LNZQB8Nvcr7x5bRERETotxAkzfvjBjBowd693jRkRw64BEAPYeXK55MCIiIs2AyeWHV+SioiLsdjuFhYVEREQ0+eeVlpYSGRlJTU0Nu3fvJjk5uck/U0RExGi8ef02Tg9MEwoNDaVPnz4ALFmieTAiIiK+ZpwAc/AgfPMNbNrUJIe/oFN3AD6evbBJji8iIiKNZ5wAM38+XHopTJzYJIff93kUAIu+ydA8GBERER8zToBpoufA1Lrt4g6YCKasKo8NGzY0yWeIiIhI4xgvwHj7Nupjrh2bwjDKAfhi/vwm+QwRERFpHOMEmCZ6kF0t828u5MqgIAC+mDu3ST5DREREGsc4AaaJh5AICODKSy8FYOn6DRQWFjbN54iIiMgpGS/ANFEPDECLq2/FRjIOl5P//vfrJvscEREROTkFGA9E3jCEYC4G4L23P22yzxEREZGTM06AGTgQnnsORo5sus+IieHpv/QCYN2PX+OsDU0iIiJyVmkpAQ9VVlYSGxtLUVER33//PRdeeOFZ/XwRERF/paUEfCgwMJDf/va3AHwwZ46PqxERETk3GSfAHDwIK1fC7t1N/lFp5W0AeOP19zWMJCIi4gPGCTD//jekp8NTTzX5R7V1dADslFQc4rvvljb554mIiEh9xgkwTf0cmF+4ZmI/UmkPwNx3/9PknyciIiL1GS/ANOFt1LUCe53Hs22PAPDRBx/gcDia/DNFRETkZ8YJME28lMCvDRkzhkggp7CQ77777qx8poiIiLgZJ8CcxSEkANuttzIYOwAvPP3Ps/KZIiIi4ma8AHOWemDo0IFK+6MALPzfPAoKCs7O54qIiMjpBZgZM2aQnJxMUFAQ6enprFy58qTt586dS5cuXQgKCiI1NZUFCxbUvVddXc1DDz1EamoqoaGhJCYmMmrUKA4ePOhZUWc7wACP/rETYdYOOJwVvPvuu2ftc0VERM51Hl/t58yZw4QJE5g8eTJr164lLS2NoUOHkpeX12D7ZcuWMXLkSO644w7WrVvH8OHDGT58OJs2bQKgrKyMtWvXMmnSJNauXcvHH3/M1q1bufbaaz0r7JJLYPJkuPJKT7/Saev/7HCenvYXAGbNmoUfPtRYRETEL3m8lEB6ejp9+/bllVdeAcDpdJKUlMQ999zDww8/fFz7ESNGUFpayueff1637YILLqBHjx7MnDmzwc9YtWoV/fr1Y+/evbRp0+a492sfRZyVlVXvUcSBgYEEBgZ68nXO2JEjR0hMTKSyspIVK1bQr1+/s/r5IiIizVVlZSWVlZV1fxcVFZGUlHT2lxKoqqpizZo1DBky5OcDmM0MGTKE5cuXN7jP8uXL67UHGDp06AnbAxQWFmIymYiMjDxpPUlJSdjt9rrXlClTGv9lvKRly5Zcd+lQAJ598pWz/vkiIiLN1ZQpU+pdp5OSkrx2bI8CTH5+Pg6Hg7i4uHrb4+LiyMnJaXCfnJwcj9pXVFTw0EMPMXLkyFOms6ysLAoLCyksLKRo2zYe/e1vITfXg2/kHcXL3L0u/134oSbzioiIHPPII4/UXacLCwvJysry2rGb1V1I1dXV3HTTTbhcLl577bVTto+IiKh7hc+aha1XL3j++bNQaX2P/DEVGynUOMuZ2Yi6RUREzgWBgYH1rtVnOmz0Sx4FmOjoaCwWC7m/6uXIzc0lPj6+wX3i4+Mb1b42vOzdu5dFixZ5/iVrp/KcpefA/NKAiYN5PWg/ANOfe47y8vKzXoOIiMi5xKMAY7PZ6N27NxkZGXXbnE4nGRkZ9O/fv8F9+vfvX689wKJFi+q1rw0v27dv5+uvvyYqKsqTstx8GGBMYaGMHD+ONkBuQQFvv/XWWa9BRETkXOLxENKECRN4/fXXefvtt9myZQt/+tOfKC0tZcyYMQCMGjWKRx55pK79uHHjWLhwIdOmTeOnn37iiSeeYPXq1dx9992AO7zccMMNrF69mnfffReHw0FOTg45OTlUVVU1vjAfBhiAgPvu476AAACeeOxpampqfFKHiIjIucDjADNixAheeOEFHn/8cXr06EFmZiYLFy6sm6i7b98+srOz69oPGDCA2bNnM2vWLNLS0vjwww+ZN28e3bt3B+DAgQN89tln7N+/nx49epCQkFD3WrZsWeMLO8trIR0nOpqowS8D0eQePcD773/omzpERETOAR4/B6Y5qH0OTL37yMeNg3/8Ax55BJ55xid1Ve4+SEKHFzjq+jvdunVn48b1mH0VqERERJqZBq/fp8k4V1cfDyEBBLZLZMu+h4iIiGDz5k288847PqtFRETEyIwTYAYPhvvvh9/8xqdlxLWO47HHHgPgscceo6yszKf1iIiIGJFxhpCakYqKCjomJXMgP5eJj/6Nvz090dcliYiI+JyGkJq5wOoabIfdizxOnTrluOfgiIiIyJkxToDJz4c9e6AZPMrfFB7GP/7QlTC6UO0o44lHH/V1SSIiIoZinAAzcSK0awcvveTrSgAY9srVfN66FIBZb77JqlWrfFyRiIiIcRgnwDSDu5B+yWQLYNDM1/g/wOlyMWbkyHpLiouIiMjpU4BpSlddxfTLhxJOMD/u3MnER5/0dUUiIiKGoADTlEwmwma+RahpKgAv/v051q1b5+OiRERE/J8CTBMLbhfPR0+lktZmME6XgzFjxni2xpOIiIgcx3gBphk+un/Aoxfz1arZREVFsX79+nqLXYqIiIjnmt/V/nQ5ne6fzawHplZsbCyvv/46AC+++CLvvTfHxxWJiIj4L+MEmMGD4c9/hl69fF3JCQ0fchkp1hsBGH3bGH788UcfVyQiIuKfjBNgbr0VZsyAyy/3dSUnZAoP4+FhN2PmYqod5Vx35ZUUFhb6uiwRERG/Y5wA4ydu/3g46y5NpA2wPSuLW2+8kZqaGl+XJSIi4leME2AKCyEvD5r76s9mM+fPm8lH7doRCHy+aBG33TwGZ+0cHhERETkl4wSYP/8Z4uJg5kxfV3Jq4eH0+fJLXgtOBCy8/9E7/Pmue/DDhcFFRER8wjgBppk+B+aEUlIY+PZ8Ivg7AP/8f6/y5JN6Uq+IiEhjKMD4UMqNPfj+Pxfx+H1TAHjyySd5/vnnfVyViIhI82f1dQFe04wfZHcyqbemkXprGlZ7NY8//jgPPvggB7JyeHH685j97LuIiIicLca5QvphD8wvTZw4kXuvuA2Al15+kZuuG6klB0RERE5AAaaZMJlM/N9VfyKcfwAWPvrsA64eMoTi4mJflyYiItLsGCfA1N6G7MfDLn3u6c+6zwbzvr0VocCi775jYO/ebN++3deliYiINCv+e7X/tUsugdGjoUsXX1dyRjpc040RG77jm/btiQU2bN9O9249mD37PV+XJiIi0myYXH748JGioiLsdjuFhYVERET4upymUVDAsqF3MnBlDi6+B+Cuu+5i+vTpBAcH+7g4ERERz3nz+m2cHhijiYyk/7I5vDr8YbpF/QGTycSsWbO44IILWL9+va+rExER8SnjBJjycigtBQOtK2SymPnjJ8P4Mf91vvzyS2JiYtiwYQM9e/Ri3J/vpbS01NclioiI+IRxAsyIERAWBm++6etKmsRll13G+sxMzgvtgwsn/3htOue1b8+CBQt8XZqIiMhZZ5wA4+e3UTdGQkIC/3hoOonmWUQRwt68PIYNG8YNV13Fjh07fF2eiIjIWaMA409MJi6ddCG7Dv0fu8f/f9xvMmEBPvriCzqldOGWG0aRlZXl6ypFRESanPECjB8/B6axAluGEv73F3l+40ZW9+xFBANw4WD2R/8hJSWFe++9l7y8PF+XKSIi0mSMc7U/F3pgfu288+ixZjWfTX2RfpH/on/6QCorK5k+fTrt2iZz+//dybZt23xdpYiIiNcpwPg7k4lBD6az4ugdLF2+hC+//JI+vXtTVlHOm+/9i86dOzNs4EC++uor/PCRPyIiIg0yToAxwFICZ8pkMnH55Zez9L0POT/4FeBqTMCCpUsZOnQoXVsl8fcXppOfn+/rUkVERM6Ica72F10EN9wAbdv6uhKfs6Uks770z/w4+ym2Xn0NfzGbCQe2Zh9gwgP3Ehcbz29/+1vmzp1LeXm5r8sVERHxmJYSOBfk5FD0j3/Qc2oFu5zfAmvr3ooIC+fS9Cv480N3MWjQRdhsNt/VKSIihubN67cCzDnEUVnDkpfXE/abQD75dDbvvPNOvduuI2w2rrjwQq4ZPZorr76aqKgoH1YrIiJGowDT0Ak4VyfxngGn08mfuv2Nf23dSxAfUMbPSxOYMJMSn8qosTdz5VWXk5aWhsVi8WG1IiLi7xRgGjoBl1wCixfD+++7lxWQRivae5TyBRnsnvc6/12yhLkVsJ2Kem3sdju9OqRx8WWXcfVNVynQiIiIxxRgGjoBF18M334Lc+bATTf5tD6/VlPDvk/X8fwzO9hcuI7gLpv57rvvKCoqqtfMbrPxm+7d6XfxJfS77DL69utLy5YtfVS0iIj4AwWYhk7AoEGwZAl88AHceKNvCzSYmpoa1n2xkBuu/x/7arYRyheU4jyuXceOHenbty/9+vWjb9++9OjRg9DQUB9ULCIizZE3A4zVSzX5nubANBmr1Urfa65mb/XVHP7pEBHrv2H9xx/x3dKl/OtAFJspBXayY8cOduzYwXvvvXdsTxNtIiLp2SGZ7j17ct5vLuL8vn3o1KkTAQEBvvxKIiLi54zTA/Ob38D338OHH8LvfufbAs8hVUUVbPp0JxWhUBJ2gJUrV/L9R5/zZeY+ILvBfSxmK93O60rXrl1JSU4mrkUbeg3sQafOnYiOjsakECoiYkgaQmroBAwcCEuXwkcfwfXX+7bAc1x1UTlr3tvG2sU76FrxDps2bWLT/v18XJFEPjlA8Qn3tQcEkBIdTUx0CnGtu9Lr4vM5r08X2rZtS1JSkp5TIyLixxRgGjoBAwbA8uXwyScwfLhP65OGHf7pEBu/yMLa0Uyh+QBbt27lh9lfMneNAzM/4eIAJ/uX0WQyEWSOIiokjoFdQmjbrh1tu3QhsWt32nTuSGJiIjExMZjP4eUkRESaMwWYhk7A/ffD5s0waRL07+/bAsUjFUfLyVmVRVxBJrtWrWL7hg1MWxFHZlEN7drtpdKaz759+6ioqDjlsSwWC/Hx8cRFtKDyiJ2EmBgG9W1BYocOJHTqRExyMrFxccTExBAcHHwWvp2IiNRSgNGTeM85LpeLDXO/4/XpP3L00G56BXzO3txc9hYV8U1NO4opxkQurpP24dQXFhZGbGwsLYPCKMkLIyYykoFpwcS2akV069aERLcmtn0ScYlxtGzZksjISD37RkTkDCjAKMDIL5RkF5O9MZ8WXaKospZw8OBBVr21iDffO0xN1X56WRZysLSU7JoathBNJYVAtcefYzKZsIeGUV0eQZgtnJ4JVbS022nZogUVpgQsQTF07JtM8nmtsNvtRNrt2CMjiYyMxG63ExgY6P0vLyLiRxRgFGDkdJSXk73pMAe3FtLi/EgcwWXk5eWR+f63zP7gKI6qbPoFfk1eSQmHKypY4WhHMdUE2fKpqCo99fFPIdBkIjIgAEwtcNCC1m0jads1jvDwcIJr4OAWJ/aIMNLPDyC8RQvCo6Iwh8QQFh1DbHIskS0jCQsLIzQ0lKCgIN2tJSJ+RwHmRHchrV3rnsQ7dKhvCxRDcFTWcHRPIcExYdjCzRw5coStX6zik7d2U1FyiPQWyzly5AiHCwuZn92BPRVBdOxUTEhMDQUFBRzam01OURVQdMrP8pwZM0HE2moItVoJDQig3BVFhbMFrZMjaZ8WT2hoKMEOF3vWVREeFkLvri7CIsIJDg+n2hSJy2YnrkMsMW2iCAkJISQkhODg4LqfelaPiHibAkxDJ6B3b3eAWbAArrzStwWKANWFZRzYeISyvCJaBe6mIDubwtxc5q+wsDnLSeeLw4npYqO4uJis7zcxb2E1LlcxA4K+paiqiuKaGjbVtKIAC2FBR6mhvFETmb3FCgSbzZhMIVS5IrCH20hoF0lwcDBBwK5NTgIsNnomFBIebCUoKIii6hYcqWhJq5QoOvVvRVBQEEEmEwfWFhIUEkxK+wDC7CEEhobisoZhCgzHHmcnvGUYgYGB9V42m029TCIGoyfxNkRP4pVmJsAeQvLAkGN/dcN+7LfzG2p8P7zU0HaXi8qiSrBYCAwLoKamhsM/7WXJ3J0UHy2kd6tdlBw9SmlREd9sDmXzQRsd0oNo1SOE0tJS8jds55PPK6lyVHJJ2BLKqqspr6lhc3UiOc4IWoYVENzCSXl5OWXFJZRVVsKxidA1QLHTCZQAJVQUQm7m8SXuKGig7o3Ax40+VSdkM5kINJkwEUiFK5yQoAASkiPcAcdsZv+WSpwE0jmyAHuQC5vVSqkjgoOlMbSMC+G8ixLr2u74JgeXKZjzkipoGWnFFhREhTOU3MJwIhMjSRnQioCAAAKsVo78mI/FFkx8go2wiEBsoaGYbcGYAoIJDg8mOMzdQ1X7stlsmuAtcpYZpwemZ0/IzISFCzWEJHIaHCXl7F6RS+GhYtpF5VJeWEB5YSFrN9ewZY+J+PMjSL6gBRUVFZTs2ssHrx2gvKqay1qtoaaqjIrKStblx7CxIJpWbWpo3y+UiooKKnLz+Pq7GpyuGrpZ1+N0VlHpcpHriuAooQRayggIrqGyspLqas8nVzcXJiAACDCZMJuslDvDsVrMxLYKISAgAKvVypHdZVQ6g0kMKqZlYBVWs5kaVxB7yhIIDg4g9TcxdaFo7/fZlFQE0yGqhER7FVarlRpXED8diiM4PIg+17bBarVitVrZv2QfxSUBtI2tJjHGhdVqxWkKZFeencDQYFKvTK5re3RzHpXlFqKjLERHW7HabJgCAikuDyIwNIjo9i2xWCxYrVZc1U4CbDYCAgMIsAXUbf/lT/WSiSc0hNTQCejRA9avhy+/hMsv92l9InJqLqeLmooanE4IDHPPt3GWl7Pxy90UHS2lbfQRqCqhsqSEHXsd/LjTTETbEFIuiaG6upqqo0f5/KUtFJU5GdB6B0GmYiorK9l2KJRVObFEJZjoMSyGqqoqqoqKWPjOAYqrzfQJ20iwqYgqh4Osigg2VyQRHlJBhx6BVFdXU11ZyU/ri6nATDz7MVNBNVBGACXYMFEFphr88D+dTcJsNmOxWHDVgMNlI9DkINRcjQWwmswcckRhxkqrdkHYbO7QU5xVTEFpCJEB5SQGF2ExmbCYzWwtTgKThfMGRBEUZMFisXD4x0NkH7IRE1xBB/tRLBYLFrOZDflJOLHSY2giYRHuHrDD67PJ2msiOqyaLnFFdUFrY3YcDlcAqVe1IaJlEBaLhaObDrJvu4MWYQ66tC7HYrVitljYlh1JjdNKp0vbYo8Jxmw2U7r7MAd2VBARaqJ9ksNdQ0AAB/KDcLoCSOwZT3hUCBaLhar8EgqyqwgONhMX525ntloprbCC2UZ4fBiBIcd67GqcOB1gDbASGOxuV3s+zWZz3av2byOERQWYhk5AWhps2ABffQWXXebbAkXEOFwucDioKa+mrLAal8mMvVUYDoeDqspKtn+1lYKjVcTbCwkyl1FVVkb+YQfrd9gIsAfQ/erW1NTUUF1dzQ+vrST3iIlucblEBRdSU11NToGZ73bGYgu3MPC2n9t+/89MDhba6N5iL/GBh6hxOMgvC+DbvA4EBMJvfteSmpoa92rx83aSXRZCh6C9xFnzqHE6Ka4JYE1ZCmazg/MvCKhru3d9Pkdrgoky5RHOUWqASpeZHKKBaux2Bw6Hg5qaGqoqqnHiBA+eryRNzUIADsyA2WSixhVIDQEEB5kIDg3AbDbjqqihoDgAs8lEvPUIZpMJM1DgiKTMGULLliYi49zhzFlczr79VsxAp9Bsd1uTidzKlhTWhBOXYCWufbi7bWEpP/3owoSZni331bU9WB7FoYpIEtvaSOrewv009KIS1i6rxGwy0zduPwEWd9jdUxzO13uWKcDUOwHnnw8bN8LXX8Pgwb4tUETEDzlrnJitPy/FUZZbTHlRNYHmKoICqqiprKSqvIqdu51U1bhIHhiP0+kOO/uW7ODg3nJahpaREFlKTVUVVZXVfL8xlGqHi763dsRkcrfd9t8f2fFTJfERxaTEHMZRU4OjpoYvNiVS43Dym7s6YbW6g9S2/25m809OEkILOC/6IA6He/sXuzpQWWNi4P+1IijEhcPhYNeibWzaEUBs4BG62/fgcDpxOJ18k9eZcoeFfpfZCY4w4XQ62f/DXrZkBdHSWkSXkF04nU4cLherS1IodwbQNc1CcIQZp9NJ/k+57M4PJdRUQlvrXpwuFw6Xi12O1lS4AkiIryEw1IzD4aA0r5hDZYEEUEkkh3ECDqCYUByYCLQ5wOR0f16NA6fLxbkYEH0WYGbMmMHzzz9PTk4OaWlpvPzyy/Tr1++E7efOncukSZPYs2cPKSkpTJ06lauuuqrufZfLxeTJk3n99dcpKCjgwgsv5LXXXiMlJaXB4zUYYG6/HXbtghdfhF69PP1KIiIiTcZZ46wLiLUhsaa0kpL8CnA6CA+pxllTg6OqirxDTspKndjbRRJkt+J0Oik9eIS9649goobkuDKcDgeOmhp2ZZk5WggJPWNpkRSC0+mkZHcuG77Jw2xykNauAEd1NU6Hgy1ZIeQdtdD2gjgSutrdx92dww+f54HLwYUpOTidTpwOBz9mhbH/aBDt+kXTpldLnE4n5Xtz+ObjfJwOJ0M67XG3dTrZlB3JniNhtOsZQYcLonE6nVQdzOOruUdxOp0MTt6GyeXA6XSyIS+CD7Yt8E2AmTNnDqNGjWLmzJmkp6czffp05s6dy9atW4mNjT2u/bJly7jooouYMmUKV199NbNnz2bq1KmsXbuW7t27AzB16lSmTJnC22+/Tbt27Zg0aRIbN25k8+bNBAUFHXdMPchORETE//h0Dkx6ejp9+/bllVdeAcDpdJKUlMQ999zDww8/fFz7ESNGUFpayueff1637YILLqBHjx7MnDkTl8tFYmIi9913H/fffz/g7lqKi4vjrbfe4uabbz7umAowIiIi/seb12/zqZv8rKqqijVr1jBkyJCfD2A2M2TIEJYvX97gPsuXL6/XHmDo0KF17Xfv3k1OTk69Nna7nfT09BMes1ZRUVG9V2VlpSdfR0RERJpQZWXlcddqb/EowOTn5+NwOIiLi6u3PS4ujpycnAb3ycnJOWn72p+eHLNWUlISdrsdu93Oj3Y7NS1bwg8/ePKVREREpIlMmTKl7jptt9tJSkry2rE9CjDNTVZWFoWFhRQWFtKvQwdCy8qgpsbXZYmIiAjwyCOP1F2nCwsLycrK8tqxPVpKIDo6GovFQm5ubr3tubm5xMfHN7hPfHz8SdvX/szNzSUhIaFemx49epy0noiIiJ/H0Gof8GOAB/2IiIgYQe3aZk3Box4Ym81G7969ycjIqNvmdDrJyMigf//+De7Tv3//eu0BFi1aVNe+Xbt2xMfH12tTVFTEihUrTnjMBmktJBERkXOGx4s5TpgwgdGjR9OnTx/69evH9OnTKS0tZcyYMQCMGjWKVq1aMWXKFADGjRvHoEGDmDZtGsOGDeP9999n9erVzJo1CwCTycT48eN56qmnSElJqbuNOjExkeHDhze+MAUYERGRc4bHAWbEiBEcOnSIxx9/nJycHHr06MHChQvrJuHu27fP/RjhYwYMGMDs2bOZOHEijz76KCkpKcybN6/uGTAADz74IKWlpdx1110UFBQwcOBAFi5c2OAzYE5IAUZEROScYZylBNq1gz17YMUKOMlTgUVERMQ3vPkcGI97YJqt88+H6GgIDfV1JSIiItLEjBNgPv3U1xWIiIjIWeLXz4ERERGRc5MCjIiIiPgd4wSYfv2gbVvYvNnXlYiIiEgTM84cmKwsyMmB6mpfVyIiIiJNzDg9MHoOjIiIyDlDAUZERET8jgKMiIiI+B0FGBEREfE7xgswZuN8JREREWmYce5C6twZjh6FwEBfVyIiIiJNzDgBZulSX1cgIiIiZ4nGW0RERMTvKMCIiIiI3zFOgOnbF7p1g337fF2JiIiINDHjzIH56ScoKYGaGl9XIiIiIk3MOD0wtbdRi4iIiOEZL8DoQXYiIiKGpwAjIiIifkcBRkRERPyOXwaYysrKej8BBZjTUFlZyRNPPFH/PMpp0bn0Hp1L79B59B6dS+9p8Pp9mkwul//Nft2/fz9JSUlkZWXRunVr98bOnaGyEpYtg8RE3xboJ4qKirDb7RQWFhIREeHrcvyazqX36Fx6h86j9+hcek+D1+/TZJzbqLdu9XUFIiIicpb45RCSiIiInNv8sgemdtSruLiYoqIiH1fjv2rPnc7hmdO59B6dS+/QefQenUvvKS4uBn6+jp8Jv5wDs2vXLjp06ODrMkREROQ07Ny5k/bt25/RMfwywDidTg4ePEh4eDgm3XUkIiLiF1wuF8XFxSQmJmI2n9ksFr8MMCIiInJu0yReERER8TsKMCIiIuJ3FGBERETE7yjAiIiIiN/xuwAzY8YMkpOTCQoKIj09nZUrV/q6JL8zZcoU+vbtS3h4OLGxsQwfPpytepLxGXv22WcxmUyMHz/e16X4pQMHDnDrrbcSFRVFcHAwqamprF692tdl+R2Hw8GkSZNo164dwcHBdOjQgb/97W9eee6G0S1ZsoRrrrmGxMRETCYT8+bNq/e+y+Xi8ccfJyEhgeDgYIYMGcL27dt9U2wzdrLzWF1dzUMPPURqaiqhoaEkJiYyatQoDh486PHn+FWAmTNnDhMmTGDy5MmsXbuWtLQ0hg4dSl5enq9L8yvffvstY8eO5YcffmDRokVUV1dz+eWXU1pa6uvS/NaqVav45z//yfnnn+/rUvzS0aNHufDCCwkICOCLL75g8+bNTJs2jRYtWvi6NL8zdepUXnvtNV555RW2bNnC1KlTee6553j55Zd9XVqzV1paSlpaGjNmzGjw/eeee45//OMfzJw5kxUrVhAaGsrQoUOpqKg4y5U2byc7j2VlZaxdu5ZJkyaxdu1aPv74Y7Zu3cq1117r+Qe5/Ei/fv1cY8eOrfvb4XC4EhMTXVOmTPFhVf4vLy/PBbi+/fZbX5fil4qLi10pKSmuRYsWuQYNGuQaN26cr0vyOw899JBr4MCBvi7DEIYNG+a6/fbb6227/vrrXbfccouPKvJPgOuTTz6p+9vpdLri4+Ndzz//fN22goICV2BgoOu9997zQYX+4dfnsSErV650Aa69e/d6dGy/6YGpqqpizZo1DBkypG6b2WxmyJAhLF++3IeV+b/CwkIAWrZs6eNK/NPYsWMZNmxYvX83xTOfffYZffr04cYbbyQ2NpaePXvy+uuv+7osvzRgwAAyMjLYtm0bAOvXr+f777/nyiuv9HFl/m337t3k5OTU+9+53W4nPT1d16AzVFhYiMlkIjIy0qP9/GYtpPz8fBwOB3FxcfW2x8XF8dNPP/moKv/ndDoZP348F154Id27d/d1OX7n/fffZ+3ataxatcrXpfi1Xbt28dprrzFhwgQeffRRVq1axV/+8hdsNhujR4/2dXl+5eGHH6aoqIguXbpgsVhwOBw8/fTT3HLLLb4uza/l5OQANHgNqn1PPFdRUcFDDz3EyJEjiYiI8Ghfvwkw0jTGjh3Lpk2b+P77731dit/Jyspi3LhxLFq0iKCgIF+X49ecTid9+vThmWeeAaBnz55s2rSJmTNnKsB46IMPPuDdd99l9uzZnHfeeWRmZjJ+/HgSExN1LqVZqa6u5qabbsLlcvHaa695vL/fDCFFR0djsVjIzc2ttz03N5f4+HgfVeXf7r77bj7//HO++eYbWrdu7ety/M6aNWvIy8ujV69eWK1WrFYr3377Lf/4xz+wWq04HA5fl+g3EhIS6NatW71tXbt2Zd++fT6qyH898MADPPzww9x8882kpqZy2223ce+99zJlyhRfl+bXaq8zugZ5R2142bt3L4sWLfK49wX8KMDYbDZ69+5NRkZG3Tan00lGRgb9+/f3YWX+x+Vycffdd/PJJ5/wv//9j3bt2vm6JL80ePBgNm7cSGZmZt2rT58+3HLLLWRmZmKxWHxdot+48MILj7uVf9u2bbRt29ZHFfmvsrKy4xbJs1gsOJ1OH1VkDO3atSM+Pr7eNaioqIgVK1boGuSh2vCyfft2vv76a6Kiok7rOH41hDRhwgRGjx5Nnz596NevH9OnT6e0tJQxY8b4ujS/MnbsWGbPns2nn35KeHh43fit3W4nODjYx9X5j/Dw8OPmDYWGhhIVFaX5RB669957GTBgAM888ww33XQTK1euZNasWcyaNcvXpfmda665hqeffpo2bdpw3nnnsW7dOl588UVuv/12X5fW7JWUlLBjx466v3fv3k1mZiYtW7akTZs2jB8/nqeeeoqUlBTatWvHpEmTSExMZPjw4b4ruhk62XlMSEjghhtuYO3atXz++ec4HI66a1DLli2x2WyN/6DTvjfKR15++WVXmzZtXDabzdWvXz/XDz/84OuS/A7Q4OvNN9/0dWl+T7dRn77//ve/ru7du7sCAwNdXbp0cc2aNcvXJfmloqIi17hx41xt2rRxBQUFudq3b+967LHHXJWVlb4urdn75ptvGvxv4+jRo10ul/tW6kmTJrni4uJcgYGBrsGDB7u2bt3q26KboZOdx927d5/wGvTNN9949Dkml0uPZxQRERH/4jdzYERERERqKcCIiIiI31GAEREREb+jACMiIiJ+RwFGRERE/I4CjIiIiPgdBRgRERHxOwowIiIi4ncUYERERMTvKMCIiIiI31GAEREREb/z/wOZX3Cn+Xm3LAAAAABJRU5ErkJggg==",
      "text/plain": [
       "<Figure size 640x480 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "styles=[\"r--\",\"b:\",\"k\"]\n",
    "labels = [\"strict LO\",\"subtracted\",\"tuned\"]\n",
    "\n",
    "plt.rcParams['ytick.right'] = True \n",
    "plt.rcParams['ytick.labelright'] = False \n",
    "plt.rcParams['ytick.left'] = plt.rcParams['ytick.labelleft'] = True\n",
    "plt.rcParams[\"xtick.top\"] =  True\n",
    "plt.rcParams['xtick.direction']='in'\n",
    "plt.rcParams['ytick.direction']='in'\n",
    "\n",
    "plt.plot(k,rateplt[1],styles[0],label=labels[0])\n",
    "plt.plot(k,rateplt[2],styles[1],label=labels[1])\n",
    "plt.plot(k,rateplt[3],styles[2],label=labels[2])\n",
    "plt.legend()\n",
    "plt.xlim(0.,12.)\n",
    "plt.ylim(-0.01,0.1)\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": ".venv (3.12.11.final.0)",
   "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.12.11"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}