{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Lab fitting from Jupyter \n",
    "\n",
    "### Import the module"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import LabModule as lab"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Load the data"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Reading csv data from ../data.csv\n",
      "[ 1.  2.  3.  4.  5.  6.  7.  8.  9. 10.]\n",
      "[0.5 0.5 0.5 0.5 0.5 0.5 0.5 0.5 0.5 0.5]\n",
      "[ 10.  15.  25.  39.  49.  54.  69.  70.  86. 114.]\n",
      "[ 3.162  3.873  5.     6.245  7.     7.348  8.307  8.367  9.274 10.677]\n"
     ]
    }
   ],
   "source": [
    "x,ex,y,ey = lab.load_data_from_csv(\"../data.csv\")\n",
    "print(x)\n",
    "print(ex)\n",
    "print(y)\n",
    "print(ey)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Initial params"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Initial parameters:  [-5.333 10.624]\n"
     ]
    }
   ],
   "source": [
    "# Give by hand or guess\n",
    "#init_params = [-3, 10]\n",
    "degree = 1\n",
    "init_params = lab.guess_poly_init_params(degree, x, y)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Run fit"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "   Iteration     Total nfev        Cost      Cost reduction    Step norm     Optimality   \n",
      "       0              1         1.6330e+00                                    2.41e-01    \n",
      "       1              2         1.5819e+00      5.11e-02       1.48e+00       7.58e-03    \n",
      "       2              3         1.5819e+00      1.43e-05       1.83e-02       1.24e-04    \n",
      "       3              4         1.5819e+00      1.03e-08       6.30e-04       3.59e-06    \n",
      "       4              5         1.5819e+00      8.45e-12       1.80e-05       9.87e-08    \n",
      "       5              6         1.5819e+00      5.33e-15       3.86e-07       6.89e-09    \n",
      "`gtol` termination condition is satisfied.\n",
      "Function evaluations 6, initial cost 1.6330e+00, final cost 1.5819e+00, first-order optimality 6.89e-09.\n",
      "\n",
      "=== Fit quality ===\n",
      "chisq per point = \n",
      " [0.344 0.075 0.084 0.039 0.022 0.202 0.004 0.778 0.079 1.536]\n",
      "chisq =  3.1637, ndf = 8, chisq/NDF = 0.39547, chisq prob = 0.92367\n",
      "\n",
      "=== Fitted parameters ===\n",
      "param 0 = -3.8631 +-  4.8552\n",
      "param 1 =  10.317 +- 0.94652\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAfkAAAGICAYAAABCyJmAAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjQuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/MnkTPAAAACXBIWXMAAAsTAAALEwEAmpwYAABizUlEQVR4nO3deVyVZf7/8dflhoKC4lZukJmmmKTikpZYllo6ozPWaGHFpFHTMk3LOC0zY30byqmmspoWo8JJSktLCystcynNBVzKPTIx1BRXFJD1+v1xw/mpLALXgfvc9/k8Hw8ewjmH+3x8n/vwOfd2XUprjRBCCCHcp57dBQghhBCidkiTF0IIIVxKmrwQQgjhUtLkhRBCCJeSJi+EEEK4lDR5IYQQwqWkyQvhx5RSJ5VSne2u42xKqd1KqavtrkMIp5MmL4QfKGmauSVNvfSrnda6qdZ6V8ljEpVS/7K71rMppR5XSs2yuw4hnEiavBD+4zclTb30a5/JwpRSDbxVmBCidkiTF8KPKaW0UqqLUioOiAGmlGzlf1rJ4+9WSv0I/Fhy22il1Eal1DGl1CqlVK/THv83pdRepdQJpdQOpdSwktvP2GuglBqqlMoo5/lGAo8C40vq2uTVAIRwOfkkLoRAaz1DKTUIyNBa//0cDx8LDABylVK9gbeB3wApwETgE6VUNyAcuAfop7Xep5QKB+pXs64vlFJPAV201hOr87tCCNmSF8KfzC/Z2j6mlJpvsJyntdZHtNa5QBzwhtZ6jda6SGs9E8gDBgJFQADQQynVUGu9W2v9k/H/QghRZdLkhfAfY7XWzUu+xhos55fTvg8DHjztw8MxoCPQTmudBvwFeBw4qJSarZRqZ/C8QohqkiYvhChV1SkpT3/cL0D8aR8emmutA7XW7wNord/TWl+O9WFAA/8u+b1sIPC05ZznhbqEEGeRJi+EKHUAqO41828CdyqlBihLkFJqlFKqmVKqm1LqKqVUAHAKyAWKS35vI3CdUipUKXUe1hZ/ZXWFK6Xk75UQ1SRvGiFEqbewjp9X+Zi91joFuB14BTgKpAGxJXcHANOAQ8CvQBvgkZL73gU2AbuBxcCcSp7mw5J/Dyul1lftvyKEAFBay54wIYQQwo1kS14IIYRwKWnyQgghhEtJkxdCCCFcSpq8EEII4VLS5IUQQgiXct3Y9a1atdLh4eF2lyGEEELUmdTU1ENa69Zn3+66Jh8aGkpKSordZTjaF198wciRI+0uw/EkR3OSoTnJ0JwTMlRKpZd3u+t21+fl5dldguPt2bPH7hJcQXI0JxmakwzNOTlD1zV5IYQQQlhc1+SbN29udwmON2rUKLtLcAXJ0ZxkaE4yNOfkDF3X5AsKCuwuwfEyMzPtLsEVJEdzkqE5ydCckzN0XZPPzs62uwTHW7t2rd0luILkaE4yNCcZmnNyhq5r8kIIIYSwuO4SuqCgoHM+Jisri4MHD8qu/QoMGDCAbdu2eXWZDRs2pE2bNgQHB3t1ub6sb9++dpfgeJKhOcnQnJMzdF2Tb9iwYaX3Z2VlceDAAdq3b0+TJk1QStVRZc5RUFBwzhyrQ2tNbm4ue/fuBfCbRt+hQwe7S3A8ydCcZGjOyRm6bnf9sWPHKr3/4MGDtG/fnsDAQGnwFThXhtWllCIwMJD27dtz8OBBry7bly1YsMDuEhxPMjQnGZpzcoaua/LnUlBQQJMmTewuwy81adJEDpEIIUQdcl2Tb9So0TkfI1vwlfPmrvrT+Vvu7dq1s7sEx5MMzUmG5pycodJa212DV0VFRenKxq7ftm0b3bt3r8OKxOkkfyGE8D6lVKrWOurs2123JX/o0CG7S3A8Jw/84EsSEhLsLsHxJENzkqGZoUOH0q1bN7vLqDHXNXm37ZmoisTERJRS5X4lJCSglGL37t2exz/++ON8/fXX9hXsJ4qLi+0uwfEkQ3OSoX9z3SV0/nbc93QffvhhmUs9zj//fL777jvOP/98z21PPPEEjz32GFdddVVdl+hX6tVz3WfoOicZmpMM/ZvrmnyrVq3sLsE2l156KV26dClze1hYWLWW07p1a2+V5NcmT55sdwmOJxmakwzNnb6R5DSu+4h3/Phxu0vwKaW78kt315fu6YiPj/fs0n/88cfP+B1vXyfvr5KTk+0uwfEkQ3OSobnDhw/bXUKNuW5LPj8/v2a/+Je/wMaN3iyl+i69FF58sca/XlRURGFhoefn8g5dfPfdd1x22WXExsZyxx13AGVHc5Jr2b1j3759dpfgeJKhOcnQXI37ig9wXZP3ZxdffPEZPw8ePLjMrrqBAwcC0L59e8/3Qggh3Ml1Tb558+Y1+0WDLWhf8fHHH5+xVd6sWTO+++67ai+nxhmKM4wZM8buEhxPMjQnGZpr2bKl3SXUmOuavD/vau7Zs2eZE+9q0uTz8/NrbdQ7f5KRkUHbtm3tLsPRJENzkqG5vLw8u0uosTo98U4p9bZS6qBSavNptz2rlNqulPpeKfWxUqr5afc9opRKU0rtUEqNqMpzZGdn10Ll/iUnJ8fuElwhNTXV7hIcTzI0JxmaO3nypN0l1Fhdn12fCIw867YvgZ5a617ATuARAKVUD2ACEFHyO68qperXXanu1ahRI3Jzc+0uQwghRC2r0yavtV4BHDnrtsVa69JTwlcDpQeVxwCztdZ5WuufgTSg/7meIygoyIsVu1OPHj1YuHAhX375JSkpKWXOvpUMvaN//3OuruIcJENzkqG5Zs2a2V1CjfnaMfnbgDkl37fHavqlMkpuK0MpFQfEAbRp04YZM2YA1srdunVrFi5cCECnTp3o1KmTZ2x2pRStWrXi6NGjnkvPWrRowalTpzxbuk2bNqVevXpkZWUBEBAQQNOmTT3XTdarV4+WLVty5MgRioqKyl1G6Qpy4sQJABo3bkxgYCBHjlifd+rXr09oaOgZywgNDSUnJ4dTp06Vu4wmTZrQuHFjjh496rkNrOs5S4exbNmypaeGw4cPc/7551NcXMyTTz7JY489xm9+8xvy8vJ46KGHmDJlCg0aNKBFixZkZ2d7Dnu0atWKrKwszyUkISEhFBYWeu4PDAykUaNGnmvrGzZsSPPmzc8Y/75169YcO3aMgoICTpw4wYEDB8jIyPDsRizvdRo+fLhnzO1GjRoRGxvL/PnzPfPRjxs3jrS0NDZt2gTAoEGDaNasGYsWLQKgc+fODBkyhMTERE+dEydOZO7cuZ7cx48fz+bNm9myZQsA0dHRNGjQgCVLlgDQtWtX+vXrR1JSEgDBwcFMmDCB2bNne9aHmJgY1q1bx86dOwEYNmwYhYWFLF++nIKCArKzs+nZsydz5szxvK7XX389s2bN8hwWiY2NZcWKFezatQuAESNGcOLECVatWgVAZGQkXbp0Yd68eYC1jo8dO5bExETP6zJ58mQWL17Mnj17ABg1ahSZmZmsXbsWgL59+9KhQwfPvNjt2rVj9OjRJCQkUFxcTL169Zg8eTLJycmeD31jxoyx/XUqKCjgggsuqNXXCSAiIsK1r1PLli3p1auX499Pdr1O+/fvp1mzZmzcuNHn30/lqfNZ6JRS4UCy1rrnWbc/BkQBv9daa6XUK8BqrfWskvvfAj7XWs+tbPlhYWE6PT29wvtlFrRzy8zMrLVR7/wp/xkzZhAXF2d3GY4mGZqTDM0MHTqU/fv3s2PHDrtLqVRFs9D5xJa8UioWGA0M0///U8deoONpD+tQcpsQQgghqsD2YW2VUiOBKcBvtdann9b9CTBBKRWglLoAuAhYe67lBQQE1E6hfqRRo0Z2l+AKnTp1srsEx5MMzUmG5pzcV+p0S14p9T4wFGillMoApmKdTR8AfFkyDOtqrfWdWustSqkPgK1AIXC31rroXM8RHBxcW+X7DcnQO4YPH253CY4nGZqTDM2FhobaXUKN1fXZ9Tdqrc/XWjfUWnfQWr+lte6ite6otb605OvO0x4fr7W+UGvdTWv9eVWe4/QTvkTNHDp0yO4SXKH0BBpRc5KhOcnQ3P79++0uocZs310vhBBCiNrhuiZfr57r/kt1rrzZ60T1ybkN5iRDc5KhOSf/TXRdR3TyRAK+olWrVnaX4AqxsbF2l+B4kqE5ydDceeedZ3cJNea6Jl86KEtdGzp0KEOHDrXlub3t6NGjdpfgCvPnz7e7BMeTDM1JhuacfJ6S65q8P89C5y2lo/8JM6WjVImakwzNSYbmnNxXXNfkhRBCCGFxXZNv0aKF3SXUucTERJRSnq+goCDCw8P53e9+xwcffEB1hy5u0aIFGzdu5PHHH/eMSy2qb9y4cXaX4HiSoTnJsOaSkpJYvXo1O3fuJDw83DPuvpO4rsnn5eXZXYJtPvzwQ7777js+++wznnzySQICArjxxhu55pprqjW17KlTp9i4cSNPPPGENHkDaWlpdpfgeJKhOcmwZpKSkoiLi/P0lPT0dOLi4hzX6H1i7HpvKp2JqC6VftrLy8sjPDyc+Ph4YmJi6ryOSy+9lC5dunh+vvnmm7nhhhu44YYbmDJlCi+//HKVliNzzXvHpk2bGDBggN1lOJpkaM4fM/TGSdClf9NPl5OTw6RJk3jzzTeNlr1s2TKj368O123J1zVf/7Q3btw4xowZw5tvvun5ADR16lT69OlDcHAwrVq14qqrrmL16v8/q+/s2bP54x//CMBFF13kOQywe/duAF555RUuu+wyQkNDad68OQMHDvRMlyiEEG5Q0V5hp+0tdt2WfNOmTWtt2eV9OqzOp726/PR2uuuuu4758+eTkpLCkCFD2Lt3L/fffz8dOnQgOzubWbNmMWTIEFJTU7nkkksYM2YMBw8e5F//+hcffvghHTp0AOD8888HYPfu3UyePJnw8HAKCwv59NNPGT16NJ9//jkjR4605f/oiwYNGmR3CY4nGZrzxwy98bc2PDyc8qYtDwsLs+1veU24rsnXr1+/Tp/PCZ/2SmehKh1/+fSxrIuKihg5ciQREREkJCQwffp02rZty4UXXgiUPQQA8Nxzz3m+Ly4uZtiwYezcuZPXXntNmvxpmjVrZncJjicZmpMMayY+Pp64uLgzDgEHBgYSHx9vY1XV57omf/z48Vpbdnmf3pzwaa/07PrSoRm/+uor4uPj+f777884se6CCy4AICsrq9LlpaamMnXqVNatW0dmZqZn+d26dauN8h1r0aJFxMXF2V2Go0mG5iTDmik9r2rSpEnk5eURFhZm2/lWJuSYvKH4+HgCAwPPuM3XPu398ssvgLW7ff369Vx33XU0bdqUt956i9WrV7Nu3ToiIyM5depUlZY1bNgwjhw5wssvv8yqVatYt24dI0eOrNLvCyGEU8TExDBw4EC6du3K7t27HdfgwYVb8gEBAXX6fE74tLdw4UIaN25M3759iY+Pp0GDBnz00Uc0bNjQ85ijR4/SvHlzoPIMv/jiC44fP84HH3zgOVYP9lzV4Os6d+5sdwmOJxmakwzNNW7c2O4Sasx1Td6O408xMTGek+x8ZRd9qXnz5vHJJ59w3333ERgYSE5ODvXr1z9jVqWvv/6aPXv2eHbXN23a1NPoz76crrSZn/4BYefOnaxcufKMpi9gyJAhdpfgeJKhOcnQXOkGkBO5bne9kycSMLVx40ZWr17NihUrmDVrFhMmTOAPf/gD11xzDU8//TQAI0eO5OTJk8TGxrJkyRJee+01Jk6cSPv27T3LOXz4MD169ADgv//9L9999x0pKSnk5+dz9dVX06BBA2655RYWL17MzJkzGT58uOfkPvH/JSYm2l2C40mG5iRDc7/++qvdJdSY65q8P7vhhhu47LLLGDFiBI899hh5eXnMnj2bL774wrO7acSIEbz00kusXLmS0aNH8/bbb/O///2vzBn0kZGRPP7443z66adcfvnl9OvXj3379hEREUFSUhLp6en89re/5ZlnnmHatGmytSCEED7Idbvr6/oSOl8QGxtbrTmj7733Xu69994zbrv66qs939erZ332mzp1KlOnTi3z+3/4wx/4wx/+cMZtEyZMqEbF/uHsEzJF9UmG5iRDc6V/E53IdU0+NDTUluf1tWPxJlq2bGl3Ca4wceJEu0twPMnQnGRorm3btnaXUGPO/XhSgaNHj9pdguPJpDTeMXfuXLtLcDzJ0JxkaC4zM9PuEmrMdU2+sLDQ7hIcr6ioyO4SXEE+LJmTDM1Jhuac3Fdc1+SFEEIIYXFdk6/KMfnSYVhF+Vq0aFEry/W33MePH293CY4nGZqTDM21bt3a7hJqzHVN/lxzoTds2FDmSz+H2hqeNjc394xBdNxu8+bNdpfgeJKhOcnQXHZ2tt0l1JjfNfk2bdqwd+9ecnJy/G7Lsqq8/SFIa01OTg579+6lTZs2Xl22L9uyZYvdJTieZGhOMjTn5GG7XXcJ3bkEBwcDsG/fPgoKCmyuxjedOHHC6yMHNmzYkLZt23ryF0IIUftc1+SrMnZ9cHCwNJtK7NixQ6aN9YLo6Gi7S3A8ydCcZGguJCTE7hJqzHW760+feEXUTIMGrvvsZwvJ0ZxkaE4yNOfkvuK6Jp+VlWV3CY63ZMkSu0twBcnRnGRoTjI0d+zYMbtLqDHXNXkhhBBCWFzX5EtnWxM117VrV7tLcAXJ0ZxkaE4yNNekSRO7S6gx1zX5oKAgu0twvH79+tldgitIjuYkQ3OSobmqnNDtq1zX5A8fPmx3CY6XlJRkdwmuIDmakwzNSYbmDh48aHcJNea6Ji+EEEIIi+uafP369e0uwfFkDAHvkBzNSYbmJENzTu4rym1Du0ZFRemUlBS7yxBCCCHqjFIqVWsddfbtrtuSl7mTzc2ePdvuElxBcjQnGZqTDM05OUPXNfmioiK7S3A8GVDIOyRHc5KhOcnQnJMzdF2TF0IIIYTFdcfk+/Tpo9evX293GY6WnZ0t4w14geRoTjI0Jxmac0KGPnFMXin1tlLqoFJq82m3hSqlvlRK/Vjyb4uS25VS6iWlVJpS6nulVJ+qPEd2dnZtle831q1bZ3cJriA5mpMMzUmG5pycYV3vrk8ERp5128PAEq31RcCSkp8BrgUuKvmKA16ryhOcOnXKK4X6s507d9pdgitIjuYkQ3OSoTknZ1inTV5rvQI4+/T3McDMku9nAmNPu/1/2rIaaK6UOr9OChVCCCFcwBcmGm6rtd5f8v2vQNuS79sDv5z2uIyS2/ZzFqVUHNbWPm3atGHGjBkA9O/fn9atW7Nw4UIAOnXqxPDhw0lISACgUaNGxMbGMn/+fM+whePGjSMtLY1NmzYBMGjQIJo1a8aiRYsA6Ny5M0OGDCExMRGAwMBAJk6cyNy5cz2X740fP57NmzezZcsWAKKjo2nQoIFnyseuXbvSr18/z3CTwcHBTJgwgdmzZ3vO4oyJiWHdunWeT5DDhg2jsLCQ5cuXAxAREUHPnj2ZM2cOAKGhoVx//fXMmjWLnJwcAGJjY1mxYgW7du0CYMSIEZw4cYJVq1YBEBkZSZcuXZg3bx6l2Y0dO5aioiJPhpMnT2bx4sXs2bMHgFGjRpGZmcnatWsB6Nu3Lx06dGDBggUAtGvXjtGjR5OQkEBxcTH16tVj8uTJJCcns2/fPgDGjBlDRkYGqamprn6d8vLyWLlyZa29TomJieTn57v6dcrLy+P48eOOfj/Z/TqFhIRQXFzs+PeTna9TZGQkGzdu9Pn3U3nq/MQ7pVQ4kKy17lny8zGtdfPT7j+qtW6hlEoGpmmtvy25fQnwN611pSPd9OzZU2/evLmyh4hz2LFjB926dbO7DMeTHM1JhuYkQwMFBfD88xzau5dWL71kdzWV8okT7ypwoHQ3fMm/pTMB7AU6nva4DiW3VerEiRNeL9DflH5qFmYkR3OSoTnJsIZSUqBfP3j4YY6vWgUOvRLNF5r8J8CtJd/fCiw47fZbSs6yHwgcP223vhBCCOF92dnwwAMwYAAcPAjz5rEkLg6UsruyGqnTY/JKqfeBoUArpVQGMBWYBnyglJoEpAN/KHn4Z8B1QBqQA/yxKs/RpEkTL1ftfyIiIuwuwRUkR3OSoTnJsBoWLYI774Tdu+GOO2DaNGjenIiVK+2urMZcNxhO79699YYNG+wuw9GOHz9OSEiI3WU4nuRoTjI0JxlWQWYm3H8/JCXBxRfDjBlwxRWeu52QoS8fk/cqmaDGXOmZq8KM5GhOMjQnGVZCa3j3XejeHT74AP75T9i48YwGD87O0BcuoRNCCCHq1q5d1q75L7+Eyy6DN98EFx7acN2WfIMG8rnFVGhoqN0luILkaE4yNCcZnqWwEJ57Dnr2hO++g1degW+/rbTBOzlD1x2Tj4qK0ikplV5KL4QQwh+tXw+33279+5vfwH//Cx07nvv3HECOyYsqmzVrlt0luILkaE4yNCcZAjk5MGUK9O8Pe/dax98XLKhyg3dyhq7bt11UVGR3CY5XOjykMCM5mpMMzfl9hl99ZV0Ot2sXTJ4MzzwDLVpUaxFOztB1W/JCCCGcYejQoQwdOrR2Fn74MMTGwjXXQP36sHSpdXJdNRu807numHzfvn116QQAomby8/Np1KiR3WU4nuRoTjI058sZljb4ZcuWeW+hWsP778Nf/gJHj1q76f/xD2jcuMaL9OUMS/nNMXkZu97cihUr7C7BFSRHc5KhOb/KMD0dRo2CmBi44AJITYX4eKMGD87O0HVNPi8vz+4SHK90ikZhRnI0Jxma84sMi4rgxRety+BWrIDp02HVKujVyyuLd3KGrjvxTgghhB/ZtMm6LG7dOrjuOnj1VQgLs7sqn+G6LXlfH1/YCUaMGGF3Ca4gOZqTDM25NsPcXHj0UYiKsiaUef99SE6ulQbv5AxdtyUvl9CZk/MavENyNCcZmnNlhkuXQlwcpKVZZ9A/9xy0bFlrT+fkDF23JX/y5Em7S3C8VatW2V2CK0iO5iRDc67K8OhRmDQJrrrKOov+q6/gnXdqtcGDszN0XZMXQgjhMlpbo9R17w4zZ8Lf/gbffw/Dhtldmc9z3e76wMBAu0twvMjISLtLcAXJ0ZxkaM7xGf7yC9x1l3W8vW9f+OILuPTSOi3ByRm6bks+ICDA7hIcr0uXLnaX4AqSoznJ0JxjMywqgpdfhh494Ouv4T//gdWr67zBg4MzxIVN/ujRo3aX4Hjz5s2zuwRXkBzNSYbmHJnh5s1w+eXw5z/D4MHWzw88ADZNJe7IDEu4rskLIYRwqFOnrCFo+/Sxzpx/9134/HNr9DpRI647Jt+wYUO7S3C8Nm3a2F2CK0iO5iRDc47JcMUK67K4HTvg5pvh+eehVSu7qwIclGE5XDdBTVRUlE5JSbG7DCGEEOcwdOhQmhYWkhwRATNmWFvsr78Ow4fbXZrj+M0ENYcPH7a7BMdLTEy0uwRXkBzNSYbmfDnDKzIzSVy3DhIS4MEH4YcffLLB+3KG5+K63fXFxcV2l+B4+fn5dpfgCpKjOcnQnE9muHcv3HMPT27dyo9Nm9Jq1Srr8jgf5ZMZVpHrtuSFEEL4qOJieO0167K4L77gtc6dubNPH59u8E4nx+RFGcXFxdSrJ5//TEmO5iRDcz6T4bZt1mxxK1daI9W98QZDJ00CYNmyZfbWdg4+k2El/OaYfFZWlt0lON7ixYvtLsEVJEdzkqE52zPMy4MnnrAGsdm2DRIT4csv4cIL7a2rGmzP0IDrmnxeXp7dJTjenj177C7BFSRHc5KhOVszXLUKeveGxx+HceOsJn/rraAUSUlJrF69muXLlxMeHk5SUpJ9dZ6Dk9dD1zV5IYQQNsvKgrvvtkaty86Gzz6D996DkuvNk5KSiIuL82yUpaenExcX59ON3qlcd0w+MjJSb9q0ye4yHG3v3r20b9/e7jIcT3I0Jxma83aGQ4cOrfT+wYcO8ZcffyQ0P5+P2rfn7QsuILd+/TMes3r16nL3ugYEBDBw4EDjGr19jN8J66HfHJMvKCiwuwTHy8zMtLsEV5AczUmG5uoqw9C8PJ7YsoX4LVvIatiQu3v35r9dupRp8FDxYVVfPdzq5PXQddfJZ2dn212C461du5ZLbZjpyW0kR3OSoTlvZ1hmK7m4GN56C/76V2vs+aee4sKHHuL1SoYYDw8PJz09vcztYWFhPnmmvZPXQ9dtyQshhKgjO3bAlVdaY8737m2NWPfII3COOUTi4+MJDAw847bAwEDi4+Nrs1q/5LomHxQUZHcJjtdXBqbwCsnRnGRorlYyzM+H+HiIjITvv7eGpf36a7jooir9ekxMDDNmzCAgIACwtuBnzJhBTEyM92v1Aievh67bXS+z0Jnr0KGD3SW4guRoTjI05/UM16yByZOtOd7/8AeYPh3OO6/ai4mJieHNN98EfH8wHCevh67bkj927JjdJTjeggUL7C7BFSRHc5KhOa9leOIE/PnPcNllcOwYfPIJzJlTowbvNE5eD13X5IUQQnhZcjJERMArr8A998DWrfCb39hdlagC1+2ub9Sokd0lOF67du3sLsEVJEdzkqE5owwPHID77rO22CMirBHsvHAdu9M4eT103WA4MkGNEEIY0toaY/7BB60R6/7xD5gyBby8EVU6sI6vH5N3Ar8ZDOfQoUN2l+B4CQkJdpfgCpKjOcnQXLUzTEuDq6+G226DSy6xzp7/+9+93uCdxMnroeuavNv2TNihuLjY7hJcQXI0Jxmaq3KGBQUwbZrV2FNT4Y03YOlS6Natdgt0ACevh647Jq+UsrsEx/P1eZOdQnI0Jxmaq1KGKSnWZXGbNlmzxb38Mpx/fu0X5xBOXg995pi8Uup+YDKggR+APwLnA7OBlkAqcLPWOr+y5cgxeSGEqKLS4+2l17r/978wdmydPb0ck/cenz4mr5RqD/wZiNJa9wTqAxOAfwMvaK27AEeBSeda1vHjx2uzVL+QnJxsdwmuIDmakwzNVZjhF19YZ8y/8ALccYd1WVwdNngncfJ66BNNvkQDoIlSqgEQCOwHrgLmltw/Exh7roXk51e6oS+qYN++fXaX4AqSoznJ0FyZDDMzYeJEuPZaCAyEb7+FV1+FkBB7CnQAJ6+HPtHktdZ7geeAPVjN/TjW7vljWuvCkodlAL49oa8QQvgqreF//4Pu3eGDD+Dxx2HDBhg82O7KRC3yiRPvlFItgDHABcAx4ENgZDV+Pw6IA2jbti0zZswAoH///rRu3ZqFCxcC0KlTJ4YPH+65HKJRo0bExsYyf/58Dh48CMC4ceNIS0tj06ZNAAwaNIhmzZqxaNEiADp37syQIUNITEwErJmTJk6cyNy5czly5AgA48ePZ/PmzWzZsgWA6OhoGjRowJIlSwDo2rUr/fr1IykpCYDg4GAmTJjA7NmzycrKAqxxndetW8fOnTsBGDZsGIWFhSxfvhyAiIgIevbsyZw5cwAIDQ3l+uuvZ9asWeTk5AAQGxvLihUr2LVrFwAjRozgxIkTrFq1CoDIyEi6dOnCvHnzAGjTpg1jS3bXlWY4efJkFi9ezJ49ewAYNWoUmZmZrF27FrAmbujQoYNn2Md27doxevRoEhISKC4upl69ekyePJnk5GTPp+ExY8aQkZFBamqqq1+nwsJCVq5cWWuvU2JiomfPlVtfp8LCQo4fP+7o95Pdr1OrVq0oTktj329/S4dt2zjQpQttV6xg/s6dHJw50yuvE9Ts/QRw8OBBz98bX32dBg4cyMaNG33+/VQenzjxTil1AzBSaz2p5OdbgMuAG4DztNaFSqnLgMe11iMqW1aPHj301q1ba71mN0tNTXX0rEu+QnI0JxkaKiwk46GH6DBjBjRoAP/+t3X83UfOFnfKiXdOWA99+sQ7rN30A5VSgcq6Bm4YsBVYClxf8phbgXPOEpCdnV1rRfqL0k+awozkaE4yNLB+PQwYQIfp02H4cOvEuj/9yWcavJM4eT30iVdba70G6wS79ViXz9UDZgB/Ax5QSqVhXUb3lm1FCiGEE+TkwF//Cv37w759fBkXBx9/DA6eLlXUnE8ckwfQWk8Fpp518y6gf3WWExQU5LWa/FX//tWKXFRAcjQnGVbTV19Zu+N37YK4OEZv2cKvqalcI4OEGXHyeugTW/Le1LBhQ7tLcLzWrVvbXYIrSI7mJMMqOnwYbr0VrrnGOva+fDm88QYnGzSQv4le4OT10HVN/tixY3aX4HilZ3sKM5KjOV/OcOjQoZ4Tx2yjNbz3Hlx8sfXv3/9uDU07ZIjnIaVnv4ua8+X18Fx8Zne9EEKIati92zqR7osvYMAAePNNa3IZIU7juiYfEBBgdwmO16lTJ7tLcAXJ0ZxkWI6iInjpJWurvV496/u77oL69ct9uC//TfT1S+dKOXk9dF2TLx1gQdTc8OHD7S7BFSRHc5LhWTZtgttvh3XrYPRoazjajh0r/ZXQ0NA6Ks69nLweuu6YfGZmpt0lOF7piEvCjORoTjIskZsLjzwCfftCejrMng2ffHLOBg+wf//+OijQ3Zy8HrpuS14IIVzl66+ty+LS0uC22+DZZ0G2zkUVuW5Lvp6M5mSsUaNGdpfgCpKjOb/O8MgRmDQJhg2zfv76a3jrrWo3eCXXyBtz8nroE2PXe1NUVJROSUmxuwwhhMvV2rjrWluzxP35z9b171OmwD/+AU2a+E6Nwuf4+tj1XiPXyZubP3++3SW4guRozu8y3LMHfvMbmDABwsIgNRWeeqpGDb7UoUOHvFigf3Lyeui6Jl9QUGB3CY5XOq2hMCM5mvObDIuK4OWXISICli6FF16A776DyEjjRcvfRHNOXg/lxDshhLDT5s0weTKsWQMjR8Jrr0F4uN1VCZdw3ZZ8ixYt7C7B8caNG2d3Ca4gOZpzdYanTlkD2vTuDT/9BElJ8NlnXm/wrVq18ury/JGT10PXNfm8vDy7S3C8tLQ0u0twBcnRnGszXLHC2hUfHw8xMbB9O9x0E9TCmfC5ubleX6a/cfJ66Lomn5OTY3cJjrdp0ya7S3AFydGc6zI8dgzi4iA6GgoKYPFiSEyEli1r7Smzs7Nrbdn+wsnroeuavBBC+BytYd486N4d3n4b/vpX61j8NdfYXZlwOdedeNe0aVO7S3C8QYMG2V2CK0iO5lyR4d69cPfdsGCBdfx94ULo06fOnl7m8zDn5PXQdVvy9SuYiUlUXbNmzewuwRUkR3OOzrC42DpTvnt3a7f8s8/C2rV11uCTkpJYvXo1KSkphIeHk5SUVCfP60ZOXg9d1+SPHz9udwmOt2jRIrtLcAXJ0ZxjM9y6FYYMsaaAHTjQ2jX/0EPQoG52niYlJREXF+c5ETk9PZ24uDhp9DXk2PUQF+6uF0KI2la6lZyXl0d4eDjx8fHExMRAXh48/bQ1Sl1wMMycCTffXK2z5kuHojVRWtvpcnJymDRpEm+++abRsmWIXGdxXZMPCAiwuwTH69y5s90luILkaM4XM6xoK5mdO4n58EPYts26LO6FF6B1a1tqrOhSYrnEuGZ8cT2sKtdNUNO3b1+dmppqdxmOlp+f7+hZl3yF5GjO2xnW1lYyQBiwOiCA57t2Za3BVLDe2FIODw8nPT29zO1hYWHs3r3bePn+xgnvZb+ZoEYmYzCXmJhodwmuIDma88UMK9oa3gP8sV8/owbvLfHx8QQGBp5xW2BgIPHx8TZV5Gy+uB5Wlet21wshREVqcyu5U1gYn3/zjfHyvSEmJgaASZMmkZeXR1hY2P8/b0D4FddtycsldObO3gIQNSM5mvO5DIuLiR82jLOr8sWt5JiYGAYOHMjFF1/M7t27pcEb8Ln1sBpcd0w+KipKp6Sk2F2GEMJtduyA22+Hb74hqUcP4tLSyMnP9+mt5NJzEOSMePfzm2PyR48etbsEx5s7d67dJbiC5GjOJzLMz4cnn4Revazr3d9+m5jNm+l32WVER0f7/FZyZmam3SU4nk+shzXkumPyhYWFdpfgeEeOHLG7BFeQHM3ZnuHq1dZc71u2wIQJ8OKL0LatvTVVk/xNNGf7emigylvySqlVSqmblVJyIboQwt1OnIB774VBgyArC5KT4f33HdfghajO7vp8YCawTyn1vFLq4lqqyUioD1y+4nTjx4+3uwRXkBzN2ZJhcjL06AH//a/V6LdsgVGj6r4OL2lt04A8buLk93KVm7zWeijQA6vR3wJsUUotU0qNV0o1rKX6qi03N9fuEhxv8+bNdpfgCpKjuTrN8NdfYfx4+M1voHlz+O47mD4dHDw5Cch88t7g5PdytU6801pv11o/ALQHYoH6wHtAhlJqmlLK9rH/pMmb27Jli90luILkaK5OMtTamuO9e3drOth//QtSU2HAgNp/7jqQk5NjdwmO5+T3co3Ortda52mt3wXuA74BWgNTgJ1KqQ+VUud5sUYhhKgdP/4Iw4bBpEkQGQnffw+PPQY+PoSpEFVV7SavlGqilLpNKbUWWAe0wWr27YA/AYMA2+YzdPK8v74iOjra7hJcQXI0V2sZFhRYs8VdcgmsXw8zZsDXX0PXrrXzfDYKCQmxuwTHc/J7ucqX0CmlLgHuAGKAIGAB8Det9dLTHvamUupX4EOvVlkNqhpTOoryNaijOa/dTnI0VysZrl1rDWrz/fdw/fXw0ktw/vnefx4fIX8TzTn5vVydLflNwFjgRSBMa33DWQ2+VBrwnXlpNZOVlWXXU7vGkiVL7C7BFSRHc17N8ORJuP9+uOwyOHzYOv7+4YeubvAAx44ds7sEx3Pye7k6H0+uBxZorYsqe5DWehtwpVFVQgjhTZ9/Dn/6E+zZY/379NMQHGx3VULUuupcQvfRuRq8L2jcuLHdJTheVxcel7SD5GjOOMODByEmBq67DoKC4Ntvrevf/ajBN2nSxO4SHM/J72XXjV0fFBRkdwmO169fP7tLcAXJ0VyNM9QaZs60LoubOxeeeMI6wW7QIO8W6AByMrI5J7+XXdfkDx8+bHcJjpeUZNvFEa4iOZqrUYY//QTDh0NsrNXkN26Ef/4TAvxzRO6DBw/aXYLjOfm97LomL4TwU4WF8Mwz1mVxa9fCa6/BihVWoxfCT/nMdQFKqeZAAtAT0MBtwA5gDhAO7Ab+oLWudC7Z+vXr12aZfiHYj45X1ibJ0czQoUM5ePAgcXFx537w+vXWbHEbNsDvfgcvvwzt29d+kT5u2bJlzJ492+4yHM/J72Vf2pKfDnyhtb4YiAS2AQ8DS7TWFwFLSn6ulExQY27ChAl2l+AKkqO5Nm3aVP6A7Gx46CHo188ae37ePPjoI2nwp5H10JyTM/SJJq+UCgGGAG8BaK3ztdbHgDFYE+JQ8u/Ycy3LyfP++gr55O8dkqO5So8nL15s7Zr/z3+swW22boXf/77Oalu2bBnLli2rs+erKVkPzTk5Q59o8sAFQCbwjlJqg1IqQSkVBLTVWu8vecyvwDkncy4q8vmr/HyeDCjkHZKjuXLfz4cOwS23wIgR1hjzK1bA669bM8eJMmQ9NOfkDH3lmHwDoA9wr9Z6jVJqOmftmtdaa6WULu+XlVJxQBxAixYtmDFjBgD9+/endevWLFy4EIBOnToxfPhwEhISAGjUqBGxsbHMnz/fs8Uwbtw40tLS2LRpEwCDBg2iWbNmLFq0CIDOnTszZMgQEhMTAQgMDGTixInMnTvXsxdh/PjxbN682TNzUXR0NA0aNPCMmtS1a1f69evnOWMzODiYCRMmMHv2bM/KFBMTw7p169i5cycAw4YNo7CwkOXLlwMQERFBz549mTNnDmAdprj++uuZNWuWZ9ap2NhYVqxYwa5duwAYMWIEJ06cYNWqVQBERkbSpUsX5s2bB1i7RseOHcvhw4c9GU6ePJnFixezZ88eAEaNGkVmZiZr164FoG/fvnTo0IEFCxYA0K5dO0aPHk1CQgLFxcXUq1ePyZMnk5yczL59+wAYM2YMGRkZpKamuvp1OnToECtXrqy11ykxMZH8/HzXvk779++noKCA48ePW6/T5s10WbOGIfPnU//kSdaPGsXGa6/lwqIi+mVn++z7ye7XKSsri+LiYse/n+x8nfLz89m4caPPv5/Ko7Qut2/WqZJZ61ZrrcNLfr4Cq8l3AYZqrfcrpc4Hlmmtu1W2rD59+uj169fXdsmulp2dLeMNeIHkaGbo0KEUFRXxzTffwM8/WyPVLVoEAwdCQgJERNhdoiPIemjOCRkqpVK11lFn3+4Tu+u11r8CvyilShv4MGAr8Alwa8ltt2JNilOp7OzsWqnRn6xbt87uElxBcjSXk5UFzz8PPXvCypXWWfPffisNvhpkPTTn5Ax9osmXuBdIUkp9D1wKPAVMA65RSv0IXF3yc6VOnTpVmzX6hdJdZcKM5Gimy8mTzNyxAx580JrzfetWuOcekMtkq0XWQ3NOztBXjsmjtd4IlNnVgLVVL4TwFzk58MQTvJGaypH69eGDD6wpYWXKVCGqzWeavLc4edACXzFsmHyu8gbJsQaWLIE77oCffuKL887jhXbt+OqGG+yuytFkPTTn5Ax9aXe9V/jCiYROV1hYaHcJriA5VsPhw/DHP8LVV0O9erB0Kc9260aW7Jo3JuuhOSdn6Lomf+LECbtLcLzSy1WEGcmxCrSG2bOt8eVnzYJHH4VNm2DoUACOHz9ub30uIOuhOSdn6Lrd9UIIh9izx7os7rPPrGFpv/oKevWyuyohXMV1W/JNmjSxuwTHi5DLk7xCcqxAURFMnw49esDy5fDii/Ddd+U2+MDAwLqvz2VkPTTn5AxdtyUvTd5cz5497S7BFXw5x6Elu8PrfOz1H36wZotbuxauuw5efRXCwip8uK8PQOIEvrweOoWTM3TdlrxMUGOudMhIYUZyPM2pU/DYY9CnjzV63XvvQXJypQ0eIDMzs44KdC9ZD805OUPXNXkhhI9ZtszaFf/UUzBxImzbBjfeWOl170lJSaxevZqdO3cSHh7uGe9cCFE9rmvyDRq47ghEnQsNDbW7BFfw+xyPHrWmgL3ySigutk6se+cdaNmy0l9LSkoiLi6OvLw8ANLT04mLi5NGX0N+vx56gZMz9IkJarwpKipKp6Sk2F2GED6tNo7Jly4TrYk+dIj7fvyRkIIC5nTsyMywMPKqeM376tWrPQ3+dAEBAQwcONCoRifM/y5ETfj0BDXeJMfkzc2aNcvuElzBH3NsnZdH/JYtPLF1K5kBAdzRpw8zOneucoMHym3wld0uKueP66G3OTlD1+3bLioqsrsExyudl1mY8asci4tZdsMN8Mgj1iVy//kP3f78ZxJqcPgsPDyc9PT0MreHhYXJlngN+NV6WEucnKHrtuSFEHVsyxa4/HJrhrjLLoPNm+GBB6CG58fEx8eXuT4+MDCQ+Ph4b1QrhF9x3TH5vn376tTUVLvLcLT8/HwaNWpkdxmO58s5euWYfF4exMfDtGkQHGwNahMT45XZ4pKSkpg0aRJ5eXmEhYURHx9PTEyM8XL9kS+vh07hhAz95pi8jF1vbsWKFXaX4AquzvHbb+HSS+HJJ2HCBNi+3bo8zkvTwcbExDBw4EB69erF7t27pcEbcPV6WEecnKHrmrycnGNu165ddpfgCq7M8fhxuPNOuOIKa4CbRYvgf/+DVq1q5elOnTpVK8v1J65cD+uYkzN0XZMXQtSSjz+2Zot780148EHr2Pvw4XZXJYSohOvOrg8JCbG7BMcbMWKE3SW4gmty3LfPOqnu44+tXfSffgp9+9bJU7do0aJOnsfNXLMe2sjJGbpuS14uoTMn5zV4h+NzLC6G11+3tt4//xz+/W9rYpk6avAg72dvcPx66AOcnKHrmvzJkyftLsHxVq1aZXcJruCrOZaOC798+fKKx4Xfvh2io6353vv1s3bNT5kCDRvWaa1ZWVl1+nxu5KvroZM4OUPXNXkhRMXOOS58fj783/9BZCRs3WqNNf/ll3DhhTZWLYSoKdcdkz97EA1RfZGRkXaX4AreztEzNryB8saFz8nJYdKkSax67jke2rmTC3Jy+KpNG1658EKOJSZCYmKVll0bo9HJfPLm5P1szskZuq7JBwQE2F2C43Xp0sXuElzBF3Os6BLT/Lw8Xt64kYMBAfytZ0/WnGOmuLrSpEkTu0twPF9cD53GyRm6rskfPXrU7hIcb968ecTFxdldhuN5O0dvbClXNC58J6DeX/7CeU8+yb+bNjV+Hm85dOiQ3SU4nryfzTk5QzkmL4QfKXdceKWIf+IJeOEF8KEGL4Qw57om37COz/51ozZt2thdgiv4Yo4xN93EjBtvpCOggLDmzZmRmEjMP/9pd2nlkvezOV9cD53GyRm6boKaqKgonZKSYncZQvienTshLg6WL2dDSAj/6dqVWWvX2l1VhbwyiY4QfsJvJqg5fPiw3SU4XmIVz6YWlfOZHPPzrdnievWCTZsgIYH7IyPJcMCVKL/++qvdJTiez6yHDubkDF3X5IuLi+0uwfHy8/PtLsEVfCLHNWsgKgr+/ncYMwa2bYNJk7w2W1xtWrZsGQ888IDdZTieT6yHDufkDF3X5IUQwIkTcN99cNllcPQofPIJzJkD551nd2VCiDokx+RFGcXFxdSrJ5//TNmW48KF1nC0GRlw993Wrvrg4DMe4pTj3bIumpMMzTkhQ785Ji9jXZtbvHix3SW4Qp3neOAA3HgjjB5tNfWVK+Hll8s0eCeRddGcZGjOyRm6bjCcikb0ElW3Z88eu0twvKFDh7J//3527NhR+0+mtTX07IMPQna2Nfb83/4GjRrV/nPXMlkXzUmG5pycoeuavBB+JS0N7rwTliyBK66AGTPg4ovtrkoI4SNct7u+efPmdpfgeKNGjbK7BFcIDQ2tvYUXFFjzu19yCaxbZ837vmyZ6xq8rIvmJENzTs7QdU2+oKDA7hIcLzMz0+4SXKHW1sWUFGuO94cfhuuusy6Lu+MO8PETg2pC1kVzkqE5J2four8K2dnZdpfgeGt9eBQ0Jzlx4oR3F5idbR13HzAAMjPh449h3jxo1867z+NDZF00Jxmac3KGckxeCCdYtMg69r57t3V53NNPQ0iI3VUJIXyc67bkg4KC7C7B8fr27Wt3Ca7Q1BszumVmwsSJMHIkNG4M33wDr77qNw1e1kVzkqE5J2fouiYvs1aZ69Chg90luEJAQEDNf1lrePdd6N4dPvgApk6FjRvh8su9Vp8TyLpoTjI05+QMXdfkjx07ZncJjrdgwQK7S3CFGk+W9PPPMGIE3HILdOsGGzbA44+DyYeGsyxbtsznR7sDWRe9QTI05+QMfarJK6XqK6U2KKWSS36+QCm1RimVppSao5Ry/ugeQlSksBCeew4iImD1avjvf63d8xERdlcmhHAon2rywH3AttN+/jfwgta6C3AUmHSuBTRywShfdmvn4rO161K11sUNG6yz5v/6Vxg+HLZuhbvucuVlcdUh66I5ydCckzP0mb8gSqkOwCggoeRnBVwFzC15yExg7LmWE+InJyTVptGjR9tdgiu0bNny3A/KyYEpU6zr3vftg7lzrUvjHHwM0JtkXTQnGZpzcoY+0+SBF4EpQOmE8C2BY1rrwpKfM4D251rIoUOHaqU4f5KQkGB3CY6WlJTE6tWrWb58OeHh4SQlJZX/wK++skase/ZZuO02a+t93DhHzPVeV2RdNCcZmnNyhj5xnbxSajRwUGudqpQaWoPfjwPiAFq0aMGMGTMA6N+/P61bt2bhwoUAdOrUieHDh3tesEaNGhEbG8v8+fM5ePAgAOPGjSMtLY1NmzYBMGjQIJo1a8aiRYsA6Ny5M0OGDCExMRGAwMBAJk6cyNy5czly5AgA48ePZ/PmzWzZsgWA6OhoGjRowJIlSwDo2rUr/fr18/zxDw4OZsKECcyePdszi15MTAzr1q1j586dAAwbNozCwkKWL18OQEREBD179mTOnDmANYTq9ddfz6xZs8jJyQEgNjaWFStWsGvXLgBGjBjBiRMnWLVqFQCRkZF06dKFefPmAdCmTRvGjh1LZmamJ8PJkyezePFizwQNo0aNIjMz0zM4RN++fenQoYPnxJR27doxevRoEhISPNMzTp48meTkZPbt2wfAmDFjyMjIIDU11XWvU3Z2No8++qhnoqT09HRuv/12vv76awYMGGC9TldeyU9jx3Lht99yrG1bgr78kuXArg8/rNbrlJiYSH5+vqtfp0OHDnH8+HFHv5/sfp2ysrIoLi525PvJV16n/Px8Nm7c6PPvp/L4xHzySqmngZuBQqAxEAx8DIwAztNaFyqlLgMe11qPqGxZ4eHhevfu3bVcsbslJCQwefJku8uoU6Xzq5tavXp1uTMhBgQEMHDAAIYdPMg9P/1Es8JC3u/YkXfDwsiv4nF3J5wN723+uC56m2RozgkZVjSfvE80+dOVbMk/pLUerZT6EJintZ6tlHod+F5r/Wplvx8VFaVTUlLqoFLhJt5q8qVbHOX5LjSUgUeOsLVZM57t2pWfqzlYjj82eSFE1VTU5H1id30l/gbMVkr9C9gAvHWuXzh+/HitF+V2ycnJjj7RpCa81UDDw8NJT08vc3uYUgzMy4Pp0+lx9928U7++V57P7fxxXfQ2ydCckzP0pRPvANBaL9Najy75fpfWur/WuovW+gatddn9oGcpPaYiaq70GJKovvj4eAIDA8+4LRCIj4y0Tqz7859BGnyVybpoTjI05+QMfa7JC+FkMTExzHjlFVrWr48CwurVY8Y99xCzfj106mR3eUIIP+Nzx+RNRUZG6tIzD0XNHDhwgLZt29pdhjMtXQpxcZCWRnLr1ozevh1CQ+2uyrFkXTQnGZpzQoYVHZN33ZZ8QUGB3SU4XkZGht0lOM/RozBpElx1FWjN/b16MbVjR2nwhmRdNCcZmnNyhq5r8tnZ2XaX4Hil13GKKtDamiWue3eYORMefhh++IENLVpw8uRJu6tzPFkXzUmG5pycoa+fXS+E7/rlF2t8+eRkiIqCL76ASy+1uyohhPBw3ZZ8UFCQ3SU4Xv/+/e0uwbcVFcHLL0OPHvD11/D889ascWc1+GbNmtlTn4vIumhOMjTn5AxdtyXfsGFDu0twvNatW9tdgu/avBluv91q6iNHwmuvQXh4uQ+VddGcrIvmJENzTs7QdVvyx44ds7sExysdS1mc5tQp+Mc/oE8fSEuDpCT47LMKGzzgGdNb1Jysi+YkQ3NOztB1W/JCeN2KFdZlcTt2wC23wH/+A61a2V2VEEKck+u25AMCAuwuwfE6yaAtlmPH4I47IDoa8vNh8WLrDPoqNnhZF83JumhOMjTn5AxdNxiOTFBjrnSqRL/20Udwzz1w4AA88AA8/jhU46TO0glvZFIZM7IumpMMzTkhQ78ZDCczM9PuEhyvdD5jv7R3L/zudzBuHJx3HqxbB88+W60GX2r//v21UKB/8et10UskQ3NOztB1TV6IGikuts6U79EDFi2CZ56BtWutE+2EEMKhXHfina/vUnGCRo0a2V1C3dq61TqxbuVKuPpqeOMN6NzZaJHLli0jMTHRO/X5Mb9bF2uBZGjOyRnKMXnhv/Ly4Omn4amnoFkzeOEFuPlmUMruyoQQolr85pi8XCdvbv78+XaXUPtWroTeveGJJ+APf4Dt263L47zY4P0ix1omGZqTDM05OUPXNXmZhc7cwYMH7S6h9hw/bo03f/nlkJNjDWgzaxbUwohWrs6xjkiG5iRDc07O0HVNXogKLVhgnVj3xhtw//3WELXXXmt3VUIIUWtcd0z+0ksv1Rs3brS7DEc7fPgwLVu2tLsM79m/H+69F+bNg169ICEB+vWr9ad1XY42kAzNSYbmnJCh3xyTz8vLs7sEx0tLS7O7BO8oLoYZM6y53hcutE6yS0mpkwYPLsrRRpKhOcnQnJMzdF2Tz8nJsbsEx9u0aZPdJVRq6NChnhHlKrRjB1x5pTUsbZ8+8P338PDDUIczw/l6jk4gGZqTDM05OUPXNXnh5/Lz4V//snbL//ADvP02LFkCF11kd2VCCFHnXDcYTtOmTe0uwfEGDRpkdwk1s3q1Ndf75s0wfjxMnw5t29pWjmNz9CGSoTnJ0JyTM3Tdlnz9+vXtLsHxmjVrZncJ1XPiBPz5zzBokHWJXHIyzJ5ta4MHB+bogyRDc5KhOSdn6Lomf/z4cbtLcLxFixbZXULVJSdDRAS88op1Bv2WLTBqlN1VAQ7L0UdJhuYkQ3NOztB1u+uFf2iRn2/tkv/gA+jZ0/p34EC7yxJCCJ/iuiYfEBBgdwmO19lwcpZapTXX7t/Pn3btgvXrrZPs/vpX8MEJJHw6R4eQDM1JhuacnKHrBsPp27evTk1NtbsMR8vPz/fNWZd+/NG6JG7pUjaGhHDpmjXQrZvdVVXIZ3N0EMnQnGRozgkZ+s1gOIcOHbK7BMfzuSlSCwqsgWx69YL163m2a1fuj4z06QYPPpijA0mG5iRDc07O0HVNXrjMunUQFQWPPgqjR5MUH88/0tNZtmIF4eHhJCUl2V2hEEL4LNcdk5dL6MwFBgbaXQKcPAn/+Ae89BKcfz7Mn0/SyZPExcV5hi5OT08nLi4OgJiYGDurLZdP5OhwkqE5ydCckzN03TH5qKgonZKSYncZ4iznHIb2NP0PH+aBH3+kTV4eC9q1480LLiCnQQNWr15d7twEAQEBDDQ8s37ZsmVGvy+EEHbym2PyR48etbsEx5s7d64tz9s8P5/Htm3jmc2bOVW/PvdeeinTL7qInAbWDqeKJh/y1UmJ7MrRTSRDc5KhOSdn6Lrd9YWFhXaX4HhHjhzx+jIr3VLWGv73P3jgAWs3/eOPE/7ww/z3rMshw8PDSU9PL/PrYWFhPrklXhs5+hvJ0JxkaM7JGbpuS144zE8/wfDhEBtrTQm7cSNMnQrljHcQHx9f5thYYGAg8fHxdVOrEEI4jOuOyffu3Vtv2LDB7jIc7fjx44SEhNTukxQWwgsvWA29YUP4978hLg7qVf65MykpiUmTJpGXl0dYWBjx8fE+edId1FGOLicZmpMMzTkhQ785Jp+bm2t3CY63efPm2n2C9euhf3+YMgVGjICtW+HOO8/Z4ME6i37gwIFER0eze/dun23wUAc5+gHJ0JxkaM7JGUqTF2Vs2bKldhack2MNQduvH/z6K8ybBx9/DO3b187z2azWcvQjkqE5ydCckzN03Yl3wkd9+aU1JO3PP1v/TpsGzZvbXZUQQria67bknTzvr6+Ijo723sIOHYJbbrFOrmvUCJYvh9df94sG79Uc/ZRkaE4yNOfkDF3X5JVSdpfgeA0aeGEHj9aQlGSdMT97tjV63caNMGSI+bIdwis5+jnJ0JxkaM7JGbquyWdlZdldguMtWbLEbAG7d8O118LEidCli3Wi3f/9HzRu7JX6nMI4RyEZeoFkaM7JGfpEk1dKdVRKLVVKbVVKbVFK3Vdye6hS6kul1I8l/7awu1ZRicJCeP55iIiAlSvh5Zfh22+hZ0+7KxNCCL/kE00eKAQe1Fr3AAYCdyulegAPA0u01hcBS0p+rlRjP9tarA1du3at/i9t2gSXXQYPPgjDhlmXxd1zD/jxhEE1ylGcQTI0Jxmac3KGPtHktdb7tdbrS74/AWwD2gNjgJklD5sJjD3XsoKCgmqpSv/Rr1+/qj84Nxcefhj69oVffoEPPoAFC6Bjx9or0CGqlaMol2RoTjI05+QMfe5sAqVUONAbWAO01VrvL7nrV6BtBb8TB8QBtGjRghkzZgDQv39/WrduzcKFCwHo1KkTw4cPJyEhAYBGjRoRGxvL/PnzOXjwIADjxo0jLS2NTZs2ATBo0CCaNWvGokWLAOjcuTNDhgwhMTERsIZVnThxInPnzvWMbzx+/Hg2b97subYyOjqaBg0aeI7rdO3alX79+nnmQg8ODmbChAnMnj3bc05BTEwM69atY+fOnQAMGzaMwsJCli9fDkBERAQ9e/Zkzpw5AISGhnL99dcza9YscnJyAIiNjWXFihXs2rULgBEjRnDixAlWrVoFQGRkJF26dGHevHkAtGnThrFjx/LSSy/RsmVLACZPnszixYvZs2cPAKNGjSIzM5O1a9fSbvt2rv7wQxpnZLB98GBWjxtHqyZNGK0UCQkJFBcXU69ePSZPnkxycjL79u0DYMyYMWRkZJCamlqj1+mmm25i3LhxrFmzxqdfp0OHDhEdHV1rr1NiYiL5+fnnfJ0A+vbtS4cOHViwYAEA7dq1Y/To0bX6OoH5++nQoUPcfffdjn4/2f06ZWVl8cADD8jfPYPXKT8/n8svv9zn30/l0lr7zBfQFEgFfl/y87Gz7j96rmV06tRJCzNvvPFG5Q84fFjrP/5Ra9D6oou0Xrq0TupymnPmKM5JMjQnGZpzQoZAii6nJ/rE7noApVRDYB6QpLX+qOTmA0qp80vuPx84eK7l1PfjY8DeEhwcXP4dWluXw3XvDu++C488Yh2Lr8Zc8f6kwhxFlUmG5iRDc07O0CcmqFHWxe0zgSNa67+cdvuzwGGt9TSl1MNAqNZ6SmXLioqK0ikpKbVar5sNLWnYZaZu3bMH7roLFi60hqVNSIBeveq8PiGEEGX5+gQ1g4GbgauUUhtLvq4DpgHXKKV+BK4u+blSTp7311eUHv8BoKgIXnoJevSAZcvgxRfhu++kwVfB7Nmz7S7B8SRDc5KhOSdn6BMn3mmtvwUqGqpuWHWWVVRUZF6Qn/Nk+MMPcPvtsGaNNbjNa69BWJi9xTmIDMxkTjI0Jxmac3KGvrIlL3xIo+Ji+PvfoU8f2LUL3nvP2k0vDV4IIRzFdU2+9NIvUTORx46xcO9eiI+3hqXdtg1uvBFkToBq8+W57p3CbRkuWbKEiy++mMDAQK688krS09NrvKzZs2fTrVs3QkJCaNOmDbfeemu5W5ylGX799df06dOH4OBgOnfu7LnUuFRmZiY33XQTISEhtGjR4ozsp0yZQseOHQkODiYsLIynnnrqjN+Ni4ujW7du1KtXz3OZnZ2GDRuGUorCwsIKH/PBBx/QvXt3mjVrRo8ePZg/f77nvsTEROrXr0/Tpk1p2rQpDzzwwBnnKYWHh9OkSRPP/cOHD/fcN3PmTPr27UtwcDAdOnRgypQpldZR68o75d7JX127djW/FsEfHTmi9eTJWoPe06iR1l9+aXdFjrdULi005qYMMzMzdXBwsP7ggw90bm6ufuihh/SAAQNqvLw9e/bozMxMrbXWJ06c0DfddJO+9957yzxu6dKlOj8/XwcHB+vXX39dFxcX67Vr1+qgoCC9ceNGz+Muv/xyff/99+tjx47p/Px8vX79es9927dv1ydPntRaa52RkaF79Oih582b57n/lVde0V999ZXu27evfuedd6r8f1i6dKmOjo6u5v+8crNmzdJXXHGFBnRBQUG5j8nIyNANGzbUn332mS4uLtbJycm6SZMm+sCBA1prrd955x09ePDgM+o8XVhYmP6ygr+Rr776ql6xYoXOy8vTGRkZuk+fPvrpp5/2zn+uEvj6JXTecurUKbtLcBatYe5c68S6d97hvY4dGd2pE1x9td2VOV7pgB6i5pyWYXh4OE8//TQ9evSgRYsW/PGPf/T8Tfroo4+IiIjghhtuoHHjxjz++ONs2rSJ7du31+i5OnbsSKtWrTw/169fn7S0tDKP27lzJ0eOHCErK4ubb74ZpRT9+vWje/fubN26FYDFixfzyy+/8OyzzxISEkLDhg3p3bu3ZxndunU7YzTRevXqnfFcd999N8OGDbN9WPHjx4/zxBNP8Mwzz1T6uIyMDJo3b861116LUopRo0YRFBTETz/9VO7jq7Me/ulPf+KKK66gUaNGtG/fnpiYGFauXFmt/4c3ua7Ji2rIyICxY+GGG6BdO5KeeILbDh7k+7Q0wsPDPSNTCSGqLikpiUWLFvHTTz+xc+dO/vWvfwGwZcsWIiMjPY8LCgriwgsv9IwQN23aNEaPHl2t5/r2228JCQmhWbNmzJs3j7/85S/lPq5t27bceOONvPPOOxQVFfHdd9+Rnp7O5ZdfDsDq1avp1q0bt956Ky1btqRfv36eUeZKTZs2jaZNm9KhQweys7O56aabqlVrXXj00Uf505/+xHnnnVfp46KioujevTuffPIJRUVFzJ8/n4CAAHqddtXQhg0baNWqFV27dmXhwoVldrnHxMTQunVrhg8f7hmBrjwrVqwgIiLC7D9morzNeyd/9ezZ02u7P1yrqEjrV17RulkzrZs00fq55/Ss//1PBwYGasDzFRgYqGfNmmV3tY6VlpZmdwmO57QMw8LC9Guvveb5eeHChbpz585aa61vu+02/be//e2Mxw8aNKhau7crkpGRoadOnap37NhR5r7SDD/55BPdpk0bXb9+fV2/fn09Y8YMz2Nuv/12DeiEhASdn5+v33//fR0SEuI5HFCquLhYr1+/Xv/zn//UWVlZZZ5r8ODBtu2uX7dunY6MjNQFBQX6559/rnR3vdZaJyQk6KCgIF2/fn3dpEkTnZyc7Lnvp59+0rt27dJFRUX6+++/1xdeeKF+6qmnPPd/++23OicnR2dnZ+unnnpKt23bVh89erTMc7z11lu6ffv2ZXKsDVSwu94nLqHzJu0Dg/vUtaHVGHEuLDubh3bu5JKsLNa2aMHz3bvz66efsnr1avLy8s54bE5ODpMmTeLNN980rrHM4Dp+wNaTbVzCiRl2PG1yprCwMM/Y5U2bNi1zYlxWVhbNmjU75zKTkpK44447ALjiiiv4/PPPz7i/ffv2jBw5kgkTJrB+/foz7issLGT79u1MmDCBjz76iGuuuYYff/yR0aNH065dO0aNGkWTJk0IDw9n0qRJAEyYMIH4+HhWrlzJmDFjPMtSStG7d28WLVrE1KlTef7556uRjGXatGlMmzbNU9upU6do3ry55/5jx46xZ88eevTo4bnt5MmTlS6zuLiYu+66i+nTp9Ogwbnb2ldffcWUKVNYtmwZffr0ITU1ld/+9rd8/vnnXHrppXTu3Nnz2EsuuYS77rqLpKQkHnnkEQAGDx7suf+RRx5h5syZfPPNN/zmN7/x3D5//nweeeQRvvrqqzMOq9Q11+2uP3HihN0l+KSGxcXE7t5NQmoqnXJyiL/4YqZccgm/NmkCUKbBl6rodnFuZ+/uFNXnxAx/+eUXz/d79uyhXbt2gDW5yum7dbOzs/npp5+qtCs3JiaGkydPcvLkyTINvlRhYWG5x5SXL1/O5s2b6dq1KyNGjKBevXp069aNUaNGeZbVq1cv1FlX0Jz9c1Weqyoefvhhjh07xrFjx0hOTubyyy/3/Hzs2DHAmqyl9P97rgYP1oellJQUxo8fz3nnneeZNa5Dhw588803ZR6/ceNGhgwZQlRUFPXq1aNfv34MGDCAr776qtzlb9u2rdINSKXUGfd/8cUX3H777Xz66adccskl56y/VpW3ee/kL5mgphzffKP1xRdrDVpPnKj1wYNlHhIWFnbGrvrSr7CwsLqv1yWcMKmFr3NahmFhYbpnz576l19+0YcPH9aDBw/WjzzyiNZa64MHD+rg4GA9d+5cnZubq6dMmWJ0dv2sWbN0enq61lrr3bt36yFDhujf/e53ZR73xhtv6LS0NB0UFKSXLFmii4uLdVpamr7wwgs9+R4+fFg3b95cJyYm6sLCQv3hhx/qFi1a6MzMTF1UVKRff/11feTIEV1cXKzXrFmjzzvvPD19+nTPc+Tl5enc3Fw9aNAgPWPGDJ2bm6uLiorO+X/w1u764uJivX//fs/X2rVrNaAzMjJ0Xl5emccvW7ZMt2zZUm/YsEFrrfX69et1aGioXrRokdZa688++0z/+uuvWmutt23bptu1a6cff/xxrbXW6enp+ttvv/X8n5955hndqlUrfejQIa211kuWLNGhoaF6+fLlxv+v6qCC3fW2N2Vvf3Xr1s1roTnesWNa33mn9TKHh2v9xRcVPnTWrFlyTN7Lvv32W7tLcDynZRgWFqafeuop3b17dx0SEqJvueUWnZ2d7bn/yy+/1N26ddONGzfW0dHR+ueff/bcFx8fr0eOHFnl53r00Ud1+/btdWBgoG7fvr2+/fbbPY1Ga61Hjhyp4+PjPRnOmTNHR0RE6KZNm+r27dvrKVOmnNGIV6xYoXv27KmDgoJ037599YoVK7TWWhcVFekRI0boFi1a6KCgIH3RRRfp+Ph4XVxc7Pnd6OjoMhsIVbn8sTYuodNal3tMvkePHmf8PXv55Zf1hRdeqJs2baovuOAC/dxzz3nue/DBB3WbNm10YGCgvuCCC3RsbKzOz8/XWmu9efNmfckll+jAwEAdGhqqr7rqKr1u3TrP7w4dOlTXr19fBwUFeb6q87rWVEVN3icmqPGm3r176w0bNthdhv0+/hjuuQd+/RX+8hf4v/+D0y6BKU9SUhKTJk0iLy+PsLAw4uPjXTcYSV06fvw4ISEhdpfhaE7LMDw8nISEBK72oUtQnZahL3JChr4+QY3X+PoENUOHDq3WiXLVtm8f/P731lebNta48//5zzkbPFjH/QYOHEjXrl3ZvXu3NHhDc+bMsbsEx5MMzUmG5pycoeuavN8qLobXX7fmev/8c5g2DdauhagyH+yEEEL4CdddQleVyydcZ/t2a7a4b7+Fq66CN96ALl1qvDi/zLAWhIaG2l2C4zktw927d9tdQhlOy9AXOTlD123Jt2jRwu4S6k5+vnWsPTIStmyBd96Br74yavAArVu39lKB/u3666+3uwTHkwzNSYbmnJyh65q8rx+T95rvvrOmgp061Tr+vm0bxMZ6Zba4AwcOmNcnmDVrlt0lOJ5kaE4yNOfkDF3X5IuKiuwuoUJJSUmsXr2a5cuX13xs+Kws66z5wYOt7xcuhPffh7ZtvVLjsmXLuP/++72yLH+Xk5NjdwmO5/QM8/Pzuf766wkPD0cpZTzy4yuvvEJUVBQBAQHExsaWub+8qWwryvDKK6+kdevWBAcHExkZyYIFCzz37d+/n9/+9re0a9cOpVSZwxARERGeaVabNm1KgwYNzhjtrbbs3r2bK6+8ksDAQC6++OIKB68B2Lt3L2PGjCE0NJQOHTrw+uuvl/u4//3vfyilSEhIOOP29evXM2TIEJo2beoZTa+Uk6aadV2T91VJSUnExcV5RpBLT08nLi6ueo3+k0+s2eJefRXuvdfaRX/ddbVUsRDCGy6//HJmzZp1zklTqqJdu3b8/e9/57bbbitz36FDh/j973/Pk08+yZEjR4iKimL8+PEVLmv69Ons37+frKwsZsyYwcSJE9m/fz9gzTI3cuRI5s2bV+7vbtmyxTMa3YkTJ+jYsSM33HBDlf4P4eHhNT534cYbb6R3794cPnyY+Ph4rr/+ejIzM8t97MSJE7ngggs4cOAACxcu5NFHH2Xp0qVnPObo0aM89dRTZUYdPHToECNHjuSOO+7g8OHD/Otf/zqjkQN8+umnngwWL17suT0nJ4cXX3yRQ4cOsWbNGpYsWcJzzz1Xo/+vN7juOvm+ffvq1NRUry3PW5e7lTc2PEBAQAADBw6s9HdD8/O5Ny2NKzMz+SkoiOe6dmVbcLDnfm+PC5+fn0+jRo28ukx/JDmac0KG4eHh3HHHHbz77rvs37+fsWPH8tprr5WZdrVDhw7MmjXLK39T/v73v5ORkUFiYqLnthkzZpCYmMiqVasAa9jcVq1asWbNmjNmVyvP2rVrGTJkCCtWrKB///6e2wsLC2nYsCE///wz4eHh5f7u8uXLGT16NL/++usZ09FWJDw8nGXLllW4vIrs3LmTSy65hEOHDnnG+7/iiiuIiYnhzjvvPOOxJ0+epFmzZhw8eNBzjlFcXBy5ubm8++67nsfdeeed9OrViw8++ICJEycyefJkwJrN7pdffvE89uz1sDrjITz//PMsXbqUTz/9tFr/3+rym+vkfXXs+hqNDa81o/bvZ+a6dQw+dIg3w8OJ69PnjAZfG1asWFGry/cXkqM5p2RY0fSy5/Lee++dswFXVUVT2Va0NQ4wevRoGjduzIABAxg6dChRNbjkdubMmYwbN65KDd7Eli1b6Ny58xkT+kRGRnqm6j1d6cbr6RuxWms2b97s+Xnt2rWkpKSU+YAA1kZZaGgogwYNok2bNkRHR7Nnz54zHuOUqWZdd62UtydU8dZWcnh4OOnp6WVuDwsLK/85du6EuDjr3+homDGD27t25XavVFO5Xbt21cGzuJ/kaM4pGd5zzz2e2ecee+wx7r333io1+ptuuslr87KfPHmyzJUxISEh5f7dKZWcnExBQQFfffUV27Zto1696m335eTkMHfuXD755JMa1VwdJ0+eLDPqXEhICHv37i3z2GbNmjF48GCefPJJnn32WbZu3cq8efM8+RQVFXHXXXfxyiuvlPt/zsjIYP369Xz55ZdccsklXHfdddx4442sXLkSsD7U9enTB60106dPZ8SIEWzfvv2M2fQA3n77bVJSUsoc769LrtuS91Xx8fEEBgaecVtgYCDx8fFnPrCgAJ56Cnr1gk2b4M034euvoWvXOqxWCFEdFU0vW5cqmsr27MMGZ2vYsCHXXnstixcvrnaz/uijjwgNDSU6OrrCx+zZs4fmzZt7vvbs2UOvXr08P7/33nsAXHvttZ4T2co7V6m6U/UmJSXx888/07FjR/70pz8xceJEOnToAMCrr75Kr169KjxU2qRJE373u9/Rr18/GjduzOjRo1m1ahXHjx8HrKlmmzRpQmBgII888gjNmzcvM9td6VSzn3/+ua1TzbpuS95XxxcuHSK20rHh16yxBrX54Qe44QZ46SXwwsk61TVixIg6f043khzNOSXDiqaXrUsRERHMnDnT83PpVLYvvvhilX6/JtPHzpw5k1tuuaXSaWk7derkmUIWKj4mX9EUuqUiIiLYtWsXJ06c8DT2TZs2VbgnJCwsjOTkZM/PN910k+d8gyVLlrB8+XI+++wzwLr0esOGDWzcuJFXXnmlzNS7w4YNq7S2iqaaXbhwoUw16+2vHj16mE7mU6uio6PLzrp04oTW992ntVJat2+v9YIFdpTm8cMPP9j6/G4hOZpzQoaVTS+rtdanTp3Subm5un379nrRokU6Nzf3jBncqqOgoEDn5ubqhx9+WE+cOFHn5uZ6ZlqraCrb8jLctm2b/uyzz3ROTo7Oz8/X7777rm7YsKFOTU31PCY3N1efPHlSA3r79u06Nzf3jGX88ssvun79+jotLa1a/4ewsLAzZt+rjgEDBugHH3xQ5+bm6o8++kiHhITog+VMna211lu3btVZWVk6Ly9Pv/vuu7ply5aexx49evSMqWkvu+wy/Z///EcfO3ZMa21NF9u8eXO9YcMGnZ+frydOnKgvv/xyrbVMNWv7l6/PJ1+myS9cqHWnTlaDv/turY8ft622Uk6bw9tXSY7mnJDhuaaXDQsLKzMNa2mTmzVrlq7OhsnUqVPLLGvq1Kme+8ubyrY0wzvuuEPfcccdWmurAfbv3183bdpUh4SE6KioKP3RRx+d8VxnP4+1Tfj/PfXUU57GVx0mTf7nn3/W0dHRunHjxrpr1676yy+/9Nx3dpYvvPCCbtWqlQ4MDNSDBw8+YzrYs0VHR+s333zzjNteffVV3a5dO928eXN9ySWX6D179mitnTfVrO1N2dtfjmnyBw5oPWGC9RL06KH1ypV2l+bhhD+sTiA5mnNChmFhYWc0G1/jhAx9nRMyrKjJu+6Y/Nknt/kcrRl54ABcfDFkZ8MTT8Df/gYBAXZX5nH6ZTii5iRHc5KhOcnQnJMzdF2TD/ChZlnG0aMsa9gQduyAyy+HGTOsqWF9TBfDCW6ERXI0JxmakwzNOTlD111Cd/ToUbtLqFhICDRqBK+9BsuX+2SDByodPENUneRozgkZ7t69u0ojn9nFCRn6Oidn6LoteZ9Wr541oYwXZooTQgghzsV1W/INGza0u4TKOaDBt2nTxu4SXEFyNCcZmpMMzTk5Q9dNUBMVFaVTUlLsLkMIIYSoM34zQc3hw4ftLsHxTp/ZStSc5GhOMjQnGZpzcoaua/LFxcV2l+B4+fn5dpfgCpKjOcnQnGRozskZuq7JCyGEEMIix+RFGcXFxdWeclKUJTmakwzNSYbmnJCh3xyTP3sqQlF9ixcvtrsEV5AczUmG5iRDc07O0HVNPi8vz+4SHG/Pnj12l+AKkqM5ydCcZGjOyRm6rskLIYQQwuK6Jt+8eXO7S3C8UaNG2V2CK0iO5iRDc5KhOSdn6LomX1BQYHcJjpeZmWl3Ca4gOZqTDM1JhuacnKHrmnx2drbdJTje2rVr7S7BFSRHc5KhOcnQnJMzdF2TF0IIIYTFddfJK6VOADvsrsPhWgGH7C7CBSRHc5KhOcnQnBMyDNNatz77RjdONbujvAEBRNUppVIkQ3OSoznJ0JxkaM7JGcrueiGEEMKlpMkLIYQQLuXGJj/D7gJcQDL0DsnRnGRoTjI059gMXXfinRBCCCEsbtySF0IIIQQua/JKqZFKqR1KqTSl1MN21+M0SqmOSqmlSqmtSqktSqn77K7JqZRS9ZVSG5RSyXbX4kRKqeZKqblKqe1KqW1KqcvsrslplFL3l7yPNyul3ldKNba7JidQSr2tlDqolNp82m2hSqkvlVI/lvzbws4aq8M1TV4pVR/4L3At0AO4USnVw96qHKcQeFBr3QMYCNwtGdbYfcA2u4twsOnAF1rri4FIJMtqUUq1B/4MRGmtewL1gQn2VuUYicDIs257GFiitb4IWFLysyO4pskD/YE0rfUurXU+MBsYY3NNjqK13q+1Xl/y/QmsP6zt7a3KeZRSHYBRQILdtTiRUioEGAK8BaC1ztdaH7O1KGdqADRRSjUAAoF9NtfjCFrrFcCRs24eA8ws+X4mMLYuazLhpibfHvjltJ8zkAZVY0qpcKA3sMbmUpzoRWAKUGxzHU51AZAJvFNyyCNBKRVkd1FOorXeCzwH7AH2A8e11ovtrcrR2mqt95d8/yvQ1s5iqsNNTV54iVKqKTAP+IvWOsvuepxEKTUaOKi1TrW7FgdrAPQBXtNa9waycdDuUV9Qcsx4DNYHpnZAkFJqor1VuYO2LklzzGVpbmrye4GOp/3coeQ2UQ1KqYZYDT5Ja/2R3fU40GDgt0qp3ViHjK5SSs2ytyTHyQAytNale5HmYjV9UXVXAz9rrTO11gXAR8Agm2tysgNKqfMBSv49aHM9VeamJr8OuEgpdYFSqhHWSSaf2FyToyilFNZx0G1a6+ftrseJtNaPaK07aK3DsdbBr7XWsgVVDVrrX4FflFLdSm4aBmy1sSQn2gMMVEoFlryvhyEnL5r4BLi15PtbgQU21lItrpmgRmtdqJS6B1iEdSbp21rrLTaX5TSDgZuBH5RSG0tue1Rr/Zl9JQk/dS+QVPKBfRfwR5vrcRSt9Rql1FxgPdZVMxtw8KhtdUkp9T4wFGillMoApgLTgA+UUpOAdOAP9lVYPTLinRBCCOFSbtpdL4QQQojTSJMXQgghXEqavBBCCOFS0uSFEEIIl5ImL4QQQriUNHkhhBDCpaTJCyGEEC4lTV4IIYRwKWnyQogaUUoFKaW2K6XWlsx5UHr7cKVUsVLqbjvrE0LIiHdCCANKqd7AauAFrfXDSqm2wCZgjdZ6jL3VCSGkyQshjCil7seau3wE8BBwCRCptT5ka2FCCGnyQggzJbOcLQSuAhoB12itl9hblRAC5Ji8EMKQtrYU3gUCgE3S4IXwHdLkhRBGlFLnAdOxpjWNVErdZ3NJQogS0uSFEDVWsqt+JpAHXA28CPxbKdXLzrqEEBY5Ji+EqDGl1IPAM8BVWuvlSqlGWGfbBwBRWutcWwsUws/JlrwQokaUUn2Ap4CntdbLAbTW+cCNQDjwvH3VCSFAtuSFEEII15IteSGEEMKlpMkLIYQQLiVNXgghhHApafJCCCGES0mTF0IIIVxKmrwQQgjhUtLkhRBCCJeSJi+EEEK4lDR5IYQQwqX+H1+PVQiPdwGxAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 576x432 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "fig, params, param_errors, reduced_chi2 = lab.fit(x, y, ex, ey, init_params)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Get help"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "\u001b[0;31mType:\u001b[0m        module\n",
       "\u001b[0;31mString form:\u001b[0m <module 'LabModule' from '/Users/gwilliam/cernbox/Teaching/PHYS105/Materials/LabModule.py'>\n",
       "\u001b[0;31mFile:\u001b[0m        ~/cernbox/Teaching/PHYS105/Materials/LabModule.py\n",
       "\u001b[0;31mDocstring:\u001b[0m   Fitting code for laboratory\n"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "lab?"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Help on module LabModule:\n",
      "\n",
      "NAME\n",
      "    LabModule - Fitting code for laboratory\n",
      "\n",
      "FUNCTIONS\n",
      "    example()\n",
      "        Example function to illustrate usage\n",
      "    \n",
      "    example_poly()\n",
      "        Example function to illustrate usage of general polynomial fitting.\n",
      "    \n",
      "    exponential(params, xdata)\n",
      "        Function for an exponential \n",
      "        - params   array of function parameters\n",
      "        - xdata    array of x data points\n",
      "    \n",
      "    exponential_diff(params, xdata)\n",
      "        Diffential of function for an exponential\n",
      "        - params   array of function parameters\n",
      "        - xdata    array of x data points\n",
      "    \n",
      "    fit(xdata, ydata, xerror, yerror, init_params, func=<function polynomial at 0x7ff3ff8de790>, diff=<function polynomial_diff at 0x7ff3ff8de820>, fmin=<function minimise at 0x7ff3ff8de9d0>, method='lsq', loss='linear', fscale=1, plot=True, plot_name='plot.png', text_loc='lower right')\n",
      "        Function to perform least-squares fit\n",
      "        - xdata         array of x data points\n",
      "        - ydata         array of y data points\n",
      "        - xerr          array of x data errors\n",
      "        - yerr          array of y data errors\n",
      "        - init_params   array of function parameters\n",
      "        - func          function describing data [default: straight_line]\n",
      "        - diff          differential of function [default: straight_line_diff]   \n",
      "        - fmin          function to minimise [default: minimise]\n",
      "        - method        Method to use: lsq or odr [default: lsq]\n",
      "        - loss          loss function used for lsq [default: linear]\n",
      "        - fscale        loss scale for lsq loss functions other than linear [default: 1]\n",
      "        - plot          Plot results or not?\n",
      "        - plot_name     Name to save plot to [default: plot.png]\n",
      "        - text_loc      Where to write the fit params on the plot [default: lower right]\n",
      "        \n",
      "        Return: (fitted parameters, parameter errors, reduced chi2)\n",
      "    \n",
      "    generate_toy_data(xarray, func)\n",
      "        Generate toy data, randomly poisson fluctuated around f(x) = func(xarray)\n",
      "    \n",
      "    guess_poly_init_params(deg, xdata, ydata, yerror=None)\n",
      "        Do a simple polynomial fit in numpy w/o taking into account the x errors, to determine starting params\n",
      "        - xdata    array of x data points\n",
      "        - ydata    array of y data points\n",
      "        - yerr     array of y data errors [optional]\n",
      "        - deg    degree of polynomial\n",
      "    \n",
      "    gui()\n",
      "        Simple TKinter gui, called if run with option --gui\n",
      "    \n",
      "    load_data_from_csv(fname, delim=',', skiprows=0, usecols=None)\n",
      "        Load data from a CSV file in form x, xerror, y, yerror\n",
      "        - fname     csv file name\n",
      "        - delim     column separator [default: ,]\n",
      "        - skiprows  skip N rows [default: 0]\n",
      "        - usecols   int or list of columns to use [default: all]\n",
      "        \n",
      "        Return: (xdata, xerror, ydata, yerror)\n",
      "    \n",
      "    load_data_from_excel(fname, sheet_name=0, skiprows=0, usecols=None)\n",
      "        Load data from a Excel file in form x, xerror, y, yerror\n",
      "        - fname       Excel file name\n",
      "        - sheet_name  Name or index of excel sheet to take [defaut: first]\n",
      "        - skiprows    skip N rows [default: 0]\n",
      "        - usecols     int or list of columns to use [default: all]\n",
      "        \n",
      "        Return: (xdata, xerror, ydata, yerror)\n",
      "    \n",
      "    main()\n",
      "        Driver function when run from command line\n",
      "    \n",
      "    minimise(params, xdata, ydata, xerr, yerr, func, diff)\n",
      "        Calculation to minimise to find the best fit using chi2: difference beteen each y point and its \n",
      "        prediction  by the function, divided by the sum in quadrature of the error on y, both from the \n",
      "        y error and from the related error in x.\n",
      "        - params   array of function parameters\n",
      "        - xdata    array of x data points\n",
      "        - ydata    array of y data points\n",
      "        - xerr     array of x data errors\n",
      "        - yerr     array of y data errors\n",
      "        - func     function describing data\n",
      "        - diff     differential of function\n",
      "    \n",
      "    polynomial(params, xdata)\n",
      "        Function for a ploynomial of order `len(nparams)`\n",
      "        - params   array of function parameters\n",
      "        - xdata    array of x data points\n",
      "    \n",
      "    polynomial_diff(params, xdata)\n",
      "        Diffential of function for a ploynomial of order `len(nparams)`\n",
      "        - params   array of function parameters\n",
      "        - xdata    array of x data points\n",
      "    \n",
      "    quadratic(params, xdata)\n",
      "        Function for a quadratic \n",
      "        - params   array of function parameters\n",
      "        - xdata    array of x data points\n",
      "    \n",
      "    quadratic_diff(params, xdata)\n",
      "        Diffential of function for a quadratic \n",
      "        - params   array of function parameters\n",
      "        - xdata    array of x data points\n",
      "    \n",
      "    straight_line(params, xdata)\n",
      "        Function for a straight line\n",
      "        - params   array of function parameters\n",
      "        - xdata    array of x data points\n",
      "    \n",
      "    straight_line_diff(params, xdata)\n",
      "        Diffential of function for a straight line\n",
      "        - params   array of function parameters\n",
      "        - xdata    array of x data points\n",
      "\n",
      "DATA\n",
      "    exp = <ufunc 'exp'>\n",
      "\n",
      "FILE\n",
      "    /Users/gwilliam/cernbox/Teaching/PHYS105/Materials/LabModule.py\n",
      "\n",
      "\n"
     ]
    }
   ],
   "source": [
    "help(lab)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "ename": "SyntaxError",
     "evalue": "unmatched '}' (<ipython-input-7-b3127725f678>, line 1)",
     "output_type": "error",
     "traceback": [
      "\u001b[0;36m  File \u001b[0;32m\"<ipython-input-7-b3127725f678>\"\u001b[0;36m, line \u001b[0;32m1\u001b[0m\n\u001b[0;31m    }\u001b[0m\n\u001b[0m    ^\u001b[0m\n\u001b[0;31mSyntaxError\u001b[0m\u001b[0;31m:\u001b[0m unmatched '}'\n"
     ]
    }
   ],
   "source": [
    "}"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.8.5"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}