WattSpec

{
 "metadata": {
  "name": "",
  "signature": "sha256:0b269c2597f8b6955f508e0f2faa72bcf9ea11b4801e967f18cb1bb8f9695295"
 },
 "nbformat": 3,
 "nbformat_minor": 0,
 "worksheets": [
  {
   "cells": [
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "from IPython.display import Math"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 99
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "import seaborn as sb\n",
      "sb.set_style(\"whitegrid\")"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 61
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "import matplotlib.pyplot as plt"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 6
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "def watt(E):\n",
      "    return 0.4865 * sinh(sqrt(2.0 * E)) * exp(-E)"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 116
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "f, (ax1, ax2) = plt.subplots(1, 2, sharey=True, figsize=(14,6))\n",
      "dE = np.linspace(eps, 10, 10000)\n",
      "ax1.plot(dE, watt(dE), label=\"W(E)\", linewidth=2)\n",
      "ax1.set_xlabel(\"Energy [MeV]\")\n",
      "ax1.set_title(\"Watt spectrum, linear scale\")\n",
      "\n",
      "\n",
      "ax2.semilogx(dE, watt(dE), 'r', linewidth=2, label=\"W(L)\")\n",
      "ax2.set_xlabel(\"Energy (log scale) [MeV]\")\n",
      "ax2.yaxis.tick_right()\n",
      "ax2.set_title(\"Watt spectrum, semilogx scale\")\n",
      "\n",
      "f.show()"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "metadata": {},
       "output_type": "display_data",
       "png": "iVBORw0KGgoAAAANSUhEUgAAA08AAAGOCAYAAABPFTgiAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJzs3XecVNXdx/HPbGMXWHqRDiIe6eqCIqgIakzsLSbGFkse\nS0yMTzSWaKImj91oihp7jDGaGMUSA1EpFsTCohSRQ5fe6y7bd54/zgysyO7O7s7MuTP7fb9evOZO\nuXN/c1n28p3TQuFwGBEREREREalbhu8CREREREREUoHCk4iIiIiISAwUnkRERERERGKg8CQiIiIi\nIhIDhScREREREZEYKDyJiIiIiIjEQOFJYmaMmWiMuabG/QONMdXGmDtrPNbFGFNmjMmv431OMsbc\nHtlua4yZktjK62aMecsY08HDcX9ojHkjsv2EMWZ8smuIF2NMX2PMTt91iEj60jUodRhjbjfGnB/Z\nrk6Vz6drmcRC4Uka4j/AMTXunwK8AZxa47HxwAfW2rp++YwEor9I20fu+3QcEPJZgLX2R9Zarxdw\nEZGA0zUoRVhrf22t/ZvvOkQSIct3AZJSJgG31bh/MnAz8KIxpp+1dhlwLPCmMaYl8GdgAO4itRP4\nAe5CdTmQaYzZDowG8owxs4AR1trq6JsbY44EHgAygTBwl7X2FWPMXyIvORDoArwF/NRaW2mMGQg8\nBHSM7PcHa+0zkfe7BPhfoArYBFwE/CbyXlOMMScBHwAfAcMin+0h4CxrbWHkPZYDZwJbgCnAZOAI\nIBu4LvLZDgJmAudaa2NahdoYMw34I1AYec83gcMj5+6X1tp/Rl73y8jxM4DlwFXW2rXGmFHAPUAL\noBvwtrX2MmNMX+B9YD7QFzjaWru+xnHPBH4JVEfOy/XW2veNMfvh/v5M5Lk/W2v/WNtx9vF59lln\nLOdCRKQWugbF8RoU+T3/10itAG9aa38Vee5S4Erc7/DNwNXWWhv57CXACGA/4J/ARlyQ3Q+4zFo7\nNfK6udbaB/Y65q3A94FKYGHkfdcbYw4Ano78/azFhcm/ASsixxgOrIt83snW2t/u9b66lknSqOVJ\nYmatXQRsMcYMM8a0x/0y+gj3beBpkZeNx/3H/zvAFmvtEdZaA3yK+yX5Me4X2YvW2luAi4ESa+2h\nNS9aEbcBv7PWjgAuAcbVeG44cDwwKPLncmNMFvAv4MbIPscA1xtjDjfGDAfuBk6w1g4HXgduttZe\nHHm/cdbaVbgL5Fxr7SBr7auR+zUvPjW3+wKvWWuH4H6h/x53URgMHAWMiunE7nnf6Hv3AyZZaw8H\nbgDuBTDGXAgMAQ6z1h4CTASejOzzU+BWa+2oyPFPNcYcEnmuB3CHtdbUDE4R9wJXWmtHArcCYyOP\nPwIssNYOxF2Y/8cY07+e4xBDnSIijaJrEBDfa9CPgCXW2oLI6wcYY9oYY8YCFwJHWWsPBe4DXtnr\ns4/CBahrgZ3W2jGR499Yo86vBTdjzMXAt3EhdTgwD/hL5OnngOettUNx15kjgLC1dirwGPAEcAvu\n7+prwSlC1zJJGrU8SUNNxF1ANgJvWWvDxph/Az82xkwAsNZawBpjlhljfgIcgLuIfBh5jxB7uijU\n1VXhn8DDxphTgHdw3yqB+4X8rLW2GMAY81fgdNy3cPsDTxtjou/RAjgYaIkLJKsjNf6+juO+X99J\niKiw1v47sr0EmG6tLYrUtAb3DVpjVFhr/xPZ/ow93UtOxnUvmRn5fJlAXuS5i4CTjDE3AQNxn7c1\nsBX3Dd+MWo71IvCqMeZN4G3cRRLct7fXAVhrdwBDI5+rruNE1VWniEhT6Bq0R1OvQROB/xhjekc+\n343W2h2RFrADgA9rfI72kcAaBt6w1lYB640xxbgWQYCl7Lle7S2EC7RPW2tLIo/9AfilMaYL7ppx\nJIC1doExZnKNfX+N+7u7Ehd09kXXMkkatTxJQ00EjgZOAqK/tKfiLg7HRR8zxlyJ+4amCHgeeIEG\n/rxZax/H/aJ7GzgBmGOMaRN5uqrGSzNxASET2GatPST6BxgDPBt5fjdjTAtjzIBaDl1UYzvM1y+u\nOTW2y/far5L4qPm+NY+fAdxd47ONwP1dgOvq8W3gS+B2YFWN/cr28Y0qAJFvXsfgunj8EJhhjAnt\n/VmMMf0i576u40TVVaeISFPoGrRHk65B1tqZuJ4Oj+NasT4xxhyBO0/P1fgMhwKjrLXRYLH3cSti\nPGTN0ErkOFlAaY37UTWvWe1wXQKrcF0l9/VZdC2TpFF4koaaChyCaxL/L4C1dhcwC7ga110C4FvA\nXyJ9vRfiBvRmRp6rwPXPhj0XnG8wxkwHDrHWPovrx90O901aCDjHGJNjjMnFdS94A7BAqTHmvMj+\nvYDZkXqnAsdF+j+D+wYr+s1UFV+/INW0kchg4kgf6W61vC4Zg33/C/zI7JlF6jbgWWNMO6AA963h\nq0BP3LeG+zyvUcaYLGPMMqCVtfYx4Me4vvLZuG8hL468ri2uS8gBMR5nn3U28jOLiNSka9C+Nfga\nZIy5G9d17TXgZ8AXuDFibwHn1qj1R5HHGnKcvV8Xxv19XRwZjwau69y7kRah6ey55vTDdb+Mdvt7\nGncNuQR4vkaAjX4OXcskqRSepEGstaW4C8QC+/XZjN7E/fKZFrl/P64P+EzgH8CrQP/Ic5Nx/Yt/\nD6wBZhlj5ke6BNT0C+COyEDeKcBt1tqvcL9Qi3BdG+bgvkH6i7W2Atfv/TJjzGzcL75brbUzrLXz\ngOuBScaYz3EX1isix3kFeN8Ys6/uADcA1xhjPgMuw32rFbV3P/R9Dsw1xrwU6faxt9r22fux6P0n\ncd+qfmSMmYfrd36RtXYbcBfuPH4AnIcbA3BAXXVZaytxF8y/G2MKcV1ULrHWluP+EzIwch4/AO60\n1s6K8Tj7rHNfNYiINISuQXG9Bj0IHGyMmYsbE7YUeMFa+xZuMoW3I5/jfOCMWo5TWw37es1TuDDz\niTFmPq618LzIcxfiAunnwJ+AZUCJMeYq3Ljd2yN1/Rc3Bmo3Xcsk2ULhcEyTgYkEhjHmGeBLa+29\nvmuJhTHmAmCTtXai71pERKRpdA2KP2PMzcDL1lobaSGaDXzbWrvAc2ki31DnhBHGmAzcTCXDgDLc\nFJRL9vG6x4HN1tqbYt1HpBkJ475tExERSbZUuAYtBP5hjKnG/d/0LgUnaahk5Zb6Zts7Hcix1o42\nxhyOW+/g9L0KuBw3leO0WPcRaYoaU7umBKuFAkVE0oauQfFnrf0Xbpp3kaZISm6pb8zTGCJTUEbW\nRhixVwGjgcNw/U9DsewjIiIiIiISZ0nJLfWFpzbAjhr3qyLNWxhjugG/wg3GC8Wyj4iIiIiISAIk\nJbfU121vB5Bf435GjfVizgY64WYo2Q9oaYxZUM8++1RYWKhZK0REAqCgoCAZ0+6nHF2nRESCoY7r\nVFJyS33haTpwCvBSZH2BOdEnrLV/BP4Iu1dqNtbaZ40xZ9a2T10KCgpieVmzU1hYqHNTC52b2unc\n1E7npnaFhYW+Swg0/dwklv5tJp7OcWLp/CZePdeppOSW+sLTBOD4yEJx4BY3Oxdoba19ItZ96itC\nRERERESkCZKSW+oMT9baMG4V7JoW7uN1z9bY3tc+IiIiIiIiCZGs3KKJHERERERERGKg8CQiIiIi\nIhIDhScREREREZEYKDyJiIiIiIjEQOFJREREREQkBgpPIiIiIiIiMVB4EhERERERiYHCk4iIiIiI\nSAwUnkRERERERGKg8CQiIiIiIhIDhScREREREZEYKDyJiIiIiIjEQOFJREREREQkBgpPIiIiIiIi\nMVB4EhERERERiYHCk4iIiIiISAwUnkRERERERGKg8CQiIiIiIhIDhScREREREZEYKDyJiIiIiIjE\nQOFJREREREQkBgpPIiIiIiIiMVB4EhERERERiYHCk4iIiIiISAwUnkRERERERGKg8CQiIiIiIhID\nhScREREREZEYKDyJiIiIiIjEQOFJREREREQkBgpPIiIiIiIiMcjyXUCiFZVU8JndwObtpbTPb8Gw\nAZ1on5/ruywRERERCZKKCigshIULYdMmyM2Fvn1hxAjo0sV3dRIQaRueqqqqeWnKIl6avIjyiqrd\nj2dkhDjm0J788KRBtG+jECUiIiLSrM2cCQ8/DC+/DDt3fvP5UAjGjoWf/ATOOMPdl2YrLcNTRWUV\n9/x1Jh9/sQ6Aof070We/fNZuLubzhRuZMnMlsxZs4IYLRzCkfyfP1YqIiIhI0i1eDNdfD6++uuex\nAw+EggLX0lRaCgsWwEcfwbRp7s/IkfDMMzB4sK+qxbO0C0/hcJiH/zWbj79YR+u8bG68cCTDD+y8\n+/k1m4p4+KXZzFm8iVsfm8EtlxxGwUFdPVYsIiIiIkkTDsMjj8DPfw5lZdCyJVxxhfszYMA3X79t\nGzz3HNx1F3z6qQtXjzwCl1yS/NrFuzrDkzEmA3gEGAaUAZdZa5fUeP4s4AYgDDxvrf1D5PFZwPbI\ny5Zaay9NQO37NLVwJZM/XUlOdia/uWI0B/Rs97Xnu3dqzR3/cwSPvTqXiR8u5/+e+YTfXD6awft3\nTFaJIiIiIuJDcTFccAFMmODun38+3HMPdO9e+z7t2rkuez/8IVx7LTz1FFx6KSxbBnfcoW58AZGs\n3FJfy9PpQI61drQx5nDggchjGGMygbuAAqAYmG+M+RuwC8BaOy72jxsfO4rLeer1LwC48syh3whO\nUZmZGVx55jDCYZg0Yzl3P/spv/vZWDq3z0titSIiIiKSNOvWwcknu0kh2raFxx+Hc86Jff/8fHjy\nSTjsMLjqKvjtbyErC37968TVLA2RlNxS31TlY4BJkTf9GBgRfcJaWwUcZK3dCXQGMoFyYDjQ0hjz\nX2PM5EjxSfHPdxayo7icYQd04tiRvet8bSgU4oozhnLwgZ3ZVlTGXc9+QmVVdZIqFREREZGkWbUK\nxoxxwWn//eGTTxoWnGr6n/+Bf/wDMjLgtttcFz4JgqTklvrCUxtgR437VZEmsWgh1caYM4HPgKm4\n9FYM3GetPQG4Ani+5j6Jsm1nGRNnLAfg0lOHEIqhCTUzM4NfXDCCLu3zWLRyG/98Z2FiixQRERGR\n5Fq7FsaPh6VL3XilGTPcxBBNcdZZrvsewDXXwAcfNL1Oaaqk5Jb6uu3tAPJr3M+w1n6tecZa+4ox\nZgLwF+BC4O/A4shzi4wxm4FuwOq6DlRYWFhPKXWbPHs75RVVHNgjl63rFlO4LvZ9Tzy0NX+ZXMKL\nb1tasZWenXKaVEu8NfXcpDOdm9rp3NRO50YaQz83iadznHjN7Rxnbt+Ouewy8pYtY5cxLLzvPqpW\nroSVK5v+5kOH0uOCC9jvueeoOP10sv7+92Z3fgMmKbmlvvA0HTgFeMkYMwqYE33CGNMGeAM43lpb\nbowpBqqAi3EDtX5sjOmOS4Fr6/u0BQUF9b2kVpVV1Tz0xlsAXHbGSAb269Cg/QuAbZXzePXdJUyZ\nV8bvrh1FZkYwBv8VFhY26dykM52b2unc1E7npna66NdNPzeJpX+bidfsznF5OXzrW25ih8GDaTlt\nGgd3ivMSNU8/DStXkj1tGr3vvpv2kydrAokEquc6lZTcUl93uglAqTFmOm7Q1bXGmHONMT+y1u4A\n/ga8Z4x5H6iO3H8KaGOMeQ94Ebh479QXb598sY5tO8vovV8+B/Vt36j3OO/bB9G5fR5L12xn0ofL\n4lyhiIiIiCRNOAxXXgnvvgvdusHEiRDv4ARuwohnn4X8fNpPnQovvBD/Y0iskpJb6mx5staGgSv3\nenhhjeefAJ7Y6/lK4IL6Pl08vf3JCgBOOLxPTGOd9iU3J4sfnTaEO//yKc9NWsCY4T1ol98inmWK\niIiISDI89JBrFcrLg9dfh169Enes3r3hwQfhssvgpz+Fb38bOjSsF5Q0XbJyS8Incki0ol3lfGY3\nkJERYuyhPZv0XqOGdOPQg7pQXFLB399aEKcKRURERCRpPvoIfvELt/3Xv8KIEXW/Ph4uuYSdBQWw\nebObgU/SVsqHp0/mr6eqOszQ/h1p27ppLUWhUIhLTxlMRgje+ugr1mwqilOVIiIiIpJwW7bA974H\nlZVuQduzz07OcUMhVl53nZu+/JFHYN685BxXki7lw9OMuWsAOGJoHStDN0Dv/dowfkRvqqrDPD9R\nrU8iIiIiKSEchosvhhUr3EK2d9+d1MOXDBgAl18OVVVwww1JPbYkT0qHp9LySmYt2ADAqCH7xe19\nzz3BkJ2VwXufr2bJqm1xe18RERERSZDHH3fjm9q1c4vY5nhYeub226FVK/jPf1z3QUk7KR2e5i/d\nQnllNf17tqVj27y4vW+X9i05aUw/AP7+Xxu39xURERGRBFi2DH7+c7f92GPQt6+fOjp3dpNGAPz6\n135qkIRK6fA0e9FGAA4e0Dnu733WuAHkZGfyyfx1LFuzPe7vLyIiIiJxUF0Nl1wCxcVwzjnuj08/\n/znk58Nbb8EHH/itReIupcPT59HwdGD8w1O7/BacMKoPAC9NXhT39xcRERGROHj4YZg2Dbp0cdu+\ndewI11zjtu+9128tEncpG562F5WxbM12srMyGNivY0KOceYxB5CVGeKD2atZtWFnQo4hIiIiIo20\nbNmeyRkeeywxC+E2xk9+Ai1awBtvgNUQkHSSsuFp3tLNhMMwsG8HWmRnJuQYndrlMX5Eb8Jh+NcU\ntT6JiIiIBEY4DFdfDSUlcO65cPrpvivao0sXuCCy9uqDD/qtReIqZcPTguVbABiyf2JanaLOHj+A\njBBMK1zFpm0lCT2WiIiIiMRowgQ3q13btvC73/mu5pv+93/d7bPPwsaNfmuRuEnZ8GS/2gqA6dMh\nocfp1qkVo4d1p6o6zJvTlyX0WCIiIiISg6KiPeOK7rwT9ovfkjVxM3AgnHgilJbCU0/5rkbiJCXD\nU0VlNYsj6y8d2Ltdwo932tH9AZg0YzmlZZUJP56IiIiI1OG222DVKhgxwi1MG1RXXeVun3zSzQoo\nKS8lw9OyNdupqKymZ5fWtG6Z+AXQDurbAdO7PUUlFUwpXJnw44mIiIhILebMgYcegowM+POfITMx\nY9/j4tvfhl69YMkSmDrVdzUSBykZnvZ02WuftGNGW59ef28J1dXhpB1XRERERCLCYbj2Wqiqcq06\nBQW+K6pbZiZcdpnbfvxxv7VIXKRkeFq4IhKeeicvPI0e1o1O7fJYvbGYwgXrk3ZcEREREYl44w2Y\nMgXat4fbb/ddTWwuucS1kk2YABs2+K5Gmiglw9PSNdsB6N8z8eOdojIzMzjlyH4AvP7e0qQdV0RE\nRESA8nK47jq3fdtt0CGxk4bFTc+ebuKIigp4/nnf1UgTpVx4Kq+oYtWGIjJC0Kdbm6Qe+1uj+tIi\nJ5PPF21k9caipB5bREREpFl75BFYtAiMgSuv9F1Nw1x0kbtVeEp5KReeVqzfSXV1mO6dWydscdza\ntM7L5uiDewAw8cPlST22iIiISLO1efOebnr33w/Z2X7raaiTT4Y2baCwEBYs8F2NNEHKhaflkS57\nfZPc6hR14mjXdW/ypysoq6jyUoOIiIhIs3LbbbBtGxx3HJx0ku9qGi43F84+222r9SmlpVx4WrZ2\nBwD9urf1cvwDerVjQK92FJVU8P5nq73UICIiItJsLFoEjz7qJl343e8gFPJdUeOcd567ff55N2ug\npKSUC0/L10TDk5+WJ4ATR/cFYOKMZd5qEBEREWkWbr3VTU3+wx/C0KG+q2m8sWOhRw9YtgxmzPBd\njTRSSoWncDjMst3d9vy0PAEceXAPWuVls3DFNhav2uatDhEREZG0NmsW/OMf0KKF67qXyjIz4dxz\n3fY//uG3Fmm0lApP24rK2Lmrgla5WXRql+utjtycLI4d2QvQxBEiIiIiCfPLX7rbq6+GXr381hIP\n0XFPr7yirnspKqXC0+oNbnrwHl1aE/Lc3/Xbo/oC8O5nqyguqfBai4iIiEjamTYNJk2C/Hy48Ubf\n1cTHyJGu696qVTBzpu9qpBFSKzxF1lbq0bm150qgV9d8hh3QibLyKt77bJXvckRERETSRzgMN93k\ntq+/Hjp18ltPvGRkwJlnuu2XX/ZbizRKSoWnVTVanoLg+MN6A/DWJys8VyIiIiKSRt54Az76CDp3\nhp/9zHc18VUzPKnrXspJqfAUbXnq2TnfcyXOEcO60yo3i8Urt+2eyEJEREREmqCqCm6+2W3fcovr\ntpdOjjrKhcLFi+GLL3xXIw2UWuEpYC1PLbIzGXtoTwDeUeuTiIiISNO98IILFX36wOWX+64m/jIz\n4bTT3La67qWclAlPFZXVrNuyi1AIunVq5buc3Y4/rA8AUwtXUlFZ5bkaERERkRRWWQm/+Y3b/vWv\n3RTl6eiMM9zt66/7rUMaLGXC07rNxVRXh+ncviUtsjN9l7Nb/55t6de9DTt3VfDRvHW+yxERERFJ\nXS++CAsXQv/+cMEFvqtJnHHjIC/PrWO1dq3vaqQBUiY87RnvFIwue1GhUGh369PbH3/luRoRERGR\nFFWz1emWWyAry289iZSXB8ce67YnTvRbizRIyoSndZt3AcHqshd1TEFPsrMy+HzRRjZs2eW7HBER\nEZHUU7PV6fzzfVeTeCed5G7ffNNvHdIgKROe1m8pBqBrh5aeK/mm/JY5HDGkG+EwTP5UE0eIiIiI\nNEhlJdxxh9tO91anqBNPdLdvvw3l5X5rkZjV+ZNpjMkAHgGGAWXAZdbaJTWePwu4AQgDz1tr/1Df\nPo21PtKiE8TwBHD84b157/PVvPPpCr53vCEjI+S7JBEREZHU8MILsGhR82l1AujdG4YMgXnz4IMP\nYPx43xWltGTllvpank4Hcqy1o4EbgQdqFJAJ3AUcCxwBXGWM6RjZp8W+9mmKDQEPT8MO6EyX9nls\n2FrC3MWbfJcjIiIikhqa01invanrXjwlJbfUF57GAJMArLUfAyOiT1hrq4CDrLU7gc5AJlAe2Wfi\nvvZprHA4HPiWp4yMEONH9AZg8kx13RMRERGJSXNsdYqKdt1TeIqHpOSW+sJTG2BHjftVkeataCHV\nxpgzgc+AqUBxffs0xo7ickrLq2iVm0XrljlNeauEGj+iFwAfzl1LSVml52pEREREAq45tzoBjB4N\nbduCtbB8ue9qUl1Sckt9oWYHkF/z9dba6povsNa+AvQAWgAXxrJPQ+1pdQreTHs1devUioF9O1BW\nXsWMuWt8lyMiIiISbM251QlcWBw3zm1Pnuy3ltSXlNxSX7yfDpwCvGSMGQXMiT5hjGkDvAEcb60t\nN8YUA1V17VOXwsLCWp+b95ULTzkZZXW+LggO6FLNl8vh1SnzaRvaGJf3DPpn9knnpnY6N7XTuZHG\n0M9N4ukcJ16gznF1NYN+9SvygOXnncfm2bN9V9RkjTm/nY2hN7Dln/9k2cEHx7+o5iMpuaW+8DQB\nON4YMz1y/2JjzLlAa2vtE8aYvwHvGWMqgNnA3yKv+9o+MXxYCgoKan1u+fZFwBbM/j0oKBgSy9t5\nYwZVMGnWJJZvKKPX/gPp0r5pY7QKCwvrPDfNmc5N7XRuaqdzU7tA/acqgPRzk1j6t5l4gTvHr7zi\nuqr16UPfW26hb3a274qapNHnNz8f7rmHDp99RodDDoGMlFlJKOnquU4lJbfUGZ6stWHgyr0eXljj\n+SeAJ/ax6977NEnQJ4uoqXVeNqOGdOP9z1czrXAV5xx3oO+SRERERIIlHIY773Tb110HKR6cmmTA\nAOjVC1auhLlzYfhw3xWlpGTllpSItuu3uvDUpX2e50piE504YsrMlYTDYc/ViIiIiATMO+9AYSF0\n6QKXXuq7Gr9CITjuOLf99tt+a5F6pUR42rytBIDOTewClyyHHNiZdvktWL2xiIUrtvouR0RERCRY\n7rrL3f7sZ5CXGl+OJ1Q0PL3zjt86pF6pEZ62lwLQsW2u50pik5mZwTGH9gRc65OIiIiIRHz0EUyd\nCm3awFVX+a4mGI491t2+9x6UlfmtReoU+PBUWl5JUUkFWZkZtGkV3DWe9hbtuvfeZ6upqKzyXI2I\niIhIQERbnX78Y7fGkUDXrjB0KJSUwIwZvquROgQ+PG2p0eoUCoU8VxO7ft3b0q97G4pKKvh0/nrf\n5YiIiIj4N28evP465Oa6LnuyR7TrntZ7CrTAh6dU67JX0/gRvQF13RMREREB4O673e2ll7rJImSP\nY45xt+++67UMqVsKhCc3WUTHtqk3mHDsoT3IyAgx88v1bC9S/1URERFpxpYtgxdfhKwsuP5639UE\nz1FHuZn3Pv4YSkt9VyO1CHx42pTCLU/t83M51HShqjrMe5+t9l2OiIiIiD/33QdVVfCDH0CfPr6r\nCZ727d24p/JyF6AkkAIfnlK55Qlqrvm0wnMlIiIiIp6sWwdPP+22b7jBby1BNnasu33vPb91SK1S\nIDylbssTwOGD96NVbhaLV21nxbodvssRERERSb4HH3RTcJ9xBgwa5Lua4Dr6aHercU+BlQLhybU8\ndUrRlqec7EyOPLgHAFMLV3muRkRERCTJtm2DRx912zfd5LeWoIuGpw8/dN33JHBSIDyldssTwLgC\n13VvWuFKqqvDnqsRERERSaKHH4adO91CsCNH+q4m2Lp0gYED3XpPhYW+q5F9CHR4qqoOs3Wnm6Wu\nfZvUDU+D+nWga4eWbNpeytwlm3yXIyIiIpIcu3bBQw+57Ztv9ltLqoi2PmncUyAFOjztKCqjujpM\nm1Y5ZGcFutQ6hUKh3a1PUwu15pOIiIg0E08+CZs2wWGHwbhxvqtJDdFJIzTuKZACnUi2RdZGapff\nwnMlTTeuoCcAH85ZQ2lZpedqRERERBKsvBzuv99t33STW8NI6hdtefrgAze1uwRKsMNTpMteu9ap\nH566d27NQX3aU1JWxUfz1vouR0RERCSxnn8eVq50s+udeqrvalJHjx7Qv78bJzZ7tu9qZC/BDk9p\n1PIEMG5EtOueZt0TERGRNFZVBffc47ZvvBEyAv1fzuAZM8bdzpjhtw75hkD/JO9ueUqT8HTUwT3I\nygzx+cJBmO97AAAgAElEQVQNbNlR6rscERERkcR49VWwFvr0ge9/33c1qeeII9zthx/6rUO+ITXC\nUxp02wPIb5nDyEH7UR2Gd2ep9UlERETSUDgMd97ptq+/HrKz/daTikaPdrcKT4ET7PAU6bbXPk1a\nnmDPmk9TZmrWPREREUlDb78Ns2a5NYsuucR3Nalp8GDIz4fly2GtxsoHSbDDU6TlqW2atDwBjBjY\nlfyW2Sxfu4Nla7b7LkdEREQkvu66y91eey3k5fmtJVVlZsKoUW5b454CJSXCU7qMeQLIzsrgqIN7\nAGp9EhERkTQzYwZMmwZt2sCVV/quJrWp614gBTs8FblJFdq1zvVcSXyNj8y69+6sVVRVVXuuRkRE\nRCROoq1OV18Nbdv6rSXVadKIQApseKquDrOtqByAdvk5nquJrwN7t6d7p1Zs3VnG7EWbfJcjIiIi\n0nTz5sEbb0BuLlxzje9qUt/hh7uFhQsLoVSzNAdFYMNTUUkF1dVhWuVlk52V6bucuAqFQrtbn9R1\nT0RERNLC3Xe728suc5NFSNO0a+cmjigvdxNwSCAENjxt2xntspc+451qOiYy696MeWvZVVrhuRoR\nERGRJli6FF54AbKy4LrrfFeTPjTuKXCCG56K0m+yiJq6dmjJ4P07Ul5RxYdzNAWliIiIpLD77oPq\najjvPLcwrsRHNDxpxr3ACGx42r7TjXdq2zq9xjvVFF3zaWqhuu6JiIhIilq7Fp55xo3PueEG39Wk\nl+ikEdOnu8WHxbvAhqcdu1x4atMqPVueAI4c3p3srAzmLtnEhq27fJcjIiIi0nAPPghlZXDGGTBw\noO9q0suAAdCxI6xfDytW+K5GCHB42hkJT/ktsz1Xkjit8rIZNaQb4bCbtlxEREQkpWzdCo8+6rZv\nuslvLekoFILDDnPbn3zitxYBghyeiqMtT+nbbQ9gXEFPwHXdC6s5VkRERFLJww9DUREcdxyMGOG7\nmvQ0cqS7VXgKhMCGpx27W57SOzwdYrrQrnULVq4vYvGqbb7LEREREYlNcTH8/vdu++ab/daSzqIt\nT59+6rcOAQIcnqItT/lp3vKUlZnB0Yf0AGBqobruiYiISIp48knYtMkt5nrMMb6rSV/RlqeZM6Gq\nym8tEuDwFJ0wIs1bngDGRRbMfe+zVVRWVXuuRkRERKQe5eVuenJwrU6hkN960lmXLtC3r2vp+/JL\n39U0e1l1PWmMyQAeAYYBZcBl1tolNZ4/F7gGqATmAldZa8PGmFnA9sjLllprL21oYTuL3cKxrdN4\nwoio/j3a0qtrPivX72SW3cBhg/bzXZKIiIhI7Z57DlavhiFD4OSTfVeT/g47DJYvd+OehgzxXU0g\nJSu31NfydDqQY60dDdwIPFCjgDzgN8Ax1tojgbbAycaYXABr7bjInwYHJ2g+Y54AQqEQ4yOtT1Nm\nas0nERERCbCqKrjnHrd9442QEdiOTOlDk0bEIim5pb6f9jHApMibfgzUnEalFDjCWlsauZ8FlADD\ngZbGmP8aYyYbYw6v/7N+XVVVNcUlFYRC0LoZhCeAYw7tSSgEn3yxjqKSCt/liIiIiOzbyy/DokXQ\nrx9873u+q2keNF15LJKSW+oLT22AHTXuV0WaxLDWhq21GwGMMT8BWllr3wGKgfustScAVwDPR/eJ\nVTQ8tMrNJjOjefSh7dQuj2EHdKKisprps1f7LkdERETkm8JhuPNOt33DDZBV5wgQiZdDD3UtfHPn\nQkmJ72qCKim5pb6f+B1Afo37Gdba3TMaRN78XuAA4KzIwwuBxZFCFxljNgPdgDoTQWFh4e7tjdtd\neMrJDH/t8XTXr2MVsxfB69O+pFP25t2PN6dz0FA6N7XTuamdzo00hn5uEk/nOPGaeo7bfPABA2bP\nprxTJ+YNHUpYf2dfk8if4YH770/LxYtZ8OKLFA8blrDjpLCk5Jb6wtN04BTgJWPMKGDOXs8/hmsG\nO8NaG13h9WLcQK0fG2O641Lg2nqOQ0FBwe7t+cs2A+vp3CH/a4+nu0FDKpk4axIrNpbTo+9B7Nex\nFYWFhc3qHDSEzk3tdG5qp3NTO/3HtW76uUks/dtMvLic42uuASDnhhs4dPToOFSVPhL+Mzx2LCxe\nzEE7dkAz/bdSz3UqKbmlvvA0ATjeGDM9eoDITBWtgZnAJcB7wBRjDMBDwFPAM8aY96L71Ex9sWgu\nazztLa9FFkcM7ca0wlVMLVzFud8yvksSERERcd5/H6ZPh/bt4fLLfVfT/IwcCU89pXFPtUtKbqkz\nPEVS2ZV7PbywxnZmLbteUNf71mfn7pn20n+a8r2NK+gVCU8r+f7xB/ouR0RERMSJjnX66U8hP7/u\n10r8adKIOiUrtwRybskdkTWemlvLE8DwAZ3p0KYFazcVY7/a6rscEREREZg1CyZNglat4Cc/8V1N\n8zRkCOTmwuLFsGWL72qarUCGp2jLU5tmMk15TZkZIcYeGlnzqVBrPomIiEgA3HWXu73iCujY0W8t\nzVV2NhxyiNueOdNvLc1YIMNTcWSq8uayxtPeogvmvv/ZaiqrwvW8WkRERCSBFixwazvl5MD//q/v\napq36EQRs2b5raMZC3R4apXX/MY8AfTt1oZ+3dtQVFLBwjWl9e8gIiIikij33uvWd/rhD6F7d9/V\nNG/R8KTZUb0JZHgqKo0uktt8F16Ltj7NWVbsuRIRERFptlasgOeecwu0/uIXvquRQw91t2p58iaQ\n4WlXM295Ahh7SE8yQrBwTSk7IlO3i4iIiCTV/fdDZSV8//vQv7/vamTQIDdpxNKlsFUTi/kQyPBU\nXKrw1L5NLgebLlRXw/uf17rIsYiIiEhibNgATzzhtm+80W8t4mRlwfDhbvuzz/zW0kwFMzxFW55y\nm294Ahhf4LruTZ2pWfdEREQkyR58EEpL4ZRTYOhQ39VIVLTrnsY9eRHM8FRaCTTvlieAw4fsR05W\nCLtiK6s3FvkuR0RERJqLLVvgT39y2zff7LcW+TrNuOdV4MJTRWU1ZeVVZGSEyM2pbSHg5iE3J4tB\nvfMAtT6JiIhIEv3+91BUBMcfD6NG+a5GalLLk1eBC0+7Svd02QuFQp6r8W94v5YATC1cSXW11nwS\nERGRBNu+3YUngF/9ym8t8k2DB7s1txYtcn9XklSBC0971nhqvtOU19SnSws6tctjw9YS5i/b7Lsc\nERERSXd//KP7T/kxx8CRR/quRvaWk7NnDNrnn/utpRkKXnjSTHtfkxEKMa6gJwBT1HVPREREEmnn\nTjdRBKjVKci0WK43wQtPmmnvG8ZFZt2bPmcNZRVVnqsRERGRtPXII26yiDFjXMuTBJMmjfAmgOFJ\nM+3trVfXfAb0aseu0ko+mbfOdzkiIiKSjoqL3aK44FqdNPY8uDRphDfBC0+lannal2jr05RCdd0T\nERGRBHjsMdi0CQ47zM2yJ8E1dKhbMNda19VSkiZ44alEY5725ehDepCZEWKW3cDWnaW+yxEREZF0\nUlIC997rttXqFHwtWsCQIRAOw+zZvqtpVhSeUkTb1i0oOKgr1dVh3v9ste9yREREJJ08+SSsX++6\ng514ou9qJBaaNMKL4IWn3d32NFX53saPcF33JmvWPREREYmXsjK45x63feutanVKFZo0wovghSe1\nPNVq5KCutMrLZunq7Sxbo0XRREREJA6eeQZWr4Zhw+DUU31XI7HSpBFeBDA8aba92uRkZzLuULfm\n01sff+W5GhEREUl55eVw111u+5ZbICNw/zWU2gwb5v6+FixwY9YkKQL3L2RXmWt5aqlue/t0/OF9\nAJhWuIpyrfkkIiIiTfHss7BiBQwaBGed5bsaaYi8PDjoIKiqgi++8F1NsxG48FRS5lqe8looPO3L\n/j3a0r9nW4pKKpgxd63vckRERCRVlZXBb3/rtm+9Va1Oqejgg93t55/7raMZCdy/kpJShaf6HH+Y\na316+xN13RMREZFGeuop1+o0eDCcc47vaqQxFJ6SLnDhqbRc4ak+Yw/tSU5WBrMXbWLd5mLf5YiI\niEiqKS2F//s/t3377Wp1SlXR8PTZZ37raEYC9y9F3fbq1zovm9HDuwPwzicrPFcjIiIiKefxx2HN\nGhg+HM44w3c10ljR8DR7NlRX+62lmQhUeKquDlNS5iZByM1ReKrLtyJd9yZ/uoKq6rDnakRERCRl\n7NoFd97ptu+4Q61OqaxzZ+jRA4qLYckS39U0C4H61xLtspebk0lGhhZoq8uQ/h3p1qkVm7aX8pnd\n4LscERERSRWPPgrr18OIEXDKKb6rkabSuKekClR4inbZy1WXvXqFQiGOP6w3oDWfREREJEZFRXDP\nPW77jjsgpC+rU57CU1IFKjyVlrsuexrvFJvxI3qREYJPvljHtp1lvssRERGRoPvTn2DjRhg1Cr79\nbd/VSDwoPCVVoMKTpilvmI5t8ygY2JWq6jBTC1f6LkdEREQCLKOoCO67z91Rq1P6UHhKqmCFJ820\n12DfOtxNHPHWx18RDmviCBEREdm3Li++CFu2wFFHwXHH+S5H4mX//aF1azd74gaNg080hacUN2Jg\nV9rlt2DVhiIWLN/quxwREREJom3b6Pq3v7lttTqll4wMN+U8uCnLJaEUnlJcVmYGx47oBWjiCBER\nEanF735HVlERjB8PxxzjuxqJN3XdS5o6U4oxJgN4BBgGlAGXWWuX1Hj+XOAaoBKYC1wFhOrapy4K\nT43zrcP78PLUxbw/ezWXnjaE1nnZvksSERGRoNiwAX73O7d9xx1+a5HEUHhKWm6pr+XpdCDHWjsa\nuBF4oEYBecBvgGOstUcCbYGTI/u02Nc+9VF4apzunVsz7IBOlJVXMU0TR4iIiEhN//d/UFzMtqOO\ngjFjfFcjiaDwBEnKLfWFpzHAJABr7cfAiBrPlQJHWGtLI/ezIo+NASbWsk+dFJ4a7zuj+wIwccZy\nTRwhIiIizvLlblHcUIjVP/6x72okUYYMgcxMWLAASkp8V+NLUnJLfeGpDbCjxv2qSJMY1tqwtXYj\ngDHmJ0Ara+3bde1TH4Wnxjt8cDfa5bdgxbqdzF+2xXc5IiIiEgS//jVUVMD551N6wAG+q5FEyc2F\ngQOhuhrmzvVdjS9JyS31pZQdQH6N+xnW2uroncib3wscAJwVyz61KSwsZOVqN1vchvWrKSzcXt8u\nzUZhYWFMrxvaO4f3vyjj+TdncdboDgmuKhhiPTfNkc5N7XRupDH0c5N4Osfxlbt4MYOee45wVhZf\nnH02oHOcaD7Pb99eveg4bx5fvfYamzIzvdXhUVJyS33haTpwCvCSMWYUMGev5x/DNXmdYa0Nx7jP\nPhUUFDDly5lAMQcN2J+Cgl6x7Jb2CgsLKSgoiOm1Pfvt4oP5b/PlylIOMENo27pFgqvzqyHnprnR\nuamdzk3t9J+quunnJrH0bzMB7rgDwmFCV1zB0FNP1TlOMO/n99hjYeJE+mzdSp80/Xuu5zqVlNxS\nX3iaABxvjJkeuX9xZKaK1sBM4BLgPWCKMQbgoX3tU18RUeq21zRdO7Sk4KCuzPxyPZM/XcmZ49Q8\nLyIi0ix9+CG8/jq0agW33OK7GkkGTRqRlNxSZ0qJpLIr93p4YY3t2toE994nJrvDU67CU2N954i+\nzPxyPZM+Ws7pY/uTkaFF8ERERJqVcBhuvNFtX3stdO3qtx5JjuhCuXPmQFWVm0CiGUlWbgnUIrml\nkfCUm6Pw1FgFA7vSqV0eazcVM2fxRt/liIiISLJNmgTvvw8dOsB11/muRpKlUyfo2ROKi2FJTEus\nSiMEKjyVlFUB6rbXFJkZIU4Y1Qdw05aLiIhIM1JdDTfd5LZvvhnatvVbjySXuu4lXKDCU1mFC08t\ncppXM2O8HX9YbzIyQnw0bx2btzfbuf5FRESanxdegNmzXQvEVVf5rkaSLRqePvvMbx1pLFjhqTwS\nnrIVnpqiY9s8Dh+8H9XVYd75ZIXvckRERCQZSkpcaxPA7bdDXp7feiT5ouFp9my/daSxYIUntTzF\nzXeO6AvApI++oqqq3mW2REREJNX94Q+wYgUMGwYXXeS7GvFB3fYSLjDhqbo6THkkPOVkKTw11fAB\nneneqRWbtpXw8RfrfJcjIiIiibRxI9x5p9u+//5mN9OaRPTrB/n5sHYtbNjgu5q0FJjwtDs4ZWdq\neu04yMgIcdKYfgD8+4NlnqsRERGRhLrjDtixA77zHTj+eN/ViC8ZGXumLFfXvYQITHiKdtnLVZe9\nuDl2ZG9yczKZu2QTy9fu8F2OiIiIJIK18Oc/u/8433uv72rEt2h4Ute9hAhOeCrXeKd4a5WXzfgR\nvQD49wdLPVcjIiIiCXHDDVBZCZdeCkOG+K5GfNO4p4QKTniq0Ex7iXDykfsDMG3WKop2lXuuRkRE\nROLq3XfhtdegVSs3w56IwlNCBSY8lZZXAmp5irdeXfM5+MDOlJVX8bamLRcREUkf1dVw3XVu+xe/\ngG7d/NYjwTB4sJswxFo3fb3EVWDCk9Z4SpxTIq1P/56+jKrqsOdqREREJC5efBFmznSh6ec/912N\nBEVeHhgDVVXwxRe+q0k7wQlP6raXMAUDu9K1Q0s2bNnFzPmatlxERCTl7doFN93ktn/7W9dtTyRK\nXfcSJjjhKdLylNsiy3Ml6SczI8TJR2rachERkbRx331uQdzhw7UgrnyTwlPCBCc8qeUpoY47rA8t\ncjL5fNFGVq7f6bscERERaayvvoK773bbf/iDFsSVb4qGJ631FHeBCU+lmqo8oVrnZTOuQNOWi4iI\npLzrr4fSUvje9+Doo31XI0FUc6Hc6mq/taSZwIQnTRiReNGue5NnrtS05SIiIqlo2jR46SU3KcB9\n9/muRoKqSxc3kcjOnbBMQzbiKTjhqUJTlSdan/3a7J62fOKM5b7LERERkYaorIRrrnHbN94IvXr5\nrUeCTeOeEiI44Und9pLijLEHAG7iiIpKNeOKiIikjCeegDlzoE8f13VPpC4KTwkRnPC0e8IIzbaX\nSIeYzvTZL58tO0p5//NVvssRERGRWGzZArfc4rYfeMB12xOpS81xTxI3wQlPanlKilAoxOlj+wMw\nYdoSwmEtmisiIhJ4v/qVC1DjxsGZZ/quRlKBWp4SInjhSRNGJNzYQ3vSPr8Fy9fuYPaijb7LERER\nkbrMnQuPPgoZGfD730Mo5LsiSQUHHAAtW8LKlbB5s+9q0kZwwlOFWp6SJTsrk5MiM+9NeHeJ52pE\nRESkVtXVcOWVe26HDvVdkaSKzEwYNsxtq+te3AQnPKnlKam+c0Q/crIzmbVgA1+t2+G7HBEREdmX\nZ5+F6dOha1f47W99VyOpRuOe4i444SnS8pSrlqekaNMqh+NGuilOX1Prk4iISPBs3rxnVr0HHoB2\n7fzWI6lH457iLjjhqdyt85SjlqekOW1sf0IhmFq4iq07Sn2XIyIiIjXddJMLUOPGwQ9+4LsaSUUK\nT3EXmPBUHllzSOEpebp3as3hg/ejsqqaN6dr9WkREZHA+Ogjt65TdjY88ogmiZDGGTrU/ezMnw9l\nZb6rSQuBCU/RBVtzsgJTUrNwxjGRRXOnL2NXaYXnakRERITKSrjiCrd9/fVw0EF+65HU1aoVDBjg\nfqa+/NJ3NWkhMEmlIjLmKTtLLU/JNKhfRwb160BxSQWTZnzluxwRERF5+GE3wL9vX/jlL31XI6lO\nXffiKjDhaU+3vcCU1Gx899gDAXjtvcWUR0KsiIiIeLBmDdx6q9v+4x/dOj0iTaHwFFeBSSoVldGW\np8CU1GwUHNSFft3bsGVHGZNnrvRdjoiISPN19dWwcyecdhqcfLLvaiQdKDzFVWCSSnmFa3lSt73k\nC4VCfHe8a316ZeoiqqqqPVckIiLSDL38MkyYAPn58Kc/+a5G0kXNtZ7CYb+1pIHAhKeq6jChEGRl\najYZH0YP7063Tq1Yt3kXH8xe47scERGR5mXrVtfqBHDPPdCzp996JH106wadO8O2bbBihe9qUl5g\nwhO4VqeQpuL0IjMjxFnj3Mx7/5qyiLC+mRAREUme66+HdetgzBi4/HLf1Ug6CYXUdS+Osup60hiT\nATwCDAPKgMustUv2ek1L4G3gEmutjTw2C9geeclSa+2lsRSjacr9Gj+iF3//r2X52h18+uV6Dhu0\nn++SRERE0t+UKfDUU5CTA08+CRn6/5DE2cEHw9tvu/B02mm+q0mIZOWWOsMTcDqQY60dbYw5HHgg\n8li0gBHAn4HuQDjyWC6AtXZcDJ/zazTTnl/ZWZmccUx/nnr9C156ZyEjB3ZVS6CIiEgi7doF//M/\nbvvWW7WmkyRGdNxTerc8JSW31JdWxgCTIm/6MTBir+dzIkXZGo8NB1oaY/5rjJkcKT4mmizCvxNG\n9SW/ZTYLvtrKvCWbfZcjIiKS3m6/HZYsgSFD4Be/8F2NpKtot73Zs/3WkVhJyS31hac2wI4a96si\nTWJECvvQWrtqr32KgfustScAVwDP19ynLmp58i+vRRanHt0fgBfesvW8WkRERBqtsBAeeMCNSXny\nSddtTyQRjIEWLWDZMjdxRHpKSm6pr9veDiC/xv0Ma21981gvBBZHilxkjNkMdANW17MfFeVlFBYW\n1veyZifZ56RX62pys0PMXbKJf/1nOv265ib1+A2hn5fa6dzUTudGGkM/N4nXnM5xqLycgeefT15V\nFevPPZdVWVkuTCVYczrHPgT5/B7Uvz+t5s/HvvQSRYce6rucREhKbqkvPE0HTgFeMsaMAubUVzVw\nMW6g1o+NMd1xKXBtDPvRrk1rCgoKYnlps1FYWOjlnKwssjw/aQGFy8OcfWIw/058nZtUoHNTO52b\n2gX5oh8E+rlJrGb3b/Omm2DpUhgwgK5PPknXli0Tfshmd46TLPDn94gjYP58TEkJBLnOOtRznUpK\nbqmvn9wEoNQYMx036OpaY8y5xpgf1bHPU0AbY8x7wIvAxTGkPgBysjXmKShOOXJ/WuVlM2/JZuYu\n3uS7HBERkfTx0Udw771uVr2//AWSEJxEmsG4p6Tkljpbnqy1YeDKvR5euI/XjauxXQlcUNf71iZb\nU5UHRqu8bM4Y25+/TVrA8/9dwF39x2jmPRERkaYqKYEf/hCqq93aTqNH+65Imos0X+spWbklUGlF\n4SlYTj5yf1rnZfPF0s3MXaLWJxERkSb75S/BWhg0CO64w3c10pwMG+Zu582Digq/taSwQKWVHE1V\nHiit8rI5faybee/v/7WEw2HPFYmIiKSw99+Hhx6CzEx49lnIDe6ETJKG2rSB/feH8nJYsMB3NSkr\nUOEpW1OVB84pR6n1SUREpMmKiuDiiyEcdpNFjNh7CRqRJEj/cU8JF6i0opan4GmZm83px7jWp79N\nXKDWJxERkcb42c/cYrjDhsGtt/quRpqrNB/3lAyBCk9qeQqmU47cn/yWOXy5fAuFCzb4LkdERCS1\nvPwyPPWU66b3979rMVzxR+GpyQKVVtTyFEwtc7M557gBADz75nyqq9X6JCIiEpNVq+BHkZmS77sP\nBg/2W480b8OHu9vPP3ddSKXBAhaeAlWO1HDi6H50apvL8rU7eO/zWhddFhERkajqarjwQti6FU48\nEX78Y98VSXPXqxe0bw+bN8Nq/X+uMQKVVrK1SG5g5WRn8oMTDgLg+UlfUlEZ07rHIiIizdf998PU\nqdClCzz9NGi9RPEtFNKkEU0UqPCklqdgGz+iFz27tGbd5l289fFXvssREREJrsJCuOUWt/3MM9C1\nq996RKI07qlJApVWtEhusGVmZnDBdwYC8I+3LaVllZ4rEhERCaDiYvjBD9xCpD/5ieuyJxIUNcc9\nSYMFKq0oPAXfEUO7MaBXO7buLOOND5b6LkdERCR4rrkGFi50k0Pcc4/vakS+Ti1PTRKotJKZGahy\nZB9CoRAXnTgIgJenLGLnrnLPFYmIiATIc8+5aclbtHDTkufl+a5I5OsGDoTsbLfu2M6dvqtJOYFK\nK1kKTylh+IGdOXhAZ4pLK3nxbeu7HBERkWCYPx+uuMJt/+lPbkFckaDJyXGtouEwzJ3ru5qUE6i0\novCUOi45dTChELz5wTJWbyzyXY6IiIhfxcXw3e/Crl1w/vlw6aW+KxKpncY9NVqg0orGPKWOft3b\nctzI3lRVh3nmjS98lyMiIuJPOAxXXulangYOhEcf1bTkEmwa99RogUorWZn6RZNKzv/OQHJzMvn4\ni3XMWbzRdzkiIiJ+PP20G+vUsiW89BK0bu27IpG6aa2nRgtUeNKEEamlQ5tczh4/AICnXvuCquqw\n54pERESSbM4cuPpqt/3II24siUjQRbvtzZkDlVp6piEClVayFZ5Szmlj+9OpbS5L12xn6syVvssR\nERFJnu3b4eyzobTUjXG66CLfFYnEpn176N3b/ewuWuS7mpQSqLSiCSNST25OFhee5KYuf27ifC2c\nKyIizUN1tZsYYtEiN6veH//ouyKRhtG4p0YJVFrJ0oQRKWnsIT05oFc7tuwo419T9e2FiIg0A7ff\nDv/+t/sGf8IEreckqUfjnholUGlFE0akpoyMEJedOgSAV6YuZu2mYs8ViYiIJNCrr8Idd0BGBvzj\nH7D//r4rEmk4tTw1SsDCU6DKkQYYvH9HxhX0pKKymide04JrIiKSpr78Ei680G3fdRccf7zfekQa\nKxqeCgvddPsSk0ClFYWn1HbxyYNpmZvFp/PX88n8db7LERERia/t2+H002HnTjjnHLj+et8ViTRe\n377QoQNs2gQrNelXrAKVVhSeUlv7Nrn84ISDAHh8wlzKKqo8VyQiIhInVVVw3nmwcCEMHerWdtJC\nuJLKQiEoKHDbhYV+a0khgUormjAi9Z08ph999stn/ZZdvDJFk0eIiEiauO46ePNN9039hAnQqpXv\nikSabsQIdztzpt86Ukig0oomjEh9mZkZXHHmMAD+NWUR6zZr8ggREUlxf/4zPPQQZGe74NS/v++K\nRPO2o1UAACAASURBVOJD4anBAhaeAlWONNKQ/p0Ye0hPyiureeLVeb7LERERaby33oKrr3bbTzwB\nRx/ttx6ReKoZnjRpREwClVYyM9TylC4uPmUQeS2y+GT+Oj6cs8Z3OSIiIg03fz5897tuvNNNN8FF\nF/muSCS+evWCzp1hyxb46ivf1aSEwISnrMwMQhp4mTY6ts3johMHAvDYhDkUlVR4rkhERKQBNm6E\nk0+GHTvg7LPht7/1XZFI/NWcNEJd92ISmPCUnaXglG6+M7ofB/Vpz5YdZfz1zfm+yxEREYlNaamb\nknzZMhg5Ep591i2IK5KONO6pQQLzm0DjndJPRkaIq885mKzMEBNnLOeLpZt9lyQiIlK3qio4/3z4\n8EPXpem116BlS99ViSSOwlODBCaxZCo8paU++7XhrHEDAPjTS59TUam1n0REJKDCYbjmGnj5ZWjb\n1k1N3q2b76pEEisangoLNWlEDAKTWNTylL7OOe5AenRuzaoNRbw0WWs/iYhIQN11Fzz8MLRo4Vqc\nhg71XZFI4nXvDl27wrZtsHSp72oCLzCJJVvhKW3lZGdy9XeHA/DS5IV8tXaH54pERET28swz8Mtf\nugH0zz8PY8f6rkgkOUIhdd1rgKy6njTGZACPAMOAMuAya+2SvV7TEngbuMRaa2PZZ5+FaMKItDak\nfydOGNWH/370FQ++OIv7f3q0WhtFRCQY/v3/7d13fFRV+sfxz2RSSYGEHkUEgUNRivSiUiyIfe1d\nV9lV11V33bWsFV12dVVc18LP3ta1Y0NFVASkiBCaIh567y2hhpT5/XFuSIiBBGRyb5Lv+/Wa15Q7\nd+4zh8kwzz3nPGckDB7sbj/1FJx7rr/xiFS2Ll3cMNVp0+DCC/2O5qBUVt5S3q/Xs4F4a20v4A7g\nsVIBdAHGA82ASEX22Rf9kK7+fntGOxqkJ7FwRbaG74mISDCMHw8XXOAKRdx9N9xwg98RiVS+kvOe\nqq5KyVvKy1h6A6MArLVTgC6ltsd7B7UHsE+ZVDCi+quVGMfNF3UC4O0vLQtXbPE5IhERqdG+/x5O\nOw127nQ9Tw884HdEIv4oWuspKwsKC/2N5eBVSt5SXsaSBpScoFLgdW/hHWSStXbFgeyzL5rzVDO0\nb1Gf0/s0o6Awwr/fmqHqeyIi4o9Zs+CUU2DbNrjkEhg+3M39EKmJGjd2hSNycmDBAr+jOViVkrfs\nd86T92KpJe7HWGvLS0cPZh927thOVtXuKoya6tYu7TMLmZQSy5LVOQx7dRwndqx90K9V3drmUFLb\n7JvaRg6GPjfRV1ltnLBkCWbwYOK2bGFz374suukmmDmzUo7tN32Oo6sqt+9RLVpQZ9UqFr/zDptO\nPdXvcA5GpeQt5SVPE4EzgHeNMT2A2eU8/2D3oU6dNDoXdRnKHllZWdWyXeo02MQdT3/LpLlbOevE\njrRumnHAr1Fd2+ZQUNvsm9pm36ryf/qVQZ+b6Kq0v81Fi+DMM2HzZjjlFNI/+ojOCQnRP24A6Psv\nuqp8+558MowfT7P162kW0PdRzv9TlZK3lDdW7gNglzFmIm4C1Z+MMRcbYwYfyD4VCJxwjLrKa5I2\nzTI4+4QWFEbgsTey2LErz++QRESkulu+HAYMgFWr4PjjYcQIt6aTiECPHu76u+/8jePgVUrest+e\nJ2ttBLi+1MPzynhev3L2KVc4RnOeaprLTm3NzPnrWbQym+EjZnPrJcE8yyEiItXAkiXQr5+77t7d\nlSevVcvvqESCo2tXN+9vxgzYtQsSE/2O6IBUVt4SmIwlHFbPU00TFxvmr5d1JiE+zNisFXyTtdzv\nkEREpDpauNAtertkCXTrBp9/Dqmp5e4mUqOkpUHbtpCX5xIoKVNgkqcYDdurkQ5vkMrgs44BYPj7\ns1i9YbvPEYmISLUyb55LnJYtg169YPRoSE/3OyqRYCoaujdlir9xBFhgkifNeaq5Tu5+BL3bZ7Iz\nt4BH35hGfkGVXV9ARESCZO5clzitXOnmOI0aBbUPvsKrSLVX9ec9RV1gkif1PNVcoVCIG8/vQP30\nJOYt28J/P5/rd0giIlLV/fCDS5zWrIH+/eGzzzRUT6Q8Sp7KFZjkST1PNVtKrXhuvaQzMSF4/5sF\nfD9njd8hiYhIVTVjhisOsX69K788ciQkJ/sdlUjwtWnjTjIsXQqrV/sdTSAFKHkKTCjik3bN63L5\noLYADHtzOms2av6TiIgcoHHjoG9f2LgRBg2Cjz6CpCS/oxKpGsJhV1QFNO9pHwKTsajnSQDO7deC\n7u0asX1nHv98dSq78wr8DklERKqKDz+EU06BnBw47zy3jlMVK7cs4jsN3dsvJU8SKKFQiFsuPpZG\ndWuxaGU2z37wg98hiYhIVfDii3DuuZCbC9ddB2+9pQVwRQ5G9+7uWslTmQKTPKlghBRJSYrjziu7\nER8bw+gpS/nq+6V+hyQiIkEVicBDD8G110JhIdx3HzzzjBt+JCIHrih5mjYN8vP9jSWAlDxJIDU/\nrDbX/aY9AMPfn828ZZt9jkhERAInPx9uugnuvBNCIXjySbj/fndbRA5OgwbQvDls3w5z5vgdTeAE\nJnnSsD0p7aTuTTmlR1N25xcy9OUpbMze6XdIIiISFFu3wllnwVNPQXw8/O9/cOONfkclUj0UzXua\nPNnfOAIoOMlTODChSID8/pz2tGtel005uQx9+XtyVUBCRERWrIDjjnNrN2VkwFdfwUUX+R2VSPXR\nq5e7njDB3zgCKDAZi3qepCxxsTHceWVXGmTUYv7yLTz59kwikYjfYYmIiF9mzHBzMmbNgpYt3aT2\n447zOyqR6qXob+rbb/2NI4ACkzxpzpPsS+2UBO75bXeSEsKMm7GC98bM9zskERHxw8iR7kfdqlXu\nevJkl0CJyKF19NFQpw4sW+Yuskdgkif1PMn+HNk4jVsv6UwoBK9/PpdJs1f5HZKIiFSWSAT+8Q84\n80w3if2yy+DLL6FuXb8jE6meYmKgd293W71Pe1HyJFVG96Mbc8WgtkQi8NgbWSxbn+t3SCIiEm1b\nt7oFb++6y91/4AF47TWt4SQSbRq6V6bAJE8xMYEJRQLs3H4tOLXnkezOL+TNcRtZsW6r3yGJiEi0\nzJ/vqn6NGAFpafDxx3DPPSpFLlIZ+vRx1yoasZfAZCzqeZKKCIVC/P6cY+jatiE7dxdy3/PfsTln\nl99hiYjIofbpp9C1K/z0E7RpA1Onwumn+x2VSM3RpYvr4Z0zBzZu9DuawAhM8qSCEVJR4XAMt13W\nhcyMONZt2sGQF79jZ65WwBYRqRYKC+Hvf4czzoDsbDjnHJgyBVq18jsykZolIcFVtgSYONHfWAIk\nMMmTep7kQCQmxHJJ33o0qluLhSuy+cfL37Nba0CJiFRt69bBaae5oXkAQ4fCe+9Baqq/cYnUVJr3\n9AtKnqTKSkkMM2RwT+qkJDBz/nr+9fo08gsK/Q5LREQOQuqUKdChA4wa5Ra+HTkS/vY3V/VLRPyh\n5OkXAvONFA4reZIDl1k/hQd+35PkpDimzFnDv9+cQWGhFtEVEaky8vLgzjtpeeONsGYNnHCCWwB3\n0CC/IxORnj3dCYysLLdMgAQneVK1PTlYzTJrM2Rwjz2L6A4fMZtIRAmUiEjgLV4Mxx8PDz3kKugN\nGQJffw2HH+53ZCICrsplp06Qnw+TJvkdTSAEJmMJq+yo/AqmaQb3/LYH8bExjJq8hJc+maMESkQk\nyN55Bzp2hO++gyZNmPfss3DvvRAO+x2ZiJTUv7+7/vprf+MIiMAkTzEatie/0jEt6nHnVd2IDYf4\ncNxCXvxYCZSISOBs2waDB8OFF0JOjqumN3Mm2zp18jsyESnLgAHuWskTEKDkSQUj5FDo0qYht1/R\nldhwiI/GL+S5D39QAiUiEhTjxkH79vDCC64M8jPPwPvvuwIRIhJMffpAXJyb97R5s9/R+E7Jk1Q7\nPY5u7PVAxTBywmKGj5itIhIiIn7asQNuuQX69nXznDp1coveXn+9m+skIsGVnAw9ekAkAmPH+h2N\n7wKUPAUmFKkGurVtxN2/7UZcbAyfT1rC0+/NUgIlIuKHSZPc3KYnnnDzme67zy16e8wxfkcmIhWl\noXt7BCZjUe4kh1rn1g2557fdiY+NYfSUpTz2RhZ5+VoHSkSkUuTkwB/+4Ib8zJ8P7dq5pOn++90Q\nIBGpOoqSpzFj/I0jAAKTsoTUbS9R0Mk04P7BPUlKiGX8zJU8+OJ37MzN9zssEZHq7aOPoG1bN6cp\nHHaL3WZlQefOfkcmIgejWzc3fG/uXFi1yu9ofBWY5ClGyZNEyTEt6vGPG3pTJyWBGfPWc9fwiWRv\ny/U7LBGR6mf1ajj/fDj7bFi50v3gysqCoUNdgQgRqZri492abFDje58Ckzwpd5JoanF4HR7+Yx8a\nZtRi/vIt3P7Ut6zbtMPvsEREqoe8PHj8cTAG3nvPnaF+4gk336l9e7+jE5FDoWjo3pdf+huHzwKU\nPCl7kujKrJfCv/54HM0y01i5fjt/+c945i1TyU0RkV9lzBhXEOLPf4atW+GMM+Cnn+Cmm7TgrUh1\ncsop7nrUKCisuXPIA5M8adieVIaMtET+eUMf2reox+atudz59AQmzFrpd1giIlXP8uVwwQXubPRP\nP8FRR8HIkfDxx3DEEX5HJyKHWrt27m973TqYPt3vaHwTu7+NxpgY4BmgPZALXGutXVhi+xnAPUA+\n8JK19gXv8elAtve0Rdbaa8oLJBSYNE6qu+SkOIb8rifD35/N6ClLefi1aawcuI0LTmylHlARkfJs\n3w7DhsFDD7n1m5KS4O67Xc9TYqLf0YlItIRCMGgQ/N//wWefQZcufke0l8rKW/abPAFnA/HW2l7G\nmO7AY95jGGPigGFAF2AHMNEY8xGwFcBa26/ib1c9T1K5YsMx3Hh+Bw5vkMLLI+fw31E/s2L9Nv54\nfkfi4zTMRETkF/Lz4ZVX4N57XWEIcMUhHn1UPU0iNUVR8vT55+67IFgqJW8pr7+nNzDKe9Ep3gGL\ntAEWWGuzrbV5wATgBKADUMsY84Ux5msv+HIpd5LKFgqFOKdvC+66qhuJ8WHGZq3g9qe+Za0KSYiI\nFItE4NNP3bymwYNd4tS5s5vr9M47SpxEapL+/V3lvSlTYMMGv6MprVLylvKSpzQgp8T9Aq9LrGhb\ndoltW4HawHbgEWvtKcB1wBsl9tknDZcSv3Q/ujH/+uNxNKpbiwUrsvnT42OZ/vM6v8MSEfHf1Knu\nx9Lpp8OcOdCsGbz5Jnz/PfQ7oAEmIlIdJCdD377upMoXX/gdTWmVkreUN2wvB0gtcT/GWltUXiO7\n1LZUYDMwD1gAYK2db4zZCDQG9jsr/+ef55K9Nr6ccGqmrKwsv0MIrEPZNlf2q8OISYXMX7WL+56f\nTL/2aRzXLrXKDinV52bf1DZyMGrS5yZxwQIav/ACGV99BUB+7dqsvuYa1p93HpH4eJgxIyrHrUlt\n7Be1cXTVhPZtcMwxNBk9mo3//S9LWrf2O5ySKiVvKS95mgicAbxrjOkBzC6x7WegpTEmHZe1HQ88\nAlyNm6j1B2NMJi7TW13OcWjXti3NMmuX97QaJysri85akb1M0WibXt0jvP3VPN4c/TPfzM5ha14S\nt1zUidopVWtxR31u9k1ts2814T/9X6NGfG5mzYIHHoARI9z9xES4+WZi77iDJnXq0CSKh9bfZvSp\njaOrxrRvWho89hh1p0yhbvv2EBdXaYcu5/+pSslbyhtO9wGwyxgzETfp6k/GmIuNMYO98YJ/Br4A\nJgEvWmtXAy8CacaY8cBbwNUlsr59B1JFz+5L9RITE+Likw33XtODlKQ4ps1dyx8f/YYZVsP4RKQa\nmzEDzjnHzWsaMQISEtw6TQsWuKp6der4HaGIBEXLltCmDWzeDOPG+R1NSZWSt+y358laGwGuL/Xw\nvBLbRwIjS+2TD1xegTe4F+VOEiRd2jTkP7f247H/ZTFn0UbufW4y5/RtweWntiEuVnX1RaSayMqC\nIUPgk0/c/cREuO46uO02aNzY39hEJLh+8xsYOhTefx9OPNHvaIDKy1sC8ytQBSMkaOqnJzH0+t5c\nNrA1MTEhPhi7gNueHM/ytVv9Dk1E5OBFIjB+vCsC0aWLS5ySktw6TYsXw+OPK3ESkf0791x3/cEH\nUFDgbyyVLDDJU0yMkicJnnBMiAtPMjz8hz40yHDV+G4eNpYR38ynoDDid3giIhWXmwuvvurKjJ9w\ngis/XqsW/OUvLml67DFo1MjvKEWkKujY0VXfXLsWJk/2O5pKFZjkSamTBFnrIzP4z5/7clK3I8jL\nL+TlkT9x25PjWbYmp/ydRUT8tHatG5rXtClcdZWb31S/PtxzDyxZAo88Ag0b+h2liFQloZAbugfF\nBWZqiOAkTxq2JwGXnBTHTRd24v7BPahXO5F5y7Zw87BxvPv1PPILyq2JIiJSuWbOhKuvdovY3n+/\nS6Lat4eXXoJly1xVvfr1/Y5SRKqqoqF777/vhgPXEAFKnvyOQKRiOrduyFN/7c/J3ZuSX1DIa5/N\n5ZZhY5mzaKPfoYlITVdQAB995Baw7dQJXnkF8vLgrLPgm2+KE6rERL8jFZGqrnt3OOwwdzLmu+/8\njqbSBCZ5UqlyqUqSk+L44wUdGfK7njSqW4ula7Zyx9MTePzN6WzZmut3eCJS0yxbBv/4B7RqBWef\nDWPHQmoq3HwzzJ8PH34IffvqTKWIHDoxMXDxxe7266/7G0slCkzypGF7UhUdaxrw1F/7c9FJhthw\nDGOmLee6h7/ms0mLKdBQPhGJpm3bXAGIAQPgyCPhrrtg0SJo3hz+/W9YscJdH3WU35GKSHV1uVfl\n++23Yfduf2OpJIFJnmICE4nIgUmIC3PpwNY8/dd+dGpVn+078xj+/mxuGjaWaXPXEqlB44BFJMoK\nCuCrr+CKK1yRh6uugjFjID4eLrzQVdCbN8/1OKWl+R2tiFR37du7y6ZN8NlnfkdTKQKTsqjnSaq6\nzPopDPldT+64oisNM2qxbM1WhrzwHfc+O5nFq7L9Dk9EqrK5c+HOO10P00knuSEyO3ZAnz7w3HOw\nZg289RYMGgThsN/RikhNctll7rqGDN2L9TuAIsqdpDoIhUL07pBJt3YNGTlhMW9/aZk5fz03DxvL\ngC5HcNHJhoYZtfwOU0Sqgg0b4M034bXXYNq04sebN3c9T5ddpiF5IuK/Sy6B22+HkSNdD1RGht8R\nRVVgkicVjJDqJC42zDl9WzCg6xG8/aXl04mL+WrqMsZOX85J3Zpy/oBW1E9P8jtMEQma3Fw39O61\n19x1fr57PC3NDcu74gro3VtnHEUkOA47DE4+Gb74ws3D/NOf/I4oqgKTPGnYnlRHacnxDD77GE7r\n04w3R1vGTV/B55OX8OX3yxjYoynnDWhJ3dpKokRqtBUr3I+OUaPgyy8h2xvmGw67YXhXXAFnnglJ\n+q4QkYC67jr3PTZ8uJtzWY2LGQQmeYpR7iTVWGa9FG69pDMXDGjFW6Mt385ayciJi/liylIGdD2C\nc044isz6KX6HKSKVITcXJkxwydKoUfDjj3tvb98errzSDYVp1MifGEVEDsTpp8Phh7ulEcaMgRNP\n9DuiqAlM8qSeJ6kJmjRM5a+Xd+GCE1vxv9E/M2n2akZNXsIX3y2h5zGNObdfS1odke53mCJyqC1c\nWJwsjRnjij0USUlx5cYHDoRTToFmzfyLU0TkYMTGwuDBcN99rvdJyVP0KXeSmqRp4zTuvLIby9du\n5YOxC/gmawWTZq9m0uzVHH1UXc4+/ii6tG1EWF2yIlXT9u1uodqihGnBgr23t28Pp57qEqZevVyp\ncRGRquzaa+HBB+Gjj9zC3Ucc4XdEURGY5EkFI6QmatIwlZsu7MSlA1vzybeLGDV5CT8u3MiPCzdS\nPz2JU3seyUndmlInNcHvUEVkfyIRV068KFkaP94NzyuSnu5KjBf1LmVm+heriEg0ZGbCeee5ZRMe\nf9xdqqHAJE8hnWGXGqxu7SSuOr0dF5zYitFTlvHZpMWs3rCd1z6by/++sPTpkMmgXs1ofWS6hriK\nBEV2Nnz9dXHCtHx58bZQCLp1c8nSwIHQtasb1iIiUp3dfrtLnp57Du6+G+rW9TuiQy4w3+TKnUSg\nVmIcZ59wFGce15yZ89bz2aTFTP1pDWOnr2Ds9BUcVj+FAV2b0L9LE1XpE6lshYUwc2ZxsjRpEhQU\nFG9v2ND1Kg0c6HqZ6tXzL1YRET907Oi+A0eNgqeecnOgqpnAJE+g7EmkSExMiGNbN+DY1g1Yu2kH\noyYv4aupy1i5fhuvfTaX1z+fS8eW9RnQ9Qh6HNOYhLiw3yGLVH+NG8O6dcX3w2E4/vji3qUOHap1\neV4RkQq54w6XPP3nP27Np7Q0vyM6pAKTPMWo60mkTA0zanHlaW25bGBrZsxbz1dTlzHlxzXMmLee\nGfPWk5QQplvbxvTpmMmxpgHxSqREomPdOmjSpLjQQ//+ULu231GJiATL8cdDnz5uSYZHH4UHHvA7\nokMqMMmTqoqJ7F84HEOXNg3p0qYhW3fsZvyMlXw9dRnzl29h3IwVjJuxgqSEWLq3a0SjlJ0c075A\niZTIofTTT9C6tcrDiojsTygEDz3kEqhhw+CGG6rVmnWBSZ5EpOJSa8VzWu9mnNa7Gas3bGfCrJVM\nmLWKRSuzGTt9BQAfTfmcjq0a0K1tQ7q0aaSKfSK/Vps2fkcgIlI19O4NZ54JH3/sep6eecbviA4Z\nJU8iVVzjesmcP6AV5w9oxaoN25g4axWjJy9gzeY8Jv+wmsk/rCYUglZN0unatiFd2zbiyMZpGior\nIiIi0fPPf8LIkfDss24NqGOP9TuiQ0LJk0g1klkvhfMHtKJ5na00adaGqXPXMvWnNcxesAG7bDN2\n2Wb+O+pnaqfE075FfTq0rE+HlvVoVDfZ79BFRESkOmnbFm6+2a339Pvfw3ffuUI7VZySJ5FqqkFG\nrT1D+3bl5jNr/nq+/2ktWT+vZWP2Lr6duZJvZ64EXFGKDi3r075FPdo0y6BBei2foxcREZEq74EH\n4N13Ydo0ePJJuOUWvyP61ZQ8idQAiQmxdD+6Md2PbkwkEmHl+m3Mmr+BWfPXM3vBBtZu2sHoKUsZ\nPWUpAPVqJ9L6yAzaNMug7ZF1aZaZRjisEswiIiJyAFJS3HpPZ5/tSpj37w/t2/sd1a+i5EmkhgmF\nQhzeIJXDG6RyWu9mFBRGWLRyC7Pmb2DOoo3MXbKJDdm7mDBrFRNmrQIgMT5MyybpHHV4bVocXoeW\nTerQqG6y5k2JiIjI/p11FlxzDbz4Ilx0EUydCslVd7qAkieRGi4cE6Jlk3RaNknnvP4tKSyMsHzd\nVuYu3sTcJZuYu3gTqzdu54eFG/hh4YY9+9VKjOWow+rsSaiOzEzjsPopxKqHSkREREp64gmYOBHm\nzoWrr4a33qqyi4oreRKRvcTEhGjaKI2mjdIY2PNIADZv3cX8ZVtYsMJdFq7Ywqac3F8kVLHhEIfV\nT6FpozSOaJzqrhul0jAjWWu5iYiI1FTJyfDee9Crl5sD1by5WwuqClLyJCLlSk9NpFu7RnRrV7zI\n3aacXS6RWr6FBSuyWbomh7WbdrB0zVaWrtkKM4v3j48Lk1kvmcz6yWTWS/Fuu+s6qQmEtOioiIhI\n9daunUucBg2Chx+GtDT429/8juqAKXkSkYOSkZZIt7aN6Na2OKHamZvP8rVbWbZmK0vX5LBszVaW\nrclhQ/YulqzOYcnqnF+8TlJCLI3rJZNZL5kG6bVokJ5E/Yxae27XSoyrzLclIiIi0XLyyfDSS3DV\nVXDXXbBjBzz4IFShk6hKnkTkkElKiKXVEem0OiJ9r8e37cxj1fptrNqwndXe9aoN21i5fjvbd+ax\naGU2i1Zml/mayUlxNEhPokF6LerXSaJenSTS0xLJSEsgIy2RjLREkpPi1HslIiJSFVxxhVvv6Yor\nYOhQNw/qlVcgNdXvyCpEyZOIRF1KUlyZSVUkEiFn+25Wb9jOqg3bWb9lB+s372Tdph2s27yT9Zt3\nsH1nHot35rF41S97rYrExcaQnpZI3bRE0tMSyEhNpE5aArWTE0hLjqd2irtOS46nsDAS7bcrIiIi\n+3PppVCnDlxyCYwYATNnwgsvQL9+fkdWLiVPIuKbUChE7ZQEaqck0PrIjF9sj0QiZG/bzfotxcnU\nxuxdbMrZxeacXDbl7GRTTi47c/NdwrVpR4WOm/rxOtKSixOqtOR4UmrFk5wUS0piHMlJcdRKiiM5\nMY6UJO9+YixJCbHq4RIRETkUTjvNlS0//3yYPdutAXXBBXD//dCmjd/R7dN+kydjTAzwDNAeyAWu\ntdYuLLH9DOAeIB94yVr7Qnn7iIhUVCgUok5qAnVSE2jZJH2fz9uZm8/mnOKkamPOLrZs3UXO9t0l\nLrnkbN/N1h15ey4r1x9YPDExIZITY71kKo6kBJdQJcaHSUqIJcG7ToyPJTEhTFJ88e3EhFh3PyFM\nYnws8XFh4mNjiIsLqxKhiIjUTK1awbRproDEgw/CO++4qnynngrXXuuuExIq9FKVlbeU1/N0NhBv\nre1ljOkOPOY9hjEmDhgGdAF2ABONMR8DfYCEsvYREYmGpIRYkuqnkFk/pdznfj91GqbNMWRvy92T\nWGVv3832nXns2JXHtp15bPcuO3blF9/flUfu7oI9idehFBsOeclUmLi4GOJjw8THxez1WEJcmLjY\n4ut47zo2XHQJuevYGMIxMcTFhgjHuPtx4RjCRdtLPbfkYyIiIpUuLg7uvtsVkRg61C2m++mn7pKc\n7Hqk+vWDzp3LW1y3UvKW8pKn3sAoAGvtFGNMlxLb2gALrLXZXlATgOOBnsDn+9hHRMRX4ZjioYIH\nKr+gcE8itX1nHrtyC9i1O59duQXs3J3Prtx8du7OJ3d3ATtzix/fc7/oObkF5OUXkJtXSF5+3c9N\nOQAADu5JREFUAfkFEfIL8tlBfhTeccXdf8nhvh5fRERqsMMPh+HDYcgQeP11eO01N5zvk0/cBVwv\n1b5VSt5SXvKUBpScpV1gjImx1hZ620qWx9oK1C5nHxGRKis2HHPQide+RCIR8gsK2Z1XyO78Aned\nV8DuvALy8gvJ9a6LHttd4nZeQSH5+REKCgvJyy8kv6CQggL3eu5S4nb+3vcLCiJu/4JCCrxrERER\n3zVoALfe6i4rV8IXX8CUKZCVVd6elZK3lJc85QAl6waWfLHsUttSgS3l7LNPWeU3SI2lttk3tc2+\nqW32raq1TQhI8C7Eepekiu4Z9i7ya1W1z01VpDaOPrVxdKl9o6BDB3cpX6XkLeUlTxOBM4B3jTE9\ngNkltv0MtDTGpAPbcV1fjwCR/exTps6dO2uwvYiIBJb+nxIRCbxKyVtCkci+1zwxxoQorkABcDXQ\nGUix1j5vjDkduBeIAV601g4vax9r7byKvWcREREREZEDU1l5y36TJxEREREREXFi/A5ARERERESk\nKlDyJCIiIiIiUgFKnkRERERERCqgvGp7UWWMiaF4klYucK21dqGfMQWFtxLyS0BTXIXiv1trP/E3\nqmAxxjQAsoABKkpSzBhzJ65yTBzwlLX2VZ9DCgTv++YFoBVQCAy21lp/o/KXt5r6Q9bafsaYFsAr\nuLb5EfiDtVaTYvfDGNMZuBFXE/42a+06n0OqVowxA4ALgVrAv6y15VbBkoNjjOkPXGytHex3LNWJ\nMaYX8Dvv7s1FC7TKoVXZn1+/e57OBuKttb2AO4DHfI4nSC4F1ltrjwcGAk/5HE+geMnls7hyk+Ix\nxvQFenp/U32B5r4GFCwnA8nW2j7AA8BQn+PxlTHmNuB5vOWjgGHA37zvnBBwll+xVSEJwC3Ap7hV\n6uXQSrLW/g54FPf3K1FgjDkK6Agk+h1LNTQYlzy9iDsRIIeYH59fv5On3sAoAGvtFKCLv+EEyru4\ncorg/p3yfYwliB4BhgOr/Q4kYE4GfjDGfAh8AnzsczxBshOo7ZUlrQ3s9jkevy0AfoNLlACOtdaO\n925/DpzoS1RViLV2EtAW+Asw0+dwqh1r7UhjTDJwE65XVKLAWrvQWjvM7ziqqbC1djfut0pjv4Op\njvz4/Po6bA9Iw63sW6TAGFPuyr41gbV2O4AxJhWXSN3lb0TBYYy5CtcrN9oboqbFK4vVB5oAp+N6\nnT4GWvsaUXBMxJ2Z+hmoixvaWGNZa0cYY44s8VDJv6NtuASzxik1lLHMoeXGmAeAlrjeumnAqcB9\nwM0+hV1lVLB9HwRa4NrzIeBea+0G34Kugg6wna+31m7xMdwqqSJtDOwwxsQDmcAa/6KtmirYxpXO\n756nHCC1xH0lTiUYY5oAY4DXrLVv+R1PgFwNnGSM+QbXVfuqMaahzzEFxQZgtLU235sHtssYU8/v\noALiNmCitdZQ/LmJ9zmmICn53ZsK1LgfU2UMZSxzaLm19l5r7cVACm5u6iPAG5UfcdVyAO17j9e+\njwANgX8aY871IeQq6UDbWYnTgatoGwPP4aYYDAZer+w4q7IDaONK53fP00Tc2d93jTE9AE0G9XjJ\nwGjgBmvtN37HEyTW2hOKbnsJ1O+ttWt9DClIJuDO1g4zxmQCycBGf0MKjGSKe7o34wpqhP0LJ3Bm\nGGNOsNaOw/WkfO13QD4oGspY9COnDyWGlhtj9hpa7n036/u54g60fa+s3PCqjQNq5yLW2ssrJ7xq\noUJtbK2djjvhKwfuQL8vKu3z63fP0we4M+MTcRnkn3yOJ0j+hhs2c68x5hvvosmcsl/W2k9xP4K/\nxw3Zu0EV0/Z4BOhhjPkWlxjcaa3d6XNMQVD0+bgVGGKMmYQ7sfaefyH5w1o7gr3nl6ZSxtDyyo2q\n+lD7Vg61c/SpjaMvyG3sa8+T96Puej9jCCpr7c1o/Hy5rLX9/I4haKy1t/sdQxB5Q1PO8TuOILHW\nLgF6ebfn4yo0SjENLY8utW/lUDtHn9o4+gLTxsqKRUREyjYRGASgoeVRofatHGrn6FMbR19g2tjv\nOU8iIiJBUzSU8QNccZqJ3n3NXTg01L6VQ+0cfWrj6AtcG4ciEU2HEBERERERKY+G7YmIiIiIiFSA\nkicREREREZEKUPIkIiIiIiJSASoYIVWeMeZIYB4wp9Sm56y1wys/IscYcxUwDLfe0n3AYi+m60o8\npyMwHbjaWvvqPl7nGuB8a+3AUo+/DMzATZpsC7S01i6LwlsREREREZQ8SfWx0lrbye8gSokAH1pr\nf+sleBuBU4wxJdcmuBBYT3E1mbK8DTxmjKlvrV0PYIypBZwG/Nla+x9jzOKovQsRERERAZQ8SQ1g\njFkNvAv0wa1WfYG1dokxpiuuZ6gWsAH4vff4WFyi0w6X3BhgCLAD10sUC7wOPGit7e0d40qgu7X2\nhhKHDnmXIttwPUXHA2O9x04Cvip6njFmoHesOFxP1WBr7SZjzAdeLE95+50NfG2t3fxr20dERERE\nKkZznqS6yDTGzCh1aedtawh8Za09FhgP3GiMiQNeAC621nbGJVHPe8+PALOsta2BVcDjQH+gC5AB\nRKy1Y4BGxphm3j5XAC9XIM53gPMAvORtNrDbu18f+CdwshfraOBhb7+XgUtKvM4VwEsVbBsRERER\nOQTU8yTVxapyhu2N8q5/xPX8tAKaA58YY4qek1ri+VO86+OASdba1QDGmFeBc7xtrwKXG2NeARpa\na6dWIM6RwFBjTAjXk/Q2cJG3rTtwBDDWiymM6wED+Bao5w3/24Wb3/RlBY4nIiIiIoeIkiepEay1\nu72bEdwQuTCwqCjhMsbEAI1K7LLTu873nluk5DC8V3BJ2S5cIlWROLYZY2bhkrJ+wO0UJ08xwARr\n7VleTIl4CZ21NuIlbpd4x3u9IscTEZF9C2rBIQBjTCrwsrX2PG84+X3W2nF+xgRgjCm01u5z5JIx\n5ljgQmvt7WVsGwscBtwBHAPcC/Sy1n5X4jn/Bm4q5xjfAk9ba98q8VgysAy4HvgHsNxa2+8A355I\nuTRsT2qaouTnZyDDGNPHu/9b4I0ynj8Z6GqMaeT1Fl2EV9zBq2y3AvdFfSDJzDvAQ8BUa22B91gE\n19vV0xjT0nvsbuBfJfZ7FTgXN+yvIkMERUSkfCuttZ1KXXxNnDz3Ac96tyPsv7BQYFhrpwNNjDFH\nl7E5AlxjrX3fu78Cbyg77DmReQLlv9eX2HsoO8BvcHOB3wGuOZjYRSpCPU9SXWQaY2aUemyctfYW\n9v4SjuDmLO02xpwPPOH18GQDV5Z+UWvtemPMTcCXuB6fJRT3SoEbdneOtXZNGTGV/vIvuj8SeBG4\nq9Sx1hpjfgu8Y4wJA8uBy0psX2GMWQfEWGuXlnE8ERE5hPwqOGSMSQNOt9b+pYyY/gZcChTg5sbe\nZq0t9P6vuhHYgjtBuNBaO6TEfnG4pKNoPvAz1toXjDFNcSfk6ntxXmut/cEYMxQ33zfDe4+/sdau\nLfF6KcDT3uuFgYdL9AS9AfwFuGo/zRsBPgLO9J6L186TgA7eMcLAI7iEKgy8Yq39N+7f5FFjTHqJ\nwkmXA495t0uOEhE5pJQ8SZVnrV0CJOxne7jE7Vfxhth5wwS6l/H8Pd38xpgM3Jd4e2/o3BO4IR4Y\nY2KBE3GFJ8qy58vbi7G5d3sbkFxi29Ulbo/EJVf7ei+n7mubiIgclLJOvl1mrZ1DccGhm4wxj+IK\nDt2J+94/zTupdQqu4NBJFBccOtcrAvQZrtjQGuA9INtaO8YY87wxppm1djGuANAdpY7fH5hV6rGQ\nMWYQcAZwLC6Zex+4zhgzAbjBezwPV9F1Qan9ewHp1tpjvf/bHvXexzPAu9ba4caYU4G7vQTNWGt7\nwp75vpfiEsYidwPTrLVXesneRGPMFO89fQu8tu8m32MDsNgY08VaOw24AHdS8npv+2DcCc/OxpgE\n4AtjzDRr7QRjzEfA+cBzxphMoJW19osKHFPkV1HyJLIfXpnwOsCPxph8IAt43hvCtxIYba39cB+7\nR4AzjTGvWGuvikZ8Xq/Zd0DjaLy+iEgNEMSCQy1xQ9pK6wf8z1qb673mS7hREwnAJ97JOYwxbwLp\npfb9wW0yo3BJXVHCdjyulwxr7efA595r3GqM+R2u96wnv0zGTgSSvBET4Hrh2gKLrbU5xpiQMSbD\nWrupjPdR0jvAecaY6bgE74+ljtHBGNPfu58MHA1MwPWi/R14DpfYVSRZE/nVlDyJlMMb+leWhuXs\nt6eXK1qstbuAjtE8hohITeZTwaECb//SYkq9Tgzut1zBfo5V9D42eUt4nAQMAqZ79/NKPt8Y0xZI\nAv6HGwb3rhdL6deMAS611s709mtEcYVYvNctZP8iwAfARNwQxHHeKI+Sx/hr0UlKrzdvq/d+Jnjz\nkQ/HJU/nlH5xkWhQwQgRERGR8lVmwaGFQNMyHh8DXGyMSfSGjl/tPfY1MMgYk2qMiccVF9orcTHG\nnA7811r7KXAzbuH2Jrj1Dy/ynnMSrkjF8cBYa+1zwFzgZPZOzopiucHbrzFuEfjDvfupQMhau6WM\n91Ak5D1nE7AUeBA3ZK/0MX5njIn15lh9C3Qrsf1V4B5gozdcUCTq1PMkIiIiNVkQCw59zd7zi/CO\n/akxpiMwDfcbbhTwpFcw4j+4hG0bbi7RzlL7j8INj5vjxfO+tfZHY8yNwAvGmBuA7cC1QA4wwmuX\nDbihfEWLwhe1yRDgGWPMD7jE6rYSCcwJwCdlvK+93k+J13oHuLdEyfKix/8PN4Rxhvd+X7TWji/x\nGq8Bi3FJZOnXFomKUCSiz5eIiIjIoeQVZbgJGFKy4JC19mmv1+h14O19zZs1xjwGjPF6iso7Vktc\nAYt/e/c/BJ6vyL7RYIx5D7cu1ZxSj38D3B/t9aqMMX2942udJznkNGxPRERE5BDzhqMVFRyahSsq\nUbLgUP5+Cg6B69mp6HpFS3FDBH8wxszGJWl+JU5dgCWlE6cSXjDG/CaKx78QV/1QvQMSFep5EhER\nERERqQD1PImIiIiIiFSAkicREREREZEKUPIkIiIiIiJSAUqeREREREREKkDJk4iIiIiISAUoeRIR\nEREREamA/wf1Jik6NHrkbAAAAABJRU5ErkJggg==\n",
       "text": [
        "<matplotlib.figure.Figure at 0x7fad3c339fd0>"
       ]
      }
     ],
     "prompt_number": 117
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "The mean $\\mu$ for a probability distribution X with probability density function p(x) is defined as:"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "%%latex\n",
      "\\begin{aligned}\n",
      "\\mu_X & := < x, p(x) > \\\\\n",
      "      & = \\int_X x \\; d\\{p(x)\\} \\\\\n",
      "      & = \\int_X x \\cdot p(x) \\; dx\n",
      "\\end{aligned}"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "latex": [
        "\\begin{aligned}\n",
        "\\mu_X & := < x, p(x) > \\\\\n",
        "      & = \\int_X x \\; d\\{p(x)\\} \\\\\n",
        "      & = \\int_X x \\cdot p(x) \\; dx\n",
        "\\end{aligned}"
       ],
       "metadata": {},
       "output_type": "display_data",
       "text": [
        "<IPython.core.display.Latex at 0x7fad3cff3bd0>"
       ]
      }
     ],
     "prompt_number": 113
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "For the Watt distribution, the probability distribution function $p(E) = 0.4865 \\sinh(\\sqrt{2E})e^{-E}$ [link](http://indico.cern.ch/event/145296/contribution/6/material/slides/8.pdf) ($\\sinh$ is the hyperbolic sine function). Technically $E$ ranges from 0 MeV to $\\infty$ MeV, but most of the time people cut it off around 10 MeV. So, using the above definition of the mean:"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "%%latex\n",
      "\\begin{aligned}\n",
      "\\mu_W & := \\int_{0}^{10}  E\\cdot 0.4865 \\cdot \\sinh(\\sqrt{2E})e^{-E}\n",
      "\\end{aligned}"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "latex": [
        "\\begin{aligned}\n",
        "\\mu_W & := \\int_{0}^{10}  E\\cdot 0.4865 \\cdot \\sinh(\\sqrt{2E})e^{-E}\n",
        "\\end{aligned}"
       ],
       "metadata": {},
       "output_type": "display_data",
       "text": [
        "<IPython.core.display.Latex at 0x7fad3c27ce90>"
       ]
      }
     ],
     "prompt_number": 129
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "It's possible to do this integral by hand, but it sucks (involves the error function). Just approximating it using numerical integration..."
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "from scipy.integrate import quadrature\n",
      "from scipy.optimize import minimize"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 123
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "mean, error = quadrature(lambda e: e * watt(e), 0.0, 10.0)"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 121
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "print(mean)"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "1.99677578646\n"
       ]
      }
     ],
     "prompt_number": 122
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "You asked about the *most probable* energy. Note in this case that this is **not the same as the mean**, as is evidenced above (this is only the case for symmetric single-mode distributions). The most probable energy is the the maximum of $W(E)$. \n",
      "(Note that we have only a *minimize* function, but finding the maximum of a function is the same as finding the minimum of the negative of a function, at least in this case)"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "res = minimize(lambda e: - watt(e), 1.0)"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 125
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "print(res.x)[0]"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "0.719636649623\n"
       ]
      }
     ],
     "prompt_number": 127
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [],
     "language": "python",
     "metadata": {},
     "outputs": []
    }
   ],
   "metadata": {}
  }
 ]
}