Blackbody_PowerVsTemperature.ipynb

{
  "nbformat": 4,
  "nbformat_minor": 0,
  "metadata": {
    "colab": {
      "name": "Blackbody_PowerVsTemperature",
      "version": "0.3.2",
      "provenance": [],
      "collapsed_sections": [],
      "include_colab_link": true
    },
    "kernelspec": {
      "name": "python3",
      "display_name": "Python 3"
    }
  },
  "cells": [
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "view-in-github",
        "colab_type": "text"
      },
      "source": [
        "[View in Colaboratory](https://colab.research.google.com/github/jrstevenjlab/ColaboratoryAtomicLab/blob/master/Blackbody_PowerVsTemperature.ipynb)"
      ]
    },
    {
      "metadata": {
        "id": "UwPPL5XqgtbL",
        "colab_type": "text"
      },
      "cell_type": "markdown",
      "source": [
        "The goal of the excercise is to fit the function \n",
        "\n",
        "$$S_{model}(T) = B + C \\cdot T^{n}$$\n",
        "\n",
        "to data for power as a function of temperature that you recorded in lab.  The parameters of the function are:\n",
        "\n",
        "*   B: a background term, the value the data would have at T=0,\n",
        "*   C: a scale factor which multiplies $T^n$, which translates the units\n",
        " from $temperature^n$ to radiated power.\n",
        "*   n: the power in the power law, which we expect to be 4\n",
        "\n",
        "By \"fitting\" the data we mean that the parameters of the model are allowed to vary until the curve defined by this function best match the data you've measured.  This is done by minimizing the $\\chi^2$ which is a metric for comparing the model function to the observed data.  The $\\chi^2$ is defined as\n",
        "\n",
        "$$\\chi^{2} = \\sum_{i=0} \\frac{(S_i - S_{model}(T_i))^2}{\\sigma_{S_{i}}^2}$$\n",
        "\n",
        "where the sum of $i$ is over your observed data points\n",
        "*  $S_i$ are the individual data points for the power (y-axis value), \n",
        "*  $\\sigma_{S_{i}}$ are the errors on the individual data points,\n",
        "*  and $S_{model}(T_i)$ is the value of your model function at the temperature $T_i$\n",
        "\n",
        "The steps below describe the process that you should follow:\n",
        "\n",
        "1.   Enter your data for power and position in the ordered lists below, replacing the example data already there.\n",
        "2.   Run the fit code by pressing the \"Play\" button, which will cause a few things to happen:\n",
        "> * Your data will be plotted as blue points in the figure below the code \n",
        "> * The `curve_fit` command in the code will find the \"best fit\" values for the paremeters of the model which gives the smallest value for the $\\chi^2$ possible (ie. it \"minimizes\" the $\\chi^2$).\n",
        "> * The function S(T) with these \"best fit\" parameters will be plotted in red overlayed with your data points\n",
        "> * The values for the parameters are output as text below the figure, along with the uncertainties on those parameters.\n",
        "\n",
        "Finally, we need evaluate the \"quality\" of our fit by looking at the $\\chi^{2}/DOF$, where DOF is the **Degrees Of Freedom** of the fit, defined as:\n",
        "\n",
        "$$DOF = # data points - # parameters$$\n",
        "\n",
        "In our case the # parameters = 3 (B, C, and n). and the # data points is given by the number of points you observed in lab.   In general, a \"good\" fit will have a $\\chi^{2}/DOF$ close to 1.  \n",
        "*  If $\\chi^{2}/DOF$ >> 1, then the errors were likely **underestimated** meaning the model function does not describe the data.\n",
        "*  If $\\chi^{2}/DOF$ << 1, then the errors were likely **overestimated** meaning the model describes the data better than it should (statistically)."
      ]
    },
    {
      "metadata": {
        "id": "DniaGSKaimDV",
        "colab_type": "code",
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 488
        },
        "cellView": "code",
        "outputId": "fcdef694-2bce-40a4-86ee-688f5216da76"
      },
      "cell_type": "code",
      "source": [
        "# import modules needed for data analysis\n",
        "import numpy as np\n",
        "import matplotlib.pyplot as plt\n",
        "from scipy.optimize import curve_fit\n",
        "\n",
        "# function for power vs distance with paramters B, C, and n\n",
        "def Smodel(x, B, C, n):\n",
        "  return B + C * np.power(T, n)\n",
        "\n",
        "##########################\n",
        "# INSERT YOUR DATA HERE! #\n",
        "##########################\n",
        "\n",
        "# enter as a list here for radiation (S) and temperature (T)\n",
        "# T = temperature (Kelvin)\n",
        "T = [800,1000,1200,1300,1500,1600,1700,1800,1800,1900,2000]\n",
        "\n",
        "# S = converted sensor reading (W/m^2)\n",
        "S = [7.95,19.32,43.18,62.50,92.05,125.0,161.36,200.0,240.91,284.09,329.55]\n",
        "\n",
        "# uncertainy (error) on converted sensor reading: YOU NEED TO DETERMINE THIS!!! \n",
        "# -choose between one of the options for uncertainties below\n",
        "\n",
        "# fixed uncertainty on power for all datapoints in W/m^2\n",
        "#Serr_value = 2.0 \n",
        "#Serr = np.repeat(Serr_value, len(S)) \n",
        "\n",
        "# scale uncertainty which is a fraction of the power for each datapoint in W/m^2\n",
        "Serr = 0.1*np.array(S) # eg. 10% error\n",
        "\n",
        "# plot the data on a new figure\n",
        "plt.figure() # create new figure for plotting\n",
        "plt.errorbar(T, S, Serr, fmt='o', label=\"data\") # plot S vs T, where Serr is the uncertainty on your power values\n",
        "plt.axis([0, 2500, 0, 500])  # set axes [Tmin, Tmax, Smin, Smax]\n",
        "plt.xlabel('Temperature (Kelvin)') # set x-axis label\n",
        "plt.ylabel('Power/Area (W/$m^2$)') # set y-axis label \n",
        "\n",
        "# fit the data\n",
        "initialParameters = [0,1e-12,4] # start with initial guess for fit parameters [B, C, n]\n",
        "finalParameters, finalParameterErrors = curve_fit(Smodel, T, S, initialParameters, Serr, True) #\n",
        "\n",
        "# show red curve for fit results and plot on top of the data\n",
        "plt.plot(T, Smodel(T, *finalParameters), 'r--', label='fit: B=%5.3f, C=%5.3f, n=%5.3f' % tuple(finalParameters))\n",
        "plt.legend()\n",
        "plt.show()\n",
        "\n",
        "# output fit parameter values and errors\n",
        "print(\"Model parameters:\")\n",
        "print(\"B   =\", finalParameters[0], \" +/- \",np.sqrt(finalParameterErrors[0,0]))\n",
        "print(\"C   =\", finalParameters[1], \" +/- \",np.sqrt(finalParameterErrors[1,1]))\n",
        "print(\"n   =\", finalParameters[2], \" +/- \",np.sqrt(finalParameterErrors[2,2]))\n",
        "\n",
        "# compare model function with obtained parameters to data to determine fit quality\n",
        "Smodel_i = Smodel(T, *finalParameters)\n",
        "deltaS = S - Smodel_i # calculate difference between data and model\n",
        "chisq = np.sum((deltaS/Serr)**2) # calculate chi^2 by summing over difference for all data points\n",
        "dof = len(S) - len(finalParameters) # Degrees of Freedom = # data points - # parameters\n",
        "print(\"\\nFit quality:\")\n",
        "print(\"chisq =\",chisq,\"dof =\",dof)\n"
      ],
      "execution_count": 1,
      "outputs": [
        {
          "output_type": "display_data",
          "data": {
            "image/png": "iVBORw0KGgoAAAANSUhEUgAAAgYAAAFcCAYAAACzwnasAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMS4yLCBo\ndHRwOi8vbWF0cGxvdGxpYi5vcmcvNQv5yAAAIABJREFUeJzs3XtgzuX/x/Hnfe9gZmN2wsw5lbOW\nrDmmVFi05JQoIkrKt4gpSVE5RoiWRiLlax0oRCclZr5ZNKYDfo6bnbDZ+XT//ljukg273ffubV6P\nf77uz/35fO73ru6vvVzX9bkug8lkMiEiIiICGO1dgIiIiJQfCgYiIiJipmAgIiIiZgoGIiIiYqZg\nICIiImYKBiIiImLmaK8PjoqKYty4cTRt2hSAG2+8kZEjRzJx4kQKCgrw8fFhzpw5ODs7s2HDBlau\nXInRaGTAgAH079/fXmWLiIhUanYLBgDt27dn4cKF5teTJ09m8ODB9OzZkzfffJOIiAhCQkJ4++23\niYiIwMnJiX79+nH33Xfj4eFhx8pFREQqp3I1lBAVFcVdd90FQLdu3YiMjGTfvn20atUKd3d3XFxc\nCAgIIDo62s6VioiIVE527TE4dOgQTzzxBKmpqYwdO5asrCycnZ0B8PLyIikpieTkZDw9Pc3XeHp6\nkpSUZK+SRUREKjW7BYOGDRsyduxYevbsyYkTJ3jkkUcoKCgwv1/SSs1Xu4KzyWTCYDBYpVYREZHr\nhd2CQa1atejVqxcA9evXx9vbm5iYGLKzs3FxcSEhIQFfX198fX1JTk42X5eYmEjbtm2veH+DwUBS\n0nmb1S/g4+OuNrYxtXHZUDvbntrY9nx83K1yH7vNMdiwYQPh4eEAJCUlkZKSQt++fdmyZQsAW7du\npXPnzrRp04aYmBjS0tLIyMggOjqadu3a2atsERGRSs1uPQZ33nknEyZM4NtvvyUvL49p06bRrFkz\nJk2axNq1a/Hz8yMkJAQnJyfGjx/PiBEjMBgMPPXUU7i7WycViYiIyMUMlXnbZXVb2Za6Bm1PbVw2\n1M62pza2vQo/lCAiIiLlj4KBiIiImCkYiIiIiJmCgYiIiJgpGIiIiIiZXZdEFhERyM/P58knR9Cg\nQUMCAtpRrZobXbt24/vvv6Fbt+4lXvfaa9P4/feDVK9eg9zcXJo2vZHx40MxGkv+N9+aNR/w/fff\nAAYee+xxgoI6XfR+YWEhYWFvs3Hjer788pvLXpeens4rr7xIeno6Vau6Mm3aDKpXr1HiZx88eIAl\nSxaSm5tLXl4enTp1Yfjwxy+7Su3ChfM4cGA/BoOBcePG06xZi4ve/9//onj33bcxGh0ICurIsGEj\nS7wuIeE006dPpbCwEC8vb1566VXzMvzWYDKZGDNmBO3aBTJixOiL3vvzzz+YN28mBgM0adKUCRMm\nA/Dzz7tZvHgBRqORvn37cd99IeZrjhw5xGOPDeGjjz6lTh0/q9V5JeoxEBGxs+TkZPLy8pgy5RV6\n9epN167dyMvLY+3aNVe8dvTosSxe/C7vvvs+J0+eIDZ2f4nnxsWd4ptvtrJkSTizZy9g0aL5Fy1F\nD7B69fvUqlX7ouXnS7ruv/9dwy233MrSpeF07dqN1atXlvjZRSHiJZ599nnCwlbw7rvv8+eff/DF\nF5+XeM0vv+zh5MkThIWtIDT0JRYsmHvJOW+9NZcZM2azdGk4u3fv4v/+70iJ14WHh9G37wCWLHkP\nf/96bNy4ocTPtsQXX3xOXl5+se8tXDiPcePGs3TpctLT04mM3EF+fj5z577B7NnzWbLkPXbvjjKf\nbzKZWLz4Lfz961m1xquhHgMRkX/xvLVlscczxzxD9ohRALiPeRynqMhLzsm7tR3n330fAJdV7+O6\nYC5n9pT8yxpg0aJ5nDp1ktdff4VatWrj4eHB0aNHOXz4EHPnzmT48JGEh4cxceKLJd4jNzeXzMxM\natb05LffDrJ48fyL3u/UqQtubu7cfnsHnJycqFmzJrVr1+Ho0f+jSZMbzOf16zcQV9dqhIe/Yz4W\nHf1zsdft2fM/Jk+eCkDHjl2YOPE/Jdb3xRdf0KVLVxo3LvosR0dHXnrpFapUcWHnzp9Ys+aDi87v\n06cvx48fpXPnOwBo2LAR58+nkZGRTrVqbgCcOnUSd/fq1KpVG4CgoI7s2bObc+fOFXvdL7/sMf9L\nvWPHznz00SoeeKBfsfW+9to0vL19+P33gyQknGbq1Bl//bK+tF0HDRrCuXPn+Prrr7j//r4kJiZc\ndE5eXh7x8XHm3o6OHTvz88+7qV69Bv7+9fD1rQXAq6++Yb5m48YNtGt3Gzt3/lRim9qKgoGIiJ2N\nHfssU6ZM4oUXXiY8PAyAwYOHEhu7nwkTQgFKDAVhYYv56KNVxMWdokuXO6hb1x+AxYvfveTcDz5Y\njodHTfPrmjVrkpKSfFEwcHWtdsl1Z86kFHtdSsrfxy8cK8mRI0e44YabLjp24bM6dOhEhw6dLrlm\n1qzXuOmmm82vPTxqkpKSYg4GxdV16tQpzp07V+x1/9zBt2ZNT1JSUkqsF4rC1ptvLubzzyP46quN\njBs3vth2BViy5C1GjRrDiRPHL3kvNfXcRSv2Fn12MqdPx+Hk5MRLL4WSnJxI374DuPvuHqSmnuOr\nrzayYMESBQMRkfLgSv/CBzi/ZNkVz8keOozsocOsUFHJRo8eS8eOnSksLGTOnDf48svPLxqnvhxL\n170t7rorLaJrMBgoLCy47DlX/tzLf0ZJbxd33dUs+tumzS0A+PjUIjb2QInn7d0bjYODA61atSk2\nGJT02SaTiYSE0yxZEk5OTjaPPTaE9u1vZ+nSRTz++JM4OtrnV7SCgYhIJWA0GunSpSvffvs1N9xw\nU7Fd3t7ePhw/fsx8LCkpEW9v7yveu6TrvL29OXMmGTc3N5KTk/D29inxHo0bN+bXXw/Qo0ew+di5\nc+fIzs7iyJHDxQ4leHt7X/Sv+uTk5Ivq9fb24cyZv9+/UJejo2Ox11Wt6kpOTjZVqrhc1c/u4OBg\n/rPJZCpxiCYpKYnffotl1KhhnDt3lry8POrW9Tf/rB4eNUlNTf1HPUVt5enpxc03N8fFxQUXFxca\nN27CqVMn2bPnfxw5chiAo0f/jxdemMBbby297MROa1IwEBEphwwG4yUTA68kNvYA9es34OabmxXb\n5X369GnWrv2QESNGk5p6jqSkJBo2bHzF+wYE3Fbsde3b3853333DsGEj2bbtWwIDg0q8R+/evXnn\nnXe59979NG/ekry8PObOfZ3bb+/AffeFFDuUEBOzj/DwMEJCHuT333/D29v7oqGOOnX8yMjIID4+\nDh8fX3bu/ImpU6eTmnqu2OvatWvPtm3fce+9vfjhh+8IDOxwlS1bpKR2/adNm74gPj7uogDk6OhI\ngwYN2bdvL23atOWHH77jwQcH0qJFK955ZzE5OTkYDAZOnDhBnTp1Wbfu70mRY8eO4sUXp5VZKAAF\nAxGRcsnb25v8/DymTJnEs88+X+LkwwtzDC48gvfCCy+XeM/atWvTu3cITz1V9IjghAlFjzbu2rWT\n+Pg4HnigH/Pnz+bw4UOkp6czduwo8+S64q7r128Q06e/xJgxI3Fzc2fq1OkAvPXWPPr3H4SfX13z\nZ1erVo15895i9uzXycnJwcHBgbvv7nHZYY9Wrdpw003NeOKJxzAYDDz33CSg6JfvhUc6J0wIZdq0\nona58867qV+/AdCg2OtGjBjNjBlTWb/+U2rXrkPPnvcB8PLLk3nhhZepUsWldP+RruCf7frMM+OZ\nM+d1TKZCmjdvyW23BQLwyCPD/2pXeOihIdSsWfMKd7U97a4oFtNuabanNi4bamfr+vTTdQQFdbzo\n2fvy3MZhYW8zYsRou43pW4t2VxQRkXLJ17dWmS7Ic61atGhV4UOBNSkYiIiIVXXq1MXeJZRKRavX\n1hQMRERExEzBQERERMwUDERELPD8kp08v2SnvcsQsToFAxGRciQzM5N+/XqX+P62bd+WYTVyPVIw\nEBEppajYBM6l55CSls3U8CiiYhOufJEVxMfH8c03W8rks+T6peczRERKISo2gbANf6+bfzIpw/w6\nsHkti+6ZkZHOiy9OJDc3l9at2wKwdetmIiLW4uBgpGHDJkya9CJvvjmLgwcPsGLFMoKD+zB9etHO\nhvn5+UyZ8op5AyWRa6EeAxGRUtgYebSE48eKPX41tmzZTOPGTViy5D2aNr0RgKysLObNW8TSpcs5\nfrxoC+aHHhpK27YBDB/+OCkpyQwf/jiLFoURHNyHTz9dZ/Hni/yTegxEREohLjmz2OPxKRkW3/Po\n0SO0bXsrALfcUvS/1atXZ/Lk8QAcO/Z/pKaeu+gaT08vFiyYS3h4GOfPp3HTTc0s/nyRf1IwEBEp\nBT9vV04mXRoC6nhVK+bsq2MygdFoAKCw0EReXh5vvjmb999fg5eXNxMn/ueSa8LDwwgMvJ2QkH58\n//037Nz5k8WfL/JPGkoQESmF4KCGJRxvYPE969dvwG+/HQQgOvpnMjMzcXBwwMvLm4SE0/z220Hy\n8/MxGv/ecfHcuXPUreuPyWTip59+IC8vz+LPF/knBQMRkVIIbF6L0X1a4PDXv/D9fdwY3aeFxRMP\nAXr0CObAgRjGjXuSEyeO4eHhwW23BTJy5COsWLGMwYOHsnDhmzRo0Ijff/+NhQvncf/9fZk/fw7j\nxz/DXXfdy9690ezevctaP6Zcx7S7olisPO+WVlmojcuGJe18YXGjOWM62KKkSkffZduz1u6KmmMg\nImIBBQKprDSUICIiImYKBiIiImKmYCAiIiJmCgYiIiJipmAgIiIiZgoGIiIiYqZgICIiImYKBiIi\nImKmYCAiIiJmCgYiIiJipmAgIiIiZgoGIiIiYqZgICIiImYKBiIiImKmYCAiIiJmCgYiIiJipmAg\nIiIiZgoGIiIiYqZgICIiImYKBiIiImKmYCAiIiJmCgYiIiJiZvdgkJ2dTffu3fn000+Jj49n6NCh\nDB48mHHjxpGbmwvAhg0bePDBB+nfvz/r1q2zc8UiIiKVl92DwdKlS6lRowYACxcuZPDgwaxZs4YG\nDRoQERFBZmYmb7/9Nu+//z6rVq1i5cqVnDt3zs5Vi4iIVE52DQaHDx/m0KFD3HHHHQBERUVx1113\nAdCtWzciIyPZt28frVq1wt3dHRcXFwICAoiOjrZj1SIiIpWXXYPBrFmzCA0NNb/OysrC2dkZAC8v\nL5KSkkhOTsbT09N8jqenJ0lJSWVeq4iIyPXA0V4f/Pnnn9O2bVvq1atX7Psmk6lUx4vj4+NuUW1y\n9dTGtqc2LhtqZ9tTG1cMdgsG27Zt48SJE2zbto3Tp0/j7OyMq6sr2dnZuLi4kJCQgK+vL76+viQn\nJ5uvS0xMpG3btlf1GUlJ521VvlD0f3K1sW2pjcuG2tn21Ma2Z63gZbdgsGDBAvOfFy1aRN26dfnl\nl1/YsmUL999/P1u3bqVz5860adOGKVOmkJaWhoODA9HR0bzwwgv2KltERKRSs1swKM7TTz/NpEmT\nWLt2LX5+foSEhODk5MT48eMZMWIEBoOBp556Cnd3dUeJiIjYgsFUmkH7CkbdVralrkHbUxuXDbWz\n7amNbc9aQwl2X8dAREREyg8FAxERETFTMBAREREzBQMRERExUzAQERERMwUDERERMVMwEBERETMF\nAxERETFTMBAREREzBQMRERExUzAQERERMwUDERERMVMwEBERETMFAxERETFTMBAREREzBQMREREx\nUzAQERERM0d7FyAiIpXb80t24uBgYOboIHuXIldBPQYiIiJipmAgIiIiZgoGIiIiYqZgICIiImYK\nBiIiImKmYCAiIiJmCgYiIiJipmAgIiIiZgoGIiIiYqZgICIiImYKBiIiImKmYCAiIiJmCgYiIiJi\npmAgIiIiZgoGIiIiYqZgICIiImYKBiIiYjNRsQmcS88h8WwWU8OjiIpNsHdJcgWO9i5AREQqp6jY\nBMI2HDC/PpmUYX4d2LyWvcqSK1CPgYiI2MTGyKMlHD9WpnVI6ZS6x+D06dMsX76c7du3ExcXB0Dd\nunXp3Lkzw4YNo06dOlYvUkREKp645Mxij8enZJRxJVIapeoxiIiIYPjw4fj7+7No0SIiIyOJjIxk\n4cKF1K1blxEjRvDJJ5/YqlYREalA/Lxdiz1ex6taGVcipVGqHoM///yTDRs24OTkdNHxG264gRtu\nuIFBgwYxb948qxYoIiIVU3BQw4vmGPx9vIEdqpGrZTCZTCZ7F2ErSUnn7V1Cpebj4642tjG1cdlQ\nO9uG8ej/sf/Rp/kocCCnPOvi7+NGcFADTTy0ER8fd6vc56qGEnJycsjNzb3oWGJiolUKEBGRysl1\n0QK6HvyRNueP41uzKq+OaK9QUAFccSjh/fff57vvvsPR0ZF69eoxefJkXFxcmDBhAh988EFZ1Cgi\nIhVNfj6OMXvJb3IDu1t00SNwFcgVg8GWLVv46KOPANixYwdPPPEE06ZNs3VdIiJSkTk6cu6r7zGe\njsf0uR5PrEiuGAwKCwvJz8/H0dGRjh070qRJEyZPnszRo0fLoDwREamwjEYK/eoCCgYVyRWDwcSJ\nEzl79iw+Pj4A1K5dm7CwML788kubFyciIhWP66zXMGRkkPl8KCb36vYuR0rpisM+t956qzkUXODs\n7Ezfvn1tVpSIiFRMhpQUXJcuosr6TzE5V7F3OWIBi/ZKWLx4MatWrcLDw4PmzZvTokULWrRoQVBQ\nkLXrExERO3l+yU4A5ozpcNXXVA17G0NmJlkvvgxVFAwqIouCQUREhHmFw9jYWA4cOMCKFSsUDERE\nrmOGc2ep+l4Yhd4+ZD38qL3LEQtZFAwaNWqEv78/AP7+/txzzz1WLUpERCqequ+FYUw/T/pzE8G1\n+OWQpfyz6NHSJ598knnz5pGTk2PtekREpCLKzqbqsqUU1qxJ1rAR9q5GroFFPQbr16/nyy+/5MMP\nP+Smm26iZcuWtGzZkvvvv/+q75GVlUVoaCgpKSnk5OQwZswYbr75ZiZOnEhBQQE+Pj7MmTMHZ2dn\nNmzYwMqVKzEajQwYMID+/ftbUraIiNiKiwup69ZjPHYM3NzsXY1cA4uCQWRkJDt37sTR0ZHffvuN\n/fv3s3v37lIFg++//56WLVvy+OOPc+rUKR577DECAgIYPHgwPXv25M033yQiIoKQkBDefvttIiIi\ncHJyol+/ftx99914eHhYUrqIiNhIfuu20LqtvcuQa2RRMGjSpAkuLi44ODjQpk0b2rRpU+p79OrV\ny/zn+Ph4atWqRVRUFK+88goA3bp1Y/ny5TRq1IhWrVrh7l60OURAQADR0dHceeedlpQuIiJW5vzN\nFgpq1aGgVWt7lyJWYFEwqFmzJhMmTCA0NJRata5tQ4xBgwZx+vRp3nnnHYYPH46zszMAXl5eJCUl\nkZycjKenp/l8T09PkpKSrure1tppSkqmNrY9tXHZUDtfysHBAFyhbbKzYcI4SEuDuDhwv/Tcq7qP\nlBsWBQNfX19iY2O5//77cXNzM69jMGrUqFLf6+OPP+bgwYM8//zz/HMH6JJ2gy7NLtHaRtW2tFWt\n7amNy4bauXgFBUV/316ubVzeD8c9Lo7Mp8aRkQ1kX3puQYEJBweD2tjGrBW8LAoGEyZMMP85Li7O\nvJZBaezfvx8vLy/q1KlDs2bNKCgooFq1amRnZ+Pi4kJCQgK+vr74+vqSnJxsvi4xMZG2bTWGJSJi\nd3l5uC6aj8nFhcwnxpZ42pwxHRS+KpBSPa7Ys2dPxo8fz2effUZiYiIAfn5+dO/enXHjxpXqg3/+\n+WeWL18OQHJyMpmZmXTo0IEtW7YAsHXrVjp37kybNm2IiYkhLS2NjIwMoqOjadeuXak+S0RErK9K\nxFocThwna8ijmK5xWFnKj1L1GGzatImYmBh27NjBc889R0ZGBoGBgXTs2JH27dtTpRTLXw4aNIgX\nX3yRwYMHk52dzdSpU2nZsiWTJk1i7dq1+Pn5ERISgpOTE+PHj2fEiBEYDAaeeuop80REERGxk/x8\nXBfMxeTkRNbY/9i7GrEig6kUg/YHDhygRYsW5tfp6ens2rWLHTt28PPPP/PFF1/YpEhLqdvKttQ1\naHtq47Khdi7eZfdKSE+n2qwZAGRMn3nFe6mNbc8ucwyee+450tPT6dChA507d6ZTp050796d7t27\nW6UYERGpINzcrioQSMVTqmCwZcsW4uPjiYyM5Mcff2TOnDn4+PjQuXNnOnfurLF/EZFKIio2gXPp\nORQUmpgaHkVwUEMCmxfNIzAkJmLy8QGDwc5Vii2UaiihOD/88APvvfceP//8MwcPHrRWXVahbivb\nUteg7amNy4ba+WJRsQmEbbj0SbPRfVoQ2MwXj7s6g9HIuY1fX/XWympj27Pb44pnzpwhMjKSHTt2\nsGfPHnx9fQkMDCz1UwkiIlI+bYw8WsLxY3Q+GY3T/l/J7tvvqkOBVCylCgZ9+vQhMzOT4OBg7rvv\nPqZOnYqLi4utahMRETuIS84s9nh8SgauEXMAyBw3odhzpOIrVTAYOHAgkZGRbN68maNHj3L8+HGC\ngoJo0KCBreoTEZEy5uftysmkjEuPu5hw2vMzOcF9KGjW3A6VSVkoVTB4+OGHefjhhyksLGT//v3s\n3LmTadOmkZycTMuWLXnjjTdsVaeIiJSR4KCGxc4x6Le36JH0zGfVW1CZWbQkstFopFGjRpw+fZrk\n5GTOnDlDdHS0tWsTERE7uPD0wXtfxlJQaMLfx437WnjQ9cON5HS/p2h7Zam0ShUMoqKi2LlzJzt3\n7uTYsWO0a9eOjh078uijj1KvXj1b1SgiImUssHktIrYdBuDVEe0BOPO/XzGeO2vPsqQMlCoYvPHG\nG3Tq1IkJEyYQEBCAk5OTreoSEZHyxsWFwtp17F2F2FipgsG0adOuuLPhvn37aNOmzTUVJSIi5Yfb\ns2MpuPFmskY9CQ4O9i5HbKxUuysuWbKE+fPnc+bMmUveO3v2LPPnz2fJkiVWK05EROyrfvwhqn74\nAVU2fQHGUv3KkAqqVD0G77zzDsuXL+e+++6jbt261KlT1KUUFxfH6dOneeyxx1i6dKlNChURkbLX\n58ePAMh4bqKWQL5OlCoYGI1GRo4cyfDhw4mJiSE+Ph6AOnXq0KpVKxzUxSQiUmn4JR6l3cEd5N0S\nQN4dd9q7HCkjFj2u6ODgQNu2ba8430BERCqu3j9+DEDmc5PUW3Ad0YCRiEgl8/ySnTy/ZOc13cPh\n0J8E7v+B47Uak3tPDytVJhWBRT0GxTl//jzu7tbZ2UlEROyroMkNLHzoZTJc3Bil3oLrisXB4NCh\nQ5w9W7TQRW5uLjNmzGDz5s1WK0xEROzIYGDvTbfbuwqxA4uCwYwZM9ixYwfJycnUr1+fEydO8Nhj\nj1m7NhERKWvZ2VSbOYOsJ8fauxKxE4vmGMTExLB582ZuvvlmPvnkE5YvX05WVpa1axMRkTJW9d0l\nuC5ZSNUwrUlzvbIoGDg7OwOQl5eHyWSiZcuW2kRJRKSCMyQk4Dp/LoVeXmSOe87e5YidWDSU0KhR\nIz788EPatWvH8OHDadSoEefPn7d2bSIiUoaqvfEqxox0zr88HVMND+aM6WDvksQOLAoGr7zyCqmp\nqVSvXp2NGzeSkpLC6NGjrV2biIiUEcdf9+Ly0WrymzUne8ij9i5H7MiiYJCXl8fGjRuJj49nwoQJ\n7N27l5o1a1q7NhERKSPVXpmKwWQiffpMcLTak+xSAVk0x2DatGkcP36cqKgoAGJjYwkNDbVqYSIi\nUnbSZ80j44Wp5HW5w96liJ1ZFAyOHDnC5MmTcXFxAWDw4MEkJiZatTARESk7BTc0JfM/E+xdhpQD\nFgUDx7+6mQx/rYaVmZlJdna29aoSEZEyUSViLY67o+xdhpQjFgWDHj168Oijj3Ly5ElmzJhBSEgI\nvXv3tnZtIiJSSlGxCZxLzyElLZup4VFExSaUeK4hIQG355+lxrCHQGvRyF8smmEyZMgQWrduze7d\nu3F2dubNN9+kZcuW1q5NRERKISo2gbANB8yvTyZlmF8HNq91yfnVZk43P55I1aplVqeUbxYFg08+\n+YQHH3yQ1q1bW7seERGx0MbIoyUcP3ZJMHCM2YfLmlV6PFEuYdFQwtdff60FjUREypm45Mxij8en\nZFx8wGSi2pRQPZ4oxbLo25Cdnc2dd95Jo0aNcHJywmQyYTAY+PDDD61dn4iIXCU/b1dOJmVccryO\nV7WLXjtv/ALnyB3k9OilxxPlEhYFgzFjxlxyLC4u7pqLERERywUHNbxojsHfxxtc9Drvjm5kjhtP\n9kMPl1VpUoFYFAzat28PQE5ODlu2bOGTTz7h8OHDhISEWLU4ERG5ehfmEbz3ZSwFhSb8fdwIDmpw\nyfwCk5s7GS++bI8SpQKwKBjs3buXTz/9lE2bNlFYWMj06dO55557rF2biIiUUmDzWkRsOwzAqyPa\nX/SeISGBKls3kz14KDg42KM8qQBKNflw2bJl9OrVi2effRZPT08++eQT6tevT3BwME5OTraqUURE\nrKDazOm4j3+GKp9F2LsUKcdK1WOwYMECbrjhBqZOncrtt98O/L36oYiIlF//fDwxJ+RBe5cj5Vip\ngsG2bdv47LPPePnllyksLOSBBx4gLy/PVrWJiIg1mExUe2ly0eOJr7yuxxPlsko1lODj48OoUaPY\nsmULr7/+OsePH+fUqVM88cQT/PDDD7aqUUREroHzxi9w3vkTOff2JO+OO+1djpRzFi1wBHDbbbcx\nc+ZMtm/fzh133MHbb79tzbpERMQacnJwmzYFk5MTGdNm2LsaqQBK1Z/Us2dPmjdvTqdOnejYsSO+\nvr64ubkxaNAgBg0aZKsaRUTEUs7OZLzwEsaTJylo0tTe1UgFUKpgsGnTJmJiYtixYwfPPfccGRkZ\nBAYG0rFjR9q3b0+VKlVsVaeIiFjCYCCnb397VyEVSKmGEgwGA61bt+bJJ59k9erVrFq1inbt2vHd\nd9/Rr18/W9UoIiKlcPZ8NmfPZ1Pl03UYzqTYuxypYK5pd0U3Nze6d+9O9+7drV2XiIhcg0YJR3B/\nczx5t3cgdf1me5cjFYh2VxR+3AwvAAAgAElEQVQRqWRqulXhie0rMJhMZD77vL3LkQrGKrsrXqDd\nFUVE7K/dwR3cfCxGjyeKRay2u6KIiNifIS2Vh74KI9/ooMcTxSIWDSW0bduWP//8kx9//NH8NEKb\nNm2sXZuIiJRStVem4pWWxMZOA/V4oljEomAwbdo0jh8/TlRUFAAHDhwgNDTUqoWJiEjpZY16ku1t\n72ZDl4fsXYpUUBYFgyNHjjB58mRcXFwAGDx4MImJiVYtTERESq/gppsJDxlPgaN2vBXLWBQMHP/a\ngOPCzoqZmZlkZ2dbryoRESkV11mv4Rizz95lSCVg0eTDHj168Oijj3Ly5ElmzJjBjz/+yODBg0t9\nn9mzZ7Nnzx7y8/MZPXo0rVq1YuLEiRQUFODj48OcOXNwdnZmw4YNrFy5EqPRyIABA+jfX6t4iYhc\n4PzFeqrNm4XzTz+yZeZKzqXnUFBoYmp4FMFBDQlsXsveJUoFYjCZTCZLLvz111/ZvXs3zs7OBAQE\n0LJly1Jdv2vXLsLDw1m2bBlnz57lgQceICgoiC5dutCzZ0/efPNNateuTUhICA888AARERE4OTnR\nr18/Vq9ejYeHxxU/IylJay3Yko+Pu9rYxtTGZaMit7MhKQnPLu0xZGSweeXXLI26dKXD0X1a2D0c\nVOQ2rih8fNytch+LhhJ69+7N5s2badasGYMGDSp1KICi3RnfeustAKpXr05WVhZRUVHcddddAHTr\n1o3IyEj27dtHq1atcHd3x8XFhYCAAKKjoy0pW0SkcjGZcH/+PxhTUsh48WW+OFL8kO7GyGNlXJhU\nZBYNJaxYsYKoqCi2bNnC7Nmz8fHxoVOnTgwbNuyq7+Hg4ICrqysAERERdOnShZ9++glnZ2cAvLy8\nSEpKIjk5GU9PT/N1np6eJCUlXdVnWCs9ScnUxranNi4bFbKdP/wQNn0BXbrg9sJE4iZ9Wexp8SkZ\n5eLnKw81yJVZFAy8vb0JDg7m1ltvZffu3WzatImwsLBSBYMLvvnmGyIiIli+fDn33HOP+XhJIxyl\nGflQt5VtqWvQ9tTGZaNCtnNhITVfn4mDazXOzF1EYUoGfl6unEzKuOTUOl7V7P7zVcg2rmCsFbws\nCgYvvPACJ06cwMfHh1tvvZVnn32Wm266qdT32b59O++88w7vvfce7u7uuLq6kp2djYuLCwkJCfj6\n+uLr60tycrL5msTERNq2bWtJ2SIilYfRyLn1m3CM+ZXCho0ACA5qSNiGA5ecGhzUoKyrkwrMojkG\nmZmZALi5ueHh4XFRV//VOn/+PLNnzyYsLMw8kbBDhw5s2bIFgK1bt9K5c2fatGlDTEwMaWlpZGRk\nEB0dTbt27SwpW0SkcsjPB8BUw4O8Tl3MhwOb12J0nxY4GIseJff3cSsXEw+lYrGox2DBggUA/P77\n7+zevZvJkydz6tQpNm+++q09N23axNmzZ/nPf/5jPjZz5kymTJnC2rVr8fPzIyQkBCcnJ8aPH8+I\nESMwGAw89dRTuLtrnEpErk/GE8fxeCCY9Jenk9s75JL3A5vXImLbYQBeHdG+rMuTSsCiYJCens6e\nPXvYvXs30dHRmEwm7r777lLdY+DAgQwcOPCS4ytWrLjkWI8ePejRo4clpYqIVB6Fhbj/5ykcjh/D\nkHHpXAIRa7AoGNx///106NCBDh068Pjjj1/VmgIiInJtXFa8h/P2H8i5tyc5A0u/qJzI1bBojsGW\nLVu49dZb2bNnD4sXL2bjxo3WrktERP7BeOQwbtOnUlizJufnLoS/lqQXsTaLegxee+01zpw5Q/v2\n7TGZTGzatIlffvmFKVOmWLs+EREpKKD6M09iyMzk/PzFmGppMqHYjkXB4M8//2T16tXm10OGDLFo\nrwQREbkyQ1oqJgcHsvs8QE7Ig/YuRyo5i4JBXl4ehYWFGI1FIxEFBQUUFBRYtTARESliqulJ6mcb\nITNTQwhicxYFg65du9KvXz9uu+02AKKioujVq5dVCxMRue7l5+N4IIb8NreA0QhubvauSK4DFu+u\nuHfvXvbt24fBYKBt27a0bt3a2rVdMy2/aVta4tT21MZlo7y2s+v8ObjOeo3zYcvJub+vvcu5JuW1\njSsTuyyJnJ6ezpIlSzhy5Ajt2rVj2LBhODpa1OkgIiKX4bA/Bte5Myn0rUVu1272LkeuI6V6XHHa\ntGlA0eJEhw8fZvHixbaoSUTkuvT8kp08v2Qn5OZSfexoDHl5nF+wGJNHTXuXJteRUv1z/9SpU8yd\nOxeALl26WLSbooiIXJ7rvJk4xu4na+gw8u4s3aqyIteqVMHgn8MGDg4OVi9GROR61yDuT1zD51NQ\nvwEZr7xm73LkOlSqYGD412My/34tIiLXJt67HtnDRpAT3AeTmzaMk7JXqmDwyy+/cMcdd5hfp6Sk\ncMcdd2AymTAYDGzbts3K5YmIXF9ynV1If2OuvcuQ61ipgsFXX31lqzpERK5rVT7+kLt2x/Dtbb3t\nXYpc50oVDOrWrWurOkRErluOe/7HL2Fr2dWuLynns5kaHkVwUEMCm2tPBCl7WoRARMSOjKfj2ffy\nfObd+x/zsZNJGYRtOACgcCBlzqJtl0VExAqys6k+/GE+ual7sW9vjDxWxgWJKBiIiNiHyYTbpOdw\n2vMzx73rF3tKfEpGGRclomAgImIfeXkY09LIa3MLft7Fb45Ux6taGRclomAgImIfzs6khX9A6rrP\nCe7YqNhTgoMalHFRIpp8KCJSpozHj+H0y56i3RKNRkweNQn0KHrvvS9jKSg04e/jRnBQA008FLtQ\nMBARKSsZGdR45CEcY/dztl598gPamd8KbF6LiG2HAXh1RHt7VSiioQQRkTJhMlH9mSeLNkcaNuKi\nUCBSnigYiIiUAdcFc6nyxefk3t6B9Bmz7F2OSIkUDEREbMx5y2ZcZ86goK4/aeGrwNnZ3iWJlEjB\nQETExlw+Wg0uLqStXIPJx8fe5YhcliYfiojYWNp7K3Hc/yv5rdvauxSRK1IwEBGxhYICHPf9UjTJ\n0NGR/LYBV7xkzpgOZVCYyOVpKEFExAaqvfYKHr264/zlBnuXIlIqCgYiIlZW5ZP/4rp4AQWNGpPX\nuYu9yxEpFQUDERErctz3C+7PjqXQvTppH3yMqYaHvUsSKRXNMRARsRJDYiLVHx0MOTmcf28lBU1v\ntHdJIqWmYCAico2eX7ITgKVH1+EQd4r0F18m956edq5KxDIKBiIiVpL+yuvk39yc7OEj7V2KiMU0\nx0BE5Bp5nz1d9IeqVcl+7HEwGOxbkMg1UDAQEbkGLsuXMXPRSNr+FmnvUkSsQkMJIiIWcvlgBb+s\n+IyPh87nhKc/v4dHERzUkMDmtexdmojFFAxERCzgsmYVvyxbx5zgCeZjJ5MyCNtwAEDhQCosDSWI\niJRSlbVrcHt2LP8NGljs+xsjj5VxRSLWo2AgIlIaBQVUXf4upho1OOFVr9hT4lMyyrgoEevRUIKI\nSGk4OJD6388xnjyJ364MTiZdGgLqeFWzQ2Ei1qEeAxGRq+C86Usco3YBYKrhQUGLlgQHNSz23OCg\nBmVYmYh1qcdAROQKnDdvpPrIRyj08ubM7n1QtSrw9wTD976MpaDQhL+PG8FBDTTxUCo0BQMRkctw\n3rqZ6iMfAecqpC1baQ4FFwQ2r0XEtsMAvDqivT1KFLEqDSWIiJTA6buvqf7YUHB0JHXNOvJvD7J3\nSSI2p2AgIlIMp59+pMajg8FoJHXVWvI6dLJ3SSJlQsFARKQYBQ0bUdCwEakrPyKvyx32LkekzGiO\ngYjIP+Xng6Mjhf71OPv9TnC88l+Tc8Z0KIPCRMqGegxERP7iuDuKmp1uw+H33/46oH87yfVHwUBE\nBHDc8z9qDOqLw7GjOBw5bO9yROzG7sHgjz/+oHv37qxevRqA+Ph4hg4dyuDBgxk3bhy5ubkAbNiw\ngQcffJD+/fuzbt06e5YsIpXNzz9TY8ADGLIySQtbTm7PYHtXJGI3dg0GmZmZTJ8+naCgvx8BWrhw\nIYMHD2bNmjU0aNCAiIgIMjMzefvtt3n//fdZtWoVK1eu5Ny5c3asXEQqC8df98Ldd2PISOf82++S\n2+cBe5ckYld2DQbOzs4sW7YMX19f87GoqCjuuusuALp160ZkZCT79u2jVatWuLu74+LiQkBAANHR\n0fYqW0Qqi4IC3J8YAampnF/0Djl9+9u7IhG7s+vMGkdHRxz/NbknKysLZ2dnALy8vEhKSiI5ORlP\nT0/zOZ6eniQlJZVprSJSCTk4MO3u56jf9igP9B9k72pEyoVyPeXWZDKV6vi/+fi4W7McKYba2PbU\nxlaWlARPPw1vvAGNGhHv15h4v8aMUjvbnL7LFUO5Cwaurq5kZ2fj4uJCQkICvr6++Pr6kpycbD4n\nMTGRtm3bXvFeSUnnbVnqdc/Hx11tbGNqY+ty2rUT99GP4RAfR6ZXLTKmzaCgwISDg0HtbGP6Ltue\ntYKX3Z9K+LcOHTqwZcsWALZu3Urnzp1p06YNMTExpKWlkZGRQXR0NO3atbNzpSJSYRQWUvWtedR4\nIBhjYgLpU14hY+qrRMUmcC49h8SzWUwNjyIqNsHelYrYnV17DPbv38+sWbM4deoUjo6ObNmyhblz\n5xIaGsratWvx8/MjJCQEJycnxo8fz4gRIzAYDDz11FO4u6tLSkSuzJCcTPWnHsf5+28pqONHWtgK\n8m8PIio2gbANB8znnUzKML/WtslyPTOYrnbAvgJSt5VtqWvQ9tTG18546iQ17+xIXkA7zi9+F5OX\nFwBTw6M4mZRxyfn+Pm7aPtkG9F22PWsNJZS7OQYiItessBDjqZMU1qtPYV1/zm3+loKGjcH49+hp\nXHJmsZfGp1waFkSuJ+VujoGIyLUwJCVRY1BfPO67B0NKCgAFjW+4KBQA+Hm7Fnt9Ha9qNq9RpDxT\nMBCRSsNp50/UvLMjztu+I79lKzCUfG5wUMMSjjewSW0iFYWGEkSk4isowPWtebjOfh0MBtJfnkHW\nk2Mv6SX4pwsTDN/7MpaCQhP+Pm4EBzXQxEO57ikYiEiF5zZ5AlXfD6egrn/RUwftA6/qusDmtYjY\ndhgHB4MmHIr8RcFARCq8rOGPY0xJ4fyc+Zg8vUp17ZwxHTRjXuQfNMdARCqeggJc58/B4beDRS+b\nNSct/INShwIRuZR6DESkQjEkJlL9yZE4b9+GU1QkqR9/au+SRCoV9RiISIXhtP2HoqcOtm8j596e\npC19z94liVQ6CgYiUv4VFOA65w1q9OuD8UwK6a++TtoHH2Oq6Xnla0WkVDSUICLlniH1HC6r3qew\nrj9py94n/9bb7F2SSKWlYCAi5ZYhLRVT9RqYPL1I/XAdhf7+6iUQsTENJYhI+VNQgOus1/C8PQBj\n3KmiQ61aKxSIlAEFAxEpV4wJp6nR/36qzZuFqWpVDGfO2LskkeuKhhJEpMxExSb8YwniagQHNbxo\nCWKnH76n+pMjMSYnkdMjmPMLl2DyqGnHikWuP+oxEJEyERWbQNiGAxQUmgA4mZRB2IYDRMUmAOCy\neiU1BoRgSD1H+oyZpK1co1AgYgcKBiJSJjZGHi3h+DEAcjt1Ib9FK859uZWsUWPAcJmtEUXEZjSU\nICJlIi45s9jj8cnpABQ2bMS5b7crEIjYmXoMRKRM+Hm7Fnu83rk4KCgoeqFQIGJ3CgYiUiaCgxoW\ne7xX5ybg4FC2xYhIiTSUICJlIvCGmji7nGD9SROnavrhn5dGz/sCaNeusb1LE5F/UDAQkbLh4MCd\nEYvpfvo0GZNfIvuxxzV0IFIOKRiIiG0UFlJl/acYzqSQPWI0ODiQ9u77FPr56TFEkXJMwUBErMtk\nwvmbLVR7fTqOB2IorOFB9qAhUK0aBc1b2Ls6EbkCTT4UEatx2rUTjz49qPHwABxi95PdbyBnv/4B\nqlWzd2kicpXUYyAiVuHw+2949OkBQE6PXmSEvqQeApEKSMFARK6opD0OHI4cwuTgSGGDhhTcdDMZ\nE0LJ7XYX+bcF2rtkEbGQhhJE5LJK2uNg3wuzqNnxNtxenWo+N3PiCwoFIhWcgoGIXFZJexysz/Gm\noFFjskMeLNN6RMS2NJQgIpdV0h4HJ3wacvbHKHDUXyMilYl6DESkWIbERFzC38XPq/g9Dur4uCsU\niFRCCgYi8reMDKpErKXGoL54tbkJ98kTuL96WrGnBgc1KOPiRKQsKO6LCIa0VNxCJ1Bl05cYMjMA\nyAu4lZwHB9DuzlsY3bKAjZHHiE/JoI5XNYKDGhDYvJadqxYRW1AwELkemUw47vuFwlq1Kazjh8nN\nHaeoSAp9fMjuN5acfgMoaNLUfHqgNwoCItcJBQOR64jx6P/h8uk6qkSsxfHQn2SMn0TmpBfBaOTc\nl1sprF1HGxuJXOcUDESuAy5rVuHy4Qc4/S8KAJOLC9khfcnr2Nl8TmEdP3uVJyLliIKBSGWUlYUx\nOYnCevUBcNr2LY4/7ya3Szey+w0gN7g3Jvfqdi5SRMojBQORyqKgAKedP1ElYi1VvlhP/q3tSF23\nHoCM0JfIeOV19QqIyBUpGIhUcA5//I7LmlVU+SwCh/g4AAr865F3azswmcBgoLBxEztXKSIVhYKB\nSAVkOJ9mHgpw3vYtrksWUljDg6yhw8jpN5C8wCAwapkSESk9BQORcigqNoGNkUeJS8nEz8uV4KCG\n3F63ClW+WF/0RMFvsaTs+x2qVCH7gf4U1K1Hbvd7oEoVe5cuIhWcgoFIOXNhN8MLLuxm6L7lLboe\n+B6A3KCOGBMTKKxXH5OPD7nBve1VrohUMgoGIuWIISGBTV/vBy5dSyAiqD+3hnQhp29/89MGIiLW\npmAgYieGtFScftqO4697cYzZh+Ov+3BIOM2p/3wCRodLzj/hUZescUPsUKmIXE8UDERsrbAQh6NH\nin7xHzxAZuhLYDDgcPgQNYYNNp9W4FeXnB698DfmcJxLdzSs41WtLKsWkeuUgoGIDTgcjMVlzQc4\n/roPx5hfMaafN7+X/dBQChs2Ir9ZCzJemEpe67bkt2qDyccHgJ7/mmNwgXYzFJGyoGAgYqmcHBx/\nP1j0y//XvRhPx5P2wccAGBNO4xq2BJPBQEHTG8lt1Yb81m3Jb92Gwlq1i653cSHzPxMuue2FzYq0\nm6GI2IOCgchfzI8IJmfi5130iKD5l3F+PjgW/d/F+csNuM6fg+NvsRjy8szXm5ycMKSlYqpeg7x2\n7Tm78Wvym7eEaqUfAghsXovA5rXw8XEnKen8lS8QEbESBQMRSn5EsOryd+ka9QUmJ2fOffdT0Zum\nQhz//J381m3Iv9AT0Ko1+Tc3/3sdATc38m8LtMNPIiJybRQMRICNkUeLPf65yY9up06R3/YWKCwE\no5Hce3uRfCTO3IMgIlKZ6G82ESAuObPY4yd8G5Fy6MTFyws7O5dRVSIiZa/CLKb++uuvM3DgQAYN\nGsSvv/5q73KkkvHzvvTxQIA63m7ac0BErisV4m+83bt3c+zYMdauXctrr73Ga6+9Zu+SpJIJDmpY\nwnE9Iigi15cKMZQQGRlJ9+7dAWjSpAmpqamkp6fj5uZm58qkstAjgiIiRSpEMEhOTqZFixbm156e\nniQlJSkYiFVdeERQROR6ViGCwb+ZTKarOs/Hx93GlYja2PbUxmVD7Wx7auOKoULMMfD19SU5Odn8\nOjExEZ+/lo8VERER66kQwaBjx45s2bIFgAMHDuDr66thBBERERuoEEMJAQEBtGjRgkGDBmEwGHj5\n5ZftXZKIiEilZDBd7YC9iIiIVHoVYihBREREyoaCgYiIiJhViDkGpfH666+zb98+DAYDL7zwAq1b\nt7Z3SRVWVFQU48aNo2nTpgDceOONjBw5kokTJ1JQUICPjw9z5szB2dmZDRs2sHLlSoxGIwMGDKB/\n//52rr78++OPPxgzZgzDhg1jyJAhxMfHX3Xb5uXlERoaSlxcHA4ODrzxxhvUq1fP3j9SufPvNg4N\nDeXAgQN4eHgAMGLECO644w618TWYPXs2e/bsIT8/n9GjR9OqVSt9j63s32383Xff2fZ7bKpEoqKi\nTKNGjTKZTCbToUOHTAMGDLBzRRXbrl27TE8//fRFx0JDQ02bNm0ymUwm07x580wffvihKSMjw3TP\nPfeY0tLSTFlZWabg4GDT2bNn7VFyhZGRkWEaMmSIacqUKaZVq1aZTKbSte2nn35qmjZtmslkMpm2\nb99uGjdunN1+lvKquDaeNGmS6bvvvrvkPLWxZSIjI00jR440mUwm05kzZ0xdu3bV99jKimtjW3+P\nK9VQQklLJ4v1REVFcddddwHQrVs3IiMj2bdvH61atcLd3R0XFxcCAgKIjo62c6Xlm7OzM8uWLcPX\n19d8rDRtGxkZyd133w1Ahw4d1N7FKK6Ni6M2ttxtt93GW2+9BUD16tXJysrS99jKimvjgoKCS86z\nZhtXqmCQnJxMzZo1za8vLJ0sljt06BBPPPEEDz30EDt27CArKwvnv7Yd9vLyIikpieTkZDw9Pc3X\nqN2vzNHRERcXl4uOlaZt/3ncaDRiMBjIzc0tux+gAiiujQFWr17NI488wrPPPsuZM2fUxtfAwcEB\nV9einUkjIiLo0qWLvsdWVlwbOzg42PR7XOnmGPyTSU9iXpOGDRsyduxYevbsyYkTJ3jkkUcuSqol\nta/a/dqVtm3V5lfn/vvvx8PDg2bNmvHuu++yePFibrnllovOURuX3jfffENERATLly/nnnvuMR/X\n99h6/tnG+/fvt+n3uFL1GGjpZOuqVasWvXr1wmAwUL9+fby9vUlNTSU7OxuAhIQEfH19i233K3Xf\nyqVcXV2vum19fX3NvTJ5eXmYTCbzv9KkZEFBQTRr1gyAO++8kz/++ENtfI22b9/OO++8w7Jly3B3\nd9f32Ab+3ca2/h5XqmCgpZOta8OGDYSHhwOQlJRESkoKffv2Nbfx1q1b6dy5M23atCEmJoa0tDQy\nMjKIjo6mXbt29iy9QurQocNVt23Hjh356quvAPj+++8JDAy0Z+kVxtNPP82JEyeAojkdTZs2VRtf\ng/PnzzN79mzCwsLMM+T1Pbau4trY1t/jSrfy4dy5c/n555/NSyfffPPN9i6pwkpPT2fChAmkpaWR\nl5fH2LFjadasGZMmTSInJwc/Pz/eeOMNnJyc+OqrrwgPD8dgMDBkyBD69Olj7/LLtf379zNr1ixO\nnTqFo6MjtWrVYu7cuYSGhl5V2xYUFDBlyhSOHj2Ks7MzM2fOpE6dOvb+scqV4tp4yJAhvPvuu1St\nWhVXV1feeOMNvLy81MYWWrt2LYsWLaJRo0bmYzNnzmTKlCn6HltJcW3ct29fVq9ebbPvcaULBiIi\nImK5SjWUICIiItdGwUBERETMFAxERETETMFAREREzBQMRERExEzBQMTOZs+ezdChQxkwYAAtW7Zk\n6NChDB06lM8//9zepZXo+++/Jy0tzSb3XrduHa+//jrHjh3jzjvvNB8/c+YMvXv35ocffrjs9V26\ndOHkyZMlvv/MM89cdsnu6dOn89lnn5W+cJFKolIviSxSEUycOBGAkydPMnjwYFatWmXniq5s+fLl\nNG3alOrVq1v1vidPnmTFihWsX7+euLg48/GMjAyeeOIJRo8eTdeuXa/pMxYuXHjZ9ydNmkTv3r0J\nCgqidu3a1/RZIhWRgoFIOZabm8u0adM4ceIEGRkZ9OnTh2HDhrFu3Tp27dpFfn4+Bw8eJCQkhKys\nLHbv3o3BYGDFihUkJiYyatQobr/9dv744w+MRiPz58/H19eXyMhI3n77bfPyqDNmzKBu3bp06dKF\nPn36cOrUKebPn8/8+fOJiooCoG7dusyaNYuPP/6YX375heeee46ZM2cybNgw1qxZg7+/Pzt37mTp\n0qWsWrWKhx56iFatWhEbG8vq1atL/Mx/Cg8PZ+DAgTg5OZmP5eXl8cwzzxASEsJ9991nPn65+xUU\nFNCxY0c+//xz87Lod911F++99x6PPvooa9asITIykv/973/k5eVx9OhR6tWrx1tvvYWzszMDBgxg\n5cqVTJo0ydb/iUXKHQ0liJRjK1asoG7duqxatYr//ve/rF+/nkOHDgFFK/vNmTOHZcuWsWjRIrp2\n7cratWsxGAzs2rULgKNHjzJo0CA++ugjAgICeP/998nIyGDatGksWbKEDz/8kEGDBjFnzhzzZzZu\n3Jj58+eTm5tLtWrVWLNmDR9//DHJycns3LmTIUOG4OnpyZtvvknjxo0vW7+7uzurV6++4mdesH37\ndjp37mx+bTKZCA0NxWQyMXjwYPPxK93PwcGBe++9l61btwJFW9J6eXldtHocwN69e5k5cyaffPIJ\nMTEx/PHHH0DR8urbt2+/qv9GIpWNegxEyrGoqCiSk5PNv+jz8/M5duwYAK1atcLZ2ZnatWtjMpkI\nCAgAoHbt2pw/fx4Ab29v82YrAQEBfPzxx/z+++8kJyfz1FNPAUX/unZ0/Puvggu7tF3YaGXw4ME4\nOTlx6NAhzp49W6r6L9zrSp95wenTpy9arjUhIYHmzZuzfv16Nm3aRK9eva76fr1792b+/Pk8/PDD\nbNq0qdhlulu3bk2VKlXM7ZaamgqAn58fp06dKtXPKlJZKBiIlGPOzs4888wzdO/e/aLj69atu+gX\nodFoxGj8uwPwwkrnhYWFFx0zGAw4Ozvj7+9f4lyGC934u3fvZsOGDaxbt46qVasyZsyYYs83GAzm\nP+fl5RV7ryt9Zklq167NiBEjuPvuu3n44Ydp3LgxN99881Xd75ZbbiEhIYHk5GS++eYb1q1bd8k5\n/w4TWiFeREMJIuVaQEAAmzdvBor+Vfzaa6+V6mmAs2fP8vvvvwOwZ88ebrrpJpo0aUJiYiKHDx8G\nYNeuXcX+0kxOTsbf35+qVaty4sQJfv31V3Jzc4GiIJKfnw+Am5sbp0+fNt+rOFf7mbVr1yY+Pv6S\n4/Xr1+fVV19l7NixnCaVFLsAAAFaSURBVD179qrv16NHDxYvXkzTpk3x9PS8fGP9Q1xc3CXzH0Su\nF+oxECnHhg4dyiuvvMLAgQPJz8+ne/fupXoSoE6dOqxdu5Y///wTgAULFlC1alXmzJnDpEmTcHFx\nwWAwMH369Euu7dKlCytXruShhx7ixhtv5Omnn2bJkiUEBgbSqVMnRo0axZw5c3jssccIDQ2lYcOG\ntG3bttg6rvYzO3fuzE8//VTs3IVu3boRExPDs88+S3h4+FXdr0+fPvTu3Zu5c+dedZsB7Ny586K5\nDiLXE+2uKFJJHTt2jOHDh/Pdd9/Zu5SrdvLkSUaNGsX69esvejKhLOXm5tK7d29WrlypxxXluqSh\nBBEpN/z9/Rk+fHixTyyUlVmzZvHEE08oFMh1Sz0GIiIiYqYeAxERETFTMBAREREzBQMRERH5/3br\nWAAAAABgkL/1vlEURRMDAGBiAABMDACABU66dySik+gRAAAAAElFTkSuQmCC\n",
            "text/plain": [
              "<matplotlib.figure.Figure at 0x7f711e41da20>"
            ]
          },
          "metadata": {
            "tags": []
          }
        },
        {
          "output_type": "stream",
          "text": [
            "Model parameters:\n",
            "B   = 0.09963406247586233  +/-  1.9283106622952189\n",
            "C   = 1.2520229787358318e-11  +/-  2.3100168937045132e-11\n",
            "n   = 4.064219058140722  +/-  0.24765692572520176\n",
            "\n",
            "Fit quality:\n",
            "chisq = 4.89900021341353 dof = 8\n"
          ],
          "name": "stdout"
        }
      ]
    },
    {
      "metadata": {
        "id": "8V3QC_0KwN9b",
        "colab_type": "code",
        "colab": {}
      },
      "cell_type": "code",
      "source": [
        ""
      ],
      "execution_count": 0,
      "outputs": []
    }
  ]
}