OLSvsODR.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "from astropy import units as u\n",
    "from scipy import odr, sqrt\n",
    "%matplotlib inline"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Up to now, we've been using the Ordinary Least Squares method to find the equation of the best fit line to our data. Unfortunately, this method only takes errors in our y-values into account. In lab 1, we'll need to fit a line to data where we have uncertainties associated with both our x and y measurements."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "def OLSfit(x, y, dy=None):\n",
    "    \"\"\"Find the best fitting parameters of a linear fit to the data through the \n",
    "    method of ordinary least squares estimation. (i.e. find m and b for \n",
    "    y = m*x + b)\n",
    "    \n",
    "    Args:\n",
    "        x: Numpy array of independent variable data\n",
    "        y: Numpy array of dependent variable data. Must have same size as x.\n",
    "        dy: Numpy array of dependent variable standard deviations. Must be same \n",
    "            size as y.\n",
    "\n",
    "    Returns: A list with four floating point values. [m, dm, b, db]\n",
    "    \"\"\"\n",
    "    if dy is None:\n",
    "        #if no error bars, weight every point the same\n",
    "        dy = np.ones(x.size)\n",
    "    denom = np.sum(1 / dy**2) * np.sum((x / dy)**2) - (np.sum(x / dy**2))**2\n",
    "    m = (np.sum(1 / dy**2) * np.sum(x * y / dy**2) - \n",
    "         np.sum(x / dy**2) * np.sum(y / dy**2)) / denom\n",
    "    b = (np.sum(x**2 / dy**2) * np.sum(y / dy**2) - \n",
    "         np.sum(x / dy**2) * np.sum(x * y / dy**2)) / denom\n",
    "    dm = np.sqrt(np.sum(1 / dy**2) / denom)\n",
    "    db = np.sqrt(np.sum(x / dy**2) / denom)\n",
    "    return([m, dm, b, db])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Orthogonal Distance Regression is works to find the best fit line that minimizes the *perpendicular* distance from the data points to the line. \n",
    "![Comparison of perpendicular and vertical distance measurements, from http://blog.rtwilson.com/orthogonal-distance-regression-in-python/](http://blog.rtwilson.com/wp-content/uploads/2015/11/Normal_vs_ODR.png)\n",
    "Unlike the OLS method, there's no formula that we can plug values into to get the answer. ODR finds the optimal line *numerically*. This means that it starts with some values of slope and intercept and then fiddles with them. Here's a function that implements this algorithm. You'll have to provide it with a 'guess' for the slope and uncertainty. You can generate this with an OLS fit.\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 56,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# Use the SciPy ODR tool to fit a line to the data, I wrote a function that takes as input your x and y values,\n",
    "# their uncertainties, and initial guesses for the slope and intercept of the best fit line. It returns the best fit\n",
    "# slope and intercept and their uncertainties in a list with [slope,slope uncert.,intercept,intercept uncert.]\n",
    "\n",
    "def ODRfit(x, y, dx, dy, m_init, b_init):\n",
    "    linear_model = odr.Model(linearEqn)\n",
    "    data = odr.RealData(x, y, sx=dx, sy=dy)\n",
    "    myodr = odr.ODR(data, linear_model, beta0=[m_init, b_init])\n",
    "    output = myodr.run()\n",
    "    return[output.beta[0], output.sd_beta[0], \n",
    "           output.beta[1], output.sd_beta[1]]\n",
    "\n",
    "def linearEqn(p, x):\n",
    "    return p[0]*x + p[1] "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let's try it out on some random data:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 53,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "xs = np.abs(np.random.randint(0,20,100))\n",
    "dxs = np.abs(np.random.randn(100))\n",
    "ys = ((2+np.random.randn(100))*xs)+3*np.random.randn(100)\n",
    "dys = 3*np.abs(np.random.randn(100))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Start with an OLS fit, which will only take the y uncertainties into account, but will provide an initial value of slope and intercept for the algorithm."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 57,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "('Slope=', 1.6449485432374951, '+/-', 0.0029388038343843453)\n",
      "('Intercept=', 5.2847513412531102, '+/-', 0.011395242763823692)\n"
     ]
    }
   ],
   "source": [
    "olsfit = OLSfit(xs, ys, dy=dys)\n",
    "print('Slope=',olsfit[0],'+/-',olsfit[1])\n",
    "print('Intercept=',olsfit[2],'+/-',olsfit[3])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now we can use the ODR fit function with the slope an intercept values found by the OLS fit."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 58,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "('Slope=', 2.4852643422524223, '+/-', 0.18446742137743272)\n",
      "('Intercept=', -3.0329852349563051, '+/-', 1.9339228436053264)\n"
     ]
    }
   ],
   "source": [
    "odrfit = ODRfit(xs,ys,dxs,dys,olsfit[0],olsfit[2])\n",
    "print('Slope=',odrfit[0],'+/-',odrfit[1])\n",
    "print('Intercept=',odrfit[2],'+/-',odrfit[3])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let's make a plot and compare the two lines."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 62,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(-2, 25)"
      ]
     },
     "execution_count": 62,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAlsAAAJPCAYAAABLghH8AAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJzs3Xt0VPW9///nTgYIILlIuCmXQMAgcomxCAbUGOqNmGrB\neo2Klvas2p5+a9sT267zPU1X17enoqfW/lbP+X5bqqCp97uOipoQQDAaDcEQJBBwCPcQCOEaIJP9\n++PDJAESMpPMzswkr8darDKbPXs+k0k7r34+7/3+WLZtIyIiIiLOiAr1AERERER6MoUtEREREQcp\nbImIiIg4SGFLRERExEEKWyIiIiIOUtgSERERcZDfYcuyrH9YlrXXsqyvWh1LsCzrQ8uyKi3LWmZZ\nVlyrf/u1ZVmbLcv62rKsG4I9cBEREZFIEMjM1jPAjWcd+xXwsW3bKUAh8GsAy7ImAXcAlwI3A/9t\nWZbV9eGKiIiIRBa/w5Zt258AdWcdvhVYevrvS4HbTv/9O8CLtm032rbtATYDV3ZtqCIiIiKRp6s1\nW0Nt294LYNv2HmDo6eMXA9tbnbfz9DERERGRXiXYBfLa+0dERESkFVcXn7/Xsqxhtm3vtSxrOFBz\n+vhOYFSr80aePnYOy7IU0ERERCRi2LYdUB16oDNb1uk/Pm8DC07//QHgrVbH77Isq69lWWOB8cDn\n7V3Utu2w//Pb3/425GPQH302kfZHn014/9HnE75/9NmE75/O8Htmy7Ks54EMYLBlWdXAb4E/Aq9Y\nlvUQsA1zByK2bW+wLOtlYANwCnjY7uwIRURERCKY32HLtu172vmnb7dz/n8C/9mZQYmIiIj0FOog\n76eMjIxQD0Haoc8mfOmzCW/6fMKXPpuexQr16p5lWVphFBERkYhgWRa2wwXyIiIiIhIAhS0RERER\nBylsiYiIiDioq01NRURExA9JSUls27Yt1MOQ8xgzZgwejyfo11WBvIiISDc4XVgd6mHIefjzGalA\nXkRERCTMKGyJiIiIOEhhS0RERMRBClsiIiIiDlLYEhEREXGQwpaIiIiIgxS2RERExBGHDx/m97//\nPVdccQWxsbEMHDiQyy67jNzcXGpqas45f9u2bURFRfHTn/60w2vv3buXX/7yl0yZMoXY2Fji4uK4\n5JJLuPvuu3njjTeceDudpqamIiIiEnSbNm3ixhtvZPv27cybN4+FCxfSp08fiouLeeqpp3jmmWd4\n5513mDlzZsDXrq6uZvr06Rw5coR7772Xhx9+GICqqiqWL1/OkiVL+O53vxvst9RpClsiIiISVMeP\nHyc7O5vdu3fz7rvvctNNNzX/28KFC3n44YeZM2cOt912G+Xl5QwZMiSg6z/++OPU1tby1ltvccst\nt5zz723NmoWSlhFFREQkqBYvXszmzZt55JFHzghaPmlpafzhD3+gpqaGxx9/PODrV1VVAZCZmdnm\nvw8dOjTgazpJYUtERESC6tVXX8WyLH7wgx+0e86CBQvo06cPr732WsDXT05OxrZt/va3v3VlmN1G\nYUtERESCqqKigkGDBjFu3Lh2z+nfvz8TJ07E4/Fw7NixgK7/i1/8gri4OH7+858zZswYcnJyeOqp\npygtLe3q0B2hmi0REZFwlJ3dfa/1zjtBvdyhQ4cYMWJEh+fFxsYCUF9fz4ABA/y+/tixY1m3bh1P\nPPEEb7/9Ni+88AIvvPACtm0zZcoUnnnmGdLS0jo9/mDTzJaIiIgEVWxsLIcOHerwPN85cXFxAb/G\n6NGj+ctf/oLH42HXrl28/PLLfOc736G8vJzs7GwOHjwY8DWdopktERGRcBTk2abuNHnyZFatWsXW\nrVvbXUo8fvw4GzduJCkpKaBZrbYMGzaM+fPnM3/+fHJycnjhhRd47733uOeee7p03WDRzJaIiIgE\n1bx587Btm8WLF7d7ztKlSzl16hTz588P6mvPnDkT27bZuXNnUK/bFQpbIiIiElQLFy5k/Pjx/OlP\nf2LZsmXn/HtpaSm/+c1vGDZsGL/85S8Dvv6KFStoaGg457ht27z99ttYlsWkSZM6NXYnWLZth3YA\nlmWHegwiIiJOsyyL3vR9t3HjRm6++WZ27NjBvHnzyMjIwOVy8dlnn5Gfn09sbCzvvPMOM2bMaH7O\ntm3bGDt2LDNmzGizWanL5eLRRx8lOzub1atXk52dTVpaGnFxcezZs4fXXnuN0tJSMjMz+eijjwIe\nsz+f0elzrICuG+oPXmFLRER6g94WtsDsjfjUU0/x+uuvU1VVhdfrZcyYMWRnZ/OLX/zinOaj27Zt\nO2+7iJiYGI4ePcrnn3/OK6+8wsqVK6murubAgQMMHDiQSy+9lDvuuIOHH36YPn36BDxehS0REZEI\n1hvDVqRxKmypZktERETEQQpbIiIiIg5S2BIRERFxkMKWiIiIiIMUtkREREQcpLAlIiIi4iCFLRER\nEREHKWyJiIiIOEhhS0RERMRBClsiIiIiDlLYEhEREXGQwpaIiIiIgxS2RERERByksCUiIiKOOHz4\nML///e+54ooriI2NZeDAgVx22WXk5uZSU1Nzzvnbtm0jKiqq+U90dDTx8fFMnDiRu+++m1dffRWv\n13vO81asWHHG86Kiohg0aBBpaWk88cQTNDY2dsfbbZdl23ZoB2BZdqjHICIi4jTLsuhN33ebNm3i\nxhtvZPv27cybN4/rrruOPn36UFxczHPPPUdsbCzvvPMOM2fObH7Otm3bGDt2LDfccAP3338/AEeO\nHGHr1q28++67bNiwgbS0NN544w1GjRrV/LwVK1Zw3XXXcc899zB37lxs22bPnj08++yzlJeXc//9\n97NkyZIOx+zPZ3T6HCuQn4XCloiISDfoTWHr+PHjpKamsm3bNt58801uuummM/69tLSUOXPm0K9f\nP8rLyxkyZAjQErZ+8pOf8Je//OWc6/7lL3/hZz/7GVOmTGHt2rVERZkFOl/YeuKJJ/j5z3/efP6x\nY8eYOHEiO3fuZNeuXQwbNuy843YqbGkZUURERIJq8eLFbN68mUceeeScoAWQlpbGH/7wB2pqanj8\n8cf9vu5Pf/pT7r33XtavX8+LL77Y4fkDBgxonjnzeDx+v06wKWyJiIhIUL366qtYlsUPfvCDds9Z\nsGABffr04bXXXgvo2gsXLsS2bdxut1/nV1VVAXDRRRcF9DrB5ArZK4uIiEiPVFFRwaBBgxg3bly7\n5/Tv35+JEyeyfv16jh07xoABA/y69tSpUwFTE3a2Y8eOsX///uaarf/5n/9h3bp13HbbbWfUeHU3\nhS0REZEwlJ3dfa/1zjvBvd6hQ4cYMWJEh+fFxsYCUF9f73fY8j3n0KFD5/zbb3/7W/7jP/7jjGM/\n//nP+cMf/uDXtZ2iZUQREREJqtjY2DbD0Nl858TFxfl9bd9zfKGrtR/+8Id8/PHHvP/++yxatIjB\ngwfzyiuvsG/fPr+v7wTNbImIiIShYM82dafJkyezatUqtm7d2u5S4vHjx9m4cSNJSUl+z2oBfPXV\nVwCkpKSc828TJkwgMzMTgBtvvJFZs2Yxe/ZsFi5cyPvvv9+JdxIcmtkSERGRoJo3bx62bbN48eJ2\nz1m6dCmnTp1i/vz5AV3773//O5Zlccstt3R47lVXXcV9993Hhx9+SGFhYUCvE0wKWyIiIhJUCxcu\nZPz48fzpT39i2bJl5/x7aWkpv/nNbxg2bBi//OUv/b7uU089xfPPP8+0adO44447/HrO//7f/5uo\nqCh+97vf+f06waZlRBEREQmqAQMG8Pbbb3PzzTdzyy23MG/ePDIyMnC5XHz22Wfk5+cTGxvLm2++\nydChQ895/qZNm/jnP/8JmDsMt2zZwrvvvsvXX3/N9OnTeeONN5obmnYkOTmZu+66i+eff56ioiIy\nMjKC+Vb9og7yIiIi3aA3dZD3OXz4ME899RSvv/46VVVVeL1exowZQ3Z2Nr/4xS/OCVrbtm07o8bL\nsiwuuOAChg8fTlpaGrfffju33XbbOUFrxYoVZGZm8sQTT/DII4+cM46NGzcyZcoUrr766vMuJ2q7\nHhERkQjWG8NWpNF2PSIiIiIRSGFLRERExEEKWyIiIiIOUtgSERERcZDCloiIiIiDFLZEREREHKSw\nJSIiIuIghS0RERERBylsiYiIiDhIYUtERETEQdqIWkREpBuMGTMGywpolxfpZmPGjHHkutobUURE\nRMRP2htRREREJMwobImIiIg4SGFLRERExEEKWyIiIiIOUtgSERERcZDCloiIiIiDFLZEREREHKSw\nJSIiIuIghS0RERERBylsiYiIiDhIYUtERETEQQpbIiIiIg5S2BIRERFxkMKWiIiIiIMUtkREREQc\npLAlIiIi4iCFLREREREHuUI9ABERkUjj8XjweDwAlJWVkZqaCkBSUhJJSUmhG5iEpaCELcuyfg3k\nAF6gHHgQGAi8BIwBPMAdtm3XB+P1REREQql1qCoqKiIjIyOk45Hw1uVlRMuyxgA/AC63bXsqJsDd\nDfwK+Ni27RSgEPh1V19LREREJNIEo2brEHASGGhZlgvoD+wEbgWWnj5nKXBbEF5LREREJKJ0OWzZ\ntl0H/BdQjQlZ9bZtfwwMs2177+lz9gBDu/paIiIiIpEmGMuI44BHMLVZF2FmuO4F7LNOPfuxiIiI\nSI8XjAL5bwGrbds+AGBZ1htAOrDXsqxhtm3vtSxrOFDT3gXy8vKa/56RkaFCQxEREQkLRUVFFBUV\ndekalm13bcLJsqxpQD4wHTgBPAOUAKOBA7ZtP2ZZ1qNAgm3bv2rj+XZXxyAiIhIqeXl5Z0waSM9m\nWRa2bVuBPKfLM1u2ba+zLOtZ4EtM64e1wN+AQcDLlmU9BGwD7ujqa4mIiIhEmqD02bJt+3Hg8bMO\nHwC+HYzri4iIiEQqbdcjIiIi4iCFLREREREHKWyJiIiIOEhhS0RERMRBClsiIiIiDlLYEhEREXGQ\nwpaIiIiIgxS2RERERByksCUiIiLiIIUtEREREQcpbImIiIg4SGFLRERExEEKWyIiIiIOUtgSEZFe\nx+PxdPkadXV15OfnU1xcTH5+PnV1dV0fGMEZm4QXhS0REel1uhpo6urqWLRoEZWVlTQ0NFBZWcmi\nRYuCErgUtnoehS0REZEAud1uXC4X0dHRAERHR+NyuXC73SEemYQjV6gHICIi0t3KysooKirq9POL\ni4tpaGgAzEyU71olJSVUVVV1aWzx8fFkZGR06RoSXhS2RESk10lNTe1SoMnPz6eyspLo6GiKiorI\nyMjA6/WSkpJCTk5Ol8bWlRAo4UnLiCIiIgHKysqisbERr9cLgNfrpbGxkaysrBCPTMKRwpaIiPQ6\nSUlJXXp+QkICubm5pKSkEBMTQ0pKCrm5uSQkJIR8bBJ+LNu2QzsAy7JDPQYREZHOysvLIy8vL9TD\nkG5iWRa2bVuBPEczWyIiIiIOUtgSERERcZDCloiIiIiDFLZEREREHKSwJSIiIuIghS0RERERByls\niYiIiDhIYUtERETEQQpbIiIiIg5S2BIRERFxkMKWiIiIiIMUtkREREQcpLAlIiIi4iCFLREREREH\nKWyJiIiIOEhhS0RERMRBClsiIiIiDlLYEhEREXGQwpaIiIiIg1yhHoCIiEik8Xg8eDweAIqLiykq\nKgIgKSmJpKSkkI1LwpNmtkREREQcpJktERGRALWewSoqKiIjIyOk45HwppktEREREQcpbImIiIg4\nSGFLRERExEGq2RIRkbDS+k6/srIyUlNTAd3pJ5FLYUtERMKKis+lp9EyooiIiIiDFLZEREREHKSw\nJSIiIuIghS0REQkpXzF8pKmrqyM/P5/i4mLy8/Opq6sL9ZA6LVI/g0ihAnkREem0YNw56PF4Iu4u\nw7q6OhYtWoTL5aKhoYHKykoqKirIzc0lISEh1MMLWCR+BpFEM1siIhKQ1rMgSUlJTJs2jR07dvDB\nBx+wY8cOpk2bdsYXd0+cNXG73bhcLqKjowGIjo7G5XLhdrtDPDIJR5rZEhGRgLSeBfFnhqejWZOy\nsjKKiora/LeioiLy8vKC+waCoLi4mIaGBsC8P9/4S0pKqKqqCuHIOic+Pl4tNhyksCUiIp3W1gyP\n73hOTo5f10hNTW33iz4vLy8sw1Z+fj6VlZVER0c39wLzer2kpKT4/b7DSXthV4JDYUtERALSeibK\nnxmenjhrkpWVRUVFRfNjr9dLY2MjWVlZIRyVhCuFLRERCUjrmSh/Zng6mjWJxMLshIQEcnNzcbvd\nlJSUkJKSQlZWVkQWx0NkfgaRRAXyIiLSaVlZWTQ2NuL1egH/ZnjOLphv64u+M20VursQPyEhgZyc\nHGbOnElOTk7EBi1Q2HKawpaIiASk9RdzfX096enpNDU1kZCQQFNTE+np6dTX17d5PnQcinxF95WV\nlc1F94sWLeowcPXEux6lZ9AyooiIBKR1ePL108rOzvbrfH8Eo+heJJwobImISLc6X6sH6HxbhZ5Y\niC89g8KWiIh0q/O1eoDOt1VQ+wIJV6rZEhGRsNKZonuRcKawJSIi3aqjGi5fW4WUlBRiYmJISUnx\na89B3VEn4UrLiCIi0q38CUW+tgpVVVV+F8UrbEm40syWiIiIiIM0syUiIhIgj8fT3NcrPj6+uTjf\n1wpDpDWFLRERkQC1DlVqNyEd0TKiiIiIiIMUtkREREQcpLAlIiIi4iCFLREREREHKWyJiIiIOEhh\nS0RERMRBClsiIiIiDlKfLRERCStqGCo9jWXbdmgHYFl2qMcgIiIi4g/LsrBt2wrkOVpGFBEREXGQ\nwpaIiIiIgxS2RERERByksCUiIiLiIIUtEREREQcpbImIiIg4KChhy7KsOMuyXrEs62vLsiosy5ph\nWVaCZVkfWpZVaVnWMsuy4oLxWiIiIiKRJFgzW08B79m2fSkwDdgI/Ar42LbtFKAQ+HWQXktEREQk\nYnS5qallWbHAWtu2k886vhG41rbtvZZlDQeKbNue2Mbz1dRUREREIkKompqOBWoty3rGsqxSy7L+\nZlnWAGCYbdt7AWzb3gMMDcJriYiIiESUYIQtF5AG/NW27TTgKGYJ8ezpKk1fiYiISK8TjI2odwDb\nbdv+4vTj1zBha69lWcNaLSPWtHeBvLy85r9nZGSQkZERhGGJiIiIdE1RUVHzZuidFZSNqC3LWgH8\nwLbtTZZl/RYYcPqfDti2/ZhlWY8CCbZt/6qN56pmS0RERCJCZ2q2ghW2pgGLgT7AVuBBIBp4GRgF\nbAPusG37YBvPVdgSERGRiBCysNUVClsiIiISKUJ1N6KIiIiItENhS0RERMRBClsiIiIiDlLYEhER\nEXGQwpaIiIiIg4LR1FRERCRoPB4PHo8HgLKyMlJTUwFISkoiKSkpdAMT6STNbImISNjweDwkJSUx\nbdo0duzYwQcffMCOHTuYNm0aSUlJzSFMJJIobImISNjweDzU1dWxaNEiKisraWhooLKykkWLFlFX\nV6ewJRFJYUtERMKK2+3G5XIRHR0NQHR0NC6XC7fbHeKRiXSOarZERCRslJWV8cEHH9DQ0ACYmS7f\nJsAlJSXcdNNNZGRkhG6AIp2gsCUiImEjNTWVxMREKisriY6OpqioiIyMDLxeLykpKYwcOTLUQ5QI\ncfw4vPYaTJ8OKSmhHYuWEUVEerhIq3PKysriyJEjlJWVsWPHDsrKyjhy5AhZWVmhHlpAIu3n3lM0\nNcGyZfDDH8JLL8HixRDqLZg1syUi0sP57vCLBL5x2md9O/oeR8r7gMj6ufcUpaXw9NOwbZt5nJIC\nDz0EVkDbRgefwpaIiISNpKQknnzySfbv3098fDwA8fHx7N+/nyVLlvDII4+EeIQSjrZtMyGrtNQ8\nHjoUFiyA2bNDH7RAYUtEpMcrKytrLjKPBMuXL6e+vh6AAwcOUFZWBkB1dXXz8UgQHx+vYn6H1dXB\nP/8JH35olgoHDIA774RbboG+fUM9uhYKWyIiPVxqamqnvvT9WQbzeDyUl5dTXFzMhg0bmDRpEjNn\nzmTKlCl+Pbetc+Li4igtLSU6Oprhw4czfPhwvF4vaWlpETWzFUkBN9KcOAFvvQWvvAINDRAVBXPn\nwt13Q1xcqEd3LoUtERFpkz9hKy4ujjVr1uByuairqyMqKoo1a9Ywe/bsTl9/wYIF1NTU4HK5SEpK\nwuv10tjYyIIFCzr3RqTHsG0oKoJnn4XaWnPsyivhwQchnG9UVdgSEenhnCzSbqsBqe94Tk5Op66Z\nkJBAbm4ubreb2tpaEhMTycrKIiEhIWjj7g4qjg+u9evhH/+AqirzeNw4+P73YerU0I7LH9bZd3x0\n+wAsyw71GERE5Fx//vOfOXjw4HnPKS4uPqMBqS9gxMTEMHPmzPM+Nz4+np/97GdBGav0XDt3wjPP\nwGefmccXXgj33w/XXWeWD7ubZVnYth1Q2b1mtkREpE3+1Hrl5+e324C0o5kt1TTJ+Rw+DC+8AO+9\nB14vxMTA/Plw223m75FEYUtERDotKyuLioqK5se++qpIa0Aq4ePUKXj3XdOQ9OhR07rh+ushJ8fM\nakUihS0REWmTPzVH9fX1pKenU1xcTEJCAk1NTaSnp1NfX99hjZVqmqQ124bVq2HpUtizxxxLTTV1\nWZH+q6KaLREREQmpykpT/P711+bxqFEmZKWlhUdT0tZUsyUiIiIRo6YGliyBVavM47g4uPdeuOEG\nOH1ja4+gsCUiIiLd6uhR05D07bdNjVafPqbw/fbbTRf4nkZhS0RERLpFYyMsWwbPPw+HDpljGRmm\nlcOQISEdmqMUtkRERMRRtg1ffGE2i96xwxybNAkWLoQJE0I7tu6gsCUiEqY8Hg8ejwcwm0mnpqYC\n5i4+3cknkWLrVlP8/tVX5vGIEbBgAVx1VfgVvztFYUtEJEy1DlW+hqEikWL/fnjuOSgsNDNbF1wA\nd90FWVng6mXpo5e9XREREXFSQwO8/rr5c+KECVZZWXDnnTBoUKhHFxoKWyIiItJlTU1QUAD5+XDg\ngDmWnm6WDEeMCOnQQk5hS0RERLqkrMzUZZ0uMWTCBNOU9LLLQjqssKGwJSIiIp2yfbu5w/CLL8zj\nxER44AG49treU/zuD4UtERERCUh9vemV9cEHZvmwf3/43vfg1luhb99Qjy78KGyJiIQpj8dDXFwc\nbreb4uJi8vPzycrKor6+Xq0fJCROnjRd319+GY4fN7NXN99sttiJiwv16MKXwpaISJgqLy9nzZo1\nuFwuGhoaqKyspKKigvT0dIUt6Va2DStXwtKlsG+fOXbFFfDQQzB6dGjHFgkUtkREwlRxcTEul4vo\n0zvy+v6zuLiY7OzsUA5NepENG2DxYti82TxOSjLF76d77IofFLZERMLUhg0bqKurA8ySYlFREQAJ\nCQkhHJX0Frt3w5IlsGaNeZyQAPfdB3PmQFRUSIcWcRS2RETC1KRJk4iKiiI6Orq5g7zX66WpqSnU\nQ5Me7MgRePFFcLvNxtF9+8K8eTB/PsTEhHp0kUlhS0QkTM2cOZM1vmkFwOv10tjYSHp6eghHJT1V\nY6MJWC++aAKXZZlZrPvug8GDQz26yKawJSISpqZMmcLs2bNxu92UlJSQkpLSfDeiSLDYNhQXwzPP\nmKVDgKlTTV3WuHGhHVtPYdm2HdoBWJYd6jGIiIS7vLw88vLyQj0M6WE2bzad3ysqzOOLLzZ3GE6f\nrqak7bEsC9u2A/rpaGZLRESkl9m3D559Fk7fc0FsLNxzD9x4o9k4WoJLP1IREZFe4vhxePVVePNN\n06DU5TJd37/3PRg4MNSj67kUtkRERHo4rxc++gjy881WOwDXXGP2MRw6NLRj6w0UtkRERHoo24bS\nUrNZdHW1OTZxIixcCCkpoR1bb6KwJSIiEiCPx4PH4wGgrKyM1NPt1JOSksJmKyWPxxS/l5WZx8OG\nwYIFMGuWit+7m8KWiEiYav2FHh8f39xBPpy+0Hur1p+Br+FsuKirM8uFH31kZrYGDoQ774RbboE+\nfUI9ut5JYUtEJEy1/kIPpy9zCU8nTsAbb8Brr0FDA0RHw9y5cPfdMGhQqEfXuylsiYiIRDDbhsJC\neO452L/fHJsxAx580PTNktBT2BIREYlQX31l6rK2bjWPk5NN5/cpU0I7LjmTwpaIiEiE2bnTbK/z\n2WfmcWIi3H8/ZGSo+D0cKWyJiIhEiEOH4IUX4P33Te+smBi4/Xa47Tbo1y/Uo5P2KGyJiEinRUIL\nhJ7g1Cl45x14+WU4etTMXt14I9x7LyQkhHp00hGFLRER6bRwboHQE9g2rF5tlgxrasyxtDRT/K4s\nGzkUtkRERMLQxo2m+H3jRvN49Gh46CG44orQjksCp7AlIiISRvbuhSVL4JNPzOO4OMjJgeuvN72z\nJPIobImIiISBo0dNTdbbb0NjI/Ttawrfb78d+vcP9eikKxS2REREQqixET74AJ5/Hg4fNseuu860\nckhMDO3YJDgUtkRERELAtuHzz03x+86d5tjkyaYua8KE0I5NgkthS0REpJtt2WKK38vLzeMRI0zI\nmjFDTUl7IoUtERERTM+wQHqD1dXV4Xa7KS4uJj8/n6ysLBI6aHpVW2v2MFy+3MxsRUUd4fvfv4C5\nc8Glb+QeSx+tiIgIgYWturo6Fi1ahMvloqGhgcrKSioqKsjNzW0zcDU0wKuvwhtvwMmTJljdcgtc\ndNFabr756iC/Ewk3UaEegIiISKRxu924XC6iT/diiI6OxuVy4Xa7zzivqQk+/BB+8AN46SUTtGbN\ngv/+b7NhdP/+3lAMX7qZZrZEREQw2w0VFRX5dW5xcTENDQ2AmRHzPa+kpISqqioAdu8eztq1l1Nf\nHw/A4MH7ufzyUhoaavl//89cJz4+Xl33ewGFLRERESA1NdXv4JOfn09lZSXR0dHN2xR5vV5SUlK4\n5pocnn4adu2CYcNgyhR44AG4+mqwrPQzruNvuJNOaGiATz81dx9MnBjSoShsiYiIBCgrK4uKiorm\nx16vlyNHXOzePY+f/MQUvw8YAHfcAdnZpkGpdAPbNrd4FhaaTSUbGsy67a9+FdJhKWyJiIhAQHci\nJiQkkJvJxz2zAAAgAElEQVSbi9vt5rPP1nLw4PVUV8+kqqovUVFw881wzz1mq51gvaacx65dJmAt\nX96yYzeYGa0w2ExSYUtEpBfxeDx4PB7A1CilpqYC5ku/t3/xB/r+4+MTuPjiHPbtu4x+/S4HYPp0\nePBBGDXKmdeUVo4ehVWroKCgZbdugCFDTAv+zEy4+OLQja8VhS0RkQ4E2n8pnCUlJREXF4fb7eaD\nDz4gMTHRr/5QkaK7Pqv16+Hpp2HzZjh2bADjxpm7C6dOdfylezevF0pLzSzWZ5/BqVPmeEyMWS7M\nzDRFcmHWGVZhS0SkAz0pbAXaHyrSOP1Z7doFS5aYumuACy+EGTM+48knU4hSMyXneDxmBquoCA4e\nNMcsC6ZNMwErPd0ErjClsCUi0ou01R/KdzwnJyeUQwtrhw/Diy+C220mV/r1g/nz4bvfhT/+cauC\nlhMOHoQVK0zI+uabluMXX2wC1nXXmSXDCKCwJSLSgUD6L4U7f/pDdVZRURF5eXldHGHXBLtv1alT\nJmC9+KIpEbIsuP56yMkxs1oSZCdPmt25Cwvhyy9NV1iACy6Aa64xIeuSS8JumbAjClsiElQ9sQA7\nkP5L4e58/aG6OrOVl5cX8rAVrFBs22ap8JlnYM8ecyw11WwWPXZsUF5CfGwbNm0yM1grV5pUCxAV\nBVdeaQLWlVdCnz6hHWcXKGyJSFC1DlW+L3MJH231h2psbCQrKyuEowovmzbBP/4BGzaYx6NGmZB1\nxRURN6ES3vbtM60aCgth586W4+PGmYCVkdFx74wIobAlItKBSJ2Ra0vr/lAlJSWkpKT0qLsRu/JZ\n1dTA0qVmcgXM9/y998INN8Dp0jbpqoYG02y0sNA0H7Vtczw+vqVdQw/675uPwpaISAd6UtgCE7hy\ncnKoqqrqcUXxnfmsjh2DV16Bt94yNVp9+sCtt8L3vme6wEsX+bq6FxTAmjUmcIH5Qc+YAXPmwOWX\n9+hEq7AlItKLtK6pi4+Pb65xiuSaus7yemHZMnj+eaivN8euvRbuvx+GDg3t2HqEnTvNDFZhIdTW\nthy/9FIzg3X11TBwYOjG140UtkREepHWoaq31tPZNnzxhSl+377dHJs0yTQlveSS0I4t4h050tLV\nvbKy5fjQoSZgZWaajaF7GYUtERHpNb75xnR+Lyszj0eMgAUL4KqrVPzeaY2NsHatCViffWYeg2ky\nOnu2WSa87LJe/QNW2BIRkS5Zt24d5eXlFBcXk5+fT1ZWFvX19WG1LHngAOTnw8cfm5mtCy6Au+6C\nrCxw6Zuwc7ZuNUuERUUt67CWZXpkZGaaBBvGXd27k37FRESk0+rq6njsscdITk4+Y/uf9PT0sAhb\nDQ3wxhvw2mtw4oQJVllZcOedMGhQ56/ba2vfDh404aqgwGyh4zNypJnBysiAxMQQDS58KWyJiEin\nud1uoqOjz9n+p7i4mOzs7JCNq6nJTLo895yZ1QKzfd6CBcEpGepVtW++ru4FBWYTaF9X90GDWrq6\nT5jQq5cJO6KwJSLSiwS7w39tbS379u1j+/btZ2z/E8q+XV99ZZqSbt1qHk+YYIrfL7ssZEOKPLYN\nGzeaxLpqVUtX9+ho064hMxOmT4/oru7dKWhhy7KsKOALYIdt29+xLCsBeAkYA3iAO2zbrg/W64mI\nSOCC3eE/MTGRIUOGkJycfMb2P02+2Y9utGOHKX4vKfGNDR54wLRz0KSLn2pqWrq679rVcjw52SwT\nXnNNj+nq3p2CObP1v4ANQOzpx78CPrZte5FlWY8Cvz59TEQk4ng8nh5Ti1NXV4fb7T6joL2zM1FZ\nWVm89957eL1eoGX7n/T09GAO+bzq602vrNdfP0ZMzAD69zcNSW+9Ffr27bZhRK7jx02z0YIC03zU\n58ILTQ1WZiaMGROy4fUEQQlblmWNBOYC/wf4+enDtwLXnv77UqAIhS0RiVA9JWzV1dWxaNEiXC7X\nGQXtubm5nQpcCQkJPProo5SXl5+x/U99vfMLGSdPwttvm+7vx47BsWPH+O53B3DPPWb3FzmPpiaz\n3lpYaILWiRPmeN++MHOmCVipqT26q3t3CtbM1pPAvwGt5xaH2ba9F8C27T2WZakfr4hIiLndblwu\n1zkF7W63u9Nb90ybNo1p06adsf2PkzVbtm3KiJYuNateYDaJvvTSzdx5p+6EO6+dO80M1vLlZ3Z1\nv+wyE7Bmzeo1Xd27U5fDlmVZWcBe27bLLMvKOM+pdldfS0R6p2AXdXdGWVlZc/F3JCsuLqbh9N50\nrQvaS0pKqKqq6tK1i4qKyMvL6+IIz2/fvkRKS9M4cGAwAHFxB0lLWwvsYffueOAqR18/Ih0+3NLV\nfdOmluPDhpmAdd11vbKre3cKxszWLOA7lmXNBfoDgyzLeg7YY1nWMNu291qWNRyoae8Crf/LmZGR\n0fNvoxWRgAS7qLszUlNTe8T/NuXn51NZWUl0dPQZBe0pKSld3pQ6Ly/PsbC1e7eZyaquNoXvEybA\nffeZmu2oqBkAPSIMB01jI3z5pVkm/Pzzlq7u/fu3dHWfNEl3DvihqKioy79bXQ5btm3/BvgNgGVZ\n1wK/sG37PsuyFgELgMeAB4C32ruG0/9PSEREjKysLCoqKpof+wras7Kyun0s56uD8xXx79xZz+bN\naezZcwWW1Ze+fWHePJg/v6U5eU+pp+sy2zb7ERUUwIoVZ3Z1T0szs1gzZ0K/fqEdZ4Q5exLod7/7\nXcDXcLLP1h+Bly3LegjYBtzh4GuJSJjr6hdiMO+g68x4esqXeUJCArm5ubjd7jMK2kPRF6u9z2Dd\nunX8538+zt69aWzceDnR0TFAOXffPZxHH734nAblvuv0lM8oYAcOmHBVUADbtrUcHzWqpav74MEh\nG54EOWzZtr0CWHH67weAbwfz+iLSIpz/33xbwagr4w32HXTQ8vPzd1zh+rPujISEBHJycs4oaA8X\ntg1vvLGbTZv+F0ePxnLy5H4mTTrC1KklTJs2lMTE9sfbkz6jDp08CcXFZpmwtNT84MB0dfe1a0hO\n1jJhmFAHeZEIFa5hq71g1JW+S07cQecTrj9Hp4TTnn5n33Swf/+FrF17OZs2xeP1nqJfv0r69n2a\nqKiDVFTAli0xbRbxx8fH94h6ug75uroXFMAnn7R0dXe5TDf3zEz41re0s3YY0iciIkHVXjB66aWX\n+PLLLzt1TSfuoOs1X9BnCac9/Xw3HdTWwrPPmnruhAQYP/4bRo/+hEsu+YaVKw92WMTf4wvja2rM\nDFZhoblTwGfCBBOwrrkGYmPbf76EnMKWSIQK11YE7QWjo0ePMn78+E5d07ZtDhw4QFRUFIcOHeLg\nwYM0NTUxogu3q5eVlZGXl9drQ1c4OHEiiueegzffNKtiLhdkZ8OkSUf4y1/e55tvoomJiWHLli14\nvV5uv/32UA+5+xw7BqtXm4C1fn3L8QsvNAErM9PUZElEUNgSiVDh2oqgvdYCTU1Nnb7zuPXS5HPP\nPUd2djaNjY1dqtnyjS0cA2tP5/XCRx/B//2/aZze5YfZs2HBAtP6CaaQkvJX3G43tbW1JCYmnreI\nv8csAzc1wbp1LV3dT540x/v2hfR0E7CmTYOoqNCOUwKmsCUiQdVea4G77rqr09dsfQedy+UKyh10\nPeYLOsJ8+aXZLLq6GuACUlJg4UKYOPHM83xF/P6I+M9y+3YTsJYvh/37W45PntzS1X3AgNCNT7pM\nYUskQoXrF4xTrQV8X76LFy8Oyh10vp9fuP4cexqPx4SstWvN42HDzEzWrFm99Ia5w4dNu4bCQti8\nueX48OGmXcN11/mm+aQHUNgSiVDhHBLCubXA2cL559gT1NVBfr5ZNrRts+3enXfCLbdAnz6hHl03\na2yEL74wAaukpKWr+4ABcPXVZhbr0kt7afrs2RS2RCTstW5XEBMTE9J2BeKfEyfgjTfgtdegoQGi\no2HuXLjrrl5245xtQ1WVWSJcsQIOHTLHLcvsnu3r6t63b2jHKY5S2BKRsBcOeyOKf2zb5Ipnn20p\nP5oxAx58EC6+OLRj61YHDpgfREGBqcnyGTPGLBNee625s1B6BYUtEREJir17h/LII7Bli3mcnAzf\n/z5MmRLacXWbEydMV/eCAigra+nqHhvb0tV93DgtE/ZCClsiItIlO3fCM89AYeEcLrnEbMN3//2m\nxrvH5wrbhg0bWrq6Hz9ujrtccOWVJmBdcYW6uvdy+vRFJKha11cVFxervqoHO3QIXngB3n/f9M5y\nuRrJyYHbboN+/UI9Ooft2dPS1X3v3pbjl1xilgmvvtrsUyiCwpaIBFnrUJWXl6f6qh7o1Cl49114\n6SWzPZ9lwY03wvbt73DnnZeFenjOOXrUdHUvKDCzWT6JiWYaLzMTRo4M3fgkbClsiYiIX2zbZI0l\nS1omcy6/HB56CJKSIC+vIZTDc0ZTk2kOVlho6rF8Xd379Wvp6j51qrq6y3kpbImISIc2boR//MP8\nJ8Do0SZkXXFFaMflmOpqM4NVVGTuLPSZMsUsE6anQ//+IRueRBaFLRERadfevbB0KaxaZR7HxUFO\nDlx/vemd1aPU18PKlSZk+W6pBBgxoqWr+9ChoRufRCyFLREROcfRo/Dyy/D226bRed++pvD99tt7\n2IROY6Pp5l5QYLq7+3bGHjiwpav7xIm94LZKcZLCloiINGtshA8+gOefN9v3gZnQue8+GDIktGML\nGl9X94ICM5Ple6NRUTB9uglYV16pru4SNApbIiKCbZsJnqefNn2zACZPNnVZEyaEdmxBU1trarAK\nC8/s6p6U1NLVvYsbpou0RWFLRKSX27rVFL9/9ZV5PGKECVkzZnS8eta6r1p8fHz49VVraGjp6r5u\nXUtX97i4M7u6izhIYUtEpBdpHY5Wr97I7t3fZu3aeGJiBjB06ADuvttsGO1vw/PWoSpseqrZNlRU\ntHR1bzjdksLlMps+Z2aanhXq6i7dRL9pIiK9SFJSEsOHJ/Haa/DPfw4lOXk8Q4fCLbfAnXfCBReE\neoRdsHt3S1f3mpqW4ykpLV3dI/oNSqRS2BKRoKurq2PJkiXs2LGD/Px8pkyZwrRp00I9rIjg8Xgc\nW35rajKTPU8/fYKNG/dRW3uIyZPX8rvfjeXSS+MdeU0fx97X0aOmL0VhIXz9dcvxxEQzg5WZCRdf\nHPzXFQmAwpaItKmzX47r1q3jscceo7a2FoBly5axZMkS/uu//qtLgauurg63201xcTH5+flkZWWR\n0IVi5tbLaWVlZaSmpgKdrzUKVphwKpSUlZm6rM2bT1JVtZXExAOMGfM3UlJG8+yzjeTm5nbp59mR\noL4vr/fMru6nTpnjMTEwa5YJWFOmqF2DhA2FLRFpU2e/HMvLy0lOTqZv376cPHmS5OTk5uOdDVt1\ndXUsWrQIl8tFQ0MDlZWVVFRUdCkgtA5VRUVFXa43ciIkBSMQbt9u7jD84gvz+PjxatLTVzN27A5W\nrNhGdPRYANxuNzk5OUEdf9B5PCZgFRVBXZ05ZllmuxxfV/eYmFCOUKRNClsiElS1tbVER0ezZ88e\nPB4PRUVFxMTENM90+QQSTtxuNy6Xi+jTLct9/9nVgBDs2bKu8ng8lJWVNd/R11pRUREHDx70+1oN\nDf0oL5/Kli3J2LaFy9XIZZdVcOLEs1RXH6G6mubPB6CkpISqqqogvZNzxcfHdy7Q1te3tGvYurXl\n+MUXmxms667rQQ3ApKeybN9tsKEagGXZoR6DiJzrz3/+c0Bf7j6ffvopu3fv5tixYxw4cIDRo0dz\n+PBhxo8fz1VXXdV8XiBhq7i4mIbTd5S1fl5MTAwzZ84MeIwAx48fZ/Xq1URFRVFdXc3o0aNpampi\n1qxZ9O9Ei/T4+Hh+9rOfdWosPr7g01YoycvLIy8vr8NrnDwJb70Fr7wCx4+bPp033QT33GO6HeTn\n51NZWUl0dHTzjJ7X6yUlJcXRma2AZg9PnYLPPzcB68svW7q6X3ABXHONCVmXXKJlQgkJy7KwbTug\nXz7NbIlIm1JTUzs1E+Fb8tuxYwerVq0iOzubLVu28Ne//vWMWaNAvnydCAj5+fn07ds3aNdsazaq\nO9k2rFgBzz4L+/aZY9Onw4MPwqhRLedlZWVRUVHR/Njr9dLY2EhWVlY3j/gstg2bNpmAtXIlHDli\njvu6us+ZY7q69+kT2nGKdILCloh0WevZpoSEBH74wx/y7//+7xw5coTa2lp+9KMfnbM8V1hY6HdA\nOXsWqrCwkKamJk6ePOnXbE9bzp4t6+pyWutrdFZ8fDy33XZbwM+rqPAVv5vHY8fC978PbZXIJSQk\nkJubi9vtpqSkhJSUlG5ZPm13FrO2tqVdg691PZhGo5mZpqt7vLN3Soo4TWFLRNoUSLF367BVV1fH\n3/72N8aNG8cFF1xAYmIi7733HpMnTz7jC/3CCy8MaNnNV1+Vn59PTk5OlwNCsGfLglEgX1RUFNA1\ndu2CJUvg00/N4wsvNHsYZmaaCaH2JCQkkJOTQ1VVVbcVxZ/xvhoaYM0aE7C++qqlq3t8vOnqPmeO\n2UJHpIdQ2BKRNnU2ODhVzB7sgBDs5bTu3Jrm8GF46SVwu83G0f36wbx55k9334zn8XiIi4vD7XZT\nW1tLYmJicxA+I4DaNpSXm4C1enVLV/c+fcy+QHPmmK7up39fRHoShS0R6bLWd9C1Xp47ePBgu8tz\nTU1N3T3MM4RqOe18OgpsjY0mYL34oilpsiz49rchJwcGD/bvNYK9l2F5eTlr1qxpDtj79+9vbsvh\n8XhI6tPHBKzly1uKyQAuvdRMwc2era7u0uMpbIlIl7Uupm+9POfxeNpdngt1QTmEZjntfNoLO7YN\n27eP4uGHzY40YOqxHnoo8D2Ug72XYXFx8TkzmX1PnmTtH/7AhOpqOHas5eQhQ1q6ul90UZdfWyRS\nKGyJSFCF7d1uZwn2DI9TNm+GxYvhk09mc8klMHKkCVnf+lZ4dD5YtWoVu3btIsq2uezECa46fpxp\nJ07Qp08fTl10EUyaFHBX92B39xcJNYUtEemy1l+ArZfnXC5Xu8tzof7SDPYMT7DV1Jg2DitWmMf9\n+p3gRz+CG24AVxj9L/d3pkxhfGws46qr6edywcCB2LbNqUsv5VBGBvzoRwEXkgW7u79IqIXRf2VF\nJFK1Dk6+WYmRI0cyfvx4Ro4cybp1686ZlQh12ApXx46ZhqRvvWV6e/bpA7feCpWV7zB37tRQD884\neLC5q3tOSQl1Bw5gRUVxODaWqjFj2DxyJD/Jy2PdunXaPkcEhS0RCTLNSgTO4/EwalQSH34I//yn\n2aEGTLP0W289yOefv0tp6aoOtxRyahNrwLSmb93V/fQNDhcMH86gO++koKmJrS4XiUOG8JPTY1Sg\nFjEUtkREQsi2Ydmy/axbl8T27ebYpZeapqRDhwa2AXfQw5ZtQ2UlFBTAqlVw9Kg5Hh1turnPmcMF\n06dDnz58p42nJyUldXpMwdq30tEAKuInhS0RkRD55ht4+mkoKBjD4MEwfDgsWADp6aaOPD/fmZ5l\nHdq3z7RqaKur+5w5pqt7XJxfl+pM2PFt+eRvyPTn9VV0L6GksCUi0s0OHID8fPj4YzN51NCwn4ED\nP2b48E189FETH31kzlu+fDn1p9cUq6urmzcGX7lyZZtbCsXHx3d+2bahwTQbLSw0zUd9Xd0TEkxX\n98zMbuvq7kRjXC1vSygpbIlIUEVKS4X2OLns1NAAb74Jr74KJ06Y1bi5c2HkyBrmzr3rnPPHjx9P\nWVkZW7dupbq6mnHjxjFu3DhSU1PbDB0B9y7zdXUvKDDb5/i6uvftCzNnmoCVmnreru4d/bxaN7wF\n0+g2voO9Dt9//31qamo4evQoR44cYcOGDQwcOJChQ4d2uG/l2dfvUgAVCRKFLREJqnBvqdARJ8KW\nbZsJo2efNbNaYLLMggVw8cVQVORt83mzZs3ir3/9K9HR0TQ1NVFTU8Pu3bv58Y9/3LUB7dzZsvlz\nbW3L8UsvNcuEs2fDwIF+Xaqjn1frhrfg36zSiBEjWLp0KS6Xi+rqakaPHk1jYyMPPPAA//Iv/3Le\n5559/XBoniuisCUi4qCvvoJ//AO2bjWPx483xe+TJ7ec015YWb16NbNmzWLLli1ERUUxZMgQkpOT\nWb16NWPHjj3n/POGxCNHYOVKE7AqK1uODx3a0tV9xIjA36ADbNvGOqv5qWVZ2L6lTZEIo7AlItLK\n2ctenXXoUCxr117Orl1mW5oBA44xdeo6YmM9vPqqWUrsSOt9Jpuamqirq+OLL75g/fr1HS6nAUQ1\nNTF8927GffMNF+/cSdTpdg2NLhfVo0fzzdix1AwZYtrUb97cqffZ0TLd2QHQn59vcXExADU1NRw6\ndIiamhoGDx7Mm2++yZ49ewIaTyQsXUvPp7AlItLK2ctegTp0CJ5/HtauNfsrT50K3/ueaUzar9/l\nAV3L330mz7F1q5nBKioyTbv694cJE8yGinPmwFVXcVm/fp1+j611FJzODjv+/Hxbv2/fsqBf77uN\n8ShsSThQ2BIRCYKTJ+Gdd+Dll00XeMuCm26Ce+4xN/R1RkD7TNbVNXd15/QNCgCMGmWWCK+7DgYP\n7txAuplT+2uuW7eO8vJyiouLefLJJ1mwYEGneneJBEphS0SklUBnQmzb9PtcutTsZwiQlmY2ix4z\npmtj6XCfyZMn4bPPTMAqLW3u6s6gQaYXVmamKRJzcMfqQH9e/pzf+n2XlJS0u79mINevq6vjscce\nIzk5mYaGBkpLS6mpqelU7y6RQClsiYi0Ekh4+PprU/zuqzcfM8YUv18e2GrheSUkJJCTk8PixYvN\nEpptmxcuLDy3q/uMGWaZcPr0bt2t2rd050+zUH9/vr73XVVVFVBvrfau73a7iY6OPqN3l8vlcr5B\nrAgKWyIiAduzx8xkffKJeRwfDzk5cP31EBUVvNdp3bNsmGWx4T/+g/jSUuKPH2fAgAHmpPHjzQzW\ntddCbGzwXtxPkdIstLa2ln379rF9+3Y8Hg8xMTEkJSVR27r1hYhDFLZERPx05IipyXrnHWhsNL0/\nv/tdmD/f1KAHW9KwYSRVVUFhISN27SJl7VqzJDhypKnBysyE0aOD/8I9UGJiYnPrjKKiIoYPH47X\n6yUxMTHUQ5NeQGFLRMJeqLvSNzbC++/DCy/A4cPmWGYm3HcfBP27uqnJNOcqLDRd3U+cAMAbHQ3X\nXGOWCVNTgzuF1gtkZWXx3nvv4fWaBrLBKroX8YfClohDtPFt5LNt+PxzeOaZlv2YJ082dVnjxwf5\nxXbsMNvmFBWd2dX9sstgzhzeWL6cSf/2b0F+0a5r/XteXFwcttszJSQk8Oijj1JeXk5JSQlpaWm6\nG1G6jcKWiEMipZalK5zcR7C1pKQk4uLicLvdfPDBByQmJvp9d1pnVVWZ4vf16+HYsWNMmDCABx+E\nK68M4s19hw+3dHXftKnl+PDhLe0ahg8H4NTq1UF60eBq/Xuel5cX1r/n06ZNY9q0aVRVVfHII4+E\nejjSiyhsiUindVfYqqurY9GiRbhcLhoaGqisrKSiosKR2/Zra80ehsuXm8eDBsFVV20hN3dKcG7w\na2yEL780Aevzz81jMEVfs2ebZcJJk8CyzKzR6ZmiSNzUu7NCvWwsEmwKWyIS9nx9plrftu87Hqzb\n9o8fN1vovPmmaV/lckF2Ntx5J5SU7O9a0LJt09W9oABWrDBt5sFMkaWlmYA1Ywac1dU90jf17qze\n+r6l51LYEpFOC9Y+gh1pvUegx+Npfs2SkhK/9ghsy8GDBzl48CBNTRYVFRdx8GAWJ08OoE8fF+PG\n7SY1dR3V1Ud4/PGO9/9r14EDpgaroACqq1uOjx5tAlZGBlx4YafGLyKRQ2FLRDqtq/sI+qsre+Wd\nT2kpPP00HDlSyaRJKaSkmOL3Sy+dDFzffF5AgfLkSSgubunqbtvm+KBBJlxlZkJysqNd3UUkvChs\niQTR2TVMdXV1uN1uiouLyc/PP6eou7tqniJdsPfK27bNhKzSUvN4wICj5OaakqlOZSBfV/eCAtPp\n9Ngxc9zlMt3c58yBK67o1q7uZ2v9u+b7vaytrT3nZoOe+jupOjAJJYUtkSA6+wuto6LuSP9i666x\n19fXk56eTnFxMQkJCTQ1NZGenk59fX1ABfJ1dfDPf8KHH5p8NGCAqclau9bN1Ventfu8dt/n3r1m\nBquw0LSV95kwwQSsa64xM1phwPe7tm7dOh577DGio6PZt28fQ4YM4b333uPRRx9l2rRpEf872R7V\ngUkoKWyJOKQ7irpDrbu+lH1flNnZ2Z16/okT8NZb8Mor0NBg+oHOnQt33w1xcfDVV94OX7/ZsWOw\nerUJWOvXtxwfPLilq/uoUZ0aZ3coLy8nOTmZ6Ohotm/fTnJyMl6vl/LycqZNm9bh89sKY76Zsh07\ndrQ5g3u+54r0BgpbIkHUumDcn6LuThdei19s29SnP/tsS5/QK6+EBx80O974rakJ1q0zy4Sffmrq\nssDs15Oebmaxpk4N667uvt/Njn4vO/qdbGup3DeD29jYeN62HApb0lspbIkEUeuCcX+KurvjTr7e\nav1605TUd7PiuHGm+H3q1AAusn27CVjLl5s7C30mTzYBa9YsZzZFDIKzg43vd/PJJ5/k008/Zffu\n3Rw5cgSv18uIESO46qqreOSRR874nfQnHPWGGVyRrlLYEnFIsIu6QVsA+WPnTliyxNwQCKazwv33\nmxU+vyaeDh0yvbCWL4fNm1uOjxjR0tV92DAnhh5U7QWlb33rW/z9738nKiqKpqYm6urq2L9/P//6\nr//q1zXObvfReqbs4MGD523LoZlc6a0UtkSCqPUXU0JCArm5ubjdbkpKSkhJSTmnliXQgNQbtgDq\nrMOHzUbR770HXi/ExMD8+XDbbebv5xPV1GTSWUEBfPFFS1f3gQPh6qtNyJo4MaLbNfh+b7Zt28bc\nuVWaUJAAACAASURBVHPZsmULDQ0NTJgwgeTkZLZt28bVV1/d4e9kamrqGb+H+fn5lJWVsWXLFqKi\nokhISCA5OZnU1NTmmS1faNNMrvRWClsiQXT2F1VCQgI5OTlUVVW1uaSi2aiuO3UK3n0XXnoJjh41\neej66yEnp4N+obZt1hgLC/nuG2/A2rXmuGWZNg2+ru59+3bL+wi29hrOtp6JamxspK6uji+++IL1\n69f7PRPVesZr1qxZ/PWvfyU6Opqmpib27dvHnj17+PGPf9zm+SK9kcKWiEQk24Y1a8ySoa/rQmoq\nPPQQjB17nifu39/S1X37dgD6njwJSUlmBuvaa3tEV/f2Gs4G0iC2rbCWlJTUvJQNsHr1ambNmsXW\nrVuJiopi6NChjBs3jtWrVzP2rA9CgUt6K4UtkdNUDxU5KitN8fvXX5vHo0aZkHXFFe2s9J040bJM\nWFbW0tU9Lg4yMvigvJyUP/85opcJ/dXVWsKkpCTefPPNNu+6bWpq4sCBAxw4cOCMmTLfDJn+eyS9\nlcKWyGmqhwp/NTWwdCmsXGkex8XBvffCDTfA6ZvgWtg2VFSYfliffGJ2mgbTxf3KK80yYVoauFzU\n5eX1uKDVXrDxp5awo2u0d9etx+PRXbcibVDYEpGwd/SoaUj69tumRqtPH7j1Vvje90wX+DPs3m0C\n1vLlpsO7T0qKWSa8+uqw6erupPaCkm8Gd+TIkdx0002MHDmSdevWtTmD689MlBN33Yr0NApbIg7R\nXmxd5/XCsmVmi51Dh8yxa681rRyGDm114tGjpqt7QQFs2NByPDGxpV1DQF1Me65gbFvT3l23Lpcr\nKHfdivQ0ClsiDmnrSy0Yd2V1tLl1Z4Tb3WK2bTowPP007Nhhjk2aZJqSXnLJ6ZO8XlN/VVho6rF8\nXd379Tuzq3sPWx4MB+3ddbt48WLddSvSBoUtkW7U1VDjz+bWoRhXMG3daorfv/rKPB4xAhYsgKuu\nOp2btm0zM1grVpzZ1X3qVDOLlZ4etl3dRaR3UtgSiSA9eWuU/fvhuefMRJVtwwUXwF13QVYWuI7W\nwzsrzD9u2dLypIsuMjNYGRlnrSuKiIQPhS2RbtReo0l/+bO5dWeEchuVhgZ4/XXz58QJc7NgVhbc\nOe8UgzaWwB8LzZqi12ue4OvqPmeOKXrXMqGIhDmFLZFu1F6jSX8F0pAyEMG8Nd/ffmVNTWY1MD+/\nZTUw/SqbBddsZUT5h/CTVWYPHjCbGk6fbpYJr7wyLLu6q0+biLRHYUskgkTCbfb+9CsrKzN1Wb5G\n5BNGHuf7E1Zy2eY34bEdLSeOHWtmsK65Brp4E4DT1KdNRNqjsCXSjbo6wxFIQ8ruHJe/tm83dxh+\n8QXQ5CXRW8MDQ97j2h1vYe04s6s7c+Z0sO+OiEhkUNgS6UbBCDUdbW7dGU6Hrfp60ytr2TKbpoOH\n6V+/h+/1f4dbL/yEvgcbTaHWzJlmmTAtrY128BKuWi+fxsTEqJ+cSBsUtkTEMV5vNK+8Aq88e5zj\nO/Zj7a/l5rhPuffiIuL6HIOJE80M1uzZ5vZDiThaPhXpmMKWiASdbcPKZcdZ+fw1lL2wAQ4f4Vvx\nm3noko8ZNdqC624xs1gXXxzScdbV1bFkyRKKi4t58sknWbBgQZeXZEVEzqawJXJa6+WQ4uJiLYd0\nhtfLhte+ZvH/d5zNVRB12EvSkC18f+pKUrMuhsxHYMqUsGjX4GsQu2PHDhoaGigtLaWmpqbLDWJF\nRM6msCVyWutQlZeXp+WQQHg87H5tDUvyXazZlQRAQt8jZIx9l0f+nE3UrN9DTExox3iWthrEulyu\nHtEgtjtpD1CRjilsiUjnHDwIK1Zw5INPeHHNaNx7v0WjHU3fgS7mZZ1kfu4U/viPlUTNuS7UIwXO\n3ZKotraW6Oho9uzZg8fjISYmhqSkJGprazv9Gk7sW+mPUG63FIyNrUV6OoUtkQgS8lmEkyehpAQK\nCmgsWYt7zxW8uOsWjnAB1tDBzLllAPf9fAiDE0O/THi2swPJ/9/evUdHVaf5/n/vqgQCBJJAuBq0\nICiChKQRuzGgYvCCIiOttJcWhVZn2unu3xn9nXXsnmOvaTxrPDPYa3ro82tP38QJindttbuD7QWI\nChjllhDkokAKCBBIQq5ASKqyf3/sVFIJlaSqUjtVST6vtVxSO5Vd36pKaj/5fp/v86SmplJZWcm4\nceNoaGhg3LhxeL1eUlNTwzq/XX0rgxFLvS1F5GI9DrYMw0gDXgTGAs3AH03T/D+GYaQArwOXAW7g\nHtM0a3r6eCLSy0wTvv7aKvf+6aeY9WcpqJ7Kf5U+xsnBk2BSKjPnJfHIPziZPDnagw2er0Cst6UN\nUE8LxPbnvpUi0jORmNnyAP+vaZqFhmEkAjsMw/gQ+AHwsWmazxqG8VPgn4GfReDxRAasXt1mX14O\nmzZZzZ+PHwfgm/rxrKl7lK+MGXDlKC5xxfPww1YnnRjIee9SoL6Uzc3NnD59mtLSUgYPHsyECRP4\n9a9/Hdb57epbGYxQeluqrZBI7+txsGWaZhlQ1vLvesMw9gFpwJ3ADS13Wwvko2BLJLY1NMDWrdYs\nVnGxNasFlCdM5EXHCvLPzYDxQxkxAr7/fbj1VqseaV/QVV/KlStXsnLlyh6d366+lcEIpbel6mKJ\n9L6IfkwahuECsoACYKxpmqfACsgMwxgTyccSkQgxTSuw2rDBCrRaZmeIj+f8rLm81fxd3i2cRGOD\nQdwIuPNO+N73YNiw6A471vSFvpUiEh0RC7ZalhDfAv6pZYbL7HCXjrdFJJqOH7eWCDdtspYMfaZN\nw3tDDh81zWfdWwnUtGRaXncdLF8OY8dGZ7g9ZfcSmV19K4Oh5T+R2BaRYMswjDisQOsl0zTfazl8\nyjCMsaZpnjIMYxxwurPv95++nz9/vqa1pVvKOwlTfT189pk1i3XgQNvxMWPgxhshJ4cdJyfwwgtw\n9Kj1pSuvhEcfhalTozPkSOmNnws7+lYGQz/zIvbJz88Paak+kEjNbL0A7DVN0z+z9M/ACmAVsBx4\nL8D3AfQ4V0IGHuWdhMDjgV27rADriy+s22AVGZ071+pNOGMG7iMGL/zOuitYM1grVlh3ifXkdxER\nu3ScBHr66adDPkckSj/MBR4Aig3D2IW1XPg/sYKsNwzDeBg4AtzT08cSkRCUlFgBVn4+rWuBhgFZ\nWVZfwmuvhYQEqqpg3W/go4+s9K1hw+Dee+GOOyA+PqrPQESkX4jEbsQtgLOTL9/U0/OLDHT+BSvd\nbjfFxcUUFBTw9ttv09TUxJw5c8jIyLDuU13N6TfeYExxMbQsswKQlmbNYM2fDy1FOy9cgHdeg7ff\ntnLinU64/Xa47z4YMaK3n2Xfp96aItKZPrJpW2Tg8g+2kpKS2Lp1a2uVcofDwReffcYNcXHw4ouw\nYweO8nIYNQoSE+H6660g6/LLW9cCTdPKiX/xRaishHPnzuFynWLhwpMcP/4lO3eGl//mP87u2tb0\nx4rn6q0pIp1RsCXSh7RWKXc4SG9qYs6OHVx65Ai1H33EiLQ0cDqpmzaN1EcesSqNdlgHLC6GNWvg\n0CHrdno6PPLIUDIyJgGTWLnyw7CDBF8AFUzbmv4YbPnPbJWWlmpmS0RaKdgSiXH+lc/3ffIJ00+d\n4lvV1cwpL2fEF19QBxxMTOTTK6/kyGWXMWTcOB7Pzm53juPH4b/+y8qPB2via/lya1Ux0snvalsD\nZ86cifYQugxoQ2mY3R8DY5HepmBLJMbNmj6d651O2LiR0upq6j0ejOHD2VNby4V58/jm0ksZ953v\ntAYy/luUa2vh1Vfh/ffB67U2IC5dCkuWwODBkR2nLygMpm1NKO1l+gr/Gazm5uaoP7/OgqRQG2Yr\n2BLpOQVbIrHINGH3btiwgWnvvdeasT4mLY3dZ89yZPJkfvnRRzwwcyYej4dHOlQpb2qCv/wF3ngD\nzp61Zq9uvRUeeADsqrHpa4cTTNuantasCVZvBgr+mxcuXLjAU0891X7zQi/p7jlr5lGk9ynYEokl\nx49b5Ro2bYKKCgCGxcXB9OmwYAGD5s7lkoMHOVNcjCM//6Iq5aYJ5eVTeewxON1SRnjWLPjBD6Cr\n630oy0rdCdS2prKyMipta3oz2PLfvDB48GAcDgdbt25l3rx5vfL4Pr7nHKjxNoTeMLs/zkKK9DYF\nW2Ibuy90kQwQgmHb86mrs6q6b9x4cVX3nByG5uTA+PGth6uqqli2bBnPP/98u5mI/fut5Pf9+637\nXnopPPwwXH111w9fVFTEqlWrcDqdHDx4kA8++ID169fz05/+lMzMzKCfhu+1CdS2Jikpqd170x+X\npWJtxqizxtuhNszurVlIkf5MwZbYxs5gK9S8k2DH4wvgSktLLwrgOvv+sJ6nxwM7d1qzWF9+2VbV\nfcgQmDfPKjp61VUXZa+73W4KCwuBth1vZ87E8+WX0/jqq5EAJCXBsmVw881W7azuFBcXk56ejtPp\n5LPPPiM9PR2v10txcXFYwRZc3Lam4wW7t4KtzmZ3OsrPz+9xJwv/GaOGhoZuZ4zs0t1MlBpmi/Q+\nBVvSJ4Uzi9BdUOQfwHk8nqACuGDO28o026q6f/JJ+6ru3/qWVQ9rzpwuM9ddLhfnzp3j3/7t3ygt\nreInP9lBY+OtGEY5TmcZU6fu57LL9lJQ4KGgoPshQfsgobq6OqJBgi+IidZSVGezOx2tXLmyx8GW\n/4xRYWFhtzNGdvEvORFIqA2z++MspEhvU7AlF4lUk+dgZxXCEWreScf7BfL5559z8uRJHA4Hp06d\n4i9/+QvNzc3s3LmTa6+9NvyAoarKapmzcSPn9u5l6NCh1vGJE60ZrBtvtGoxBHWqKrZt28XRozNp\nanqYEycuITkZHn54Io89NpTU1OnAXSENb926dXzxxReMGjUKt9sd0SDBF8QMhKWoWJsx6up3NZSG\n2Qq2RHpOwZZcJFJNnoOdVQhHqHknAKtXr6a6urrTcxqGwciR1lKcaZokJye3HgcreAw0+xEwCGts\ntIpabdhgLReaJgB1wNA77rCCrClTQipyZZrw7/++iU8+WYrHMw6oJy5uP8nJH5GaeiOpqY8FfS5/\nixYtYv369a3PN9pBQiT1ZqDgP2PkcDi6nTGyi4IjkdijYEv6pHBmEboL/rpbBuos8GydtTFNK0t9\nwwbYvNmquQAQFwezZ8OCBew/e5axCxaE/HwPHbKS319+2UVz81kMoxT4DcOHnyIxcQKbN2/mscfC\nC7ZSUlK4//77qampIS4uzpYgIVoBQDQet7y8vN3/e5uCrYvZlT/qy/GsqKggNTXV1uBa9c76NgVb\nYhs7PxhCzTsJRncBXGfPJ76qCl57zdpNePJk2xcuv9yawbr++tY6WbtWr2bTZ58FPaZz54ZSVDQT\nt3sSAPX1ZcTHv4XD8QFnz57i7NlUjhw5QlVVVY9yjpKTk3n88ccv2uEYKQPhIuHL+SstLW1dfj59\n+nS3OX9iPzsCFf8cT6fTSWVlZVA5nuFSsNW3KdgS29j9wRBK3kkw4+luGajd958/D1u2wIYNfOvL\nL8GXhzVypJWDlZNj1V7oINil1YYGeOsteOcdGDTIKrO1eDHU1W3m9df3U18/ipqaSsaMGUNiYiLL\nly/nhz/8Ybfn7cxAyKmyW25uLqWlpZSXl+NwOFpntnJzc3niiSeiPDqJtFgr9SGxTcGWDBjBBH++\nAO7nP//5xR+Yzc2tVd3ZutXKywKGJifDtddauwkzM8HhCHuMzc3w8cewbp2VVw9WJYjly2HcOCgq\nms1nn71PSUkJDoeD+Ph4JkyYwJw5c8J+TIkMwzBIT0+nvLyc5uZmANLT01tz/iR67Niss2nTJmpa\ndhTX19eTmJgIwKeffmpLqQ8Vl+3bFGyJdOfYMWuJcNMmqKxsOz5jhjWDNXdu28xWN7oK+Hbtghde\ngJaNoEydCo88AtOmtd0nKSmJ5cuXU1BQwJEjR1i4cCFz5swhKSkp9OcV5LgkOKmpqVRWVjJu3DiS\nk5MZN24cXq+X1NTUaA9twLNjs05SUhI7d+7E6XRSVlbW+n7PmjXLlplMzT73bQq2RAIY1twMeXnW\nLNY337R9Ydy4tnIN48aFfN5AQc3Ro1aQtWOHdXvMGGsm67rrLt6s6NspunjxYl5++WWeeeaZkMfQ\nka/Uh9vtJiEhoV2dpv4chPmXOElOTu7x8/bl/Hm9XqB/7eqUi61YsYLTp08TFxeHy+Vqfb9XrFgR\n7aFJDFKwJbaLVN0uu7Qmnno81G7axKE//pFVZWUcfvJJ0iZOZFBSEpXTpjHqnnusaaYILQtVV8PL\nL8MHH1gbGYcOhXvusXKzBg3qfKy+Zsf19fURaXYcqVIfgcYay/yfdySesy/nLzc3l88//5xZs2ax\nYsUKJcfHADs+Z/xzPHtjN2IsfFZK+BRsie3suphHhGlyassWXB9+yIWPPuL4rl0McThwmiaFTicf\njRrFPb/6FcX79jF/+vSIPGRjI7z7Lrz5ppUI73DAbbfB979vtdrpin+z48TExKg1Ow5GrAdbdkhJ\nSeGJJ57gvffeU1J8DLErUPHlePYGBVt9m4KtPszu3oO90eQ5atuZz5yxqrpv2ED6rl0wahTlpaXU\npKTgnjyZn5WVMX/WLCZOnEjeRx+RlpbW44c0TeshX3wRKiqsY9dcAz/4gVVMPhjaASU9FenlUzvE\n+my4SKgUbPVhdhbqC7XJc7h6Ndi6cKGtqvuuXa1V3SsuXOBz4E/Dh1NiGHDqFGWNjWzevBmPx0Nc\nXBxTpkzpUYLq6dOj2bVrFmfOWBXqk5OrmDVrF83Np1izJvjz2NnHECLTkNnHV5FeYkukl0/tENOz\n4SJhULDVz4UTzOTl5dHU1MT+/fspLS1lz549pKenhzx7EhNF+EwT9u61AqwtW+DcOet4XBx8+9uQ\nk8Op+noWL1hATYcK8vPmzWPixIlMnTqVtLS0sD7wT5yA3FxrQ2NqKlxxBTz0kJVf73CEXq7Bv8p9\npPsYQmQaMvvk5+drB5WICAq2IiYagUUwtWPCGdenn37KgQMHaGho4Ny5c1RVVfHll1+yc+fOkGZP\nAj12x5mTcGvHdLsUUlZmlWrYuNH6t88VV3BqxgzGLl0Kw4cDYLZ8b1cV5IuKikIaX12dVVQ+Lw+8\nXhg8GJYuhSVLICEh5KfbKtaaHYuISPcUbEVINIKtYGrHhDMFv3TpUi5cuEBcnPXjERcXx4ULFxg5\ncmRIsx6BHrvjzEm4Mx8Bl0LOnrWKjf7ud+AXkDBqVFtV94kT2Zefz9iWQMt3Lri4grz/brJg31uP\nxwqwXn3VGo5hwM03w7JlVnH5nqqpqSE7O5uCggISEhJobm4mOzubmpqamNv1FvVZTRGRGKFgq58L\np3Ly9u3bqa+vB6CxsbG1SvL27dtDCrZ6peJxczMUFlrLhAUFrVXdGTwYsrOtAGvmzC6rurtcroDB\n8ujRo4MehmnC559bS4a+9ohZWfDwwzBpUojPqRvDhw/n5ptvZvTo0a2Jw7FIwZaIiEXBVoTY0Q6i\nO/4J0p3xBRHV1dVBJywPHz4cj8dDQ0MD1dXVJCYmkpCQwHC/2aBgFBYWsnLlSqqrqzl16hSlpaXs\n27ePDRs2kJaWxtixY4HwZreSamqYdPgwriNHGHL+fOvx02PGcHjSJI5NnIinttaqsfDuu+2+Nzk5\n+aJdTW63m6SkpNaNAc3Nze02BnQ1c/n117BmjZUaBtbOwocfhquvjlhJrlZ9IbnZjt1u2p0mIn2Z\ngq0IsaMdRCT4lvJCWU4cN24cL730Ek6nk927dzN9+nS8Xi8PPvggjz32WMiP7dvdOGXKFI4dO8b8\n+fPxeDyh726srYVPPrHysHy5YxMnwvjx1gxWTg5Tx4zhuiDGFSh46qqsQqDSD6dPw9q18Omn1u2k\nJHjgAbjlFnA6B26AYEdAaMfutM2bN1NVVUVBQQF79+5l+vTpTJo0iZtuuqlfvz8i0vsUbPVz4Vw0\n7r33Xvbv34/b7cbhcDBy5EhcLhf33ntvWI/do9pQHg9s22YtE27fbmWbAwwbZvWzycmBK6/s8RRS\nYWEhf/vb31rLKjQ0NLQrq7Bw4cLWC/zZs/DWW/Dee9DUBPHxcOed8L3vtW+RqO3rsW3//v1s3LgR\np9NJcXExiYmJfPHFF1xzzTVhnc8/uB5IbY9EpHsKtiIkVj9MfeMKZZnTt+xXWVmJw+GgsrKSIUOG\n8PTTT4dVOynk2lCmycgzZ5hcUsJlR44wqCUPyzQMTo4bR8nkyRy/5BK85eXw+ushj6e6uhqXy9Uu\n+MnKyiI1NbVd6Qf/sgppaWl4vVZrnVdegZY0Nm64wSrlMGbMxY/jf/EtKCjQxTfGlJSUkJ6ejtPp\n5LPPPiM9PR2A4uJiMjMzQz6fgmsR6YyCrQiJ9YtnuMuckai7FHRtqIoKq8T6xo1WYSrDAJfL+i8n\nB+bP58qUFG6k50t0HQPPsrKyTssq3H77It544zBvvGENC6wWiY8+atXN6ozL5SIpKYm8vDwOHjxI\naWmprb3TYqKumU3sCFz37t1LVVUV0PZHQEJCAhW+8v7djCfQ44bSeaE/v18i0p6CLbFdRkYG69ev\nx+l0EhcXx6FDh/B6vSxdutSq6v7559YyYVFRa1V3kpJg/nwryJo8+aJz2jGL0LH0w9SpU5kxYzG/\n+lUSn38+jaFDrfSwFSvg2mu7X7n0r8Tv8Xh6XIm/u4tzTy7esd6eyf/9XrlyZUTe7+nTp+NwONr9\nEXDo0CFSU1O7/d5AzyPUzgsKtkQGDgVbvSiaCdPR/FDPzMzkueeeIy8vj/Lycm695RYWT55MUn4+\n/Ou/Wt2Ywarq/p3vwIIF8K1vWbdt0vH1GDduHNDWWPaf//mXnDmzjH/5Fyv+GzNmKPfdB4sWBT+s\nSPcxVHumyJozZw5bt25tve31evF6vWEXiFXfShHpjIKtXhTNnI5oX+BSUlJYtmABQ/70J+7Oz4c3\n3mj74tSp1gzWdde1VnW3W8fXw5fT5vHEsW/flZSVPc1zzx3AMEyuuOJrXK497NrVyK5dwT9GQUEB\ndXV1VFZWcurUKV577TVGjRoVdh9Dt9vdunwWqJRHuHXN+nqQEG4Al5KSQkZGBq+99hrnzp2jqKiI\nBQsWBFUgtmMOZHV1Nfv372/NTfR/rzp7v3ulDp2IxAQFW2Kvs2dh82ZrmXDfPmZ8/bWV6JSa2lqu\ngUsuifYomTkzixMnruSll6CuLg6ncxfXXRfPwoVlzJ49E5fr70I+5+9//3vWrl3LmDFjKCsrY8yY\nMXg8Hr773e/ywx/+MOTzrV69murqaiBwsOWraxYq/w0MwQQJPWFHgBFusHXVVVeRl5dHZmYmtbW1\nZGZmUlZWRlJSUrff2zEHMj8/n9LS0tbcRN8fU131rVTfyM756t51XNquqamJ+h+OIuFQsCWR5/W2\nr+re1GQdT0igxOVi6jPPQEZG5Ct+hqmoCH7zmyk0NY1j0CCYMwcqK3/H88+/DVycLxYs0zQxOjxH\nwzAwfXlpIfK/wAeaGQ13ttR/A0MwQUJPxFKAEekZPfWtjJzi4mK2bt160dJ2dna2gi3pkxRs9XF2\n5IGFXQHc7bZ2EubnQ8suLwzDapeTk8ORtDQKXn2VhTNnhjWujnqa1F1aCi+8YJXxOnduJJdeCsuX\nW+Uc/vznIh5//PGwSl34FBQUAHD69GkuXLjA6dOnGTVqFO+++y5l/s2xg+Q/6xRohijc97uvBwnh\ndm/oyYxex+4NvvfDt8Fi27ZtTJ06tcufSQUNnSsoKAgYCBcUFLB48eJoDk0kLAq2+jg78sD8z9nt\n+WpqrKruGzbA4cNtxydMsBLdb7wRWnoMlkRwVqMnSd01NVatrL/9zWqtOGQIPPTQUO68EwYNsu6T\nmJjIkiVLevR6+s8Y5ebmct999/Voxsh/uSxQcBHuxdt/F2YwQUJP2BFghFvWJJIzer73w7fB4uDB\ng92eQ8FW5/zLcvgHwrHWbF0kWAq2JHRNTW1V3XfsaF/V/frrrSDriitsXSYMZwmosRH+/Gd48004\nd84a3sKFVoudHkxgdSrSM0Z2XpxDCRJ6IpYCjL4+o9ef+Zfl8A+Em5uboz00kbAo2JLgmKbVcXnj\nRqsZYH29ddzhgGuusQKsa65pmxoKwLfc09MiqRDaEpBpwtGjl1FYmMm5c8MAGD/+JN/61i7Ky2tY\nvfri81dUVLRWkQ9XTU0N2dnZFBQUkJCQQHNzM9nZ2UHtduuMb4m3rKyM1atXR3zZ2M5K9z0t++Bb\nNi4tLW23bBzuOWtqaiI2oxdLQWR/EKgsh8fjITs7O4qjEgmfgi3pWkUFbNpkzWIdP952fPJkayfh\nDTcEPS2UlZVFdXV1RIKtYJeA9u2DNWusp5GWZhWjf+QRyMqaCszv9Pz5+fmtgUy4fMHK4sWLiY+P\nj8jztqO4px3nDKSnRVc7KxDbkyBz/vz5EZnRU7AVWRkZGcybN++iQLjG1ydLpI9RsNXH2TIr0dDQ\nVtV99+62qu7JyVZV9wULrKglirpbAjp5EtauhS1brK+npMCDD1pDdziiMWLpib5eC0xC4/vs6hgI\nK2dL+ioFW32cy+WipqaG4uLi1v57GRkZoQdapgnFxdYy4ZYtbVXd4+PbV3VvuciFO9ZI+eqrrxg+\nfDh79+5l7969pKWlMX36dLZt28exY9n89a/g8VirmnfdBXffDQkJvTvesHd1xpBIVnwPd9cgRL5A\nLES25lcsVcYXkdijYKuPq6qqYtWqVaSnp7cur6xfv57nnnsuuL8CT5ywAqxNm+D0aQDOnDnD4cGD\nOXb55bxx4gTfGTwYNm9mdkMD8+bNs/kZBWfevHmtRSnXr1/PzTffxqBBd5Kbm0hdnZX8vmABLFtm\n1U+NhpB2dcaoSAYR/rsGQz3v73//e9asWdPuWFNTU8gFYrvb0RkuBVsi0hUFW31cbm4uFRUV59oi\nUgAAIABJREFUDBo0iNraWo4dO0ZFRQW5ubk88cQTgb/p7Fn47DNrmXD//rbjo0dDTg4jc3IYOWEC\ns4Ffz5/P448/HpGx+mZ5grlfdxeuqqoqnn76aUpK3Jw9m8EvfpGKw+Hmiiuu4OqrB/HIIwH7Vwd9\n/lDuJ6EL9bWtr6/n+PHjOJ1OTNPk3Llz1NXVUe/bqGHT44qIRIKCrT6uqqqK6upqTp8+TW1tLYWF\nhQwaNKi1Rk0rrxd27rRmsb74ol1Vd+bOtaaBZsyIiaruwVwQX3/9dT799ATl5UvxeFKoq0skPv4g\nmZlf86//eleXTyNWL7jBjivcYq7hPO+eLP115F8I1H8HaTBee+01Lly4QFNTE42NjdTU1BAfH88f\n//hH6urqgj6P/+MePXqU3NxcAMrKylixYgVgLS+GWsxWfQ5FpCsKtvq4lJQUkpOTmTBhAgcOHCAr\nK4sTJ060XXxLStqqurf01cMwIDPTCrCuvTb0ZKYQ+fqcPf/88xQUFFBUVMQVV1zBkCFDAvb46+5C\nfO7cUF56yaC6+h8wDDDNcpzOtxk0KJ/164cxZMjubscT6PzV1dWtvQdLS0t56623gN7LsQp2Ri/c\nYq7hBFvhFgztjn+fx2DEx8czbNgwDMPg7NmzpKSkYJom8fHxIT1uVlZWxGZq/dnVhsiODhEi0vsU\nbPVxK1asYNu2bUycOJHBgwczceJE4urr+fsxY+C//Tcr2PK55JK2qu69mMjk63NmGAZpaWnMmDED\nj8fDP/3TP1FUVHTRxbyzC3FTUxx7907nwIErqasro7m5gYSED7lw4QUcjqFcuGDS2Nj9xTeYC66v\nBEKszVZ0tStv3rx5PSqtEKiGVTiCCepCDeLi4+N55ZVXcDqdxMfHM2LECLxeL/fddx9PPfVU0OeJ\npd6MwfAPqiLVIUJEep+CrT4uJSWFn/70p3y1axfXNjezdO9eXGfOkPjaa9YdEhOtqu45ObZXde9M\nZ33O8vLySEtLu+j+HS/EXi989BG8/LK1wzA9Ha68spajR/+DoUPr2L3bwRVXXIHH42H58uWtCdOd\nXfS7uuD6L9E1NTWRmZkZctAR7jJlMEt2XRVz/fjjj7t83M6WuiJdw8qOZdof/ehHVFRUUFJSQmlp\nKZMnT2bSpEn86Ec/iujjhEuzTCLSFQVbfZlpwoEDZG7ZQuannzLG4yGzsRFGjIDZs9uquoe41BJp\nvj5nDQ0NlJWVtQsQFi5c2OVf6zt2WM2ijx61bk+dCo8+CmPHTuV//a8rKCkpweFwMHLkSCZNmsQ9\n99zT+r2hXvQ7LtHV1tby7LPPBrVE5y/cYCOY2Z6uirmmpaV1+f2dBXLRqGEV6uuTkpLCv/zLv5CX\nl8eePXtYunRpWLNvdgVFCrZEpCsKtmJIZxfpi46fPm3lYG3cCMePc+7cOQ4eOsSBxka+aWritlWr\ncGVm9vxxIqSzPme+AKEjl8uF220FWbt2WcfGjoUVK6xcfmtyLvyLb2fPMTc3l9LSUpxOJ6Wlpbhc\nLqqrq7ve2dnLuirmunbt2i5nxjrLVfOfLfNPYg+3hlUwyeLh/Jz5+jc+//zzYQeBCopEJBoUbMWQ\nLoOgceOsYqMbN1pV3VvUx8Xxa7ebXamp/LWsjNn19bz42GO88sorTJo0KfTHseFi1Fmfs0WLFlFU\nVNTuvlVV8Ne/uvjwQ2vibtgwuPdeuOOO9hN0vsThtLQ0pkyZQlpaGkVFRe1yXELdSecfdFRUVLQm\nJpeXl4fUJiTcnWnBvPYpKSmd9vPrbmass/fXf7bM18ImUOujYPW1vCgREbsp2IplLVXd0958E557\nrn1V9zlzYMECnvzDH9g9fjxxcdZb6fv/s88+y29/+9tojbydlJSU1obMZWVl7Roy+y7+Fy7Au+/C\nW29ZT9PphNtvh/vus1ZFOwomcTjUJOxIBR3hBhvBBrq+GZ5Q+/l1dv7uWh+JiEjPKNiKIb6ZmOG1\ntUwqKWGS283Qc+eoqqritNfLkaFD2ZWSwu6kJBrWr4f169mzZw8ejweAxsZGCgsLAdi/fz/79u0L\n+DgJCQkBAwK7agWlpaWxv6V4akVFBQCnTp0iIyODyy5zsXEjvPSS1SwarO5AP/iBtXmyM3a0wsnI\nyGD9+vU4nU7i4uI4dOgQXq+XpUuXhnW+3hbu8/afLYuLi2s3W9ab4xAR6a8UbMWK+npyGhqYWV8P\n33xjHUtLgzFj2DtuHNN/8hOyx4/n/g7f9o//+I8UFxcTFxdHQUEBWVlZeDweMjIyOp3ZCjQT5B+8\nQPelAEJZckxKSmLjxo2tLYUcDgcbN25k2rTvsXo1HDpk3W/yZHjkEZg5s/tz+gdVnQWIoV70MzMz\nee6558jLy6O8vJxbb701ppKw7XzcSORDRWIcIiL9kYKtaPJ4rAzwDRvgiy+YcOoUjBplFRmdN8/a\nTXjVVZz+5BOmjx8f8BRPPvkk3//+9/1O6cHr9fLkk0+GNJSOgVZnpQB8gYd/sNVd4JWXl4fH42HP\nnj1cuHCBbdtOUFr6d/zkJ3WkpSUxahQ89JBV/qtjZYqe5JH1JAk71CW6nj6uiD//n/uuOgbEajeE\nnuoPTdxF/CnYiobDh5m1c6e1vc6XeG0YxF19Ndx110VV3bv6cJk0aRKvvPIKzz77LDt27CAjI4Mn\nn3yy0+T4rs7nO+5fCsDj8XRbCqC7D/yjR4/yxRdfUFcHDQ0PsXHjnTidg5gypY4HHoDvfhcGDw78\nvcFcTGL1gqPq3xIu3890dx0DYvVnv6eCmbkW6UsUbPWWqir45BNrFsvtZurXX1tFRtPSrBms+fNJ\n7qSqe3cfppMmTeK3v/0t77//flBJ8YHO579zL5hSAKHkdx04cJgzZ26guXkp4CAuLh7D+JiMjHru\nu+83QZ2jK7F6wQkmiT9U+ot/YIlGDTQRiTwFW3ZqbIQvv7QCrJ07obnZOj58ON9cfjlTf/UrmDIl\nJpo/++/c89+Vt3///oC78vwT7DsrsWCacOzYpfzpTzdx9uzwlqe5E693HaZZwtatE1i5suu2QWrw\n295A/YvfP8j03+DRX4NM3+9UVx0DDh48qN8PkT5CwVakmSbs32/Vw/rsMzh71jrudFrb7BYsgNmz\n2f7MM9xx+eUReciioiKKi4upra1l3bp1ZGRkkNlS1DScWZ9FixaxY8cOSkpKqK2t5dVXXyUlJYXG\nxkZWrlzZel7fB3+gx6ioGMXnn19Fff0lmOYJBg06TFzci5w9+xkJCcNxOhNJCKIBdmFhYetjdiZW\nLzj+AUJBQUG/DxB6y8KFC6M9BNv5/vjpqmPAsmXLVNNMpI9QsBUpp0/Dpk1WkHXiRNvxKVOsvoTX\nXw9JSRF/2KqqKlatWkV6ejrNzc0cOHCA9evX89xzz4WU09HxPqZpAtayxbRp03C5XPzsZz8LmJjr\nv0R26hSsXQtHjsDgwRVMm5bKLbccZMeOXOLiDHbvHsHMmTPxer08+OCDPPbYY12OK5jlt1i84Phe\nn6SkJPLy8jh48CClpaUDIrnZDrE2oxfovXO73RQXF1NQUMDbb79NU1MTc+bMISMjI+D7HMz7rxpo\nIv2Dgq2eOH8etm61lgmLi9uOjxwJ8+dbQdZll7Ue3rx5M1VVVRd9GKekpDBv3rywhpCXl4fT6WyX\n0+F0OkPO6fD/0M/Ly2P48OFkZWW1lpPwer3tzul/f5fLxdmz8MYb8Oc/W5ssBw2C+fPL+fnPU7lw\nIZunn74at9vd2sfQ5XJx7733hjSuvsTtdpOUlNTlrs5YCbY2b97M9u3bAWs35urVqwGYPXt22D+X\n/V2g9y4pKYmtW7e2JrM7HA62bt3a6WvY1fvvO95VxwD/+4lIbFOwFarmZqtdzsaNVqB14YJ1fNCg\n1qruZGZay4YdVFVVBfwwzs7ODns4FRUVlJeXc+zYMRoaGsjPzychIaG1eGioLWugfYL8uXPnuuyV\n19zs4ODBKRQX19HYOAgAl8vNlVcWMWbMIIYOncbQoSn84he/iGgfw1DvEw19Jbl53rx5rQHB448/\nHuXR9F2RfL/9f6a7KkcSqz/7ItKegq1glZZaAdamTW2lzgGuusoKsLKzrUZ+XSgoKAj4YVxQUMDi\nxYvDGlZqaiqjR48mPT2dwsJC5s+fz6FDh0ht2dkYassaCK5tjWlauf//9V9QVwcuF8yYAQ8/DJdf\nPhW4tV2QF8mimR3F4gWnsLCQv/3tb13u6ozVXDPpXqA/YkJt6K33X2TgULDVlbo6K8l9wwb4+uu2\n42PHWkuEOTkwblzQp9u7dy9VVVVA+w/jcNuigJXTsX79erxeL2DldHi93h7ldHSXJ3L4MKxZ09YP\ne/x4K8j6zndiYmNlTMjKyiI1NbXLoDUWc80iYSDUFwv0R0yovTX76/svIhdTsNWRxwM7dlizWF9+\nad0GGDKkrar79OlhRRXTp0/H4XBc9GHc7CsJEYaamhpycnIoKSkhLi6O5uZmcnJyqKmpISUlJayL\nW2e98pqbU1i92nppTBOGD7caRd9+O8QF+EnqLxfWcA3U5GY76ov1BQP1/RaR7inYAityOHzYiiI+\n+aRdVXdmzbICrO98p/My50GaM2cOW7dubb3t+zDuSc6Wy+Xi0UcfBWDLli0888wzF309HP7LfkuX\nLuPtt+FPf7JKh8XFwR13wL33QmJi12MbqFwuV7cNnmPl9RkIM1GRFuh1CbWht15bkYFjYAdbZ860\nVXU/cqTt+MSJVoB1443WzsIwBNpplJGRwbx58y76MK7xBXcxpKqqir/8JY+vv76Mm246TGJiGvHx\ng8jOtroMddKqMaBoF6S0c9dfZ+fuSWJ/b+9SHKgzUT0R6P3x39WZkZFBRUUFa9eu7XRXp4ItkYFj\n4AVbjY1QUGDNYu3cac1qgbUm5ivXkJ7e4+SjQBdM3+2OieI9ydmyQ1FRET/72WuUlNxIRcX9OBw1\nJCbu43//73TuuuvKkM8X7Yt5sMFLOK1wujp3MA29ezLecHR27q6aHffGuPoD7eoUkc4MjGDLNGHf\nPivA2ry5rap7XBxcc401i3X11YETj3og0MXJd1ErLS2NyYva0aPw9NMGhw4txzAM4Etuvvkkl1xy\nmHPnqoDQg62+ItKFM2Ox9ENnP5NdNTvu6ntFRKR7/TvYOnWqrar7yZNtxy+/vK2q+4gRtjx0YWEh\nubm57S5O58+fZ8uWLTgcDioqKvjDH/7A7373O+bOncuQIUNa7xfuRc3tdnfb2qYzDQ2DKS6eycGD\nU6iqOgdUMHr0JxjGrygpuZqSEti+PfAW9lDk5+eHPcZw2bnFvqs6ZqGWAvDp7ZIAsRgUikWNx0X6\nh/4XbJ07B1u2WAHWnj1tx0eObCvXMHGi7cPIysq6aHv4unXrGDRoEE6nk6NHj5KTkxNwa3goS20d\n++/5vi/YD+PGRnjvPXjzTXA44MorIT5+A07n6wwZ0sjJk/G4XC68Xi+zZs3iiSeeCPIVCGzlypW9\nHmzZucW+qzpm//mf/8nOnTtxOp2UlpYG/TraOd7u6kMFanbso7pQvS/Ss60iEh39I9hqboaiIivR\n/fPPrQgCrKru2dnWMuHMmVY0EYKe7tLqeGELdqYjlItauPlQpmntDXjxRSgvt45dcw384Adw5kwq\nq1bV4/U6GTlyZGvtrpycnKDOLZYVK1Zw+vRp4uLiSEtLY+LEiXg8HlasWBG1MXVXHypQs2Mf/5/l\n7mZfg8kB60jLlCLSX/XtYOvYsbaq7pWVbcdnzLBmsObOhaFDwz59TxK7A/016j/TsXv37k5nOuwu\ndvjVV1ZR0m++sW5PmgSPPGJ1GQKYODGT5557jry8PMrLy7n11luDbq/TW0IJhO28gHd17u762oVz\nTjuEUx8q2I0BXeWAhXJOEZG+rO8FW3V11pTMxo1t0QJYtQhycqxyDWPHRuSh3G43H3/8MSUlJbz9\n9tuUl5ezcOFCMjIycLlc3TaS9QUDPitWrODo0aOUlJRw/vx5qqurmTRpUsgzHeFelE6ehNxcq6Uj\nWCurDz5ovWwdJ/266scWC0IJhKMVbEF4r2NvjzfYoDDYcSkHTESkvb4RbHk8sH27tUy4fXtbVfeh\nQ+G666xoYdq0iPeKSUpKYuPGjaSnp9PQ0EBtbS1bt25t3d7dXdAT6Gumr9REJ7c7+z5/oQZbdXXw\n+uuQl2e9dIMHw113Wf8lJAR9GukHAv3c+GYJ09LSWLhwIWlpaRQVFV00S+j/72A3BnSVA9aRcsJE\npL+K3WDLNOHgwbaq7nV11nHDsMo05OTAnDlWXpZNcnNzqaioYNCgQdTW1lLektyUm5sbVKJ4xwtb\nXl4ew4cPJysri/z8fLKysvB6vRf9xR+pmQ2PxwqwXnsN6uutl+6mm2DZMhg1KiIPIf1AOEnYXW0M\nCDYHrCP1ChSR/ir2gq0zZ6wcrA0brJwsn8susxLdb7gh7Kru/oKZHaqqqqK8vJz9+/dTXV3N119/\nTUpKSmsz6c7+uq+uriY5Ofmi4/5/8Tc0NAT9F39H3c0AmKa1TyA3t63iRWam1Sx68uSgHyYqQpm1\nCycJO5zHkdCoR6CISHuxEWxduGBVdd+wAQoL26q6jxjRVtV98uSILhMGc7GdMGECjY2NuFwuTpw4\nQWJiIkePHuWGG24AOv/rvrMcIv+/+AsLC4P+iz/Q+TvzzTfw/POwd691Oy3NCrJmz474Kqstgg2C\nwk3CDvVxJLD+sDFARKS3xEaw9eCDcP689e+4OPj2t60Ay4aq7qEwAkQnpmkGPB4MO//iLy+HrVuz\n2bHDuj1iBDzwANxyS1RfQtsoCTu6+trGABGRaIqNy/D583DFFdYy4XXXWX0KbdZVgq/Ppk2baGpq\nwu124/V6qaiooLm5mVdeeYWTJ0+2q5Xlr7NlPv+/+B0OR9B/8Xfkf1Hat+8Ia9eeY+vWVCorpzBq\nVDnZ2RU89NBQpk+/LKTz2iWUKtjBvC8QfhK2j5KxRUSkt8RGsPXb31rrXb3IfwmwsyWlpKQkPv/8\nc06ePElpaSnp6emMHz+ea6+99qIEef9zdBYs+O/6GjlyZKe7vrpj1eeCDz+El1++jJoaSEqCxYtH\ns3w5jBkzOuhzdTVOiEyLkFASsLtKvPYXbhK2j5KxRUSkt8RGsNXLgVZHnQVbS5Ys4Y033sDpdOJw\nODBNk2PHjrFq1aqgz+HPP+hIS0sLa2bFNGHHDnjhhbb9A9OmWUVJp04N+XTdjjNWZ3/6QhK2+tqJ\niAjESrAVBcFc7LZs2cLcuXM5dOgQAKNHjyY9PZ0tW7YwadKkLs8dKPjavHkz27dvB6CsrIzVq1cD\nMHv27NbaXV0pKbGCrMJC6/a4cbBiBUybVsX69Xm8/34FqampvVrtPdKJ5sGeK9wk7FAfpyf6QtAq\nIiL2U7BF53lCmzZtoqampvV2SUkJJSUlAfOC/HOAXC4X+fn5F13Q582b1xpUvfvuuzz++ONBjfXM\nGVi3Dj7+2JrZGjYM7r0X7rgD6uvbduU5nU4qKytD2pXXU9EKtqBnVe41syQiIr3F9mDLMIyFwGrA\nAawxTfPiNbgo6yxPaMqUKa15QW63myVLlnSaF2RHDlBDA7zzDrz9tlUdw+mE22+H++9v20OgXXki\nIiKxzdZgyzAMB/AbYAFwAthmGMZ7pmnut/NxI6UneUHd7apzu92sXLky4NdME0pKJrF7dybnzw8B\n4JJLSsnKKuTEiTr+4z/a7vvOO+9QUVEBwIULFxg8eDAAqampIRVKDZd29UlHylUTEWnP7pmtbwPf\nmKZ5BMAwjNeAO4GYCrY6uwDU1NSQnZ1NQUEBcXFxNDc3k52dTU1NTbdNervbVZefnx8w2Nq9G9as\ngaoqmDgRpkyxkt9nzJiKFbO25z/7Fs6uvJ6K1q4+XdBjl3LVRETaszvYugTw67lDKVYAFlM6uzj7\nji9evJgtW7bwzDPPhHyOjnwtZkpLS9u1mCkttZLft22z7peaCg89ZBXQ76qGaqDZt8rKypjalWcH\nXdBFRKSvGLAJ8sFwu90UFxezefPm1uAoIyODzMzMbr83UPDl32LG4/Fw4MABduw4SFra/+DTT4fR\n3AwJCfC978GSJcH12A60Ky8pKanXdiNqFklERKRrdgdbx4FL/W6ntRxrx385bf78+TEzU1FcXMyr\nr75Kenp6a3C0fv16nnvuubD6vPkns5tmPAcOzOCrr6YzaFAVEycOY+FCq8VOgB7WXeq4K683l/YU\nbImISH+Wn5/f4+uqYfqaPtvAMAwncAAr2egk8CVwv2ma+/zuY9o5hp546qmn2LZtG42NjRQWFpKV\nlcWZM2dIT08Panaro4KCAs6fb6CmZiZ7917D8OHpAKSkHOaOO8pJTq7p5gxd8+VsJScnB11WQkRE\nRIJnGAamaYbUJNnWmS3TNL2GYfwE+JC20g/7uvm2mLF3715Onz5Nst9U04QJE7jhhhvCCmZWrfoz\n7747kvPnRxMXd4z09DgyM7dzww0jWLbsie5P0I2VK1eycuVKtaIRERGJIbbnbJmm+TcgQo1ketf0\n6dNJTEwkPT0dt9vN/PnzOXToEKmpqSGdp6wMcnMhP38h5eWHGDLkPMOHv8gtt7jweptYtOhJe56A\niIiIRJ0j2gOIZZMmTaK2tpZDhw61/v/kyZNkZGQE9f319dYOw3/8R9iyBQyjib/7u/Pcdts7eL3r\nMU1vaymJSFIelYiISOzQbsQu3HTTTdx9993k5eWRl5fHrbfeGtRuRI8H3n8fXn0V6uqsYzk58OCD\nw0hNnQXM4o031nZZSqInFGyJiIjEDgVbXfAFLcuWLePnP/95t0VCTRO+/NKazTpxwjqWkWEVJU1P\nt0pJ7NnjBiAuLk6FOEVERAYABVstOmuo7CtCWltb264Iacfv9XhcrFkDe/ZYxy65BG699QRLlkxo\nLUrqcrlISkoiLy8PgNLS0k7Pp+BLRESkf1DOVgtf6xd/viKkBw4coLm5mQMHDvDss89SVVXVep+K\nCnj2WQ9PPGEFWsOHwz/8A/zmN5CS8nW76u9FRUX8+Mc/5oMPPgDggw8+4Mc//jFFRUXdjkVERET6\nJs1sdcG/CCnQ+v+8vDzuvnsZb78N77wDJ08mM3YsLF4M994Lw4YFPl9xcTHp6ek4nU6OHTtGeno6\nXq+X4uLisOp2iYiISOxTsNWisLDwovpUBQUFNDQ0ANDQ0EB+fj6mafD++8388pdFNDQkAHDFFeX8\n9rfzGDeu63P6n8/tdrd+bdu2bRw8eLD1fsnJyTFTRV9ERER6pl8HW6HkPmVlZV0U4Kxbt47CwkIO\nHz6MYRg4HLOpqvoucBmXXprG1KlW8vupUycvCrQCnXPdunV88cUXjBo1qrXau9frZerUqe2S71WU\nVEREpP/o1zlbPc19mjt3Llu2bOHYMYOmpqcoLPwubrfJ5Zcn8eST8MtfwrRpwZ9v0aJFVFZW4vV6\nAfB6vXg8HhYtWtSjcYqIiEjs6tczW6EINAP24YdfMnTo/2Dv3jRMs5yhQ02uvvogt99+lOuu+36X\n3xvoeEpKCvfffz81NTVs27aNqVOnBtyNqJ2IIiIi/Ue/DrYC5WEFw+t1sn//lXz22QSampwYRiMO\nx/u4XCXU1JzltdcScLu/DmtMvibRBw8e7LRul4ItERGR/qNfB1uB8rC6YpqQnw8vvggNDTB+fClJ\nSQeYNauIF15Yx003PR4wxyoUyscSEREZWPp1sBWKPXtgzRrwbQqcPBl++tMR5OV9TFyc9TIpx0pE\nRERC1a+DrWCW444fh9xcKCiwbo8cCQ89BDfeCA7HCK666kny8vJwOByd5lgFy+12U1ZWRn5+PsnJ\nyWrXIyIiMgAYpmlGdwCGYUZjDHV1VqPo9evB64WEBLj7bliyxPp3Ry6XK+Yru69cuZKVK1dGexgi\nIiL9lmEYmKZpdH/PNv16ZiuQpib461/h9dfh7FkwDLj5Zli2zJrV8ud2u1sDLDWOFhERkXAMmGDL\nNGHrVmvJsKzMOpaVZRUl7Sxu8g+q0tLSVNVdREREQjYggq0DB6zk9337rNsTJ8LDD8PVV9OuUbSI\niIhIpPXrYOv0aVi7Fj791LqdlAQPPAC33AItPaVFREREbNUv2/WcPWstFz72mBVoxcfD974Hf/gD\n3HZbzwKtWEySr6qqYt26dRQUFLBu3TqqqqqiPaSLxOLrJiIi0hv61cyWxwMffggvvwy1tdax+fOt\nUg6jR4d+Pv8E+YSEBPLz8yksLGTJkiUxkyBfVVXFs88+S1xcHA0NDRw4cICvvvqKJ598MuwSFXYI\npSm4iIhIf9Ivgi3ThO3b4YUXoLTUOjZ9Ojz6KFx+efjn9U+Qz8/Pb02Qj6WgIS8vj7i4OJwt03W+\n/+fl5YVd5V5EREQip88HW4cPW8nvu3dbt8ePhxUr4Npr7Ul+D7ffol0KCgpoaGgArNkj39i2bdvG\nQV85/BiQnJys3ZwiIjIg9dlgq7ISXnoJNm60ZrYSE+G++2DRIoiz8VmF2m/RbuvWrePAgQM4nc7W\n2bee9m+0QywFqCIiIr2pzwVbDQ3w9tvwzjtw4YIVWC1aBPfeC8OHR3t0vW/RokV89dVXrbfVv1FE\nRCS29Jlgq7kZNmyAdevgzBnrWHa2tWQ4fnzvjSOW8rUAUlJSePJJq3/jtm3bety/0S6x9rqJiIj0\nlj4RbBUWWnlZvuoBl19uVX6/6qreH0ssBg0pKSksW7aMgwcPxtTSob9YfN1ERER6Q0wHW8eOWTsM\nt2+3bo8eDcuXw/XXq/K7iIiI9A0xGWzV1Fi1sj74wFo+HDLEKkp6550waFC0RyciIiISvJgKthob\n4b334M034fx5a/bqttusFjtJSb0/Hv+ipsnJya076vzrb4mIiIh0JSaCLdO02uqsXQvl5dax2bOt\nZtETJ0ZvXP5BVSyVexAREZG+IyaCrf/+3+Gbb6x/u1xW8ntWVlSHJCIiIhIRMRFsffNM4Xp3AAAG\nRklEQVQNpKTAgw/CggXg6JftsUVERGQgiolg67774O67ISEh2iMRERERiayYCLYeeCDaIxARERGx\nhxbsRERERGykYEtERETERgq2RERERGykYEtERETERgq2RERERGykYEtERETERgq2RERERGykYEtE\nRETERgq2RERERGykYEtERETERgq2RERERGykYEtERETERgq2RERERGykYEtERETERgq2RERERGyk\nYEtERETERgq2RERERGwUF+0BSM+43W7cbjcAycnJ5OfnA+ByuXC5XFEbl4iIiFgM0zSjOwDDMKM9\nBhEREZFgGIaBaZpGKN+jZUQRERERGynYEhEREbGRgi0RERERGynYEhEREbGRgi0RERERGynYEhER\nEbGRgi0RERERGynYEhEREbGRgi0RERERGynYEhEREbGRgi0RERERGynYEhEREbGRgi0RERERGynY\nEhEREbGRgi0RERERGynYEhEREbGRgi0RERERGynYEhEREbGRgi0RERERGynYEhEREbGRgi0RERER\nGynYEhEREbGRgi0RERERGynYEhEREbGRgi0RERERGynYEhEREbGRgi0RERERGynYEhEREbGRgi0R\nERERGynYEhEREbFRj4ItwzCeNQxjn2EYhYZhvG0Yxgi/r/2zYRjftHz9lp4PVURERKTv6enM1ofA\nVaZpZgHfAP8MYBjGdOAeYBpwG/B/DcMwevhYUZWfnx/tIUgn9N7ELr03sU3vT+zSe9O/9CjYMk3z\nY9M0m1tuFgBpLf/+O+A10zQ9pmm6sQKxb/fksaJNP/ixS+9N7NJ7E9v0/sQuvTf9SyRzth4G1rf8\n+xLgmN/XjrccExERERlQ4rq7g2EYHwFj/Q8BJvCUaZp/abnPU0CTaZqv2jJKERERkT7KME2zZycw\njBXA3wM5pmleaDn2M8A0TXNVy+2/Ab8wTfOLAN/fswGIiIiI9CLTNEPKQ+9RsGUYxkLgP4DrTdOs\n9Ds+HXgZ+A7W8uFHwOVmTyM7ERERkT6m22XEbvx/wCDgo5bNhgWmaf7INM29hmG8AewFmoAfKdAS\nERGRgajHy4giIiIi0jlVkO+GYRgLDcPYbxjG14Zh/DTa45H2DMNwG4ZRZBjGLsMwvoz2eAYywzDW\nGIZxyjCM3X7HUgzD+NAwjAOGYXxgGEZSNMc4UHXy3vzCMIxSwzB2tvy3MJpjHKgMw0gzDGOjYRhf\nGYZRbBjGf2s5rt+dGBDg/fl/Wo6H9Pujma0uGIbhAL4GFgAngG3AfaZp7o/qwKSVYRiHgatN06yK\n9lgGOsMw5gH1wIumac5sObYKqDRN89mWP1ZSTNP8WTTHORB18t78AqgzTfNXUR3cAGcYxjhgnGma\nhYZhJAI7gDuBH6Dfnajr4v25lxB+fzSz1bVvA9+YpnnENM0m4DWsF1lih4F+jmOCaZqbgY5B753A\n2pZ/rwWW9OqgBOj0vQHr90eiyDTNMtM0C1v+XQ/swyoQrt+dGNDJ++OrGxr0748uUl3rWJy1FBVn\njTUm1gaNbYZh/H20ByMXGWOa5imwPrSAMVEej7T3k5bets9rmSr6DMNwAVlYHVnG6ncntvi9P74y\nVkH//ijYkr5urmmas4DbgR+3LJdI7FLeQuz4v8Dklt62ZYCWE6OoZYnqLeCfWmZQOv6u6HcnigK8\nPyH9/ijY6tpx4FK/22ktxyRGmKZ5suX/5cA79PEenP3QKcMwxkJr7sPpKI9HWpimWe5XkuePwDXR\nHM9AZhhGHNaF/CXTNN9rOazfnRgR6P0J9fdHwVbXtgFTDMO4zDCMQcB9wJ+jPCZpYRjG0Ja/NjAM\nYxhwC7AnuqMa8Aza5zH8GVjR8u/lwHsdv0F6Tbv3puUC7nMX+t2JpheAvaZp/trvmH53YsdF70+o\nvz/ajdiNlu2cv8YKTNeYpvnvUR6StDAMYxLWbJaJVaD3Zb0/0WMYxivAfGAUcAr4BfAu8CYwETgC\n3GOaZnW0xjhQdfLe3IiVf9IMuIEf+nKEpPcYhjEX+BQoxvosM4H/CXwJvIF+d6Kqi/fn+4Tw+6Ng\nS0RERMRGWkYUERERsZGCLREREREbKdgSERERsZGCLREREREbKdgSERERsZGCLREREREbKdgSERER\nsZGCLREREREb/f9dURXvyfIkPAAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x107a88090>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "plt.errorbar(xs,ys,xerr=dxs,yerr=dys,fmt='o',color='black',alpha=0.5)\n",
    "plt.plot(np.arange(-2,25,1),olsfit[0]*np.arange(-2,25,1)+olsfit[2],color='red',lw=2,alpha=0.7,label='OLS')\n",
    "plt.plot(np.arange(-2,25,1),odrfit[0]*np.arange(-2,25,1)+odrfit[2],color='blue',lw=2,alpha=0.7,label='ODR')\n",
    "plt.gcf().set_size_inches(10,10)\n",
    "plt.legend(fontsize=18)\n",
    "plt.xlim(-2,25)"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 2",
   "language": "python",
   "name": "python2"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 2
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython2",
   "version": "2.7.12"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 0
}