analysis.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# SIMULATIONS\n",
    "## Import and helper functions"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "import workload_simulation\n",
    "import forgetting_curve as fc\n",
    "import matplotlib\n",
    "import importlib\n",
    "import numpy\n",
    "from operator import itemgetter"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "importlib.reload(workload_simulation)\n",
    "import workload_simulation as ws"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "import matplotlib.pyplot as plt\n",
    "import numpy as np\n",
    "import pickle\n",
    "SimulationResult = ws.SimulationResult\n",
    "%matplotlib inline"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "def adapt_simulations_to_newcardcost(sims, cost_new_card):\n",
    "    new_sims = []\n",
    "    for sim in sims:\n",
    "        pbt = sim.pbt\n",
    "        factor = sim.factor\n",
    "        workload = sim.w_load + cost_new_card\n",
    "        efficiency = sim.r_rate / workload\n",
    "        new_sims.append(SimulationResult(pbt, sim.r_rate, workload, efficiency, factor))\n",
    "    return new_sims"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [],
   "source": [
    "def get_simdata(fp, cardcost=0):\n",
    "    sims = ws.get_simdata(fp)\n",
    "    return adapt_simulations_to_newcardcost(sims, cardcost)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [],
   "source": [
    "def get_closest_factor(new_factor, total_output):\n",
    "    for r in total_output:\n",
    "        if r.factor >= new_factor:\n",
    "            return r"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [],
   "source": [
    "def printsim(sim, name, forcomparison=False):\n",
    "    name = name.center(10)\n",
    "    print(f\"{name} | {sim.factor:0.2f} | {sim.pbt:0.2f} | {sim.r_rate:0.2f} | {sim.w_load:0.2f} | {sim.efficiency:0.2f}\")\n",
    "\n",
    "def print_header():\n",
    "    print(\" \"*10+\" |factor| pbt  | ret  | load  | effi\")\n",
    "\n",
    "def analysis(total_output, new_factor=None):\n",
    "    factors = [output.factor for output in total_output]\n",
    "    r_rates = [output.r_rate for output in total_output]\n",
    "    w_loads = [output.w_load for output in total_output]\n",
    "    effs = [output.efficiency for output in total_output]\n",
    "    \n",
    "\n",
    "    \n",
    "    plt.scatter(factors, r_rates,  color='blue')\n",
    "    #t = np.polyfit(factors,r_rates,1)\n",
    "    #new_rates = [t[0]*f + t[1] for f in factors]\n",
    "    #plt.plot(factors,new_rates, color=\"red\")\n",
    "    plt.xlabel(\"Factors\")\n",
    "    plt.ylabel(\"Retention rate\")\n",
    "    if new_factor is not None:\n",
    "        plt.plot([new_factor, new_factor], [min(r_rates), max(r_rates)], color=\"green\", lw=1)\n",
    "    plt.plot([2.5, 2.5], [min(r_rates), max(r_rates)], 'k-', lw=1)\n",
    "    plt.show()\n",
    "    plt.scatter(factors, w_loads,  color='red')\n",
    "    plt.xlabel(\"Factors\")\n",
    "    plt.ylabel(\"workload\")\n",
    "    if new_factor is not None:\n",
    "        plt.plot([new_factor, new_factor], [min(w_loads), max(w_loads)], color=\"green\", lw=1)   \n",
    "    plt.plot([2.5, 2.5], [min(w_loads), max(w_loads)], 'k-', lw=1)\n",
    "    plt.show()\n",
    "    plt.scatter(factors, effs,  color='green')\n",
    "    if new_factor is not None:\n",
    "        plt.plot([new_factor, new_factor], [min(effs), max(effs)], color=\"green\", lw=1)  \n",
    "    plt.plot([2.5, 2.5], [min(effs), max(effs)], 'k-', lw=1)\n",
    "    t = np.polyfit(factors,effs, 3)\n",
    "    #new_effs = [t[0]*f*f*f + t[1]*f*f + t[2]*f + t[3] for f in factors]\n",
    "    #plt.plot(factors,new_effs, color=\"blue\")\n",
    "    plt.xlabel(\"Factors\")\n",
    "    plt.ylabel(\"efficiency\")\n",
    "    plt.show()\n",
    "    \n",
    "    \n",
    "    # PRINT INSIDE TABLE\n",
    "\n",
    "    best_productivity = max(total_output, key=lambda x:x.efficiency)\n",
    "    least_work = min(total_output, key=lambda x:x.w_load)\n",
    "    print_header()\n",
    "    if new_factor is not None:\n",
    "        sm2 = get_SM2_Point(total_output)\n",
    "        printsim(sm2, name=\"SM2\")\n",
    "        print(\"---\")\n",
    "        closest_new_factor = get_closest_factor(new_factor, total_output)\n",
    "        print(\"User Defined Factor\")\n",
    "        printsim(closest_new_factor, \"Absolute\")\n",
    "        printsim(compare(closest_new_factor, sm2), \"VS SM2\")\n",
    "        print(\"---\")\n",
    "    print(\"Max Eff\")\n",
    "    printsim(best_productivity, \"Absolute\")\n",
    "    printsim(compare(best_productivity, sm2), \"VS SM2\")\n",
    "    print(\"---\")\n",
    "    print(\"Least Work\")\n",
    "    printsim(least_work, \"Absolute\")\n",
    "    printsim(compare(least_work, sm2), \"VS SM2\")\n",
    "#analysis(total_output)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [],
   "source": [
    "def get_SM2_Point(results):\n",
    "    for r in results:\n",
    "        if r.factor == 2.5:\n",
    "            return r\n",
    "\n",
    "def compare(source, target):\n",
    "    return SimulationResult(source.pbt - target.pbt, \n",
    "                            source.r_rate - target.r_rate, \n",
    "                            (source.w_load)/(target.w_load), \n",
    "                            source.efficiency/target.efficiency, \n",
    "                            source.factor/target.factor)\n",
    "\n",
    "def compare_best_to_sm2(results):\n",
    "    o = get_best_prod(results)\n",
    "    s = get_SM2_Point(results)\n",
    "    return compare (o, s)\n",
    "\n",
    "def get_best_prod(results):\n",
    "    return max(results, key=lambda x: x.efficiency)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 23,
   "metadata": {},
   "outputs": [],
   "source": [
    "def plot_comparison(difs, c_wlds, c_rets, c_pbts, c_effs):\n",
    "    fig, ax1 = plt.subplots()\n",
    "    ax1.scatter(difs, c_wlds, color=\"red\")\n",
    "    ax1.set_xlabel('Initial Success Rate')\n",
    "    ax1.set_ylabel('Workload Fraction', color='red')\n",
    "    ax1.tick_params('y', colors='red')\n",
    "\n",
    "    ax2 = ax1.twinx()\n",
    "    ax2.scatter(difs, c_rets, color=\"blue\")\n",
    "    ax2.set_ylabel('Retention Rate Loss', color='blue')\n",
    "    ax2.tick_params('y', colors='blue')\n",
    "\n",
    "    fig.tight_layout()\n",
    "    plt.show()\n",
    "    \n",
    "    # copy-pasta\n",
    "    fig, ax1 = plt.subplots()\n",
    "    ax1.scatter(difs, c_wlds, color=\"red\")\n",
    "    ax1.set_xlabel('Initial Success Rate')\n",
    "    ax1.set_ylabel('Workload Fraction', color='red')\n",
    "    ax1.tick_params('y', colors='red')\n",
    "\n",
    "    ax2 = ax1.twinx()\n",
    "    ax2.scatter(difs, c_pbts, color=\"blue\")\n",
    "    ax2.set_ylabel('Success Rate Loss', color='blue')\n",
    "    ax2.tick_params('y', colors='blue')\n",
    "    plt.show()\n",
    "    \n",
    "    plt.scatter(difs, c_effs,  color='green')\n",
    "    plt.xlabel(\"Initial Success Rate\")\n",
    "    plt.ylabel(\"Efficiency Gain\")\n",
    "\n",
    "    fig.tight_layout()\n",
    "    plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 35,
   "metadata": {},
   "outputs": [],
   "source": [
    "def display_overall(results, sm2s, diffs, display_SM2=True):\n",
    "    #results are post-optimized selection for a range of factor\n",
    "    factors = [output.factor for output in results]\n",
    "    r_rates = [output.r_rate for output in results]\n",
    "    w_loads = [output.w_load for output in results]\n",
    "    effs = [output.efficiency for output in results]\n",
    "    pbt = [output.pbt for output in results]\n",
    "    c_wlds=[]\n",
    "    c_rets=[]\n",
    "    c_effs = []\n",
    "    c_pbts = []\n",
    "    for sim, sm2 in zip(results, sm2s):\n",
    "        comparesim = compare (sim, sm2)\n",
    "        c_wlds.append(comparesim.w_load)\n",
    "        c_rets.append(comparesim.r_rate)\n",
    "        c_effs.append(comparesim.efficiency)\n",
    "        c_pbts.append(comparesim.pbt)\n",
    "    \n",
    "    plt.scatter(diffs, factors, color=\"blue\")\n",
    "    #t = np.polyfit(np.asarray(pbt, dtype=float),np.asarray(factors, dtype=float),2)\n",
    "    #new_factors = [t[0]*p*p + t[1]*p + t[2] for p in pbt]\n",
    "    #plt.plot(pbt,new_factors, color=\"red\")\n",
    "    if display_SM2:\n",
    "        plt.plot([min(diffs), max(diffs)], [2.5, 2.5], color=\"black\", lw=1)\n",
    "    plt.ylabel(\"Optimal Factors\")\n",
    "    plt.xlabel('Initial Success Rate')\n",
    "    plt.show()\n",
    "    \n",
    "    plt.scatter(diffs, effs, color=\"green\")\n",
    "    plt.ylabel(\"Efficiency\")\n",
    "    plt.xlabel('Initial Success Rate')\n",
    "    plt.show()\n",
    "\n",
    "    fig, ax1 = plt.subplots()\n",
    "    ax1.scatter(diffs, w_loads, color=\"red\")\n",
    "    ax1.set_xlabel('Initial Success Rate')\n",
    "    ax1.set_ylabel('Workload After Opti', color='red')\n",
    "    ax1.tick_params('y', colors='red')\n",
    "\n",
    "    ax2 = ax1.twinx()\n",
    "    ax2.scatter(diffs, r_rates, color=\"blue\")\n",
    "    ax2.set_ylabel('Retention Rate After Opti', color='blue')\n",
    "    ax2.tick_params('y', colors='blue')\n",
    "\n",
    "    fig.tight_layout()\n",
    "    plt.show()\n",
    "    plot_comparison(diffs, c_wlds, c_rets, c_pbts, c_effs)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 36,
   "metadata": {},
   "outputs": [],
   "source": [
    "from ipywidgets import interact, interactive, fixed, interact_manual\n",
    "import ipywidgets as widgets\n",
    "import os"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 37,
   "metadata": {},
   "outputs": [],
   "source": [
    "length_to_path = {\"1year\": \"simunight\", \"10years\": \"simunight10years\", \"30years\": \"simunight30years\"}\n",
    "def analyse(success_rate, cost_card, new_factor, length_simulation):\n",
    "    path = length_to_path[length_simulation]\n",
    "    try:\n",
    "        fpath = f\"./{path}/{success_rate}.pkl\"\n",
    "        sims = get_simdata(fpath, cost_card)\n",
    "    except:\n",
    "        fpath = f\"./{path}/{success_rate}.0.pkl\"\n",
    "        sims = get_simdata(fpath, cost_card)\n",
    "        #TODO : code get_sim_data_from_settings\n",
    "    print(\"IM :\"+ str(new_factor/2.5))\n",
    "    return analysis(sims, new_factor)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 38,
   "metadata": {},
   "outputs": [],
   "source": [
    "def range_analysis(cost_card, length_simulation):\n",
    "    path = length_to_path[length_simulation]\n",
    "    bests = []\n",
    "    sm2s= []\n",
    "    diffs = []\n",
    "    for file in os.listdir(f\"./{path}/\"):\n",
    "        sims = get_simdata(f\"./{path}/{file}\", cost_card)\n",
    "        diffs.append(float(file.replace(\".pkl\",\"\")))\n",
    "        #print(file)\n",
    "        #print(max(sims, key=lambda x:x.efficiency))\n",
    "        bests.append(max(sims, key=lambda x:x.efficiency))\n",
    "        sm2s.append(get_SM2_Point(sims))\n",
    "    display_overall(bests, sm2s, diffs)\n",
    "        \n",
    "    \n",
    "    #display_overall(results)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Analyse results of past simulations\n",
    "### Range analysis\n",
    "\n",
    "For a broad range of initial success rate, find the factor that maximises efficiency and display the consequences on workload, retention and efficiency. \n",
    "\n",
    "The cost of creating new card (how many reviews could have been done instead) and the average length of the simulations are passed as extra hyperparameters.\n",
    "\n",
    "Black line represents SM2 point (with factor == 2.5)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 39,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "application/vnd.jupyter.widget-view+json": {
       "model_id": "ee11595101e74fa08f60b9a453d7cd69",
       "version_major": 2,
       "version_minor": 0
      },
      "text/plain": [
       "interactive(children=(IntSlider(value=5, description='cost_card', max=15), Dropdown(description='length_simula…"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "interact(range_analysis,\n",
    "         cost_card=widgets.IntSlider(min=0,max=15,step=1,value=5),\n",
    "         length_simulation = widgets.Dropdown(options=[\"1year\", \"10years\", \"30years\"],value='1year',description='length_simulation:')\n",
    "        );"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Analyse result of simulation with specific succes rate"
   ]
  },
  {
   "cell_type": "raw",
   "metadata": {},
   "source": [
    "TODO : \n",
    "- drop the new factor line (add a section where this could be chosen)\n",
    "- show comparison results in a readable HTML table"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 29,
   "metadata": {
    "scrolled": true
   },
   "outputs": [
    {
     "data": {
      "application/vnd.jupyter.widget-view+json": {
       "model_id": "b3cdb007041240bfabad5cf268978549",
       "version_major": 2,
       "version_minor": 0
      },
      "text/plain": [
       "interactive(children=(IntSlider(value=90, description='success_rate', max=97, min=70), IntSlider(value=5, desc…"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "interact(analyse, success_rate=widgets.IntSlider(min=70,max=97,step=1,value=90), \n",
    "         cost_card=widgets.IntSlider(min=0,max=15,step=1,value=5),\n",
    "         new_factor=widgets.FloatSlider(min=1.2,max=20,step=0.1,value=4),\n",
    "         length_simulation = widgets.Dropdown(options=[\"1year\", \"10years\", \"30years\"],value='1year',description='length_simulation:')\n",
    "        );"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## TRASH"
   ]
  },
  {
   "cell_type": "raw",
   "metadata": {},
   "source": [
    "## Run an analysis for specific options (need rework)\n",
    "\n",
    "nsims = 4000\n",
    "dif = 0.9\n",
    "ndays = 365\n",
    "onfail = \"reset\"\n",
    "\n",
    "list_of_factors = [i / 120 for i in range(3 * 60, 10 * 300) if (i % 30 == 0) or (i % 40 == 0) or (i % 12 == 0)]\n",
    "cost_new_card = 5\n",
    "\n",
    "vals_to_print = []\n",
    "total_output = []\n",
    "# trick to get this kind of factor 1.0, 1.1, 1.2, 1.25, 1.3, 1.3333333333333333\n",
    "# Why ? These factors produce sensible enough different simulations because intervals are integers \n",
    "# Intervals get rounded if they become float, so for the sake of the simulation, let's pick factors that\n",
    "# would result in neat intervals\n",
    "\n",
    "# How ? : take floats that would be close to integers when multiplied by common integers\n",
    "focus_output = []\n",
    "focus_output_too_high = []\n",
    "focus_output_too_low = []\n",
    "for factor in list_of_factors:\n",
    "    pbt = fc.get_pb_success_from_interval_modifier(factor / 2.5, dif)\n",
    "    l = ws.sim(pbt, nsims=nsims, onfail=onfail, ndays=ndays\n",
    "            difficulty=dif, factor=factor)\n",
    "    quotient = l[0] / l[1]\n",
    "    row = ws.SimulationResult(pbt, l[0], l[1], quotient, factor)\n",
    "    total_output.append(row)\n",
    "    \n",
    "## ANALYSIS\n",
    "total_output = adapt_simulations_to_newcardcost(total_output, cost_new_card)\n",
    "analysis(total_output)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Retention rate point loss is distributed between -3pt and -6pts whereas workload diminish fast with Initial Success rate.\n",
    "\n",
    "Two practical implications :\n",
    "1. Whatever your initial success rate, if you use the suggested IM, the retention loss is not large.\n",
    "2. If you have a large Initial Success Rate, you have more interests to use optimized IM's."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Partial reset mode (failed cards are partially reset), need rework))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## diff 90%"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYsAAAEKCAYAAADjDHn2AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvOIA7rQAAIABJREFUeJztnX2wJUd5n593VytgJbDRrsrG4L0LkSvEISBLaxwcRagMIpikDAoEimxAAgpZWyYmTgiKrTJyxVECmLiIy8bK2pEFddfCBoz5cBUBEz7sRHFYfSAJJEMZtIpAQULrBGHJMaDOH2dOdvZoZvpjuufjnN9T1XXvnXvOTE9Pd/+63/ftHnPOIYQQQnSxY+wMCCGEmD4SCyGEEF4kFkIIIbxILIQQQniRWAghhPAisRBCCOFFYiGEEMKLxEIIIYQXiYUQQggvp4ydgVzs3bvX7d+/f+xsCCHErLjhhhu+7pw70/e5tRGL/fv3c/To0bGzIYQQs8LMjoV8TmYoIYQQXiQWQgghvEgshBBCeJFYCCGE8FJMLMzsGjO718xu83zuh83s22b2ktqxi83si1W6uFQehRBChFFyZnEt8PyuD5jZTuAtwEdrx84ArgR+BHgmcKWZPb5cNoUQQvgoJhbOuU8Dxz0f+6fA+4B7a8f+HvAx59xx59yfAx/DIzpCCCHKMprPwsyeCFwE/PrKv54I/M/a33dXx4QQQozEmA7utwOXO+ceTj2BmV1qZkfN7Oh9992XMWtCCCHqjLmC+wDwbjMD2Au8wMy+DXwFuKD2uScBn2w6gXPuMHAY4MCBA65gXoUQYqMZTSycc09e/m5m1wIfds79fuXg/rc1p/bzgJ8dIYtCCCEqiomFmV3HYoaw18zuZhHhtAvAOXd12/ecc8fN7BeBz1SH/rVzzucoF0IIUZBiYuGce3nEZy9Z+fsa4JrceRJCCJGGVnALIYTwIrEQQgjhRWIhhBDCi8RCCCGEF4mFEEIILxILIYQQXiQWQgghvEgshBBCeJFYCCGE8CKxEEII4UViIYQQwovEQgghhBeJhRBCCC8SCyGEEF4kFkIIIbxILIQQQniRWAghhPAisRBCCOFFYiGEEMKLxEIIIYQXiYUQQggvEgshhBBeJBYDc+QI7N8PO3Ysfh45MnaOhBDCzyljZ2CTOHIELr0UHnxw8fexY4u/AQ4eHC9fQgjhQzOLAbniihNCseTBBxfHhRBiykgsemBmUZ+/666440IIMRUkFgOyb1/ccSGEmAoSiwG56irYvfvkY7t3L44LIcSUkVgMyMGDcPgwbG2B2eLn4cNybgshpo+ioQbm4EGJgxBifmhmIYQQwovEQgghhBeJBVpVLURO1J7Wk433WWhVtRD5UHtaX4rNLMzsGjO718xua/n/C83sFjO72cyOmtl5tf+91cw+Z2a3m9mvWOzqtwg2dVW1Rn+iBJvanjaBkmaoa4Hnd/z/48AznHNnA68GfhPAzH4U+DvA04GnAT8MPLtUJjdxVfVy9HfsGDh3YvQnwRB92cT2FMOcB2nFxMI592ngeMf/v+mcc9WfpwHL3x3waOBU4FHALuBrpfIZs6p69UHPFY3+RCm0S0E7cx+kjergNrOLzOwO4A9YzC5wzl0PfAK4p0r/2Tl3e6k8hK6qbnrQy+NzQ6M/UQrtUtDO3Adpo4qFc+79zrmnAi8CfhHAzM4C/gbwJOCJwI+Z2d9t+r6ZXVr5O47ed999SXkIXVXd9KCXx+eGRn+iFKV2KZiz+WbJ7AdpzrliCdgP3Bb42S8Be4F/Cfx87fibgDf6vn/uuee6kpg5t5hT1BPOrOhli7C97dzu3Sffy+7di+NCjM32tnNbW4s2t7Xl3KFD61Fft7aa+pDF8TEBjrqAPnq0mYWZnbWMcjKzc1j4J+4H7gKebWanmNkuFs7tYmaoUNZpNK49qsRUaTL3Xn31vM03S+Zuoiu2zsLMrgMuAPaa2d3AlSyc1TjnrgZeDLzSzL4FPAS8zDnnzOy9wI8Bt7Jwdn/EOfehUvkM5aqrTo4frx+fI9qjSkyRJnPv/w+DWWE25puKZXu74opF3vftW/Qfc2mH5tqexMw4cOCAO3r0aNFrHDly8oM+dsxYl/ITYgrs2NEuDqtsbcGddxbNzkZgZjc45w74PqftPiI4eHBROR9+WJVUiBK0mXVXl+XOyXyzLkgsKtYh2kKIudNm17/sMvnYxkZiwfwXy4yNhFbkoi344h3vOHlWP0WhWPt2EBIyNYfUJ3Q2NaRtUXybjcJw+7EaJuort9jPi2GYczsgMHR29E4+V+ojFs1rKJx3DYXEIk5o1dGdTGwHM+cOad2Z6hqKEELFQmYo1msNxdC0hS8eO3byNFymvkfi2/5h1azx+tevx3qDdWT2q7MDkFgw/8UyY9IlqHUxmPu+OCXo6mCaxPX+++POI4ZjEwacEgu0orkPTUK7pC4GmzDyiqWrg2nbiyzmPCn0ddKmfn/uzuGNGHCG2KrmkEruDdVma0c+C+fcojya7LV1v8+cbbql6PJBtPnRVlNOn0Vfn0jq99fFFzNXnxxycOehqyJLLE7gE4N16RBys73t3J49J8pkz54TnU6b+K5+Nhd9BT30+6udav3+NZAYnlCxkBnKg2ztYfim4XM19Q1hHnnooRO/33//wlfxghc0m/eca/5eDvqaCkO+H+uLmbt5aq0IUZQ5pFIzi66wWjSzOIm5TsPbGGI21DUar5fnzp3lR99DzCzaPtOU9uzRbHQIkBkqnqbOrqsBSCzWmyH8LKFrfFLXAsUwhM8ixhcj89QwhIqFzFAVbesAmswBaxflIBoZIoIrNORyiNDMvqbCkO+35XfPnkd+7/jx5s/2Kf9Ys5bMYDVCFGUOqe/MItQcoGiozWGImUXoaH5dAgRi7iN3+WvFfDPIDBVHyjRfYjFt+vpQ2jqLQ4fy+mZC87kuPqGY+83ZWceKz6aEe2cVC2ALeG71+2OAx4Z8b8hUcmbRRohYrEsDnxu5Opp1fR/0XMjZfmIHhEP4iaZANrEAXgt8Bviz6u8fAD4ecvIhU1+xSOlcfGKxKdPYKVJqVJjrvBpEdFOifDSzaCanWNwMnArcVDt2a8jJh0yloqG68InFXCrbOnZcpUaFMedtK9ehBxEln2+Jc5cqH/ksmskpFn9S/byp+nkKcEvIyYdMJbf7aMMnFlOaxk6l4+pDTMc09syiq1yHHET0eb5tq8tznLuLkuWj94c8kpxi8Vbg54A7gAuB9wNXhZx8yJRLLGIqx1xmFkN2XKUaV9M9NHVgXZ8fcnTaVa5DDiJSn+/2tnO7dj3ye6eeeuJeS9XvKQ2yNoGcYrGj8lu8B3gv8NqQEw+dcpmhYjqYUj6L3B3uUB1XyVlK18rftmuUFC7febvKdchBROrz7SrvZT5LdepTGWRtCjnF4vUhx8ZOOcQitpKWiIYqITBDdVwlG7lv5e/UOhLfup2hTH+pz6SrvEvvJDyGaXQTzE1t5BSLGxuO3RRy8iFTDrGIHSmVWGeR0gB96wG6OticDTOk/FIbZdd95BjN5sZXrkN1TqnPN2RmUbJTH7LznpPfrgS9xQJ4OfAh4M+BD9bSJ9YxdNa5MjOLWFKm9m357hod+jqulMbqK7++ztYmn8VUZxbOTWe0mpKPEJ9F6rlzkPO6m272yiEWW8AFwPXAs2vpHOCUkJMPmabos0ghpeKGbs62OqNoo48prOt7fRvlanTOMu3atTheXziX2pFMpYOfAr5oqDHzlXMmkOp7WZe6ou0+EskZDZV6/diG4DPRxJpr+nTqqb6TGOrX2LNnMdrtuue2vZaaZlNTNkesS+fUl9wzgZym39RZeds1hnjeOX0Wf7tawf1N4K+A7wDfCDn5kGmK6yxSyeEUb+uYQxrUnKJcQoWyfo22hh6yJfaYZpcpC9mQ5K6fOQdoud7BMeTzzikWR4GzgJuAncCrgH8XcvIh0zqJRQqrnVifPYzmFOUSaoKrdyQxM7H6d8fssDfdrl6nRFnEDgJSTL8xDPm8s4pF9fOW2rG1jIbqoqkyTUksmkgdBZecYg+1hqSrkaU29DE77K48b5pJKka0S80EUwccoQy5MDGnWHy62hvqXdVq7p8BPhty8iFTSbFoq5xTF4s+TN2mv8zfsgF1NdTVPKaaEMZcWezrnMY0SYV2yDk77pBz9Y2+6zp/jOl3Wbdi7n2uM4st4NHA44ArgV8Gzgo5+ZCppFi0N9T1FYsmpmIK6WqoIdFQqTOn1PvP0Un6QofHeA5t+VrtkNui2EoLXJ/nFfpCqpABy65djwzC8N377HwWlY/iSMiJxk4lxaK9IkxXLEpMv2NH1kObAGI6yyYfT4mRas5GX++cxprhrOJ7Fj6RKylwubc62bkzbhCxvL/Ud4mHzG5ytK+cM4s/Bk4NOdmYKVUsQgp8iJlF7in62LuBlhwZjRkNE/ucho4AG9p/4XsWIbb9UnlOLXtfflfrh68MSpgvc7avnGLxrip09ueBf75MAd+7BrgXuK3l/y8EbmHxvoyjwHm1/+0DPgrcDnwe2O+7XopYxEw3S/oscnesQ0YzrS6IW+a5a3TW12YcK1qpg4Eco94u80QqvtH6kP4LX9mFBhOUyHNqu9q5MyzPoYEPIfVrzEFITrG4sikFfO/8arV3m1icDlj1+9OBO2r/+yRwYe1zu33XSxGLvp1OLrHI3VmVdMTWy+G009obfehWI23X8O2r1EfkVz9XYrHgso60PVuz9M6xzQeQW+hC85Kyar8pxTqBQ/MXe87Q/IaGVOeqz3VytvFJrOAG9reJxcrnngXcXv3+g8Afx14rRSz6FnguscjduQ/hiO4ShK0tfyfRlZdcI7HQcsjlA2nbzLHPAsnV+13Wi9CObAi6nkWIY74tjRXhFSpwMfUxd+DEJGcWfZJPLICLWLxU6TjwrOrYi4APA79XLQT8JWCn71qlZxZNTHVm4YsWytEAfY7WkE6izbmcqwMMFeEcZsCuZ9hVBkuTR9dzSelwQ2bHffbQimH1ur5ZUY420De/vvLOKWQpg8VJ+iz6pIiZxfnAH1a/vwT4P8BTWLzC9X3Aa1q+d2nl7zi6b9++6ELq69icqs+int+6UOSs6L5FYss8hNp/S3QcuX0bKeXhE8CQ5xJjymk7z9AdYBcx4jdGhNcyj6HC2rfupA4WJxcN1SeFikX12S8Be1nsRfWp2vFXAL/m+37JaKjSDu7QfKRQyiQVaovvY4aop507mx3pXZSMyAotj2VeQ8ug7uuIEZq2sE5f/nLXi1BCzWpd+UltM1OLPByynjaR08F9Jot3cB+uIpyuAa4JOnmHWLDYb2rp4D4H+ApgLNZ2fBY4s/rfbwE/5bvWWIvyYjb7W1bQPXvSOr6UCh4zxY25Rpup69Ch7vP2FY3YxlRKhJuu43NihnbYTe+R6FsWKXtolWZ7u/teu+4rtYNt+l5bRF8IuQZjQ9XTJnKKxX8D3gK8FHjxMgV87zrgHuBbwN3Aa4DLgMuq/18OfI5F6Oz1nBw6eyGLsNpbgWsJWOcx1qK8kArnG1n2iQ7yVbLQypzS+FIqeKxJpSuNYc/uoo/TPTbF+qCmNrPw5aktxNonur78h5RDTL2fguj2JadY3BxyorHTWDOLkAoXUkFTIh9CtkMOFYFS5qpVcpml5tYgl/S5/z4jzyn5LJbEdrahQRNN34ud2bbV+9DnN7WBTBc5xeLfAC8IOdmYaeiNBH1iUa8wfUMdY803TbMGX0czZIfc5TxseplRjtDTEvS1mUNcAEDu/A4VDdVG7AwhZdCVKs5t9T7HzGRq5BSLB4CHgb+sfn+ADXz5Ub2hnWjg3WIRu+VBG7GmC18H39RpTKlDbspfCQdgaue5vZ13Y7yQwcCePen3mUIJG3pTeYe853uJr5xiZsyp7dE3GxpDdPsyiWioIdOQLz86MVoJm1mU8lmkbFDW5phuq/xjVfyUjryvg74pnX76yUEJXeVVL/8cTlJf55l67z5KROe0nfPQofD3fPtmIrEC0zaT7bpXX9TbmDO1VLKKBfATwNuq9A9CvjN0GvpNedvbzi2joUIqXIloqJRGHTvSGoNUZ3vMd3I62lMGAL68xziwQ+49piNL6ZR95PCJtUVPdYlpyHX7DjKWojdEgMjq91L6klVymqHeDHwceHWVPsYavlY15cHV11mUsgf78hWb7xxOvtKkdCyx38kZxuvrWENNW6n1xXfvsUKaYu7xkWtLm9jZdKlZ0uqzCnkGOUyrfa0UTeQUi1uAHbW/d1J7xepUUh+xSK1QbYvyclXQEhW9rVKvNua2NROxtHWCXZ1jW0NICQLo46jMkYZYbOW797Z7bQtP7etjayJXtF3q1hg5TZix+Yox/frKo8SzyS0WZ9T+PmPdxCK1IreJRa6Gkes8ddoE6DnPSdsWpKuRpUzZt7fTGlJsWaVGycSktkin3DM23713zRRStwaJnRE0BVGkCOeQbSLX7C6n6bdvZGUTOcXi5cCxanHcO4EvAy8LOfmQqY9YpE6RgcbOMvR8vtFMrqn7KinT6LbzdHX6bZ1lVyfaNfPJ6bNoKocmwfSlPXvaxS+1Qwh5Xqv/X7Xl79p14nO+zmq5NXj92dQd+m3PsK9PJmXmOuRsu+89LvMVU6d89XzSM4vFuXhC5eT+CeB7Q74zdBprZtHUCE4/3X++kAbUlq8SYZQpwtSVv9TY9q6G5SNHNEqb83DPnu73dzRdu0ssY/IT4rxeDbCoO31TZ1Fts8DYTjr3bCB31FHKzCvEtLr0W8aWe+wALeWZ1OktFsBTq5/nNKWQkw+ZxvJZtD00X3RUyAg6JfojlZQGneokTplZjOVsrxOat66ZWqjwxVwzNOInZQdgn/CFPJdSM+Rc+J6Vb5CXYsrrM/OcXDQUcLj6+YmG9F9CTj5kGisaqu1h+9761dXR1itn6sveY+8rRTBTRk0hPovUKJHQeP1UQjq9kFF8jPCFXDPG7Bk7w1ieo0+HP+UBgHP+cqnfY+i9dLWNZd1MNbnmJqfP4tEhx8ZOQ6+zcM51ioWvEXVVpvp3U0cfzsV3vG3T6C5becxiwboN2OcYX/2f7/O+GVgO00WfUX6M8OW6ZlNnvFoOoa9n7dPh5/QzlIpsCp09hYpmiGkrNZgjNznF4saQY2OnscTC97C77JshFaWP3XuIBt7WsTeN1FJH+768dHXQyzz1Ec0u2//qeXwzxtj7D/VZpHbGXaPq+jn6Oqlz+ZOGLotQ83HMzKL++bb/D2miy2GG+l7gXOB24Idq/ooLgDtCTj5kGkssukICmypffSvzJofp6lbnXRXOx5img1WzUJ9O2meK6+qgu8qxbeTd1WH4Or0SJpdSo+nV70L3a15zhb/GUL+vkJF/jnqbMpuO9Vn41sDMamYBXFz5Jx5Y8Vd8EPiHIScfMo0lFs61V7AQm35dHGJ2XA0J6+uqiL5GEdIB97n+Km2N8NCh9Aa3vFbMaD9HZzPEQrwxGLpjC/WxpPhv+uYr1A/oE7gp1JecZijvi46mkMYUizZCo4V8duG284SsI0jdxyZE6HzX9zXcttnHqij6ys3ns4jxI/QJ3a2Xe1+TyxQZOqopNIAiZGYRs2YihJyBI2PXl5xi8SjgH7N4teqblink5EOmUmLR9SCBJJNEW2NL6TR9o7qm/PtGiCGdeMj1fTOb2NeHrqbVBtcWDdXlH1q9j6lEqAxJaGc19MwiZLDV1PG21aumz6ZE0OXygU2FnGLxEeB3gDcC/2KZQk4+ZCohFl2VYnvbudVoqBBnW1dj63Jm5xzVdZ0rNM8h1+8qv1AhbSuT2IWJofdROkIlxCY+ZMcS0/ENbTLpmiV0ld+qKbfpGfpmoyn5Cg0caXu+MbOVnHUkp1jcFnKisVMJseiqFIv/PTJ0drXC+DrG+ojVJyhtx2Mrj/++mhtoSgNpy1voqDHXi49Co1O6nkOIMPs6A5/zfGj7dWyHPKSYpZRHyCDE1yZ9dTp14OYbfKY6zmOi0prIKRaHgb8VcrIxUwmx6KoUi9+b11k0EfKQfYKQ6n8IyUvKPjZ9OrKQzrveiS47qNR3gXSZ1UL9Nb6Flr4G7xuRjhEZk2LqGZJYcQq9n67/p66R8j2nlEHa6jnb6nEfM2lOsfg88FfAn1Y70N66brvOttH1ABcj7UeKRdf6hxATRJOgLK9Z32tmOdLPPeIPNQ/1XSGdYgZIGWm2mdV27DhRTiHPYdeu9C1cls/CNyLNaWrs64dou4c2O/8QM44QE07IvXS1G1/baRt4hAhq1/P1RewtB0ohzyiWnGKx1ZRCTj5kGtpnsfi7e2YR24BCRsBdG7v17WDa7rmtsfXtEGIdjCkjutRR4OqzC9l2xdfZDzWz6OuHaLuHNoHfuTPu1aQpxJpwQlKOwUpIvV2SMrOImemntvncu86eB7yq+v1M4Mkh3xsyDR0N5fNZpERMhHbSIRWnj+mifs9d1xjaPJEy8s41Wg85T0iU2RA+i1jRqT/vlE0eS9TBmHvqylfX/QwxWKkTK3gxQtGnvHPOLK4EPgR8ofr7+4D/GnLyIdNY7+Bua9ixFSu2IXal1VXgfTr0GMdwV1nlMFEMObNIOU9IZx9iiuxbVn0EMpc/K/R6Oe6pK1+5AiR8eVil7TnGREPFlHWfgVtOsbgZMOCm2rGN8Fn46FpnEdtgQxuiL7y2aRV4n4rkm/H0iQDJkZc+zvwS1x4yWqiNHCPgPv6sVEFOvSdfZNPQg5XSM8TV1Nd/mFMs/kf188bq52kSiwVdK7hLzCxCop9KRNRsb6c70nPnJ7Th942g6nPtXN9NvV5OgV49b5vPYrV+1N/UlwOfCSdkMNNXNELLNafvySfGOco4p1i8AfiPwJeA1wLXAz8dcvIh09TEIofPos2c1FXpS23JkNoBDb1FRJ+8Tik/fe+hlLA12fkPHep+U1+uPDZ9fnmsq0NdfjfXDNOX55x1vu/7bELI7eC+EPgl4G3AhSHfGTrlFouQSuHbGypHY4ilxMyiT/5K5ifHNYcwG43pa4klpVNNyWuOzttnIi09424j57WGGPjknFm8JeTY2CmnWIQ+IJ9YjEFbVEWfFZ6581N6hBw6shtqBjJmFFcsKR1dSl5zdKghvoo+eUwld70qPaDJKRZNLz9aa59FaEWeolg4N857B7oY2vYe+vyGGm3OaWaRs+Pv2r8rR+cdc46hy3MKgQ6h5HifxaFqtfZfVCu3l+nLwHbIyYdMOcXCVwlP2EkZtSLERq2U7mhyk3ofoSIz1GhzzCiuWFJNSrEL3Nps8TEbRMaaG6fkx5oSOcTiu4D9wHUrq7fPCDnx0GmomcXJlY7RKl1X5R/LhJGbvusFfCO7nFErvmuViIYqcd3UTjXWEZtDLFKCSOYy2h+Skiu49677Cu6uSnhyB0NnoyhZObs6uk2fWYRSysk6xOAh5Lrb28NEKTmXb11R7IBGAtAfreDuSVslPLmS01rJS3ciXY1tXabcQ9xH386mlKD58hVy3Ryj91C6nM1N7/VelwHNOjD6Cm7gGuBeWt6HAbyw8oHcDBwFzlv5/+OAu4FfDbmRodZZhM4sSjcG3/nXZcQ19fsoYfILEcmQ67Z13nDytXKUry+MdfU+5jCgGaLuTaF+j76CGzgfOKdDLE4HrPr96cAdK///D8BvT00sQn0Wpf0Gc2hsm0CJQUHIOUM+4xOLUiGePsGYw4BmqFntFNrwJFZwVw5y75v2gGcBt9f+Phd4N3DJ1MTCOX801PZ2+vYYodfOuY1F7HW3tk68W2OKDX1ISjT4kIFG23Xrz6XtPEszVKnZr2+fszkEWgxhJpuKKW4SK7h9YgFcBNwBHAeeVR3bAXwSeNJUxWIJ9fm8WzTgmDeytdE24pqSMzX13oZkqJFr7uuEdiJNAu57TvU9m0rNfn2zizn4JYaIKJxK1GJWsTjpC4vO/GDgZ0NnFucDf1j9/jrgjdXvnWIBXFr5O47u27evRDl2UhcLX4ca+rKg8EisE8n3us++hJgWxu4EQjrOKQpaE6mDgrbn1PY+7ZLO+ba2MJdnoJlFhFhUDuafBX4VeF7l5H4dcCfwgaCTB4pF9dkvsQjLPQLcVV3n68A3gDf7vj/2zMLXoYaOFroqUOg25svP5RKOmOuGkHsk3tQ5teW5REMcaw3FKrEj1T4z1dD8rz6LvttpD8VYPouc76IJJYdYfAC4FvhJ4Hcr09CngLNDTuw8YgGcVXNwnwN8Zfl37TOzMUP5OtTQTqqrwYeO8PtU8KZOIOfMokQjjCmXnFtWp95PqY4oZYfSVKELzf9UnLglBD3Hd1Z9kE3vointH8whFrfWft/JIgz20SEnrb5zHXAP8C0WIbCvAS4DLqv+fznwORahs9ezEjrrZiYWXR1WTOPomlmEhiemClWX0zSXz6LE1Dt05lO/Tq5OLOV+SpTB9nb8lhupxOR/CqaWoQTLd53UtTN9+pMQcojFjV1/Ty2NLRZtHXnstDu2wnU51Jcphwmsy3kfc48lnHpt+e7aTLHL/xPDVHaUzXU/IcTkfwpO3KEEK3agF7p2pnTec4jFdyp/wTeAB4Bv137/RsjJh0xji4Vz+WzxMedpG1GmVKzQkM0+91hqVO0LI13Na1fDjLmnqcwshuyU5zazGKpsUkzI9XKINTPnml0Ui4aaapqCWIxFjnBd54Zp2KVMArEi1tUwY+53Kj6L2EioPszNZzGFmUXq2pmulKscJRYDMBWxyDVKHtK2O8T6B18e2spsiM3shogIK9W5LK8XM/sd83lPwWcRKljb2+0LeksJnsRiAPqIRc4GlHPkNHbDHpIh3m88JPVnV2oHgbkyVL1uu07sbMxnWk4d2DQhsRiAVLHIPdKZwlR/jqxzuU3BsSxOJnY2Vh/M7NhRTvwlFgOQKhalnLybMiPIybqW2xQcyyIfJQc2oWKxXBQ3ew4cOOCOHj066DXNjJTy27Fj8bgfeT54+OEMGRMbz5EjcOml8OCDJ47t3g2HD8PBg+PlS6Rz5AhccQXcdRfs2wdXXZXnWZrZDc65A77P7eh/KRHLvn1xx0U7R47A/v0LAd6/f/G3WHQihw/D1tZiELK1NT+h0LM9mYMH4c7HzB3VAAAJ20lEQVQ7FwPKO+8c/lmeMuzlBCxGBE2jvquuGi9Pc2R19Hzs2OJvmFenWIqDB+dbDnq200MzixFYh1HfFLjiipMFFxZ/X3HFOPkR+Uh9tpqNlENiMRJjTynXgbvuijsuTjD1TjXl2S5nI8eOLXyCy9lIqXubehnmRmLRwKZUgineZ0ye5PtJY+hONYWUZ5trphlSB+dQhtkJCZmaQ8oVOhsTosZEVnCnUGKtx9Dbfa/zOomSzCGsNuXZ5lhbEnrdOZRhKGidRRoxlWDOYpF71fdY232v6zqJksxlwV6u/b5i6nToOeZShiFILBKJqQRzFoume1ymWHIJzzo1wCkzh1Fx6n5bfQctoXVwDmUYSqhYyGexwqbYwXfujDveRS5H86aU/dhcddUiVLvOlEK3U/0BOaIMQ+vg1MuwCCGKMockn0UcU5xZyAcxHKt7D03p3dhjjtpjN/xbBxMoMkOlE1oJ5iwWU/RZLM+1Dg1w6kxZmMc2R25aHZRYDMAcxaK+t37Xa0dTz7spDWyqhD6HKdvcp5a3da/bEosBmJtYNI0ml4Kxjo1g04iZLYw9eu9iSrOeKeWlFBKLAZibWExtxCbyMrd3Y3cxldF8qXKayv05Fy4W2qK8B6lblI+FtkZfb2Ker7YwD6NEm5la2WuLcvEIFJq63sQ8X21mucC3tUeJNjPXDTAlFhvERsaGF2ZK+2vFPt9N38wyZD1HiTYz2w0wQ2xVc0jyWYQxJVvp3Jmi81PPN5xQf0TuMp2avwj5LMozN5+FyMv+/YvR6CpbW4uRupg2Y/nw5LMQYsOYrTlBAOP58ObqL5JYCJGIAgbmzZg+vDn6iyQWEaw6M8Vmo4CBeTPXEf5YSCwCaYqcWB4Xm4k6m2mREpk2xxH+WMjBHUizM9PY2nJyZgoxMlNzGs8JObgzI2emECeY0voSmO9CtzkhsQhEzkwhFqS+nKgkGsyVR2IRSJMzc3lciE1iiqN4DebKU0wszOwaM7vXzG5r+f8LzewWM7vZzI6a2XnV8bPN7Hoz+1z1/5eVymMMTc7M5XEhNokpjuIVmVaekjOLa4Hnd/z/48AznHNnA68GfrM6/iDwSufc36y+/3Yz++6C+QxmNXJiXYIDhIhhiqN4RaaVp5hYOOc+DRzv+P833Yne9jTAVce/4Jz7YvX7V4F7gTNL5VMIEcdUR/EKgy3LqD4LM7vIzO4A/oDF7GL1/88ETgX+bOi8CSGa0Sh+Mym6zsLM9gMfds49zfO584E3OeeeWzv2BOCTwMXOuf/e8r1LgUsB9u3bd+6xpl3dhBBCtDKrdRaVyeopZrYXwMwex2K2cUWbUFTfO+ycO+CcO3DmmbJUCSFEKUYTCzM7y8ys+v0c4FHA/WZ2KvB+4F3OufeOlT8hhBAnOKXUic3sOuACYK+Z3Q1cCewCcM5dDbwYeKWZfQt4CHiZc86Z2UuB84E9ZnZJdbpLnHM3l8qrEEKIbrQ3lBBCbDCz8lkIIYSYNhILIYQQXiQWQgghvEgshBBCeJFYCCGE8CKxEEII4UViIYQQwovEQgghhBeJhRBCCC8SCyGEEF4kFkIIIbxILIQQQniRWAghhPAisRBCCOFFYiGEEMKLxEIIIYQXiYUQQggvEgshhBBeJBZCCCG8SCyEEEJ4kVgIIYTwIrEQQgjhRWIhhBDCi8RCCCGEF4mFEEIILxILIYQQXiQWQgghvEgshBBCeJFYCCGE8CKxEEII4UViIYQQwos558bOQxbM7D7g2Nj5GJi9wNfHzsTIqAxUBqAygPQy2HLOnen70NqIxSZiZkedcwfGzseYqAxUBqAygPJlIDOUEEIILxILIYQQXiQW8+bw2BmYACoDlQGoDKBwGchnIYQQwotmFkIIIbxILGaKmd1pZrea2c1mdnTs/AyBmV1jZvea2W21Y2eY2cfM7IvVz8ePmcfStJTBL5jZV6q6cLOZvWDMPJbGzL7fzD5hZp83s8+Z2eur4xtTFzrKoFhdkBlqppjZncAB59zGxJab2fnAN4F3OeeeVh17K3DcOfdmM/tXwOOdc5ePmc+StJTBLwDfdM69bcy8DYWZPQF4gnPuRjN7LHAD8CLgEjakLnSUwUspVBc0sxCzwTn3aeD4yuEXAu+sfn8niwaztrSUwUbhnLvHOXdj9fsDwO3AE9mgutBRBsWQWMwXB3zUzG4ws0vHzsyIfI9z7p7q9/8FfM+YmRmR15nZLZWZam3NL6uY2X7gh4A/YUPrwkoZQKG6ILGYL+c5584Bfhz4qco8sdG4hU11E+2qvw78NeBs4B7g34+bnWEws9OB9wH/zDn3jfr/NqUuNJRBsbogsZgpzrmvVD/vBd4PPHPcHI3G1yr77dKOe+/I+Rkc59zXnHPfcc49DPwGG1AXzGwXi07yiHPu96rDG1UXmsqgZF2QWMwQMzutcmphZqcBzwNu6/7W2vJB4OLq94uBD4yYl1FYdpAVF7HmdcHMDPhPwO3OuV+u/Wtj6kJbGZSsC4qGmiFm9hQWswmAU4Dfds5dNWKWBsHMrgMuYLG75teAK4HfB34X2Mdi1+GXOufW1gHcUgYXsDA7OOBO4Cdrtvu1w8zOA/4IuBV4uDr8cyxs9htRFzrK4OUUqgsSCyGEEF5khhJCCOFFYiGEEMKLxEIIIYQXiYUQQggvEgshhBBeJBZCeDCz79R28by52l4h9hyXmNn35c+dEMNwytgZEGIGPOScO7vnOS5hsUDqq6FfMLOdzrnv9LyuEFnQzEKIBMxsv5n9kZndWKUfrf3v8updI581szeb2UuAA8CRambyGDN7jpndVH3uGjN7VPXdO83sLWZ2I/CPzOynq3cW3GJm7x7pdoXQzEKIAB5jZjdXv3/ZOXcRi32HLnTO/aWZ/QBwHXDAzH6cxVbZP+Kce9DMznDOHTez1wFvcM4dNbNHA9cCz3HOfcHM3gUcAt5eXeP+apNIzOyrwJOdc//XzL57sDsWYgXNLITw85Bz7uwqXVQd2wX8hpndCrwH+MHq+HOB33LOPQjQst3EX2chOl+o/n4nUN81+Hdqv9/CYkbyT4Bv57kdIeKRWAiRxs+w2JvpGSxMTKdmPPdf1H7/+8CvAecAnzEzWQPEKEgshEjju4B7qq2gXwHsrI5/DHiVme2GxXuhq+MPAI+tfv9TYL+ZnVX9/QrgU6sXMLMdwPc75z4BXF5d8/QC9yKEF4mFEGm8A7jYzD4LPJVqNuCc+wiLrbKPVn6ON1Sfvxa4ujpmwKuA91RmrIeBqxuusRPYrj5zE/Arzrn/Xe6WhGhHu84KIYTwopmFEEIILxILIYQQXiQWQgghvEgshBBCeJFYCCGE8CKxEEII4UViIYQQwovEQgghhJf/B3lPO+GavGteAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYUAAAEKCAYAAAD9xUlFAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvOIA7rQAAIABJREFUeJzt3X2QJGd9H/Dvd1d7tvZOIGl1IbLQ7uIEQzCFhFgFIhFK5OwYyTFIhSDlbIR0UWqtK+FIiXGUcFUxLtelIGWrgAKknNFFkm8tbErCJgRjCBFIBiNn7zj0doUsh7uzxJlbncA6cSRId7/80T1zvb3dTz/98kz3zHw/VVO7M9PT/fQzPf3r57VpZhAREQGAibYTICIi3aGgICIifQoKIiLSp6AgIiJ9CgoiItKnoCAiIn0KCiIi0qegICIifQoKIiLSd1rbCSjrnHPOsfn5+baTISIyVPbs2fOMmW0uWm7ogsL8/DxWVlbaToaIyFAhedBnOVUfiYhIn4KCiIj0KSiIiEifgoKIiPQpKIiISJ+CgoiI9CkoiIhI33gHheVlYH4emJiI/i4vt50iEZFWBQsKJM8neT/Jx0k+RvImx7IXk3yR5NWh0rPO8jKwtAQcPAiYRX+XlhQYRGSshSwpvAjg18zsNQDeBOBGkq9JL0RyEsCHAHwxYFrW274dOH587WvHj0evi4iMqWBBwcwOm9ne+P9jAPYDOC9j0V8FcC+AI6HSkunQIe/XSQZOjIhINwykTYHkPIDXA3go9fp5AK4CcNsg0rHG7Gy510VExkDwoEByE6KSwM1m9lzq7Q8DuMXMThasY4nkCsmV1dXVZhK2YwcwPb32tenp6HURkTEVNCiQnEIUEJbN7L6MRRYAfIrkAQBXA/gEySvTC5nZTjNbMLOFzZsLZ371s7gI7NwJzM0BZPR3587odRGRMRVs6mxGFfF3ANhvZrdmLWNmr0gsfyeAz5nZH4VK0zqLiwoCIiIJIe+ncCmAawA8QnJf/Nr7AcwCgJndHnDbIiJSQbCgYGZ/BsC7246ZXRcqLQCi8Qfbt0fjESYngRMnoiqjHTtUWhARiQ3dndcq6Q1U641LOHEi+tsbsAYoMIiIYFymucgaqNajAWsiIn3jERTyBqrlvZ+eE0lEZEyMR1AoGpCWfD9rTqTe6yIiI248gkLWQLUeErjiilPPb7opu6rp2msVGERk5I1HQ3OvEbnX+yjJDLjrLuDSS6PnR49mr+PECTVKi8jIo5m1nYZSFhYWbGVlpfoK5ufXBwYg6p76/POZQYEALLncgQPVty8i0gKSe8xsoWi58SgpJOU1OmcFijKfFxEZAePRppCU1+g8OVnv8yIiI2D8gkJWozN5akCbi2ZRFZERN35BAQBOP33tc1e7ykScRZpFVUTGwHgFhd4YhLweRmnT08Ddd0f/HziggCAiI2+8goJruou0yUmVDERk7IxXUPDtYQT4tTGIiIyY8emSurwcNSiXGZfRG6wmIjImxqeksH17uYAAaAZVERk74xMUXIPOdu+u9jkRkREzPkEhb9DZ3FzUmDw3V+5zIiIjaHyCQtagteRgtKL3RUTGQLCgQPJ8kveTfJzkYyRvylhmkeTDJB8h+XWSF4RKDxYXoy6mc3NRg3N6MFrR+yIiYyDYLKkkzwVwrpntJXkGgD0ArjSzxxPLXAJgv5l9n+TlAD5gZm90rbf2LKkVkMSwzSYrIpLkO0tqsJKCmR02s73x/8cA7AdwXmqZr5vZ9+On3wDw8lDpWSN9u03dPEdEBMCAximQnAfwegAPORa7HsCfBE9Mb6qL3sjmgwd18xwRkVjwhmaSmwDcC+BmM3suZ5m3IgoKt+S8v0RyheTK6upqvQRlTXWh8QgiIgACBwWSU4gCwrKZ3ZezzOsAfBLAO8wsc6Y6M9tpZgtmtrB58+Z6icobd9B7XVVLIjLGQvY+IoA7EDUk35qzzCyA+wBcY2ZPhErLGnnjDmZnT1UtHTwYjX5OVi2JiIyBkCWFSwFcA+CfkNwXP64geQPJG+Jl/hOAGQCfiN8P363INR4hr2pJRGRMBOuSGkojXVKXl6MAcOhQVELYsSNqZJ6YyJwfiYC6pIrIUPPtkjo+s6QmLS5m9zSanS03vbaIyIgZn2kufORVLYmIjAkFhaS8qS5ERMbEeLYplKRpLkRk2LU+zYWIiAwfBQUREelTUBARkT4FBRER6VNQEBGRPgUFERHpU1AQEZE+BQUREelTUBARkT4FBRER6VNQEBGRPgUFERHpU1AQEZE+BQUREelTUBARkT4FBRER6QsWFEieT/J+ko+TfIzkTRnLkORHST5J8mGSF4VKj4iIFDst4LpfBPBrZraX5BkA9pD8kpk9nljmcgCvjB9vBHBb/FdERFoQrKRgZofNbG/8/zEA+wGcl1rsHQDutsg3AJxJ8txQaRIREbeBtCmQnAfwegAPpd46D8BfJ54/hfWBAySXSK6QXFldXQ2VTBGRsRc8KJDcBOBeADeb2XNV1mFmO81swcwWNm/e3GwCRUSkL2hQIDmFKCAsm9l9GYs8DeD8xPOXx6+JiEgLQvY+IoA7AOw3s1tzFvssgPfEvZDeBOBvzexwqDSJiIhbyN5HlwK4BsAjJPfFr70fwCwAmNntAD4P4AoATwI4DmBrwPSIiEiBYEHBzP4MAAuWMQA3hkqDiIiUoxHNIiLSp6AgIiJ9CgoiItKnoCAiIn0KCiIi0qegICIifQoKIiLSp6AgIiJ9CgoiItKnoCAiIn0KCiIi0qegICIifQoKIiLS55wlleS/c73vuE+CiIgMoaKps8+I/74KwMWIbooDAL8E4C9CJUpERNrhDApm9psAQPIBABeZ2bH4+QcA/I/gqRMRkYHybVN4GYAfJ57/OH5NRERGiO+d1+4G8BckPxM/vxLAXWGSJCIibfEKCma2g+QXALw5fmmrmX0zXLJERKQN3l1SzWwPgHsAfAbAUZKzruVJ7iJ5hOSjOe+/lOR/J/ktko+R3Foq5SIi0jivoEDy7ST/EsB3AHw1/vsnBR+7E8DbHO/fCOBxM7sAwGUAfofkBp/0iIhIGL4lhd8C8CYAT5jZKwD8HIBvuD5gZg8AeNa1CIAzSBLApnjZFz3T0x3Ly8D8PDAxEf1dXm47RSIilfkGhRfM7CiACZITZnY/gIWa2/4YgH8A4LsAHgFwk5mdzFqQ5BLJFZIrq6urNTfboOVlYGkJOHgQMIv+Li0pMIjI0PINCj8guQnAgwCWSX4EwA9rbvsXAOwD8FMALgTwMZIvyVrQzHaa2YKZLWzevLnmZhu0fTtw/Pja144fj14XERlCvkHhHQB+BOBmAF8A8FeIRjXXsRXAfRZ5ElE7xatrrnNwlpejkkGWQ4cGmxYRkYZ4BQUz+yGAcxA1HB8F8Km4OqmOQwC2AADJlyGaSuP/1FznYPSqjfKcffbg0iIi0iDf3kfvRjTX0bsAvBvAQySvLvjMPQD+HMCrSD5F8nqSN5C8IV7ktwBcQvIRAF8GcIuZPVN1RwYqq9oo6dgxtSuIyFCimRUvRH4LwM+b2ZH4+WYA/zPuTjpQCwsLtrKyMtBtksSafJqYiBqWXWZmgGeGI8aJyOgjucfMCjsI+bYpTPQCQuxoic+OnlnnuL3I0aMqLYjI0PE9sX+B5J+SvI7kdYhmSP18uGR12PIy8PzzfsuqF5KIDBnfuY9+neQ7AVwav7TTzD7j+sxI6jUwu9oTktQLSUSGjO8sqTCzewHcGzAt3ZfXwDwxAZzMGHfnU80kItIhRbfjPIZoOgrGf/tvATAzyxxsNpLm5/PHJZw8CUxPrw0Y09PAjh0DSZqISFOcbQpmdoaZvSTx9yXJ54NKZKt6jcV5AQEA5uaAnTujv+Sp54uLg0mjiEhDfMcpXJ/x2gebT04HFTUW90oEi4vAgQNRqeHAAQUEERlKvr2P3kmyf5Yj+XEAfydMkjrG1VisEoGIjBjfhuZ3AvgsyZOIprr4gZn9q3DJ6pDZ2eyqo7m5qEQgIjJCnCUFkmeTPBvA6QD+NYB/D+AYgN+MXx99WY3FakQWkRFVVH20B8BK/LgfwJkAfjHx+ujrVQ2pEVlExoCz+sjMXkFyAsA/MrOvDShN3aSqIhEZA4UNzfHd0D42gLSIiEjLfHsffZnkO+P7KYuIyIjyDQq/AuDTAH5M8jmSx0g+FzBd3TM/H01nMT+/fvbT5WX3+yIiQ8J3QrwzQieks9Ijmg8ePHXXtcXF9ZPkpd8XERkiXjfZAQCSbwfwlvjpV8zsc8FS5TDwm+zMz4MHD2JdLvXGKeTNiaRxDCLSIb432fEqKcRTWlwMoFcvchPJS83sP9ZI43DIG9F88KD7DmyaNltEhpDviOYrAFwY90QCybsAfBPA6AeFvBHNgPuWnJo2W0SGUJlbap6Z+P+lRQuT3EXyCMlHHctcRnIfycdIfrVEWganyshljXgWkSHlGxT+M4C9JO+MSwl7ABSd9e5ENE9SJpJnAvgEgLeb2c8CeJdnWgYrPaLZRSOeRWTI+VYf/TMAuwB8H8ABALeY2d+4PmBmD5CcdyzyLwDcZ2aH4uWPeKalHb1GYzUsi8gI8y0p3BH/fTuAjwD4OMmbam77ZwCcRfIrJPeQfE/N9Q3Gjh1R9VCSqotEZET4jlO4n+QDiHogvRXADQB+FlGAqLPtNwDYgmgW1j8n+Q0zeyK9IMklAEsAMNt2A26vWmj79qiH0ezsqZvsiIgMOd87r30ZwNcA/HMA3wZwsZm9uua2nwLwp2b2QzN7BsADAC7IWtDMdprZgpktbN68ueZmK0qOWt6+PQoEde6yplHQItJBvtVHDwP4MYDXAngdgNeSPL3mtv8YwJtJnkZyGsAbAeyvuc5wtm6N2hLMor9bt1Y/kfdGQSfXt7SkwCAirfMe0QwAJM8AcB2A9wH4u2b2E45l7wFwGYBzAHwPwG8AmAIAM7s9XubXAWwFcBLAJ83sw0VpGPiIZgAk149oBoCZGeCZZ8qvUI3VIjJgviOavYICyfcC+MeI2gAOAHgQwINm9r9qprO0TgUFwD2ALU/eSGgyqpISEWlYo9NcAPhJALcC2GNmL9ZK2bAJUaWTN0q67UZ0ERl7Xm0KZvbbZvbQWAaE3oynWWZmqq1X3VpFpKPKTHMxfrZvPzUldtqGDcBHPlKtF9HiYjTqWfd9FpGOUVBwcc10umtX9NfVi8gVMBYXo0blOt1aRUQa5tumMJ7y6v7n5qKT+Pz8+pLE8eNRCQPQzXdEZOiopOBSVPefV5I4dCi76ikZMEREOkhBwaVX9w9k1/3n9RaanXUHDBGRjlJQKNILAFl1/66ShCtgiIh0lIJCHa5eROp2KiJDSEGhrrxeRKG7nWpCPREJQL2PQlpcDNPTqDeoTj2bRKRhKikMI/VsEpFAFBSGkXo2iUggCgrDSD2bRCQQBYWm+TQA120kVs8mEQlEDc1N8mkAbqKRWPeJFpFAVFJokk8DsGuZMiWIXlfY3/u96Pk116hrqojUVup2nF3Q2p3XfPLJ545qecsAURVQMmBMT7vHNqRLHT6fEZGx5HvnNZUUmuTTAJy3zORk+W6m6poqIg0LFhRI7iJ5hOSjBctdTPJFkleHSsvAZDUAk1G7Qa9qZ8eO6AY9SRs2ACdOZK/T1c00VNdUjZYWGVshSwp3AnibawGSkwA+BOCLAdMxOMmpLYAoIPSqinoNyl/72vrqI7P8W3u6upmG6Jraq5LKu3GQiIy0YEHBzB4A8GzBYr8K4F4AR0KlY+B6DcAzM+tP/sePR0HjhRfWvv7CC8CzGVlV1M00RNfUQVVJqTQi0kmttSmQPA/AVQBuaysNwSwvA0ePZr+XV02UDiAzM8UNxiEm3RvEaGmVRkQ6q82G5g8DuMXMThYtSHKJ5ArJldXV1QEkrSbXVfXkpN86Nm3yO7k3fa/nQYyWVgN581wlr3EplXV5P7uctjQzC/YAMA/g0Zz3vgPgQPx4HlEV0pVF63zDG95ggxZlUwmkWXQNvP6xbZvZ9HT++70HGWZniuzevT5909PR603Jy5+29rln926zubkoHXNzze5zSK7vbBDfZxd0eT87kjYAK+Zz3vZZqOrDFRRSy90J4GqfdQ5FUJibyz7pzcxE7ydPPpOT2cvOzTW8FyWEPjnm5U/b+9yBH24lrvzsYl6H0OX97EjaWg8KAO4BcBjACwCeAnA9gBsA3JCx7GgFhTInmGE+GZXVCza9UkFT+9xEEAvxwy1KV9V0pz/nKm22USpro8TV1dJnh9LWelAI9RiKoGBW7oeRXHZmJnoMWxVGkazg1/ux1NnPpoJq0z/conRVTbcrH7tQUmjrIqcjV+OZOpI2BYUGVQoKVYxyqSHUD6Op9TadvqL1Vd1e3ufySl6DPqbaOgF2+bfTkbQpKDQoeEmhpyNXFEGEKkI3td6mf7hF6aqablcnhrzjzXUsVjlOXZ9ps6qkrY4CPtttqqqwxj4pKDSokTYFn6qSjtQ9Zqp7cHa9pGDW7Ell0CWFqvtbNhAWfWaUL2yy+ORhnYDQ4IWKgkKDGut9VPTFFv2gBnEllLWNOgfntm35PayaKEJ3pGheOl1Ntikkq4rKHB9VTuA+x2hT30ed431QpYam8yNgz0QFhQY1Ok7B9cW23d88bxszM9UOzm3b3PvfVNq7Or5gUL2Pqh4fVUqmPp9p4vuoc7wP8kKhKD/KBN6sdJf9fhwUFBrUeEmh90jbvXvtCXhm5tSBXPXEnF6/68fqm27fgzPvSmdy0j/Nkq/uVWWIkkJT6mxnkFVYRdsqE3h9f38qKQxhUPCJ+OT6ekdXKaFoXT7dXouunsoEhOQBunv32mqiyUl3KaFsfib3oYslApdQaW7iqrLomEtfiGzcaLZhw/ptlK2+8lmuTvtayLa5dNqzZiggo9fN8k/0k5Pr99+nhkFtCkMaFMxOHTy+Ed91xeF7BeE6YHyunvKuNkn3CWhiIv9zeT+IKvnZxbYDl5BpbuqqMq8qamrKb/1AdAL03Vff5bpYUshL+5Yt5boE5/1uywSQChQUGlR7nELewZC8cvG5SvB55B1APldPPutuIo2bNpU/yKv+0NssXYSsxvA5XqoODixbjei6cEnvq281aJ2G+KxtNBGMXSdt1z75VvMFvvBRUGhQ7aCQdzD15kJyLVPn4dNVMHlAhkjDxo2nShKk2WmnVTvoXSfBPFlXvFNTfieWokASugokb/15J730CSq97akp/5HyZS9QfKfTcFWD9qpT01UzZQK666o82T5X5kKhzu8i63suyqeAFzEKCg2qHRSKiuN5dZNFB5zPcslJ+Fzrn54unwbfR1HxuE6VQLptJinv5JkMxlnflU+/c58ruirbd61/2zb3cVSUz0XpLcpr1zHm09DtWu/MTP2r5KJ05x3jZeYly3qUaeQfZCN4ioJCgxqZ5sJ1dZc8YH1/kGVO4OkrpLzlkiWGpqqzkuuue+XsmuMn+4vLf+Tx+dHWrSopCgplqyl675VpsHTlW1HHhvR2sxqfe8e0TxVmL0/KpDGLz343cQJP72PdQDOgtjEFhQZV6n2ULgKWGbtQtGyv3aAo0KSX90lLssjqqvYqqiNNP3rbr/PDd607K89d6cn7vnw+41rO9yToqiLxyc+s9fWUudLPO3Z9PpusGsw65tKlq7xjZeNGv/1yca2/TB4mJ6T0+UyvWqqo6iv5vMyklw1WJykoNKhUUMi7EvA5gRcNeMn7TJkD3yctZQbM+Raxm2hIcwWprPVu2pS/fNn0J6upXN+PbztO3qNMr590/tbZpzKf86lqTLcl5C07NZVf0kjvV5nf3CAfVX4TdUevV6Cg0KBSQaHMiSvvB1C2S2DZR9m0+DS6ugJNenqHqlc+ZQPuzMz6E86GDWu3WebEncwTV/75LtfUI71PyXz23aeyeVF0Ve6z3slJ/+PGJUQHibKPKvlYFPAabn9QUGhQqaDgqjdP/lCLeods27b2B9N0HX9RaSF9pZdePtmbI7lcb/+SA9l6B7JvLx5X4ChTNdfL87xR4mWratLfYdW8a/p7rNOTKES36HTaXN+Pb1WcS1F1aNP7lfeoko95io7NCqUFBYUGNVJSyOqHnaxjzCpCJ+srmzx4fQ7aXruBz/qyTvg+Rd9t27KDYzoviq4YXXnu6slT9ip+40b/qjJXPjTx8L1aLPoOQ5UUeo8NG9yNyGWuhvMuFnzWkZeGvDaRsr8X3+rFvM+k97OJ6qcUBYUGNdKmUOWk1vtsyCvNvB+8q57X5yAt+qH6Nr67ThI+eV6lJ0/ew+cE4psPdR5lTgiuE4xPXXhekC4TVF3dTX1/L67lfNaR9xsqCvRVjtEyFwJVez0V/SYyKCg0qJHeR67l6pwgtmxp/krU5+TnOlB9ejhVCXR5VUvpKppe+gdd15xVlVOlG2Z6mbKDuFzHWbI6L6tHTHoZV3Web68tYG0Po4mJ6LjNqkrNqw7zucjIqyY0y09XuufazIy7N1TRI9mtPJmPPr+ZsueDknM5KSg0KMjtOJuoVkgOTBtk3WnRI1TpJuuqtaj3imv+pbw8rZqPZebw8Wnsz1vnII6/pgeO1clDM7+RwK6ShO8Yl6ar+3wHFBYdxz5pL9B6UACwC8ARAI/mvL8I4GEAjwD4OoALfNY7MkGhzI+ozGjPJn6cdUoKvUedq62yP44qP9S8AUd1g1l6TEidLspZP/6qvbeySgs+x13WerKuyJts9yLXHj/JUkxenud1P3a1WwCnSiy9/AxxMZPXtlX3GB62NgUAbwFwkSMoXALgrPj/ywE85LPekQkKvlejyasd3yqpOgff1FRzE99VORn0Ztxset1Zg4XK9GSq88Otux2fq2GXOsdEus0hq6t0r0tsyGOjSjtXL+/aLkGn5zWqsx6g8iC21oNClAbM5wWF1HJnAXjaZ50jExRc1Qp1RzBWPfhc/cabKD2U2dem2wPq9tapGijztltlfETVfut18tK3d1LRFXmTx0vZ9DedrrLHQvr7aer7KGnYgsL7AHzS8f4SgBUAK7Ozs5UzpaqBtSk0XX9cdKObso+munYmP1vUaNzkI92omjcOosliflZjYN4+urrjFl2JZzUSl5lLyyf9RdN2hOp+m05LmeW3bAlzPPnuZ6/0W/Tbr3M8eRqaoADgrQD2A5jxWefIlBTMqtcP+64768CrWpTupa+oF4pPvWzvhJXV/TFU1ZVrv7O6ZpY9kfhOtJb3veTNoZPXLjGoR5mSQm//fL7DjRurXfU3faGT9SgqkSR/C8nvfm7OfcOd9HFQJWAXTaboMBRBAcDrAPwVgJ/xXedIBYWQXNVTZafQKFuCKap7dqWv6FElqJWZjLAo/1wnCZ+SX5UqoEF3rc1Lv8/32luuqO9/2dJZbxuh27vSU7I0NeYo7/stuz+jHBQAzAJ4EsAlZdapoJCSd/VeZuoHnxNemW333nP1Gx9U45/vCTVdLK/aEJx19ei733l5PciG0qLZO4u+1+Ryru241rlx49reRxs31usq7HqkeznlXdH7luZdx0h6XVVKPcNcfQTgHgCHAbwA4CkA1wO4AcAN8fufBPB9APvih1eCFRQSXFcyZa5Yqly91m0TGcTVb5mePnVLCr75UrTeMqWLqnlStcG6LJ/tFOVXyKozn2qeMlzHmu+4lLLHaQmtB4VQDwWFBNePrsxJexDFZJ9tNvnwHTSUt6++6fOddTVZ5+47Q23TeVV2aom6fLZTlF+hLh58B7OVkZfWXueGumms+R0pKDSos0HBZ5Snb9G3bDG5aNtFQvcyKqpGSz+KZn3NWk/6M64qk/TMqa6TXV5VQ97zvO7CvcZc39lmqyhaT9H7RflVpcooK58A/2knesuUnV6kqHdW1WO5oY4oCgoN6mxQGFQ1QNPbzrqCLDsdhW9vnzJXmln3JfCV1wjryhdXZ4AyV/KDuvJveruu6paikkKZ22pmqVICqVN6rrK9hn/HCgoN6mxQaOtkUHfbrmJ2en1501H43he3TACq80P0ufL3zcMq9ysO2b05T92LEtdxUNQDqMx9kbNUrZKr2s5WdnsBfscKCg3qbFAwa+dkUHfbPjdFSa4vbzu+208v5/PDL5uPRevMk7UPrvUM6vv1ydu6VYiu4yCZBiC7N1fdY78or6vsW1GPvLwqwDL3ba5IQaFBnQ4Kw6jNai/X9utcdTZV+ihaV5eqhUKVFMp2lKijbLXOoI7RABQUGqSg0LCiYvYgGkd9BvA1UT9dtp3CtyRTtH918sn3ZN9Em0LZLtXp0kJdWXf/KwrIbZbOa1BQaJCCQgBZP6yydbJ1rhx9ej/5VoO41lE2fT4nKFe6msinMtVCTVThZH2+KB+aKDW46vnzeh+12Y5Xk29QYLTs8FhYWLCVlZWBbpMkhi2fhtL8PHDw4PrX5+aiv3nvHTgQbrs+6z7tNODEifWvT04CL77YTFp801V3X5paR11186HONkLnb0tI7jGzhaLlJgaRGOmQ5eXowJ6YiP4uL7edosjycv5J4NCh6JH3XtntzM8DZHQyJ4Hnnwc2bFi73PQ0sGOH3zqzAkLy9TJ5vmNHtO08RelqIp+y0lAmP5pwxRXRd+NS9rv3/bxrvWU+09XfWhGf4kSXHqo+qqGrRd+i7nquft5lG3HztjM1Vb33R9MNpskqlbK9UppqxG+7V5tP9826jb5V8mpQ7S0BQG0KzRmZoNB2r5+y6Wq6TaGoEbdqPlRpMA2V5x08GZXWdO+wPFUD9iB6ZgWgoNCgkQkKdfuVDzpdQPXeR1nLFjVe1smHsg2mIfN8SHvH9BVNF9HkPlXJq0GM4QjANyioodnDyDQ0d7WRrOl0LS8DS0vA8eOnXpueBk4/HTh6NP9zIfKhq3neZaOQZx3cBzU0y3pdaEDM0nS6tm9fGxCAU8/zGnFD5UNX87zLRiHPhnkffIoTXXqo+qimrlYtNJmuohsMuaZOCKGNbQ67rh6nPjr6fUPVR80ZmeqjcdHBontuldbOncDiYjtpkuZ1+HtW9ZGMry4W3fOqtLZvbyc9EsYIfM8KCjJ6FhejK7O5uWgA1Nxc+1dqTQ2+k24bge/5tLYTMAxUdTSEFhdbL66vMTubXaWpZ8u2AAAFnElEQVQ1Ozv4tEg4I/A9ByspkNxF8gjJR3PeJ8mPknyS5MMkLwqVFpHWdbFKS5o3At9zyOqjOwG8zfH+5QBeGT+WANwWMC0i7epilZY0bwS+56C9j0jOA/icmb02473/CuArZnZP/PzbAC4zs8OudbbR+0hEZNgNQ++j8wD8deL5U/FrIiLSkqHofURyieQKyZXV1dW2kyMiMrLaDApPAzg/8fzl8WvrmNlOM1sws4XNmzcPJHEiIuOozaDwWQDviXshvQnA3xa1J4iISFjBximQvAfAZQDOIfkUgN8AMAUAZnY7gM8DuALAkwCOA9gaKi0iIuInWFAws18ueN8A3Bhq+yIiUt5QNDSLiMhgDN0sqSRXAeTc4X1knQPgmbYT0TLlgfIAUB4A1fNgzswKe+oMXVAYRyRXfAadjDLlgfIAUB4A4fNA1UciItKnoCAiIn0KCsNhZ9sJ6ADlgfIAUB4AgfNAbQoiItKnkoKIiPQpKHQcyQMkHyG5j+RYzBmedYMmkmeT/BLJv4z/ntVmGkPLyYMPkHw6Phb2kbyizTSGRPJ8kveTfJzkYyRvil8fm+PAkQdBjwNVH3UcyQMAFsxsbPpmk3wLgOcB3N27FwfJ/wLgWTP7IMn/AOAsM7ulzXSGlJMHHwDwvJn9dptpGwSS5wI418z2kjwDwB4AVwK4DmNyHDjy4N0IeByopCCdY2YPAHg29fI7ANwV/38Xoh/HyMrJg7FhZofNbG/8/zEA+xHdb2VsjgNHHgSloNB9BuCLJPeQXGo7MS16WWIW3b8B8LI2E9Oi98b3NN81ylUnSfEdHF8P4CGM6XGQygMg4HGgoNB9bzazixDd0/rGuFphrMWTKY5jvedtAP4egAsBHAbwO+0mJzySmwDcC+BmM3su+d64HAcZeRD0OFBQ6Dgzezr+ewTAZwD8w3ZT1JrvxXWsvbrWIy2nZ+DM7HtmdsLMTgL4XYz4sUByCtHJcNnM7otfHqvjICsPQh8HCgodRnJj3MAEkhsB/FMAj7o/NbI+C+Da+P9rAfxxi2lpRe9kGLsKI3wskCSAOwDsN7NbE2+NzXGQlwehjwP1Puowkj+NqHQARPe++H0z29FikgYieYMmAN9DdIOmPwLwhwBmEc2S+24zG9mG2Jw8uAxRlYEBOADgV0b1boUk3wzgQQCPADgZv/x+RHXqY3EcOPLglxHwOFBQEBGRPlUfiYhIn4KCiIj0KSiIiEifgoKIiPQpKIiISJ+CgggAkicSs07ui6cVKLuO60j+VPOpExmc09pOgEhH/MjMLqy5jusQDST6ru8HSE6a2Yma2xVpjEoKIjlIzpN8kOTe+HFJ4r1b4vtcfIvkB0leDWABwHJc0jid5BaS34yX20XyJ+LPHiD5IZJ7AbyL5L+J58x/mOSnWtpdEQAqKYj0nE5yX/z/d8zsKkTz6vy8mf1fkq8EcA+ABZKXI5rC+Y1mdpzk2Wb2LMn3Anifma2Q/EkAdwLYYmZPkLwbwDYAH463cTSe6BAkvwvgFWb2/0ieObA9FsmgkoJI5EdmdmH8uCp+bQrA75J8BMCnAbwmfv3nAPw3MzsOADnTLLwKUXB5In5+F4DkDLd/kPj/YUQljH8J4MVmdkekGgUFkXz/FtG8Qxcgqhra0OC6f5j4/xcBfBzARQD+N0mV4KU1Cgoi+V4K4HA8RfE1ACbj178EYCvJaSC6b3D8+jEAZ8T/fxvAPMm/Hz+/BsBX0xsgOQHgfDO7H8At8TY3BdgXES8KCiL5PgHgWpLfAvBqxFf3ZvYFRFM4r8TtEO+Ll78TwO3xawSwFcCn4+qnkwBuz9jGJIDd8TLfBPBRM/tBuF0ScdMsqSIi0qeSgoiI9CkoiIhIn4KCiIj0KSiIiEifgoKIiPQpKIiISJ+CgoiI9CkoiIhI3/8H6BVlp6abmnIAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYUAAAEKCAYAAAD9xUlFAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvOIA7rQAAIABJREFUeJztnX94HVd557+vFClINhFYFi0lkZRuXWiKGwpa2i2UsrjsExtMEnc3TVY2IqFVHTdgd8smsNrNr33UXw+7a/fZxKlKTYysJg3UAYc4DV2XJbCFrpWQYJKUxE+wjANLbDkRdiyQLb37x71zMxrNOXNm7sydufd+P8+jR/fOnR/vnDlz3nPe9z3vEVUFIYQQAgAteQtACCGkOFApEEIIqUClQAghpAKVAiGEkApUCoQQQipQKRBCCKlApUAIIaQClQIhhJAKVAqEEEIqnJe3AHFZuXKl9vf35y0GIYTUFY8++ugJVe2J2q/ulEJ/fz8mJyfzFoMQQuoKEZly2Y/mI0IIIRWoFAghhFSgUiCEEFKBSoEQQkgFKgVCCCEVqBQIIYRUoFIghBBSgUqBEEJIBSoFQgghFagUCKkRE4cm0L+9Hy23taB/ez8mDk3kLRIhS6i7NBeE1CMThyYw/MAwzpw9AwCYmpnC8APDAIDB1YN5ikbIIjhSIKQGjBwYqSgEjzNnz2DkwEhOEhESDpWCAyKStwikzjk6czTWdkLyIjOlICK7ROQFEfm24ffLReRbIvK4iEyKyDuzkoWQvOnt6o21nZC8yHKkcDeAyyy/HwBwqaq+BcB1AD6VoSyE5MromlF0tnUu2tbZ1onRNaM5SUSyot4DCjJTCqr6CICTlt9Pq6qWvy4DoKZ9Cal3BlcPYmz9GPq6+iAQ9HX1YWz9GJ3MDYYXUDA1MwWFVgIK6kkxyCvtcgYnF+kH8EVVfbPh9ysB/DGA1wF4n6p+3bDfMIBhAOjt7X3b1JTTWhGpISLIspwIIY1B//Z+TM0sbZ/6uvpwZNuR2gvkQ0QeVdWBqP1ydTSr6v2q+iYAVwD4r5b9xlR1QFUHenoiV5MjhJBcaISAgkJEH5VNTT8rIivzloUQQpLSCAEFuSkFEfk5Kcd6ishbAZwPYDoveQghpFoaIaAgy5DUewB8HcAbReSYiHxYRDaLyObyLr8F4Nsi8jiAOwD8ttJwT1KiSBEgRZKFZEsjBBRk6mjOgoGBAZ2cnKzpNelori+CKSWAUm8tj5ezSLKQ5qYuHM2EZEGRUkoUSRZCXKBSIA2HLQKk1qacRohGIc0FlUKT0Qz2bVOkx4qOFTWfWNQI0SikuaBSaCIaYbalC6YIEAChppyh+4cyK4NGiEYhzQWVQhPRLPZtUwTIydnwrCvzOp+ZcmyEaJQ8aIYRbVFh9JEDjRJ91HJbCzQkxZRAsHDLQg4S1RZTCgKPIqQiIIzYygpGH5ElFNm+7dIzrLb3GGbK8UPnbzFolhFtUaFSaCKKat928XWk4Q/xTDmt0hr6exGUI6ltxBbNVEuhUmgiimrfdukZptV7HFw9iN1X7i6kckyLem/oajWiTdrRqPfyjYJKockYXD2II9uOYOGWBRzZdiR3hQCYe4BTM1OVFy7N3mMWyrEoDUUjRJjVakSbpKMRp3yLUifiQqVAcsfWA/ReuLR7j0HlCCDxC5x2Q1xNY9II9vhajWiTdDRcy7eelTOVAskdmwPYe+Gy7D1W+wKn2RBXK0ujzKBOU2mbSNLRcC3felbOVAokd7yeoYmjM0cz7T1W+wKn2RCbZNm4d6NTY1jkCLOkZNXrTtLRcC3felbOVAqkEAyuHkRfV1/ob94Ll5U/pNoXOM2G2HZNl8awqBFm1ZBVrztJR8O1fOtZOVMpkMKQV4NW7Qs8umYUbS1ti7a1tbQlkjvqmlGNYVEjzJLg+VZMEw7T6HXH7Wi4lm89K+fz8haAEA/vxRo5MIKjM0fR29WL0TWjmTdoo2tGQ2fQxnmBy4sIGr9XI0uQqMZwcPVg1WU2cWii5s8heP2ocsir1+1SvnnV5TTILM2FiOwC8H4AL6jqm0N+HwRwEwABcArA9ar6RNR5meaCZEE1jaCpN5s0bcaWB7fgrsm7QlOSVHNeV2qZZsJU7lEpSdKUJ28FWCtc01xkOVK4G8D/BPAZw+/fBfAbqvqiiKwFMAbgVzKUh5AKYQ1B0oY2bafi/mf3GxVCLUwQNht+mo1lUPl4PhPAXnZ9XX2pNdw2GRpRMbiQmU9BVR8BEJ6WsvT7P6rqi+Wv3wBwYVayEOIn7WgWkxmjRVoSndPWINbCP1CryBmb8jGVqTdKSqsM6jl0NCuK4mj+MICHTD+KyLCITIrI5PHjx2soFsmKPGd7pt0QmOZZJE3JbWsQa9F7rVXkjE351MpRW6RV+opC7kpBRP41SkrhJtM+qjqmqgOqOtDT01M74Ugm5D3b05ZWIwm2RHsuaROCDU/ekSu1ur5N+dQqiqpIq/QVhVyVgoj8EoBPAbhcVafzlIXUjryH7KaGQCCJX/rB1YOY1/nQ30xKyKQcAeQaVlqrBjlK+dQiT1fcVfqawayUm1IQkV4AewFsUtVn8pKD1J60bNZJh/eja0YhWBoyqtDEL/3EoYnQcwJmJRTl0M0zcaHr9asxsRRhTkXcVfrqYUZytWQZknoPgHcDWAnghwBuAdAGAKp6l4h8CsBvAfDG7OdcwqUYklr/pBHCWW3YpNwW3oAnXYXOdE8CwfiG8VCZ6n0lvEZeIS3tMGMg/9DX3FdeU9VrVPX1qtqmqheq6l+p6l2qelf5999R1deq6lvKf5HCkngU1VGWhs06jgkqrByiUmrExdSDVKjxxa/nVAhA7cyAedTjpHXUJGvefrQ45O5oJtlQi0qY9GU1DdkB90yYriYoUzmsW7UuVWeqLWLIRJ4O5TQa2lqYAWvZmPrlGDkwgqFLh2KZtmyy5u1HiwOVQoOSdSWs9mUNS40c53yuvWxTOex/dn+kPTtOw2kKSz09d9p4XF429bQa2jRGOlGy2OpxmiOIMDl2P7Ebo2tGnf06NlmjFGiRRvWZ+RSygj6FxZjslFnbq9O2ucY9n6s9O2k5JLGXTxyawNaHtmJ6dnEgXa3s7K4267SeXRo+hShZTM/Pu1Za/ow0ysRW13q7eo3nX7dq3ZK0JlnUmdx9CiR7bL2srO3Vac96jXs+11520nIw9fq2PrTV2KMbXD2I5e3Ll5wrjRFaVE8yTu8/rWfn8gyi5I6SxfScWqU11ZFwGmViq2thI8m2ljacOHMCOyd3LlEmeZqWqBTqGNtwNWt7dVqmA6/BaJHwqriiY8WShmXLg1vQv70fm/ZuAgCMbxg3Du9NL+PpudPWobqpMZienbY2vFmkiHBp8OOYC9PsMNhCV13kjpLFVI/jzgkJw6X+xSkT2zsXVKDdHd0QEbx89mXj+fIKf6VSqGNsDVDW9upqlU6wwTC95C/9+CVc94XrFjUsOyd3OtvDTS/j9Oy09XjXxiDY8GYxQnNp8OOsOVALB/fEoQkM3T8UKbfLBLawemxy4K/oWOEsX1T9i1smUe+cX4Eub1+Oufk56/nyikKjUqhjohqgOBOg4jq6qlU6YQ1dGPM6H/nyuJh1bC9jWI/atm50EH/D69Lgxi1rFydlnIlzWXcYvAbXpTfvydLd0V3Z1nFexxJ5g0EJp+dOh5771NypJdFLYWVtqn+t0lpVmbi+c1GjAIFg3ap1uTifuchOHZPG4jBA8vTBtsVGopyeaQ+Np2enKw5em/yu5h3vuI17N0Ze29/wRi2ukqSsTU5K77ojB0aMDk5TXUhjIR4TUQo/TFHNnputfJ6enTaWSdTiO3Pzc5UZ4UlSc8/rPPSWV8oyqwlnpmcKlJ7bey5+D3Y/sTuXlN4cKdQxafX40g5frcaWnBZB+b0eoymSxdSjts0zABYrYe8aNl9HkrKOGn2YGhfbxLkssSn8sE5LnDJxGWFOzUyh5bYWq/nKJf9V0tBdl5GgaSTaIi3YPLAZh08ezm1eA0NSHShySGoSgr0fW48lzZQP/vC+sB5fW0sbzi6cjX09E578Ub1LW/if7ViBYPPAZtz5vjtrEiI7cmAEUzNTaJVWzOu8MZzRo7ujG8vbly/q5QKLRzHrVq3D/mf34+jM0Yo9/uTsyap6xabn3yqt2H3l7kov3pPDpKiBV8I5o0Kt08Srp0nCVOOE6dpCmE11tZqQcteQVCoFBxpJKUQ1kEGSrHLl2vCFDc3DXpIgy9qWYWXnyspxp+dOG4/p6+qL/D3q/iYOTeCD938QC7r0ZYzbgFQTDx/32bW3ti/yn7S1tEFEIn00Hklj5aMaxrj34T/eU4xZ4tXTJAo87vO15cwyKfoTN55wuIulcJ4CCcXVweuRZLarawROmFNux9odkQ7eM2fPVI6L8p9MzUwZFYJAnDOQhikE4BVTiauvoprIn7jPLtj4n10466wQgHBzhYtpJCziq+O8Dmzauwn92/ux9aGtse7DL8vomlG0ZNxsefU0SSRZ3JBkW86svKBSaDJs9t6wRWKA+AvFVJNXKCwaJYj3Uno9zqiRRdR5orDde9wGpBo/UB5x6/5rhtnYN+3dhC0PbllynKfwxzeMY/bc7KIQ4KTPy5NlAeEK2lR/4+Cvp0kUuCkkNq1QZVNK7zShUmgybJVwQReMoY1hieb6t/dDbhNs2rtpSc6YoUuHrGGGNkwzg4HFETVxe85+4kRp2RrjJA2Iv8EEUOlBTxyasPbEXRuQzrZOq1KNg/+aYeWtUNw1eZdxJBnnGfV19Vkd+71dvVYFPa/zS+qvS1ixKQw1rgKfODSBH/3kR0u2t7e2G+uaaW0PE7WYu0CfggON5lPYtHdT6PC0VVrxmle9JrQn193RjR1rd1Rsuiabp3//2XOziXPT2ByKfV19kQ7KMHmCTlfXcFqTTyLo2/A7bV2uEeZoD9r849riPR8JgNh2+zD2bNhTuYeoZ3Jk2xHnIIYg3n2a6qYni+33YJ30AgH2P7vfKEdUkEGccFSTfyDKD2Ba2yOOrC7Q0ZwitVQKWcVF+9ny4BZjxEpbSxsWdGHJxKNWaUVrS2ssm7QJf8Nlute4Djgb7a3t2HX5LqdydG2soxrwKEz3F0Ywass16Z4/YikJwcbMJrO3mFCw7GzPq1VasaALi559VMMap9yAV+ratZ+/dklkm9fRiRMVZHu+SSPLbPfs2pFxIXelICK7ALwfwAuq+uaQ398E4NMA3gpgRFU/6XLeRlYKaWSdDJ7PNolq6P6h0FmnLdJidKymha1RBRD6UiZRCEC8iA3bCwqgIpOpjFyzasYJrQxrVOJ0Hkz3tKxtGX587sfGFA9hSsbUU/fMPnEabNN9RUUu2UYLYdcY3zCO675w3aK6ZuooRI3GgiMif4iw9990jIm033sTRVAK7wJwGsBnDErhdQD6AFwB4MVmVAquZook6ahdKlpWMd+dbZ3oOK8jkUNxWdsyKHTJS9nd0Z3YQRkntjtOqmbXa4U14HF68NUsAeld32X04yldW5iuaZSZ5PnY0qLbFJ6rucW7BhCurMKuHzUSMSkZG2EjkuA9rlu1Dvc9eV+lDE2jmGrIPSRVVR8BYHSVq+oLqnoQQHqzleqIsEgO00t1dOao1QEZ9pvL4iS2YX1cPGeZ54y76hevin0OAHj57MuhDe/y9uVGJ2SUoy6Oc860b4u0ONnng8ebZsWGRWiFkUaiujCH6QXnX7CkUfMUgi1M98733YnxDeNLHNlxFUJUCo6w/EFevbWd049XdqZAgbDGPyrCq7erF5u/uNmqEIJyeGk7bDOlP/XYp3Bq7pT1mFrlQcrUpyAi/QC+GDZS8O1zK4DTtpGCiAwDGAaA3t7et01NZTt5JeT6qY8U4thGbU5bYKlDMapHG5zU5KezrRNDlw4tyrsSRVjPMq7tN4okNmsgviMRQKj92YWwa9kmM/lHDGH3kUVv0aPaRZjSeL56izqbwVxmpQ9dOhTq6Lf5p8Y3jDvXW+8aOyd3Jrpfz0eQxCeShnkp95FCmqjqmKoOqOpAT09P3uKkgmvMudebDOv1D90/FDoR6MzZM9bevkkheL38O993pzU9sZ89G/aE9iyj7q+tpS3y3H56u3pDe7w2hWBb+9nUgweAC86/IJZsHmFht1HpzY9sO4Luju7Q+1jevjyz3EVx4ubDeqnVzpno6+qLlVvIFtrqr7dhIwxT2KdCnbPjdnd0Y2z9GPY/uz/J7QJ4ZS2OOEzNTDmlIE+TulAKjYirSWNs/Zhxwsq8zhuH7aa0xSaCs3u9BitKMZgqpu3++rr68OkrPm2MpTeZAfxyeS++Sb7ujm7r2s8mZTpyYCTxBKHgkB+IbnwnDk0Yn+HUzFRmpgLXeRWmhtt13QLA/DzjJMIzKSGv3gJLFb/H4OpBY+chrJH2K/fujm7s2bAHJ248gcHVg5HKMM05Ih5pLCgUByqFnFi3al2kLbyvqw+DqwdrMmHFdI3RNaPWXr2pYpoaHf/IIiylRWdbJzYPbHaeMGSSz8urb2p4bA3xsvZlob+5TDJyWURGIJUGf+tDW63nS5JmxAXXiVmm8gPcJoa1t7Zj88Dm0ImMNlu/68S93q5epxGHaeTs3x42Q96f0tsmB/DKaMIlVUuQtpY2tLe2xzomq3Yhy+ijewC8G8BKAD8EcAuANgBQ1btE5KcBTAK4AMACSpFKl6jq0imBPuoh+ijKTuoyCanaBGJxaGtpwwXnX2DMjumSEM6PKaNnmL04jXkZK/9spTFqK+4kNxumyKggpsR/LpP+TFQbgeRKnAy64xvGK/u2SEtoj9bzi4TZxKMi1KLegagkef4ys0UseesnuMwXWNGxAqfmTi2J2vKy5XoEy3F6dtq4MFBw3o6LiSlLnwInrzkQRym4hILaUgsHJ/P4z2tqmOM6sASCFR0rcHL2ZGglN8Wouzi7ahVz7cfmNDU1bGHO+yhcG8K4GTFdr500ZbJHks6KSYkF7zGrZxCcuBcmv4vT3CV7qWuIdlQnKsjEoYnQ4AXTXAmXeuKfae5KQzma6wmXUFDTA1/QBesyfmFD4PbWduxYuwNHth1xzqGyomMFTtx4ItbylNWaG7JcHMRmWjCZsXas3VFxQse5jt+nsfvK3bESplVjA26RlqrCEV1MLKbcRjYfj4ftGZjue3p2GkOXDln9VsGlO8McyS5Oc5sPJSpEO8jZhbNY3r7caZlboFSuYdFsr25/deixUUvBemblrKBSSBmbndR7KU3YbIQuFcvVxuh3pMZJ6euy/mzcFMFpYHvhbcrMtrKaS0MYN2Ga6/M5r2XpKrnzOh9r9a8gLsralsY56h5tz8B237uf2I3RNaPG5+BSZqZr+9c4HjkwUlFA/vsAEPlehhFnf1O5mgIabJmC05i3EgWVQsqYKnGrtEb6EGwP21ax/CMQl9GCX8YkKX2TZPJMyykWdm1b4xxlMjE1KK7ObhdFabtWGOe3nm99jmfOnsHGvRtjjRpMjZi/XpmekWdiMd2j36HvjWb9ZWa77zNnz2DrQ1tD7e1x0637n5c31yaYvXd0zeii7LRh4Z7++7ZFyLmWfZJ3YnD1IE7ceAJ7NuypernduNCn4EAaPgWXzJY2J6zNkWeyyy5rW7ZkYZWgfT+uD8C2PxCesygtn0KaskalHIjKdhpn0pXp3Gk4wF3KdsuDW4wTrqKWSI1KfRGnjDfu3eh8X9VO3LM5jV18SZ4/wrXsbOThZwuDjuYUSSP6yCVCIngOlwglWwSHbZZnlLxxk6vZlg9Ma1ZuWksdVrvGblpO97RmfUdFgJkIm9Fri5QKc7CaruFfj9kjaWbYJFSb18s0kziIq8O3FtmPo6BSSJE00lzE7S2YwiyDEUpRGSNN0RtJF2mP+7IFX+5qXo64qRmyWmM3rfWYTXUiSTJBL6wSiLeWs/84Py4NeNQI2EWZmqg22qoahRtnPeg8evxJYfRRwYjjlLTNcp3X+UW23ShbvWczDkafTM9OL1oi0dV5Gdc3ELWcYxynaVzbrGnWrW02rouj3NWZbmpQvO2mOhE2+cnmYwjat11XO3ON+jERlU4l6MgeXD2IoUuHYvu94mDzr9lmG4etvuZaBv7IwlokrMsaKoUa4uqUtIVv+huAiUMTxgkxHt7LFdVQmMJGg5XdNbtn8PomGeKEq1az6L0rLorHVTm5zKINqxNhymJ8wzj2bNjjlMPHpTGzZSkF7IrTz7zOW+tDUJb9z+6PHGkmfab+TgeARdexzTbubOvE7it3O4e6BvFHFrp0doquQKgUMiL44Lc8uCWyIkTNYwBeaQDCpuQH8b9cLg1F2DrMwcrurb/sNVhRabb9L3e14apxQ0BNIX+23EYuisdVOZly1rjkpTIpC1OD6hJF5OHNwDWV28Sh8LWGw/AaW1M9CMpie9ZhkWNxGk9bx8dLVxGnDrlGi4VFFto6WdWMlmsBlUIGhD34nZM7rRUh2MuxMTUzFZrQDTAvQu7S6wnuY+rZ7392/6IJXDb8L1u14apx/RFJQwGjGg3XhsVknnHJPmvCdGyLtFTqkynfknf8+IbxRSkZgpjmxIThrQHgOpHPFvbqV4BJGk+bwvE30q4j9uBz7u7oXpKfqLOtM1bCujwmd8aFjmYH4jqaXZ1cLs7LuJgcdC756IMNm6ujdvkfLcfLZ19esp8XeeR3bruk1AgjSVhf3qGAWVzf9hz90UFJAwmA+MEE3nKnLkrbtUySRI65rJqWRaqQOJGF1a5jUQ10NOeIqznExXkZF1NPLKzX093RvaSnu+XBLTjv9vMgt5kTtwWv8Rfr/2JJptK2ljZc9YtX4drPX7vIuX1u/lzodaNI0sOKa25Kmyyu750zzFxzduFsJXhgenYas+dmMb5h3NobDjPRxHXyTs9OVyYRRvXAXcskiakxytyTxgTKsHuMMif6y7hFwpvcWmRCdoUjBQfyHCnEWfu22l6obaJOUKaoNWdH14yGTmTzjj9x44lYsuXZwyoirr35JHMykqwulnYW1yQjBaB0T1lOoLRdN2yUFDcjcpakOlIQkXRXjWhwXBxU7a3tTs7LHWt3GO3Iy9qWpdoLHXvULUFc2GIyYT0okzKLG4MPZJ8+o95wvW9bz9rmM4q7UEzaua2SRprllR7CNEoyOb9Nvr8iRCa5mo++ISKfFZF1IuKWirOJ8Q+RTQQzJNqG1aaFZM4unMXomlHnbI1RxFmtrdbOMZcsl0UN8csC18gYm/KwmWh2rN3hNJ/A5TpJqNb0FicnVbXY6p+pjMMyIhclMsnJfFRWBL8J4DoA/xLAfQDuVtVnLMfsAvB+AC+o6psN59wBYB2AMwA+pKqPRclSD+YjP2mZPWwLyaQ1bD/v9vNiKYaoezDJnMR8BIQP0QEUIq9MHgRnqMd14keZaGwL0/hplvIOI2kqk7D3Nsk6K3FI1XykJf5eVa8B8LsAhgD8XxH5ioj8K8NhdwO4zHLatQBWlf+GAcQzYtYJcc0epl6HKbY+zWH78NuGQ7cvb18euj2qd7hj7Y7QJQanZ6cT9ejDen9FDvHLegTjL48TN57Arst3xepZR5lobOtf5+W8LxpR9S+OGcz0LlebNj0uzj4FEdkqIpMAPgbgIygts/mHAP467BhVfQSAbQX0ywF8pqxwvgHgNSLy+ljS1wFxKoVt+FgLm/qd77sT1w9cX4lsaZVWXD9wPe56/12J7bteQwUsTtWQVgXPY/0GF1xMAWkrjbgmkygTjc3PFZVKu1nMeVH1L44ZzOVdrkWHx9V89AyAcQCfVtVjgd9uUtU/NRzXD+CLBvPRFwH8iap+rfz9AICbVNVqG6o38xHgPunKNtQMy9hY65j7arI8Jo0myeu81ZI0GV7Ret1xn3u93FdapFn/XBMGJo24SzVLqoiIJmgV01IKIjKMkokJvb29b5uaqn6SVxzSyJLqQpT/oQjpd5OSVUhpmo1QmuUbdb9FVWbV0qj3ZSJtJeiyhkrSskx78tqXROQ1vpO/VkQeji3VYp4HcJHv+4XlbUtQ1TFVHVDVgZ6eniovW1yiTES1jKhIm6zMX2lNEEs78iPqfqsxe8U1z9TSnFNUc15WpD1B0f+Ox10DPC1clUKPqr7kfVHVFwG8rspr7wPwQSnxqwBmVPUHVZ6zrqlFFtC8yPLe0lCWaTuso+43qZKMq7xqHebYjPNJsuqs5TUj31UpzItI5amKSB9gn04pIvcA+DqAN4rIMRH5sIhsFpHN5V32A3gOwGEAfwlgS2zpG4y80zJkSdHvLe0eblInbpSSjKu8ah2dlWfHphEd3HlYB1x9CpcBGAPwFQAC4NcBDKtqtSak2NSLo7me7f/NSB628CR1xJbeoq+rb8m58kgPkkfdL7qDuwjtQerLcYrISgC/Wv76DVWNP/soBepBKcStoEWoMM1O0RsVD9c1sqOWlGw0x2+RHdxFqVtZZEk9H6V5Bz8CcImIvCupcI2Oaci+9aGtS/YtytT2Zqfo5i0P01oJwdGAZyJqZD+VnyI7uIs8wTKM81x2EpE/BfDbAJ4E4I05FcAjGclV15gqoj/FsIetwhStQWp0vNXNiow/0Zo3sjRl5D06czR0/0YciZrKoQgO7iIrrDCclAKAKwC8UVV/kqUwjYLtRfV6B95L6rK8YiNCk1lygsrLZDrxhzI3etmaJncWYURUZIUVhqv56DkAS9N0klBsFTG4yLeJOEtU1lvEBU1m6dIsJiIbRTb/1dvzcY0++lsAlwI4AKAyWlDVj2YnWjj14GgGzBlCW6U1MhNplktUupB1L77ITsF6hSOvYlOE55N2mouhsO2qal+1PQOKrhS8hz81MxUaERKV18RTGl6+I1vFyaJxrUWkBFdRI6T2pJ06ezdKayh8Q1V3e3/VCtlo+M0iAKDQSmZQbzhrS0fc2dZZGUW4mFSycGDVIlKiGWe9ElIvuKbOXg/gcQB/V/7+FhHZl6Vg9UhYg6rQSs/dtIqa9z1uY5xF41qLSIl6s7ES0ky4OppvBfB2AC8BgKo+DuBnM5KpbnFtUIMrmoqIcd1iW2OcReNai158kZ2ChDQ7riGpZ1V1JtCY0fgbwCX0bOTAyKIlEwFgbn7O6IC2NcZZxKBq8BeEAAASdklEQVTXKrSvGcIkCalHXJXCkyLy7wG0isgqAB8F8I/ZiVWfuDSotiX3go5o19XN0mxcm2WyEyEkHFfz0UcA/CJK4aj3oJTqYltWQtUrLmYRU8/f74jO26RSz+s2EEKqwzkhXlEoekhqFEVJjkUIaS5cQ1Kt5iMR2a6q20TkAYSsn6CqH6hCxqaE5hlCSJGJ8imMl/9/MmtBmgk6WQkhRcXqU1DVR8sfJwF8VVW/oqpfAfA1AAejTi4il4nId0TksIh8POT3PhE5ICLfEpH/LSIXJrmJZqQecx4RQoqPq6P5AAB/QHwHgP9lO0BEWgHcAWAtgEsAXCMilwR2+ySAz6jqLwG4HcAfO8rT1DChHCEkK1yVwqtU9bT3pfy507I/UJrsdlhVn1PVOQD3Arg8sM8lAP6h/PnLIb+TEOpt0Q5CSP3gqhReFpG3el9E5G0AZiOOeQOA7/m+Hytv8/MEgA3lz1cCeLWIdDvK1LTU26IdhJD6wVUpbAPwWRH5qoh8DcDfALghhet/DMBviMg3AfwGgOcBLJnWKyLDIjIpIpPHjx9P4bL1DRPKEUKywjVL6kEAbwJwPYDNAH7B54Q28TyAi3zfLyxv85/3+6q6QVV/GcBIedtLIdcfU9UBVR3o6elxEblm5OHwZUI5QkhWWJWCiLyn/H8DgPUAfr78t768zcZBAKtE5GIRaQdwNYBFmVVFZKWIeDJ8AsCu+LeQLV4jH9boRzl8s1IYTChHCMmKqHkK70LJEbw+5DcFsNd0oKqeE5EbADwMoBXALlV9UkRuBzCpqvsAvBvAH4uIAngEwO/Hv4Xs8Bp9AIsafaDUMEc5fP0zl4PHVgvnOhBCssCa5kJEtqrqDhF5p6p+rYZyGallmovKyma3ovRXxlsfwbaCmCljKpecJITkQVorr11b/v/n1YtUf5iieaZmpjBxaMLq8GWEECGkHolSCk+LyLMA3lSedez9HRKRb9VCwDyxRfMMPzCMdavWGR2+jBAihNQjUWkurgHw6wCeRcmv4P29H+F+hoYiLMrH48zZM9j/7H6jwzcqQohpKgghRSQqS+oBVV0jIg+r6lIDeYPjOXI3jmwM/X1qZsro8LVlQw2mz07bCU0IIUmJcjQ/BeB3APwVgCWtlao+lp1o4eS1noLf0ezRKq04d/O52OerOLAD0AlNCMmKVNZTAHAzgP+C0sQzL322t1CzAnhPYgkbgLA1lV2gE5oQUlSsSkFVPwfgcyJyM4DvArhYVW8XkV4AP10LAYtMX1dfouNM4ap0QhNC8sY199HPAPhVANeUv59CKS1205BmWgmmqSCEFBVXpfB2Vf19AD8GAFV9EUB7ZlIVEFtaibiRRExTQQgpKlZHc2UnkX8C8GsADqrqW0WkB8CXyonsakpejmZTOQUjiYBSr5+NPCGkSKQ1o9njzwHcD+B1IjKK0nKcf1SFfHWNf2QwdP8QF7whhDQMUdFHAABVnRCRRwGsQSn66ApVfTpTyQpG//Z+HJ05ihUdK3Bq7hTm5ucAmCOQGElECKlHnJQCAKjqPwP45wxlKSSef8CLFpqenXY6jpFEhJB6xNV81LQkMQMxkogQUq9QKUTgagZqlVZGEhFC6h5n81Gz0tvViynY0z4x2ogQ0ihkOlIQkctE5DsiclhEPh7ye6+IfFlEvllOyb0uS3mSEGYGamtpQ3dHN0cGhJCGI7ORgoi0ojTr+b0AjgE4KCL7VPUp327/GcB9qrpTRC4BsB9Af1YyJWFw9SA2YiP6uvqWZDslhJBGI0vz0dsBHFbV5wBARO4FcDkAv1JQABeUP3cB+H6G8lQFs5cSQpqBLJXCGwB8z/f9GIBfCexzK4AvichHACwD8JsZykMIISSCvKOPrgFwt6peCGAdgHERWSKTiAyLyKSITB4/frzmQgLgKmmEkKYgS6XwPICLfN8vLG/z82EA9wGAqn4dwKsArAyeSFXHVHVAVQd6enoyEjcc/+Q1hWJqZgqb9m6C3CZUEISQhiNLpXAQwCoRuVhE2gFcDWBfYJ+jKKXOgIj8AkpKIZ+hgIGwyWuKUnI8bxlNKgZCSKOQmVJQ1XMAbgDwMICnUYoyelJEbheRD5R3+0MAvysiTwC4B8CH1CVtaw2JmrzG5HeEkEYi08lrqrofpTBT/7abfZ+fAvCOLGWoFpfJa0dnjmLi0ARGDowwbJUQUtfk7WguPC45jFZ0rMDwA8OL/A40KxFC6hEqhQi83r63HrNAFv3e2daJn8z/hGsqEEIaAioFR45sOwK9RTG+YXzRMppDlw7h9Nzp0GO4pgIhpN6gUojJ4OpBHNl2BOMbxgEAOyd3GvflmgqEkHqDWVITELYucxhcU4EQUm9wpJCAkQMjkQqhu6Ob0UeEkLqDSiEBUb4CgWDH2h01koYQQtKDSiEBUb4ChXKUQAipS6gUEjC6ZhSdbZ3G373wVUIIqTeoFBIwuHoQY+vH0N3RveS3zrZOOpgJIXULlUJCBlcP4sSNJ7Bnw55F8xa4NCchpJ6RguWfi2RgYEAnJydrek0RQb2VEyGE+BGRR1V1IGo/jhQIIYRUoFIghBBSgUqBEEJIBSoFQgghFTJVCiJymYh8R0QOi8jHQ37/HyLyePnvGRF5KUt5CCGE2MksIZ6ItAK4A8B7ARwDcFBE9pVXWwMAqOof+Pb/CIBfzkoeQggh0WQ5Ung7gMOq+pyqzgG4F8Dllv2vQWmdZkIIITmRpVJ4A4Dv+b4fK29bgoj0AbgYwD9kKA8hhJAIiuJovhrA51R1PuxHERkWkUkRmTx+/HiNRSOEkOYhS6XwPICLfN8vLG8L42pYTEeqOqaqA6o60NPTk6KIhBBC/GSpFA4CWCUiF4tIO0oN/77gTiLyJgCvBfD1DGUhhBDiQGZKQVXPAbgBwMMAngZwn6o+KSK3i8gHfLteDeBeZXIhQgjJnUzXaFbV/QD2B7bdHPh+a5YyEEIIcacojmZCCCEFgEqBEEJIBSoFQgghFagUCCGEVKBSIIQQUoFKgRBCSAUqBUIIIRWoFAghhFSgUggwcWgC/dv70XJbC/q392Pi0ETeIhFCSM3IdEZzvTFxaALDDwzjzNkzAICpmSkMPzCcs1SEEFI7mnKkYBoNjBwYqSgEj+B3QghpZJpupGAbDRydOZqnaIQQkjtNN1IwjQZGDoygt6s3J6kIIaQYNJVSmDg0gamZqdDfjs4cxeiaUXS2dS7aHvxOCCGNTNMoBc9sZKK3qxeDqwcxtn4MfV19EAj6uvowtn6shlISQki+SL2tbTMwMKCTk5Oxj+vf3m8cJXS2dWJs/RgGVw+G/i4iqLdyIoQQPyLyqKoORO2X6UhBRC4Tke+IyGER+bhhn6tE5CkReVJE/jorWWxOZE8hcI4CIaTZySz6SERaAdwB4L0AjgE4KCL7VPUp3z6rAHwCwDtU9UUReV1W8vR29YaOFLo7uisK4bovXIe5+TkApaik675wXVbiEEJIIclypPB2AIdV9TlVnQNwL4DLA/v8LoA7VPVFAFDVF7ISZnTNKNpa2pZsPzV3Clse3IIP3v/BikLwmJufw9aHtmYlEiGEFI4slcIbAHzP9/1YeZufnwfw8yLyf0TkGyJyWVbCDK4exAXnX7Bk+9z8HHZO7sSCLoQeNz07nZVIhBBSOPKevHYegFUA3g3gQgCPiMhqVX3Jv5OIDAMYBoDe3uRzCU7Onkx8LCGENANZjhSeB3CR7/uF5W1+jgHYp6pnVfW7AJ5BSUksQlXHVHVAVQd6enoSC7SiY0XiY+l4JoQ0A1kqhYMAVonIxSLSDuBqAPsC+3wepVECRGQlSuak57IQZuLQRFUjBS8dBhUDIaSRyUwpqOo5ADcAeBjA0wDuU9UnReR2EflAebeHAUyLyFMAvgzgP6pq6kZ8b+Kaorq5Bl46DEIIaVQy9Smo6n4A+wPbbvZ9VgD/ofyXGWH5jpLCpHmEkEamKdJcRDXky9uXV9JaXD9wPbo7uo37MmkeIaSRaQqlYGrIW6UVezbswalPnMLCLQs4su0I3tH7Dsyemw3dv7OtE6NrRrMUlRBCcqUplIIp++nuK3cvyXdkMjW1Sqs1PxIhhDQCec9TqAleQz5yYARHZ46it6u30uPv396/aJvJ1LSgC1QIhJCGpymUAlBSDP5G3bQC24qOFaGzmOlLIIQ0A01hPgrDth5z2MI69CUQQpqBplUKJjPRydmTSxbaAUDTESGkKWga81EQUyptbwW2RUpgWw0FI4SQHGnakYIpIolmIkJIM9O0SsG0HjPNRISQZqZp1mgmhJBmphBrNBNCCKkvqBQIIYRUoFIghBBSgUqBEEJIBSoFQgghFTJVCiJymYh8R0QOi8jHQ37/kIgcF5HHy3+/k6U8hBBC7GQ2o1lEWgHcAeC9AI4BOCgi+1T1qcCuf6OqN2QlByGEEHeyHCm8HcBhVX1OVecA3Avg8gyvRwghpEqyVApvAPA93/dj5W1BfktEviUinxORi8JOJCLDIjIpIpPHjx/PQlZCCCHIPyHeAwDuUdWfiMjvAdgN4D3BnVR1DMAYAJR9EEsz2TU2KwGcyFuInGEZsAwAlgGQvAz6XHbKUik8D8Df87+wvK2CqvpXs/kUgD+LOqmq9qQiXR0hIpMu09MbGZYBywBgGQDZl0GW5qODAFaJyMUi0g7gagD7/DuIyOt9Xz8A4OkM5SGEEBJBZiMFVT0nIjcAeBhAK4BdqvqkiNwOYFJV9wH4qIh8AMA5ACcBfCgreQghhESTqU9BVfcD2B/YdrPv8ycAfCJLGRqEsbwFKAAsA5YBwDIAMi6DukudTQghJDuY5oIQQkgFKoWCIyJHRORQOQ1IU6wuJCK7ROQFEfm2b9sKEfl7EXm2/P+1ecqYNYYyuFVEnvelhVmXp4xZIiIXiciXReQpEXlSRLaWtzdNPbCUQab1gOajgiMiRwAMqGrTxGaLyLsAnAbwGVV9c3nbnwE4qap/Us6j9VpVvSlPObPEUAa3Ajitqp/MU7ZaUI5MfL2qPiYirwbwKIArUApGaYp6YCmDq5BhPeBIgRQOVX0EpWg0P5ejNLkR5f9X1FSoGmMog6ZBVX+gqo+VP59CKVz9DWiiemApg0yhUig+CuBLIvKoiAznLUyO/JSq/qD8+f8B+Kk8hcmRG8ppYXY1sunEj4j0A/hlAP+EJq0HgTIAMqwHVArF552q+lYAawH8ftms0NRoyebZjHbPnQD+BYC3APgBgP+WrzjZIyLLAfwtgG2q+iP/b81SD0LKINN6QKVQcFT1+fL/FwDcj1L22Wbkh94M+PL/F3KWp+ao6g9VdV5VFwD8JRq8LohIG0qN4YSq7i1vbqp6EFYGWdcDKoUCIyLLyg4miMgyAP8GwLftRzUs+wAMlT8PAfhCjrLkQiAtzJVo4LogIgLgrwA8rar/3fdT09QDUxlkXQ8YfVRgRORnURodAKXZ53+tqqM5ilQTROQeAO9GKRvkDwHcAuDzAO4D0AtgCsBVqtqwjlhDGbwbJZOBAjgC4Pd89vWGQkTeCeCrAA4BWChv/k8o2dSboh5YyuAaZFgPqBQIIYRUoPmIEEJIBSoFQgghFagUCCGEVKBSIIQQUoFKgRBCSAUqBUIAiMi8L+vk4+W0AnHP8SER+Zn0pSOkdmS68hohdcSsqr6lynN8CKWJRN93PUBEWlV1vsrrEpIaHCkQYkBE+kXkqyLyWPnv13y/3VRe5+IJEfkTEfm3AAYATJRHGh0iskZEvlneb5eInF8+9oiI/KmIPAbg34nIR8s5878lIvfmdLuEAOBIgRCPDhF5vPz5u6p6JUp5dd6rqj8WkVUA7gEwICJrUUrh/CuqekZEVqjqSRG5AcDHVHVSRF4F4G4Aa1T1GRH5DIDrAWwvX2O6nOgQIvJ9ABer6k9E5DU1u2NCQuBIgZASs6r6lvLfleVtbQD+UkQOAfgsgEvK238TwKdV9QwAGNIsvBEl5fJM+ftuAP4Mt3/j+/wtlEYYGwGcS+d2CEkGlQIhZv4ApbxDl6JkGmpP8dwv+z6/D8AdAN4K4KCIcARPcoNKgRAzXQB+UE5RvAlAa3n73wO4VkQ6gdK6weXtpwC8uvz5OwD6ReTnyt83AfhK8AIi0gLgIlX9MoCbytdcnsG9EOIElQIhZu4EMCQiTwB4E8q9e1X9O5RSOE+W/RAfK+9/N4C7ytsEwLUAPls2Py0AuCvkGq0A9pT3+SaAP1fVl7K7JULsMEsqIYSQChwpEEIIqUClQAghpAKVAiGEkApUCoQQQipQKRBCCKlApUAIIaQClQIhhJAKVAqEEEIq/H+Cm82nYAsuigAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "           |factor| pbt  | ret  | load  | effi\n",
      "Max Eff\n",
      " Absolute  | 13.90 | 0.56 | 1.37 | 1.04 | 1.32\n"
     ]
    },
    {
     "ename": "UnboundLocalError",
     "evalue": "local variable 'sm2' referenced before assignment",
     "output_type": "error",
     "traceback": [
      "\u001b[0;31m---------------------------------------------------------------------------\u001b[0m",
      "\u001b[0;31mUnboundLocalError\u001b[0m                         Traceback (most recent call last)",
      "\u001b[0;32m<ipython-input-16-df20d169dae9>\u001b[0m in \u001b[0;36m<module>\u001b[0;34m\u001b[0m\n\u001b[1;32m     12\u001b[0m     \u001b[0mrow\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mws\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mSimulationResult\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mpbt\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0ml\u001b[0m\u001b[0;34m[\u001b[0m\u001b[0;36m0\u001b[0m\u001b[0;34m]\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0ml\u001b[0m\u001b[0;34m[\u001b[0m\u001b[0;36m1\u001b[0m\u001b[0;34m]\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mquotient\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mfactor\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m     13\u001b[0m     \u001b[0mtotal_output\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mappend\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mrow\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m---> 14\u001b[0;31m \u001b[0manalysis\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mtotal_output\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m     15\u001b[0m \u001b[0mcompare_best_to_sm2\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mtotal_output\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n",
      "\u001b[0;32m<ipython-input-7-fd08ff823e83>\u001b[0m in \u001b[0;36manalysis\u001b[0;34m(total_output, new_factor)\u001b[0m\n\u001b[1;32m     59\u001b[0m     \u001b[0mprint\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m\"Max Eff\"\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m     60\u001b[0m     \u001b[0mprintsim\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mbest_productivity\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0;34m\"Absolute\"\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m---> 61\u001b[0;31m     \u001b[0mprintsim\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mcompare\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mbest_productivity\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0msm2\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0;34m\"VS SM2\"\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m     62\u001b[0m     \u001b[0mprint\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m\"---\"\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m     63\u001b[0m     \u001b[0mprint\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m\"Least Work\"\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n",
      "\u001b[0;31mUnboundLocalError\u001b[0m: local variable 'sm2' referenced before assignment"
     ]
    }
   ],
   "source": [
    "list_of_factors = [i / 120 for i in range(3 * 60, 10 * 300) if (i % 30 == 0) or (i % 40 == 0) or (i % 12 == 0)]\n",
    "focus_output = []\n",
    "focus_output_too_high = []\n",
    "focus_output_too_low = []\n",
    "dif = 0.90\n",
    "total_output = []\n",
    "for factor in list_of_factors:\n",
    "    pbt = fc.get_pb_success_from_interval_modifier(factor / 2.5, dif)\n",
    "    l = ws.sim(pbt, nsims=100, onfail=\"stable\", ndays=1,\n",
    "            difficulty=dif, factor=factor)\n",
    "    quotient = l[0] / l[1]\n",
    "    row = ws.SimulationResult(pbt, l[0], l[1], quotient, factor)\n",
    "    total_output.append(row)\n",
    "analysis(total_output)\n",
    "compare_best_to_sm2(total_output)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##  70% success rate"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "list_of_factors = [i / 120 for i in range(3 * 60, 10 * 300) if (i % 30 == 0) or (i % 40 == 0) or (i % 12 == 0)]\n",
    "focus_output = []\n",
    "focus_output_too_high = []\n",
    "focus_output_too_low = []\n",
    "dif = 0.70\n",
    "total_output = []\n",
    "for factor in list_of_factors:\n",
    "    pbt = fc.get_pb_success_from_interval_modifier(factor / 2.5, dif)\n",
    "    l = ws.sim(pbt, nsims=100, onfail=\"stable\", ndays=1,\n",
    "            difficulty=dif, factor=factor)\n",
    "    quotient = l[0] / l[1]\n",
    "    row = ws.SimulationResult(pbt, l[0], l[1], quotient, factor)\n",
    "    total_output.append(row)\n",
    "analysis(total_output)\n",
    "compare_best_to_sm2(total_output)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##  95% success rate"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "list_of_factors = [i / 120 for i in range(3 * 60, 10 * 300) if (i % 30 == 0) or (i % 40 == 0) or (i % 12 == 0)]\n",
    "focus_output = []\n",
    "focus_output_too_high = []\n",
    "focus_output_too_low = []\n",
    "dif = 0.95\n",
    "total_output = []\n",
    "for factor in list_of_factors:\n",
    "    pbt = fc.get_pb_success_from_interval_modifier(factor / 2.5, dif)\n",
    "    l = ws.sim(pbt, nsims=100, onfail=\"stable\", ndays=1,\n",
    "            difficulty=dif, factor=factor)\n",
    "    quotient = l[0] / l[1]\n",
    "    row = ws.SimulationResult(pbt, l[0], l[1], quotient, factor)\n",
    "    total_output.append(row)\n",
    "analysis(total_output)\n",
    "compare_best_to_sm2(total_output)"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.6.7"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}