{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#AUTHOR: Lisa Blum Moyse\n",
    "#        lisa.blum-moyse@inria.fr\n",
    "#\n",
    "# REFERENCE: Blum Moyse & Berry. Modelling the modulation of cortical Up-Down state switching by astrocytes\n",
    "#\n",
    "# LICENSE: CC0 1.0 Universal"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "import matplotlib.pyplot as plt\n",
    "import numpy as np\n",
    "import math as math\n",
    "from ipykernel import kernelapp as app\n",
    "from scipy.optimize import curve_fit\n",
    "from mpl_toolkits.mplot3d import Axes3D\n",
    "#%matplotlib qt\n",
    "%matplotlib inline\n",
    "import warnings\n",
    "warnings.filterwarnings('ignore')\n",
    "from scipy import interpolate\n",
    "import pandas as pd"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "# with astrocytes Jai>0 3x\n",
    "\n",
    "taue = 20; taui = 10; tau_astro = 160\n",
    "M = 10\n",
    "Ne=400*M; Ni=100*M; Na = 200*M\n",
    "th = 20 ; tha = 13; taub = 1; tauba = 1;\n",
    "Cae = 0.5; Cai = 0.5\n",
    "Jee = 280/Ne*taue; Jei=-70/Ni*taue; Jie=500/Ne*taui; Jii=-100/Ni*taui;\n",
    "\n",
    "Jia = 0.4*11/Na*10*20*taui; Jaa = 2/Na*tau_astro; Jea = 1*11/Na*10*20*taue; \n",
    "Jai = 0.2/22*100/Ni/1.5*(1/Cai)/10*3*tau_astro; Jae = 0.2/8*400/Ne/1.5*(1/Cae)/10*tau_astro\n",
    "   \n",
    "siga = 3\n",
    "Ka = 600; \n",
    "Vli=6.5; Vle=7.6; Vla = 7; Vr=14; Vra = 9; Nmax = 2\n",
    "disp = 1\n",
    "LE = []; LI = []; LA = []; E01=[]; E02=[]; E1=[]; I01=[]; I02=[]; I1=[]\n",
    "beta = 1\n",
    "Lsig = [1.5]\n",
    "epsilon = 1e-5; epsilona = 1e-5\n",
    "delta = 8e-5\n",
    "leng = 0.007#0.01\n",
    "ini = -0.01\n",
    "du = 0.025\n",
    "\n",
    "RE = np.arange(ini,leng,delta)\n",
    "\n",
    "re = np.array([np.array([RE,]*math.ceil((leng-ini)/delta)).transpose(),]*math.ceil((leng-ini)/delta)).transpose()\n",
    "ri = np.array([np.array([RE,]*math.ceil((leng-ini)/delta)),]*math.ceil((leng-ini)/delta))\n",
    "ra = np.array([np.array([RE,]*math.ceil((leng-ini)/delta)),]*math.ceil((leng-ini)/delta)).transpose()\n",
    "\n",
    "\n",
    "for sig in Lsig:\n",
    "    sige = sig\n",
    "    sigi = sig\n",
    "    print('/Ce'+str(int(round(sig*1000)))+'.csv')\n",
    "    Se=0; Si = 0; Sa = 0\n",
    "    for u in np.arange(1e-10,Nmax,du):\n",
    "            print(u)\n",
    "            mue = Vle+Ne*Jee*re*taub+Ni*Jei*ri*taub+0.1*Na*Jea*ra*tauba - Ka*beta*re\n",
    "            Se += du*np.exp(-u**2)/u*(np.exp((th-mue)*2*u/sige)-np.exp((Vr-mue)*2*u/sige))\n",
    "            mui = Vli+Ne*Jie*re*taub+Ni*Jii*ri*taub+0.1*Na*Jia*ra*tauba \n",
    "            Si += du*np.exp(-u**2)/u*(np.exp((th-mui)*2*u/sigi)-np.exp((Vr-mui)*2*u/sigi))\n",
    "            mua = Vla+Cae*Ne*Jae*re*taub+Cai*Ni*Jai*ri*taub+Na*Jaa*ra*tauba\n",
    "            Sa += du*np.exp(-u**2)/u*(np.exp((tha-mua)*2*u/siga)-np.exp((Vra-mua)*2*u/siga))\n",
    "    Ge = 1/(taue*Se)\n",
    "    Gi = 1/(taui*Si)\n",
    "    Ga = 1/(tau_astro*Sa)\n",
    "    E = np.abs(re-Ge)\n",
    "    I = np.abs(ri-Gi)\n",
    "    A = np.abs(ra-Ga)\n",
    "    idxe = np.nonzero(E<epsilon)\n",
    "    idxi = np.nonzero(I<epsilon)\n",
    "    idxa = np.nonzero(A<epsilona)\n",
    "\n",
    "    \n",
    "    if disp == 1:\n",
    "        fig = plt.figure(figsize=(8,8))\n",
    "        ax = fig.add_subplot(111, projection='3d')\n",
    "        ax.scatter(re[idxe],ri[idxe],ra[idxe],color='r')\n",
    "        ax.scatter(re[idxi],ri[idxi],ra[idxi],color='b')\n",
    "        ax.scatter(re[idxa],ri[idxa],ra[idxa],color='g')\n",
    "        ax.set_xlabel('re')\n",
    "        ax.set_ylabel('ri')\n",
    "        ax.set_zlabel('ra')\n",
    "        plt.show()\n",
    "        \n",
    "\n",
    "    folder = './data/idx_FP_Jai_sup_0_2000A'\n",
    "    Ce = np.concatenate((re[idxe].reshape(len(re[idxe]),1),ri[idxe].reshape(len(re[idxe]),1),ra[idxe].reshape(len(re[idxe]),1)),axis=1)\n",
    "    np.savetxt(folder+'/Ce'+str(int(round(sig*1000)))+'.csv', Ce, delimiter=',')\n",
    "    Ci = np.concatenate((re[idxi].reshape(len(re[idxi]),1),ri[idxi].reshape(len(re[idxi]),1),ra[idxi].reshape(len(re[idxi]),1)),axis=1)\n",
    "    np.savetxt(folder+'/Ci'+str(int(round(sig*1000)))+'.csv', Ci, delimiter=',')\n",
    "    Ca = np.concatenate((re[idxa].reshape(len(re[idxa]),1),ri[idxa].reshape(len(re[idxa]),1),ra[idxa].reshape(len(re[idxa]),1)),axis=1)\n",
    "    np.savetxt(folder+'/Ca'+str(int(round(sig*1000)))+'.csv', Ca, delimiter=',')\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": 66,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAA64AAAFSCAYAAADyyYI5AAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/Il7ecAAAACXBIWXMAAAsTAAALEwEAmpwYAABoGUlEQVR4nO3dd5xcVf3/8ddntqR30gsJHZLQklCloxSRYqEJSA2KdBVpSsQfSBMpolIEQREQkd5Umgh+hUCQGjqEFhLSy/b5/P44dzKzu7Ozs3dmtr6fD+7jtnPvOTt3L9nPnGbujoiIiIiIiEhnlejoAoiIiIiIiIjkosBVREREREREOjUFriIiIiIiItKpKXAVERERERGRTk2Bq4iIiIiIiHRqClxFRERERESkU1PgKiIiIiIiIp2aAlcRERERERHp1MoLudjMNgC2BcYAw4HewCJgIfAG8Iy7ry60kCIiIiIiItJzmbu37QKzbYGZwJ7AiFaS1wMvArcCf3T3ZXEKKSIiIiIiIj1X3oGrmR0GnAFMBizj1EpCLetioAoYGi1rkW6K7NG524Dz3f2jYhReREREREREur9WA1cz2xm4DNiCELAuBu4C/gX8193faeG6/sB0YGtgX0KTYoBq4ErgQndfUfBPICIiIiIiIt1aPoFrMtp8FPgd8JC717U5I7OJwBHASYQa2Vnu/vO23kdERERERER6lnwC14cJQeZ/i5KhWV/gRGCVu19TjHuKiIiIiIhI99XmwZlERERERERE2pPmcRUREREREZFOTYGriIiIiIiIdGqxAlczqzSzCWY2Ksu5/mZ2mZn9z8xeNLOfm1mfwosqIiIiIiIiPVGsPq5mdgJwNXCzux/d5NxTwJdIz/XqwNPALq4OtSIiIiIiItJGcZsK7xGt/5x50Mz2BXYgBKu3AjcAddGxw2PmJSIiIiIiIj1Y3MB142j9QpPjhxKC1ovd/XB3nwmcSqh9PTRmXiIiIiIiItKDxW0qvAQod/cBTY5/DqwFrOfu70fH+gErgM/dfXThRRYREREREZGeJG6Naz8gmXnAzCYCw4GPUkErgLuvApYCQ2PmJSIiIiIiIj1Y3MB1MdDfzAZnHNs1Wj+bJX05sDJmXiIiIiIiItKDxQ1cX4zWxwCYWSLaduCJzIRmNhzoD8yPmZeIiIiIiIj0YHED1z8QBly6yMweBp4DtiXUqt7ZJO0O0fqNmHmJiIiIiIhIDxYrcHX3vwA3AWWEqXG2BKqB77r70ibJDyJLTayIiIiIiIhIPmKNKrzmYrPtge0Igy/9M3NQpuh8JfBroAI4v+l5ERERERERkdYUFLiKiIiIiIiIlFpeTYXNbBczKyt1YURERERERESayqvG1cyShClwHgLuAx5xd01vIyIiIiIiIiWXb+D6ErBptOtALfAYIYi9z9011Y2IiIiIiIiURN59XM1sArAfsD9hiptyQhDrwGzgHuBed9e0NyIiIiIiIlI0sQZnMrPBwFcJQeweQH9CAAvwLiGIvQ94xjX6k4iIiIiIiBSg4FGFoylvdgUOAPYBRkenHPgCuJ8QxP7d3asLykxERERERER6nKJPh2NmWxFqYvcDNo4OO1AF/AO4l9CkeElRMxYREREREZFuqaTzuJrZuoSa2P2AbQnT7zjwM3c/v2QZi4iIiIiISLdR0sC1UUZmawH7Rsu/3f2ydslYREREREREurR2C1xFRERERERE4kjEucjMKs1sgpmNynKuv5ldZmb/M7MXzeznZtan8KKKiIiIiIhITxR3OpwTgKuBm9396CbnngK+BFh0yIGngV00NY6IiIiIiIi0VawaV8LcrQB/zjxoZvsCOxCC1VuBG4C66NjhMfMSERERERGRHixu4Jqa5uaFJscPJQStF7v74e4+EziVUPt6aMy8ADCzvma2l5mda2Z/M7MPzcyjZVbMe/4u4x4fFFI+ERERERERKY3ymNcNB1ZnmYt1l2h9Q8axPwK/ATaLmVfKVsBDBd5jDTPbGZhZrPuJiIiIiIhIacStce0HJDMPmNlEQkD7kbu/nzru7quApcDQmHllWgI8BlwKHALMj3MTM+tLCK7rgdlFKJeIiIiIiIiUSNwa18XAcDMb7O5Lo2O7RutnW8hnZcy8Up5290bBr5ldFPNeFwDrRutxwPQCyyYiIiIiIiIlErfG9cVofQyAmSWibQeeyExoZsOB/sSsHU1x94ZCrs8ozzbAycBbwP8rxj1FRERERESkdOIGrn8gDLh0kZk9DDwHbEuoVb2zSdodovUbMfMqGjPrBdxIKPvx7l7dwUUSERERERGRVsQKXN39L8BNQBlhapwtgWrguxlNh1MOIktNbAf5KWFE5N+7+5MdXBYRERERERHJQ9w+rrj7MWZ2I7AdYfClf2YOygRgZpXAMuAWijgicBxmtgVwBvB5tBYREREREZEuIHbgCuDuzwDP5DhfSzTljJkNKiSvQphZOaGJcDlwcpZpfPK9z0yin6dfv37TNtpoo+IVUkSkE3lz9WoANuzbt4NLIiIiIp3VCy+88IW7D2+PvGIFrmb2/9z93DakHwT8gzAXa0c4E9gceCBq5hyLu18HXAcwffp0nz1bM+mISPe085w5ADy5xRYdXBIRERHprMzsw/bKK+7gTGeb2ffySWhmA4BHgWkx8yqImW0C/IQwcNQJHVEGERERERERiS9uU+EG4Coz+9zd/9ZSIjPrDzxCqGn9NGZehboGqATOA5ZEZcqU+gws41yNu9e1VwFFRERERESkZXFrXI8ljCj8JzPbIVsCM+tHGJBpW8IcrrvFzKtQk6L1L4AVWZZvR+cnZBz7fjuXUURERERERFoQdzqcm4GzgN7APWY2JfO8mfUFHgS+RBjFdzd3f7PAsoqIiIiIiEgPFLfGFXe/GLgaGAI8bGbjAcysD3A/sCOwENjd3d8oQlnjlnOiu1tLC3BzlPTDjONXdFR5RUREREREpLFCp8M5xcxGAgcCj5rZl4E/ALsAi4Avu/trBZcyYmZDCE2UU1KBd18zWyvjeLW7ryxWviIiIiIiItJxYte4ZjgceBLYCHiL0Jd1CSFofbkI9880h1CLm1rGR8d/1OT4r4ucr4iIiIiIiHSQggPXaPTd/YCXgD7AUkLQ+lKh9xYRERERERFptamwmf00z3s9B2wOPA3sY2b7NE3g7ue3qXTNr59YyPUt3PNI4Mhi31dERERERESKI58+rrMAb8M994mWbAoKXEVERERERKTnySdw/RdtC1xFREREREREiqbVwNXdd26HcoiIiIiIiIhkVYxRhUVERERERERKRoGriIiIiIiIdGr59HFtlZlNBqYDI6JDC4DZ7v5aMe4vIiIiIiIiPVdBgWs05c0vgE1aOP86cI6731dIPiIiIiIiItJzxW4qHM3vei8wGTCggVDTuiDatujc3WY2q+CSioiIiIiISI8UK3A1sz0J87saYbqcrwD93X20u48G+kfHnozS/MTM9ihGgUVERERERKRniVvjenq0vhPYxd3/6e61qZPuXuvu/wR2i9JYxjUiIiIiIiIieYsbuE4HHDjd3b2lRNG5H0S7M2LmJSIiIiIiIj1Y3MC1Eljq7p+0ltDdPwaWABUx8wLAzPqa2V5mdq6Z/c3MPjQzj5ZZrVw71sxOMLM7zewdM6uKlvfN7DYz27WQsomIiIiIiEjpxB1V+D1gQzOrzGwinI2Z9SL0eZ0bM6+UrYCH2nqRmY0HPiQ0V05ZHe1PjJaDzexGYKa7NxRYThERERERESmiuDWufybUoB6RR9rDo7R/jplXpiXAY8ClwCHA/DyuKSMEqY8B3wHGuns/QjA9mTAyMsDRhAGnREREREREpBOJW+P6S2Bv4Cozq3P3m7MlMrMjgKuAp6NrCvG0uw9tcv+L8rhuCTDN3V/MPOjuSeB1MzuAUJO7J3CqmV3g7tUFllVERERERESKJG7gehZhGpypwI1m9jPC1DefEAZtGgfsBEwAlhEC17PNrNmN3P38fDKM24TX3ZcBL+Y471Ez4T0JtbAbA3Pi5CUiIiIiIiLFFzdwnUUIUCE0w51AaBKcKRWlDiIEui3JK3Atscwa1rIOK4WIiIiIiIg0Ezdw/RfpwLU72Dla1wJvdWA5REREREREpIlYgau771zkcnQYM5sEfDfavcPdl3dkeURERERERKSxuKMKdwtm1ge4E+gLLCJ3k2bMbKaZzTaz2QsXLmyPIoqIiIiIiPR4PTZwNbNywhQ904A64FB3/yTXNe5+nbtPd/fpw4cPb49iioiIiIiI9HixAlczqzSzCWY2Ksu5/mZ2mZn9z8xeNLOfRzWbnYaZlQF/AvYH6glB6987tFAiIiIiIiKSVdzBmY4FrgZuBo5ucu5B4EukRxXeDNjBzHZx9w4f0CkjaD0IaAAOc/e/dmypREREREREOjd3+OQTeP55mD27ffOOG7juEa3/nHnQzPYFdgCS0bkq4Ijo2OHALTHzK4ooaL2VxkHrHR1ZJhERERERkc5o4cIQoD7/fDpYnT8/nCuPG0nGFDe7jaP1C02OH0qYJudidz8HwMxeBH4TneuwwLWFoPX2jiqPiIiIiIhIZ7FsGbzwQjpAff55+PDDcM4MNtoIvvIVmDEDpk+HzTaDvn3br3xxA9fhwGp3X9Lk+C7R+oaMY38kBK6bxcyrYFHQ+mfgQEKfVtW0ioiIiIhIj7R6NcyZ07g29a230ufXWQe23hpOPDEEqltsAQMHdlx5IX7g2o/QDHgNM5tICGjnufv7qePuvsrMlgJDY+aVmccQoCzjUGpwqb5mtlbG8Wp3XxldU0YInlNB66HufmehZREREREREensamvhlVcaN/d97TVoaAjnx4wJwekRR4Sa1OnTYdiwji1zNnED18XAcDMb7O5Lo2O7RutnW8hnZcy8Ms0B1s5y/EfRknIzcGS0vT1wSLTtwNVmdnWOPE5RbayIiIiIiHQ1DQ3wxhuNm/v+738heIUQkM6YAfvum27yO2ZMx5Y5X3ED1xcJAzQdA/zSzBLRtgNPZCY0s+FAf+CNAspZiMwpfyqAka2k71RT94iIiIiIiDTlDu+807i574svhmbAAAMGwLRpcMop6SB14sTQX7Urihu4/gHYE7jIzHYnNBHeElgBNG2Gu0O0LjhwdfeJMa55kvTUPCIiIiIiIl2KO3z8cePmvrNnw9Kl4Xzv3qEf6jHHhCB1xgzYYANIJHLetkuJFbi6+1/MbA/gKNJT41QD381oOpxyEFlqYkVERERERKS5BQsaN/edPRs+/zycKy+HqVPhwAPTNamTJ0NFRceWudRiz77j7seY2Y3AdsBS4J+ZgzIBmFklsIwwDc5DBZRTRERERESk21m6tPk0NPPmhXNmsPHGsOeeIUCdMSNMQ9O7d4cWuUMUNG2suz8DPJPjfC0ws5A8REREREREuoNVq5pPQ/P22+nz664L224LJ5+cnoZmwICOK29nUlDgKiIiIiIiIs3V1DSehub55+H11yGZDOfHjg3B6ZFHpqehGVrwBKLdlwJXERERERGRAtTXN5+G5uWXm09Dc8AB6X6po0d3bJm7mlYDVzP7IXCNu1cVK1MzmwGs5e4PF+ueIiIiIiIipZZMNp+GZs6c9DQ0AweGaWhOPTXdL3XttbvuNDSdRT41rpcAPzCzS4CbsowanDcz+xJwJrAX8DNAgauIiIiIiHRK7vDRR82noVm2LJzv0yf0Qz322PQ0NOuv372moeks8glcLwROAy4DLjSzB4HbgKfdfUGuC82sAtgc2Bc4FJhImFP1OeCeuIUWEREREREpts8/bz4NzYIo4ikvh003hYMPbjwNTbk6X7aLVj9mdz/XzH5LCGAPBb4OHABgZh8B/wMWAouBGmAIMBRYB9gMqIxuZcC7wE/c/fbi/hgiIiIiIiL5W7IkPQ1NKkj96KNwzgw22QT23jvd3HfTTXvmNDSdRV7fD7j7J8B3zOwswvQ2RwPjgAnR4lkuS7XirgceBK4FHnX3bGlFRERERERKYtUqePHFxrWp77yTPr/eerD99unmvltsAf37d1x5pbk2VWy7+6fALGCWmU0BdgS2BsYAw4HewCJCDezrwL+AZ9x9RRHLLCIiIiIiklVNDfzvf40HT3rjjfQ0NOPGheD06KPT09AMGdKxZZbWxW6R7e6vAq8CvylecVpmZn2BnYBpwJbRekJ0+mfuPiuPe4wEzgD2ia6tAl4DbgZ+r9pgEREREZGuo74+zI3adBqaurpwfq21QpD6jW+k+6WOGtWxZZZ4ulJX4q2Ah+JebGbTgEeBYdGhlcAA4EvR8i0z29fdawotqIiIiIiIFFcyCW+/3Xwamqpo0s6BA0Ngevrp6X6pEyZoGpruoisFrgBLgBczll8BrX5nYmaDgAcIQetc4HB3n21mlcBx0X2+Eq1PKE3RRUREREQkH+4wb17jgZNeeKHxNDRbbgkzZ6b7pa63nqah6c66UuD6tLsPzTxgZhflee0PCQFuFbC3u78P4O61wDVmNpAwavJMM7vC3d8qYrlFRERERKQFX3wBr77afEkFqRUVYUTfQw5JN/fdZBNNQ9PTdJnH7e4NBVx+RLS+PRW0NnE1cDbQH/g2cF4BeYmIiIiISBMrVoT+qJnB6SuvhLlTUwYPhqlT4dBDw3r69BC09urVYcWWTqLLBK5xmdmGpAdxejhbGndfaWZPA3sRmgy3W+D6l7/8pb2yEhHJW6rPxF/efrtDyyHSXRx44IEdXQSRdlNTA2++mQ5MU0HqBx+k0/TpA5Mnw157wZQpYZk6FUaPVp9Uya7bB67AlIztV3Oke5UQuG5S2uKIiIiIiHR9DQ3w7rvNm/i+9VY4B6E574YbwjbbwLHHpoPUSZPUH1XapicErmMytj/JkS51bqCZ9Xf3lSUsk4iIiIhIl+AOH3/cPEB9/XWorg5pzGCddUJQ+vWvpwPUDTaAysqOLb90Dz0hcB2Qsb06R7rMcwMI0+U0YmYzgZkAEyZMaHpaRERERKRLazpQUqqp7/Ll6TRjxoSg9IQT0gHqJptAv34dV27p/ooWuPaEWkp3vw64DmD69OnewcUREREREYllxQp47bXmtaiZAyUNGRKC0m9/O/Q/nTIl9EsdOrTl+4qUSjFrXJ81s6+6+0dFvGcxrMjY7gssbyFd3xauERERERHpkmpqYO7c5gFq5kBJffuGgHTvvdM1qFOmaKAk6VyKGbhOBP5rZvu7+3OtJTazfu6+qoj5t+TTjO2xtBy4jo3Wy9uz5lijDIpIZ7TznDkAPLnFFh1cEhERyUe2gZJeeQXefrvxQEkbbaSBkqRrKmbguj3wAPCEmR3p7ndmS2RmawEnE2ZbWKuI+bckcyThKcAbLaRLjT78emmLIyIiIiIST+ZASZlTzbzxRvaBkr7xjXQz3/XX10BJ0nUVLXB191fMbAZwH3Cbma3v7hemzpvZ2sAPgaMIzXKfKlberZTrTTObR5jLdU+gWUBtZv2AHaLdv7dHuUREREREclm4sHkT32wDJU2dCrvumq5B3XhjDZQk3U9RRxV29wVmthNwC/BzM9sQuBo4FTgQKAPuBS7KpzlxEd0CnAscbGY/d/cPmpz/PtAfaABubcdyiYiIiEgPl+9ASVOnwmGHpQNUDZQkPUnRp8Nx9xrgIDO7GzgsWuqAPwEXu/ubce9tZkMIwW9KqjV+36gJckp1k36qlwHHAqOAB83sCHd/wcwqgWOAn0fprnP3t+KWT0RERESkJZkDJWU28/3ww3SabAMlTZ0Ko0ZpoCTp2YoeuJrZPsCPCX1eV0Z5fAr8wt3fLvD2c4C1sxz/UbSk3Awcmdpx92VRuR4FNgFmm9kKoDdQESX7O3BageUTERERkR4uNVBSZnD66qvZB0rabjuYOTMdpE6cqIGSRLIp5jyuhwNnEALDRcB5wK+BjYG7CdPlfN3dny5Wnm0R1bBOJgTV+wDjgVWEwZtuBm5092RHlE1EREREuh53+Oij5k18X3891K5CqCVdd90QlH7zm+kAVQMlibRNMWtcbwY+BE4Bfu/uVdHxZ81sW+BB4B9mNtPdb4mTgbtPLKSA7v45cHq0iIiIiIjkpelASa+8EvqlZg6UNHZsCEozB0raZJPQ/FdEClPMwPUI4DZ3b2h6wt3fM7NtgLuAm8xsA3c/t4h5i4iIiIgUbPnyUGPatJnvggXpNJkDJaWmmpk8ORwXkdIo5nQ4f2rl/DIz2xP4HXC2ma3n7gcXK38RERERkXxVV6cHSspcmg6UNGUK7LNPugZ1yhQNlCTSEWIFrmZ2JmGgpH9lNAlulbvXA8ea2TvA/4uTt4iIiIhIvurrw0BJTZv5vv02JKPRTSoqNFCSSGcXt8b1QsKcpyOAvAPXFHe/yMwKHWFYRERERARoPlBSqqnvG29kHyjpW99KN/Ndf/0QvIpI5xU3cF0M4O5L4mbs7nfFvVZEREREeq4FC5o38X31VVixIp0mNVDS7runa1A33lgDJYl0VXED1/8BO5hZn7Y0FRYRERERydfy5WHk3qYBauZASUOHhprTI45IB6gaKEmk+4kbuP4e2AU4BLixeMURERERkZ6m6UBJqWa+8+al0/TrFwLSzIGSpk6FkSM1UJJITxArcHX3P5vZd4BLzey/7v5akcslIiIiIt1MaqCkplPNZBso6UtfajyS79pra6AkkZ4s7qjC9wGvA9sCz5nZj4HfRaMGi4iIiEgP5h5qS5s28W06UNJ664Wg9MAD0wGqBkoSkWziNhXeB/Bo24ArgfPN7EHgCeBF4FUFsiIiIiLdWz4DJY0bp4GSRKQwcQPXXwGbAZsDQwnB62Dg0GgBqDWz1whB7IvAi+7+XCGFFREREZGOkRooqWkz34UL02kyB0pKTTUzeTIMHtxhxRaRbiJuH9cfpLbNbAKwRZNlHNAL2DLaP4ZQQxs3UC4aM/sycBywNTCSUK7PgP8A17n7Ux1YPBEREZEOVV0dmvQ2rUFtOlDSlCmw776N+6FqoCQRKZWCA0l3nwfMA+5NHTOzoTQOZLcE1i80r0KYmQG/BY7POFxNCFwnRcuhZvYrdz+9A4ooIiIi0m7q6+Gdd5oHqE0HStp4Yw2UJCIdryQ1oO6+GHgsWgAwsz6lyKsNjiQdtP4VONvd3wYwsw2Bi4H9gNPM7Gl3v7tDSikiIiJSRE0HSko19X3jDaitDWmyDZQ0dWo4poGSRKQzaLemu+5e1V55teCIaP0OcEjmwFHu/qaZfQuYC6wDHAgocBUREZEuZcGC5n1QX3ut+UBJU6fCV76SrkHdaCMNlCQinVuH9zltR6Oj9f+yjXbs7nVm9hIhcO3fngUTERERaYtly0JA2rSZb+ZAScOGhQD1O99JB6gaKElEuqqeFLi+B2wIbGZm5U2DVzOrIIySDDC7ncsmIiIi0kxVFcyd2zg4feUV+OijdJqmAyWlRvMdMUIDJYlI99GTAtffAnsB6wG3mdlZ7v4OrOnjehGhtvVdwnQ/IiIiIu0iNVBS02a+77yTHiipsjI06d1hh3RwOmUKTJiggZJEpPvrMYGru99vZqcRBmH6JvBNM0v1u+0DLCUEt+e6+/KOKaWIiIh0Z+7w4YfNm/hmDpSUSKQHSjr44HSAqoGSRKQn6zGBK4C7X2FmbwM3AiMIAWtKL2AAMAhYnO16M5sJzASYMGFCaQsrIiIiXZZ7GCipaYD66quwcmU63fjxISjNHChp442hT0fPxSAi0sn0mMDVzPoCNxFGDJ4NHAa8CBhhrtkLo2N7mtlu7v5y03u4+3XAdQDTp0/3diq6iIiIdCK1tfDZZ/DJJ/Dppy2vMwPU1EBJRx6ZbuY7eTIMGtRhP4aISJfSYwJX4FJC0PoWsGOT6Xn+YWb/Bl4CNgCuAXZo9xKKiIhIh0km4YsvWg9IM0fuTamshDFjwrLZZrDXXjBxogZKEhEplh4RuJrZAKImvsCvs80p6+5VZvZr4CrgS2Y2wt0XtGc5RUREpDRWrGg9IP3sM6ira3ydWQg6x4wJzXq33hrGjg37methwxSYioiUUo8IXAm1qKmf9d0c6d7O2J4EKHAVERHpxFLNdnMFpJ980rjZbsrAgengc+edGwejqe1RozQgkohIZ9BTAtdkxvbaOdKNzNheUaKyiIiISCtSzXZbC0hba7a76aaw557Na0nHjIH+/dv/5xIRkXh6SuA6F6gijCJ8rJld7+71mQnMrIx0c+IlwJvtW0QREZGeIVez3dR2a812x42DrbZSs10RkZ6iRwSuUf/VG4CTgC2B+83sDOC1KMkUwuBN20X7V7h7Q/uXVEREpOsqtNluKvDcaafsAama7YqI9Fw9InCN/BhYH9gzY6mJzvXKSHcbcEH7Fk1ERKTzKqTZbkVFOvCcOjU0220akKrZroiItKbHBK5RrevewDcI87VOA0YADnwEPAfc5O4PdlwpRURE2teKFa0HpNma7UJotjt2bFi22ip7QDpsGCQS7f9ziYhI99JjAlcAd3fgr9EiIiLSbdXWwvz5rU8BsyLLUIQDBqQDz512yh6QjhoVBkESERFpD0UJXM3MgGFAX3efV4x7ioiISHPJJCxa1HpAuiDLhG6pZrtjxoRmu3vskX203QED2v/nEhERyaWgwNXMtgTOBXYH+hGa3ZZnnB8CXBQdP9ndawvJT0REpDtbubL1gPTTT1tutpsKQGfMaHm0XTXbFRGRrih24GpmhwM3AC2O7+fuS8xsErAbcD+g/qMiItLj1NXlN9pua812d9ih5dF21WxXRES6s1iBq5ltDFxPCFqvAm4BHiE0F27qFkKN7H4ocBURkW7EPf/Rdt0bX5vZbHfKFPjKV9RsV0REpCVxa1xPByqBa9z9VAAza2ne08ej9bYx8xIREWl3K1fmN9pubZZOMJnNdqdPV7NdERGRQsUNXHcl9Fu9uLWE7v6pma0GJsTMS0REpGgym+3mCkyXL29+7YAB6cBTzXZFRETaT9zAdQywyt0/zjN9FTAoZl4iIiKtKrTZ7ujRIfCcPDk02802BYya7YqIiHSMuIFrDdDbzCyaG7VFZtYHGAwsjZmXiIj0cIU02x0+PB14Tp+ePSBday012xUREenM4gauHwBTgfWBt1pJuzdQBrweMy8REemm6upg/vzWp4DJ1my3f//Go+1mC0hHj1azXRERke4gbuD6CLApcArw/ZYSmdkw4BJCf1iNKCwi0kO4w6JFrQekCxY0b7ZbXp4eUXfyZPjyl7P3JVWzXRERkZ4jbuD6K+AE4LtmtgC4PPNk1Dz4AOACYG3gC+B3BZRTREQ6iVWrWg9IP/205Wa7qcBz2rTsAama7YqIiEhTsQJXd//czA4F/gr8FDgLSACY2WvAOoTpcozQH/YQd8/S0KtjmNlA4HuEuWXXBwYCC4G3gaeAK9x9aYcVUESknbiHQHTZstAcN7VeCNTWwVl/aR6YZmu2269fCDrHjoXtt295tN1evdr9RxQREZFuIG6NK+7+gJntCPwamJ5xauOM7TnAd939+bj5FJuZ7QLcBoyMDtUDK4Gx0bIzcA/wUvuXTkQkfzU1zQPOtq6XL4dkMsvNfxVWl12Wbra7ySaw++7Ng9IxY2DgwHb90UVERKSHiR24Arj7c8BWZrYp8CXCNDllwHzgGXefXXgRi8fMtif0te0D/BOYBfzH3ZNR8+ZNCE2cl3VYIUWk22toSAeNhQSe2ZriNtW7dwgqBw1Kr9ddt/F+tvWpHgY1erpGzXZFRESk4xUUuKa4+8vAy8W4V6mYWV/gFkLQehdwoLuvqWdw9yrghWgREWkm1ay20IBz1arW8yorax5Qjh4NG22UO+Bsuo47om6/OWGtoFVEREQ6g1iBq5lNABrc/ZM8048Byt19Xpz8iuRwQt/bKkLz5WyN40Skm6qpKTzgbLFZbRP9+zcOIAcNggkT2hZw9u0LZqX/XERERES6gkLmcf2M0Cc0H88A4wvIrxiOiNb3uvsXHVgOEWmDhgZYsaLwgLOmpvW8evVqHkCuu27bAs4BA0JtqYiIiIgUTyGBZFvrAjqs7sDMepEeQOopM1sHOAfYAxgOLAH+C/zO3R/umFKKdC/usHp14QHnypWt55VINA8gR4+GDTfMP+AcOFAj3oqIiIh0Vu1VA9qbMHpvR5lImJ4HYByhP24/oBZYTRhheF9gXzP7nbt/ryMKKdJZ1NYWZ7TahobW8+rfv3kgOX58/gHnoEFqVisiIiLS3ZU8cI36tw4HFpQ6rxyGZGyfBSwHDgHucvc6MxsPXAIcDHzXzOa6+5VNb2JmM4GZABMmTCh9qUXaKNWsttC+nPk0q62sTPffTAWQkya1LeBUs1oRERERyUdegWs0X+vOTQ73N7Of5roMGAzsHW3/N0b5iiXRZPu77n5H6oC7f2Rm3wY2BLYAzjWza9y9US2xu18HXAcwffp0L32xpadINastNODMt1ltqmlsKogcORLWX79tQaea1YqIiIhIe8m3xnUX4DwgM1jrFx1rjQHVwC/aVrSiWpGx/VFm0JoSzeX6S+BPwFrANDo22JYuora2OKPV5tOstl+/5gHk2LFtCzj79VOzWhERERHpWvINXD8AnsrY3wmoA/6T45okoUnuq8DN7v5OnAIWSea0PXNzpHsjY3ttFLh2a8lkcUarra5uPa+KiubNaidObHuz2vKOHJdbRERERKSD5PVnsLvfDNyc2jezJLDY3XcpVcGKyd0Xm9knhOl7cjXxzayHytkUeN48+F6TIZya1mIVut9V7tkZy+0ems1mCzZT2ytW0Cqz5gHkiBGhWW1bpkjp3bv1vEREREREJLu49TdHAVXFLEg7+Duh3Bubmbl7tsB044zt93PdbMkS+Nvf0vtN71bofle5Z0eVOx99+zYPIMeMaVvA2b+/mtWKiIiIiHS0WIFrVAPb1dxECFzHAwcBt2eeNLMEcHq0+wnwYq6bbbYZzJ5dglJKLNmC30Qie1oREREREelaesyf9u7+NPDXaPe3ZnaQmVUARNPh3EoYURjgHHdPdkAxJSazxouCVhERERGR7qOgoV7MbDPg+8CXgHGEkYZb4u7e0UPLHAmMAHYk1LjWmNlqGs/zen4XrVEWERERERHplmIHkmZ2InA5UEbjQY06LXdfZWa7AEcDhwNTgAGEpsFPA1e7+7MdWEQRERERERFpIlbgamZbA1dGu78BHgQeAhYDBwKjgN2BQwlT4pwMfFZoYYshagJ8Q7SIiIiIiIhIJxe3xvVkQi3rFe5+OoCFoVdr3f3xKM2fzewq4FHg58CWBZZVREREREREeqC4Q9hsT5jn9Momxxs1GXb3l4CTgHWBH8XMS0RERERERHqwuIHrSKDG3T/MOJYEemdJezdQB3w9Zl4iIiIiIiLSg8UNXFcTgtFMK4CBZtYr86C710Xp146Zl4iIiIiIiPRgcQPXT4D+ZjYw49i70XpGZkIzGwMMoouMPCwiIiIiIiKdS9zA9eVovWHGsScJwelPzaw3gJlVAldF51+JmZeIiIiIiIj0YHED1wcIQepBGceuAWqA3YCPzewZQs3sAYSBnH5dQDlFRERERESkh4obuD4E/Ax4O3XA3d8nzNu6AhgKbAsMIwStl7j7rYUVVURERERERHqiWPO4uvtyQuDa9PjdZvYUsDcwHlgG/N3d3ymolCIiIiIiItJjxQpczWzHaPNld1+aec7dFwN/KrBcIiIiIiIiIkDMwJUwEFMSGF68ooiIiIiIiIg0F7eP6zJgqbsvKWZh2puZnWlmnlo6ujwiIiIiIiLSXNzA9R1ggJn1KmZh2pOZbQic19HlEBERERERkdziBq63AxXAgUUsS7sxswTwe6A38J8OLo6IiIiIiIjkEDdwvZIQ8P3azPYuYnnay0nA9sCtwN87uCwiIiIiIiKSQ9zBmc4G/gVMBe43s9eAZ4AFQENLF7n7+THzKxozmwRcACwCTgO+37ElEhERERERkVziBq6zAAcs2p8CTM7jug4PXIHrgX7ACe6+0MxaSy8iIiIiIiIdKG7g+i9C4NqlmNlxwG7AP939lo4uj4iIiIiIiLQuVuDq7jsXuRwlZ2ZjgUuBKuD4Di6OiIiIiIiI5Cnu4Exd0bXAIGCWu78X5wZmNtPMZpvZ7IULFxa3dCIiIiIiIpJVjwhczeww4KvAS8Dlce/j7te5+3R3nz58+PBiFU9ERERERERy6PaBq5mNAK4gjHZ8nLvXd2yJREREREREpC3iDs7UlVwMDAN+C8w1s/5NzlemNjLO1bp7bTuVT0RERERERHLo9jWuwKRo/T1gRZblrIy0qWOXtGcBRUREREREpGU9IXAVERERERGRLqzbNxVubeoeM5sFnBeltXYokoiIiIiIiLSBalxFRERERESkU+v2Na4iIiIiIiLSivp6qK2FmpqwpLZzHWtHClxFRERERETaQzLZejDYUceSyY7+dHIqSuBqZkaYcqavu88rxj3bi7vPAma1+cJ334VDD4XKyvTSq1fb9ttyTVkZmLrgioiIiIjk5N55gsGm2w0Nxf1ZKypCvJCKGTLXmdt9++aXLtuxXOc33bS4P08OBQWuZrYlcC6wO9AP8Mx7mtkQ4KLo+Mndam7U6mqYPTv8Imb+gqaWYjMrPPgtRgCd61h5uYJrERERka7KPdS61dWFZqN1dY23cx0rVfp8gsGm5+vqivu5lJfnF8QNHFicYDDfY5WVPepv79iBq5kdDtwAVLSUxt2XmNkkYDfgfuDBuPl1OpMnh8A1G/fwwmQGsk0D26b7+aTJ9x7Ll+d3TbGlguv2Dphz1VJnLolE82OZSw968UVERCSmhoaODeRKnb69lZWFwLCiIiyp7dS6acDWty8MGVK6YLDpsYqKUEbpcLECVzPbGLieELReBdwCPEJoLtzULYQa2f3oToFrLpkBXGflHv4HVcoAOtf+ihX53aM9meUObFsLfPNZCr1HVy2DvhQQEek5mtbYdbbArND07u37eSYS2QO6XMFeRUX4O7Rfv9zpWrtHPscKSZ/6u0IkD3FrXE8HKoFr3P1UADNrqcH249F625h5SSmYpf+n0Vmlgus4AXSqD0Hmkkw2P9bWJc49MstSaBlS17f3P5qFavqlQE8N4EtRBv2DL9L1pFpmddbArND07T3AS+pvmjhBVp8+pQ/UCr2H/j8vAsQPXHcl9Fu9uLWE7v6pma0GJsTMS3qqzOC6X7+OLk3n4t6+wXdn+AIg2/V1daG/eTHz72pfCkBpAuejjgr3/uEPG+eV7fPJ91ip0nb169szr46+vj3z6ujrIfw/JVtw1xEjd8YNmnr37hzBW65zZWXt/3mKSLuLG7iOAVa5+8d5pq8CBsXMS0SaMgv/WJdrRquiy/dLgY7+AqDUXyCk/hDP1h8+W9PvthzLN22heXWl69szr46+vj3z6ujrW2ri2Z7BXkVFKIe6bIhIFxf3r94aoLeZmXvu6gkz6wMMBpbGzEtEpP3oS4Fgzpywfvrpji2HiIiICBC30fwHhIGZ1s8j7d5AGfB6zLxERERERESkB4sbuD4CGHBKrkRmNgy4hNAftmeMKCwiIiIiIiJFFTdw/RWwEviumZ1nZgMyT5pZHzM7FJgNTAIWAb8rqKQiIiIiIiLSI8UKXN39c+BQoA74KbCQaA5XM3sNWAz8EVib0B/2EHdfXowCi4iIiIiISM8Se2Iod38A2BF4gTCnazmh+fDGQK9oew6wo7s/VnhRC2Nmw8zsKDP7k5m9bmarzKzGzD42s3vM7ICOLqOIiIiIiIg0V9Cwme7+HLCVmW0KfIkwTU4ZMB94xt1nF17EoplP45+3mlBjPDZa9jOzh4FvuvvqDiifiIiIiIiIZFGU+R7c/WXg5WLcq4TKgeeAPwCPuvt7AGY2ETgXOAbYC7gWOLxjiigiIiIiIiJN9aSJCnd19yeaHnT3D4BjzaweOB44zMzOdveP2ruAIiIiIiIi0lzsPq5dTbagtYnfZ2xPL2VZREREREREJH+xalzN7MY2XlINLAVeAx5z9/lx8i2x6oztsg4rhYiIiIiIiDQSt6nwkYC3Ib1lpK83s5uB0919Zcz8S2HnjO1XOqoQIiIiIiIi0ljcwPUWQiC6LzAEWE2YFucTQpA6BpgG9CXM6foAMBjYEhhHGAhpfTPb3d0bCih/UZjZYOCsaPdpd3+zA4sjIiIiIiIiGWL1cXX3Iwlztw4GZgGj3H0ndz/U3Q9x952AkcB5URrcfX93n0Cora0jzAF7WGHFL5yZJYA/AqOBGuCkHGlnmtlsM5u9cOHC9iqiiIiIiIhIjxYrcDWz44CDgbPd/fxsTX7dfZW7/xw4BzjczI6Mjt9CCHYNOCRmuYvpSmCfaPsEd/9fSwnd/Tp3n+7u04cPH94+pRMREREREenh4jYVPgZIAtfkkfYa4P8BxxHmUAW4EbgQ2Dxm/kVhZpcBJ0a7p7l7/oNOVVfD/fdDr17ppbKy8X6vXrDWWpDoMYM3i4iIiIiIFF3cwHUjYFk+gyu5+0ozWw5Mzji2wMyWETUj7ghmdgnwg2j3R+5+RZtusGgR7Ltv6+kWLIDhw+Gcc+CKK2DVqnD8Bz+Au+9uHuhmLv36wS23hPS33Qbvvw9nnx32//pX+OST3NcPHAibbx7Sf/45mMGIEWG/vh7KysIxERERERGRTixu4JoABpvZEHdfkiuhmQ0BBhEGcMpUAXTIqMJmdinww2j3DHe/rM03GTECHngAampyLwMHhvQ77hgCxZSNNoLtt2+eftUqWLw4bGfW1D7+OPznP+nA9Te/gSdamZp2nXXg3XfD9re/HWqJ//3vsL/ppjB3bvZa4tQyYwZce21If+KJMGlSCLgBfvhDaGjIHTivu274GQH+9S8YOzYcSybhnXeyX1OmmYhERERERKSxuIHrq8DWhJF4z2gl7ZmEQPe11AEzG0YYcXhezPxji5oHp2paz3D3S2PdqKICpk3LP/0ee4Ql5bjjwpKv669vvP/AA1BVlTtorqhIp//BD0KgmXLiiTB/fn5BN8BHH0Hfvun9O++EJUtCutra7GU+5JB04Lr33nD88fDLX4bgfMMNs19TVpYOYk86CX72s5B+m23gjDPg8MPh44/h2GNzB829esGee4brli4NNdS77hqC+UWL4PnnW7++f38oj/uKiIiIiIhIscT9q/x6YBvgB2Y2CLjQ3T/MTGBmE4CzCX1bHbgu4/TO0XpOzPxjaRK0/tDdf9me+RdV376NA8nW7LVX4/0TTmhbfvfe23j/w4zH7R6C16aBb58+6TQPPwyjRoXtXr3g1ltbr61ONXN2h/XXh0GDwn5dXbpWuqWlvh6GDAmB6yefhC8J7rgjBK5z5jT/PLK5/XY46KBQs/21r8Hf/w7bbQf33ANnntl64Hv22bDBBiG/u+6C00+HoUPhhRfyC5w33jisV64MX1KstZaadouIiIhIjxQrcHX3m8xsD+BA4FjgWDObB3xKCFLHAGtHyQ24091vyrjFN4BlwMNxC95WZnYx6aD1dHf/VXvl3e2ZpYOtluywQ3q7shIOPTT/+/fvD3/7W3p/0iR47rnc1ySTIeCFEDzOmxeCRghNoJ99tvXAeYstQvoxY0Jt8ZgxYX/IENhss8Zpq6pCzW7msZOimZVefhl+8YsQPA8dGoL4n/yk9Z/7/fdh4kS46qrQR7q6OnzGJ58cmnC3Fvj+85+h1v3GG+G//003+77pJnjjjdzXDhgA+0SDbc+dG74smDo17H/2WfhsM9OXlyuoFhEREZGSMU/9cd/WC8P8p2cQmgIPbCHZcuBi4BJ3b2ghTclFtb+pKsIk0NokrJe11u91+vTpPnv27GIUT3oK9xDcrV4Ny5a1HjjvuWeotX7hhdC/+fvfD9ffd19+gffjj4d+0rNmwaOPhnsAHHxwqEGvqUkH902NHBmakgPst18I/OdEDSS22irUGGfK/PKiV68Q9D/ySDj3ne/AsGFw+eVh/5hjYMWK3IHzhhvC178e0t91F4wfH/KFUAOeq2+2+ksXxc7R834y9QWOiIiISBNm9oK7T2+PvGJ34HP3JHCRmV0J7AFsCQwn1LAuAF4E/u7uTQdl6giJJtsjW0nfv4RlkZ4qVSPZ1mbe06Y17k+97775jWidMmtWWFJuvz2s3UOT6mxBbzLZ+PrVGa/xOee03j86VTsNoYn3gAHp/Q8+gE8/zX5dXV1I87WvpQPXE06AAw4IgWsyGfoqt2bmzFDDnEyGUb3PPBN+9KMwyvfuu7ce9O6/f+iXvXx5GA38a18LwfjCheGLg9auHzcu/Nz19aGPdv/+CqRFREREClDwyDPuXgXcEy2dkrt/QAioRSTFLDQlrqgIgVVLmta47bdf2/K56qrG+4891nLaZDL0l84MnP/zn8aB/hNP5A6aa2vDqNmp+x16KEyOZuMyCyNbZ6ZfvrzxtTU1oX/x3nuHAcjOOy/U+G6xRRgN+9hjW/+Zf/97OPpomD0btt0WHnoo9Ku+774wwne2YDezFvkXvwhN2mfPhmuugfPPD2X4z3/gwQdbD5x33DEMrrZgQfiSYMqU0Jx7xYrw86XSVVSoibeIiIh0CRoyVUQ6j0QCevdufGyddRqf33nn/O9XXg5XX53eHz48zJ+cr7XXDrWmKdOmhYHJWmumPWNGSD9hQhhJe+ONw/7EiaGvc65rV69ON+H+/PPQ5Puss8L+7Nlw4YUtN/FOefXVEKzfdhucemoYSXvoULjgArj44sZpswW+qW4Q8+eH2u6//S0EuPfeC2+/DYMHh2XQoObrXH3dRURERGKK3ce1p1MfVxHpELmaeKdqjKdMCf2j3303DA62zz6hdvXZZ0Of6dYC7xtuYOe5c+GTT3jykkvCPMwQRtn+y19yl2/YMPjii7A9a1YI9G+Kxua77rpwLhXoZga9qe3+/VULLCIi0kV0iT6uAGa2GfB94EvAOKBfjuTu7qrhFREpRL5NvCE0i1533fT+dtuFJV9jx6aDVgj9o2+4IYygvXRpGGSs6XbmfNHJZONm33/4Q3qQsJZsuin8739h+6ijwjRQl0bTbZ9/frh/rhrfQYMaz2EtIiIi3ULsQNLMTgQuB8pQ/1ERke7PLAy0NWBA6HPbmvPPb7z/7LNhWqdly9KBbtPgd2DGIPV9+zaeD/qmm8LgXrnsu2963umddgp9lX/84xBAf//72Wt5M7eHDm3eXF1EREQ6XKzA1cy2Bq6Mdn8DPAg8BCwmzO06CtgdOJQwJc7JwGeFFlZERLq43r3DMrK1wd0JA1Nlev/9UOO6YkXLNb4TJqTTjxuXnr959Wr4619Dmsx+y02dckoYSbqqKtT+nnMOHHlk6G989tnZg97M9ciRbRs1XERERPISt8b1ZEIt6xXufjqAhT5Jte7+eJTmz2Z2FfAo8HPCdDkiIiLxlZWlA8fW3Hprert//zCdkXsISluq8U2NQF1XFwbZGjEi7C9ZEuZDXrYMVq5sOc8rr4STT4a5c+ErXwn9evfcMzR/vvTS1mt8x49X4CsiIpJF3MB1e8BJ17qmNGoy7O4vmdlJwO3Aj4CfxsxPRESkcGbpuZRHj2453cCB8Oc/p/c32gg+/jhs19c3bu6cGQBvvXVI06sX7LZbumb5iy9CU+lUusy+v5nuvz8MpvXww6GP7z/+AVOnwiOPwI03tl7ju8EGauosIiLdUtzAdSRQ4+4fZhxLAtn+tbwbqAO+jgJXERHp6srLw+jJw4a1nGbSpPRoyhCC2PfeC9vusGpV86B32TLYMmqcNHJkmDN5rbXC/qJF8Mor6XRVVdnzfe012GQT+PWv4Sc/CXkOGRIG1br//tZrfDffPPx87hrdWUREOpW4gevqLNeuAAaaWS93r0kddPc6M1sNrB0zLxERke7DLDRd7t8/9MPNZsst4dpr0/vf/nZYUmprc/fxnTIFDj88DKQFoV/wBx+k0y1fnn0+4NWrQ+B62mlwxx3wWTQ8xU9/Cs8/33JN7+DBIZBP1TjX1obRnRX8iogUnbvT4A00JBuarZOebPFca+s417anuIHrJ8BkMxvo7sujY+8CWwAzgH+nEprZGGAQIdgVERGRQlVWwvDhYclm553DknLaaWFJSSabD3K1bFm6mfGuuzauUa6vD7W+776bTl9b2zjPUaPSge63vhWaVr/wQtg/4oiwn6vGd8wY2GqrkH7FitCcu6ysrZ+MiLQTd8fxEOxEQU+xlwZvWHP/OMFYya4tQQDYlnI7Wb547AHiBq4vA5OBDYHno2NPEgZg+qmZ7evu1WZWCVwVnX+lkIKKiHQH7mFg3Lq6xkttbXp7zJjQxXLZstA6dOrU8Pf9Rx+FOCAzbbbrDz88jPEze3bopnnOOSEGefDBMBVr02ubXp+aTWb+5zBtGjz9dIghrrsutDatrAxdOFtan3VW2H7mGXj77TAoL4QKuwULsl+X2u7VK90tNJkMFXaqtCuBRCI97202++4blpQLLwxLint6aqNU0JsZyB58cONBrPr2DefffjudfsWKxnluu23oBwxhvuH114e//S3s77BDyDNXje/666ebWi9YEI736tXWT0baKDN4yRmAlCiw6axLKvjozkt3C54MoyxRRpmVUZYoI2GJNdtx1glLUJYooyJRQe/y3vldW4S82/va8bPymB6vSOIGrg8Qpro5iHTgeg3wfWA34GMzexPYABhKGMjp14UVtTjMbADwA+AbwCSgAXiLMIDU1e5em+NyEelg1dXZg63MZdAgWGedkP6xx2Ds2DC2Tm0t3Hln9msy77f99rDXXqEb4WmnwTe+AV/+MnzySZgtJVfedXXwox+FYO3tt8O9rr0WDjgAHn88dHVsze23w0EHhSB1t93gqadgxx3D+vDDW79+++1D4PrOO6Fr4wknhMD1449DMFlR0XyprAyVbZmtO8vKQhBdURH2V6wIn0FtLdTUNF6ntuvr4cwz0z/HrbemA9eLL4a77spd9sGDwwC+ED6DN96AV18N+1/7GsyZkztoXndduPrqkP6ii8I0sKecEvZ/9avQtTSVvqLCqezlVFaGdUWlM2qUM2NG+GP85VecgQOd8RMcd+eTT52KCqc8uq6iwjELf7g5vuYP+NbWQN5p415T6jI1un/C8cHRubn3hvXmvXF64a/fFe5/8u6479b4Pg31+OrV+OpV+KpVeMLw/90S7j9zK7xfP3zOTSHtJjX4yhX46o/wxVX4x6vxqtXhXtHvq8+Yjtd/Jxz74Q/wL30J329fvKYa/+Uv8T698T598N69oE+fsN87Yxk9Ch85Ek824MuX4X374uXlBX+m3T246W7BSxwJS+RcUkFMKZbyRHnJ7h13KebPa2ZFCR7zSWv6lrTTs9T/jNt0kdlA4DRgvrtfm3H8AOAmIGMGeZLApe5+VoFlLZiZrU2oGZ4YHVoNlAGpr2TnALu5+5LW7rXO5HX8gtsvSN03nUfGwMqp49mOZR7v6td31rI2/WMsczvbHyFtOd9Z7pVv2oZk+g/s6mqnusapb3AaGqC+Pr3d0ODU1ztJd9ZdN9zrw3nOqlXORhuFe732urNkqZOMrm1IhuuTDdCQdBoanN59YJddQhmefMpJOuy8UyjXPfc6S5eGtMmGcH0yGV0bXT9mNBzw9XD9jTfB8BHOV/cO119xBdTVOZhD6g+m1Hb0M26wgbPHnmH7d79zNpnsbL+9U1cH19/QOG16O32vKVOdGTOc+nq4805n8y2cDTeEFSudxx5zEmVOwsASTiLhYV1G2DZn4kRnzFioqgqf18SJzpAhsHKV88GH6WvNomub7A8ZAr37OLW1ztJlzqBBUF7u1NY5VVWEdBbKnUiEMqf2Qw1l4b97b477LgDrf/SbNv2x3uj3rsFJJsPn5Th1ddEx9+gP+tQ1UTmi/d59wnZ1Tbi+V++wX1Ud9t1T77c3W2MhoEzlhzllZVF/oIboWa959iJtZx7+rTGitYU9SyTCgoU0iQTAmlqMjgpsOuPSnX5eBTvS05nZC+4+vV3yihO45ryh2VBgb2A8sAz4u7u/U9RMYjCzMkJgOhX4DDjC3f9pZgngW8D1wADgYXffu9X7jTHn+FKWWKQH8cwvHsIfgolE2E8mw9GysnC8oT6dNvzBEH1dkfrjkXBteXnYbmiARMIoS0TXN0AildYsCvSiPK3xvTO/MGm6ba2kbe18IfcqJN987/V/ww4GYPvFdza7Pp91Xte0NX2OnzHf8njSaGgIv1fJhrA0ZKx7VRrDh4f077xj9O9njB0Tfs+efz78/jU0NFnXG/UN0FBvTJpobLNNSH/t74xp04xttjZWrjSuvsqorzfq64z6eqhPXRsda6g39t/f+PoBxtKlxsknGUcfbey+m/H+B3DO2Ra9Ky2vv3u8sddexnvvGqedDj8519h6K+Oll4zLLzcqyqOlwigvNyobreHww42NNzLee8+49x7jsMOM0aOMN980/u8/4brKitSasK4Mxyorjc02Mwb2N5YsNRYuMNZbz6gsN1avNmpqoLIyfMaVFRbe6UJ+h5KO1dVhdfVY1ATCamux2rpwvKZ2zbmwXddoTW0NVlWNLV+OLV+BAVZRia1cic19Exs8BOvbF1u0CHv9jXTzgro6ogfYeFqjPn1Cc4Xq6vTxqVPh5Zez/m9PRKQ7aM/ANVZTYTPbMdp82d2XZp5z98XAnwosVykcSQhaAb7h7v8BcPckcEcUwP4Z2MvMdnP3x3LdbINhG3Du/ucCoYmIYYT/oj+qHDDYaK2NqEhUsHDVQhauXshGa22EmTF/5XyWVS9LX5Plj1qA8QPHY2YsqV5CTX0NI/uHzl9LqpZQ01CTvja6JHM/YQkGVIYRJavrq0l6kj4Vfdbsr/nSokkZwh914VhlZUiyYgW4OQP6h/1PP3Ma6kNfvVBLll4nk6HmbsBAmLh2qOl44UUYMsSZNClc86+nQw1dMkmorUs2uT4Z0m6xeajVu/0O2HRTZ8qUUJY77khdQ1SrA8lk42M7fMnYdltjyRK4/HLjG183ttwS3n/fuPyX6T821/zxWGdRDWTY/sEPjK98GV59NfwBedllxowZ8OQTxnnnpdqmZf7RCJl/QF55JWy7jfH448aZPzZuvx3WX9+4/Tbj0ksbX19eDhXl4Q/IsvLQhPH31ydYb90E9z/g3HBjPX+6uZLBgxLcdkc9Dz1SR3l5GPwzNF0M9ygvg/IKp6wczjtjMAP6lfH3J1fxwv+qOOfUEZgZDz2+lHc+qFpzfXlFkvIyKCsPtXpl5VBZAQd+eT0A5ry1gKWrV7H7FmH/pffnsaJm5Zr05dE6kUjX1lWWVbL5qM0xM+YunEu917PpyE0xjDnz57CyZiVJT4aPK0uN3YBeA9hu/HYAPP3h0/Qu782MsTMAePjth8PvL6EJXrbrxw4cy45rh/9N/fX1vzJu4Di2GbcN7s5NL93UYm1h6n6TR0xm54k74+5c9d+r2GbcNmw9bmtW1Kzg2heubbXWcce1d2THtXdkec1yLv/P5Xxtg68xbcw0Pl7+Mb99/rc583acb23yLbYdvy3zls3jsmcv4/hpxzN5xGRe/vzlFq9fs8Y5detT2WzUZrz42Yv86v9+xc93+TkTB0/kn+/9kxtevKHFfFPry79yOcd8uJzFVUtgMfxh/z8wtM9Q7nj1Dv7y+l9aDRR/vdevGdJnCHe/cTf3vXUfN+0XpoW5+aWbeXre0xi2pqai6bUJS3DVXmFohDtevYM3vniDWTvPAuCmOTfxxhdvZM03db+BvQbyw+1+CMDtr97OytqVHLvlsQDcOOdGFqxa0GK5E5ZgZL+RHDL1kJD/kDvoX9mfr27wVQBuGXALVXVVLZbdMMYPGs+uk3YFYOyAuxk3cBwzxs7A3Zkw/M6cZceNCYMmsvnoqdTXO0Ouf4SpY9Zns/HrsWhCNYMu/XcIelNLfXpdFx3bc7sJ7L7hBN7vU8srB85hz6mT2G6DEQxdvornt5hLXZ2ll1qjvjZB7WqjKtpfu3I0O08cwvLXq/n3I/M4+8gxbD2qktl3ruSGXy0jV9CMG//3s1fYYvwyrnpwFKfcuTEf/+RexvZewdkPbs4v/r15RtoEZZak0hqotAZ6WT2VVs/srU5npC/kNx/uwS0L9+DZtQ8nUVfH1QsP5vHV29DLq6hMVtOroYpKr6YXNVRSSy9q6EUNZ/MLAP7FDnzGaA7iLwA8xwyWMGRN2mzr3ok6hvZe3bjDdWp7+PDmx1Kds6OaVQYMCEtdXeiQvuGGYV5dEREpirh9XJ8kNAFuYTjDTuk70fqJVNDaxO3ABYR+r0cAOQPXJatXcMQ9R7Sa6bsnfMo6w0dz5B/O4+YPz8fPC8HiqbdfyWNLbmj1+lT67S86jpdWP8iq8z8FYP3z9uGdxIM5r62sGU3NhSH9uDP2YSXzWXrJbAD6/2BLVg2ck/P6vsu2ZNXlYUTIvqfMYHDZWD69/B4Ayq6aSLLv/Ox/wHgCMMZV7cVHV9wGwB5XTGFKnz146eJfUlvr7HntmObXNVlvufBQzvrmBdTWOl9btCFLXzmJMw44iTlvLeCe1TtHgXbGt/OWwMpTIbgxftHxXLDBTF59fwFvjfsyq2p/yt7rf4P7F7/Kgi0OXNO0snEzyyQV5lTi9B95IbtMOoSGlS8z6OivMnLMTey49u58sP7fGXTcdyCziaI5TjK9jzN2vduZMXZP5m32EH2/fyCTJj7NlqO34P+2u57K75/YKFho8CT1NG79MGbCbDYePo2npv6Ot3f/HuOHf8roAaNJbnwezy08v5XfHLhiQEg/b+B5/NbP5zcDw/3/a+dxw9LWf/e+u29If/uHP+bBtx/k023D79Ilc87kobcfynnt6P6j+fQHIf1Zj5/F5ys/Z/bM8Ls38/6ZzJmf+3dv2uhpa9Kf8sgpjBs4jvsOuQ+Ao+49is9XfZ7z+q9t8LU1geuJD53IfhvuxzbjtgHgmPuOyXktwMwtZ4bAFefUR09l1k6z2Hrc1iyvWc6P/vGjVq+ftdMsdlx7R1bUrOBnT/2MMQPGMG3MNOavnM/Fz1zcYuCSOrbpyE3Zdvy2LFq9iD++/Ee+uv5XmTxiMp+t+Iy73rir1drFw6YeBsDiqsU8+9GzrKwNg+R8sfoL5syf0yhQy3Z9bUPo6t/gDXy0/KPwJUN0/ZtfvNlqU+G6ZB0A7y99nyfef2LN5/Lqgld5+J2HcwbNwJrA9ZmPnuHRdx9dE7je/9b9ja7PFvyPGTBmTeD6x5f/yOcrP18TuF793NW8NP+lnM9u2uhpawLXi565iPEDx68JXM/4xxl5/e6lAtfvPfg99t9o/xC44hz014NyXgtw/JQj+d32vyBRXcXhj+zNrI1PYLMFR1O9/BNO/O9+rV4/4OFd2f2+7amo/YIb+/yWbW7Zie2WrUc5n3Dv+Edavf6TCwbBgeWMHlLF4sNW88XBwNswaUPg1FYv5+Pf9WLrt2voveEIOGMBH11VztjP6klOXxvO/LBR2gagKlpSXnt6SwbU9+aZTf/Nf7c4i9fn7sawsv7cyov8d72fp/+NgKxfHI4echL9KgZxxitvMG/EAywZfz51iQQ/ffEJlg59pvH1mTx8cXzIpB9Sl6znwdcfo6bPh3x5zP7UJ+v510ePU9fr3ejSjP9XW/QtNU6ZVzK2YlMavJ7P+74DiSfZ+OPd+R9Htv7BiYhIq+L2cV0CNLj7WsUvUvGZWV/CPLMJ4Ax3v7SFdL8Bvkfouzs61z1HjJnqC3tfyZr+cZZMb2esP3hiN9Ye25vvnjuX6/42l+Tr+wNwwPdf4IF/vx/6xJWFfnJl0TqR2i93Fvwz/AF61E+fYfbbH/PKbeEPn6+f8SgvfzCPRJljieSa61J97hIJZ3C/vjz5q/BH+lEX382KmlX89afhfgdeeDOfLpsf8o761WHJRn32xg8ezU0nzgTgmN/9hiH9BnBZNDLMwb+5gFV1y0Owl+qbl7HtlmTTUZP52b7h+qNuPYvNR0/llF0PJZl0DvrT8RnBYpMAEoeEs+vEXZi51VEkk86hfzuUb27yDb65yTdZXLWY4x84Pusfr6mBItydQ6ceyqFTD2Vx1WKOue8Yvjvtu+yx3h68u/hdznzszKy1LZnHjtr8KHaZtAsfLP2A8586n1O2PoXNRm3Gy5+/zDXPXdNq8HDctOOYMmIKry98nRvn3MgpW5/C+EHj+e/H/+XuuXe3GjwcN+04xgwYw4ufvcjDbz/MyVufzIBeA/j3vH/z7EfPtlrrdcwWx9Cvsh/Pf/I8sz+dzfdmfA8INZhzv5ibM28z47BNw+/Kc588x2crPmO/jcIfzc/Me4b5K+e3WGNlGL3Le7PbOmEUohc/e5Hahto1geNznzzHqtpVLZY9YQn6V/Zn6sjQQOK1Ba9RWVbJ+sPWB0LwU5+sz1n2/pX9GTcwzI/5wdIP6FfRj+H9huPuzFs2L2fZzYw+5X0Y0GsA7s6S6iX0Ke9Dn4o+JD2Zs+ypdapPm7uvadba1ew8J3y58OQWW3RwSeJLTRlQngjf0dbWVpGsqSZZXRUG7KmtbbyuqSFRV8/AhnKorWXJ6kUkauoY1FAONTV8Xv0FDTXVeG0NybpavKYGr6vDa2vCPepq6VvdwKjVCait5V1bQr/qJKNWOF5Tzeu9V+J1teGaulq8tpZkfS319XXUJaA+AYOrYewKqEvAc2NhrVUwdHUZK8oTvDTaqEmUUZcwahMJ6hKJsG0J6hMJas0YvLqCQat7sbpXGW+OqWVIVX/61fdleW9n3rDV1CeM+rIE9WUJGhLQkEhQn3AaEkZDwuhfN5DKZD+qy2tZ2Gch/euGkfBeVJWvZnnlIhrMcXOSJEmuWYdtJ0lFcgB4OXWJKurKllLW0A8M6hNVNCRW4yRxGsCSYduSQDL8G9pZXpWG8ig2NajvC8lyKKuBREOjL2ebfeHaUAFfbAReDoM+gESSEcv25vMbr+nIn0ZEpKTas6lw3MD1eWBTYKC71xS9VEVmZtOA2dHu3u7+cAvpTiCMjgwwLGr2nNVmm033P/5xNmVlYeTN8nKybg8dGloR1dWFJrAamV9EuoLYgWsy2XiY4abDD0fHkjXV1Fevpra6mpqqGmqqqqmprqGmupqa6lpqamqpramjtq6Wmto6amtrqa2rp6aunrqGOmrrGqirr6O2oYG6hvqwTjZQ29BAvTt1yXrq3KnzBuo8Sb0nqUs4DWZRkGjUJ5z6hNGQcOoT0BAFj2Hfo8DOSZrTkAhLMuE0JJLR2vGE01AW9j2RJJlIkkw0rNn2RENYrCG9nagPgVuiPgRDnUEyEQK0ZDl4WXq72VIGXo4lyzHPWCjDvJwE5c3WCcpIWGq7nDLLXJdRZuWNl0RYl1sZ5YmwX54op9zKo/0yKsqiY2XlVCTKqSgrj46VrdmuKA/ryrJyKsvLKS8ro7K8PPtSkTqXoKLCGv1bnm3d0rlEInRzFRHpKTp9H1dCs9ppwIHAH4tXnJIZk7H9SY50mefGAC0GrnM/fot9L/tyoeUSEemU5u91GOCMveI0ktZAAw0hICNJ0hpCLVsiOhYFZUlLBWohKPNEAyTSARuJBkjUh6W1kX0ro6XUGloO1EJQVpaxXY5RjnkvElGw1jhAK6fcy0k0lJFoKF8TrGUGamGevsZBWvma7RCorQnSUsFZIh2QlWcEaRVlIRhLBWeZwVqv8nIqysO5XhXla9aV5dF+eVnGfhmVFYm8grVUd04REZH2FjdwvRL4OvBrM1vk7rk7vHW8ARnbq3Okyzw3oOlJM5sJzASwEZV8OrjDB0sWESmJ+vIqwFjU/zPMy0h4GeYJEslyEl5GmfeiIlkWgjQvC7VqlFEWBXNlVh72UoFaKkBLVEQ1aBWUlZVRXlZBeVllWJdXUFFeSXl5JeXlFZRXVFJRUZmuUSuPateigKwyVaNWHtWqVZRTmRGspYK0zBq13hUhXVgnKC+3RsGZAjMREZHOKW7gejbwL8Iovfeb2WvAM8ACwngLWbl76yPKdGLufh1wHcD06dN99lWzW7lCRKRrWtNU+FtvdnBJREREROIHrrNYM+ELAFOAyXlc11GB64qM7b450mWeW9FiKhEREREREWk3cQPXf9FoPPhO79OM7bFAS7OBj23hGhEREREREekgsQJXd9+5yOUotTcI884mCLXDWUcVjs5BmA6nxYGZREREREREpP30iGEo3H01oQ8uwJ7Z0liYbHGPaPfv7VEuERERERERaV2PCFwjN0frXcxs6yznvwWsE23f0j5FEhERERERkdYUJXC1YC0zm1CM+5XIzcArhAGl7jKz3QDMLGFm3wKuj9I97O6PdVAZRUREREREpAlzjz/GkpltCZwL7A70A9zdyzPODwEuIgzkdLK71xZW3MKY2UTgCWBidGg1IXjvHe3PAXZz9yV53GsFoHkiOre1gC86uhDSIj2fzk/PqHPT8+n89Iw6Pz2jzk3Pp/Pb0N0HtEdGcUcVxswOB24AKlpK4+5LzGwSsBtwP/Bg3PyKwd0/MLNNgR8CXwcmAXXAa8BtwNVtCK7fdPfppSmpFIOZzdYz6rz0fDo/PaPOTc+n89Mz6vz0jDo3PZ/Oz8xmt1desZoKm9nGhKa1FcBVwHRa/jbkFkLz3P3i5FVs7r7C3c9z96nu3t/dB7r7dHf/ZUfXCIuIiIiIiEhzcWtcTwcqgWvc/VQAM2toIe3j0XrbmHmJiIiIiIhIDxZ3cKZdCf1WL24tobt/SuhL2pkHborjuo4ugLRKz6hz0/Pp/PSMOjc9n85Pz6jz0zPq3PR8Or92e0axBmcysyqgzt0HZhz7DBjh7mVZ0i8EBrl7ZSGFFRERERERkZ4nbo1rDVBpZtZaQjPrAwwGlsXMS0RERERERHqwuIHrB4SBmdbPI+3eQBnwesy8REREREREpAeLG7g+Qhgp+JRcicxsGHAJoT9sh06Fk8nM+prZXmZ2rpn9zcw+NDOPllkF3ntWxr1yLesV6cfplsxsmJkdZWZ/MrPXzWyVmdWY2cdmdo+ZHVCEPEaa2S/N7E0zqzKzxWb2tJkdm09rgp6slM9H71BxmNmWZnaemd1nZnPNbJGZ1UXrZ8zsHDMbWmAeeodiKuXz0TtUOmZ2ZuZnWOC99P6UQDGekd6h4jCzI/P8HHcvIA+9RwUo5TMqxXsUd1ThXwEnAN81swXA5U0K2gc4ALgAWJswVc7vYuZVClsBD5U4jzpgcY7z9SXOv6ubT+Pfz2rCZzo2WvYzs4eBb7r76rbe3MymAY8Cw6JDK4EBwJei5Vtmtq+718T/Ebq1kj6fiN6hwhwNfD9jvxqoAoYC20XLqdHv+X/aenO9QwUr6fOJ6B0qIjPbEDivSPfS+1MCxXxGEb1DxZEEFuY4H+v3XO9RUZXkGUWK9h7FqnF198+BQ6OC/JTwgw4DMLPXosL9kRC01gCHuPvyOHmV0BLgMeBS4BDCH+LF9Ky7j8qxfFDk/LqbcuA5whck67p7H3fvD0wCfh+l2Qu4tq03NrNBwAOE39m5wAx3HwD0A04k/F5/hfAFjWRXsueTQe9QYZ4DfkSYimxI9IwGEv5RP5Lw/+21gHuidyJveoeKomTPJ4PeoSIxswTh/229gbhfJKTupfenBIr5jDLoHSqOj1r5HJ9u6w31HhVd0Z9RhqK9R3GbCuPuDwA7Ai8Q5nQtJzQf3hjoFW3PAXZ098fi5lMiT7v7UHff3d3PcPfbKeybBCm+Xd19a3f/rbu/lzro7h+4+7GkA6LDzGx8G+/9Q2AUoXZjb3efHd271t2vIf1t7Uwz26CwH6PbKuXzkSJw91vc/TJ3/z93X5pxfKW73wwcFh0aAezTxtvrHSpQiZ+PFN9JwPbArcDfC7yX3p/SKOYzks5P71EPFDtwBXD359x9K2BzwrcbFxLmdj0N2Mrdp7n78wWXssjcvaGjyyC5ufsTrST5fcb29Dbe/ohofbu7v5/l/NWE5iZlwLfbeO8eocTPR9rH/2Vsj2vjtXqHSq+Q5yNFZGaTCF2fFhH+vimU3p8iK8Ezks5P71EPVFDgmuLuL7v7b9z9XHc/y92vTH3zIVIi1RnbzeYObknU/2VCtPtwtjTuvhJINYn4SqzSSaznI+1qh4ztd/O9SO9Qu4n1fKQkric0Pzzd3XP1AWuV3p+SKdozks5P71HPFStwNbO3LIzIO7HI5elOJpvZq9EIZyuj0c6uN7MtOrpg3cTOGduvtOG6KRnbr+ZIlzq3SRvuLWk7Z2y35flk0jtUZGbWy8wmmtmJhHEIAN4B7m/DbfQOlUiRnk8mvUMFMrPjgN2Af7r7LUW4pd6fIivBM8qkd6g4hpvZC9FnWGVm71mYlWDnmPfTe1R8xX5GmYr2HsWtcV0P+Bnwrpk9aWbHmNnAmPfqrtYi9PddTejzuwFwLPCCmf2/jixYV2dmg4Gzot2n3f3NNlw+JmP7kxzpUucGmln/Nty/xyvw+WTSO1QkZlZtYVqIauB9QhOqIcAzwG7ethEX9Q4VWZGfTya9QwUws7GEARyrgOOLdFu9P0VUomeUSe9QcfQFtgRqCbHHJELz3SfM7EYza+ssJ3qPiq/YzyhT0d6juIHrBcAHhAGYdgSuA+ab2W1mtnc0sltP9TZwBrAh0NvdhxGar+xBGMjKgHPM7AcdV8SuK/rd+iMwmjCg1kltvMWAjO1c07RknhvQYipppAjPB/QOlcJ84HNgVcaxJ4BT3X1eG++ld6j4ivl8QO9QsVwLDAJmZQ5CVyC9P8VVimcEeoeK5VNCRddmhM9xKCFA2h74Z5TmKNo+8q/eo+Ip1TOCUrxH7h57iX6oawnT3ySBhmiZT5jbdYtC7t+eCyEQd8L//EqVR2/CFAgOrAAGdfTP3dUWQk2ER8vRMa4/O+P68hzpjstIN7qjf+6ushT6fPK4v96hwj/DEcAPMv6/fX4br9c71ImfTx731zuU3+d0WPQZzWn6ew7MSv1ux7iv3p9O/ozyyFfvUHE+xwRwT/Q5NgDrt+FavUed/Bnlce9Y71Ghowo/4+7HE4aj/iahH0494R/eU4DZUZvmH5nZmBy36hHcvZrwsgH0J/TJkDyZ2WWE0asBTnP3G2PcZkXGdt8c6TLPrWgxlaxRpOeTk96hwrn7Anf/JbAn4R+Mn5hZW6Zb0TtUQkV4Pq3dX+9QK8xsBHAF4Q+149y9voi31/tTBCV+RjnpHSoOd08SprSBECB9rQ2X6z1qBwU+o9buHes9KtaowrXu/jd335/QRPBEQhRthA7RFxFqNKXxpNjrdFgpuhgzu4RQCwHwI3e/IuatPs3YHpsjXerccg8j00kORXw++dA7VATu/hzw72h3Zhsu1TvUDgp4PvnQO5TbxcAwQjeouWbWP3MhzF0PQMbxypZu1oTen+Io5TPKh96hInD3d4Avot22fI56j9pJAc8oH21+j4reF9XdF3uYGmdbYDIwmxDAakoMicXMLgV+FO2e4e6XFXC7zNHnprSYKn3u9QLy6hGK/HykfaUGrlivDdfoHWo/cZ6PFG5StP4eoZam6XJWRtrUsUvyvLfen+Io5TOSzk/vUQ9VkkGUzGwrM7saeAqYVoo8urBtMrazTZgsGaLmp6lmCme4+6WF3M/DCLepwU72bCHPfqTnUPx7Ifl1d8V+PnnSO1Q8qW84825CpXeoXbX5+eRJ71AH0fvTbegdKgIzW5cw4iy04XPUe9R+4j6jPLX5PSpa4Gpma5vZOWY2l1D1ewLhB60D/gYcUKy8Oiszs1bO9yKMyAxh9MjHSl6oLiwKilLNT39YxKAoNc/bwZZ9LuLvE9rbNwC3FinPbqcUz0fvUHGYWVken+VuwFbR7pNtzELvUAFK+Xz0DhXO3Xd2d2tpIYzAmUqbOn5qG7LQ+1OgUj4jvUPFkcfnaISpjCAMRPdAG7PQe1SgUj6jkr1HBY4INZAwD89ThEGZGqIfLEmYf+54YHCxR7kq0mhWQwiBdWqZRxgI45Imx/s3uW4W6RHKJjY5txNh6OjDgHEZxysInY6fy7j2jI7+DDrzQui/kvqsTmvjtS0+o+j8IOCz6PxrwLToeCWh2VFNdO43Hf05dNalVM9H71DRns9E4KXo/8HrAJZxbjxwJrAy+hwXAaPyfUbReb1DnfT56B1ql+e35vNv7bzen673jPQOFe0ZTIw+q0b/nyNUmm0DPJLxOTb7Xdd71LWfUaneo7g/6D7AHYT5kTKD1XcJ33Kt29EPI4+f4YOMDyzX8oc2PKSdm1y7GlhImMw3dawBuKCjf/7OvAATmnxe81tZfpjvM8pIM43Q2TyVbnmT5/Qo0KujP4vOuJTy+egdKtozmtjkc6yJPseVTY6/R5Zpy/QOdd3no3eoXZ7fms+/tfN6f7reM9I7VLRn0PT/c9XR51jd5PiNZJnSRu9R135GpXqPyonnvihDA5YBdwK3uPu/c17V/b1C6O+3LTCVUGM7mPCwXgeeBq5z91c6qoBdRKLJ9shW0vdvawbu/oKZTQZ+TPgiZjyhqcKrwM3AjR6GAZfmSvl89A4Vx6fAgYR/OLYmjPa+FuEfiXnA/4B7gT+7e1WcDPQOFaSUz0fvUBeg96dT0ztUHJ8DJxE+x82B4YTWjtWE/ozPEn7Pn4mbgd6jgpXyGZXkPUpVCbeJmdURqo9vAe5z95o230REREREREQkD3FrXMe4+8KilkREREREREQki1ijCrclaDWzhJl9zczuiZOXiIiIiIiI9Gxxa1xbZWYbAEcDR9B6HzgRERERERGRrIoauJpZX8KAE0cD26cOR+s3ipmXiIiIiIiI9AxFCVzNbBtCsHoQ6RFEDZhLGHH4Tnd/tRh5iYiIiIiISM8SO3A1s+HA4cAxwEapw9HagRnu/kJhxRMREREREZGerk3T4ZiZAXsRgtV9CIGvAVXAPYQ5kx4hBK4D3H11kcsrIiIiIiIiPUxeNa5mti6hKfB3CBOlGyE4/TdhLte/uPuKKG1pSioiIiIiIiI9Ur5Nhd8mBKoGvAf8EbjF3d8vVcFEREREREREoO19XK8CznD32lIURkRERERERKSpRJ7pagm1rScBn5rZNdFIwiIiIlJCZvb/zMzN7IyOLkscZvZIVP5dO7osIiLSdeUbuI4CTgZeBoYC3wOeMbM3zexsM5tQqgKKiIj0VGY2DjgdWAhc0475Xh8Fm4vMrFcbrnsnuu6+jMOzovVlZpbv3x0iIiKN5PUPiLsvdfdfu/sWwDTgt8AyYH3g58B7Zva4mR1VuqKKiIj0OBcAfYBL3H1VO+b7+2g9FNgvnwvMbCdg3SbX4+7/BzwKbAEcVsQyiohID9Km6XAaXRi+gf0mYWqcnUiPNJxafwN4wN3ri1NUERGRnsPMxgIfAg3AaHdf3M75vwZsAjzi7nvlkf4PhNkHPgfGZf77b2b7A3cDr7n7lJIUWEREurXYTXbcvcbdb3X3XYH1gAuBT6LTBtwFLDCzm8xsbzNr60BQIiIiPdlxQBnwUHsHrZFUrelXoibLLTKzAYQvsyHMOtD0S+uHgMXAZDPbvrjFFBGRnqAofU3c/X13PxdYG9gb+BtQDwwGjgDuJ3wDKyIi0qmZ2T5mdqeZfWBmq6M+my0tx5SoDEZo0QTw5xzpnozKMcvMyszsNDObY2YrzWyBmd1jZptlpO9rZuea2atmtirqw3pHNF97U38E6gh/K3ynlSIfBPSLtm9sejKajeCuaHdmK/cSERFpJnZT4VZvbLYWIWg9CpgMuLuXlSQzERGRAplZH+B2YN/okBPGc+hPy9PHzXD32SUoy1TCgIgQmgnPbyHdk4TuOhcCWwG7E2YCqCMdSK4EdgHeB/5B6GtaTfj5+kRpFkQ/y7wm9/8roevPO+6+fo7yPgNsBzzj7l9qIc1hhGB4gbuPbOleIiIi2ZRsdD93/8LdL3f3qYR/zH7f2jUiIiId6M+EoLWBMPDgKHcfAgwATstI91/Cl7JHAa+UqCw7RuuPWgpamziBEJB+ixBoDyAEsu9F+1cC1wNDgD0IQW1/QqC7EBhBCH6bSv3bvZ6Z7ZjlPGa2IeHfechS25rhv9F6hJltlMfPJCIiskbJalxFRES6CjM7HLgl2j3W3Zt92WpmtxOaxC5y97VKXJ5bgMMJgxx+LUe6Jwk1rgA7uPu/m5zfFXgs2q0CNnX3d5qkOZoQoFYBg9y9LuNcgjBA1DjgZnc/MksZLgbOINTsjnb3lTnKu4IQMB/t7je1lE5ERKQpzacmIiI9WtSf9Nxo955sQWvk/mg9zMxGl7hYY6L1wjzT/7tp0Bp5CqiJtv/aNGiNPBqt+xCmuVvD3ZPAzdHuN82sf+Z5MysjBNgAd+QKWiOLovWYnKlERESaUOAqIiI93bbABtH2ZTnSLcrYrkhtmNnh0QBJzQYwMrOpZrY8GgxpYBvKNDxa5zua8HPZDrp7A/BFtPt8C9dmDp44JMv5Gwn9YfsRapwz7QWMzkjXmtTPMzxnKhERkSYUuIqISE+3W7SeD/wnR7oR0ToZpU25jdCc9oyo9hYAMxsJPEBogruPuy9vQ5l6R+uanKnSVuQ4V58rTZOpayqynH8PeDLaPbrJ6dT+XHd/tvViUhWte+dMJSIi0oQCVxER6emmRuvnoqaxLUn1JX05mt4FWBP4/RLYBPgagJn1Bu4lBLv7u/sHbSxTqnY3Ww1oR0g1n94uGowpNXvAPk3Ot2ZotF6UM5WIiEgTClxFRKSnGxutv2gpQdSXMzVNzkNZkvw+uv7HUa3rzYRRfY9y91y1uC1J9W0dmjNV+7kLWBptHxWtDyfU0NYTprnJR+rnybfvroiICKDAVUREJNVUdliONMcDaxGmyrm+6Ul3Xw1cTZgW5m7gQOA8d789Zplej9brxLy+qNy9mjBdEMARUSCfCmAfcPfPs1+ZZmYDCJ8hwBvFL6WIiHRnClxFRKSnmxutdzGzwU1PRnOOXhLtXpWj2e+vCVPC7Afc6u4/L6BM/4rWm5lZrwLuU0yp5sCjgZ+QbmKdbzPh6YS/O+qBZ4pbNBER6e4UuIqISE93a7QeCNxpZhMg9FON5nd9mjCi7jPAOTnusxFQHm1fVWCZniEEeJXA5gXeqyjc/UXgpWj3J9H6M+DhPG+xdbR+MY9pc0RERBpR4CoiIj2au/+LUFsKsDvwoZktIYzCewuheevdwN7uXpXtHma2DnAPIZCrA04tsEzLgQej3X1zpW1nqdrV1N8PN0dT7uQj9XP8OWcqERGRLBS4iohIj+fuJwHfBp4iDELUlzDlzV+APdz96y1NZxM1L36QUNu6F6EG90Azm1Rgsa6N1odmTrPTwW4FqjP285m7leiz2JYwHc4tJSiXiIh0c+buHV0GERGRLsnMKghNZXcEvuLuT0Z9Yl8Hfuvu3y/g3gngLWBdYKeoZrhLMrOfAj8DbnL3pnPBioiItEqBq4iISExmdgNwDHC0u9+UcfxuYA9gorsvKOD+hxCa1j7s7nsXWt6OYGb9gA+AAcCG7v5hx5ZIRES6IjUVFhERicHMziQErRdlBq2Ri4A+wMkFZnM78Bywl5lt3VriTupEQj/hqxS0iohIXKpxFRERaSMz+yah/+tdwIGe5R9TM3sS2BRY291XFJDX5sD+wPPu/mDu1J2PmZ1ImCP3Cndf1tHlERGRrkmBq4iIiIiIiHRqaiosIiIiIiIinZoCVxEREREREenUFLiKiIiIiIhIp6bAVURERERERDo1Ba4iIiIiIiLSqSlwFRERERERkU5NgauIiIiIiIh0agpcRUREREREpFP7/4GgDUZaFjclAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 1080x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "csvfile = './data/FP_Jai_sup_0_2000A.csv'\n",
    "df = pd.read_csv(csvfile)\n",
    "\n",
    "Lsig = df['sig']\n",
    "re0 = df['x0']; ri0 = df['y0']; ra0 = df['z0']\n",
    "rei = df['xi']; rii = df['yi']; rai = df['zi']\n",
    "re1 = df['x1']; ri1 = df['y1']; ra1 = df['z1']\n",
    "\n",
    "ini = 0; end0 = -4\n",
    "plt.figure(figsize=(15,4))\n",
    "plt.plot(np.array(Lsig[ini:end0]),re0[ini:end0]*1000,'r')\n",
    "plt.plot(np.array(Lsig[ini:end0]),rei[ini:end0]*1000,'-.r')\n",
    "plt.plot(np.array(Lsig[ini:end]),re1[ini:end]*1000,'r')\n",
    "\n",
    "plt.plot(np.array(Lsig[ini:end0]),ri0[ini:end0]*1000,'b',label='stable')\n",
    "plt.plot(np.array(Lsig[ini:end0]),rii[ini:end0]*1000,'-.b',label='unstable')\n",
    "plt.plot(np.array(Lsig[ini:end]),ri1[ini:end]*1000,'b')\n",
    "\n",
    "plt.plot(np.array(Lsig[ini:end0]),ra0[ini:end0]*1000,'g')\n",
    "plt.plot(np.array(Lsig[ini:end0]),rai[ini:end0]*1000,'-.g')\n",
    "plt.plot(np.array(Lsig[ini:end]),ra1[ini:end]*1000,'g')\n",
    "\n",
    "plt.axvline([3],color='c')\n",
    "plt.plot([1.6,4],[10,10],color='darkgray',linewidth=6)\n",
    "plt.ylim(-0.5,14)\n",
    "plt.xticks(fontsize=25)\n",
    "plt.yticks(fontsize=25)\n",
    "plt.xlim([1.5,5.5])\n",
    "plt.xlabel(r'$\\sigma_X$ (mV)', fontsize=25)\n",
    "plt.ylabel('Average rate  $r_X$ (spks/s)', fontsize=25)\n",
    "\n",
    "plt.savefig('./figures2000A/FP_Jai_sup_0_2000A.pdf.pdf',bbox_inches='tight')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.8.3"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}