{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "from numpy import log10, sin, sqrt\n",
    "from scipy.integrate import quad\n",
    "%matplotlib inline\n",
    "plt.rcParams['figure.dpi']=150"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {},
   "outputs": [],
   "source": [
    "import ROOT\n",
    "def setup_minuit(fn, params):\n",
    "    from ctypes import c_int\n",
    "    PRINT_LEVEL = 0\n",
    "    UP = 1.0\n",
    "    npar = len(params)\n",
    "    minuit = ROOT.TMinuit(npar)\n",
    "    minuit.SetFCN(fn)\n",
    "    minuit.SetErrorDef(UP)\n",
    "    minuit.SetPrintLevel(PRINT_LEVEL)\n",
    "    \n",
    "    status = c_int()\n",
    "    print(\"{:<10s}{:>10s}{:>10s}{:>10s}{:>10s}\".format('param', 'guess', 'step', 'min', 'max'))\n",
    "    for ii, t in enumerate(params):\n",
    "        name, guess, step, pmin, pmax = t\n",
    "        print(f\"{name:<10s}{guess:10.3e}{step:10.3e}{pmin:10.3e}{pmax:10.3e}\")\n",
    "        minuit.mnparm(ii, name, guess, step, pmin, pmax, status)\n",
    "        if status.value != 0:\n",
    "            raise ValueError(f\"** mnparm({name}) status = {status.value}\")\n",
    "    return minuit\n",
    "\n",
    "def perform_fit(minuit, params):\n",
    "    from array import array\n",
    "    from ctypes import c_int, c_double\n",
    "    \n",
    "    LINE = '='*79\n",
    "    print(LINE)\n",
    "    print('Running Minuit')\n",
    "    \n",
    "    # define fit control parameters\n",
    "    MAXITER = 5000\n",
    "    TOLERANCE = 1e-5\n",
    "    args = array('d')\n",
    "    args.append(MAXITER)\n",
    "    args.append(TOLERANCE)\n",
    "    \n",
    "    swatch = ROOT.TStopwatch()\n",
    "    swatch.Start()\n",
    "    \n",
    "    # run MIGRAD minimizer\n",
    "    status = c_int()\n",
    "    minuit.mnexcm('MIGRAD', args, 2, status)\n",
    "    if status.value != 0:\n",
    "        raise ValueError(f'** mnexcm status = {status.value}')\n",
    "    print(f'real time: {swatch.RealTime():10.3f} s')\n",
    "    \n",
    "    # print results\n",
    "    out = open('fit.txt', 'w')\n",
    "    value = c_double()\n",
    "    error = c_double()\n",
    "    results = []\n",
    "    print('{:<10s}{:>11s}{:>13s}'.format('name', 'value', 'uncertainty'))\n",
    "    for ii, t in enumerate(params):\n",
    "        name = t[0]\n",
    "        minuit.GetParameter(ii, value, error)\n",
    "        record = f'{name:10s}{value.value:11.3f}{error.value:13.3f}'\n",
    "        out.write(record + '\\n')\n",
    "        print(record)\n",
    "        results.append((value.value, error.value))\n",
    "    out.close()\n",
    "    print(LINE)\n",
    "    return results"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {},
   "outputs": [],
   "source": [
    "def parse(fname):\n",
    "    with open(fname) as f:\n",
    "        rows = [r.split() for r in f.readlines()[1:]]\n",
    "    bins = [tuple(map(float, r[1:])) for r in rows]\n",
    "    zs, xs, sigmas = map(np.array, zip(*bins))\n",
    "    return zs, xs, sigmas\n",
    "    \n",
    "zs, xs, sigmas = parse('SCPUnion2.txt')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {},
   "outputs": [],
   "source": [
    "def mu(z, O_p, O_m):\n",
    "    def integrand(a):\n",
    "        O_a = 0.5*(O_p+O_m)/a**3 + (1-O_p)/a**2 + 0.5*(O_p-O_m)\n",
    "        return 1/(a**2 * sqrt(O_a))\n",
    "    u = quad(integrand, 1/(1+z), 1)[0]\n",
    "    \n",
    "    O_k = sqrt(O_p - 1)\n",
    "    return 5*log10((1+z)*sin(O_k*u)/O_k)\n",
    "\n",
    "def chi2(O_p, O_m):\n",
    "    return sum((x-mu(z, O_p, O_m))/sigma for z, x, sigma in zip(zs, xs, sigmas))\n",
    "\n",
    "def chi2_fit(npar, grad, fval, xval, flag):\n",
    "    O_p = xval[0]\n",
    "    O_m = xval[1]\n",
    "    fval[0] = chi2(O_p, O_m)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "/usr/lib/python3.6/site-packages/ipykernel/__main__.py:8: RuntimeWarning: invalid value encountered in double_scalars\n"
     ]
    }
   ],
   "source": [
    "N = 10\n",
    "\n",
    "delta = 0.2\n",
    "Ops = np.linspace(1, 1+2*delta, N)\n",
    "Oms = np.linspace(-0.4, -0.4+2*delta, N)\n",
    "\n",
    "vals = np.zeros(shape=(N,N))\n",
    "\n",
    "for i, O_p in enumerate(Ops):\n",
    "    for j, O_m in enumerate(Oms):\n",
    "        vals[i][j] = chi2(O_p, O_m)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<matplotlib.colorbar.Colorbar at 0x7fe07a4b6208>"
      ]
     },
     "execution_count": 20,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAA1IAAAItCAYAAAAzJUunAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAAXEQAAFxEByibzPwAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMS4wLCBo\ndHRwOi8vbWF0cGxvdGxpYi5vcmcvpW3flQAAIABJREFUeJzs3Xm4ZGV57/3vDwQRUMABYiRIDqKE\nMSiRIaAYUYOAhsGQowkq8ZUTPSpg9GRAPQZioggIjlGP4IkkRgEHBD0KCgZlUhARBIUwiCKDyDx0\n0/t+/1irQqWo2l21d9Ueen8/Xutau9Zaz1AN0n33/az7SVUhSZIkSRreavM9AUmSJElabAykJEmS\nJGlEBlKSJEmSNCIDKUmSJEkakYGUJEmSJI3IQEqSJEmSRmQgJUmSJEkjMpCSJEmSpBEZSEmSJEnS\niAykJEmSJGlEBlKSJEmSNCIDKUmSJEkakYGUJEmSJI3IQEqSJEmSRmQg1SPJWkneneQnSR5M8osk\nn0qy8Qz6Wj/JB5LckOSh9nx8kvWnabNakkOTXJ7kgSS3Jfl8ki1n980kSZIkjUuqar7nsGAkWQs4\nG9gFuBn4d2BT4LnAbcDOVXXtkH09CTgf2Bz4D+B7wFbtcQ2wU1X9qqdNgM8BBwB3tnN5MvA84EHg\nBVV14ay+pCRJkqRZMyP1X/0NTRB1PvDMqjqwqnYE3go8BfjUCH0dRxNEnQY8q+1ra+CDwDOAY/u0\neS1NEPVTYIuqOqCqdgdeATwOODnJY2b0zSRJkiSNjRmpVpI1gFuB9YFnV9WlPfcvA7YFdqiq76+k\nr98Afg6sAH6rqm7puvdY4GfAE4Gn9dy7AtgS2LeqvtjT55eAlwEHVNWpM/6ikiRJkmbNjNQjdqUJ\noq7tDaJap7TnfYboa0+aX9tvdwdKAFX1EHA6sHr7HABJfpsmiHoAOGOW40uSJEmaIAOpR2zXni8Z\ncP+SnufG3Vfn5x9V1fJZji9JkiRpgnzf5hGbtOebBty/qee5cfc1zvE7ywT7eSZN1utnw/QjSZI0\nj34LuL+qfmO+J5Lky8BmExzi2qp62QT715gZSD1i3fZ8/4D79/U8N+6+xjn+dFZbjdUev/bq6y25\ncur12DXmewrzYsWame8pzIuppfmPm1pjab73uvpjpuZ7CvNi7ccsm+8pzIt1V39ovqcwL9bN0vv/\n97XXL+ehZQvme2+25prZ8hmbjv83mGuuX86yhfM9NSQDqUd0/rQ56N/iUf40OpO+VtZmJFW1Vd+B\nkyvWXn29LXd98p+MY5hFZdkWI28Ftkq4a7O15nsK8+LeTZZmAPnAxv1WBq/6NnranfM9hXmx80bX\nzfcU5sVL1r98vqcwL/7wcUsvcN7m+Tdw5U+WLZhVNM/YdA0uP/fpY++3/Z5j71eT5TtSj7inPa8z\n4P7a7fneCfW1sjad68OML0mSJGmCzEg94sb2PChtsXHPc+Pua5zjS5IkaeyKKSaxlNhlfYuRGalH\nXNaenz3gfuf6DyfUV6fN1u2eVrMZX5IkSdIEGUg94jvAXcBmSbbvc/+A9vyVIfr6GjAF7JZkw+4b\n7Ya8+7T3v9q5XlXXAT8GHgfsNcvxJUmSNGYFrKipsR/moxYnA6lWVS0DPtR+/FCS/3xXKcnhwLbA\neVV1cdf1/5nkqiT/0NPXzcC/AmsCH0nSvYTyfcBTgH+pql/2TOPYzjPdAViS/YCXAdcBX5zF15Qk\nSdIsTFFjP7Q4+Y7Uf3UUsAewC/DTJP8OPB3YEfgV8Nqe558MPAt4ap++DgV2AvYHrkryPWArYGvg\nWuCwPm0+BbwU2Ldtc3Y7xvOBB4E/HbBZryRJkqQ5ZEaqS1U9CLwAOJJmP6c/AjYFPg1sX1XXjNDX\n7cDvAR+kyUztC6xHk/V6bnu/t80U8ArgrcAvgL2BbYAvADtU1Xdn+t0kSZI0OwVMTeB/5qQWJzNS\nParqAeCd7bGyZ/838L+nuf9r4M3tMez4K2iW+B27smclSZIkzQ8DKUmSJGkoxYqaRP7InNRi5NI+\nSZIkSRqRGSlJkiRpSL7RpA4DKUmSJGkIBayYQCBlaLY4ubRPkiRJkkZkRkqSJEkakkv71GFGSpIk\nSZJGZEZKkiRJGkLBRMqfm+NanMxISZIkSdKIzEhJkiRJQ5qa7wlowTCQkiRJkoZg+XN1c2mfJEmS\nJI3IjJQkSZI0pBWmj9QyIyVJkiRJIzIjJUmSJA3JYhPqMCMlSZIkSSMyIyVJkiQNoanal4n0q8XH\njJQkSZK0CkjynCR/leS0JD9PUkkenOb5lyX5dJLLk9yeZHmSW5OcmWSvIcY7IMnX27YPJrmxHXvX\nAc+vn+QDSW5I8lB7Pj7J+tOMsVqSQ9s5PpDktiSfT7LlSua2d5Jzk9yV5O72571X9p1GYUZKkiRJ\nGkbB1CTSR+Pr8x3Ay0d4/iBgP+AK4ELgHmBTYE9gzyRHVtU7exslWR34v8ArgfuA84A7gU2AvYBL\n2mvdbZ4EnA9sDvwH8EVgK+DNwEuT7FRVv+ppE+DfgAPa/s8AngzsD+yV5AVVdWGf+b0ZOB54GDgL\neAh4MXB6krdU1Qkj/BoNZCAlSZIkDWERLO07H7gMuLg9frmS5/8eOKRPALMjTQByRJLPVtWVPe3e\nQxNEnQn8WVXd0dV2A5pgp9dxNEHUacCBVfVw+/wJwJuAY4FX97R5LU0Q9VNgt6q6pW2zP3AKcHKS\nLTp9tfeeCRxDEzy9oKrO77r+XeCYJF+tqp+u5NdmpVzaJ0mSJK0Cquq9VfWuqvpKJ+hYyfOX9gZR\n7fULgc8CAXbvvpdkc+CtwI3AK7qDqLbtr3uDlCS/AbwKWA68oTvwAd4G3Aa8KslGPVN5a3t+e/f3\nqapTgS8Dm/HoDNxbaJJFH+sEUW2bn9AEjo+hyYLNmoGUJEmSNKQVZOzHArWiPS/ruf46YHWaQOX+\nIfvakybu+HZvgFdVDwGnt33u2bme5LeBLYEHaJb09TqlPe/Tc33vnvvdPj+gzYy4tE+SJEnSf0qy\nLXAgTQbp7J7bL2zP32iDnf8OPB24Azi7qs7q0+V27fmSAUNeAhzc9Vx3mx9V1fIBbbqfoy1asUn7\n8dLeBlV1U5LbgacnWa+q7hown6EYSEmSJElDKGCqFvQ7UjOSZB+aAg5r0AQiu9AEUa+vqut6Ht+q\nPe9I8y7SY7vu/VWSs4D9q+ruruud4OamAVO4qee52bb5dVXdN027J7fPXj7gmaEYSEmSJEnzb7Mk\nV/S7UVVb9bs+RtvxXws9PAgcSlOZ7z8lWQtYq/34AeCbwF8C19MEVh8H9mjPf9LVdN32PGgp4H09\nz02qzaB2M+I7UpIkSdKQVsV3pKrqqKoK8DhgG+BTwEeBLyVZs+vR1bt+/jmwT1VdXlX3tEv6Xg5M\nAX/cFqXo6HzJQcm3fr8IK2szXT/TtRnbL7gZKUmSJGkIRVgxgTxENX+2v3YOMk/Tz6PqQeBHwBuT\nPExT3e5NNEv4qKr7kkzRJGM+U1XLetpfnuR7wHOB59OULYdmfyqAdQYMvXZ7vrfr2sradK6P0mbQ\nWDNiRkqSJElSr8+0597y4jf0nHtd35437Lp2Y3veeECbjXuem22bDZIMCqb6tZsRAylJkiRpSFOV\nsR8L1O3t+Sk91zvV8J44oN2T2nN3xuey9vzsAW0613/Yp83WSdYYpk1V3ckjAdL2vQ2SbExTaOLG\n2VbsAwMpSZIkSY/2/PZ8bc/1L7fnF/Q2SPJ4Hglwukudf43m3andkmzY0+axNPs6TQFf7VxvqwX+\nmOa9rb36zO+A9vyVnutn9Nzv9ooBbWbEQEqSJEkaQjGZYhPzUf48yYZJ3truvdR770XA+9qPJ/bc\n/izN8r2XJHl1V5vHAMcDG9C8Z/Wdzr2quhn4V2BN4CPtsx3vo8l6/UtV/bJnrGM7z3QHYEn2A14G\nXAd8safN8TSbCf+PJDt1tdkc+Nv23gm933kmLDYhSZIkrQKS7AW8o+fymkku6Pp8ZFWdQVN04f3A\nkW2BiJtoijQ8E9iiffa4qjq1u7OqeijJfwfOAk5K8maawOrZwKbAr4BXVlVvfHgosBPNflVXtWNu\nBWxNk/U6rM9X+hTwUmDfts3ZNEvznk9Tov1Pezfrraqrk7yNJgj79yTfAJYBL6bJbh1eVVf3GWtk\nBlKSJEnSkFbUgl7Q9RSa/Zy6peda552nW4G3A7vTBDQ70KxWu5km6/RPVXVOv0Gq6oIk2wPvotk3\nahvgFuATwFFV9ahCDlV1e5LfA94N/BFNcHQL8CHgXVV1R582U0leAbwFOBjYm2YfqC8A76yqQftu\nHZfkGuBtwG7t5e8DR1fVl/u1mQkDKUmSJGkoYWoib8aMp+BEVZ0EnDTks/cDR7fHTMb6KfCnI7b5\nNU1J9TeP0GYFTXbp2JU929PudOD0UdqMakGH1JIkSZK0EJmRkiRJkobQKTYxiX61+JiRkiRJkqQR\nmZGSJEmShrTAi01oDvlvgiRJkiSNyIyUJEmSNIQCpnxHSi0zUl2S7JLkzCR3JLk3yUXdOzbPoL+9\nk5yb5K4kd7c/7z3g2U2T1DRH707PkiRJkuaJGalWkn2Bz9MEl98GbgdeSLNj83ZVdfiI/b0ZOB54\nmGbn54dodlQ+PclbquqEAU1vAb7W5/pdo4y/UiumxtrdYpCppfn3PVl6/6gbS/R7Z2r8f1O6GKyo\npfm9l9fq8z2FebFU31GZWpL/YVtov3eHFQt4HynNLQMpIMkGwInA6sD+VXVae30j4DzgsCSnV9W3\nhuzvmcAxNMHTC6rq/K7r3wWOSfLVdiOzXldV1Wtm+50kSZI0XsVkAvmFFi5qOEvzr3Qe7XXAesCX\nOkEUQFXdAry9/ThKRuotNEHqxzpBVNvfT4C/b+8NvaOzJEmSpIXFQKrReW/plD73zgAeBPZIstYY\n+vt8e95n+OlJkiRpIZhitbEfWpxc2tfYtj1f0nujqpYl+RGwA/As4LLpOkqyPrBJ+/HSPv3dlOR2\n4OlJ1quq3nefNkrybuCpNO9FXQh8uaqWjfKFJEmSJE3Okg+kkjwBWL/9eNOAx26iCaQ2YSWBFI8E\nUb+uqvum6e/J7bOX99zbAnhnz7Ubk/xxVV24krH/U5IrBtzabNg+JEmS9IiqTKS4TS3RgjmLnblE\nWLfr5/sHPNMJiNYdcL9ff4P6GtTfQ8BHgd2BjWje2doZOJMm4Ppakk2HGF+SJEnShK0SGakkpwBb\nj9jsoKq6iOHqTY7y1wSdZ6crwPKo/qrqZuANPZcvAPZKcjLwSuBvgNcPM4mq2qrvwE2masth+pAk\nSdJ/NZny51qMVolACtiU5v2lUazdnu/puXb3NM/eO0S/nf7WGWLsYfoDeA9NIPWSIZ+XJEnSBEwt\n0X3M9GirRCBVVTvMou3dSe6iWUq3MXBln8c2bs83DtFl55kNkqwz4D2pUfoD6Ow39dQhn5ckSZI0\nQYbUjU4BiWf33kiyBs2ywYeAq1fWUVXdySMB0vZ9+tuYptDEjX0q9g2yQXseNoMlSZKkMSvCClYb\n+1EjvUWihcJAqnFGez6gz729gbWAs6vqwTH094r2/JXhp8f+7fn7I7SRJEmSNCEGUo1P0rwb9fIk\n+3UuJtkQeF/78djeRkmuao+n9dw6HlgB/I8kO3U9vznwt+29E3r6OqjNVvWOsR/wj+3Hj4z6xSRJ\nkjQ+K9oS6OM8tDitEu9IzVZV3ZHkYOBzwClJzgVuB/ag2WPqhKo6u0/TToGLNXr6uzrJ22iCr39P\n8g1gGfBi4HHA4VXVu0zwYODEJFcB1wMP0lTX26K9f3RVfWF231SSJEnSOBhItarq1CTPA44AdgLW\nBH4MfLiqTpxBf8cluQZ4G7Bbe/n7NAHRl/s0+QRwG/C7wK40AddtwGnAR6vqrFHnIEmSpPEpYGoC\nC7qm2zNHC5eBVJeq+g6w5wjPT5uLrarTgdOH7Otk4ORhx5YkSdLcW2H5c7X8N0GSJEmSRmRGSpIk\nSRpCEaYmUKrc8ueLkxkpSZIkSRqRGSlJkiRpSL4jpQ7/TZAkSZKkEZmRkiRJkoZQwArLn6tlRkqS\nJEmSRmRGSpIkSRpGwdT024jOuF8tPgZSkiRJ0lAykaV9WP58UXJpnyRJkiSNyIyUJEmSNIQCpiZQ\n/tyVfYuTGSlJkiRJGpEZKUmSJGlIK3yfSS0zUpIkSZI0IjNSkiRJ0hCKTOgdKbNci5GBlCRJkjQk\nl/apw6V9kiRJkjQiM1KSJEnSkCaxtE+Lk/8mSJIkSdKIzEhJkiRJQyhghRvyqmVGSpIkSZJGZEZK\nkiRJGkqYmkjVPisBLkYGUpIkSdIQXNqnbi7tkyRJkqQRmZGSJEmShlEwVRNYhmdKalEyIyVJkiRJ\nIzIjJUmSJA2hgBUTyEOYkFqczEhJkiRJ0ojMSEmSJElDyWTekbL8+aJkRkqSJEmSRmRGSpIkSRrS\nlHkItQykJEmSpCE0G/KOfxmexSYWJ0NqSZIkSRqRGSlJkiRpSJMpNqHFyIyUJEmSJI3IjNQSVFXU\nihXzPY05lxVLcwVypuZ7BvNjqX5vppbm35SumFqafy/48NTq8z2FebG8lur3Xnq/dy+837nDVE3i\nvzdL87/di93S/J1HkiRJkmbBjJQkSZI0hAJWTCB7tPAybxqGgZQkSZI0JItNqMOlfV2S7JLkzCR3\nJLk3yUVJXj2Dfp6c5HVJPp7kB0keTlJJ/mSItr+T5PNJbkvyQJLLkxyWxH9WkiRJ0gJhRqqVZF/g\n8zTB5beB24EXAicl2a6qDh+hu12BT8xgDjsBZwNrAxcB1wPPA44Ffj/JK6rK7K8kSdI8qGIixSb8\n093iZJYDSLIBcCKwOnBAVe1eVQcAWwDXAIclecEIXd4CfAR4LbA18M9DzOExwGdogqjDq2rHqjoQ\n2Bw4H9gfeM0Ic5AkSZI0IQZSjdcB6wFfqqrTOher6hbg7e3HoTNSVXV+Vb2xqk6qqiuAYQox7wts\nBlxWVcd19XUv8MZR5yBJkqTxmyJjP7Q4GUg19m7Pp/S5dwbwILBHkrXmYw5VdSnwH8DWSTad4Bwk\nSZIkDcFAqrFte76k90ZVLQN+BKwFPGuCc9hu0Bx6rm834L4kSZImKqyo8R9uyLs4LfliE0meAKzf\nfrxpwGM3ATsAmwCXTWgqmwwxh+7nJEmSNIeKCRWbGHuPmgtLPpAC1u36+f4Bz9zX59lJzWMsc0hy\nxYBbm40yKUmSJEmPtkoEUklOoamON4qDquoihsulzmW+ddBfSpjzlSRJmmduyKuOVSKQAjZl9PeX\n1m7P9/Rcu3uaZ+8dcYxR3AtsAKwz4P5Ic6iqrfpdbzNVW448O0mSJEn/aZUIpKpqh1m0vTvJXTTl\nzzcGruzz2Mbt+caZjjOEG2kCqY2BH87THCRJkjQNy5Wrw6p9jU4BiWf33kiyBs2ywYeAq+djDj3X\n+wVZkiRJkuaQgVTjjPZ8QJ97e9OUPj+7qh6cjzkk2R74b8CVVXXdBOcgSZKkAZqqfRn7YdW+xclA\nqvFJmnejXp5kv87FJBsC72s/HtvbKMlV7fG0MczhC8B1wHZJDusaYx3gw4PmIEmSJGnurRLvSM1W\nVd2R5GDgc8ApSc4Fbgf2oNlj6oSqOrtP006BizV6byS5oOtjp+T4kUkObX++pKre0DWH5Un+FDgL\nODbJgcANwG7AU4EvAifO9DtKkiRptjKRfaQszrw4GUi1qurUJM8DjgB2AtYEfgx8uKpmEsDs2Ofa\nM9oD4FHLBKvqu0l+D3g3sDvwu8C1wDHAB6pqagbzkCRJ0jjUhMqfu7ZvUTKQ6lJV3wH2HOH5gf9P\nmu7eSvq8gv7vakmSJElaIHxHSpIkSRpC0ZQ/H/cxroRUkuck+askpyX5eZJKMrBYWpKXJfl0ksuT\n3J5keZJbk5yZZK8BbXZv+x10XNCvXdv2VUm+k+SeJPcmuTjJ66Z5/pyVjPWH07Q9KMlF7Th3tN9p\nl0HPz4QZKUmSJGnV8A7g5SM8fxCwH3AFcCFwD7ApzQqtPZMcWVXvHND2WuC8AdcfJclHgL+g2VLo\nfOB+4PeBTyTZuar+fJp5ngrc2+f6zweMdSxwGPAA8HWaCtwvAl6c5BVV9YVpxhqagZQkSZI0pIm8\nIzU+59PsTXpxe/xyJc//PXBIVf2q+2KSHWkKoB2R5LNVdWWftudV1WuGmVSS/WmCqF8De1TVJe31\npwHfAA5O8vWq+rcBXfxlVV0/5Fh/QBNE/QrYuap+2l7fGTgHODHJOVX162H6m45L+yRJkqRVQFW9\nt6reVVVfqapbhnj+0t4gqr1+IfBZmnKCu49han/Rnt/fCaLacX4OvK39+L/GMA7AW9vzUZ0gqh3r\nfOBjwHrAweMYyEBKkiRJGtIkNuRdoFa052Vj6Os57fmcPvc617ZP8luzGSTJWsAL24+n9Hmkc22f\n2YzT4dI+SZIkaQjFZAKfWmD7SCXZFjgQWA7020sVYPMk/wA8iWb/1fOArw3Yrmed9vyo5XRVdV+S\nh4DHAtsBP+vT/s+TPAmYAn4CfLGqbuzz3BZtP7dV1U197neyYdsO+E4jMZCSJEmS5t9mSa7od6Oq\ntprkwEn2AfYH1gA2AXahCaJeX1XXDWi2S3t0uzzJ/t1L6lq3Ab8JPJ1mn9busX+DJviBptBFP0f0\nfH5/WwjjyJ7rm7TnfkFUJ2i7E9ggyeOr6p4B4w3FpX2SJEnSkFbRpX3bAa8GXgnsSrOc7y3A/+3z\n7F3A0cBONNmoJ9Esp7sA2Ab4RpL1etqc255f06e/13b9/Piee98G/gzYDFgbeBbwt8DDwN8leUvP\n8+u25/v7jNNxX8+zM2ZGSpIkSZp/10468zRIVR0FHNW+Y/QMmuIQHwX2bjNMy7qevRS4tKeLbybZ\nFfgWsBvwRuA9XfePBg4ADkxyA/Ah4MH22jtpAqPH0Czd655Xb+n1nwDvSfI94P8B707y8ap6oL3f\niUqn25prbJGrGSlJkiRpSJPYkHehqKoHq+pHVfVG4IPA3sCbhmy7Anhv+/ElPfcupcl4PQC8HbgR\nuBX4CPBd4Mvto0OVJK+qrwPfo6nAt1PXrc5SvXUe1egRa7fnfvtSjcRASpIkSVKvz7TnUTb47bwb\n9dTeG1X1rzRL9A6jKUN+PM1mwC8CNm4f6/uO2AhjdQpQbEwfSdYB1gfunO37UeDSPkmSJGkoxWQ2\n5J1uHdo8ur09P2WENhu0577Znqq6GfhA97UkTwS2p8kmXdKv3QhjXQ08BDwlycZ9Kvc9uz3/cIRx\nBjIjJUmSJA2jJlRsYmFGUs9vz9eO0Gb/9vz9Edq8iaZa4D93ves0rSRPoXkXC7qCr7b9N9uPB/Rp\n2rn2lRHmN5CBlCRJkrTEJNkwyVuTrN/n3ouA97UfT+y5d0i7p1P3tSQ5hGbZXtEs3evtc4c+1/6U\npgrf7cC7eu7tlOQFSdJzfVPgCzTvQX25T9bp2PZ8RJLNu9rtDBwC3A38n965zIRL+yRJkqQhLZBy\n5X0l2Qt4R8/lNZNc0PX5yKo6g6bowvuBI9sqeDfRBCfPpNnYFuC4qjq1p7+/Bj6Y5ErghvbaNsBv\n01Tde0tV9ctIXZzkGpp9pO4HfpemnPmvgD2r6vae57egCeJuTvIT4Jc07z49B1iL5n2q/693kKo6\nK8nxNOXbf5DkG8CaNO9irQa8qqru6DO/kRlISZIkSauGpwA79lxLz7XOO0+30lTQ2x3YCtiBJtC4\nGfgs8E9VdU6fMY4BXty2eSHNsrybaYpTnFBVFw+Y23HtWLsBj6MpDHEs8N6qurXP8xfSlGDfEdgS\n+H2aPaB+AHwe+OigpYBVdWiSHwD/kyaAWg6cDRxVVecNmN/IDKQkSZKkIRST2UC3xlQCvapOAk4a\n8tn7afZ3OnrEMT5IUxp91LkdPuLzPwbeMOo4Xe1PYshfi5nyHSlJkiRJGpEZKUmSJGlItYDfkdLc\nMiMlSZIkSSMyIyVJkiQNaWpM7zNp8TOQkiRJkoZQTKb8+cLcj1cr49I+SZIkSRqRGSlJkiRpSBab\nUIcZKUmSJEkakRkpSZIkaRg1mXekfElqcTIjJUmSJEkjMiMlSZIkDSUTekfK964WIwMpSZIkaQiW\nP1c3A6ml6uGH53sGcy5TS/M/U5ma7xnMj6X6vVmi33tqamn+be7yWpor9JfX0vzjy/JaMd9TmHNl\niKEFbGn+l0iSJEmagTK2U2tp/lWWJEmSJM2CGSlJkiRpSFMWhlDLjJQkSZIkjciMlCRJkjSkyZQ/\n12JkICVJkiQNwfLn6ubSPkmSJEkakRkpSZIkaRg1ofLnpqQWJTNSkiRJkjQiM1KSJEnSkCw2oQ4z\nUl2S7JLkzCR3JLk3yUVJXj2Dfp6c5HVJPp7kB0keTlJJ/mSaNpu2zww6fjm7bydJkiRpXMxItZLs\nC3yeJrj8NnA78ELgpCTbVdXhI3S3K/CJGU7lFuBrfa7fNcP+JEmSNBaZUEbKLNdiZCAFJNkAOBFY\nHdi/qk5rr28EnAccluT0qvrWkF3eAnwEuLg9/hfwZ0O2vaqqXjPC9CVJkiTNMQOpxuuA9YAvdYIo\ngKq6JcnbgdOAw4GhAqmqOh84v/M5ydR4pytJkqT5MIl9pLQ4GUg19m7Pp/S5dwbwILBHkrWq6sG5\nm5YkSZIWimIy5c+tfr44GUg1tm3Pl/TeqKplSX4E7AA8C7hswnPZKMm7gafSvBd1IfDlqlo24XEl\nSZIkDWnJB1JJngCs3368acBjN9EEUpsw+UBqC+CdPdduTPLHVXXhsJ0kuWLArc1mPDNJkqQlzvLn\n6rD8Oazb9fP9A565r8+z4/YQ8FFgd2Ajmne2dgbOpAngvpZk0wmOL0mSJGlIq0RGKskpwNYjNjuo\nqi5iuHqTE/+rh6q6GXhDz+ULgL2SnAy8Evgb4PVD9rdVv+ttpmrLWUxVkiRpaaoJZaR8SWpRWiUC\nKWBTmveXRrF2e76n59rd0zx774hjjMt7aAKpl8zT+JIkSZK6rBKBVFXtMIu2dye5i2Yp3cbAlX0e\n27g93zjTcWbpp+35qfM0viRJkjB5pEf4jlSjU0Di2b03kqxBs2zwIeDquZxUlw3a83xlxCRJkkSz\ntG/chxYnA6nGGe35gD739gYjIDTdAAAgAElEQVTWAs6exz2k9m/P35+n8SVJkiR1MZBqfJLm3aiX\nJ9mvczHJhsD72o/H9jZKclV7PG22E0hyUJKN+1zfD/jH9uNHZjuOJEmSZqEmcGhRWiXekZqtqroj\nycHA54BTkpwL3A7sQbPH1AlVdXafpp0CF2v03khyQdfHzt5NRyY5tP35kqrqrtJ3MHBikquA64EH\naarrbdHeP7qqvjDyl5MkSZI0dgZSrao6NcnzgCOAnYA1gR8DH66qE2fQ5Y59rj2jPaAJlLp9ArgN\n+F1gV+Bx7efTgI9W1VkzmIMkSZLGyHea1GEg1aWqvgPsOcLzA/+fNN29Ac+fDJw8ShtJkiRJ88NA\nSpIkSRpCATWBd5p8TWpxMpCSJEmShjKpcuUuF1yMrNonSZIkSSMyIyVJkiQNo4BJZKRc27comZGS\nJEmSpBGZkZIkSZKGNIliE1qczEhJkiRJ0ojMSEmSJEnDMiOl1oLISCXZLYmb0UqSJElaFOYtI5Xk\nCcCrgUOALYEr5msukiRJ0jAms4+UFqM5z0gl+b0k/wf4BfD3wAXA71fVNnM9F0mSJGkkNYFDi9Kc\nZKSSrA28CvgfwPbAhcChwL9W1X1zMQdJkiRJGpe5Wtr3C2AZ8BngoKpyGZ8kSZIWHZf2qWOulvat\nSRNIPdAekiRJkrRozVUg9ZvA+4F9gZ8m+VaSVyZ57ByNL0mSJM3OJN6P8j2pRWtOAqmqurOqPlBV\nWwJ7ALcCnwJuTvKhJNvPxTwkSZIkaRzmvGpfVX2rqg4Efgs4GtgT+F6SS+Z6LpIkSdJoMoFDi9G8\nbchbVbdV1T8AzwD2Bm6Yr7lIkiRJQ3FZn1rztiFvR1UV8NX2kCRJkqQFb94CqSRPBNatqhvnaw6S\nJEnSSMwgqTXWQCrJasD+NO89bQmsD9wJ3Ah8DTi1qu5qHz8G+LNxz0GSJEmSJm1sQUyS59JsuLsZ\nj35r7rk0AdZ7k/x1VX2y02xc40uSJEkT54a8ao0lkEryQuAMmo13rwA+B1wO3Eezh9RmwIHA5sA/\nJXk281joYsmrolasmO9ZzL0VSzMXn6ml+r2X6G90S/R7r5hamr+lTNXS/N7LavX5nsK8mHJNmbSg\nzDqQSrIBTeC0OvCGqvrYgEffmWRP4JPAIcDy2Y4tSZIkzaUynlVrHH+V9RZgA+DQaYIoAKrqq8B2\nwPdpsleSJEnS4mH5c7XGEUjtA9xcVR8e5uGquh3YF7htDGNLkiRJ0pwbxztSmwHfHKVBVf0c2GgM\nY0uSJElzozKZYhMWsFiUxpGRWgN4aAz9SJIkSdKiMI6M1C+AbUZpkOTxwNeBDatqszHMQZIkSZq4\n+E6TWuPISH0T+J0kfzDMw0kCnATsCGw6hvElSZIkaU6NI5A6HpgCTk6y7XQPJnkqcA7wR2MYV5Ik\nSZpbVu1Ta9aBVFVdCbyDpnjEhUlOSLJrkieksX6S5yU5DrgK2BW4EDh9tmNLkiRJ0nwYxztSVNU/\nJlkOvAd4Y3v06pQj+XR7f6hy6ZIkSdKCYYU9tcaxtA+AqjoG2Ar4GHADTeDUOX4NfAZ4XlW9tqru\n55HASpIkSVocXNqn1lgyUh1VdQ1tNirJ44D1gbvawKnXPwAnjnN8SZIkSZoLYw2kulXVA8AD09y/\nGrh6UuNLkiRJY2cGSa2xLe2TJEmSpKViYhkpSZIkaZUyqXeazHItSmakJEmSpFVAkuck+askpyX5\neZJK8uA0z78syaeTXJ7k9iTLk9ya5Mwkew1os3vb76DjgmnGe1WS7yS5J8m9SS5O8rqVfKe1krw7\nyU+SPJjkF0k+lWTjlbQ7KMlF7Th3tN9pl+najMqMlCRJkjSshV3+/B3Ay0d4/iBgP+AKmn1e7wE2\nBfYE9kxyZFW9c0Dba4HzBlx/lCQfAf4CeAg4H7gf+H3gE0l2rqo/79NmLeBsYBfgZuBL7fxeC+zd\ntnvUeEmOBQ6jqdfwdWAt4EXAi5O8oqq+MOA7jWTigVSSHYC1q+rbkx5rttoo9QhgJ2BN4Ergw1X1\n6RH7eQ6wD/BCYDPgScCtwLnA+6rqh9O0/R3g74DdgXWBa4BPAcdX1dSIX0mSJEljlIW9DO984DLg\n4vb45Uqe/3vgkKr6VffFJDsCZwFHJPlsVV3Zp+15VfWaYSaVZH+aIOrXwB5VdUl7/WnAN4CDk3y9\nqv6tp+nf0ARR5wMvrqp723aHA8fQ/Bn5+T1j/QFNEPUrYOeq+ml7fWfgHODEJOdU1a+Hmft05mJp\n3z8D35yDcWYlyb7At4E/BH4IfA3YHDipjWqH7ecxwPeAdwFbAJcCX6aJvl8FfC/JAQPa7tS2PQD4\nj7bdk4Fjgc8lWdB/BSJJkqT5U1Xvrap3VdVXquqWIZ6/tDeIaq9fCHyWZt/X3ccwtb9oz+/vBFHt\nOD8H3tZ+/F/dDZKsAbyp/fjGThDVtjuW5s/rz2sTGN3e2p6P6gRRbZvzafa7XQ84eHZfpzFX70gt\n6AAgyQY0e1qtDhxQVbtX1QE0gdA1wGFJXjBClxcCewMbVdVeVfUK4Jk0Uf8awKeSPLlnDo+h2bR4\nbeDwqtqxqg6kCebOB/YHXjOLrylJkqTZWjob8q5oz8vG0Fcn2Dmnz73Ote2T/FbX9V1p9qS9tqou\n7dPulPa8T+dCuxTwhT33p20zGxabaLyOJjr9UlWd1rnYRvJvbz8ePkxHVfVwVe1UVWd0L8Vrf34H\ncBXweKD3Bb59aZYBXlZVx3W1u5d2k+Nh5yBJkiTNVJJtgQOB5TTvKPWzeZJ/SPLxJO9J8tIkg2KL\nddrzo5bTVdV9NCu3ALbrutX5+RL6u6TnOWiSII8Fbquqm6Zps+2APkdisYnG3u25X+R6BvAgsEeS\ntapqYOWTlamqSnI5zT/k3xx2DlV1aZL/ALZOsmlVXT/TOUiSJGlB2izJFf1uVNVWkxw4yT40q5/W\nADaheS9pOfD6qrpuQLNd2qPb5Un2715S17qN5s++Twd+3DP2b9AEP9AUkujYpD33C4i6r2/SdW3a\nNlV1X5I7gQ2SPL6q7hnQ91DMSDU6UemjIt6qWgb8iKbax7PGMNZ/a8+9L//NJOqWJEmSZms74NXA\nK2mW1C0D3gL83z7P3gUcTVOc7Unt8ULgAmAb4BtJ1utpc257fk2f/l7b9fPju35etz3fP2DO9/U8\nN0ybQe1mZMlnpJI8gWb9JUwf8e5AE+VeNouxdqVZI7qMpphFt5lE3dON1fdvNGiWD0qSJGkGJli1\n79pJZ54GqaqjgKPad4yeQVMc4qM0Jcb3bxMLnWcvpSmm1u2b7Z9zvwXsRvNaynu67h9NU0ztwCQ3\nAB+iWfF1APBO4GGauKS7QnWnxsKgX/F+NRhW1mZQuxkxI/Vfo9FRIt6RtAHbp9qPx1XVzQPmMbE5\nSJIkaZYq4z8WiKp6sKp+VFVvBD5I8+rJm1bSrNN2BfDe9uNLeu5dSpPxeoCm/sCNNFsDfQT4Lk2l\naviv71B1lt2tQ39rt+d7u66trM2gdjOySmSkkpwCbD1is4Oq6iKGi0pn9W94ktWBf6GpwHcRTeQ9\nyChR9+BOBvyNRpup2nKUviRJkrTkfAZ4M80Gv8cM2abzbtRTe29U1b8mOYemiMWzaApMnEuzye75\n7WPdK6pubM8bDxhr457nVtomyTo0K9HunO37UbCKBFI0L6aN+v5SJxq9p+fa3dM8O9PI9eM0Vfqu\nBvbqTo92uRfYgNGibkmSJM2VSZUrX5gl0G9vz08Zoc0G7bnvn1fbFVkf6L6W5InA9jR/Ju+uFdB5\nnebZA8bqXP9h17WraQK0pyTZuE/lvn5tZmyVWNpXVTtUVUY8zmnb3k3z0hyMFvEOJcnRNJt+/Qx4\nUVXdPuDRmUTdkiRJ0iQ8vz1fO0Kb/dvz90do8yaaaoH/XFUPdF3/Ds2f0TdLsn2fdge05690LrTt\nv9lzf9o2s7FKBFJjMDDibXdV3pomur16lE6T/DXwlzRrQF9UVT+byRx6ro8lgpYkSdIMrCIb8ibZ\nMMlbk6zf596LgPe1H0/suXdIkif1XEuSQ4DDaL7Rx/r0uUOfa38K/C1N9utd3ffaFVwfaj9+qF2W\n12l3OE3V7fOq6uKebo9tz0ck2byrzc7AITSrz/5P71xmYlVZ2jdbZwDPo4lSP9Nzb2+a0udnjrKH\nVJLX01QruRN4SVWtLAg7AzioncNRPX1tT1M2/cppavlLkiRpCUuyF/COnstrJrmg6/ORVXUGzWsj\n7weOTPI9mgrR6wDPpNnzFJoCaaf29PfXwAeTXAnc0F7bBvhtmqp7b6mqfhmpi5NcQ7OP1P3A79K8\nmvMrYM8Bq7aOAvag2a/qp0n+nWYvqh3bdq/tbVBVZyU5nqZ8+w+SfANYE3gRTRLpVVV1R5+xRjYX\nGakrGLw30kLxSZro9OVJ9utcTLIhj0Tjx/Y2SnJVezyt5/oBNCUj7wVeWlU/GGIOXwCuA7ZLclhX\nX+sAHx40B0mSJM2d1PiPMXoKTZDROaApWNZ9rfPO0600FfS+RbO9zh/RVNtbC/gs8IKqOrzPGMcA\n/w94As3+UX9IE1N8Btipqj44YG7H0bwHtVs71mo0f7bdsqq+169Bm8R4AXAkTfD1RzS1ET4NbF9V\n1wxodyhNkPVjmgBqF+Bs4Pl9AsMZm3hGqqr6rU9cUKrqjiQHA58DTklyLk2KcQ+ayh4nVNXZfZp2\nClys0bnQBl8n0/zLcR1wSJvq7PXFqvpi1xyWt+nNs4BjkxxIE+XvRlP55Iv0pFYlSZKkjqo6CThp\nyGfvp9nf6egRx/ggTWn0UefWLygbpt0DNBWvp6t63a/dSQz5azFTLu1rVdWpSZ4HHEGzU/OaNFHs\nh6tqlABm7bYtNGnObQY8dz1NcNQ9h+8m+T3g3cDuNCnPa2ki/w9U1RSSJEmaPwuzwp7mgYFUl6r6\nDrDnCM8/am+nqrqeWew7VVVX0L/KiCRJkuabgZRaVu2TJEmSpBGZkZIkSZKGNObiEFrEzEhJkiRJ\n0ojmJCOV5ElV9au5GEuSJEmaiAIe/Yr8ePrVojNXS/tuTXIDcHHX8f2quneOxpckSZKksZmrQCrt\n8UTgb2j2ZppK8hMeCay+V1UXDO5CkiRJmmdmj9Saq0DqFTS7GT8G2BVYBjwX2AH4PWB/4HHA6nM0\nH0mSJGkkYTLFJiawWFBzYE6KTVTVqcDvAN8HLgL+HPhCVR1eVbsBT6DZfFaSJEmSFrw5q9pXVfdV\n1V8CO9NkpX6c5OXtvamqunyu5iJJkiTNSE3g0KI05+XPq+ryNgv1d8DHk5yeZI25nockSZIkzdSc\nbsib5InAdl3H3cBLgccDd8zlXCRJkqRRuSGvOuZqH6mv0AROv0kTPP2A5n2pdwGXVJVBlCRJkqRF\nY64yUi8FrgOOAD5bVdfN0biSJEnS+JiRUmuu3pE6mabk+VHANUluS3Jmkr9Lsk+S35ijeUiSJEkz\nM4lCExacWLTmJCNVVX8GkOTxNPtGdY5X02SpKskvquq35mI+kiRJkjQbc1psoqruAb7ZHgAk2ZBH\nNueVJEmSFiyLTahjTgOpbm0Fv3Wr6kbgK+2hubJixXzPYM5lamq+pzAvsjS/NizR771U/3lPTWW+\npzAvHq4538VkQZia+91bFoTlS/A/bMYsWsjGGkglWQ3YH9gT2BJYH7gTuBH4GnBqVd3VPn4McBCw\n+jjnIEmSJEmTNrZAKslzgc8AmwG9fzX4XJoA671J/rqqPjmucSVJkiRpro0lkEryQuAMYE3gCuBz\nwOXAfTR7R20GHAhsDvxTkmczdxUDJUmSpPFwvaFasw6kkmxAEzitDryhqj424NF3JtkT+CRwCLB8\ntmNLkiRJ0nwYR1boLcAGwKHTBFEAVNVXge2A79NkryRJkqRFIzX+Q4vTOAKpfYCbq+rDwzxcVbcD\n+wK3jWFsSZIkae64Ga9a4wikNgMuHKVBVf28qjaqKiv2SZIkSVp0xlFsYg3goTH0I0mSJC1sZpDU\nGkdG6hfANqM0SPL4JOcnuXYM40uSJEnSnBpHIPVN4HeS/MEwDycJcBKwI7DpGMaXJEmSJm8ChSbi\ne1KL1jgCqeOBKeDkJNtO92CSpwLnAH80hnElSZIkaV7MOpCqqiuBdwAbARcmOSHJrkmekMb6SZ6X\n5DjgKmBXmuIUp892bEmSJGlOWbVPrXEUm6Cq/jHJcuA9wBvbo1fa86fb+0OVS5ckSZIWCvd9Usc4\nlvYBUFXHAFsBHwNuoAmcOsevgc8Az6uq11bV/TwSWEmSJEnSojKWjFRHVV1Dm41K8jhgfeCuNnDq\n9Q/AieMcX5IkSZooM1JqjTWQ6lZVDwAPTHP/auDqSY0vSZIkSZMysUBKkiRJWuWYkVJrbO9ISZIk\nSdJSYUZKkiRJGkKYTNU+K7AtTgZSkiRJ0jAmte+TywUXJZf2SZIkSdKIzEhJkiRJwzJ7pJYZqS5J\ndklyZpI7ktyb5KIkr55BP89J8r+T/HuSXyR5KMnPknwmybYD2myapKY5fjn7byhJkiRpHMxItZLs\nC3yeJrj8NnA78ELgpCTbVdXhQ/bzGOB77cfbgYuA+4HtgVcBf5zklVV1yoAubgG+1uf6XcN+F0mS\nJE3GJIpNaHEykAKSbACcCKwO7F9Vp7XXNwLOAw5LcnpVfWvILi8EjgS+WlVTbV+rAX8H/C3wqSTn\nVNXtfdpeVVWvmdUXkiRJkjRRLu1rvA5YD/hSJ4gCqKpbgLe3H4fKSFXVw1W1U1Wd0Qmi2utTwDuA\nq4DHA3uNa/KSJEmaIzWBQ4uSgVRj7/bcb7ndGcCDwB5J1prNIFVVwOXtx9+cTV+SJEmS5o9L+xqd\nAhCX9N6oqmVJfgTsADwLuGyWY/239jyoeMRGSd4NPJXmvagLgS9X1bJZjitJkqRZ8h0pdSz5QCrJ\nE4D12483DXjsJppAahNmEUgl2RV4DrCM/gUlALYA3tlz7cYkf1xVF850bEmSJI2BgZRaSz6QAtbt\n+vn+Ac/c1+fZkbQB26faj8dV1c09jzwEfBT4N+DHNMsJt6R5r+qlwNeSbF9V1w853hUDbm024tQl\nSZIk9VglAqkkpwBbj9jsoKq6CMgwQ4w+q67GyerAvwCb05RD78040QZWb+i5fAGwV5KTgVcCfwO8\nfjZzkSRJ0gxNqjiEWa5FaZUIpIBNad5fGsXa7fmenmt3T/PsvSOO0fFxmip9VwN7zeB9p/fQBFIv\nGbZBVW3V73qbqdpyxPElSZIkdVklAqmq2mEWbe9OchdN+fONgSv7PLZxe75x1P6THA0cDPwMeNGA\nvaNW5qft+akzaCtJkqQxmdUyJa1SLH/e6BSQeHbvjSRr0CwbfIgmozS0JH8N/CVwK00Q9bMZzm+D\n9jzTjJgkSZKkMTKQapzRng/oc29vYC3g7Kp6cNgOk7yeZknencBLqmqkIKzH/u35+7PoQ5IkSbPl\nhrxqGUg1PknzbtTLk+zXuZhkQ+B97cdjexsluao9ntZz/QCaCnz3Ai+tqh+sbAJJDkqycZ/r+wH/\n2H78yJDfR5IkSROQGv+hxWmVeEdqtqrqjiQHA58DTklyLnA7sAfNHlMnVNXZfZp2Clys0bnQBl8n\n0wSp1wGHJDmkT9svVtUXuz4fDJyY5Crgeh4pf75Fe//oqvrCDL+iJEmSpDEykGpV1alJngccAewE\nrEmzn9OHq+rEEbpau20LsE179HM90B1IfQK4DfhdYFfgce3n04CPVtVZI8xBkiRJk2AGSS0DqS5V\n9R1gzxGef1ThlnbD3JELulTVyTSZLEmSJEkLnIGUJEmSNCwzUmpZbEKSJEmSRmRGSpIkSRqSVfbU\nYSAlSZIkDWNS+z4ZnC1KLu2TJEmSpBGZkZIkSZKGECaztG/kcs9aEMxISZIkSdKIzEhJkiRJw/J9\nJrXMSEmSJEnSiMxISZIkSUOy/Lk6DKSWoipqxYr5nsWcy4ql+V++LL1/1ABkar5noLlUtTRf1X54\navX5nsK8WF5L9Xsvvd/Hlt431mJiICVJkiQNy+hOLQMpSZIkaVgGUmpZbEKSJEmSRmRGSpIkSRpG\nTajYhFmuRcmMlCRJkiSNyIyUJEmSNCyzR2qZkZIkSZKkEZmRkiRJkoaUJbifl/ozkJIkSZKGZRyl\nlkv7JEmSpFVAkuck+askpyX5eZJK8uA0z78syaeTXJ7k9iTLk9ya5Mwkew1os3vb76DjgmnG++9J\nzk1yZzvWL5J8LslzBjx//UrG2mJAu9WSHNp+rweS3Jbk80m2XNmv4SjMSEmSJElDmkj58/F5B/Dy\nEZ4/CNgPuAK4ELgH2BTYE9gzyZFV9c4Bba8Fzhtw/VGSnAC8CVgOfBv4FbAl8Apg3yT7VtVXBoz1\n6QHX7+ozToB/Aw4A7gTOAJ4M7A/sleQFVXXhgP5GYiAlSZIkrRrOBy4DLm6PX67k+b8HDqmqX3Vf\nTLIjcBZwRJLPVtWVfdqeV1WvGWZSSbajCaJ+DexcVVd33TsMOBb4ENA3kBp2nNZraYKonwK7VdUt\n7Tj7A6cAJyfZoqoeHqHPvlzaJ0mSJA2rJnCMa2pV762qd1XVVzoBxEqev7Q3iGqvXwh8Fgiw+xim\ntlt7/rfuIKod6/9v796jJSvqQ49/f/IQgWtABVERIYRIGGQUQdHIw6BRHj5gMKxEBS5xiYoG0cRr\nvCoaNYlGQRQlXHnpxURlEFAxGEXRcMMjCkGYcRSG10IQwSHAADMI87t/VHVsN91nep/TZ87p7u9n\nrb1qunZV7eqamsfv1N61T6CsHD0jIrYcwrXeWdN3dY9BZp4DfA3Ynnardn0ZSEmSJElqeqSmDw2h\nrdUDXu9Rt+q1ERHbUW4XfJByS1/T4pq+YibX6fDWPkmSJGkQOUvPSM2z564iYhfgUMrzTBf1KbZD\nRPwd8ETgLsrzUhdm5poeZb8HPAwcGhGf7HFr32bAP2Vmz4ArIv6KspK0mvI817mZeWePogtrem1m\n/rrH+Ssb5WbEQEqSJEka1DwLeoYhIl5B2YxhA2Ab4IWUIOqNmXljn2ovrEe3ayJiUWZe152ZmddH\nxF8CJ9Qy3ZtN7Ah8EXjTFF38WOPzCRHxF5l5WiN/m5re2qedWxvlZsRASpIkSZp720fEkl4nMnPB\nLF97IXB41+dVwNuBL/Qoew/wD8A5lA0dAJ5N2bhiD+DbEbEwM3/rNr3MPDEi7gBOA/btOnUD8N3M\nXNnjWl+jrGb9CLgT+F3gSOAY4NSI+FVmntdVftOaPtDne97fKDcjPiMlSZIkDSAot/YN/Zjj75WZ\nH87MAB4HPAs4HTgZOD8iNmyUvSoz35WZl2fminp8F3gR8G/AM4Cju+tEcQLwz5Tg7PcpwczzgduA\n0yLi4z369ReZeW5m3pKZD2bmksx8J/CWWuSjjSqdoVwn64YGUpIkSdLcW56ZC3od66oDmbkqM6/N\nzKOBTwMHUrYtH6TuI/wmsHlZ4/ThlBWu8zLzzZl5XWben5lXAPsDPweO7feC3R5OBX4J/H7dYKLj\nvppu0qdeJ7/X6ldrBlKSJEnSoObx9udDdlZN22wV3rnV7ymN/NfXdHEjn8y8D7iQEpfsNchF6oYW\nnRf/dl/rlppu3afq1o1yM2IgJUmSJKnprppu0aLO5jVtrvh0Aph7+9Tr5D9hhte6uqY7R8QGPers\nWtMft7hOXwZSkiRJ0oBm4xmpeWrvmi6fstRvW1TTHzXyf1HT3frU272mNw1ykYhYADyTsqnEsk5+\n3WHwJ5RnvQ7oUfWQmn5jkOusjYGUJEmSNGEiYsuIeGdEbNbj3Ev5zZbjZzTOHRURT2zkRUQcBRxL\nuVnxHxtNdnbWe0dEPK9R962UjSruA/61K/9lEfHcHn3bBTibsrHEqZnZfGHw8TX9WERs2VXvYOCV\nwI1d/ZkRtz+XJEmSBpXzdwkpIg4A3tfI3jAiLuv6/KHMvADYGPg48KGI+CHlHUubUHbU62z6cEJm\nntNo76+BT0fEUuDmmvcsYDtgDXBMZjZXpE4GDqYETJdGxKWU3foWUN4l9QhwdGau6KrzAuC4iLiZ\nsip2Z73GrpQY5vu1L02nUzawOAhYFhEXAU+irLCtAl7X52W9rRlISZIkSYOYrVvxhtfmFpQtxbtF\nI6/zzNMvgXcB+1ACmt0od6vdDnwJOCUzL+5xjU8Af1zr7Et5ie/tlM0pPpWZ/9GskJmrImJf4K3A\nocAutU93Ujag+ERmXtao9i3g6ZTb/hYCv0N5luoSygt8z6g7BTavtSYiXkN519SRlJ0H7wfOBd6f\nmT3f1TUdBlKSJEnSGMjMM4EzByz7AOXFuv/Q8hqfpmyN3rZvD1Fuuzt+bWVr+UuBS9tep9Z9pM21\npstASpIkSRrU/L2zT+uYm01IkiRJUkuuSEmSJEkDijVz3QPNF65IdYmIF0bENyNiRUSsjIgrIuLw\nabSzV0R8LiKujIg7IuKh2ub3IuJ1a6n7BxFxdkTcGREPRsQ1EXFsRPh7JUmSJM0TrkhVEXEQZU/6\nxwA/oLzNeV/gzIhYmJnvaNHcK4E3AD8DrgLuBp4G7AnsExF/nJmH9ejDHsBFlO0or6C8lGwvyoNy\nfxgRr8mcx3tuSpIkjTv/J6bKVQ4gIjanvGxsPeCQzNwnMw+h7KF/PXBsRLy4RZOnA0/LzGdm5ssz\n808zc6/a3m3A6yPi5Y0+rE/ZNnJj4B2Z+fzMPBTYgbJjySLgiBl9UUmSJM1I5PAPjSYDqeINlL3p\nz8/Mr3YyM/MOyv76AAOvSGXm0sy8rUf+9cBn68c/apw+CNgeuDozT+iqsxI4um0fJEmSJM0eA6ni\nwJou7nHuAspbkF8SERsN4VqdF4c9NGgfMvMq4AZg54jYdgh9kCRJUlsJZM7CMddfTNNhIFXsUtMr\nmyfqy8OuBTYCnjmTi0TE04Gj6scLG6cX9utDI39hn/OSJEmS1pGJ32wiIh4PbFY/3tqn2K3AbsA2\nwNUt2n4BJXBaD3gq8L2D2CYAACAASURBVCLKmL83My9pFN9mgD50l1vbtZf0ObX9IPUlSZL0aD7T\npI6JD6SATbt+/UCfMvf3KDuI7YHu7dPXAMcBH5+iH8PugyRJkqQhG4tAKiIWAzu3rHZYZl4BxCCX\naN8ryMyzgLMiYkNgW+Aw4H3AgRGxX2be3avaMPqQmQt6NlJWqnZq05YkSZIqV6RUjUUgRQlS2j6/\ntHFN72vk3TtF2ZUtrwH893NWPwPeGxG/orwX6m+At3UVWwlsDmyylv5Oqw+SJEmShmcsNpvIzN0y\nM1oeF9e69wL31Ka27nOJTv4tQ+juWTV9VSO/0/a66IMkSZJaCmbnPVLTuvVJc24sAqkh6GwgsWvz\nRERsQLltcDXw0yFcawXlWaktBu1DI//HQ+iDJEmSpmM2tj/XSDKQKi6o6SE9zh1I2fr8osxcNYRr\n7UkZ9+WD9iEingP8LrA0M28cQh8kSZIkzYCBVHEq5dmoV0XEwZ3MiNgS+Fj9eHyzUkQsq8fTGvkf\niIitepTfDfhc/XhG4/S5wI3Awog4tqvOJsBn+vVBkiRJ685s3Nqn0TQum03MSGauiIgjga8AiyPi\n+8BdwEso75j6VGZe1KNqZ4OLDRr5xwHviYgrgZuADYHtgGfX818BTmz04dcR8TrgO8DxEXEocDNl\nBespwHk8OviSJEmSNAcMpKrMPCci9gLeC+xBCX5+AnwmM9sGMG8DXkwJnHamBFp3AucDZ2bmeX36\n8O8RsTvwQWCfWn858Angk5m5pu33kiRJ0pAks7P9uatSI8lAqktm/j9gvxble26ykpknASdNsw9L\n6P2sliRJkqR5wkBKkiRJGpDPNKnDQEqSJEkaSMIa7+1T4a59kiRJktSSK1ITKh95ZK67sM49ZgK/\nM0DMyk/O5r+Y0K1ZYk3PRzfH3prej6yOvYdzMn8e+lBO5n9fHpnAVYt5+Y3nZac0Fybzb2BJkiRJ\nmoHJ/JGOJEmSNA1uNqEOV6QkSZIkqSVXpCRJkqRBJJCzsCTlKtdIMpCSJEmSBuStferw1j5JkiRJ\naskVKUmSJGlQrkipckVKkiRJklpyRUqSJEkaUMzGZhMaSa5ISZIkSVJLrkhJkiRJg1oz1x3QfOGK\nlCRJkiS15IqUJEmSNIDI2XlGyndTjSYDKUmSJGlQBj2qvLVPkiRJklpyRUqSJEkalNufq3JFSpIk\nSZJackVKkiRJGpAbQ6jDFSlJkiRJaskVKUmSJGlQPiOlykBKkiRJGkRCrJmddjV6vLVPkiRJklpy\nRUqSJEkaSM7SrX0uSY0iV6QkSZIkqSVXpCRJkqRBuXikyhUpSZIkSWrJFSlJkiRpQOH256oMpCRJ\nkqRBGUip8tY+SZIkSWrJFSlJkiRpEAn4Ql5VrkhJkiRJUkuuSEmSJEkDcrMJdbgiJUmSJEktuSIl\nSZIkDcoVKVWuSEmSJElSSwZSXSLihRHxzYhYERErI+KKiDh8Gu3sFRGfi4grI+KOiHiotvm9iHhd\nnzrbRkROcfxi5t9QkiRJM5I5/EMjyVv7qog4CDibElz+ALgL2Bc4MyIWZuY7WjT3SuANwM+Aq4C7\ngacBewL7RMQfZ+ZhfereAVzYI/+eFteXJEnSbJiN7c81kgykgIjYHDgDWA9YlJlfrflPBi4Bjo2I\nr2fm9wZs8nTg+My8rXGd3wO+D7w+Iv4pM3sFTMsy84hpfhVJkiRJ64C39hVvAH4HOL8TRAFk5h3A\nu+rHgVekMnNpM4iq+dcDn60f/2j63ZUkSdI6l0nMwuHtfaPJQKo4sKaLe5y7AFgFvCQiNhrCtR6p\n6UNDaEuSJEnSHPDWvmKXml7ZPJGZD0XEtcBuwDOBq6d7kYh4OnBU/djrtj6AJ0fEB4GnUJ6Luhz4\nWmYaeEmSJM01V49UTXwgFRGPBzarH2/tU+xWSiC1DS0CqYh4ASVwWg94KvAiypi/NzMv6VNtR+D9\njbxbIuJPMvPyQa8tSZIkafZMfCAFbNr16wf6lLm/R9lBbA90b5++BjgO+HiPsquBk4EvAz+h3E64\nE/A+YH/gwoh4TmbeNMiFI2LJFH2SJEnSdLgipWosAqmIWAzs3LLaYZl5BRCDXKJ9ryAzzwLOiogN\ngW2BwyiB0YERsV9m3t1V9nbgLY0mLgMOiIgvAn8GvAd443T6IkmSpCEwkFI1FoEUJUh5Zss6G9f0\nvkbevVOUXdnyGkB5zoryTqn3RsSvgOOBvwHeNmATf0sJpF7W4poLeuXXlaqdBm1HkiRJ0qONxa59\nmblbZkbL4+Ja915+87LbrftcopN/yxC6e1ZNX9WiznU1fcoQri9JkqTpWjMLh0bSWARSQ9DZQGLX\n5omI2IBy2+Bq4KdDuNYKyh+ZLVrU2bym01oRkyRJkjRcBlLFBTU9pMe5A4GNgIsyc9UQrrUnZdyX\nt6izqKY/GsL1JUmSNA2RzMoLecPHrkaSgVRxKuXZqFdFxMGdzIjYEvhY/Xh8s1JELKvH0xr5H4iI\nrXqU3w34XP14RuPcYRHxqFsLa3/+vn787OBfSZIkSdJsGZfNJmYkM1dExJHAV4DFEfF94C7gJZR3\nTH0qMy/qUbWzwcUGjfzjgPdExJXATcCGwHbAs+v5rwAnNuocCZwREctqnc725zvW8/+QmedO6wtK\nkiRpONy1T5WBVJWZ50TEXsB7gT0owc9PgM9k5hlTVn60twEvpgROO1MCrTuB84EzM/O8HnU+V8s8\nm/Li3sfVz18FTs7M77T+Un08Y6et+faSs4fVnCRJ0qzYYP0FwNK57kaXhDWzEUgZnI0iA6kumfn/\ngP1alO/5fqnMPAk4qeW1vwh8sU0dSZIkSXPDQEqSJEkalLf2qXKzCUmSJElqyUBKkiRJGlTm8I8h\niYjnRsS7I+KrEfHziMiI6Pv6noh4ZUR8PiKuiYi7IuLXEfHLiPhmRBzQp84+td1+x2VTXO9PI+L7\nEfFf9Vq3RcRXIuK5U9TZLCI+GRE3R8Tqmp4YEZtNUecxEfH2+r0ejIg7I+LsiNipX53p8NY+SZIk\naTy8D3hVi/KHAQcDS4DLgfuAbSl7BuwXER/KzPf3qbscuKRP/qNExKcoG7L9GvgB8CvKDtWvAQ6K\niIMy8xuNOk8ELgV2AG4AzgMWAH8B7B8Re2Tmrxp1Avgy5f2w/0V5X+yTKO9lPSAiXpyZl/cbkDYM\npCRJkqRBJLPzjNTwmrwUuBr4j3r8Yi3lPwIc1SMYeT7wHeC9EfGlzOy1deIlmXnEIJ2KiIWUIOpu\n4AWZ+dOuc8dS3td6EvCNRtUTKEHUV4FDM/PhWqcTlB0PHN6o8z8pQdR1wJ6ZeUetswhYDHwxInbs\ntDUT3tonSZIkjYHM/GhmHpeZ3+gEEGspf1UziKr5lwNfAgLYZwhd27OmX+4Oouq1TqCsHD0jIrbs\n5EfEVsBrKStYb2kEPn9FeU3QayPiyY1rvbOm7+oeg8w8B/gasD3tVu36MpCSJEmSBrUmh3/MT4/U\n9KEhtLV6wOvd0/V5P0qs8oNmUJiZq4GvA+vR9eqiiNiOcrvgg5Rb+poW1/QVA/d8CgZSkiRJ0kAS\ncs3wj3n2Qt6I2AU4lLIadFGfYjtExN9FxP+JiL+NiP0jol9s8T3gYeDQiHhm41rHAptRVqu6A66F\nNb2yT5tXNsp1//razPz1gHWmzWekJEmSpAkWEa+gbMawAbAN8EJKEPXGzLyxT7UX1qPbNRGxKDOv\n687MzOsj4i8pzzxdExHdm03sCHwReFOjrW1qemuf69/aKDfdOtNmICVJkiQNavZeyLt9RCzpfclc\nMFsXrRby25s2rALeDnyhR9l7gH8AzqFs6ADwbMrGFXsA346IhZnZfZsemXliRNwBnAbs23XqBuC7\nmbmycZ1Na/pAnz7f3yg33TrT5q19kiRJ0gTLzA9nZgCPA54FnA6cDJwfERs2yl6Vme/KzMszc0U9\nvgu8CPg34BnA0d11ojgB+GdKcPb7lGDm+cBtwGkR8fFGt6JzyT7djiny1sm9kgZSkiRJ0iCS2dls\novy3f3lmLuh1rLOvl7kqM6/NzKOBTwMHUrYZH6TuI8BH68eXNU4fTlnhOi8z35yZ12Xm/Zl5BbA/\n8HPg2IjYsavOfTXdpM8lN65p90rW2up08purX9NiICVJkiSp6ayattkqvHOr31Ma+a+v6eJGPpl5\nH3AhJS7Zq+vULTXdus+1tm6Um26dafMZKUmSJGlQs/eM1HxzV023aFFn85o2V3w6Acy9fep18p/Q\nlXd1TXftU6eT/+MedXaOiA167NzXq860uSIlSZIkDSpz+Mf8tHdNl7eos6imP2rk/6Kmu/Wpt3tN\nb+rKuxBYA+zZ/aJegIh4LOVdUGuAf+nk1x0Gf0J51uuAHtc5pKbf6PsNWjCQkiRJkiZMRGwZEe+M\niM16nHsp8LH68YzGuaMi4omNvIiIo4BjKU98/WOjyfNq+o6IeF6j7lspG1XcB/xrJz8zb6dsTrEh\n8NmI6L6T7mOUlbJ/ysxf8NuO75TpDsAi4mDglcCNXf2ZEW/tkyRJkgY1f1eQiIgDgPc1sjeMiMu6\nPn8oMy+gbNbwceBDEfFDyjuWNqHsqNfZ9OGEzDyn0d5fA5+OiKXAzTXvWcB2lBWiYzKzuSJ1MnAw\nJWC6NCIupezWt4DyLqlHgKMzc0Wj3tspW6ovApbVfi4AdqaslB3bYxhOp2xgcVCtcxHwJMoK2yrg\ndX1e1tuagZQkSZI0HragbCneLRp5nWeefgm8C9iHEpzsRrlb7XbgS8ApmXlxj2t8AvjjWmdfykt8\nb6dsTvGpzPyPZoXMXBUR+wJvBQ4Fdql9upOyAcUnMvOyHvXuiojdgQ8Cr6YER3cAJwHH9Qi8yMw1\nEfEa4BjgSMrOg/cD5wLvz8ye7+qajsh5HFVr+CJiyU477bTTkiVDm0OSJEmzYsGCBSxdunTputwC\nvJ+IWLLp+pvv9KIn/enQ277krn9m5cN3z4vvqcH5jJQkSZIkteStfZIkSdIgktl5RsobxEaSgZQk\nSZI0KB+LUeWtfZIkSZLUkitSkiRJ0qDWuCKlwhUpSZIkSWrJFSlJkiRpAEmSuWZW2tXocUVKkiRJ\nklpyRUqSJEkaRDI7z0i5IDWSXJGSJEmSpJZckZIkSZIG5XukVBlISZIkSYNaM/zNJjSavLVPkiRJ\nklpyRUqSJEkalLf2qXJFSpIkSZJackVKkiRJGkQmORvPSLnKNZJckZIkSZKkllyRkiRJkgbl6pEq\nAylJkiRpUGsMpFR4a58kSZIkteSKlCRJkjSo9IW8KlyR6hIRL4yIb0bEiohYGRFXRMThQ2r79IjI\neuwxRbk/iIizI+LOiHgwIq6JiGMjwt8rSZIkaZ5wRaqKiIOAsynB5Q+Au4B9gTMjYmFmvmMGbb8Y\n+J9AAjFFuT2Ai4CNgSuAm4C9gOOBP4yI12T6hKMkSdKcyCRn4xkp/3s3klzlACJic+AMYD3gkMzc\nJzMPAXYErgeOrcHQdNreCDgFWAJcOkW59YGzKEHUOzLz+Zl5KLBDrbcIOGI6fZAkSZI0XAZSxRuA\n3wHOz8yvdjIz8w7gXfXjdFek3gf8HvAm4NdTlDsI2B64OjNP6OrDSuDoGfZBkiRJw5Brhn9oJHlr\nX3FgTRf3OHcBsAp4SURslJmrBm00InYG/go4PTMvieh7V9+UfcjMqyLiBmDniNg2M28atA+SJEka\nnlm5tU8jyRWpYpeaXtk8kZkPAdcCGwHPHLTBujnE54B7+M2q1lQW9utDI39hn/OSJEmS1pGJX5GK\niMcDm9WPt/YpdiuwG7ANcPWATR8N7AEcnpkrBii/zQB96C43pYhY0ufUjsuXL2fBggWDNCNJkjRn\nli9fDvD0ue5Hx/3cy6VrLpyVdjV6Jj6QAjbt+vUDfcrc36NsXxGxNfAR4OLM/ELLfgylD1N4zOrV\nq9csXbp02QzbGTfb13T5nPZi/nFcenNcHs0x6c1x6c1x6c1xebQdgcfNdSeq5TCrQY+/7yNmLAKp\niFgM7Nyy2mGZeQVTbEfefYmWbX8GeCzw5pb1oGyRPuM+ZGbPJafOSlW/85PKcenNcenNcXk0x6Q3\nx6U3x6U3x+XRprjDZp3LzFfOdR80v4xFIAVsS4vnl6qNa3pfI6/Xjxk6ZVeurdGIWAS8EvhQZrZZ\n9VkJbA5s0uf8wH2QJEmSNLvGIpDKzN1mUPfeiLiHsv351sDSHsW2ruktAzT5ipq+NCL2apx7dk0/\nGxH3AidlZmeXvlsogdTWwI9n2AdJkiRJs2gsAqkhuBrYC9iVRiAVERtQbhtcDfy0RZt7THHuOTU9\nr9GHhbUP3+xRZ9ea9gqyJEmSJK1Dbn9eXFDTQ3qcO5Cy9flFg7xDKjOPyMzodQDfr8VeUPM+OUgf\nIuI5wO8CSzPzxkG/lCRJkqTZYSBVnEp5NupVEXFwJzMitgQ+Vj8e36wUEcvq8bQh9OFc4EZgYUQc\n23WNTSibV/TsgyRJkqR1LzJ9OzP89yYRX6Hsjvd94C7gJZR3TH0qM4/pUaczeNtl5k0DXONiYG/K\nitRlPc6/EPgOZZvPy4GbgT2Bp1BuA1yUmWvafjdJkiRJw+WKVJWZ51Cek/oWZVOI/Sn7+R/ZK4ia\npT78O7A7cA7we8CrgLuBvwQOMYiSJEmS5gdXpCRJkiSpJVekJEmSJKklAylJkiRJaslASpIkSZJa\nMpCSJEmSpJYMpCRJkiSpJQMpSZIkSWrJQGoERMRzI+LdEfHViPh5RGRErJpBe5tFxCcj4uaIWF3T\nEyNisynqPCYi3h4R10TEgxFxZ0ScHRE7TbcfMzHXYxIRZ9Zr9jveNP1vN33DHJeI2DsijouIC+rv\nd0bEsgHqzau5Uvs0p+My7vOl/vn5s4j4p4hYGhH3R8R9EXF5RBwTERtMUXdezZe5HpMJmCvrR8QH\n6p+fG+qYrIqI6yLiMxGxzRR159VcqX2a03EZ9/nSp+0d6u9/RsSFU5Sbd/NF48f3SI2AiDiP8nLe\nbqszc6NptPVE4FJgB+AG4IfAgnpcD+yRmb9q1AngK8AhwH8BFwFPorzAeBXw4sy8vG1fZmIejMmZ\nwOGUFzj/okezn8/M77Xty0wNeVz+E1jYyP5pZu44RZ15N1dqv+Z6XM5kjOdLRHwY+N/AGuAqyp+b\nLYA/BB4LXAK8LDMfaNSbd/NlHozJmYz3XNkUuA9YCfwYuA3YEHg2sA1wD/BHmXllo968myu1X3M9\nLmcyxvOlT9vfBfYBAvhWZr68R5l5OV80hjLTY54fwP8CPggcCDwZSGDVNNv6Qq1/DrB+V/6nav7n\ne9Q5sp77GfDkrvxFNf/67rYmZEzOrOf2mev5MYvj8jHgPcBLgefUtpatpc68myvzZFzGer4A7wY+\nAjytkb8DcHNt929HYb7MgzEZ97myPiWYXL+Rvx7wt7Xdy0ZhrsyTcRnr+dKj3T+vbZ1S0wv7lJuX\n88Vj/I4574DHNH7Tpv8X9VbAI8BD3X+x1HOPBX4JPNzj3JJ6zVf3aPP8em7RhI3JvPzHa1jj0qOd\nbRksYJj3c2WOxmWi5kujzT+t7d44ivNlDsZkkufK+sCDte1NRm2uzNG4TMx8AbYEVgDfpqxITRVI\njcR88Rj9w2ekJst+lOfifpCZd3SfyMzVwNcpP/3ar5MfEdsBO1H+Er+gR5uLa/qK2ejwOtB6TNTb\nBMwVTc/VNX1qd+aEz5eeYyKScivkGsoPsICJnyvQZ1wm0KeAxwFvnqqQ80Xr0vpz3QGtU53nOq7s\nc/5KynJ49/MfnV9fm5m/7lOnu9yomc6YdDs4IhZRgq0bga9n5lo3ZBhT4z5XhmES58vv1rT5/MYk\nz5d+Y9JtouZKfabl3cDGwHfqD7I6JnaurGVcuo31fImI/YFDgfdn5vURsfUUxSd2vmjdM5CaLJ1d\nf27tc/7WRrnp1hklM/1+b2t8/mhEnAwck5mT9pPDcZ8rwzCJ8+WYmp7fyJ/k+dJvTLqN/VyJiI9S\nnp95PLALsD2wDHhjo+hEzZUW49JtbOdLRGwCfBb4KfDRAapM1HzR3DKQmiyb1vSBPufvb5Sbbp1R\nMt3vdxVlp7/vUv5S3opy+9+HgbdQnrk6dqg9nf/Gfa7MxETOl7r18ksou2b9feP0RM6XtYwJTNZc\nWUQJEjquBV6bmTc2yk3aXBl0XGAy5suHgWdQdi18aIDykzZfNId8RmqyRE1zLefb1Bl10xkTMvPE\nzDwlM6/LzAcz88bM/Cxla9WHgLdFxNNnob/z2bjPlWmbxPkSEXsDJ1Lmw5GZeVuzSE0nZr4MMCYT\nNVcy8/cyMyjbwr8cWA38KCIObxSdqLnSYlzGfr5ExG6U1bYv5ODbuE/UfNHcMpCaLPfVdJM+5zeu\n6coWdTr5K/ucn++mMyZ9Zea1wNco96m/ZGZdGznjPleGblznS0TsApxHeQ/OMZl5bo9iEzVfBhyT\nvsZ1rgBk5l2Z+S1gX8r7k05u/Od/ouZKxwDjMlXdkZ8vEbE+8DnKO7T+skXViZwvmhsGUpPllpr2\ne0hz60a56dYZJbPx/a6r6VOm1aPRNe5zZbaM1XyJiO0pLwfdDPhAZn66T9GJmS8txmRtxmquNGXm\nPcA3KDuzvbTr1MTMlV6mGJe1GfX5sjXlhcQPAWdHxMWdA/hkLfO8mveNrnoTPV+0bvmM1GTpbLm7\na5/znfwf96izc0Rs0GMHnF51Rsl0xmRtNq/ppP20a9znymwZm/kSEU+lvONlK+DEzPzgFMUnYr60\nHJO1GZu5MoW7arpFV95EzJW16DUuazMu82WrevSyObA3ZdWqw/midcYVqclyIeU9FHtGxJbdJyLi\nsZR3KqwB/qWTXx9u/QnlJ2EH9GjzkJp+o8e5UdB6TKZS63TG6UdD7Oe8NwFzZejGab5ExOaUVZft\ngDNYywPukzBf2o7JWtoam7myFnvXdHknYxLmygAeNS5TGYf5kpk3ZWb0OoAX12LfqnmbddVzvmjd\nmes3Anu0P1jLG8KBt1K2Sv27HufOqvUXA+t35XcegP6/Peq8oZ77GbBlV/7BNf8GYINJGRPgmcCr\ngPUa+VsA59Y6/wnEKM+VRrlta1vL1lJu3s+VdT0ukzBfKM8SXlrrf7n5XUd5vqzLMZmQufJKyq5y\n0cjfGPhIbfd2YNNRmyvrelwmYb5MUX6f2uaFfc6PxHzxGP3DW/tGQEQcALyvkb1hRFzW9flDmdl5\ng/eTKH/B9rov+u3AHpTtVZdFxA+BBcDOlJ909fqJ6enA/sBBtc5F9Rp7A6uA12Xvl97Nmjkek6dQ\nHhr/VUQsA34ObAk8F/gflC1o/yQz1/mOQcMcl4h4A+UfI4DH1vQZjbbekpndLzOed3MF5nxcJmG+\nfITyZ+gR4GHgtPIe0d+WmUc0subdfJnjMZmEubIrcBxwW0RcRbklayvKszBPqJ//JDObt6PNu7kC\ncz4ukzBfpmtezheNHwOp0bAF8PxGXjTyBrpvOjPviojdgQ8Cr6b8JXMHcBJwXGau6FFnTUS8hvIS\nySOBAynvYTiX8pbxJe2+zlDM5Zj8jPKg6x6Ud308j7I97c+Ar1Oeg7i71bcZnqGNC+WB3GZbGzXy\nHt99cp7OFZjbcZmE+dJ5FmM94M+mKHdE94d5Ol/mckwmYa58lfKf/D2B3SlBwoPA9cApwKcz8/Zm\npXk6V2Bux2US5su0zOP5ojETc/CDCkmSJEkaaW42IUmSJEktGUhJkiRJUksGUpIkSZLUkoGUJEmS\nJLVkICVJkiRJLRlISZIkSVJLBlKSJEmS1JKBlCRJkiS1ZCAlSZIkSS0ZSEmSJElSSwZSkiRJktSS\ngZQkSZIktWQgJUmSJEktGUhJ0oiKiAURcUpEXB8RD0TEPRFxZUR8MCKeMNf9m65x/V6SpPESmTnX\nfZAktRQR7wY+DKwH/AS4Bngs8HxgK+Bu4E8y8ztz1slpGNfvJUkaPwZSkjRiIuKdwMeBO4HXZua3\nu849BjgaOAFYA+yVmZfNSUdbGtfvJUkaTwZSkjRCImIHYAmQwB6ZeVWfcm8GPgssB3bMzIfXQd+2\nBW4EPp+ZR7SsO2+/lyRJvfiMlCSNlrcDGwCn9gs2ADLzZEpgsj2waB31bSbG9XtJksaUgZQkjZb9\na/rlAcp+pab7zVJfhmlcv5ckaUwZSEnSiIiIjYFt68cfD1DlmpruNCsdGpJx/V6SpPFmICVJo2Oz\nrl/fN0D5+3vUm4/G9XtJksbY+nPdAUnSwO7p+vXvACvWUv6pNf2vYXekblO+YyN705q+KCLO7FHt\n7zNzWY/8efO9JEkalIGUJI2IzLw/IlYATwB2AS5eS5Xta3pLd2ZE/F/gdcC/ZOb+jXOLgMX145aZ\neWeftl8O7D3FdbfvkX8m8KhAaljfqyMiNgd+DjwMPDUzV66lPUmSWvPWPkkaLd+s6UEDlD20phc2\n8ncHbgIWdmdGxKbAJykBys1TBFFk5j6ZGd0HsF09/fnmuXpcPMvfq+MIymrVQ8BrB2hPkqTWDKQk\nabT8n5oeGRFP61coIvYDdqDcJndOV/7jgd8HTgWeGhFP6qr2AeAG4DrgP4bb7bWa0fdqeCPweeAs\n4E3D7KQkSR0GUpI0QjLz34AvUJ5H+lJErNcsU4OjU+vH92Tm3V2ndwMC+BLl2aSFtc6zgLcARwPP\nZR0HUkP4Xp0yL6Y8u3V6LfvsiHj+rHVckjSxDKQkafS8BbgMeBGwuK4yARARW1Nuk3sqcEpmntKo\n+zxgRWYuB/4TWBgRAZwMnASspuyGt65XpGBm36vjTcAPMvO6zLwWuBxXpSRJs8BASpJGRETsWHfD\n+wzQeX7p1cC/dhX7EeUZqF8DG0fEmfXo3MK3O78Jkn5EWZE6EtgG+GA9n/XcOjGk70VEbFnrndpV\n7zTg0Ihwq3RJ0lC5a58kjY6tgMP75HdsWdMNgNd35X8AuIsSjHyh5l1JCZ72B46qu+ftDvwsM+8d\nYr/XZhjfC+DPSnGAmgAAAT1JREFUgQf5za6DAP8MHF/bP3EIfZUkCYDIzLnugyRpHYiIJwO/AF6d\nmedHxB8AS+naBj0iLgFuzMzXT9HUvBMRjwGuB76VmW9unDsNeEFm7jQnnZMkjSVv7ZOkybF7TTu3\n9i0DtqBuOV43eHgOc/N81Ey9jLL9+mk9zp0G/EFE9HvvlSRJrRlISdLk2B24PTNvA8jirsxcXc8v\nADYGfjhXHZyBo4AfZ+aj+p6Z/05ZeXPTCUnS0HhrnyRJkiS15IqUJEmSJLVkICVJkiRJLRlISZIk\nSVJLBlKSJEmS1JKBlCRJkiS1ZCAlSZIkSS0ZSEmSJElSSwZSkiRJktSSgZQkSZIktWQgJUmSJEkt\nGUhJkiRJUksGUpIkSZLUkoGUJEmSJLVkICVJkiRJLRlISZIkSVJLBlKSJEmS1NL/B6MoRIOObJHP\nAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x7fe07a5332e8>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "plt.pcolormesh(Ops, Oms, vals)\n",
    "plt.xlabel(r'$\\Omega_M + \\Omega_\\Lambda$')\n",
    "plt.ylabel(r'$\\Omega_M - \\Omega_\\Lambda$')\n",
    "plt.colorbar()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "metadata": {},
   "outputs": [],
   "source": [
    "params = [\n",
    "    ('O_m', 0.3,  1e-3, 0, 1),\n",
    "    ('O_l', 0.7,  1e-3, 0, 1),\n",
    "          ]\n",
    "# minuit = setup_minuit(chi2_fit, params)\n",
    "\n",
    "    \n",
    "# results = perform_fit(minuit, params)"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.6.2"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}