Skip to content
Snippets Groups Projects
check_timeseries_stationarity.ipynb 29 KiB
Newer Older
many's avatar
many committed
{
 "cells": [
many's avatar
many committed
  {
   "cell_type": "markdown",
   "source": [
    "Import the necessary modules. We need stationarity from ntrfc and some other modules for the definition of a signal and for rendering a plot"
   ],
   "metadata": {
    "collapsed": false,
    "pycharm": {
     "name": "#%% md\n"
    }
   }
  },
many's avatar
many committed
  {
   "cell_type": "code",
many's avatar
many committed
   "execution_count": 1,
many's avatar
many committed
   "outputs": [],
   "source": [
    "from ntrfc.timeseries.stationarity import stationarity\n",
    "import numpy as np\n",
    "import matplotlib.pyplot as plt"
   ],
   "metadata": {
    "collapsed": false,
    "pycharm": {
     "name": "#%%\n"
    }
   }
  },
many's avatar
many committed
  {
   "cell_type": "markdown",
   "source": [
    "Lets define a signal generator. We want to generate a signal with some noise, an initial transient and a constant deterministic fluctuation."
   ],
   "metadata": {
    "collapsed": false,
    "pycharm": {
     "name": "#%% md\n"
    }
   }
  },
many's avatar
many committed
  {
   "cell_type": "code",
many's avatar
many committed
   "execution_count": 2,
many's avatar
many committed
   "outputs": [],
   "source": [
    "def signalgen_abatingsine(amplitude, noiseamplitude, frequency, mean, abate, time):\n",
    "    resolution = 2048\n",
    "    step = (resolution * frequency ** -1) ** -1\n",
    "\n",
    "    times = np.arange(0, time, step)\n",
    "    noise = np.random.normal(-1, 1, len(times)) * noiseamplitude\n",
    "\n",
    "    values = amplitude * np.sin(frequency * (2 * np.pi) * times) + mean + np.e ** -(times * abate) + noise\n",
    "    return times, values"
   ],
   "metadata": {
    "collapsed": false,
    "pycharm": {
     "name": "#%%\n"
    }
   }
  },
many's avatar
many committed
  {
   "cell_type": "markdown",
   "source": [
    "Lets define the input arguments for the signal generator and lets generate a signal."
   ],
   "metadata": {
    "collapsed": false,
    "pycharm": {
     "name": "#%% md\n"
    }
   }
  },
many's avatar
many committed
  {
   "cell_type": "code",
many's avatar
many committed
   "execution_count": 3,
many's avatar
many committed
   "outputs": [],
   "source": [
    "test_amplitudes = 0.1\n",
    "test_noiseamplitude = 0.01\n",
    "test_frequencies = 6\n",
    "test_times = 40\n",
    "test_mean = -1\n",
    "test_abate = 1\n",
    "\n",
    "\n",
    "timesteps, values = signalgen_abatingsine(amplitude=test_amplitudes, noiseamplitude=test_noiseamplitude,\n",
    "                                          frequency=test_frequencies, mean=test_mean, time=test_times,\n",
    "                                          abate=test_abate)\n"
   ],
   "metadata": {
    "collapsed": false,
    "pycharm": {
     "name": "#%%\n"
    }
   }
  },
many's avatar
many committed
  {
   "cell_type": "markdown",
   "source": [
    "Lets compute the index of the stationary timestep."
   ],
   "metadata": {
    "collapsed": false,
    "pycharm": {
     "name": "#%% md\n"
    }
   }
  },
many's avatar
many committed
  {
   "cell_type": "code",
many's avatar
many committed
   "execution_count": 5,
many's avatar
many committed
   "outputs": [],
   "source": [
    "stationary_timestep = stationarity(timesteps, values)\n",
    "\n",
many's avatar
many committed
    "well_computed_stationarity_limit = -np.log(0.05) / test_abate\n",
many's avatar
many committed
    "well_computed_stationary_time = timesteps[-1] - well_computed_stationarity_limit\n",
    "stationary_time = timesteps[-1] - stationary_timestep\n"
   ],
   "metadata": {
    "collapsed": false,
    "pycharm": {
many's avatar
many committed
     "name": "#%%\n"
many's avatar
many committed
    }
   }
  },
many's avatar
many committed
  {
   "cell_type": "markdown",
   "source": [
    "Lets plot the result"
   ],
   "metadata": {
    "collapsed": false,
    "pycharm": {
     "name": "#%% md\n"
    }
   }
  },
many's avatar
many committed
  {
   "cell_type": "code",
many's avatar
many committed
   "execution_count": 6,
   "outputs": [
    {
     "data": {
      "text/plain": "<Figure size 640x480 with 1 Axes>",
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAi8AAAGdCAYAAADaPpOnAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjYuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/av/WaAAAACXBIWXMAAA9hAAAPYQGoP6dpAABJ3ElEQVR4nO3deVQT5/4G8CdhCYsQRHZFATfcF6wUq5YKVdRrtXaztb9qa7W12tXbW+1160r3xeqt9Xaxi1tt61KvpcW9KoIiqKCgIgqCgIoQFtmS+f2BRAIJEMhkCDyfc3IOTGYm38lkefK+78zIBEEQQERERGQh5FIXQERERGQMhhciIiKyKAwvREREZFEYXoiIiMiiMLwQERGRRWF4ISIiIovC8EJEREQWheGFiIiILIq11AWYmkajQXZ2NpycnCCTyaQuh4iIiJpAEAQUFRXBx8cHcnnDbSttLrxkZ2fD19dX6jKIiIioGTIzM9GlS5cG52lz4cXJyQlA9cY7OztLXI35lVSUwOdjHwBA9oJsONo6SlwRERFR41QqFXx9fbXf4w1pc+GlpqvI2dm5XYYXqworwK76b2dnZ4YXIiKyKE0Z8sEBu0RERGRRGF6IiIjIojC8EBERkUVheCEiIiKLwvBCREREFoXhhYiIiCwKwwsRERFZFIYXIiIisigML0RERGRRGF6IiIjIojC8EBERkUVheCEiIiKLwvBihPRrJfhqfxpuVqilLoWIiKjdanNXlRbTPR/tAwBcvF6KyKkDpC2GiIionWLLSzNsiMuQugQiIqJ2i+GFiIiILArDCxEREVkUhhciIiKyKAwvREREZFEYXoiIiMiiMLwQERGRRWF4ISIiIovC8EJEREQWheGliQRB0P5tLZdJWAkREVH7xvDSRGlXi7V/jx/gLWElRERE7RvDSxPJZLdbW+ys+bQRERFJxSzfwqtWrYKfnx/s7OwQHByMuLi4BuffvHkzAgMDYWdnhwEDBmDnzp3mKLNBtla3nyoHWysJKyEiImrfRA8vmzZtwiuvvIJly5bh+PHjGDRoEMaNG4e8vDy98x8+fBiPPvooZs2ahYSEBEyZMgVTpkxBUlKS2KU2qNaQF2gEw/MRERGRuEQPL5988glmz56NJ598En379sXq1avh4OCAb7/9Vu/8n3/+OSIiIvDqq6+iT58+eOuttzB06FCsXLlS7FKbTC0wvRAREUlF1PBSUVGB+Ph4hIeH335AuRzh4eGIiYnRu0xMTIzO/AAwbtw4g/OXl5dDpVLp3MQg4HZgERheiIiIJCNqeLl27RrUajU8PT11pnt6eiInJ0fvMjk5OUbNHxkZCaVSqb35+vqapvg6aucVNfuNiIiIJGPxh80sWrQIhYWF2ltmZqYoj1M7rjC7EBERScdazJW7ubnBysoKubm5OtNzc3Ph5eWldxkvLy+j5lcoFFAoFKYpuAG1u4o0TC9ERESSEbXlxdbWFkFBQdi9e7d2mkajwe7duxESEqJ3mZCQEJ35ASA6Otrg/FLQcMwLERGRZERteQGAV155BTNmzMCwYcMwfPhwfPbZZygpKcGTTz4JAHjiiSfQuXNnREZGAgBefPFF3H333fj4448xceJEbNy4EceOHcOaNWvELrVBteOKmtmFiIhIMqKHl0ceeQRXr17F0qVLkZOTg8GDByMqKko7KDcjIwNy+e0GoBEjRmD9+vVYvHgxXn/9dfTs2RNbt25F//79xS61QV7OdhjV0w1/n7vGlhciIiIJyYQ2dtyvSqWCUqlEYWEhnJ2dTbruH2MuYsm2ZIzv74UvHw8y6bpNpaSiBB0iOwAAihcVw9HWUeKKiIiIGmfM97fFH21kTvJbV5NmywsREZF0GF6MIL91cUa1RuJCiIiI2jGGFyNYydjyQkREJDWGFyPcyi4ML0RERBJieDGClbym24jhhYiISCoML0aoCS9seCEiIpIOw4sRZDK2vBAREUmN4cUIHLBLREQkPYYXI8g5YJeIiEhyDC9GuH2SOokLISIiascYXowg55gXIiIiyTG8GMHq1rPVxi4HRUREZFEYXoygPdqI4YWIiEgyDC9GsOK1jYiIiCTH8GKEmjEv7DYiIiKSDsOLEeS3ni0O2CUiIpIOw4sReJI6IiIi6TG8GKHmPC9pV0vw/IYEtsAQERFJgOHFCDVjXgDg9xPZ2JOSJ2E1RERE7RPDixHkMt3/yyrV0hRCRETUjjG8GMGqTnq5cLVEokqIiIjaL4YXI9TuNgKAT3edlagSIiKi9ovhxQh1w4t13X4kIiIiEh3DixHkdZ4tGys+fURERObGb18jWNVpeXnyLj9pCiEiImrHGF6MIKsTXtydFBJVQkRE1H4xvBih7tFGHPFCRERkfgwvRqjbbcTz6xIREZkfw4sR6mQX8OoARERE5sfwYgR5nW4jgRdoJCIiMjuGFyPUPa0LL8xIRERkfgwvRqg75kXNlhciIiKzY3gxQt1DpTVseSEiIjI7hhcj1O02YnYhIiIyP4YXI9S9thHHvBAREZkfw4sR6oaXvKJyiSohIiJqvxhejCCr82xtiMuQphAiIqJ2jOHFCHVbXoiIiMj8GF6MUPdQaSIiIjI/hhcj1M0u4/p5SlMIERFRO8bwYoS63UYdFDYSVUJERNR+MbwYoe55XnhtIyIiIvNjeDFC3ZYXDcMLERGR2TG8GKHumBeeo46IiMj8GF6MUO/aRmx5ISIiMjtRw0t+fj6mT58OZ2dnuLi4YNasWSguLm5w/ueffx69e/eGvb09unbtihdeeAGFhYViltlszC5ERETmJ2p4mT59OpKTkxEdHY0dO3bgwIEDmDNnjsH5s7OzkZ2djY8++ghJSUlYu3YtoqKiMGvWLDHLbDa2vBAREZmftVgrPnPmDKKionD06FEMGzYMAPDFF19gwoQJ+Oijj+Dj41Nvmf79++PXX3/V/t+9e3e88847ePzxx1FVVQVra9HKbRaGFyIiIvMTreUlJiYGLi4u2uACAOHh4ZDL5YiNjW3yegoLC+Hs7GwwuJSXl0OlUunczIUDdomIiMxPtPCSk5MDDw8PnWnW1tZwdXVFTk5Ok9Zx7do1vPXWWw12NUVGRkKpVGpvvr6+LarbGDzPCxERkfkZHV4WLlwImUzW4C0lJaXFhalUKkycOBF9+/bF8uXLDc63aNEiFBYWam+ZmZktfuymYssLERGR+Rk9iGTBggWYOXNmg/MEBATAy8sLeXl5OtOrqqqQn58PLy+vBpcvKipCREQEnJycsGXLFtjYGD4Nv0KhgEKhaHL9psQxL0REROZndHhxd3eHu7t7o/OFhISgoKAA8fHxCAoKAgDs2bMHGo0GwcHBBpdTqVQYN24cFAoFtm/fDjs7O2NLNBu2vBAREZmfaGNe+vTpg4iICMyePRtxcXE4dOgQ5s+fj2nTpmmPNMrKykJgYCDi4uIAVAeXsWPHoqSkBN988w1UKhVycnKQk5MDtVotVqnNduDsVfx2/LLUZRAREbUroh57vG7dOsyfPx9hYWGQy+V44IEHsGLFCu39lZWVSE1NRWlpKQDg+PHj2iORevToobOu9PR0+Pn5iVlus7zy8wmMCfSAi4Ot1KUQERG1C6KGF1dXV6xfv97g/X5+fjpH7ISGhlrkETxllRqpSyAiImo3eG0jE3jrf6elLoGIiKjdYHgxgT9OXZG6BCIionaD4cUELK+ji4iIyHIxvJiABQ7TISIislgMLyagtDd8Ej0iIiIyLYYXE+jp0UHqEoiIiNoNhhcTsLe1kroEIiKidoPhxQTuH9JZ6hKIiIjaDYYXE7CSy6QugYiIqN1geDHSkn/0rTdNzSs0EhERmQ3Di5FmjfTHD08N15nG8EJERGQ+DC/NYG2l202k4YleiIiIzIbhpRnksrrhRaJCiIiI2iGGl2aoG17YbURERGQ+DC/NUPfgInYbERERmQ/DSzPI2PJCREQkGYaXZqjb8sLwQkREZD4ML81gLdd92thtREREZD4ML81Q/1BpiQohIiJqhxhemsHGimNeiIiIpMLw0gxWdbqNBHYbERERmQ3DSzNYy9ltREREJBWGl2awsdJ92thtREREZD4ML81gVaflhd1GRERE5sPw0gx1B+zml1ZAVVYpUTVERETtC8NLM1jX6Tb66UgGBi7/Cxp2HxEREYmO4aUZ6p5ht0alRmPeQoiIiNohhpdmqHtV6RpVara8EBERiY3hxYQYXoiIiMTH8NIMhlpe2G1EREQkPoaXZjCQXVCpZnghIiISG8NLM3DMCxERkXQYXprBQMMLW16IiIjMgOGlGQx1G1XxPC9ERESiY3hpBhm7jYiIiCTD8GJCvEAjERGR+BheTKiKh0oTERGJjuHFhNjyQkREJD6GFxPigF0iIiLxMbyYEFteiIiIxMfwYkJseSEiIhIfw4sJaQSGFyIiIrExvJiQwPBCREQkOoaXZjr42j2YPcpfZxqvDkBERCQ+hpdm6tLRAf18lDrT2G1EREQkPlHDS35+PqZPnw5nZ2e4uLhg1qxZKC4ubtKygiBg/PjxkMlk2Lp1q5hlNptcrnuZAA0H7BIREYlO1PAyffp0JCcnIzo6Gjt27MCBAwcwZ86cJi372WefGbyGUGtRJ7uA2YWIiEh81mKt+MyZM4iKisLRo0cxbNgwAMAXX3yBCRMm4KOPPoKPj4/BZRMTE/Hxxx/j2LFj8Pb2FqvEFpPXCVdqdhsRERGJTrSWl5iYGLi4uGiDCwCEh4dDLpcjNjbW4HKlpaV47LHHsGrVKnh5eTX6OOXl5VCpVDo3c6nb8sKjjYiIiMQnWnjJycmBh4eHzjRra2u4uroiJyfH4HIvv/wyRowYgcmTJzfpcSIjI6FUKrU3X1/fFtVtjLrdWpeul6K8Sm22xyciImqPjA4vCxcuhEwma/CWkpLSrGK2b9+OPXv24LPPPmvyMosWLUJhYaH2lpmZ2azHbo663UafRJ/F5JWHzPb4RERE7ZHRY14WLFiAmTNnNjhPQEAAvLy8kJeXpzO9qqoK+fn5BruD9uzZg7S0NLi4uOhMf+CBBzBq1Cjs27ev3jIKhQIKhcKYTTAZKz3RLyWnyPyFEBERtSNGhxd3d3e4u7s3Ol9ISAgKCgoQHx+PoKAgANXhRKPRIDg4WO8yCxcuxNNPP60zbcCAAfj0008xadIkY0sVXWs/GoqIiKgtEu1ooz59+iAiIgKzZ8/G6tWrUVlZifnz52PatGnaI42ysrIQFhaGH374AcOHD4eXl5feVpmuXbvC39+/3nSp1e02qiEIAoMNERGRSEQ9z8u6desQGBiIsLAwTJgwASNHjsSaNWu091dWViI1NRWlpaViliGaukcb1bhaXG7eQoiIiNoR0VpeAMDV1RXr1683eL+fn1+jhxe35sOPDbW8WMt51QUiIiKx8Fu2BdgzREREZH4MLy1gZSC9qHmdACIiItEwvLRA3Qsz1mjNXV1ERESWjuGlBQwN2GXDCxERkXgYXlrAysDAXF6gkYiISDwMLy1gaMyLhk0vREREomF4aQErg2NezFwIERFRO8Lw0gKGwgu7jYiIiMTD8NIChsKLhuGFiIhINAwvLWC424jhhYiISCwMLy1g+CR1Zi6EiIioHWF4aQErK3YbERERmRvDSwsYPFSa4YWIiEg0DC8tYPAMu+w2IiIiEg3DSwvI2PJCRERkdgwvLWCo5eW345dxPq/YvMUQERG1EwwvLWCo5eX7mEsI/2S/mashIiJqHxheWsDOhk8fERGRufHbtwUcbK3xwYMDEejlJHUpRERE7QbDSws9PMwX04O7Sl0GERFRu8HwYgKGxr4QERGR6TG8mICc4YWIiMhsGF5MwIrPIhERkdnwa9cE2G1ERERkPgwvJsBuIyIiIvNheDEBQ2faJSIiItNjeDEBK6YXIiIis2F4MQGOeSEiIjIfhhcTMNTwIvDq0kRERCbH8GICVgZaXtQahhciIiJTY3gxAUPdRmq2vBAREZkcw4sJGO42Mm8dRERE7QHDiwkYOtqI3UZERESmx/BiAoZOUsduIyIiItNjeDEBhbX+p/HLfWnIuF5q5mqIiIjaNoYXE7haXK53+pf70jD6w71mroaIiKhtY3gxgSo1u4eIiIjMheHFBOxtraQugYiIqN1geDEBQ2NeiIiIyPT4rWsCho42IiIiItNjeDEBZhciIiLzYXgxAba8EBERmQ/DiwkYOsMuERERmR7Diwmw4YWIiMh8GF5MoLFuoz+Tc1BWqTZTNURERG2baOElPz8f06dPh7OzM1xcXDBr1iwUFxc3ulxMTAzGjBkDR0dHODs7Y/To0bh586ZYZZpEY+HlmR/jsWxbspmqISIiattECy/Tp09HcnIyoqOjsWPHDhw4cABz5sxpcJmYmBhERERg7NixiIuLw9GjRzF//nzI5a27gcivk0Oj82w6lmmGSoiIiNo+mSCY/tLHZ86cQd++fXH06FEMGzYMABAVFYUJEybg8uXL8PHx0bvcnXfeiXvvvRdvvfVWsx9bpVJBqVSisLAQzs7OzV6PsU5eLsCHf6bi73PXDM5z8b2JotdRUlGCDpEdAADFi4rhaOso+mMSERG1lDHf36I0acTExMDFxUUbXAAgPDwccrkcsbGxepfJy8tDbGwsPDw8MGLECHh6euLuu+/GwYMHG3ys8vJyqFQqnZsUBnZxQXf3DpI8NhERUXsiSnjJycmBh4eHzjRra2u4uroiJydH7zIXLlwAACxfvhyzZ89GVFQUhg4dirCwMJw7d87gY0VGRkKpVGpvvr6+ptsQI/GQaSIiIvEZFV4WLlwImUzW4C0lJaVZhWg0GgDAM888gyeffBJDhgzBp59+it69e+Pbb781uNyiRYtQWFiovWVmSje2hEcUERERic/amJkXLFiAmTNnNjhPQEAAvLy8kJeXpzO9qqoK+fn58PLy0ruct7c3AKBv37460/v06YOMjAyDj6dQKKBQKJpQvfjWxRquk4iIiEzDqPDi7u4Od3f3RucLCQlBQUEB4uPjERQUBADYs2cPNBoNgoOD9S7j5+cHHx8fpKam6kw/e/Ysxo8fb0yZrRJPZEdERGQaoox56dOnDyIiIjB79mzExcXh0KFDmD9/PqZNm6Y90igrKwuBgYGIi4sDAMhkMrz66qtYsWIFfvnlF5w/fx5LlixBSkoKZs2aJUaZREREZIGMankxxrp16zB//nyEhYVBLpfjgQcewIoVK7T3V1ZWIjU1FaWlpdppL730EsrKyvDyyy8jPz8fgwYNQnR0NLp37y5WmWbDhhciIiLTEC28uLq6Yv369Qbv9/Pzg75TzCxcuBALFy4UqyzJyNhvREREZBKt+9S1bQijCxERkWkwvJhJlcbkJzImIiJqlxheiIiIyKIwvJiQNc+wS0REJDqGFxN6Iayn1CUQERG1eQwvJmRvYyV1CURERG0ew4sJHblwvcH7H/4qBhnXSxuch4iIiBrG8GJCaj3nraktLj0f/9x8wkzVEBERtU0MLyZ03yCfRue5XlJuhkqIiIjaLoYXE7JqwtFGaVdLzFAJERFR28XwYkJyXgKAiIhIdAwvJtS/s1LqEoiIiNo8hhcT8ndzxLZ5d0ldBhERUZvG8GJig3xd8OjwrlKXQURE1GYxvIiAQ1+IiIjEw/AiAl7iiIiISDwMLyJo6KgjV0dbFJdXmbEaIiKitoXhRQQNNbzkl1Sg/7I/8dORS2arh4iIqC1heBFBwxcJqLZ4a5LodRAREbVFDC9ERERkURheRNDI9RmJiIioBRheRCA0qeOIiIiImoPhRQRseSEiIhIPw4sINh+7LHUJREREbRbDiwgq1BqpSyAiImqzGF6IiIjIojC8SGjGt3GoZCsNERGRURheJLT/7FVsOpopdRlEREQWheFFYku38Uy7RERExmB4kZiGh1UTEREZheGFiIiILArDCxEREVkUhhcR7HxhlFHzbz7GQbtERERNxfAigr4+zrj43kT4uto3af5XfzkJNQe/EBERNQnDSytxs1ItdQlEREQWgeFFRMZcoDHlikq8QoiIiNoQhhcRXS0qb/K89rZWIlZCRETUdjC8iKi8qumn/jemlYaIiKg9Y3hpJSavOoTl25OlLoOIiKjVY3hpJdQaAWsPX5S6DCIiolaP4YWIiIgsCsMLERERWRSGFyIiIrIoDC+tzPYT2TiRWSB1GURERK0Ww4uIHJpx7pYXNiRg8qpDIlRDRETUNogWXvLz8zF9+nQ4OzvDxcUFs2bNQnFxcYPL5OTk4P/+7//g5eUFR0dHDB06FL/++qtYJRIREZEFEi28TJ8+HcnJyYiOjsaOHTtw4MABzJkzp8FlnnjiCaSmpmL79u04deoUpk6diocffhgJCQlilSkqmdQFEBERtUGihJczZ84gKioKX3/9NYKDgzFy5Eh88cUX2LhxI7Kzsw0ud/jwYTz//PMYPnw4AgICsHjxYri4uCA+Pl6MMomIiMgCiRJeYmJi4OLigmHDhmmnhYeHQy6XIzY21uByI0aMwKZNm5Cfnw+NRoONGzeirKwMoaGhBpcpLy+HSqXSubUWLTnjf/TpXOSpykxWCxERUVshSnjJycmBh4eHzjRra2u4uroiJyfH4HI///wzKisr0alTJygUCjzzzDPYsmULevToYXCZyMhIKJVK7c3X19dk29FSM0f4NXvZ2T8cQ/gn+01XDBERURthVHhZuHAhZDJZg7eUlJRmF7NkyRIUFBRg165dOHbsGF555RU8/PDDOHXqlMFlFi1ahMLCQu0tMzOz2Y9vai+F98KPs4Y3e3lVWZUJqyEiImobrI2ZecGCBZg5c2aD8wQEBMDLywt5eXk606uqqpCfnw8vLy+9y6WlpWHlypVISkpCv379AACDBg3C33//jVWrVmH16tV6l1MoFFAoFMZshtnYWssxqqe71GUQERG1KUaFF3d3d7i7N/5lHBISgoKCAsTHxyMoKAgAsGfPHmg0GgQHB+tdprS0FAAgl+s2BllZWUGj0RhTZpsy7tMD+M/jQ9HdvYPUpRAREbUKoox56dOnDyIiIjB79mzExcXh0KFDmD9/PqZNmwYfHx8AQFZWFgIDAxEXFwcACAwMRI8ePfDMM88gLi4OaWlp+PjjjxEdHY0pU6aIUaZFSM0twtyfeLQVERFRDdHO87Ju3ToEBgYiLCwMEyZMwMiRI7FmzRrt/ZWVlUhNTdW2uNjY2GDnzp1wd3fHpEmTMHDgQPzwww/4/vvvMWHCBLHKtAi5qnKpSyAiImo1jOo2MoarqyvWr19v8H4/Pz8Igu7BxD179uQZdfW4WamWugQiIqJWg9c2sgAVVe13zA8REVFdDC8Wpm5rFRERUXsjWrcRmdZj/z0CT2c7/H3uGn6bOwJdOzlIXRIREZEk2PJiIQ6nXceWhCxcKy7Hh3+lSl0OERGRZBheLBC7joiIqD1jeLFA53KLpS6BiIhIMgwvFsjHxU7qEoiIiCTD8GIGYwI9Gp/JCFeLy/G/k1dQUl6FwtJKk66biIiotePRRmbw+bTB2H0mDy9tSjTJ+pKyVJi3/rj2/+Q3xsFRwV1JRETtA1tezMDJzgZThnQWbf3p10pEWzcREVFrw/DSBhxOu4ajF/OlLoOIiMgs2NfQBry7MwUAkPbuBMSmX5e4GiIiInExvLQhk744iKQreYC91JUQtT1llWpYy2WwtmKDNZHU+C5sQ05fUdWbVqXW4NL11jsmpqJKw5PuUatXUl6FwCVRuPfTA1KXQkRgeGnznv7hGO7+cB/+OHVF7/3/PXABaw6kYfznf+PHI5fMWlt+SQX6L/8Ts384ZtbHJTJWQkYBAA6OJ8tTXqXGpqMZyCq4KXUpJsXw0obd/59D2Jd6FQDw6a6zOHLhOsoq1SgorQAApF0txjs7z+DdnSk4c0WFJVuTzFrf7yeyUVGlwa4zeXrvr1RrcDzjBirVGrPWBZj+EgxRSTl4cWMCSiuqTLpeMcSkXcestUdx+Uap1KW0GgKa9nrIzC/FzlNXWk1r4r7UPCRk3JC6DJLQf/am4bVfT+HeT/ZLXYpJMby0YSk5Rdq/z+YWY9qaIwhcEoXBb0bj2MV8rNp7XsLqAE0jH/Bv/J6Mqf85jLd2nG5wvioTh5vD569h0Bt/YfuJbJ3pKTkqHD5/rVnrfPaneGxLzMbq/RdMUaJJXb5Rind3ntH+Mnv0v0ewOyUPr/x8Qu/8FVUaSQIlAFwvLseYj/Zh5Z5zDc63cs85jPl4H/JLKpr9WFFJV/BJ9Fm9QeRmhRpf7D6HM3W6akd9sBfPrTte77UjhayCm5j53VHc/5/DZn/smxVqiwjq7cGBc9U/YEsr1EYtV1RWCY2mdYRwfRhezOi350YgvI9pz7bbXA+ujsFvx7OMWkYQhCZ/ICVmFuC5dfHIzL/9671KrUFmfili0q7fWl/D6/jpSAYA4IeY6u6spKxCrD2UrvOGyi+pwJA3o/HKrRMAllWqsS72UoOtBhqNgDk/HDMYip74Ng6qsiq8sCEB62MztNsc8dnfeOzrWFxsQdfBjpO6X2r5JRUmD18AsDUhC498FYNrxeUAqved2sAH0RPfxGHNgQt48rs4nenZt8LM94cv4q0dpyEIAqrUGgx7Oxohkbub9MG289QVfHcovYVbc9tXBy7gwrUSfPTX2Qbn++ivs7hwtQRf7U9r8rqTswvx5b40VFRV749nfzqOFbvPYf/Zqzrz/X4iG4u3JuHj6LMY//nfAIBcVRl+PpqpnSc2XfpTF2TdaLibIL+kAp/vOtfiFra64U6tETBg+Z/ou/RPLPz1ZIt+JKnKKnHfyoNYc0D/fhQEAYU3K+tN++fmE3h35xmjH6+5LWbXi8txNreo8RmboKxSjfIq44KGITWvZWOlXS3GgOV/Ydb3R01ShxgYXsxoaNeOWPnYUKnLaNDvJ7Jx9GI+3vz9NM5cUel8QS34+QT6Lv0T5269SZdvT8biradQWFpZ700yZdUh7DyVg/m3zgT8z80nEPT2Loz6YC8e/e8RnMgsaGJDfLUjF67jH18cxPLfT2Nr4u3QtfFoBorKq/BbQvW0z3efw7+3JCHis78NruvE5QL8dToX3xzU/VItLq/C139fQFWtbX59yyks2ZqMXFWZdtrKvefxY8xFnWUFQcDJywUoKa/SfvDo+3K/cLU6+Px05BI+/DMFQ9+KxrQ1RwzWWvvDNDm7UPvc11VYWon/+yYWn+06i+8OpeOlTYmITc/Hh1GpAID56xMw8v09OuFTEASczS3ChVth7GxuMe7/z6F66162PRnfHExHYmYBcovKoSqrwrXiClwvqdCp7/qtoFSl1iAlR4X3o1Lw3LrjeOP300jJUWkfs/Yy6ddKGv2A1WgEbUuPsS0+VU0IWGWVavyZnIOJKw7i/agUBL+7C1NrPQ/5JRU6Qfv5DQn49fhlnXXct/Ig/vXrSe3/ctnt+wRBwJKtSfj674Zb3b47lI6QyN1Iu1p94dXfjl/GtkT9PzDKq9RYsfscTl0uNLg+faFifWwGzudVr3/SFwfx6a6zGPn+3gbrasjGuAzcGblbpwWqpKJK+7xvPJqJD/9MNVhTY749mI6Tlwu1p4OoodYISMoqxBu/n8agN/7Ce3+kYPX+NAiCgLSrxfgl/jLWHNB9vjUaAX+fu6q9pEpM2nVMWxOD83nV76m3d5zG3R/uqxeGmiLo7V0Y++kB7XNbW35JBZZtS0JSlu6+EgRB24Vfo1KtwcDlfyHorV0Gf2w01cVrJQhc8od2vJYx1sdW/3Dcm3q1kTmlw0Olzcy69qdaK/T8hgTt398eSoedjRxfPh6Eu7q7aQPCvZ8egI/SDtmF1V/oPx3JgI/SDocXhdVbX1K2CgkZN/BLvO6H/YnLBQY/yFJzivDsT/E602p/wX/811nIZMCYQE+U1WkK/ftWE2lxeRVe2ZSIZ0O7o5enk848lWr9j/v2jtPYWOvXc41fj1/W+bL6Jf4yfom/jNDeHvB1dQAA7DyVo3PJhrt6dMLJy4X4+olh8HDWvZBm/KV8LK41vujYJf1jEtYcSMPXf6dj87Mh6NRBgYkrDgIAVj02FDZWMgzu6oJOjgpYyWX4Ys85/H3uGv4+p9utVXCzOmD879aA7V1n8jBpoDdkMhn+SMrBc+uO68zf0Afd6Ssq5BTeDnF3vLMLEwd6w0omQ66qDLHp+RjYRYnkbFW9D97rxRXQaARMW3MEChs5fnhqOHadycPsH45hSFcXPHWXP0b3dIfSwQaCICBXVQ4vZfXzdv9/DiGroAyHFt5jsDZDat5t+1LzUFapQUR/r3rzLN+erLPfb5RW4kat58GqkffstDUxyFWV60yzkt1eJjGzQDsY/ulRAdrpF6+VYOn2ZNwZ4IqD567h8K0WyTd+P40V0wZru+3+Ss5F5472WHPgAmaEdMMbk/vj67/T8Un0WXwSfRYX35tYr6brxeU4V+uL9J3/nUZcej5O3Ao7F9+baHAAZ03olt/a7o/+TIUAAb29nNHTowP6eDsjLj0fb+5IRlJWdWh57deT2D5/JM7mFmHHyfoHB3z0ZyquFpVj07FMBHo5Ieql0dr7opJy8POxTHz44EB06qDQWa5cT7AtLK3Ef/adx1e1wsnqWy1sfp0cdeYVBAGyW/tiXewlLNmWjAA3R+z5Zyge/W/1Z8q0NUdwbPG9+PrWj5lNRzPw1F3+eObHeAR6O+Fo+g2EBrrjudAeSL9Wgv/sPY+5od0R4N6hXm3HL91AD4/b049dzMf7USk4evEGvo+5hDcn90NYH09k3biJrYlZWB+bgbem9Eewvyt6eTohV1WGCrUGFWoNgt6Oxvj+XoicOlC7vpQcFdbHZuD5MT3h7qTAtwfT8c3BdIR074Tsgpv47sk7YCOX480dp/HTkUuom3/+9csJDPbtiH2peVjyj76o0gjwd7v9nAmCgE1HM3Eis0A7bfLKg5g9OgD/GOiD+Ev5iD6dh5fCe8LOxqre9psTw4uZWdo5IsoqNXjyu/pNh9m1vsRq/t+akIU7/F3R2eX2iWbUGkFvn/vSbcnQ952Qfq0E89cfb/CojqyCm3h5k/7xGLX9lpCFP5NzkPxmBK4WlSO74CYG+bpAVudxswpu4vNdZ/Hzscv6V2TA61tO4dHhXTFhgLdOaxAAHDpf/UX0iJ5WlQe+jNG7Po1GwKmsQvTxdoattVz7azNyZwqWTuqrna92SHJxsEFBAxfn/DM5F6M+uP3L+oUNCfg1/jK+f2o4NsRlNLh9dVtE/r2l/oDu/9X5ojppoCVAIwjIKriJuFtngi6tUGu/0BMyCvB8RgIGdlFi6T/64tfjl7EhLhMR/bww0Fep/cK9893duFFrW5dvT8boXm7wd+uAbq4O2i/b2irVGly+UYqZt17DxxaHw62DAuVValjL5dh1JldvYK3txY2JDd5/5EL9LqKD569hQ1wGHh3eFXlFt4NN+rUSXLpegqXbkpFxq0v1QJ1uqQNnr+LNWl2a/6t1pOD3MZdQWqHG5nj9r9WVe85hQ1xmvWDy37+b1nUnCAImrTyIKrWAP14chaLyKqys0+1z/p3xePgr3ddwzX4fa+BQ8trrqD0WD4D2h0rQ27swI6QbnO1tcDa3CF9OD0LtPZpxvRR7U/OwbHuywfovXS9BSfnt1sXjGQXo7eWEDgpr7Y+vC3U+W64VV+h8WctlMkSfzsXulDzsTqk+mCDuYj6eC+2Bf24+gfhLN/D7yWz8NCsYQd06asMRAJ3X4LGL+Xhwte7ztHRbMpZu062/5kCJKYN9kJx9uwWroLQSG+IyEX/pBh6/sxvu6e2BR746gsKblUi7Wox1T9+pfZ3U/DjceeoKlPY2WHv4ot7n5+djl7Wfc3+dzgUArJ8djBHd3bAu9hKiknLq/QA6cbkQ89cnYNXeNG0Lm72NFV4M76n3McyF4YVMpubCk+F9PJs0f+1fBV/uS0NeURm+O3Sx2Y9fVqnW/hKsUXKrZeaOd3YBADY/G6Jz//tRKfhyX9PHRdRW09Ix/54eSKz14dccvRf/gQD3DjhzRYW+3s4I7e2uvU+AAGsr/b/+GwouNS7XGfuw/+xV+C38X6PL1f7Sban/+0Z3PE2/ZX/Wm+fk5UKdD/uo5BxEJedo/79RZ1vXHr6o/ZDu7GIPhY0cHz44UNsaBlR/2X8fc/sUAKezVVCVVWL++gSIKe1qCRb9dgpnc4t0XtP3fLSvScs3NB5NX3C5XlyOLQlZjY4FqlF3HNKCn0/gvsE+2JuSp/0CXbwtCWo9rZQ9/v1Hkx6jIfd8tA9jAj1ws1K35bT2vloXewlFZbeDyOgPG+/e0gi6XYUPfFn9w+nrJ4Y12Ko4edXtbsK3/2d4rEz8rVbSskoNHlwdg4eCuuDDhwZp76/dsl7TktZUWxP1D/I+m1t8K/DcDj2JGQUoq6w/LqaiSoP8EuO6veatO44xgZ71ukLrqt01WNMNLCWZ0FqO6TMRlUoFpVKJwsJCODs7S12OXk354mguDcqQaf8gAMD35i+Qw66RJdq+Lh3t632BW5IJA7ywbFI/BL+72+yPPa6fJ/5MzjX741LTvTCmB1bskfbIQQBwUlijqFzaI4w8nRX1uvD0mT3Kv8mtUY35z/ShOt2vI3u4ITm7sF7YNjVba3m91lEXBxsEejnpbQ00pfA+nvh6xjCTr9eY72+2vFCbZ8nBBageT6O0t5HksRlcWr/WEFwASB5cADQpuABN70Zrirrjxg4283QKxtI30L2gtFL04AIAVRppTpVQm2UNwCBqpzbENTwug4jIXKoMHPRgTgwvRERE1GQtPYzbFBheJDCqp5vUJRARETVLY2dHNweGFwn8OCsYT4/0l7oMIiIiozG8tGNj+9U/WRYREVFrx26jdiy7jV2enIiI2odWkF0YXqTS2CnHiYiIWqOWnpTTFBheJNLXp3WeQI+IiKi1Y3iRSJeO9o3PRERERPUwvEjERs6nnoiIqDn4DSoRfVfAJSIiau0cbK2kLoHhhYiIiJquUwdbqUtgeCEiIiLLwvAioXfu74+RPXipACIiImMwvEhoenA3vD6hj9RlEBERNVkruDoAw4vUeLI6IiIi4zC8SKy3l5PUJRARETUZW14IAPD7/JFwbAWHnhEREVkChpdWYEAXJT6bNkTqMoiIiCwCw0srwZEvRERkCYRW0G/E8NJKSP9SICIisgyihZd33nkHI0aMgIODA1xcXJq0jCAIWLp0Kby9vWFvb4/w8HCcO3dOrBJbldaQZImIiBrTGr6tRAsvFRUVeOihhzB37twmL/PBBx9gxYoVWL16NWJjY+Ho6Ihx48ahrKxMrDJbDU1reDUQERFZAGuxVvzGG28AANauXduk+QVBwGeffYbFixdj8uTJAIAffvgBnp6e2Lp1K6ZNmyZWqa2Cn5uD1CUQERE1qjV0FLSaMS/p6enIyclBeHi4dppSqURwcDBiYmIMLldeXg6VSqVzs0SBXs5YPqmv1GUQERE1SNMK0kurCS85OTkAAE9PT53pnp6e2vv0iYyMhFKp1N58fX1FrVNMM+/yl7oEIiKiBllceFm4cCFkMlmDt5SUFLFq1WvRokUoLCzU3jIzM836+KY23N9V6hKIiIgMulZcgec3JEhag1FjXhYsWICZM2c2OE9AQECzCvHy8gIA5ObmwtvbWzs9NzcXgwcPNricQqGAQqFo1mO2RoFeTohLz5e6DCIiIoNOZxdK+vhGhRd3d3e4u7uLUoi/vz+8vLywe/dubVhRqVSIjY016oglS5eZXyp1CURERA2S+qLCoo15ycjIQGJiIjIyMqBWq5GYmIjExEQUFxdr5wkMDMSWLVsAADKZDC+99BLefvttbN++HadOncITTzwBHx8fTJkyRawyWx1VWZXUJRARETVILpM2vIh2qPTSpUvx/fffa/8fMqT62j179+5FaGgoACA1NRWFhbebnv71r3+hpKQEc+bMQUFBAUaOHImoqCjY2dmJVWar09HBVuoSiIiIGiSTOLzIhDZ2aleVSgWlUonCwkI4OztLXY7Rjl7Mx0OrDR8a3hgNypBp/yAAwPfmL5Cj/QQ/IiIyj34+zvjfC6NMuk5jvr9bzaHSVO0OP1ccX3Kv1GUQEREZJHHDC8NLa+TqaIuXw3tJXQYREZFeUvfZMLy0Us/d0x1vTu4ndRlERET1MLyQXjZWcjwR4id1GURERPVIPViW4YWIiIiMIvWxPgwvrVxHBxupSyAiItLBbiMiIiKyKILEHUcML62c1CcCIiIiqostL9SgTx4eJHUJREREOjhglxoU2tsDqW9HYET3TgCAXp4dJK6IiIjaO6kH7Ip2bSMyHYW1FdbPvhOCIEAjAFk3bmL0h3ulLouIiNopqVteGF4siEwmg5UM8HBWSF0KERGRZNhtZIGs5BzES0REEuKAXTKWFY9AIiIiCUndbcTwYoHkbHkhIiIJST1gl+GFiIiIjHK1qFzSx2d4sVCLxgfiqbv8pS6DiIjaoZIKtaSPz6ONLNQzd3cHAMwN7Y7Pdp3FutgMiSsiIiIyD7a8WDh3JwV6eTpJXQYREZHZMLy0ASG3zr4LALbW3KVERCSujg42kj4+u43agF6eToh+eTTcnRSQycrh8sHt+w4tHIMnv4vD2dxi6QokIqI2JcBd2kvV8Gd6G9HT0wkuDrawttLdpZ1d7PHjrGAsm9RXosqIiKitkfqEHWx5aQc8ne3w5F3+GNHdDefzipGjKsPR9HxEJefozDftDl9sPJopUZVERGQppD5JHcNLO9Lbywm9vaoH984a6Y+9KXl4cu1R7f0LxvbGc6E9eNFHolaqs4s9AtwdoREEHDp/XepyiCTD8NKO3RPogZWPDcH89QkAALkM6NrJAalvR0CjAQ6ev4YR3TvhZqUaw97eJXG14glwd8SFqyVSl9HqPHt3d6zenyZ1GVTL/ldDtV3Dfgv/J3E17VcPjw44n9e+xxFK3W3EMS9t2EPDujQ6z3B/V+3f8lvXTFJYW8He1gr39vWEo8Iabh1a31Wse3s64bsn7zDJuvYsCDXJelqjxRP7NHvZheMDTVhJ4yYO8BZlvV9OH4rPpw02ermEJfeavhg9Pn1kEMICPZo0b+2Lso7v7wUAGNatI+L+HdbiOrq7O7Z4HbXNHOFnsnW1pqMoh3R1wa5X7sbOF0bhtYiWv0eeC+3e6Dxxr4ch7vUwvDqud4sfr61oPa8IMrk3JvVrdJ7al6eQN+OCj0r7hg+XG9ato8H7nrrLXycYXXxvYoPr6uxir/27n48zenrcHu3u62qPp0f64+XwXo2VjJfDe+GThwfpTHPRc9jfoC5Knf/fmzrA4Drv8NPdzmfuDtA7X0Q/L+3fjw73bfQL28HWCsP9XBucxxC/Tg4NnoV51kh/2NkY9xHwy7MhWD87uMnzOykMN+6+MKYH3J2q9/+onm5YNX2owXnXPnkH5oy+/ZyufGwIdjw/Ehtm34mVjw2BrZXh7dA0s3O+o6Nto/P07+zc6Dz/ntAHvTzrH5kxcYA3fp07AvcP6YL/PjEMhxaOwaPDuwIAts67C9vm3YWL703E908N1y4jq/Uefe+BgYicOgDfzLgDHk52+HbmsAafh7l1viRr1z6kqwt21wnxz4/p0ei26fPocF/seH4klt/XDy+F92zWOuo6uWwsNj8bove+0N7uSH07AjGLxuCXWvP88NRw2NtY1Zt/Rki3etNCAjrh5PKxjdax+vEgrJ1ZvT/6+jhjbmh3DOnq0uAyXs52GNBZafD+V8f1bvSwYw9nO3g422HePU3fJ4N8G66rOWq/jqW+PjC7jdowWRNeXa63PqCt5DJ0sGv6yyGoW0cs/UdfDOyixPeHL2L576cBAO/ePwCDfJXIunETwf6doHSw0WneHtLVBQkZBQCAp0f5Y+H4QLy+5VSDvzxdHW3x06xg9PZyQvfXdwIA/N0c0aWjA36cNRwdHWzRv9aHw6e7zuos37+zM5KyVNr/HRVWGNnTDQC0H27OdjYoKK0EUB08ZoT4wcfFXqf2h4f5YuFvp/TWuPnZETrzLhrfB1/tv6D9v5+PMzY/G4KScrV2oPSr4wLh6miLjyvVqFRr8POxy3hrx2ntMnv/GQq/Tg6QyWTada96bCi2JmYh+nSuweerxmPBXRv8gFnyj754Iawn1h66WO85++DBgTr/9/LsgAkDvBHUrSNi0uqPtejsYo+sgpv1ph9fei96/vsPvY8/f0xPzA3tgSMXruucq0jfOvt4O8NKLsOaA9XP6T8G+ujM+4+BPtrnyFtph32vhqL34qhbyzrhxOUCvTUY4uFUv7Wxu7sj3rivP7yUdgAEuDvZQWlvg6tF5fgj6QqWbksGAAS4OWLTMyH4/UQ2LlwrxtOj/FFaocbZXN3nuHZYk8tl6Oxij8ipAxBZJyRbG7gQq9LeRht2AGBMoCfOvjPeYHfSi2E9sfPUFVy6XgoA2PH8KO28/p10W12GdeuIBWN7Y7CvC2ys5Bjdyx3Xi8sx/etYPBjUBU+PCsDLmxKxJSEL4X080NvLCav2Vncxhvb20L4fnxzhj+2J2bhw7Xa37KePDMIdtwL5yPdvj6+bOrQzfjuepf3frYMCSntrdOqggMJajjv8XHFi6Visi7uED6JStfMN6KyEwtoK3kp7eCvtkR45AUD159+2+Xdh7KcHMHVoZ4zs4YY7AzrBx8Ue38dc0tneAV2UcLazwek3x+Hk5UJMW3NE73MY0d+r3rTNz4TgRmkl7nhHf9f63n+Gwt7WCoU3K6GwliM5uxBXiyrw7E/x2DD7TshkMsQvvhfbTmTh5U0n6i1fN3Tq49ZBgXH9PKERgA1x1WdbXzvzDjz1/VHt5+1X/xeEZ36Mb3A9vz03Amdzigx+zv318t3o8fpOVGkEhATUf8+aE8NLO2djJcfpN8dBLpPpNEnX9eq43lh35BLcnBQ4ebkQ794/QDv4t5vb7Q++x4KrP0z7+dT/pfHe1AGYNrwrMvNLUXizEj63WlI+euh2K8jf/7oH9608iH9FBMLORo4xvT3hbG9dL4jVNCOP6une4PadeTMC9rZWKK2owgdRqdidkouH7/CFs50NYl8PQwc9LQOLxuvvatF3Ne/wPh4Yeqt1yd1JoXOxspkj/LD28EUsGh+ovZxDaa3rgdRcldXOxgp2NlaYNdJfJ7z413pe1z0djMs3SjFxoDcKb1bWCy8fPTQIw7p1RGmFGl1c7ZF0uRDBAZ10nre7enTCnf6d8HzY7V/DSnsbvBjeUye8xL0eBg9nO531v3v/AAy79YVTuyHjnfv7w9/NESO6u0EQBJzLK8bYTw8AqA6d1nIZOjna4npJRb3nTiYD7G2tcE+t4Lru6WB8d+gi3p7SH57OCqRdLUZxuRqeznbwcFLgi0eHINBL/xmlXR1tkV9SgZE93KCwtsLf/7oH14rLEeDeAYmZBdr5RvZww8Hz17T/D+qixInLhTrrqmlR/GbGMFwvqcDUIZ0hM/AecXdS4IkQP7h3UODj6LNY+dgQuDsp8NTI261eU4d2rhcQm8qmgdaUpujkaIu7erjBzsaq3rp+nRuCzccua7s/wvt4YNeZPMy6VXtYH8/b6+mgQNRLo7X/f/zQIMy7pwe6uztCJpNpw0vtFlylgw32/DMUX+5Lw/tRKQCA+4fc7s52tLVCSYUajw7visipA7ThZf3TwbjD3xVWMhlksts/xJQONph7d3dMGugDVVkldp/J02mRA3R/tPXydMKZNyNgZyPXme7WwRbXiivQy7MD7u3rqW3RcLC1xp0BnfDl9KHYlphd74hMfayt5NoWRKB6PExfb2eczS2Cp7Md7G2rfyDVvKaCulW/j2q3NMvlMvRwr/+6TnkrAnZ1Wo/eub8//r0lSft588uzIQjq1hEymQx5RWXa8AIAa2cOx/t/pmDqkM4Y5ueKFY8OwXeH0vHOlAEor1IjMbMAl66XYu3hiwAAK5kMDw/zhdLeBv07K7H/7FUs3pqk8/h7FoRi/9k8PDTMt9HnRkwMLwQH28ZfBvPu6YHnQrtDEIDiiio4291u5gzt5Y7FE/ugr4/+JvSUtyJwo7QC3srqsOLr6gBDL3tfVwckLG28+XZoA91RQ7u64HhGAXyUtz84HGytsfy+flg2qa/2Q8yzzhd0Q2qCxJTBPtiamK2d/vWM2+Nuvvq/ILy8KRFLJlafU2fZpL6YOcIP3To5aOep3c3mrKfLbeOcO7FkaxLentJfZ/pdPdy0f4f38cDrW3SXezBId3zTiFrz13g5vJc2gNT1+J1d8dORDNw3yEcnuMS9HoZL+aU6y93h54oeHh0Q4OaI6cG3m+BlMhl6eTrh3r6eiD6di3+N6w2ZTIaYRWHQCALKKzWwtZbj+5iLsLOW6/1SvquHm8629vBw0ln/pEE+9ZapsX3+XfjjVA4evRWgfV0d4Ota/dzf29cTnRxtEdStI9Y8MQy5qjJ8/FcqngjxQ//OSvx45BIid57RhsuafVb7y7sx4wd4Y7yBbkBfVwckvTEO6VdL8OKmBKPGSgR164hh3Tqia63XkTF+mTtC+/r97JHBmPldnHbsRFA3V+2XKVDdLZJdUNakx5LLZejhUb87TN9PoJkj/HAqqwBj++q2XPzx4mjsTLqCx++sfh39Pn8k0q+X6H39atcvk2n3q74fSXXVfAbUtnXeXdh56goeHd4VTnb134fjB3jj8o2bTQovNQ68eg+yC2/izma2SPSs1SWzeGIfDPZ1qRdcAGB6cDc8PMwX1nIZisp1P4vrUjrY4N37b7fk3TfIB/fVeg8N6Vr9OWpjJUNm/k0M7KKETCbTvo4fv7MbvjqQhsz8262qXTs54P9C/Jq1jaYkEwRB6sO1TUqlUkGpVKKwsBDOzo33R7c1JRUl6BBZ/SYoXlQMR1vTDsKTWmZ+KS5dL9V2++iTpyrDd4cv4rHhXbUfco2ZtiYGRy7kA9D9RfT417E4eP4adi+4G91vnVGyUq3BvHXHMdzfFU+P0j+2pSE3Siog4HaXXXOoyiohaIAZ38VhbD9PPBdquC+8pmtgz4K7DZ4Vs7xKjaPpNzDMr6PeD8y6NBpB5xdxbZVqDS5dL0F39w5N6ro0lyq1BlZymcGa1BoBMWnX8d2hdLw1pb+2ZdDSvPF7Mv5MykF2YRkA6Lx2geoWPzH2yx3v7MLVonIkLLm3SeOFWrsbJRUI/WgfRvdyx8DOSozu5a5tbRbLot9O4WZFFT59ZHCz9tHVonJt95Wp9sOoD/Zow0tj4xJbypjvb4aXNqathxexZBXcxLJtSXhqZPXJ/GpoNAKKyqsaHZjcmv2ZnINrxeU6rSTUtgmCgNEf7kVRWRWO/ju8xV1PTVFepUZZpcai3yt1Vak19c5a3pqJEV5Gvr8Hl2+0vvDCbiMiVA8Ord0FVEMul1n8h/G4fvUHGVLbJpPJsHdBKDRCy8fMNJXC2goK68Zb7SyJJQUXoProxBpNaUFtio8fGoRH/3sE/57Yui4xw/BCRNQGWdoXL7Wco8Iaa2+d/0rfWJ/mCA7ohLNvj291ryeGFyIiojYitHfTTnhojNYWXACepI6IiIgsDMMLERERWRSGFyIiIrIoDC9ERERkURheiIiIyKIwvBAREZFFYXghIiIii8LwQkRERBaF4YWIiIgsCsMLERERWRTRwss777yDESNGwMHBAS4uLo3OX1lZiddeew0DBgyAo6MjfHx88MQTTyA7O1usEomIiMgCiRZeKioq8NBDD2Hu3LlNmr+0tBTHjx/HkiVLcPz4cfz2229ITU3FfffdJ1aJREREZIFEuzDjG2+8AQBYu3Ztk+ZXKpWIjo7WmbZy5UoMHz4cGRkZ6Nq1q6lLJCIiIgvUqq8qXVhYCJlM1mC3U3l5OcrLy3WWAQCVSiV2ea1SSUUJUFb9t0qlgtpWLW1BRERETVDzvS0IQqPzttrwUlZWhtdeew2PPvoonJ2dDc4XGRmpbeWpzdfXV8zyLILPez5Sl0BERGSUoqIiKJXKBucxKrwsXLgQ77//foPznDlzBoGBgcastp7Kyko8/PDDEAQBX375ZYPzLlq0CK+88or2f41Gg/z8fHTq1AkymaxFddSlUqng6+uLzMzMBgOVpWrr2we0/W3k9lm+tr6N3D7LJ9Y2CoKAoqIi+Pg0/sPbqPCyYMECzJw5s8F5AgICjFllPTXB5dKlS9izZ0+jT4xCoYBCodCZ1pSjm1rC2dm5zb4ogba/fUDb30Zun+Vr69vI7bN8YmxjYy0uNYwKL+7u7nB3d29WQU1RE1zOnTuHvXv3olOnTqI9FhEREVkm0Q6VzsjIQGJiIjIyMqBWq5GYmIjExEQUFxdr5wkMDMSWLVsAVAeXBx98EMeOHcO6deugVquRk5ODnJwcVFRUiFUmERERWRjRBuwuXboU33//vfb/IUOGAAD27t2L0NBQAEBqaqr26KCsrCxs374dADB48GCdddVeRkoKhQLLli2r103VVrT17QPa/jZy+yxfW99Gbp/law3bKBOackwSERERUSvBaxsRERGRRWF4ISIiIovC8EJEREQWheGFiIiILArDSxOtWrUKfn5+sLOzQ3BwMOLi4qQuyWSWL18OmUymc2vpWZKldODAAUyaNAk+Pj6QyWTYunWrzv2CIGDp0qXw9vaGvb09wsPDce7cOWmKbabGtnHmzJn19mlERIQ0xTZDZGQk7rjjDjg5OcHDwwNTpkxBamqqzjxlZWWYN28eOnXqhA4dOuCBBx5Abm6uRBUbpynbFxoaWm8fPvvssxJVbJwvv/wSAwcO1J7ELCQkBH/88Yf2fkvedzUa20ZL3n/6vPfee5DJZHjppZe006TcjwwvTbBp0ya88sorWLZsGY4fP45BgwZh3LhxyMvLk7o0k+nXrx+uXLmivR08eFDqkpqtpKQEgwYNwqpVq/Te/8EHH2DFihVYvXo1YmNj4ejoiHHjxqGsrMzMlTZfY9sIABERETr7dMOGDWassGX279+PefPm4ciRI4iOjkZlZSXGjh2LkpIS7Twvv/wyfv/9d2zevBn79+9HdnY2pk6dKmHVTdeU7QOA2bNn6+zDDz74QKKKjdOlSxe89957iI+Px7FjxzBmzBhMnjwZycnJACx739VobBsBy91/dR09ehRfffUVBg4cqDNd0v0oUKOGDx8uzJs3T/u/Wq0WfHx8hMjISAmrMp1ly5YJgwYNkroMUQAQtmzZov1fo9EIXl5ewocffqidVlBQICgUCmHDhg0SVNhydbdREARhxowZwuTJkyWpRwx5eXkCAGH//v2CIFTvMxsbG2Hz5s3aec6cOSMAEGJiYqQqs9nqbp8gCMLdd98tvPjii9IVZWIdO3YUvv766za372qr2UZBaDv7r6ioSOjZs6cQHR2ts01S70e2vDSioqIC8fHxCA8P106Ty+UIDw9HTEyMhJWZ1rlz5+Dj44OAgABMnz4dGRkZUpckivT0dOTk5OjsT6VSieDg4Da1PwFg37598PDwQO/evTF37lxcv35d6pKareZklq6urgCA+Ph4VFZW6uzHwMBAdO3a1SL3Y93tq7Fu3Tq4ubmhf//+WLRoEUpLS6Uor0XUajU2btyIkpIShISEtLl9B9TfxhptYf/NmzcPEydO1NlfgPTvQdHOsNtWXLt2DWq1Gp6enjrTPT09kZKSIlFVphUcHIy1a9eid+/euHLlCt544w2MGjUKSUlJcHJykro8k8rJyQEAvfuz5r62ICIiAlOnToW/vz/S0tLw+uuvY/z48YiJiYGVlZXU5RlFo9HgpZdewl133YX+/fsDqN6Ptra29S7Caon7Ud/2AcBjjz2Gbt26wcfHBydPnsRrr72G1NRU/PbbbxJW23SnTp1CSEgIysrK0KFDB2zZsgV9+/ZFYmJim9l3hrYRsPz9BwAbN27E8ePHcfTo0Xr3Sf0eZHghjB8/Xvv3wIEDERwcjG7duuHnn3/GrFmzJKyMmmvatGnavwcMGICBAweie/fu2LdvH8LCwiSszHjz5s1DUlKSRY/Daoih7ZszZ4727wEDBsDb2xthYWFIS0tD9+7dzV2m0Xr37o3ExEQUFhbil19+wYwZM7B//36pyzIpQ9vYt29fi99/mZmZePHFFxEdHQ07Ozupy6mH3UaNcHNzg5WVVb0R1Lm5ufDy8pKoKnG5uLigV69eOH/+vNSlmFzNPmtP+xMAAgIC4ObmZnH7dP78+dixYwf27t2LLl26aKd7eXmhoqICBQUFOvNb2n40tH36BAcHA4DF7ENbW1v06NEDQUFBiIyMxKBBg/D555+3mX0HGN5GfSxt/8XHxyMvLw9Dhw6FtbU1rK2tsX//fqxYsQLW1tbw9PSUdD8yvDTC1tYWQUFB2L17t3aaRqPB7t27dfo225Li4mKkpaXB29tb6lJMzt/fH15eXjr7U6VSITY2ts3uTwC4fPkyrl+/bjH7VBAEzJ8/H1u2bMGePXvg7++vc39QUBBsbGx09mNqaioyMjIsYj82tn36JCYmAoDF7MO6NBoNysvLLX7fNaRmG/WxtP0XFhaGU6dOITExUXsbNmwYpk+frv1b0v0o+pDgNmDjxo2CQqEQ1q5dK5w+fVqYM2eO4OLiIuTk5EhdmkksWLBA2Ldvn5Ceni4cOnRICA8PF9zc3IS8vDypS2uWoqIiISEhQUhISBAACJ988omQkJAgXLp0SRAEQXjvvfcEFxcXYdu2bcLJkyeFyZMnC/7+/sLNmzclrrzpGtrGoqIi4Z///KcQExMjpKenC7t27RKGDh0q9OzZUygrK5O69CaZO3euoFQqhX379glXrlzR3kpLS7XzPPvss0LXrl2FPXv2CMeOHRNCQkKEkJAQCatuusa27/z588Kbb74pHDt2TEhPTxe2bdsmBAQECKNHj5a48qZZuHChsH//fiE9PV04efKksHDhQkEmkwl//fWXIAiWve9qNLSNlr7/DKl7BJWU+5HhpYm++OILoWvXroKtra0wfPhw4ciRI1KXZDKPPPKI4O3tLdja2gqdO3cWHnnkEeH8+fNSl9Vse/fuFQDUu82YMUMQhOrDpZcsWSJ4enoKCoVCCAsLE1JTU6Ut2kgNbWNpaakwduxYwd3dXbCxsRG6desmzJ4926LCtr5tAyB899132nlu3rwpPPfcc0LHjh0FBwcH4f777xeuXLkiXdFGaGz7MjIyhNGjRwuurq6CQqEQevToIbz66qtCYWGhtIU30VNPPSV069ZNsLW1Fdzd3YWwsDBtcBEEy953NRraRkvff4bUDS9S7keZIAiC+O07RERERKbBMS9ERERkURheiIiIyKIwvBAREZFFYXghIiIii8LwQkRERBaF4YWIiIgsCsMLERERWRSGFyIiIrIoDC9ERERkURheiIiIyKIwvBAREZFFYXghIiIii/L/tQ5isQvHCb8AAAAASUVORK5CYII=\n"
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
many's avatar
many committed
   "source": [
    "plt.figure()\n",
    "plt.plot(timesteps, values)\n",
    "plt.axvline(stationary_timestep, color=\"green\")\n",
    "plt.show()\n"
   ],
   "metadata": {
    "collapsed": false,
    "pycharm": {
     "name": "#%%\n"
    }
   }
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 2
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython2",
   "version": "2.7.6"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 0
}