Skip to content
Snippets Groups Projects
plots_fig3.ipynb 71.3 KiB
Newer Older
{
 "cells": [
  {
   "cell_type": "markdown",
   "source": [
    "# Create plots for panels used in figure 3"
   ],
   "metadata": {
    "collapsed": false,
    "pycharm": {
     "name": "#%% md\n"
    }
   }
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "collapsed": true,
    "pycharm": {
     "name": "#%%\n"
    }
   },
   "outputs": [],
   "source": [
    "import os\n",
    "import numpy as np\n",
    "import pandas as pd\n",
    "import matplotlib.pyplot as plt\n",
Felix Reichel's avatar
Felix Reichel committed
    "import matplotlib.patches as patches\n",
    "import seaborn as sns\n",
    "from scipy.optimize import curve_fit\n",
    "\n",
    "import warnings\n",
    "warnings.filterwarnings('ignore')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "outputs": [],
   "source": [
    "# folder to save all panels for figure\n",
    "savefolder = r\"plots\\fig3\"\n",
    "\n",
    "# file containing the data for the controls\n",
    "results_ctrl_file = r\"data\\shape_analysis\\histograms_HealthyControl_deformed_undeformed.txt\""
   ],
   "metadata": {
    "collapsed": false,
    "pycharm": {
     "name": "#%%\n"
    }
   }
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "outputs": [],
   "source": [
    "#define a color seed for each patient\n",
    "color_dict = {'VS': 'C0', 'VL': 'C1', 'RS': 'C2',\n",
    "              'KM': 'C3', 'LM': 'C4'}"
   ],
   "metadata": {
    "collapsed": false,
    "pycharm": {
     "name": "#%%\n"
    }
   }
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "outputs": [],
   "source": [
Felix Reichel's avatar
Felix Reichel committed
    "def asymptotic_exponential_growth(x, lambda_):\n",
    "    \"\"\"(Inverted) exponential growth function with maximum at 1 for x->infinity:\n",
    "    f(x) = 1 - exp(-lambda * x)\"\"\"\n",
    "    return 1 - np.exp(-lambda_ * x)\n",
    "\n",
    "def deformed_probability_curve(df, v_min=0, v_max=3, binsize=.25):\n",
    "    \"\"\"Compute the values for the shape probability diagram to find a cell\n",
    "    in a deformed state for velocities between v_min and v_max in the DataFrame df\n",
    "\n",
    "    returns: *tuple* (deformed_bins, deformed_hist_normal)\n",
    "        - deformed_bins: *array* limits for the bin ranges of the histogram\n",
    "        - normalized counts for each velocity range\n",
    "    \"\"\"\n",
    "\n",
    "    bins = int(v_max/binsize)   #number of Bins in histogram\n",
    "    #find index of cells in a deformed state. Class definitions are:\n",
    "    #1-parachute, 2-slipper, 3-asym. parachute, 5-multilobe, 7-undefined deformed\n",
    "    #4-discocyte/undeformed, 6-tumbler\n",
    "    deformed_index = ((df['shape'] == 1)\n",
    "                      | (df['shape'] == 2)\n",
    "                      | (df['shape'] == 3)\n",
    "                      | (df['shape'] == 5)\n",
    "                      | (df['shape'] == 7))\n",
    "\n",
    "    #create new column in df that is True for deformed state\n",
    "    df['deformed'] = False\n",
    "    df['deformed'][deformed_index] = True\n",
    "\n",
    "    df_deformed = df[deformed_index]\n",
    "\n",
    "    deformed_hist, deformed_bins = np.histogram(np.array(df_deformed['velocity']),\n",
    "                                                range = (v_min,v_max),\n",
    "                                                bins = bins)\n",
    "    #get the counts for all events to use for normalization\n",
    "    all_hist, all_bins = np.histogram(np.array(df['velocity']),\n",
    "                                      range = (v_min,v_max),\n",
    "                                      bins = bins)\n",
    "\n",
    "    #normalize the deformed histogram\n",
    "    deformed_hist_normal = deformed_hist/all_hist\n",
    "\n",
    "    return deformed_bins, deformed_hist_normal"
   ],
   "metadata": {
    "collapsed": false,
    "pycharm": {
     "name": "#%%\n"
    }
   }
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "outputs": [],
   "source": [
    "#define dict to store fit values\n",
    "dict_fitvalues = {}\n",
    "\n",
    "def dict_fit_values_patient(patient, dict_fitvalues):\n",
    "    result_summary_folder = r\"data\\shape_analysis\\result_summaries\"\n",
    "\n",
    "    v_min = 0.\n",
    "    v_max = 3.\n",
    "    binsize = 0.25\n",
    "\n",
Felix Reichel's avatar
Felix Reichel committed
    "    # bounds of the parameters in the exponential growth function\n",
    "    fit_bounds=(0, np.inf)\n",
    "\n",
    "    result_file = os.path.join(result_summary_folder, patient + \"_results_MCFM.tsv\")\n",
    "    df_results = pd.read_csv(result_file, sep='\\t')\n",
    "\n",
    "    dates = np.unique(df_results['date'])\n",
    "    dates = np.sort(dates)\n",
    "    day0 = pd.to_datetime(dates[0])\n",
    "\n",
    "    #create dataframes to save fit parameters\n",
    "    df_fit_all = pd.DataFrame()\n",
    "    df_fit_healthy = pd.DataFrame()\n",
    "    df_fit_unhealthy = pd.DataFrame()\n",
    "\n",
    "    for num, date in enumerate(dates):\n",
    "        df_date = df_results[df_results['date']==date]\n",
    "        #create new Dataframe to work with, leave out skipped cells\n",
    "        df = df_date[df_date['shape'] != 0]\n",
    "\n",
    "        healthy_index = df['health'] == 0\n",
    "        df_healthy = df[healthy_index]\n",
    "        unhealthy_index = df['health'] == 1\n",
    "        df_unhealthy = df[unhealthy_index]\n",
    "\n",
    "        #calculate percentage of healthy cells in sample\n",
    "        percentage_healthy = len(df_healthy)/len(df)\n",
    "\n",
    "        bins, deformed_curve = deformed_probability_curve(df, v_min=v_min, v_max=v_max, binsize=binsize)\n",
    "        bins_healthy, deformed_curve_healthy =  deformed_probability_curve(df_healthy,\n",
    "                                                                           v_min=v_min, v_max=v_max, binsize=binsize)\n",
    "        bins_unhealthy, deformed_curve_unhealthy =  deformed_probability_curve(df_unhealthy,\n",
Felix Reichel's avatar
Felix Reichel committed
    "                                                                               v_min=v_min, v_max=v_max, binsize=binsize)\n",
    "\n",
    "        bins_plot = bins[:-1]+binsize/2\n",
    "\n",
    "        #exclude nan values before fitting\n",
    "        ind_nonnan_all = ~np.isnan(deformed_curve)\n",
    "        ind_nonnan_healthy = ~np.isnan(deformed_curve_healthy)\n",
    "        ind_nonnan_unhealthy = ~np.isnan(deformed_curve_unhealthy)\n",
    "\n",
Felix Reichel's avatar
Felix Reichel committed
    "        x_all = bins_plot[ind_nonnan_all]\n",
    "        y_all = deformed_curve[ind_nonnan_all]\n",
    "        x_healthy = bins_plot[ind_nonnan_healthy]\n",
    "        y_healthy = deformed_curve_healthy[ind_nonnan_healthy]\n",
    "        x_unhealthy = bins_plot[ind_nonnan_unhealthy]\n",
    "        y_unhealthy = deformed_curve_unhealthy[ind_nonnan_unhealthy]\n",
    "\n",
    "        popt_all_exp, pcov_all_exp = curve_fit(asymptotic_exponential_growth,\n",
    "                                               x_all, y_all,\n",
    "                                               bounds=fit_bounds\n",
    "                                               )\n",
    "        popt_healthy_exp, pcov_healthy_exp = curve_fit(asymptotic_exponential_growth,\n",
    "                                                       x_healthy, y_healthy,\n",
    "                                                       bounds=fit_bounds\n",
    "                                                       )\n",
    "        popt_unhealthy_exp, pcov_unhealthy_exp = curve_fit(asymptotic_exponential_growth,\n",
    "                                                           x_unhealthy, y_unhealthy,\n",
    "                                                           bounds=fit_bounds\n",
    "                                                           )\n",
    "        #days since treatment start\n",
    "        treatment_days = (pd.to_datetime(date) - day0).days\n",
    "\n",
Felix Reichel's avatar
Felix Reichel committed
    "        df_fit_all = df_fit_all.append({'lambda': popt_all_exp[0], 'lambda_err': np.sqrt(pcov_all_exp[0,0]),\n",
    "                                        'days': treatment_days,\n",
    "                                        'percent healthy': percentage_healthy\n",
    "                                        },\n",
    "                                       ignore_index=True)\n",
Felix Reichel's avatar
Felix Reichel committed
    "        df_fit_healthy = df_fit_healthy.append({'lambda': popt_healthy_exp[0], 'lambda_err': np.sqrt(pcov_healthy_exp[0,0]),\n",
    "                                                'days': treatment_days\n",
    "                                                },\n",
    "                                               ignore_index=True)\n",
    "        df_fit_unhealthy = df_fit_unhealthy.append({'lambda': popt_unhealthy_exp[0], 'lambda_err': np.sqrt(pcov_unhealthy_exp[0,0]),\n",
    "                                                    'days': treatment_days\n",
    "                                                    },\n",
    "                                                   ignore_index=True)\n",
    "\n",
    "    dict_fitvalues[patient] = {'all': df_fit_all, 'healthy': df_fit_healthy, 'unhealthy': df_fit_unhealthy}\n",
    "\n",
    "    return dict_fitvalues"
   ],
   "metadata": {
    "collapsed": false,
    "pycharm": {
     "name": "#%%\n"
    }
   }
  },
  {
   "cell_type": "markdown",
   "source": [
    "Fill dictionary with patient data"
   ],
   "metadata": {
    "collapsed": false
   }
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "outputs": [],
   "source": [
    "patients = ['LM', 'KM']\n",
    "labels = [\"P4\", \"P5\"]\n",
    "\n",
    "for patient in patients:\n",
    "    dict_fitvalues = dict_fit_values_patient(patient, dict_fitvalues)"
   ],
   "metadata": {
    "collapsed": false,
    "pycharm": {
     "name": "#%%\n"
    }
   }
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "outputs": [
    {
     "data": {
Felix Reichel's avatar
Felix Reichel committed
      "text/plain": "<Figure size 792x432 with 2 Axes>",
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAwsAAAGoCAYAAAAJsXnTAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjQuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/Z1A+gAAAACXBIWXMAAAsTAAALEwEAmpwYAACj3klEQVR4nOzdd3hTZfsH8O/Japu26aALWqCF0pa9yxDkZYkgKDhAhihukCHOV1T8qS/gVqbgwMGSjaAIMlRAtsheZdOW0kF3mn1+f4QGQpKONE3T9vu5Li7b5zznnPtp6mnuPEsQRVEEERERERHRHSRVHQAREREREXkmJgtERERERGQXkwUiIiIiIrKLyQIREREREdnFZIGIiIiIiOxiskBERERERHYxWSAiIiIiIruYLBARERERkV1MFshjrVixAvHx8ejUqROmT58Ok8lU1SEREdV4fPZWruTkZMTHx+O///1vVYdCVCZMFshjRUZG4sknn4TBYMAPP/yAv//+2233jo+PR3x8PHr27AmtVmu3Tq9evRAfHw+DweC2uGoS/sEk8kxV+ey93Zdffml5Fl+4cKFKYnBGTX621eS2kWNMFshj3XXXXXj99dfx7LPPAgCOHDni9hhSU1Pxww8/uP2+RERVxROevaIoYuXKlRAEAQCwcuVKt8dARGZMFsjjtW7dGgBw+vRpt943ICAAgYGB+Oqrr3Djxg233puIqKpV1bMXAHbt2oWUlBQMGTIEoaGhWLt2LXQ6ndvjICImC1QNREVFAQDOnDnj1vt6e3tj7NixyM/Px9y5c8t17saNGzFy5Ei0b98erVq1wqBBg7BgwQKbP3a3d+levHgRL774Irp06YKEhATs27fPps6VK1cwceJEdOrUCW3btsWTTz6Js2fPAgBu3LiBt99+G926dUPLli3x0EMPYe/evRWKr9jRo0fx4osvonv37mjRogW6deuGJ598Ehs3bgQAnD9/HvHx8Xjssccc/kwGDRqE5s2bIz09HbNnz0bv3r0BAGvXrrUMNYiPj8eaNWuszjty5AgmTpyIu+66Cy1atECPHj0wdepUXL9+3eYe27Ztw+OPP45u3bpZ4hw1ahSWLFniMC4isq+qnr3ArZ6ERx55BIMGDUJ2dja2bt1a4jmlPaeKrVmzBhMmTEDv3r3RqlUrtGvXDo8++ih+/vlnm2ve/vxNTk7G5MmT0alTJ7Rs2RIPPvgg/vjjD6v65Xm2FV+/tGsWK+9zuyw/k/I8u6dOncrndi0liKIoVnUQRCWZMmUKVq9eDUEQ8M8//8DX17fS7xkfH4/w8HBs27YNAwYMwLVr1/DLL78gOjraUqdXr15ISUnBiRMnIJPJLOWfffYZFixYgKCgIPTr1w9KpRI7d+5EUlISEhMT8e2330KhUAAw/6Ho3bs32rdvj6SkJERHR6Nt27bQaDQYNmwYmjdvbqmTmJiIs2fPonHjxmjVqhVSUlKwZcsWBAQEYPny5Xj66afh5+eHjh07Ijc3Fxs3boQgCNi8eTPq1avnVHyAebLj//3f/0EikaBXr16Ijo5GVlYWjh8/DpVKhUWLFgEARo8ejX379mHTpk2IiYmx+nkeOnQIw4cPR79+/TBr1izs27cPW7duxY8//oiEhAT06dPHUrdPnz5o2rQpAGDVqlWYOnUqFAoFevXqhYiICFy+fBnbt29HnTp1sGLFCkvbli9fjqlTpyI0NBQ9e/ZEUFAQsrKycObMGYiiiNWrV7vot4OodqiKZy8AZGZm4j//+Q8iIyOxefNmnD17FoMGDULnzp0dDgst63MKAFq1aoXY2FjExcUhNDQUOTk5+Ouvv3D9+nWMHTsWL774oqXu7c/fc+fOoX79+mjTpo3lGWswGPDdd9+hc+fOAFCmZ1t5rwmU/7ldnp9JWZ/dI0eO5HO7thKJPNjOnTvFuLg4y79Dhw655b5xcXFi9+7dRVEUxd9++02Mi4sTX3jhBas6PXv2FOPi4kS9Xm8pO3TokBgXFyf26NFDTE9Pt5Tr9XrxueeeE+Pi4sQvv/zSUn716lVL2z799FO7sdxeZ968eVbH5syZI8bFxYkdO3YU3377bdFoNFqOrV27VoyLixOnTZvmdHxJSUlis2bNxI4dO4pnz561ie3atWuWr4t/Th988IFNvddff12Mi4sTd+3aZdOu119/3W67L1y4IDZv3lzs06ePmJaWZnVs9+7dYkJCgjhu3DhL2ZAhQ8TmzZuLmZmZNtfKysqyew8isq+qnr2iKIoLFiwQ4+LixPnz51vKhgwZIsbHx4uXLl2yqV+e55QoiuLly5dt6mi1WnH06NFis2bNrJ43tz9/Z8+ebXXOjh07xLi4OPHpp5+2Ki/t2Vbea5b3uV3en0l5nt18btdOHIZEHqugoABvv/02VCoV7r//fgBV0x1+7733om3bttiyZQsOHjxYYt3iT0HGjh2L0NBQS7lMJsPrr78OiURid6JeSEgIxo8fX+K1IyMjLRMOiw0ZMgQAoNPp8Nprr0EiufW/9KBBgyCTyXDq1Cmn41u2bBkMBgPGjRuHJk2a2MQUERFh+bpPnz4IDQ3FmjVrrLrF8/Ly8Ntvv6FBgwbo2rVriW283bJly6DX6/Hmm28iPDzc6liXLl3Qq1cv/PHHHygoKLBqx+29PMWCg4PLfF+i2q4qn73izYnNEokEgwcPtpQ/+OCDEEURK1assDmnPM8pAGjQoIFNHYVCgZEjR8JgMGDPnj02xyMjIzF27Firsu7du6NevXo4evRoWZvn1DWd+btSVc9uPrdrJttXh8hDfPzxx0hNTcUHH3wAuVyO9evXl3ui3ezZs23KhgwZYhmLW1avv/46Hn30UXz00Ud2/1gVO3nyJABYdR8Xi4mJQUREBJKTk5Gfnw9/f3/LsYSEBJsu5Ds1bdoUUqnUqiwsLAwAEB0dDT8/P6tjUqkUderUsRojWt74Dh8+DMD8B6w0MpkMQ4cOxdy5c7F582YMGjQIAPDzzz9Do9Fg6NChlpVNyqL43vv378exY8dsjmdlZcFoNOLSpUto0aIFBg0ahA8++AD33XcfBgwYgMTERLRr145/cIjKqSqfvXv37sWVK1fQrVs3qzebAwcOxAcffIC1a9fixRdfhFwutxwrz3MKMK9y9/XXX2PPnj24du0aNBqN1XF74+oTEhJsnr+A+U138f3Lq6zXdObvSlU9u/ncrpmYLJBH2rt3L5YvX47//Oc/GDJkCM6fPw+g/J9uzZkzx6YsMTGx3MlC27Zt0a9fP2zevBkbN27EgAED7NbLz88HAKtPf24XGhqK1NRU5OXlWSULISEhpcZwe/1ixZ/G2DtWfPz2fSDKG19x/Ts/IXJk2LBhmD9/PpYvX275g7NixQrI5XI89NBDZbpGsZycHADAt99+W2I9tVoNABgzZgyCgoKwdOlSLFq0CD/88AMEQUDHjh3x2muvoWXLluW6P1FtVNXP3uXLlwMw9yTcLjAwEL169cLmzZuxbds23HvvvZZj5XlOXb16FQ8//DDy8vLQoUMHdOvWDX5+fpBKpUhJSXG46pJKpbJ7PZlM5vSmdWW9pjN/V6rq2c3nds3EZIE8jlqtxltvvQV/f3+89957AMyfnHt7e+Ps2bMQRbHMn3K4suv85Zdfxvbt2/Hpp59aTey6XfGDOjMz025Xd0ZGhlW9YuX5xL0iyhtf8X+vX79u03NhT3h4OHr16oUtW7bg/PnzyM3NxdmzZzFgwIByf1JUfL9//vmnTPcGgMGDB2Pw4MHIy8vDv//+iy1btmD16tV4+umn8dtvv/HTKqISVPWz98aNG5YVj1566SW89NJLduutWLHCKlkoz3Pqu+++Q05ODmbMmGGTkPzyyy9Yu3ZtueOubM78XamqZzef2zUT5yyQx/nss89w9epVqzGPUqkUTZo0QUFBAZKTk6skroYNG2L48OFITk7G4sWL7dYpXg2ieNnT212+fBlpaWmIiopy+IlSZStvfG3atAEA7Ny5s8z3GDFiBADzJ4TFQ7aGDRtmU6+4+91oNNq9TvG9S5snYo9KpUKPHj3wv//9D0OGDEFOTg4OHDhQ7usQ1SZV/exdu3Yt9Ho9mjdvjocfftjuv+DgYOzevRtXr161nFee59Tly5cBAPfcc4/Nsf3797ukHaU928rLmb8rlfXs5nO7dmKyQB7l4MGDWLJkCXr27Gk1uQ0wj+8EqmaSc7EXXngBKpUK8+fPR2Fhoc3x4u7aL7/80mojN6PRiA8//BAmkwkPP/yw2+K9U3njGz58OGQyGebNm4dz587ZXC8tLc2mrEuXLoiOjsa6devw22+/ISYmxu5YW5VKBUEQcO3aNbuxjhw5EnK5HDNmzMDFixdtjut0Oqs/SHv37oVoZyXo4nZ6e3vbvQ8Recazt/gN6v/93/9h2rRpdv8NGzYMoihi1apVlvPK85yKjIwEYJsY7Ny50+qaFVHas628nPm7UlnPbj63aycOQyKPodFo8Oabb1p1gd+u+NOV06dPOxwGVNkCAwPx3HPP4eOPP7Z7vF27dnj66afxzTffYODAgejXrx98fHywc+dOnD17Fu3bt8dTTz3l5qidjy82NhbvvPMO3nnnHQwePBi9e/dGdHQ0srOzcfz4cfj6+lqtXw6Yh1QNHz4cM2bMAGC/VwEAfH190bp1axw8eBAvv/wyYmJiLOuBJyQkoHHjxpg2bRrefPNNDBw4EN27d0d0dDQMBgNSU1Pxzz//ICgoCJs2bQIAjB8/HkqlEm3atEFkZCREUcTBgwdx7NgxNG/evFwrMRHVJp7w7N23bx8uXbqEuLg4tGrVymG9hx9+GPPnz8fq1asxYcIEyGSycj2nRowYgTVr1mDSpEno168fwsLCkJSUhJ07d6J///42G7g5o7RnW3k583elsp7dfG7XTkwWyGPMnDkTly5dwocffmhZ5ed28fHxAGDZsbiqjB49GkuXLkVKSord46+++iqaNWuGxYsXY926dTAYDGjQoAFefPFFPPnkk6WuelTZyhvf0KFD0aRJEyxcuBD79+/Htm3bEBgYiPj4eDzyyCN27zFkyBB8+OGHkMvlNp9S3u6jjz7CjBkzsGvXLvz6668QRRERERGWP6gPPPAAEhIS8N1332Hfvn3YtWsXlEolwsLC0K9fP/Tv399yrZdffhm7du3CiRMn8Ndff8HLywv16tXDK6+8guHDh1utnkJEt3jCs7e4V8HRM6VYVFQUunbtir///ht//PEH+vbtC6Dsz6mEhAT8+OOP+OKLL/DXX3/BYDAgISEBc+bMgb+/v0uSBaD0Z1t5OfN3pbKe3Xxu1z7cwZmIXG7fvn0YPXo07r//foe9MERE5Fn47CZ7OGeBiFzum2++AQCMGjWqiiMhIqKy4rOb7OEwJCJyiTNnzuDPP//EiRMnsGPHDvTs2ROtW7eu6rCIiKgEfHZTaZgsEJFLnDhxAp999hn8/Pxw77334p133qnqkIiIqBR8dlNpOGeBiIiIiIjsYs/CTSaTCUaj83mTVCpU6Pzqgu2sWdjOmqW6tFMul5a5Lp/NZVMb2lkb2giwnTVNdWqno2czk4WbjEYROTlqp88PDFRW6Pzqgu2sWdjOmqW6tDM01L/MdflsLpva0M7a0EaA7axpqlM7HT2buRoSERERERHZxWSBiIiIiIjsYrJARERERER2MVkgIiIiIiK7mCwQEREREZFdTBaIiIiIiMguLp1KRERERJWmqKgQBQW5MBr1lrLr1wXUhn2Bq7qdEokUXl4+8PVVQSaTO3UNJgtEREREVCn0eh3y87MRGBgCudwLgiAAAKRSCYxGUxVHV/mqsp2iKMJoNEKjKcSNG9cRHBzuVMLAYUhEREREVCny83Pg5xcAhcLbkiiQewiCAJlMBj+/ACiV/igszHPqOkwWiIiIiKhSGAw6eHn5VHUYtZ63ty+02iKnzmWyQERERESVwmQyQiKRVnUYtZ5UKoXJZHTqXCYLRERERFRpOPyo6lXkNeAE5xrkzJmTNmXx8c2qIBIiIiIiqgmYLNQgV65csiljskBEREREzmKyQERERERudSklGYVqdZXd39vLG9GRURW6xvHjR7F8+VIcPXoYubk5UCp9ERcXj379BuDrr79Eevr1Uq+xa9dBbNy4AdOnv2spk8lkCAsLR+/e9+CJJ56Gl5eX5dihQwcxceLz+PzzuejYsVOF4i8rJgtERERE5FZqrQZ+St8qu3+BurBC569YsRSzZ3+Odu06YuzYCYiIqIv8/Dzs378Pn3wyAyNHPo7ExC6QSgUYjSI+++wDGI0mvPrqFIfXfP/9DxAaGo6iokL89defWLToO6jVhZg8+bUKxVpRTBaIiIiIiMro8OFDmD37czz00FC8+OKrVse6d/8PHn10JIqKihAb28SyKZtS6Quj0YgWLVo6vG6TJvGIiqoPAOjYsTOSk6/g11/XY9KkVyCRVN2aRFwNiYiIiIiojJYs+QEqlQpjx060ezwyMgqxsU0qfJ+4uARoNBrk5uZU+FoVwWSBiIiIiKgMjEYjDh06iI4dO1vNJagMaWnX4OfnB5UqoFLvUxqPGoa0adMm/Prrrzh+/DiysrJQt25d3HPPPXjuuefg5+dX4rnx8fF2y9etW4emTZtWRrhEREREVIvk5uZAq9UiIqKuy69tMhlhMBigVquxY8cf+Ouv7Zg48SVIpVW7qZ1HJQsLFy5E3bp1MXnyZERERODkyZOYM2cO9u3bh59++qnU8VoPPvgghg0bZlUWHR1diRETEREREVXciBEPW30/ZMgjeOihYQ5qu49HJQvz589HcHCw5fvExEQEBgbi9ddfx759+9ClS5cSzw8LC0ObNm0qOUoiIiIiqo1UqgB4eXkhLe2ay689ffonCAsLQ3Z2NpYvX4K1a1eiWbPm6N9/oMvvVR4eNWfh9kShWMuW5lnj16+XvlYtEREREVFlkclkaNu2PQ4e3AedTufSazdq1BgJCc3Qpctd+OijL1C/fgPMmzcLRUVFLr1PeXlUsmDP/v37AQCNGzcute5PP/2EFi1aoHXr1hg9ejQOHjxY2eERERERUS0yatQTyM3Nxbx5s+weT01NwblzSRW6h0KhwAsvTEJ29g2sXbuyQteqKI8ahnSn69evY9asWejataulh8GR+++/Hz179kRYWBhSUlLw7bff4vHHH8fChQvRqVPpO9xJpQICA5VOxyqVSip0visoFLYvp6tj8oR2ugPbWbOwndVXTXg2u0NtaGdtaCNQ89p5/boAqdT+Z9MSB+XuIEglDuMqTfv2HTBx4kuYNeszXL58EffdNwjh4RHIz8/HwYP7sX79Wrz77nTL4jtSqQSCIEAQ7P8sJBLB8t/bj/fo0RNNmzbHTz8twSOPPApvb2/L8WPHDkN9x8ZyUqkUPXr0dNxmwbnnqccmC4WFhRg7diykUilmzJhRav2PP/7Y8nWHDh3Qu3dvDBo0CF988QWWLVtW6vlGo4icHOe3HQ8MVFbofFfQ6Qw2Za6OyRPa6Q5sZ83CdnqW0FD/MtetCc9md6gN7awNbQRqXjtFUYTRaLIpV3p5Iy8/vwoiMvP28rYbV1k98shwJCQ0w/Ll5p2cc3JyoFT6IiGhKV59dQq6dOkGo9Fk2ZRNFEWHPwuTSbT8987jzzwzFi+9NB5r1qzEsGEjLccXLvza5jo+Pj7YsmWnw5hFseTnqaNns0cmCxqNBs8//zySk5OxaNEiRERElPsafn5+6NGjB1atWlUJERIRERGRs6Ijoyr0Zt0TtGzZGi1bti5T3TlzvnJ4bMCAQRgwYJDdY4mJnbFr161h9e3adbD63h08LlnQ6/WYOHEijh8/ju+++87h/gllJQiCiyIjIiIiIqpdPCpZMJlMeOWVV7B3714sWLCgQsugFhQU4M8//0SrVq1cFyARERERUS3iUcnCu+++i02bNuH555+Hj48PDh8+bDkWERGBiIgIpKSkoG/fvhg3bhzGjx8PAPj2229x8eJFdOrUCWFhYUhNTcXChQuRmZmJTz75pIpaQ0RERERUvXlUsrBzp3lSxvz58zF//nyrY+PHj8eECRNuTg4xQhRFy7GYmBhs2bIFW7ZsQUFBAfz8/NC2bVtMmzaNPQtERERERE7yqGRh+/btpdaJiorCmTNnrMp69eqFXr16VVZYRERERES1ksdvykalO/N3GlZOPYCiZB+r8qJkH6ycegBn/k6rosiIiIiIqDrzqJ4FKr8zf6fhyKarAADNBV8AgE9UEYqSfSzfFx+Pv6v8S9ASERERUe3FZKEauz1RMBOgueALQ7YchmwFgFvLxjJhICIiououPf26TVlYWHgVRFJ7MFmoxo5uvmqnVLBJFG6vz2SBiIiIqquMDCYL7sY5C9VYq371HRyxvxFdq3sd1SciIiLyXOkX8vDLp0dQeN161+fC6yb88ukRpF/Iq6LIaj72LFRjxb0E1kOR7Gvdvz7iu97qVfDkbrwzZ07alMXHN6uCSIiIiKiqpV/Iw87FSTDqTbi8A2h4txy+4RIUXjfh8g49RCOwc3ESuo9qgrBGqkqPZ+PGDZg+/V3L9z4+StSrF4n77x+MBx54CDKZ9dvrlJRkjB49DFqtFj/9tBZRUdXrw1smC9Vc/F0RuH4+D2lJuQ7rRMQFWCUKgGd34125csmmjMkCERFR7XN7ogAAohG4vEOPsBZSpB83QjSa6xn1JrcmDADw/vsfIDQ0HGp1If74Yys+//xjZGdn4+mnn7eq9+mnH8LPzw9ardYtcbkahyFVc2f+TisxUQCAtLO5OLOby6cSERFR9bJ/7UVLolBMNMIqUShm1Juwf+1Ft8XWpEk8WrRoicTEznj99bfQtm17rFy5zKrO5s2/ISnpDEaOfNxtcbkaexaqMdvVkBw78tvN1ZC6coIzERERVQ1NgR7711xE+sU8mAyi09e5M1Eops7RYcXbB0o9XyITEBajQuKDMfD2kzsdx+2aNm2Gf//9B9nZNxAUFIy8vDzMmvUZXnhhEkwmU+kX8FDsWajG7K+GVEL9MiYWRERERJVh/5qLSEvKrVCi4Aomg4i0pFzsX+O6nojU1FRIpVL4+CgBAF9+OQsNG0bj3nvvc9k9qgKThWrM8WpIrqlPRERE5EpZVwuqOgQrFYnHZDLCYDAgLy8P69atxo4df6Br127w9vbGkSP/YtOmX/Hqq2+4MNqqwWFI1Vj8XRHIz9LgwoGMO46IuHP51LBG/txjgYiIiKpUnfp+pc61dKc69f2cPnfEiIctX0skEvTtey8mTXoZer0eH388HUOHjkBMTCMYjdV3CBLAZKHau3PSDwBIfIwwFVm/tDlpRTDoTZDJ2ZlEREREVSPxwRiXzFmoqNvnLDhr+vRPEBYWBqVSifDwuvDy8gIALFnyA/Lz8/DII48iPz8fRqMJGo0GAKBWF0KtLoRS6euSdrgDk4VqTK81IuVktk25d1QR1Of8APFW74JObcCVI1lo1CHUnSESERERWXj7yXH36DhIpZIyfeL+y6dHoM7R2ZQLUvuTnJWBCgx8ubUrQi1Vo0aN7e6ZcOnSRWRlZWHw4P42x558chRiY+Pw/fdL3RGiSzBZqMZSTmbDoLP+H02QmaAI18CQK4cu3dvqWNKe64hpHwJBsL/DMxEREZEnSRwSY7XPAmBOFO7cZwEApHIJEoc431PgKqNGPYH+/QcCgCUp2rdvD5Ys+QFTp76PBg0aVnGE5cNkoRq7dDjTpkweqoUgAbwii2yShdz0IqRfzEe4mzYrISIiIqqIsEYqdB/VxJIwCNJbOzj7BEssOzhL5RK3bshWkoYNo9GwYTSAW8lCWto1AECzZi2q3Q7OHMBeTRXmaJF+Md+m3CvcPCZO5m+ATKW3OZ6023bnZiIiIiJPVZwwKAMVlkQBAHzDJWh4txzKQIXHJAo1EXsWqqnLh7PMix7dRuJjgNTfYPneK1INQ16AVZ3UszkouKFxR4hERERELhHWSIWBL7fGiRNHrcp9wyUY+HIrt8YyYMAgDBgwqNLP8RTsWaiGRFG0OwRJEabF7dMR5CE6KAMUVnUCb5xFwWMPQn7mtFW5/Mxp3Hj4fugOHayUmMtCvXwJMnt0Ruie3VbloXt2I7NHZ6iXL6miyIiIiMgThIaG2/yjysVkoRrKulqIgiytTbki3LrHQBCA2E5hlu8Ds8+i9bEvIcvLROCcWZaEQX7mNILmzoLpehryXnupShIG9fIlUM+ZCZhMiNzyO0L37gFgThQit/wOmExQz5nJhIGIiKgWCwsLt/lHlYvJQjVkr1chLMYfUm/bJchi2odCKpdYEgWpybz8mESvQ9DcWVBu/g1Bc2dB0N1clkyrcXvCYEkUbhIARP6+GY2WLkbklt+ttpdjwkBERETkPpyzUM0Y9SZcPXbDpjy6bQiSMi/YlHspZWjYpg7CZy2yJArFBJ0Ofht+hmAwWJ+k1SD/nTfh99oUwGgETEaIRiNgNCLPSwpNfpHlexiN1l+bbn1tc+y2r2+/pm7HnzZxCwBU587B3iKv6nmzoRw2shw/NSIiIiJyBpOFaib1TA70GutdSGQKCSKbBSFph/1zmnQOx77fH7PqWShmkyjcJOZkI3/KqzblBc6F7RRHu0EoX5jkxiiIiIiIai8OQ6pmLv1rOwQpslkQ5F5Sh+cEhPnA965OyBz6BnBzK/LqyuuhoVAOHV7VYRARERHVCuxZqEY0BXqkncu1KY9uE1LquV0fbQxBiEWB6Sw0K5ZVRnhuoV29AmL6dSifHw9ZNdsBkYiIiKi6YbJQjVw+kgXxjjnMygAFwmL8Sz1XEATo/jkAzcrllRSd++h2/gXd7l3wvn8IlGOehiQouKpDIiIiIjdIT7fdXJYrIlUuJgvViL0hSA1b14EgcTS6/xbdoYPIe/VF2GQb9kgkkMY3g7ROMCCVAlIpBIkEkEqh8PGCzgjL95BKIUhlN7+W3PxeCtwsK0s97e5d0G3cUL4fhtEIzdpV0G7aCJ+Rj8Fn2EgI3t7luwYRERFVKxkZTBbcjclCNZF9TY3c60U25Q3b1CnT+QXT3gX0eptyEXYmEptMEG9kQvXVQpv6gYFK5OSoy3TPssqf+obdcrux3VmnSA31NwugWbcGyqeehVf/geYkhIiIiGoM3aGDKJj+HuQjRkEfn2Apl585jRvvvAW/KVOhaNehCiOsuTjBuZqw16tQp74vVKE+ZTpf1qad3fKiLndBVFjv8gwvb/hNmVruGJ2lHDfBpkwEkBfbBGIZr2HKzEDBh9OQ8+Qo6PbuhiiW9UwiIiLyZLpDB5H32kswXU9D0Nyq31R248YN6Natg+Vf37534/HHh2P16uUw3Fxlctq0/0O3bh3QpUs7q7ozZ37qlhhdiT0L1YDJaMKVo1k25Q3LMLEZAIypKdD9uc2mXBfTCHmPPYEbsZ1Q76c5kOh1gJc3VB995tbsvHjPhOKN2UQAKff0Q0bnLpYdnIt7GOSdusJw7AhEdaHdaxkvnEfeqy9C3q4DfMdNhOy2Tx+IiIioeilOFKDVADDvERU0dxYK7hsEv1832Gwq6873MO+//wFCQ8OhVhfijz+24vPPP0Z2djaefvp5AEBgYBA+/vhzGI23PsAMCSnbezdPwmShGkg7lwdtofV+CBKZgAYtyzaxt3DuTEBnvb+CKJHgco/HkfarDvrCxsjrMhYtziyvsm48S8IwbzZS+vRFRucuAICMLl0BAFHbtkI5bgKUw0bClH0D6u+/hebnNeYN3uzQHzqInKdHw+/1N+E98AH3NIKIiIhcqmD6e5ZEoZig01knCsW0GhRMfw/Bq9a7JbYmTeIRFVUfAJCY2BnJyVexcuUyS7Igl8vRokUrGI1lmC/qwZgsVAP2hiDViw+Ewqf0l093YJ/dHZKv1u+Nq5duZbdpkji0mLkMikjfCsVaEcphI6EcNhL/btloVZ7RpSvaTP2f5XtJUDD8Jr8Kn4eGovCredD99Yf9C3p7Q965a2WGTEREROVgyr6B/OnvQX/ooM0HmeVhkygUX/96GjK7J5Z+AYUC8nYd4D9lqstWVWzatBn+/fcfZGffcMn1PAXnLHg4rdqA1NM5NuXRbW27sRo0iLb+Vy8KhTM/s6kn1AlBaov+NuVJe2xXGPBk0gYNofrfhwiY+zVkLVraHPcZNhLSkNAqiIyIiIjsyZ/+HvR7d1coUXAJnQ76vbuRP/09l10yNTUVUqkUPj5KAEB29g3ce28v9OjRCY8++iAWL/4eRgcjIjwZexZuyi8sxOylP6JD85Y4ffECCtSF6NGxE/46sA9REXXhrVDg3JXL6NSqDY6cOQWDwYjOrdtg16GDiI6Mwv29emD20iXo1q4D9h45DJlMitbxTbHv6GHENmgIjU6H5LRrlmv6KX2RENMIB08cQ0JMY+Tk5yEtM8NyPNBfhejIKJzZcw0qY12rWI1SA/46txMReSE4nnQWbRKa4VJKMnLy8yznR4SEosnhf2G8fNGmrbvbtEF+nVzICsOsyi8dzYSksRp7T/2D6Mgoc1lKslWbxo14FLOX/lChNh0+fRItmsQhLTMTmdk3LMdDgoIRERKCM5cuIFgVgIIiNXR6PcLrhGD20h8QERKKQH8VTl88b/s6deuO1s1boNHW3+GdlQWDvz9W+PpAs3yJ1etkr032XqdXnhqD2Ut/cFmbSnqdHLapjL97ZW2Tvddp0uOPYeXvm2pUm+y9TuNGPordh47WqDbZe50e6tcHp89f9vg23d+nR5mfzXkF+fjsh+/QJr4pkq5cQmFREe5q0x5/H/4H9ULD4KVQ4GJKMto3a4Hj587CYDSiQ7OW2HfsMOpH1MXD/e/B5z/8gE4t2+DgyWOQSaVoERuHf04eR0xkFLQ6HVIz0i3X9PXxQZMG0Th85hSaNIhGbkE+0m9kWY6r/PzQIKIejp87i4SYxki/kYUbuTmW48EBgQgLroPTF8+jRWwcrqSlIq+gwHI8LLgOAvz8kXTlktNtAoCrades2jTx8cfw2Q/f1qg23fk6DenXB4ePn6lRbbL3Ot3VsS1+37G7xrTpnk4dkZaZAW+FF4wmI0zHjpS60qE76Y8dRVpmBvx8lCjUFEEA4OPtjcKiIngpFBBFETq9Hn5KXxSoCyGVSFCkMQ+LysnLhSLLB3n5efhn/z7s2PEHOnTsBK1Bj/C69fDU08+jadOmSM/IxKF/DmDBgrk4d/4cJr/0Goq0Wvh4e0On18NoNFquL5fJIJVIodFpofT2gVanhdFkshxXyOUQBAFanQ6+Pj4o0mggAvD19kFBkRoKuRwAzDHf0aa8ggLs2brJ4e/esEH32v0ZCSKXjQEAFBZpcCU53enzo+uH49JV138yf3DpVeRd01qV1W8fgCb/KeUT85wc4LkxgPqOZU6bNgM+/AxJZ88i5XcJTDrr/2WjuwShUVfHy7FWVjtvt3f3XzZlnbuW8c2FwQBs+hVQ+gK9+pRc98/tQGgo0Ny2V8Id7fQEbGfNUl3a2bRJgzLXzS9Q49zFZKfvFd+4Ac6cv+L0+dVFbWhnbWgjUAPbaSxCWFiU5Vv1f1+BYf/eKgzImiyxM5QffFKuczZt+hUffTjNqkwikaB373vwwvgXoVKpLOXeXl7QaM3v4+bOnYnVq5bjx0XLLXMd3Ck9PRmQOl5Fs23LOLvl7FnwYIU3dDaJAgBENFPZqX2HH7+zTRQEAXh2HCAIkEgBv4Yi8pKsk4WUI3mITgyGROZJeX85yGRAWSY05+YCX842/4w6dQEefxKoX/Y3MERERFR+3q+/Cc2H02A4fKhqhyLJFZC1bQfv1990+hLvvT8DoaFh8PFRIiIiAgqFV4n1e/Xqi9WrluPMmVNVkiw4i8mCB0s7kW9T5heqgH9Yyb+MOHsG2LrZtvye/kBsE8u3/jEi8s6JgHgrMdCrjbh+Jh91m5chIanOViy9lUzt2wMc2Af06w8MfwwICqra2IiIiGooSVAQlB98YvWJe0nyhz8E8bptj6lJoYDETrIhhIfDf9lql8RampiYxoiMjCq94h0EoXp9IMsJzh5KNIlIO2mbLESU9ibeZAK+mgfcObrM1w947HGrIpkPoKxnOwot+VBuzd7U7FoqsPEX6zKTCfjtV/PQrZ+WwFRku1s2ERERuZfPa28CXtYfkpoUCuQPGAiTzaayXub6Hmrb1s0QBAHx8U2rOpRyYbLgobKvFkFbYL23giAAEQl+JZ/4xzbg5s6GVkaOBgICbYpVjWyTgvx0LXJTNDblNcamjea5DfYUFQFLfkTK/fcDm39zuI8DERERVT5Z23ZQTv/YkjCYFArcGPsCCvvegxtjX7iVMHh5QTn9Y8jatqvCaM3S0q5h0qSx+Hndauzbtwe7d+/CRx9Ow9q1qzBw0GCneiOqEocheSh7vQrBMUoofEt4ydSFwPff2pY3aAgMGGj3FK9gQBEoQpdj3SV29d8cBEY5ngRTrY0eY/6ZLP4eyLTdwwIAjBkZwJwvgPVrgTFPA+07mrM1IiIicqvihKHoo2m48egI6OLiAQC6uHjcGPsCQn9aCp/X3vSIRAEAlEpfqPxVWPbTYuRkZ0MQBNRv0BDjJ0zGAw88WNXhlRuTBQ9k0JmQfrbAprxuaRObf1oK5GTblj87DpBKHZ6maiwi8x/rN8KZSYXQ5OnhrZKXKeZqRSoFevcFut1tTgZWLgeK1PbrXrkMvPs20KqNOWm4bc4HERERuYesbTv4L1uNlFPHrcp1cfFum6NQ7N5778O9997n8LhKpcL7//sQAMo8N8OTcRiSB8o4WwCTwXp4kMxLgjqNlY5PunrF/Mb3Tnd1B1q3KfF+ynoipF7W9xNFIPlwbllDrp68vIBHHgW+/s68glIJCRWOHgYmjwc+/RBI9/ylKYmIiGqikJAwm39UuZgseKBrdoYghSf4QSpz8HKJIvDNfNvx9QoF8OQzpd5PkAB+MbZzF1KP5cGoN5Up5motIBB4bhww72tzclWSP7cDzz0FLPwaKLB9nYiIiKjyhISG2fyjysVhSB6mKFePnKu2K/GUuLfCgX3AoX9syx8eBoSFl+m+/tEi8pMEmIy3kgaDxoS0k/mIbB1Qpmu4SkTdSLfez6JeJPDft4DTJ+G1+HtojxyxX8+gB9auArIygVffcG+MRERERG7EZMHD2JvYrAySQ1XXwd4KOh3w9Xzb8rBw4MFHHN4nKNh2l2ZjghHX7tjbIfnfXNRrpXLrmsDRMbFuu5ddCc0Q8f33uLxiDfDDQiA1xbaORAIMG+H+2IiIiIjcqMzJwpw5cyp0o8GDByMqqnotFeVuouhobwV/x2/Wf14DpF2zLX/qWZt1iW8XHBxiUyZvp7VJFgqzdMi+XITg6BLmS9RAgiAAXbsBiZ2BzRuBZYvNuz4X69vPvKISERERlUgUxWq3EVlNU5H9s8qVLAiC4NTNBEFA+/btmSyUIjdVg6IcvU15RFN/+ydkZgDLl9qWt24DdLmr3Pf3D/NCYJQ3cpJv7bHgG6KAiBq8QVtpZDLgvvuBnr2B1SuBdWsAAcCIx6o6MiIiompAgMGgh1yuKL0qVRqDQQdnpyqXaxjSG2+8gd69e5frBrm5uXjwweq3pmxVsNerENTAx/Hypd99A9y5HJdEAjwz1uk9Aeq3C0ROchpCGvsiql0Agur78NMAAFD6Ao89AfS/DziXBNgZxmXlj23AyePmpCIo2C0hEhEReRpRkCE3JxMBgSGQyeR8T+FGoijCZDJBp9WgoDAXoiCFMz/9ciULQUFBiIws3+RTpbJ2DV9xllFvQvoZ270VIpo56FU4cQzY8adt+cD7gYbRTscR0tgXnZ9qCGVgDdxfwRVCQs3/SqLVAj8uNG/49ud289yRwQ8BPjV0kzsiIiIHBIkMRhOQnZ0B3DZSQS6TQW8wVF1gblL17RQACBAlcghCJfcsLF++HA0aNCj3DVQqFZYvX47GjRuX+9zaJPN8IQxa62VKpXIBoU38bCsbjcCCebblqgBgeMWGxwgSgYlCRa1fd2tnaI0GWLoI+O1XYORjQJ9+Je/n4AFu3LDd1dreHBciIqKyECQy3PmWs1GjBjhz/krVBORGntLOivTnlDlZaN26tVM3kEqlZT5306ZN+PXXX3H8+HFkZWWhbt26uOeee/Dcc8/Bz8/Om+bbaLVafPHFF9iwYQPy8vLQtGlTvPLKK+jYsaNTcbubvSFIoXF+kCnsZIG//wZcvGBbPnoMUMrPiSpZbg6w6ifb8uwbwJyZwM9rgSeeAjp2cnqoWGXLvpFlU8ZkgYiIqHbyqE3ZFi5cCIlEgsmTJ+Obb77B8OHDsWzZMjz55JMwmUreHGzKlClYuXIlJk6ciAULFiA0NBRPPfUUTp065abonactMCDrktqmvK69IUj5ecCi723LY5sAfe5xfXBUPvaSuNtdvQK8/w7w5mtA0ln3xERERETkJI/aZ2H+/PkIDr41GTQxMRGBgYF4/fXXsW/fPnTp0sXueadPn8Yvv/yC6dOn46GHHgIAdOzYEffddx9mzpyJ+fPt7EPgQa6fysedCw55q2QIrG9njPuSH4F8OzsHPzvW44e31Apt2gFffQ+sWAps/AVwNE7x2FHgpQnA3f8x9wiFR7gzSiIiIqIyqVDPgk6nK/fqSCW5PVEo1rJlSwDA9evXHZ63bds2yOVyDBgwwFImk8lw3333YdeuXdDpdC6L0dVEUcQ1e3srNLOzt8LFC+ax73fq2Rto2rySIrSWm2G7uzTdISDAvCLVvK+Bu7qXXHfHn8DzTwPffgUU2EkCiYiIiKpQhYchpaTY2d3Whfbv3w8AJU6QPnfuHCIjI+Fzx2ozsbGx0Ov1uHz5cqXGWBEF6VoUZtomMzarIIki8NWXwJ3DsXx8zGPgK5HJJCL9bAH++SkZq6cdQkGGtvSTCKhbD/jvW8DHX5SczBn0wLrVwDNPAGtXAXrPTW6JiIiodil1GFJJPQeVvSPf9evXMWvWLHTt2tXSw2BPbm4uAgICbMoDAwMtx0ujkMsQXT/c6Vi9FHKnzt97wHaMe1iMP5q1qm9VVrh5MzKOH7WpG/Tccwho3azc9y2rpP3pOPTbVRRm30oQcs5q0aJd+VfGqk6cfT3tqh8OsU8PqP/4A9kzZ8LgKHktKAAWfg3Zpl8ROGECfPv1gyCp3GlF9tqZk2W7I7jLfhZVxKWvpwerie309pIjvrHzzxsvL0WFzq8uakM7a0MbAbazpqkJ7Sw1Wbhx4wZeeukl1K1b1+aYXq/HSy+9VCmBFRYWYuzYsZBKpZgxY0al3ON2Or0BV5LTnT4/un44Ll11PFTKHpNRRNIB23sGNfa2vpZGA3z8se0FIqOQfXcfZJfzvuVx/XquVaIAAOcOZiCsrS8Uypo7R8KZ17NUTZoDX8wzr2a1bDGQk2O3miE1FZnffofMZm3Mm+xVInvtzMrOs6nn8p+Fm1XK6+mBqks7mzYp+x9OjVaPcxeTnb5XfGPPWLawstWGdtaGNgJsZ01TndrZtmWc3fJSk4WmTZuiTp066NOnj80xnU4HURTtnFUxGo0Gzz//PJKTk7Fo0SJERJQ8+VOlUtkdDpVz882YvV4HT5B1UQ19kdGqTCIVEBZ/x/Knq5bfWrf/dk8/D8grd0+Eus39cWFXltUeECaDiNRjuYjuxJ2Jy00mAwYMAv7TG1izEli7GtDZGdb15DOVnigQERERlabUdyOPPfYYgoKC7B6TyWQu/9Rfr9dj4sSJOH78OL766ivEx8eXek5sbCxSUlJQVGQ9+fb8+fOQy+Vo2LChS2N0lbSTtp/ghjT2hdz7tk/s066Z31TeqWMnoEPl7yEhlUtQr6XKpjzlcC5MRtcnirWGUgmMehz4aiHQ917rxCCxM9CyVdXFRkRERHRTqclC//79HS5ZKpFIMGTIEJcFYzKZ8Morr2Dv3r2YN28e2rRpU6bzevXqBb1ej02bNlnKDAYDNm7ciG7dukGhULgsRlfRFxmReb7Qpjyi+R0Tm7/9CtDrrctkcuDp5yoxOmuRbQNstv7TFhiRkVTgthhqrDohwMTJwMx5QPuO5qShtAnromgemkZERERUyco9zuHAgQMwGo2lV3TCu+++i02bNmHMmDHw8fHB4cOHLf/S0tIAmFdfatasGebMmWM5r1mzZhgwYACmT5+OlStXYs+ePXjppZeQnJyMCRMmVEqsFXX9dD7EOxY2UvhKERytvFVw6B9g727bkwc/CNSLrNwAb+OjkiM01tem/OqhHLfFUONFxwD/9z9g7ldA/VLGc//7j3nlpN9+BSrp/0UiIiIiwIlN2UaPHo09e/ZYVhpypZ07dwIwb85250Zq48ePx4QJEyCKIoxGo81ciRkzZuDzzz/HF198gby8PCQkJOCbb75B8+bu2X+gvOztrRCe4A+J5OZH+AYD8PWXticG1wGGDq/k6GxFtQtERpJ1T0jeNS3yrmmgquvt9nhqrKj6JR83GoHvvgFysoF5s4AN68w9ER07AZW4MhkRERHVTuVOFipjQnOx7du3l1onKioKZ86csSn39vbGG2+8gTfeeKMyQnOpwiwd8tNsJ7XWvX0I0i/rgeSrtiePedq8t4KbBUZ6wy9UgYIM6z0Arh7KQfP7uPuw2/yxDbh08db3V68A778DtGhl/t2IK32ODxEREVFZcbmVKnDthO3EZr9QBfxCvczfZGcDyxbZnti0GdCjZyVHZ58gCKjfLtCmPP1sAbQFBvcHVBtpNMDiH+wfO34UeHki8PEM4OaQPSIiIqKKYrLgZqJJxPVTtkOQ6ja/bcWhRd8BarV1BUEAnh1XpUNNwhL84O1nvVSraAJSjpS+6R25gFYDNCtlWN2OP4GxTwPfLgDybZNSIiIiovJgsuBm2VeKoC2wnpQqSIDwpjf3Vjh7Gtiy2fbEfv2B2CZuiNAxqUyChK62u8OmHMmF0WCycwa5VEAg8NoU4JOZQPMWjusZ9MC6NcAzY8zL7up0jusSERERlaDcyYLASZQVcs3O3gp1opVQKGWAyQQsmGd7kq8fMOqJyg+uDBLuqgvhjt8afZEJ109xGVW3iU8AZnwCvPV/JU+ILiwwT4Ye+zTw53bz7xcRERFROZQ7WajMCc41nUFrsllRCAAiiocgbd8KnLWdvI1RowEP2YVaGaBAWLy/TXnyvzn83XAnQQA6dQHmLADGTQBKWp0s/Trw6YfmOQ1HDzusln1Fjd1fX4Imw7pckwHs/voSsq+o7Z9IRERENVa5k4UBAwZ45CZn1UH62QKYDNZvqGXeEoQ08gUKC4EfFtqe1DAa6D/QPQGWUf12tolLQYYOOclFdmpTpZJKzb8fC74DHh0JeHk5rnsuCXjzdeDdt82T6G+TfUWNI2uvQZNnQPo+CTSZ5nJNBszf5xlwZO01JgxERES1TLmXTv3ss88qI44qp9XqsGzjenRo3hKnL15AgboQPTp2wl8H9iEqoi68FQqcu3IZnVq1wZEzp2AwGNG5dRvsOnQQ0ZFRCK3TA8s2rke3dh2w98hhyGRStI5vin1HDyO2QUNodDoU7lNAAevNzfJ8MnEkqQh11qxEg5xsm7i2tW6NsAtJSMvMRGb2DUtMIUHBiAgJwfGks2iT0AyXUpKRk59nOR4REopAfxVOXzzvdJsA4FJKslWbxo14FL8e+h11A+NhypFaxXrg93OI6eOPgyeOISGmMXLy85CWmWG5Z6C/CtGRUTh8+iRaNInzmDbd+Tolp13DK0+NwbKN6+Gn9EVCTKPq0abefXDY2xstD+5DwxMnIDjo6SlMOoNrWRk4sGcnJj3+GLb/fhDicRUE0fzZgWgUcH2PgILAa/DLqQuYzEMPTQYRh1Ylw9QsBxGNVR7xOhVfs7TXadzIR3HqwjnPeJ1c1CZ7v3sPqfogNeO6x7epaZNSNh68jU6vw4rNv6JNfFMkXbmEwqIi3NWmPf4+/A/qhYbBS6HAxZRktG/WAsfPnYXBaESHZi2x79hh1I+oi8h6IVix+Vd0atkGB08eg0wqRYvYOPxz8jhiIqOg1emQmpFuuaavjw+aNIjG4TOn0KRBNHIL8pF+I8tyXOXnhwYR9XD83FkkxDRG+o0s3MjNsRwPDghEWHAdnL54Hi1i43AlLRV5BQWW42HBdRDg54+kK5ecbhMAXE27ZtWmiQ0ew4rNv9aoNt35OoWGBiDlelqNapO91yk0NAD/njpRo9pk73UKDQ3Apr//qlFtsvc61Y0IxorNv1aLNrVtGWf3OSyIHDsCANDrjcjJcf5T08BAZYnnF2RrsfGzozblvZ9rigBjJnIeH26zG6+iZ2+o3pvhdEyVobidV4/dwJ4V562OefvL0X9SS8i9pA7Orj5Kez09neHieRR+OQf6PX/bHPP/v2nw6t0XgLmdi9/eA3WO7SRoQQqIdjaIVgYqMPDl1i6PuTJV99ezrKpLO0NDbYcyOlLZz+aaoja0sza0EWA7a5rq1E5Hz+Zy9yzc7uLFizh37hyysrIgCAKCg4PRpEkTREdHV+SyNdLlfzNtyvxDvBFUT4n8Vz61SRTg5QXfcRPdFF35RTYLhI9KjqI8PQIjlIjrGo76LYIhlXOBLU8gi2mMgI8+h+7QQRTOnQXj2dPm8qbNoW/XHdnncpGXoYE2Tw+5t/3kzl6iIJVLkPhgTGWGTkRERB6k3MnC+fPnsWzZMmzevBmZmeY3wMWdE8UrJdWpUwf9+/fHo48+isaNG7sw3OpJFEVcOpxlUx7dpg70u3dBf2CfzTHlyMchvdmV5IkkUgna3dcQCqUMIQ39uEqWhzEaTCjI0iBP0Qj5I6fDtHs7gnYtwxHffsj67FiJ5/rnX0G+X32bPT2kcgm6P9YEYTEqB2cSERFRTVPmZOHKlSv45JNPsGXLFnh7e6N9+/YYNmwYGjRogMDAQIiiiNzcXFy5cgWHDx/GqlWrsHjxYvTt2xevvvoq6tcvYYnHGi7zSgEKs7XWhQLQoJk/Csd/blNfElEXPiNGuSk650U2Cyq1Tnr6dZuysDDbvRrItf5eeg5pSbdvltcUQvt3IEpKHiLml38VHf/5EDkBjZHU+EHkq6IBmBOF5r3qMVEgIiKqZcqcLAwYMABxcXGYMWMG7rnnHiiVyhLrq9VqbN68GT/++CMGDBiAY8dK/jSzJrtkZwhSeCMVhN9WwZSaYnPMd/yLELy83RFapcvIYLJQUaJJhDpPh/wMDfIyipCXqUH95sEIb+z4jbsq1PuOZAGlJgoQRcSeXwsACMw9j46HPsb10PY43+h+aBCCE9tTERzpy4SBiIioFilzsjBz5kz07t27zBdWKpUYMmQIhgwZgq1btzoVXE1g0JuQfNx2laOYhiLU07+zKZd3SITi7v+4ITLyNEa9CflZGuRnapCXoUH+zcQgP1MDo956QzVvX3mJyYJ/SPmTzeAbJxGcY73PR3jGPwjNPIzkyLtxqeG92LkoiUORiIiIapEyJwvlSRTu1KdPH6fP9XTFw2yKiryRn68BYP3JeeqpbOi11jNFZQoJAv5YBL1GY30xqRS+E1/i+P8aTqs2ID+z6GZCoEFeZhHyMzQozNairGuT5WWUvKeFKtTH4TGZQgJVqA/8Q7wR1kCFo9uuQluoQ+yFdXbrS0QjGiT/gbrX9uJyw344sFLEfa91KFugREREVK1VaDUkujXMprBQAbXavPzk7cmCvSFIsYFp0K/YYlPu/eAjkMU0qqRI3Y+r8traPPs4ctMrvnldfqamxOP+od7w8ZfDP9TbkhioQr3hH+oDH3+5JSENDFTCN0SBnYvO4nzM/Yi9sBa+atuhYwAgNxYh9sI6NM7dDc3mF+DVtx8ECVe/IiIiqslcmiwYDAZs27YNubm56NmzJ0JDQ115+WqnKE+H6+fzrMoEkxF19yyyqSsEBkE55hl3hVapspILkLTnOrx85VA0rOpoKpdl1aEM87+Y9iFQqhzvcK5QuuZ/ufwsDUwmERKJ/V4ob185Br3WpkzXCmukQvfH4rBzsYD9wc1Q7/puNE7ZCFlBnt36QlY6Cv73DoqWL4XvuAlQdEh0thlERETk4Zx+5/LRRx9h3759WL16NQDzp8hjxozBwYMHIYoiAgMDsWLFCjRoUPadOmuay0eybIaVROfug3D1gk1d3+fGQeJf9o2KPI3JaELyyWwk7bmOrKuFAMwr6DSpK4NUUf2HVZV16FBghE+JyYJ/qDcyLuWX+/4yLwlUIT7mnoIQcw+B+cau+dmGNVKh+6gm2L/2IuS9e+FGQDcot2yG75bNEHS2G7YBgDHpDPImj4e8c1f4jp0AWSMuk0xERFTTOJ0s7Ny5E127drV8v337dhw4cABPP/00mjZtivfffx9fffUV/ve//7kk0OpGFEWbIUgyfQEanv3Zpq4soSm8BgxyV2iVwqgXcfDnSzBoTbeVmZB9wYiQhOox2q141aG8DA2uFtzA9at5lhWItIWGMl0jL0ODyKaOj6tKmXjso5LDP8Tn5pAhb0uCcPvQocoS1kiFgS+3xokTRyHCG4WDHkDR3T3gu2E9lLt3ASaT3fP0e3cjZ/9eePUfCOVTz0IaGlapcRIREZH7OP0uLi0tDQ0b3hpj8scffyAqKgqvvPIKACApKQkbNmyoeITVVHaqGnkZ1uPKG138FZKiApu6vpNervZjv+XeUsS0DUHS3nSr8htnjagTJ4XgYLiMJ9AVGfDnd2fsrjpUXvmZpU88FiQC/IK9LHMIVMXzCUJ8HO6mXFVMAYHIHzUa9Z4dh8L5c6DfvctBRRO0v66HdutmqKZ/DEViZ/cGSkRERJXC6WRBr9dDJrt1+r59+6x6GurXr4+MjIyKRVeNXTps3avgV5CMyGs7bep53TsA8hat3BVWpYrtHI6kfenAbcNy9GogP8UEVX33vwnWqg2W5Udj2oU4/GRe7i1FQVbFEwUANgnincIa+eOhqe0gkVav5FAW0wgBH34G3b//QD13JgxnTtutJ3h7Q9a0uZujIyIiosridLIQERGBf//9F0OHDkVSUhKuXr2KiRMnWo5nZWWVunFbTZCfb54EqtfLodHoAZgnvV45euNWJVFEk6SVEO6YwCD4KOH73Hi3xVrZ/Ot4o25cAK6dsd4MLOussdKSBdEkQp2rQ16mebhQ8bCh/EyN1dChiFgVlAFedq8hCAL8Q72RnaJ2KgYfldyy4lBQvZJ/56tbknAnRdv2kH/1PXTbt6Dwq3kwXbtmdVz5+FPVeu4NERERWXM6Wbjvvvswb9483LhxA0lJSfDz80OPHj0sx0+dOlUrJjcXJwtarQw6nfnNadrZXOjUt96ohmUcQlDuOZtzfZ54CpKQEPcE6iZxXcJtkgV1hoiibBN8giQovG7CL5uOIHFIDMIalX1jL8uGZbftYpyfoUF+GXsE8jI0DpMFAFCF+JSYLAgSAf51vG4uQepjtSSp3Muzhg5VNkEigVefflDc3ROaNSuh/vE7iPl5kERGwXvwQ1UdHhEREbmQ08nCc889h2vXrmHbtm3w8/PDhx9+CJXK/OYvPz8f27dvxxNPPOGqOKuV24cgSYxaxJ5fY1NHWr8BfB551J1huUVYIxWUgQqoc6xX0LlxxojAGODyDj1EI7BzcRK6j2pSYsKQm16Eo5uvIi9Dg8IcrdXwpvLKy9AgIjbA4XFVqHnisdxLCv9Qb4TU84NXgPzmXAJv+AV7VfteAVcTFAr4PDoSXvcNQtGP30PWoiUEubzEczRbNkHetDmkUfXdFCURERFVhNPJgkKhwPTp0+0e8/X1xa5du+DtXfLKLzWRSS/g2tlbn6xHX/4d3tocm3q+k14u9Y1VdZRxMR+afL1Nec5lE3KvmiDe3MzaqDdhx6KzuPuxOIcJg0Rq/bOsiNImHjfqEIrodiHw9jOvOhQYqEROjnPDkmobib8Kvi9MLLWe8VoqCma8D4givAc/ZB6yFBhY+QESERGR08r1UWmPHj3w/vvvY8+ePTAajY4vKpHA398f8hr4Zrg0unQvmIzmj8B9ijLQ4OpWmzqKu7pD0amLu0Nzi/1rL1rab0WEJVEoZjKI2L/2osNr+QZ6QSJ1fhUliVSAKtQbkc2CEFTPt8S6Xr5y+PgrKn150tpM/fV8QK8HDAZoVi1H9vAHoV7yA0RtyZPCiYiIqOqUq2ehd+/e2Lp1K5YsWYKAgADcfffduOeee9CtWzf4+PhUVozViu76rd6U2HNrIBHvWJ9fLofv+BfdG5QbJQ6Jwc7FSWVeWajtAMfzWiRS8xKjpa0wJPeWWpYe9b+5DKkq1Ae+QRVLNsh1DGdOQ7tlk1WZWFAA9fy50KxZBeUzz8Prnv7VfglhIiKimqZcycLUqVMxdepUHD16FFu2bMHWrVuxYcMGeHt7o0uXLujbty969uyJoKCgyorXoxkLpTAWmHtTgm+cRGjWUZs6PsNG1ujx2sU7Ae9cdBZGQ+mTDLyUJf8KqkJ9LMmCMkBhSQaK9yfwD/WBt5+MPQIeTr10kcNjpvTrKJj2LopWLIPvuIlQdEh0Y2RERERUEqfmLLRq1QqtWrXCyy+/jPPnz2Pr1q3YunUr3nzzTUgkErRr1w59+/ZFnz59UK9ePVfH7LG0N3sVBJMBTZJW2RyXhIRC+dgTbo7K/cIaqdC8VySObUmGWEq+kJehQUhDx0ttNvtPPSTcXbdWrjpUk/i//iaKoqOhXrYYKLI/f8SYdBZ5k8dD3qkLfJ8fD1lsEzdHSURERHeqcJ9/48aN8dxzz2HlypX4888/8cYbb0AqleKjjz5C7969MWTIEOzYscMVsXo0UTTPVwCAqJS/4Ft03aaOctwECLVg74n0C3k48UdqqYmCf4g3ZIqSfwUD6yoRHOnLRKGaE5RKKMc8g+Blq+H9wBBA6vj11O/bg5wnRyF/xvswZqQ7rEdERESVz6UDhMPDwzFq1Ch8//332LVrF2bMmIHIyEgkJSW58jYeSXdDBlEnhUKbi5hLG22Oy1q2hlefflUQmXulX8gr85wFda4O3v61bxJ8bSapEwK/V95A4PdLoeh2t+OKogjtxg3IHv6QefO3wgL3BUlEREQWTi+dWprAwEAMHjwYgwcPrqxbeBRNqrlXodHF9ZAZ75iQKwjwe/GVWjGufv/ai3YTBUFquxqSUW/C/jUXMfDl1m6KjjyFLDoGqhmfQP/vIRR+OQuGUyftV9RqUbToe2g2/AzlmKfhff8QCLJKe2xRKdLTbXtMw8LCqyASIiJylwr91d2wYQOWLFmCy5cvIycnx+a4IAg4edLBm4AaRDQI0GTI4Z93CfXS9toc975/MGRx8VUQmfvZWw1JkAJhLaRIP260ShikcgkSh8RUQZTkKeRt2yFg/kLotm819yBcS7VbT8zJRuHnH8N07VqZ9nSgypGRwWSBiKi2cTpZmDdvHmbPno06deqgbdu2CAhwvDtuTZR+IQ/7116EpIEcpiIpYBQRn7TCpp7gr4Ly6bFVEGHVsKyGdDNhEKRAw7vl8A2XwCdYYtnBWSqXlLqDM9UOgkQCrz73QHH3f6BZuwrqHxZCzM+zrajwgvfDw9wfIBERUS3mdLKwdOlSJCYm4ptvvql1m6+lX8jDsVm/oPXxH3Gq2Sjk1W2Mumn7oMq/bFNX+fRztW6X2uKEYf/aiwhrK8I33Dw1xjdcgoZ3y5H+r4DEITFMFMiKoFDAZ9gIeA0YiKJF36No9QpAp7Mc93nkUUjD+Sk2ERGROzk9wbmwsBD9+/evtYlCi3/nwUd7A62PzEfglVNofOFnm7pF3iHIb9mnCqKsemGNVBj4cmtLolDMN1yCgS+3ZqLgwUJDw23+uZPEXwXfcRMRtGQlvPreCwAQAgLhM+pxt8ZBREREFehZaNq0Ka5du+bKWKqF01/9hhb/zoPUZP7EU2rSodWJryARbSf1KrQ5OP3t7wj7gEMnqPrwlDHo0oi68J/6HnyGjYAx/Tokfn4l1i9atxpifj58hj4Kwcu7xLpERERUNk73LLz44ov46aefasUE5ts1Pb3IkigUs5coAIBUNKDpacc71xJR6WTxCfDq3qPEOqb8PKi/+hLqr+Yhe8TD0Pz2C0SjscRziIiIqHRO9ywkJiZi2rRpGDp0KNq0aYPIyEhIJNa5hyAImD59eoWD9CQBU/8Pua9OhqDTllpXVHgh4J133RAVUe1W9ON3lknRpvR0FEx/D0UrlsF33EQoOnaq4uiIiIiqL6eThSNHjuC///0vDAYDDh48iIMHD9rUqYnJgqJdBwR8/DlyX5kMQe84YRDlXgj45Aso2rZ3Y3REtY/xWqp5MvSd5eeSkPfSBMg7doLP668C4Q2qIDoiIqLqzelkYdq0aZDL5Zg3bx46dOgAlar2TFjNCYzDpegBiD63ARLRYHPcKJHjUvQAxAQ2QVgVxEdUm+gPHwJMjncM1x/Yh6uPPAKvfgOgfOZ5SD1kTgYREVF14PSchTNnzuDJJ59Er169alWiULwaUsPzv9pNFABAatIj+vyvODrzV6RftLNePBG5jHf/gQj8YRkU3e52XEkUod30K7KHP4zCBfNgKihwX4BERETVmNPJQp06dWrdsqmA7WpIjkhNOrT8dy5OL9jkpsiIai9Zw2ioZnyCgDkLIGva3HFFnRZFi79H9qMPomj1CogG+wk/ERERmTmdLDz44INYv349DLXsj6291ZAAwCi1TZykJh2anv7RHWEREQB567YIWLAQ/u9Og6RepMN6Ym4OCr/4BNmPDYP2z+0QRdGNURIREVUfTicL7du3h0QiwdChQ7Fq1Srs3bsXBw4csPlX0wRM/T+ICi+rMqNMjow+vWCUWScMosILAVP/z43REZEgCPDq1RdBi5bDd+JLEEoYJmlKvor8t/+L3HHPQH/sqBujJCIiqh6cnuA8ZswYy9dvvfUWBEGwOi6KIgRBwKlTp5yPzgNZVkO6uXyqUSbHhREjoItrgvywCDRathRSg96cKHz8ORTtOlR1yFXK3bv/EhUTFAr4PPIovO69D6aVS5CzZDGgsz980HD8KNTff4OAT2e5OUoiIiLP5nSyMGPGDFfGUa0UJwwF099DUr9+KIiOgQJAQUwMLgwfgfjfN8NvytRanygAnrMbMNVeEn9/BL/0EoT7BkP99ZfQbv7NtpIgwHfsBPcHR0RE5OGcThaGDBniyjiqHUW7DghetR4FWzZalRfExCB41foqioqIHJGGR8D/rXfhM3QECufNgv6fW8Mkve69D7LYJlUYHRERkWdyes4CEVF1JIuLh+rzOVB9MhPSRo0BLy8on36uxHNEkwmi0eimCImIiDyH0z0LRETVlSAIUHTqAnmHRBjPnS11ozbttt9RtORH+I6dAEWnLm6KkoiIqOpVKFnYsGEDlixZgsuXLyMnJ8fmuCAIOHnyZEVuQURUaQSpFLL4piXWEbVaqL/6Eqa0a8h7ZRLkHTvBd+wEyJrEuSlKIiKiquN0sjBv3jzMnj0bderUQdu2bREQEODKuIiIPELRmpUwpV2zfK8/sA85B/fDq19/KJ8eC2k4J/ETEVHN5XSysHTpUiQmJuKbb75x6U7OaWlp+Prrr3H8+HGcPn0aGo0G27ZtQ1RUVKnn9urVCykpKTblc+fORZ8+fVwWIxHVDqa8XBT9+J3tAVGEdtNGaLdvg88jj8Jn1OOQ+Pm5P0AiIqJK5nSyUFhYiP79+7s0UQCAy5cv47fffkPz5s3RoUMH7Nq1q1znd+vWDRMmWC+BGBMT48oQiaiWEHNzIW3YEIYTx+1X0GlRtOQHaH5ZB+UTT8P7gQchuPiZSEREVJWcThaaNm2Ka9eulV6xnDp27Ijdu3cDAFauXFnuZCEoKAht2rRxeVxEVPtI6zdAwJffQvfndhQumAtTSrLdemJuLgpnfoqi1Svg++w4KP7Ty2ajSiJPkZ5+3aaMe+IQkSNOL5364osv4qeffnL5BGaJhKu5EpHnEAQBXj17I2jRcvhOehlCCfOzTMlXkT/1DeSOfQr6o0fcGCVR2WVkXLf5R0TkiNM9C4mJiZg2bRqGDh2KNm3aIDIy0uaNviAImD59eoWDLI8//vgDrVu3htFoRLNmzfDss89yvgIRVZggl8Pn4WHw6jcARUt+QNHKnwCdzm5dw4njyH3hGSju/g98n3sB0gYN3RwtERGRazidLBw5cgT//e9/YTAYcPDgQRw8eNCmjruThZ49e6Jly5aIiopCZmYmlixZghdeeAEfffQRHnjggRLPlUoFBAYqy33P0NA6AACJRIDJJAKAU9epLqRSSY1uXzG2s2ZxaTsDlcB/X4P+8cdwY85s5G/YAIii3aq6HX9C9/dOBDz8COq8/DIkPj6uicGByn49lUqFTVll//44+2y+dT5/x+9UFa+jK/C1rFnYzurD6WRh2rRpkMvlmDdvHjp06ACVSuXKuJzy9ttvW33ft29fDB06FJ999lmpyYLRKCInR13ue8rl5j/+SqUCarX5U0ZnrlNdBAYqa3T7irGdNUultNMnAF6vvgXp4KEo/HI29Af22a9nNKLg0L+QFRkhaCv3Z13Zr2fxM+52ztwvNNS/zHWdfTYX4++4LVe9ju7G17JmYTs9j6Nns9MTBM6cOYMnn3wSvXr18ohEwR6pVIp7770XaWlpSE9Pr+pwiKgGkjWJQ8Bns6H6ZCakjWLt1vEdNwEC52MREVE15PRfrzp16rh82dTKxJVJiKgyKTp1QeDCRfD779uQhIZZyuWdukDRIbEKIyMiInKe08nCgw8+iPXr18NgMLgyHpcyGAz47bffUK9ePYSGhlZ1OERUwwlSKbzvG4SgpaugfGYsBD9/+I6dUOp5xuSrboiOiIio/Jyes9C+fXv8+eefGDp0KEaMGIGoqChIpVKbeh07diz3tTdt2gQAOH7cvBHSjh07EBwcjODgYCQmmj+ha9asGQYPHmyZQP3LL79g27Zt6NGjByIiIpCVlYUlS5bgxIkT+Oyzz5xtJhFRuQne3lCOHgOfRx6FUMqkZt3B/cibPB5e9/SH8tmxkIZHuClKIiKi0jmdLIwZM8by9VtvvWUzzEcURQiCgFOnTpX72pMmTbL6/t133wVgXq510aJFAACj0QiTyWSpExUVhaysLHz00UfIzc2Fj48PWrRogW+++Qbdu3cvdwxERBVVWqIgmkwonDcbAKD9/Tdo/9wGn4eHwWfUE5D4l30SMBERUWVxOlmYMWOGK+OwcubMmXLXadOmDX788cfKComIyOW0W3+HMem2Z5lOh6Kli6D55WcoH38K3kMehlCN5oYRUeXhzttUVZxOFoYMGeLKOIiIahVRq4X6q3n2j+XloXD25yhavQK+z42DomcfLtJAVMvZ22mbyQK5A9fyIyKqClIpfIaPghAQ6LCKKTUF+e+8idznn4L+6GG3hUY1U/qFPPzy6REUXjdZlRdeN+GXT48g/UJeFUVGRJ6szMnCnj17nL7J7t27nT7X04WGhiM0NBwREXUtXxMRlUaQyeDz0FAE/bQGPqOeABReDusaTh5H7gvPIu+NV2C4ctl9QVKNkX4hDzsXJ0Gdo8PlHXpLwlB43YTLO/RQ5+iwc3ESEwYislHmZOHpp5/G6NGj8ccff8BoNJZaX6/XY8uWLRg1ahSeffbZCgXpycLCwhEWFo66detaviYiKiuJnx98nxuHoGWr4NV/IFDCcCPdrh3IGf0oCj79EKYbWW6Mkqqz4kTBqDcnCKIRuLxDj8xTBlzeoYd480+6UW9iwkBENso8Z2Ht2rX44IMPMHbsWAQHB6NLly5o1aoVGjRogICAAIiiiNzcXFy+fBmHDx/G3r17kZeXh7vuugvr1q2rxCYQEVV/0rBw+E+ZCp+hj6Jw3mzoD+yzX9FohGbdamg3/wafkY/BZ+iIUlddotpt/9qLlkShmGgE0o8bLYlCMaPehP1rL2Lgy63dGCERebIyJwtxcXFYuHAh/v33XyxduhTbtm3Dr7/+anfJVD8/P/Tt2xfDhw9Hq1atXB40EVFNJYuNQ8Bns6HbvxeF82bBeP6c3XpikRrqbxZAs24NVB9/AVlsEzdHStVF4pAYq56FYncmCgAglQmo08AP5/alw0clh49KAR+VHF6+ckgknGRPVBuVezWktm3bom3btjAajThx4gTOnTuHGzduQBAEBAcHo0mTJmjWrBkkEs6dJiJyliKxM+TtO0L7+29QfzMfpvR0+xVlMkjrN3BvcFSthDVSofuoJnYThttJ5RK0HdAAB3++hKtHb1gdEySAt19x8qCAUnXra0tS4S+HTGG7OSsRVW9OL50qlUrRqlUr9hwQEVUSQSqFd/+B8OrVB0Url6No8fcQCwut6iifGQvBy/HkaCIACGnoh+BIJW4kF8JoEG2OS+USNO9VDz4q+/t6iCagKE+Pojw9gEK7dQBA7i21JBCt+kYhqJ6vq5pARFXE6WSBiIjcQ/DyhnLU4/C+736of/gWmnWrAaMR0rgEePW5p6rDIw+nKzLgz+9OI+dakcM6Rr0JJ7anonFiaIXupdcYodcUIS+9CKbetknJ7Q79ehnaAv2tHgp/OXwCbn0tlXGEApEnKHOycOjQIcTExCAoKKhcNzAajThy5Aji4uLg5+dX7gCJiMhMEhQEvxdfgc/Dw1C4YC68Bz8EoZQhn9rtWyDv2AkSf5WboiRPkp+lwZ8LT9/sESiZUW9C0h7bjb+c5aNSlHg87WwuCm5oHR73UspumzehQFCoEoKXcCuxUCmg8JFyw0KiSlbmZGHkyJH46KOPMGjQoHLdIC8vDyNHjsTChQvRpUuXcgdIRETWpFH1oXr/g1LrGc4lIf//3oLg5w/l40/Ce8jDEBQlv4GjmiP9Yh52LzsHXZH95c4Fqe0kZ9EEyLwkiGoejKI83c2hRzroNaUvmW51bQHw9rU/pAkwL4ZSlF9yAqNVG6BVG5CT5rhHRCoT0KhjGNoO4LwdospS5mRBFEXk5OQgNTW1XDfIzc2FKJbcFUlERK5XOH8OIIoQ8/NQOOcLFK1eAd9nx0HRuy8/ja3hLvyTgX/WX4ZocvD3VwKEtZDaLJ8qlUvQbUQThDWy7oky6IwoyjcnDkW5Oqjz9CjKv5VMFOXpoCnQQ7w5f9rbTw6J1PHvmF5jLHGydVkZDWKpqzSd25+OE9tToAxQwMdfYbXKk7mXQgGfADnkXuylILKnXHMWpk+fjunTp5f7Jvyfj4jIvXQH90O/b49VmelaKvLffQuyFUvhO24i5G3aVVF0VFlEk4ijW5JxZleawzqCBGjYQw7fcAl8giWWjdmkcgm6j7JNFABAppDCv44U/nW8HV7XZBKhLdSjKFcHQymJQFGeruyNKoWjSdnF1Lk6aAsN0BYakA21w3oyhcRBMnHzvwHmr/mehmqbMicL48ePr9CN6tevX6HziYio7NTffuXwmOHUSeROeB6KbndD+dwLkEXHuDEyqix6rRH7Vl1A6ukcu8e9/eSAAER0AHzDzXNdfMMlaHi3HOn/CkgcEmM3USgriUQwv9n2L32omzLQCz2eiENRnh7qPOseiqI8PTSFeqCMgxJKmxtR1sTEoDMhP0uD/CyN3eMSqYCH3mnv8HzRJCI/SwMflQJyL9ctIZt+IQ/7115EWFvR8roBQOF1E37ZdKTCrxtRadyWLBARkfv4vzsN6m8WQLvpV8DBUFDdrh3Q7fkb3gPvh3LMM5DUCXFzlLXDmTMnbcri45u59B4F2Vr88c1p5KTZ/+Q8pKEfug6PhbevHCdOHLU65hsuwcCX3bsMutxLivDGAQ6Pm4wmaAr0liRCnaeHSWtCTob6VmKRr4dRbypDslD65O6yKK1XQas2YNOs4wDM7fNxsBdFcZm3rwxCKUOo0i/kWfbHuLwDaHi3uUeo8LrJ0iO0c3GSwx4hIlfg0qlERDWQNCwc/lOmwmfocBR+ORv6/XvtVzQaofl5LTS/b4Jy+GPweXQkBB8f9wZbw125csmmzJXJwo2UQuxeeg5qB5+gN2xTBx0eiK5WS5FKpBIoA7ygDLi1h0hgoBI5ObeSIVEUodcYIZWX3C5NvmuGPJWnB0OvNUKfYURehv1eCgAQJMLN5ME6qQiJ8IMyzAsFmRqrjfREI3B5h95mrolRb2LCQJWKyQIRUQ0mi22CgE9nQXdgHwrnzYLxXJL9ikVFUC/8CkXrVsP3qefgNWAgBBn/RHi65BM3sG/1RYeThVv2iUTC3XVr5Dh7QRCg8Cn9d7TvuObQFE/OztNDna+zWump6OZkbZOdzepuV9rciNJWd7qTaBKhztVBnavDnRvd3TUiFv9uvGLzuopG2ExKB8wJw/61FzHw5dblioGoLPiXgIioFlB07AT5Nz9Cu2UT1F9/CVN6ut164o0sFHw8HUUrl8F37ATIu9xVI99oVneiKOL0jms4tjXF7nGpXILEB2NQv0WwmyPzPFKZBL5BXvANcrzTuSiK0KkNlsRBfccciqI8HfxDSu5xc+2kbQUSh8RY9SxYYrWzim3x601UGZgsEBHVEoJUCu9774NXz94oWrkcRYu/h1hYaLeu8dJF5L3+EnwnvwppdAwKpr8H+YhR0McnWOrIz5zGjXfegt+UqVC06+CuZtR6RoMJB3++hMuHs+we9/aTo9vIWARHcSPUshIEAV6+cnj5yhFYV+n0dZQBChTl6x0vWVtGPio5fPx90X1UE7sJw+2kcgm6P9YEYTEcgkSVg8kCEVEtI3h5QznqcXgPfADqH76FZu0qwGj7caXg5w9JaBjyXnsJ0GoQNHcWsl+YCH18AuRnTiNo7iyYdDrkvfYSVB99xoTBTdS5OqSeyrF7LDBCiW6jYq3G+pN7NO4YhsYdwyCaRGjVBvPE7Fzdrf0pbu+tyNc73OhOkABeNze0C2ukQvOe9XDij1S7CYNULkHzXvWYKFClYrJARFRLSQID4TfpZfg8NBSFC+ZC9+d2q+OK3n2Q/+7bgNY8SVPQ6RA0dxYK7hsEv183QNDdHHah1TBhcCP/Ot7oMqwxdi4+a9kEDQDqJQSi08ONXLpsJ5WfIBHg7SeHt58cQfV8HdbTa43Q5BcvHXsrmZBKJJaN5tIv5DlMFADzXIUT21MRHOnLhIEqjUuWRtDpdLh+/Tp0OteN1yMiIveQRtWH6v0PEDD/W8hamidISiLqQr9ntyVRKCbodNaJQjGtBgXT33NXyLVeRJMAtB3Q0PJ9q95R6Do8lolCNSL3ksI/xBvhjVSIbhOCpnfXRbuBDdFjRDwA62VTS2LUm7BzURLSL+a5I2yqhSqULJw4cQKjR49Gu3bt8J///Af//PMPACArKwuPP/44du/e7ZIgiYio8smbt0TA3K/gP+1D+L34CvzefAfwst2x1yZRAAAvb3N9cpvYTmGI6xqODg9Eo/PgxpZPo6lm2L/W/ipXgp180Kg3Yf+ai26Iimojp5OFU6dOYeTIkbh69SoeeOABq2N16tSBVqvF2rVrKxwgERG5jyAI8Lq7JxR3dYeiXQeoPvrMbsJgxcsbqo8/h6Kt491tqXK06d8AjTqEVnUYVAkSh8TY7CEhSIGwFlKbhEEqlyBxCFdDosrhdLIwc+ZMhIWF4ZdffsHLL78M8Y4dQjt37oyjR486OJuIiKoDRbsOUI552nHC4OUN5ZPPMFFwMb3GiDN/p9n8baXaI6yRCt1HNbEkDILUvINzSFMZGt4ttyQMUrmEG7JRpXI6Wfjnn3/wyCOPwNfX1+4a3PXq1UO6g3W8iYioetAdOgj1d9/YzF2w0GqgXvg1dP/+497AarDCbC22f30KRzZdxck/U6s6HKpCxQmDMlCBhnfL4RtuftvmGy5Bw7vlUAYqmChQpXN6NSStVgt/f3+HxwsKCpy9NBEReQDdoYOWZVNLpNUg79XJHIrkAplXCvD30iRoCw0AgBPbU+FfxxsNWtVx2T1CQ8Nddi2qfGGNVBj4cmucOGE9WsM3XIKBL7eqoqioNnE6WWjQoAFOnDjh8PjevXsRGxvr7OWJiKiKFUx/z26iICoU9ldDmvYugletd1N0Nc+Vo1nYv/YiTAbroUf7116Eb5AX6tR3zSZrYWFMFoio7JwehjRw4ED8/PPPViseFQ9HWrhwIXbu3Gkz8ZmIiKoPvylTbeYqiAoFCgbeD1GhsK7s5W2uT+UmiiKOb0/B3pUXbBIFAJBIBGjVhiqIjIioAj0LTz75JP7++2889dRTaNSoEQRBwIwZM3Djxg1kZmaia9euGDFihCtjJSIiNypeDal4KJKoUFh2cNY3jEbQ3FnmHgYvb27I5iSD3oQDay/i6rEbdo/7qOToNioOQXWVbo6MiMjM6Z4FhUKB7777Dq+//jq8vLzg5eWFS5cuISgoCK+++ioWLFgAicQle74REVEVKU4YJOERlkQBAPTxCch+YSIk4RFMFJykKdDjr+9OO0wUgiN90ef5ZkwUiKhKOd2zAAAymQxPPPEEnnjiCReFQ0REnkbRrgOCV63HtTsmWOrjEzhHwUk5aWrsWpwEda6dDe4ARDUPQuKDMZApuCMzEVUtpz/6f+ONN3DkyBGHx48ePYo33njD2csTERFVa2f+TsPKqQdQlOxjVV5wxg+/zz3hMFFo2qMuugxtzESBiDyC08nC2rVrceXKFYfHk5OTsW7dOmcvT0REVG2d+TsNRzZdhSgCmgu+KEr2gSgC+Sf8ob9uf4M7iVRA4kMxaNknCoLEdv8iIqKqUKFhSCVRq9WQySrt8kRERB6pOFG4RYDmgi+0yT4QdfZ7C7yUMnQdEYvQho73LyIiqgrlejefmpqKlJQUy/cXLlzAgQMHbOrl5uZi2bJlaNiwYcUjJCIiqkaObr5qp1RwmCgAQO/nmsIv2H6PAxFRVSpXsrBmzRrMmTMHgiBAEATMnz8f8+fPt6kniiIkEgmmT5/uskCJiIiqg1b96t/Rs1CyFn3qMVEgIo9VrmShT58+iIyMhCiKmDJlCoYOHYq2bdta1REEAUqlEi1btkTdunVdGiwREZGni78rAgDKlDC06heFhG78W0lEnqtcyUJCQgISEsxrbKempuKee+5BXFxcpQRGRERUXUU1D8Lxrckw2tmRuVhEXAATBSLyeE7PQB4/frwr4yAiIqoRCrO1+H3e8RITBQBIO5uLM7vTEN81wk2RUXUWGhpe1SFQLVXh5YqOHTuGo0ePIjc3FyaTyeqYIAh44YUXKnoLIiKiaqHghgZbvjwBvcZUemUAR34zD1ViwkClCQtjskBVw+lkQaPRYPz48fj7778hiiIEQYAomj9FKf6ayQIREdUW+Vka/LnwdAmJggjAdv+Eo5uuMlkgIo/l9KZsc+fOxd9//43nn38eP/74I0RRxAcffICvv/4aHTp0QMuWLfHrr7+6MlYiIiKPlJdRhD++PY2iPL2DGiJkQTqYEwZrrfrVr9TYiIgqwulkYfPmzbj33nsxadIkNGnSBAAQHh6O7t2747vvvoNer8fatWtdFigREZEnyk0vwp8LT0OT7zhR8G5UCP+WefBuVIjbE4bW99a3rJ5EROSJnE4Wrl27ho4dOwIApFLzRjN6vflBKZPJcN9997FngYiIarScNLU5USgw2ByTKiSAAHg3KoRPVBEAwCeqCN6NCiEITBSIqHpwes6Cr68vjEaj5WuJRIL09HTLcX9/f2RmZlY8QiIiIg+UfU2Nv74/A53aNlHwUcnxnycT4F/HG1u2bLQ+FlWE+8f0dFeYREQV4nTPQoMGDXDp0iUA5p6F2NhYbN68GYB5B+ctW7ZwUzYiIqqRslML8dd3p+0mCsoABXo+ZU4UiIiqO6eThS5dumDz5s2W3oVhw4Zh586d6NOnD+655x7s3r0bDz30kMsCJSIi8gQ3kgvw53dnoCsy2hzzDTQnCn7BTBSIqGZwehjSs88+iwceeMCyXOrIkSOh0+mwfv16SCQSTJ48Gc8884zLAiUiIqpqmVcKsPPHs9BrbRMFv2Av9BgTD99AryqIjIiocjiVLGg0GmzZsgUxMTGQyW5dYsyYMRgzZozLgiMiIvIUGZfzsfPHszDobPdR8Kvjhf+MSYAyQFEFkRERVR6nhiEpFAq89dZbOHnypKvjQVpaGt5//30MGzYMrVu3Rnx8PJKTk8t0rslkwoIFC9CrVy+0bNkS999/v2UeBRERkbPSL+Y5TBT8Q7zR80kmCkRUMzmVLEgkEtStWxcFBQWujgeXL1/Gb7/9BpVKhQ4dOpTr3JkzZ2L27NkYOXIkvv76a7Rp0waTJk3CX3/95fI4iYiodrh+Pg87FyXZTRRUod7o+VQCfFQVTxTS06/b/CMiqmpOz1kYPHgw1q9fj8cffxwKhes+TenYsSN2794NAFi5ciV27dpVpvOysrLw7bff4tlnn8VTTz0FAOjcuTMuX76MTz75BD169HBZjEREVDukncvF30uSYDTY7rwcEO6DHmPi4e0rd8m9MjJsk4OwsHCXXJuIyFlOJwvt2rXDli1b8MADD2DEiBFo2LAhfHx8bOoVb9xWVhKJcws07dy5E3q9Hvfff79V+f33348pU6bg6tWrqF+/vlPXJiKi2ufa2Rz8vewcTHYShcAIJXo8EQcvFyUKRESeyulk4faJzNOmTYMgCFbHRVGEIAg4deqU89GVw7lz56BQKNCwYUOr8iZNmgAAzp8/z2SBiIjKJPV0Dnb/dA4mo22iEFRPibsfj4eX0uk/oVZ0hw6iYPp7kI8YBX18gqVcfuY0brzzFvymTIWiXfmG5VaGM2ds5ynGxzergkiIyJ2cftLNmDHDlXFUWG5uLlQqlU3SEhAQAADIycmpgqiIiKi6ST6ZjT3Lz0M02SYKwVG+uHt0HBQ+rksU8l57CdBqEDR3FrJfmAh9fALkZ04jaO4smHQ65L32ElQffVblCcOVK5dsypgsENV8Tj/thgwZ4so4qpxUKiAwUFmB8yUVOr+6YDtrFrazZqnsdiqVtvPTKvvn6uyz+dq1awAAnS4Pxpu9A3Xr1i31vAv/ZmDP8nMQbecyIzxGhf5jW5Y7UVAobOsHBiqh3r8fWa+bEwUAEHQ6BM2dBd2QIVCsXQtBpzNX1mqQ//pLqDt3HpSJiXbv4Y7fcUftcBf+f1yzsJ3Vh2s+GvEAKpUKeXl5luFPxXJzcwEAgYGBJZ5vNIrIyVE7ff/AQGWFzq8u2M6ahe2sWSq7nWq1zqbMmfuFhvqXua6zz+YLFy4DMCc4xXH7+ASUeM6Vo1nYt/qC3UQhpKEfuo6MhVqrg1pr+3MoiU5nsCnLyVHjxpQpEDUaq3JBp4NixUoIdwQhajRImzIFwavW272HO37HHbXDXfj/cc3CdrqOvZXTnFkcwdGzucYkC02aNIFOp8OVK1es5i2cO3cOANC4ceOqCo2IiDzc5cOZ2L/mIkTbkUcIjfFHt5FNIPeSuvSeflOmWoYg3e7ORKGYqSAfOc8+AUndSEjr1oO0Xj1I6taDtF4kRN8Yl8ZGRNVHZa+kVmOShe7du0Mul2PDhg0YP368pXz9+vWIi4vj5GYiIrLr4qEMHFh3CbCTKIQ1UqHbyFjIFK5NFABA0a4DfJ8Zi8K5X8BulnKnwkIYTp0ETtlONM6WSCAJDbMkD9K6NxOJm98LderYzOkrK/XyJVDPm43Q3n2Q0aWrpTx0z25k/u89KMdNgHLYSKeuTUSezyOThU2bNgEAjh8/DgDYsWMHgoODERwcjMSb4zWbNWuGwYMHY/r06QCAOnXq4IknnsCCBQvg6+uLZs2aYePGjdi7dy++/PLLqmkIERF5tPMHM/DPz5fsHouIVaHriCaQyZ1b0rskxrRrKJw3C7o/trnmgiYTTNfTYLqeBsPhQ7bHFV6Q1q1rSSYsicTNpELib3/4gXr5EqjnzAQARG75HRAEZHTugtA9u83fA5bjTBiIaiaPTBYmTZpk9f27774LAEhMTMSiRYsAAEajESaTdVft5MmToVQq8eOPPyIjIwMxMTH44osv0LNnT/cETkRE1ca5/ek4tOGy3WN14wLQ9dFYSF2cKAh6HcL//hvZH8wAdFqXXrtEOi2Mly/BePkS9Pbi8ldZJxD16sGQdBbaDetu1QEQ+ftm+F84D9W5c7i9n4IJA5H7uWvZZY9MFs6cOeNUHalUinHjxmHcuHGVERYREdUQSXuv499fr9g9Vi8hEF2GNYZU5ppEoUGDaEAUoTiwH8pVyyG9ccMl13UlMT8Pxvw8GM+eLrGeANgkCsXUc2cBajUgSACJBJBKAEECQSIBpFJAECz/FSRSyLt0hTQ0zOG99CeOQczNNZ8jkUCtUkKn1pmvJ5ECEgGQSG9+f9s/qdQ85OrmeYKXFyRBwY7bbjIBJpO5rpMbwxK5mzuXXfbIZIGIiKiynPk7DUc2XbV7LLJZEDo/0shliQIANJbIUDDzUxiO/Fum+qJCcWvZ1NsIIaHwf+v/YEpNgfFaKkzXUmG8+U90YwLicOaDKEK98OsyX0f1+ZwSkwX1V19Cf+ig5fu8Ml/ZmqxlawTOcxyX4dgR5I5/7laBJfEwJySCRHrze+FmmQSCJfm5WU8QIEil8PvvW5C3aOXwXoVzZ8Jw9owlkbEkKLfdT+OtgN4gWpIt6/uZz5NGRsHn4WEO72O8fh3arZvMydRt59ne0147revJmjWHJCDQ8b2Sr0LUFFl+NrbXl9y89s2f3812mRSAqNXcSi6ZrJXZ7YkCcGvZ5YL7BsHv1w1Wyy67ImGoULLw77//YvHixbh8+TJycnIg3jFBSxAEbN26tSK3ICIicplTO67h2JZku8eimpsTBYnUNW9YTNnZUH87H5oNP5s/uS4DUaFAwcD74ffLeuuEwcsb/m+/a/6D376jzXkqBXDjzHmYUm8lELeSimsQ1YUuaZNLlfLGUCzjz8zl9ynuaYB5qVh7U88dTUcX7SR5tzOcOQ39v/+UWKcsC/PKWrcpMVkwpaZAPX9uGa5UOtUXc6Gw8ztXrOCj6aW2yZ4701tZq9YInOs4qdMfP4q8Ka/Z9iRZkpPihERys8fJOjm5PYnxHTcRsrh4h/dSL/4BxssXbyWHdu9pmwAJt/WoQSKBNDwCeHiww/uYbmRBd2A/BOnN6woSSwLqOMmTIv/dt2xXUdPprBOFYloNCqa/53DZ5bJwOllYt24d3njjDchkMkRHR5dpsxsiIqKqcvLPVBzflmL3WINWwUh8sBEkUudWDLqdaDBAs3YV1Au/hliQ77iiTAbF3f+B7u+dgFYLUaGwDCXQN4xG0NxZ5j/8Xt6lfjIoUSohi2kMxNguEy6KIsS8XBhTU2G6dluvxM3EwpR2DTDY7qFQ2Ur9FNlkdNN9XJSUABCE0hITV7Wp5NW5RAfL71bKvVyW1JVyH70eYvYNh4laeZgKC0o8rtu3x/5CAeUka9OuxGTBcPkSCv73ToXvU8xejyS8vOH3ZsXu4XSy8OWXXyImJgbfffcdwsNdt5YrERGRK4miiIwTRmQct58oNGxdBx0fjIFEUvFEQXdgHwpnfQbjpYsl1pN3uQt+41+EtEFDyyTFzNsmKerjE5D9wkSELF1c4UmKgiBACAg0DyVp2szmuGg0wpSZaU4kihOIm0mF4dw5oLJ6Jdz1Jt6NyQJK65UyueKtLtzbptL+v3DRvUpN6oyuSbQAlL6MsMuSulLu48I22eXlDdXHn0PRtn2FLuN0spCamorXXnuNiQIREXmM9At52L/2IsLaivANl0AURaQfMyLzpP0/ytHtQtDhgegKJwrGlGQUzvkCul07Sqwnrd8AvhNfgqLzrf0KFO06IHjVelw7cdSqrj4+oUJDB8pKkEohDQ+HNDwc8jbtLOXq5UtgOHqkzNeRJ3aGrFkL8xstowkQzcN5RKPRvI+EyWR+cySKEI1GSILrlHg9WfOWEPz8zeeYTJBJBei1+pvDhIzmN94m460Jyrf9E4vrGE2QBDue3AzAHJtU6po3bqX0LLjqTWipyYIr34SW8om/3S3PnbqPOxOgUtrkpt6SMu2v4iwvbyiffKbCiQJQgWQhIiICulLG5hEREblL+oU87FycBKPehMs7gAZ3y5CVKTpMFBp1CEX7QQ1L//SvBKJaDfWi71G0fAmgt7coqZng6wufJ56Gz0NDIcjlTt/PndTzZtstF2F/krP+4H4EfDrLZff3G/+i1feBgUrk5Khddv1iisTOCPlzDwBzL1RxcgLRBBhvSzwcJSdGo+U8ad16Jd7L98VXIBYU3LqHJbG5dQ+ltwyF+UWWZAtGk3lYUXGyZTJBEh5R4n0kdevBZ/go61hvxmgdv/med7bx9rgElarEe0nrRUFUF9mcV/zzs9fG4uOi0WQ5r9RkwZVvrN3VA1RqD4YLE6A7aTVQL/wasqbNqq5n4dFHH8WGDRvwxBNPQCp1/c6WREREZXV7ogAAohG4/KfB4UzUegmBFUoURJMJ2t83QT1/DkxZmY4rCgK8BgyC77NjS/0kPTTUs3rqleMmWPZPKCYCyIttAtW5JJuEQTlugttiqyyCIAAy67dGFR+cdos8wXYY2J38A5UwVjApkjVoCNm4iRW6Rln5v/2uU+eVN/mTt22P4HUbbyYet5Ia62TETu+SnSRP2jCmxHspn3gKppxsh8mavXuKdyRbMJogrd+gxPtI6tSBV99+5qTp9nuIonUidUd7TPl5MF25XHoCpdUg79XJFR6KVOZk4cCBA1bft2jRAr///jseeeQRjBgxAlFRUXaTho4dHc+gJyIicoX9ay9aEgWLEv6O5lxTO50o6E+dROHMT2A4cbzEerIWLeE76eUyvUEEgLAwD0sWbm6wVpwwiABS7ulntYNz8U9QOX4SN2SjSiUoFBDqhLjlXooud7nlPrLYOPhPfb/c5914+H67iYLdZZe1GhRMe9c9qyE99thjNhNCipdKfeutt+weEwQBp06dcjo4IiKiskgcEmPVs1ASqVyCxIdK/mTRHlNWJgq/+hLajRtKrCcJCYVy7AR49e1X+kRKD2dJGObNRkqfvsjo3AUAkNHFPOciattWKMdNYKJA5EZ+U6Za7bMAlLzsst+UqRW6X5mThRkzZlToRkRERJUlrJEK3Uc1KTVhkMol6P5YE4TFlDwO+3aiXo+ilT+h6IeFJe9XoFDAZ9gIKEc9AUGpLE/4Hk05bCSUw0bi3y0brcozunRFm6n/q6KoiGovRbsOUH30mSVhqMiyy2VR5mRhyJAhFboRERFRZQprpELznvVwdEuy3SFIUrmA5r3qlStR0O3ehYLZn8OUbH/H52KK7j3gO/5FSOtFljdsIqJyK04YKnPZ5WJOb1P5xhtv4MgRx0uqHT16FG+88YazlyciIiqX9At5OPFHqsO5Cka9iBPbU5F+Ma/UaxmuXEbuqy8i7/WXSkwUpDGNoPp8DlTTP2aiQERuVbzscnGiUKx42WVXJApABZKFtWvX4sqVKw6PJycnY926dc5enoiIqMzuXA3JEaPehJ2LkhwmDKaCAhTM+QI5ox+Ffu9uh9cR/PzhO+llBC5cDEWHxArFTkTkyZxeOrU0arUaMlmlXZ6IiMjC7mpIMO+JdOc+WEa9CfvXXMTAl1tbykSTCdqNG1D41ZcQs284vpFEAu9Bg6F8+nlIAgNdFH310KBBdFWHQERVoFzv5lNTU5GSkmL5/sKFCzZLqgJAbm4uli1bhoYNG1Y8QiIiolLYWw1JkAL12nkh5ZAW4m0Jg1QuQeKQW6sh6Y8dNS+FeuZ0ifeQtW4Lvxdfhiw2zuXxVwfx8WVbApaI3Kuy92gpV7KwZs0azJkzB4IgQBAEzJ8/H/Pnz7epJ4oiJBIJpk+f7rJAiYiIHLlzNSRBKiK0ox7KSBNCJQZkHJRDNArm1ZBGNUFYIxWMGelQfzkH2i2bSry2JCwcvi9MhKJnn2q/FCoR1TyVvUdLuZKFPn36IDIyEqIoYsqUKRg6dCjatm1rVUcQBCiVSrRs2RJ169Z1abBERESOFCcM+9dehKRBOvQKPXJzZdB7GeDbXA7TlTAkDolBaKQX1D8uhHrR94BG4/iCXl5QjhwNn+GPQfD2dls7iIg8SbmShYSEBCQkmGdcp6am4p577kFcXO3sjiUiIs8TmHMWXfe8hzP+/VAQeGuoUVD2WcTvmQ2vOvcje+MvMF1LLfE6il594DtuIqThEZUdMhGRR3N6BvL48eNdGQcREXm4yh4XW1G6QwctmxQ1+mkpLgwfAV1cE/hdvIjGy5bAZDCg6NuvSryGNLYJ/Ca+DHnbdm6KmojIszmdLNx3333o0qULOnXqhE6dOkGlKvsmN0REVP1U9rjYirg9UQAAqV6PRsuWIrN7d4T9+ScEU8lLqgoBAVA+/Ty8Bw2GIJW6I2QiomrB6WRBqVRi2bJlWLx4MSQSCRISEtCpUyd07twZHTt2hLIGbXVPRESerWD6e5ZEoZhUr0fY9u0ocUqyVArvIQ9D+eQzkPjzQy8iojs5nSysXLkSBQUF2L9/P/bt24d9+/bh+++/x3fffQeZTIYWLVqgS5cumDRpkivjJSIisuE3ZapVz0KxkhIFeYdE+E6cDFlM48oNjoioGqvQrml+fn7o1asXevXqBQDIycnBjh078NVXX+Hw4cM4cuQIkwUiIqp0inYdoProM7sJw52E4Drwe/l1KLr34FKoRESlqPAWyyaTCceOHcPevXuxZ88e/Pvvv9BqtQgJCUHnzp1dESMREVGpFO06QDnmaagXzAVE0W4deWJnqKZ/DMHLy83RERFVT04nCz/88AP27t2LAwcOoKCgAAEBAejYsSNee+01dO7cGY0bs1uXiIjcR3foINTffeMwUQAA/ZHD0J88DkXb9m6MjIio+nI6WZgxYwakUikGDhyI0aNHo1mzZuzOJSKiKnHnakgOaTXIe3UyVB9/zoSBiKgMJM6eeNddd0GhUODnn3/Gc889h1deeQUrV67E1atXXRkfERFRqeythgQAJrnctrJWg4Jp77ohKiKi6s/pnoVvv/0Wer0eR44cwZ49e7Bv3z689957MBgMqFu3Ljp37owuXbpg0KBBroyXiIjIhr3VkExyOfIGDIBq40ZI9Ppblb284TdlahVESURU/TjdswAAcrkcHTp0wIQJE7B48WIcOHAAn376KXx8fLB27Vq89tprroqTiIjIoeLVkODlDQAQFQrkjJ8E6UMPIWf8JIgKhbmilzdUH30GRbsOVRgtEVH1UeHVkLRaLQ4ePIi9e/di7969OHnyJIxGo2WjNiIiIncoThgKpr+HzBGjoI9PgByAPj4B2S9MRMjSxfCbMpWJAhFROTidLMyZMwd79+7FkSNHYDAYIIoiYmNjMXz4cHTu3BmJiYlQqbgbJhERuY+iXQcEr1qPayeOWpXr4xMQvGp9FUVFRFR9VShZqF+/PgYPHozOnTujc+fOqFOnjitjIyIiIiKiKuR0srB9+3bUq1fPlbEQERG5RGhoOADA398b+fmlLKdKREQOOZ0sMFEgIiJPFRZmThYCA5Xw8VFXcTRERNVXhSY4q9VqfPPNN9iyZQuSk5MBAFFRUbjnnnvw1FNPQalUuiRIIiIiIiJyP6eThZycHIwcORLnz59HcHAwmjZtCgC4dOkS5s6di02bNmHJkiUIDAx0VaxERERERORGTicLs2bNwoULF/D222/j0UcfhVQqBQAYjUYsX74c//vf/zBnzhy89dZbLguWiIiIiIjcx+lN2bZv345HHnkEI0eOtCQKACCVSjFixAg89NBD2Lp1q0uCJCIiIiIi93M6WcjMzLQMPbKnWbNmyMzMdPbyRERERERUxZxOFkJCQnDq1CmHx0+dOoWQkBBnL09ERERERFXM6WShZ8+eWLVqFX766SeYTCZLuclkwvLly7F69Wr06tXLJUESEREREZH7OT3BeeLEidi9ezfeffddzJ49GzExMQCAixcv4saNG2jQoAEmTJjgskCJiIiIiMi9nO5ZCAoKwurVq/Hss88iMDAQx44dw7FjxxAUFIRnn30Wq1evRlBQkCtjJSIiIiIiN6rQpmx+fn6YPHkyJk+e7Kp4iIiIiIjIQzjVs6DRaLBu3TocOXLE1fEQEREREZGHcCpZUCgUeOutt3Dy5ElXx0NERERERB7CqWRBIpGgbt26KCgocHU8RERERETkIZye4Dx48GCsX78eOp3OlfEQEREREZGHcHqCc7t27bBlyxY88MADGDFiBBo2bAgfHx+beh07dqxQgEREREREVDWcThbGjBlj+XratGkQBMHquCiKEAShxF2eiYiIiIjIczmdLMyYMcOVcRARERERkYdxOlkYMmSIK+OwuHbtGmbMmIG///4boiiia9eumDJlCurVq1fqufHx8XbL161bh6ZNm7o6VCIiIiKiGq3cyUJ6ejoEQUBoaCgAQKvVYunSpTb1IiIi0L9//3Jdu6ioCI8//jgUCgU+/PBDAMDMmTMxevRorF+/HkqlstRrPPjggxg2bJhVWXR0dLniICIiIiKiciYLFy5cwKBBgzBp0iQ8++yzAAC1Wo0PP/wQgiBAFMVbF5bJ0LRp03K9UV+xYgWuXr2KTZs2oWHDhgDMvQX9+vXD8uXLreZJOBIWFoY2bdqUp1lERERERGRHuZZOXbNmDQICAvDEE0/YHHv99dfx448/4scff8T3338PX19frF69ulzBbN++Ha1bt7YkCgBQv359tGvXDtu2bSvXtYiIiIiIqGLKlSzs3bsXvXr1gkKhsDkWHx+PxMREJCYmonPnzujfvz/27NlTrmDOnTuHuLg4m/LY2FicO3euTNf46aef0KJFC7Ru3RqjR4/GwYMHyxUDERERERGZlWsY0uXLl8s8sblRo0bYuHFjuYLJzc2FSqWyKQ8ICEBeXl6p599///3o2bMnwsLCkJKSgm+//RaPP/44Fi5ciE6dOpV4rlQqIDCw9DkRjs+XVOj86oLtrFnYzpqlJraTz+ayqQ3trA1tBNjOmqYmtLNcyUJRUZHNJOOAgABs2LABUVFRVuV+fn5Qq9UVj7AcPv74Y8vXHTp0QO/evTFo0CB88cUXWLZsWYnnGo0icnKcjzcwUFmh86sLtrNmYTtrlurSztBQ/zLX5bO5bGpDO2tDGwG2s6apTu109Gwu1zAkf39/ZGRkWF9AIkGTJk1sdm/OzMyEv3/Z/yAAgEqlstuD4KjHoTR+fn7o0aMHjh07Vu5ziYiIiIhqu3IlC3Fxcdi1a1eZ6u7atcvu/IOSxMbGIikpyab8/PnziI2NLde1bnfn7tJERERERFS6ciUL99xzDw4cOFDqykRbt27FgQMH0K9fv3IF06tXLxw5cgRXr161lCUnJ+PQoUPo1atXua4FAAUFBfjzzz/RqlWrcp9LRERERFTblStZeOSRRxATE4MXX3wRM2fOREpKitXxlJQUfPHFF5g8eTIaN26Mhx9+uFzBDB06FJGRkRg3bhy2bt2Kbdu2Ydy4cYiIiLDaaC0lJQXNmjXDnDlzLGXffvst3nrrLWzYsAH79u3D2rVrMXz4cGRmZmLy5MnlioOIiIiIiMo5wVmhUGDBggV49tln8eWXX2L+/Pnw8/ODn58fCgoKUFBQAFEU0ahRI8yfP9/uEqslUSqV+OGHHzBjxgy89tprEEURXbp0wZQpU+Dr62upJ4oijEaj1SZwMTEx2LJlC7Zs2YKCggL4+fmhbdu2mDZtGnsWiIiIiIicIIi3v+MuI61WixUrVmDz5s1ISkpCYWEhfH19ERcXh379+uGRRx6Bl5dXZcRbafR6I1fcKAO2s2ZhO2uW6tLO8qyGxGdz2dSGdtaGNgJsZ01Tndrp6Nlcrp6FYl5eXnjsscfw2GOPVSgoIiIiIiLyXOWas0BERERERLUHkwUiIiIiIrKLyQIREREREdnFZIGIiIiIiOxiskBERERERHYxWSAiIiIiIruYLBARERERkV1MFoiIiIiIyC4mC0REREREZBeTBSIiIiIisovJAhERERER2cVkgYiIiIiI7GKyQEREREREdjFZICIiIiIiu5gsEBERERGRXUwWiIiIiIjILiYLRERERERkF5MFIiIiIiKyi8kCERERERHZxWSBiIiIiIjsYrJARERERER2MVkgIiIiIiK7mCwQEREREZFdTBaIiIiIiMguJgtERERERGQXkwUiIiIiIrKLyQIREREREdnFZIGIiIiIiOxiskBERERERHYxWSAiIiIiIruYLBARERERkV1MFoiIiIiIyC4mC0REREREZBeTBSIiIiIisovJAhERERER2cVkgYiIiIiI7GKyQEREREREdjFZICIiIiIiu5gsEBERERGRXUwWiIiIiIjILiYLRERERERkF5MFIiIiIiKyi8kCERERERHZxWSBiIiIiIjsYrJARERERER2MVkgIiIiIiK7mCwQEREREZFdTBaIiIiIiMguJgtERERERGQXkwUiIiIiIrLL45KFa9euYeLEiWjfvj3atWuH8ePHIzU1tUznarVafPjhh+jWrRtatWqFYcOG4cCBA5UcMRERERFRzeRRyUJRUREef/xxXLhwAR9++CE++ugjXL58GaNHj4ZarS71/ClTpmDlypWYOHEiFixYgNDQUDz11FM4deqUG6InIiIiIqpZZFUdwO1WrFiBq1evYtOmTWjYsCEAID4+Hv369cPy5csxZswYh+eePn0av/zyC6ZPn46HHnoIANCxY0fcd999mDlzJubPn++WNhARERER1RQe1bOwfft2tG7d2pIoAED9+vXRrl07bNu2rcRzt23bBrlcjgEDBljKZDIZ7rvvPuzatQs6na7S4iYiIiIiqok8Klk4d+4c4uLibMpjY2Nx7ty5Us+NjIyEj4+Pzbl6vR6XL192aaxERERERDWdRw1Dys3NhUqlsikPCAhAXl5eqecGBATYlAcGBlqOl0QulyI01L/swdpR0fOrC7azZmE7a5aa1k4+m8uuNrSzNrQRYDtrmureTo/qWSAiIiIiIs/hUcmCSqWy24PgqMfhznPt9R7k5OQAgN1eByIiIiIicsyjkoXY2FgkJSXZlJ8/fx6xsbGlnpuSkoKioiKbc+VyudWkaSIiIiIiKp1HJQu9evXCkSNHcPXqVUtZcnIyDh06hF69epV6rl6vx6ZNmyxlBoMBGzduRLdu3aBQKCotbiIiIiKimkgQRVGs6iCKqdVqPPDAA/D29sakSZMgCAJmzpyJwsJCrF+/Hr6+vgCAlJQU9O3bF+PGjcP48eMt50+ePBm7du3Ca6+9hqioKCxbtgx//PEHfvrpJzRv3ryqmkVEREREVC151GpISqXy/9u797ia0v0P4J9KbinKRIfkVmun2y6pSBqFaST3qZxDDtHIcU25zmuOg/MSI4UMSkxicpBL6hBOZhBdxjg4aBgdql2jQaNUql09vz/67XVse5ddU23t/X2/Xr28evaz1vOsZ60+y7puHDp0CKGhoVi1ahUYYxgxYgTWrVvHHygAAGMMtbW1ePc4JzQ0FBEREdixYwdKS0thbm6OmJgYOlAghBBCCCGkGT6oKwuEEEIIIYSQD8cH9cxCe/PLL79g6dKlsLe3x9ChQ7F48WIUFhYqu1sKe/bsGTZt2gRfX18IhUIIBAKIRCKZelVVVdi6dStcXFxgY2MDX19f/PDDDzL16urqEBUVBXd3d1hbW2PSpEm4cOFCWyxKg1JSUrBkyRK4ubnBxsYGHh4e2L59O8rKyqTqlZSU4IsvvoCTkxNsbW0xZ84cPHz4UGZ+io5FW7t27Rpmz56NkSNHwsrKCq6urli2bJnMlxkqus0qOh7KNm/ePAgEAkREREiVt/f1mZmZCYFAIPMzbNgwqXrtfTlbS3vOZsrl/2nv27e65jKgmtms1rnMSLNUVFSwcePGsQkTJrBLly6xS5cuMS8vLzZmzBhWXl6u7O4pJCMjg40YMYLNnz+f+fv7M47jWH5+vky9FStWMHt7e3bs2DF248YNtmjRImZtbc0ePHggVS88PJxZWlqymJgYlp6ezr788ksmEAjY999/31aLJMPb25stXbqUJSYmsszMTPbNN98we3t75u3tzWpraxljjNXV1bEZM2awUaNGsaSkJHblyhU2c+ZM5ujoyH755Rep+Sk6Fm0tKSmJbdmyhZ0/f55lZmay06dPM09PT2ZnZ8dEIhFjTPFttinjoUxJSUls5MiRjOM4Fh4ezperwvrMyMhgHMexuLg49u9//5v/uXv3Ll9HFZazNbT3bKZcplxuz7nMmOpmszrnMh0sNFNsbCwzNzdnT58+5cvy8vLYkCFD2MGDB5XYM8VJQpkxxo4fPy53p5Sdnc04jmMJCQl8mVgsZp988glbsGABX/bixQtmaWnJdu7cKTX97NmzmZeXVystwfu9fPlSpuz06dOM4zh248YNxhhjly5dYhzHsfT0dL5OaWkpc3BwYJs2beLLFB2LD0VOTg7jOI4dOHCAMab4NqvoeCjTq1evmLOzM0tKSpLZIanC+pTslK5fv95gHVVYztbQ3rOZcplyub3mMmOqnc3qnMt0G1IzXb58GUKhUOr7G/r164ehQ4ciNTVViT1TnKbm+1d/amoqtLW14enpyZd16NABEyZMQFpaGqqrqwHUX24Vi8WYNGmS1PSTJk3Co0ePpF6H25YMDAxkyqytrQEARUVFAOrXZa9evTB8+HC+jq6uLtzc3KTWpaJj8aHo0aMHAEBLSwuA4tusouOhTGFhYTAzM4OXl5fMZ6q6Pt+lLsvZVO09mymXKZfbay4DlM2quox0sNBMjx8/BsdxMuWmpqYy9yO2Z48fP0bfvn3RpUsXqXJTU1OIxWLk5uby9Tp27Cjz5XdmZmYA6r8c70ORlZUFABg8eDCAxtdlYWEhysvL+XqKjIUy1dbWorq6Gk+fPsX69ethaGjIh7ai26yi46EsN2/exJkzZ/DXv/5V7ueqtD5DQkIwZMgQODk5ITg4WOo+ZlVazpakDtlMudy+tm91yGVAfbJZHXP5g3p1antSUlICPT09mfLu3bujtLRUCT1qHSUlJejevbtMueTsSElJCf+vnp4eNDQ0pOpJpn316lWr9lNRRUVF2LVrF5ydnfkzWSUlJejbt69MXckylpaWQkdHR+GxUCZvb2/cv38fANC/f38cOnQIPXv2BKD4NqvoeChDdXU11q9fD39/fwwaNEhuHVVYn7q6uvD394eDgwO6deuGBw8eICoqCllZWThz5gx69uypEsvZGtQhmymX29f2req5DKhHNqtzLtPBAlEb5eXlWLhwIbS0tBAaGqrs7rSKbdu2oaysDPn5+Th48CDmzp2L+Ph4GBsbK7trLSImJgaVlZVYuHChsrvSqiwsLGBhYcH/7ujoCAcHB3h7eyMuLg5BQUFK7B0hLYdyWTWoQzarcy7TbUjNpKenJ/csVUNnCdorPT09uUe4kjNSkiNjyXiwd762QzKt5GhZWSorKxEYGAiRSIQDBw7AyMiI/6yhdSlZRsn6VHQslGnw4MEQCoXw8vJCbGwsKioqEB0dDUDxbVbR8WhrhYWF2LdvH5YtW4bq6mqUlpby/ZT8Xltbq1Lr822WlpYYMGAA7t27B0C1ttuWpA7ZTLncvrZvVc5lQL2zWV1ymQ4WmsnU1BQ///yzTHlOTg5MTU2V0KPWYWpqioKCArx580aqPCcnB9ra2vy9sGZmZqiurkZeXp5UPck9l5L7UJVBLBZj6dKluHfvHqKjoyEQCKQ+b2xd9unTh7+0q+hYfCj09PRgYmLCrxNFt1lFx6Ot5efno6qqCitXroSDgwP/AwAHDx6Eg4MDHj16pLLr813qspxNpQ7ZTLncfrdvVctlgLL5baq6jHSw0Ezu7u64c+eO1NskRCIRbt26BXd3dyX2rGW5u7tDLBYjJSWFL6upqcG5c+fg4uKCjh07AgBGjRoFbW1tJCUlSU1/9uxZcByHfv36tWm/Jerq6hASEoKMjAzs2bMHtra2MnXGjBmDoqIi/gE7ACgrK8N3330ntS4VHYsPxYsXL/DkyROYmJgAUHybVXQ82tqQIUMQFxcn8wPUv90lLi4OJiYmKrs+//Of/+DJkyewsbEBoLrb7e+lDtlMudx+t29Vy2VAvbNZXXKZnlloJh8fH3z77bf4y1/+gmXLlkFDQwM7d+6EkZERfH19ld09hUk2VMkltKtXr8LAwAAGBgZwdHSEhYUFPD09sXnzZtTU1MDY2BhHjx6FSCRCWFgYP5+ePXtizpw5iIqKgo6ODiwsLHDu3DlkZGRg7969Slk2ANiwYQNSUlIQGBiILl264Pbt2/xnRkZGMDIygru7O+zs7LBy5UqsWrUKenp6iI6OBmMM8+fP5+srOhbKsGjRIlhYWEAgEKBbt254+vQpYmNjoaWlhblz5wJQfJtVdDzamp6eHpycnOR+1qdPH/4zVVifwcHBMDY2hqWlJXR1dZGdnY2oqCj07t0bfn5+AFRjOVuDKmQz5TLlcnvJZUB9slmdc1mDvXszI1FYYWEhQkNDcf36dTDGMGLECKxbt65dPbT07qVfCUdHRxw+fBhA/X2lERERSE5ORmlpKczNzRESEiITDrW1tYiKisKJEyfw/PlzDBw4EIsWLcKnn37a6svREHd3dxQUFMj9bPHixViyZAmA+vsEt27ditTUVFRVVcHW1hZr166Fubm51DSKjkVbi46ORkpKCvLy8iAWi2FkZAQnJyd8/vnnUtujotusouPxIRAIBAgMDJR6uKy9r8+oqCgkJyejsLAQlZWV+Oijj+Dq6oolS5agV69efL32vpytpb1nM+Uy5XJ7z2VA9bJZnXOZDhYIIYQQQgghctEzC4QQQgghhBC56GCBEEIIIYQQIhcdLBBCCCGEEELkooMFQgghhBBCiFx0sEAIIYQQQgiRiw4WCCGEEEIIIXLRwQJp1yIjIyEQCCASiZTdFUIIIf+PspkQ1UEHC2ooMzMTAoGA/xkyZAgcHBzg5eWF1atX4+rVq6Cv32h52dnZiIyM/KB2npGRkfjXv/6l7G40S2v1PTY2FqdOnWrx+b5NJBIhMjIS2dnZrdoOaV8om5WDsrllUTarHjpYUGNeXl746quvEBoaiuXLl2P48OHIyspCQEAA/P39UVpaquwuvtfChQtx9+5d9O3bV9ldea/s7Gzs3r27wW8uVYbdu3e32x1Sa/U9Li4Op0+fbvH5vq2goAC7d++mHRKRi7K5bVE2tyzKZtXTQdkdIMpjYWGByZMnS5WtXbsW27ZtwzfffIMVK1YgJiZGSb1TTIcOHdChg2puxrW1taiurkaXLl2U3RXSQsrKytCtWzdld4N84CibP2yUzaqHsrlxdGWBSNHS0sKaNWtgb2+Pa9eu4ebNm/xnRUVF2LJlCyZPngwHBwdYW1vD09MT0dHRqK2t5etdunQJAoEAx48fl9vGhAkTMG7cuPdeTv/+++8xa9YsODk5wcbGBqNHj8bixYvx5MkTvo68+2IlZf/9738RHh4OV1dXWFlZYdKkSbhy5Yrcti5cuAA/Pz8MGzYMQqEQHh4e+Pvf/47q6mq+DmMM8fHxmDZtGoRCIezs7ODn54eMjIzGB/X/+7R27VoAwOzZs/nbDNasWQMAOHXqFAQCAW7cuIGvv/4aY8eOhY2NDc6fP9/ktr/99lv4+/tj1KhRsLKygouLC0JCQqTGSCQSQSAQAABOnz4tdeuDhKR/6enp8PX1hVAohKurK6KjowEAJSUlWLduHUaMGAGhUIgFCxagqKhIpj+vX7/Gtm3bMG7cOFhZWWH48OFYsWIF8vPzpepJxiA9PR0HDhzA2LFjYWVlBQ8PD6mzSYr0XZ66ujrExsZi4sSJsLOzw9ChQ+Hh4YF169ZBLBbzy1xQUICsrCyp+UrGLi0tDcuXL8eYMWNgY2ODYcOGwd/fH1lZWTLt+fn5wd3dHfn5+Vi6dCkcHR1hb2+PU6dOYfbs2QDq/wMoacPPz6/R/hP1RtlM2SxB2UzZ3NZU87Cf/G6fffYZfvzxR1y5cgXDhg0DADx8+BAXL17EuHHjYGJiArFYjGvXrmH79u0QiUTYuHEjAMDNzQ2GhoY4efIkfHx8pOZ7+/ZtPH78GEFBQdDQ0Giw/aysLCxcuBBmZmZYsGABdHV18euvvyI9PR15eXkYOHDge5dhzZo16NChA/z9/SEWi3Ho0CEsWrQIKSkpMDY25utFRERg3759MDU1xZw5c2BoaIi8vDxcvHgRS5cuRceOHQEAK1euxD//+U94eHhg2rRpqK6uRlJSEvz9/REZGYkxY8Y02Jdx48bh+fPnOHbsGAIDAzFo0CAAgImJiVS9rVu3oqamBj4+PtDR0eGXsyltHzx4ELa2tvDz80OPHj3w6NEjJCQkICMjA0lJSdDX14eBgQG++uorrFq1CsOGDZNZTxIPHjzAd999Bx8fH0yePBnnz5/H9u3b0alTJ5w5cwZ9+/bF4sWLkZeXh8OHD2P16tWIjY3lp3/9+jVmzJiBwsJCTJ8+HWZmZnj+/Dni4+Ph7e2NkydPytymEBERgcrKSvj6+qJjx444evQo1qxZAxMTE9jb2yvc93ft3bsXu3btgpubG2bMmAEtLS2IRCJcvnwZ1dXV0NbW5m/90NfXR2BgID+tgYEBgPodYElJCaZMmQIjIyMUFRXhxIkTmDNnDuLi4vi/FYny8nLMmjULQ4cOxfLly1FcXAwHBwcEBgZi37598PX1hb29PQDgo48+Umg5iHqjbKZsBiibJSib2wgjaicjI4NxHMdiYmIarHPv3j3GcRxbvHgxX/bmzRtWV1cnUzckJISZm5uzoqIivmz79u2M4zj2888/S9X94osv2JAhQ9izZ88a7ePmzZsZx3HsxYsXjdbbtWsX4ziO5efny5R9/vnnUv29c+cO4ziOhYWFyZT5+fmxyspKqXnX1dXx01+8eJFxHMf+8Y9/SNURi8Vs6tSpzM3NTe7YvO3kyZOM4ziWkZHR4GeffPIJq6iokPqsqW2Xl5fLzP/GjRuM4zgWHR0tVc5xHFu9erXc/nIcxwQCAbt9+zZfVlVVxUaOHMkEAgHbtGmTVH3JOsvJyeHLNm3axKytrVl2drZUXZFIxOzs7KTalozB5MmTWVVVFV/+7NkzZmlpyYKCghTuuzxTpkxh48ePf289Nzc3NmvWLLmfyRvb58+fM0dHRzZ//nyp8lmzZjGO41h4eLjMNJK/wZMnTyrYe6IOKJspmyUom2VRNisP3YZE5JLcu1dWVsaXde7cmT/jVF1djVevXqG4uBguLi6oq6vDvXv3+Lre3t7Q0NBAQkICX1ZRUYFz587B1dUVvXv3brR9XV1dAPWXoGtqapq1DLNnz5Y6Q2ZjY4OuXbsiNzeXLzt79iwAIDg4GJ06dZKaXkNDg5/+7Nmz0NHRwdixY1FcXMz/lJaWwt3dHQUFBXj69Gmz+vm2P/7xjzL3wTa17a5duwKov7T7+vVrFBcXQyAQQFdXF3fv3m1Sf2xtbSEUCvnfO3bsCGtrazDGZC7NSs7cSMaXMYakpCQ4ODigV69eUn3v0qULbG1tkZaWJtPmn/70J/6MIQD07t0bAwcO/N3j261bNxQVFUndvtFUkrEF6s9M/fbbb9DU1IRQKGxwbOfNm9fs9gh5F2UzZTNA2fwuyubWRbchEbkkO6K3H/ipqalBdHQ0EhMTkZubK3Nf69tv6OjXrx+cnZ2RmJiI4OBgaGtr4/z58ygvL8dnn3323vZnzpyJ1NRUbNiwAWFhYbC3t8eoUaPg5eXFX3Z8n379+smU6evr47fffuN/z83NhYaGBszNzRudV05ODsrLy+Hs7NxgnZcvXyp0Cb4x8qZvatvp6enYs2cP7ty5g6qqKql6JSUlTeqPvDHs3r07AEjdLgAAenp6AIBXr14BAIqLi/Hq1SukpaVhxIgRcuevqSl7vkJemz169PjdbypZsWIFFi1ahJkzZ6JXr15wdHTE6NGj4eHhIbUDbExeXh4iIiKQlpYm80YaebduGBgY8ONCSEugbJZG2fw/lM2Uza2FDhaIXA8fPgQgHZBbtmzB4cOH4enpicDAQBgYGEBbWxv3799HWFgY6urqpObh4+ODZcuW4fLly/Dw8EBCQgIMDQ0xevTo97avr6+PhIQE3Lx5Ezdu3MAPP/yA0NBQREZGIjo6GnZ2du+dh7ywk+fts1QNYYzBwMAA27dvb7COmZmZQu01pnPnzr+r7bt372LevHkwMTFBcHAwjI2N+bOOQUFBTX5Hu5aWVpM/k7Qh+dfZ2RkBAQEKt6noemsqOzs7XLp0CWlpacjMzERmZiaSk5Oxd+9exMfHo0ePHo1OX15ejpkzZ+LNmzf485//DI7joKOjA01NTURFRcl9oJHelkJaGmWzNMpmxT+jbP4fyuamoYMFIpfkEvXHH3/MlyUmJsLBwQERERFSdd++dPy2MWPGoGfPnkhISICZmRlu3bqFgIAAhV+np6WlBScnJzg5OQEAfvrpJ0yfPh179+7l3/rwew0YMABXr17FTz/9BBsbmwbr9e/fH0+fPoVQKISOjk6z2nrfTq8l2k5OTkZtbS32798vdRaooqKizd/NLjlzU1ZW1uiZt7ako6MDDw8PeHh4AKh/O8nGjRuRkJCA+fPnNzpteno6fv31V2zevBnTp0+X+mzHjh1N6kdztwVCKJulUTY3HWVzwyib5aNnFoiU2tpabN26FT/++CM+/vhj/m0AQP1ZhXfPflRUVEi9YeFt2tramDp1KtLS0vD1118DgEKXuYH6y6TvGjRoEDp16tTky7WNmThxIgAgPDxc6lV8EpLlnTJlCurq6hAeHi53Pi9evHhvW5J7Kpva/6a03dAZpaioKJmzi5I+SS5NtzRNTU1MnDgRd+/eRUpKitw6L1++bPb8m9p3eduUpaUlAOl1oqOjI3e+krF9928gLS0Nd+7cUbgfQPO3BaK+KJulUTY3H2Vzwyib5aMrC2rswYMHSExMBFB/Ge/JkydITU1FQUEBXFxcZC6tenh44NixY1i+fDmcnZ3x4sULnDx5stFLhD4+Pjhw4ACSk5Ph6OiIAQMGKNS3L7/8Es+ePYOLiwv69OmDyspK/r7ad7+s6PewsbFBQEAA9u/fj2nTpmH8+PEwNDSESCTChQsXcOLECejp6eHTTz/FtGnTcOTIEdy/fx9ubm7Q19fHs2fPcPv2beTm5iI1NbXRtqytraGpqYl9+/ahpKQEXbt2hbGxsdRDavI0pe2xY8ciNjYWAQEB8PX1hba2Nq5fv46HDx9CX19fZt62trZIT09HdHQ0+vTpAw0NDUyYMKH5A/qOoKAg3Lp1C8uXL8f48eMhFAqhra2NwsJCXL16FZaWltiyZUuz5t3Uvnt6esLW1hY2Njbo1asXnj9/juPHj0NbW1tqOqFQiISEBOzYsQODBw+GpqYm3NzcYG9vD0NDQ2zduhUFBQUwMjJCdnY2EhMTwXEcHj16pHDfTU1NoaOjg/j4eHTu3Bl6enowMDBo8P5hol4omymbKZspmz8kdLCgxpKTk5GcnAxNTU107doVRkZGcHBwwN/+9je4urrK1F+7di10dHSQkpKC1NRU/OEPf4Cvry+sra0xZ84cuW30798fTk5OyMjIkLk82JjJkyfj1KlTOH36NIqLi9GtWzeYmppi165d/GXKlhISEgJzc3McOXIEMTExYIzByMgIrq6uUvephoaGwsnJCcePH0dUVBTEYjEMDQ1hYWGB4ODg97bTp08fbN68Gfv378eGDRsgFosxderU9+6QmtK2vb09IiMjsWfPHuzcuROdOnWCs7Mzjhw5glmzZsnMd/369di4cSP27duH8vJyAGjRHZKuri6OHj2KgwcP8tuNlpYWjIyMYG9vD29v72bPu6l99/f3x5UrV3D48GG8fv0aPXv25L+w6O2HKIOCglBSUoL4+HiUlpaCMYbU1FQYGxsjJiYG27Ztw5EjR1BTUwMrKyvs378fCQkJTdohde7cGREREdixYwc2b96M6upqODo6qv0OidSjbK5H2UzZTNn8YdBgTX2qhpAmCggIwO3bt3Ht2jW5D4kRQghpe5TNhBBF0DMLpFXl5uYiLS0NkyZNop0RIYR8ICibCSGKotuQSKu4c+cOcnJycPjwYWhra2Pu3LnK7hIhhKg9ymZCSFPRwQJpFUePHsWZM2fQr18/hIWFyXxJDCGEkLZH2UwIaSp6ZoEQQgghhBAiFz2zQAghhBBCCJGLDhYIIYQQQgghctHBAiGEEEIIIUQuOlgghBBCCCGEyEUHC4QQQgghhBC5/g9fT5uoDfAQdQAAAABJRU5ErkJggg==\n"
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
Felix Reichel's avatar
Felix Reichel committed
    "# plot variables\n",
    "fontsize = 18\n",
    "linewidth = 5\n",
    "markersize = 12\n",
    "errbar_width = 5\n",
    "xlabel = 'Day since treatment start'\n",
    "\n",
Felix Reichel's avatar
Felix Reichel committed
    "# color for the control interval\n",
    "ctrl_clr = 'darkslategray'\n",
    "\n",
Felix Reichel's avatar
Felix Reichel committed
    "# compute control fit values\n",
    "results_ctrl = np.loadtxt(results_ctrl_file)\n",
Felix Reichel's avatar
Felix Reichel committed
    "\n",
    "v_ctrl = results_ctrl[:,0]\n",
    "probs_ctrl = results_ctrl[:,3]\n",
    "probs_ctrl_err = results_ctrl[:,4]\n",
    "\n",
Felix Reichel's avatar
Felix Reichel committed
    "v_min = 0.\n",
    "v_max = 3.\n",
    "binsize = 0.25\n",
    "bins = int(v_max / binsize)\n",
    "\n",
    "ind_vmax = v_ctrl <= v_max\n",
    "v_ctrl = v_ctrl[ind_vmax]\n",
    "probs_ctrl = probs_ctrl[ind_vmax]\n",
    "probs_ctrl_err = probs_ctrl_err[ind_vmax]\n",
    "\n",
Felix Reichel's avatar
Felix Reichel committed
    "fit_bounds = [0, np.inf]\n",
    "popt_ctrl, pcov_ctrl = curve_fit(asymptotic_exponential_growth, v_ctrl, probs_ctrl,\n",
    "                                 sigma = probs_ctrl_err, absolute_sigma=False,\n",
Felix Reichel's avatar
Felix Reichel committed
    "                                 bounds=fit_bounds\n",
    "                                 )\n",
    "perr_ctrl = np.sqrt(np.diag(pcov_ctrl))\n",
    "\n",
Felix Reichel's avatar
Felix Reichel committed
    "# limits of the 95% confidence interval\n",
    "ci_lower = float(popt_ctrl - perr_ctrl)\n",
    "ci_upper = float(popt_ctrl + perr_ctrl)\n",
    "\n",
Felix Reichel's avatar
Felix Reichel committed
    "with sns.axes_style('darkgrid'):\n",
    "\n",
Felix Reichel's avatar
Felix Reichel committed
    "    fig = plt.figure(0,(11,6))\n",
    "\n",
    "    plot_titles = ['Normocytes', 'Acanthocytes']\n",
    "    para = 'lambda'\n",
    "    para_label = r'$\\lambda$'\n",
Felix Reichel's avatar
Felix Reichel committed
    "    ylim = [0, 2.7]\n",
    "\n",
    "    for jj, patient in enumerate(patients):\n",
    "        data = dict_fitvalues[patient]\n",
    "        color = color_dict[patient]\n",
    "\n",
Felix Reichel's avatar
Felix Reichel committed
    "        for n, health in enumerate(['healthy', 'unhealthy']):\n",
    "            ax=plt.subplot(1,2,n+1)\n",
    "\n",
    "            df_plot = data[health]\n",
    "            xdata = df_plot['days']\n",
    "            ydata = df_plot[para]\n",
    "            yerr = df_plot[para + \"_err\"]\n",
    "\n",
    "            if patient=='LM':\n",
    "                # plot data on treatment\n",
    "                plt.errorbar(xdata[:-2], ydata[:-2], yerr=yerr[:-2],\n",
    "                             c=color, label=labels[jj],\n",
    "                             ls='-', lw=linewidth, marker='X', markersize=markersize,\n",
    "                             ecolor='gray', elinewidth=errbar_width)\n",
    "                # plot data off treatment\n",
    "                plt.errorbar(xdata[-3:], ydata[-3:], yerr=yerr[-3:],\n",
    "                             c=color, #label=labels[jj],\n",
    "                             ls='--', lw=linewidth, marker='X', markersize=markersize,\n",
    "                             ecolor='gray', elinewidth=errbar_width)\n",
    "            else:\n",
    "                # plot data on treatment\n",
    "                plt.errorbar(xdata[:-1], ydata[:-1], yerr=yerr[:-1],\n",
    "                             c=color, label=labels[jj],\n",
    "                             ls='-', lw=linewidth, marker='X', markersize=markersize,\n",
    "                             ecolor='gray', elinewidth=errbar_width)\n",
    "                # plot data off treatment\n",
    "                plt.errorbar(xdata[-2:], ydata[-2:], yerr=yerr[-2:],\n",
    "                             c=color, #label=labels[jj],\n",
    "                             ls='--', lw=linewidth, marker='X', markersize=markersize,\n",
    "                             ecolor='gray', elinewidth=errbar_width)\n",
    "\n",
    "            plt.ylim(ylim)\n",
    "            plt.xlabel(xlabel, fontsize=fontsize)\n",
    "            plt.tick_params(axis='both', which='both', labelsize=fontsize-2)\n",
Felix Reichel's avatar
Felix Reichel committed
    "            plt.xticks([0, 100, 200, 300, 400, 500])\n",
Felix Reichel's avatar
Felix Reichel committed
    "            plt.title(r'{} - {}'.format(para_label, plot_titles[n]), fontsize=fontsize+2)\n",
    "\n",
    "            # plot control region at end only\n",
    "            if patient==patients[-1]:\n",
    "                if health=='unhealthy':\n",
Felix Reichel's avatar
Felix Reichel committed
    "                    ax.axhline(ci_lower, ls='--', lw=.5, c=ctrl_clr, zorder=0)\n",
    "                    ax.axhline(ci_upper, ls='--', lw=.5, c=ctrl_clr, zorder=0)\n",
    "                    axis_limits = ax.get_xlim()\n",
    "                    ax.add_patch(patches.Rectangle((axis_limits[0], ci_lower),\n",
    "                                                   np.diff(axis_limits), ci_upper-ci_lower,\n",
    "                                                   color=ctrl_clr, alpha=0.15, zorder=0,\n",
    "                                                   label = 'CTRL'\n",
    "                                                   )\n",
    "                                 )\n",
    "                    ax.get_yaxis().set_ticklabels([])\n",
    "\n",
Felix Reichel's avatar
Felix Reichel committed
    "                else:\n",
    "                    ax.axhline(ci_lower, ls='--', lw=.5, c=ctrl_clr, zorder=0)\n",
    "                    ax.axhline(ci_upper, ls='--', lw=.5, c=ctrl_clr, zorder=0)\n",
    "                    axis_limits = ax.get_xlim()\n",
    "                    ax.add_patch(patches.Rectangle((axis_limits[0], ci_lower),\n",
    "                                                   np.diff(axis_limits), ci_upper-ci_lower,\n",
    "                                                   color=ctrl_clr, alpha=0.1, zorder=0,\n",
    "                                                   )\n",
    "                                 )\n",
    "            # set alpha of errorbars\n",
    "            for collection in ax.collections:\n",
    "                collection.set_alpha(.4)\n",
    "\n",
Felix Reichel's avatar
Felix Reichel committed
    "    fig.supylabel(\"Growth rate [(mm/s)$^{-1}$]\", fontsize=fontsize)\n",
    "    plt.legend(loc='upper right', ncol=1, fontsize=fontsize-2, title_fontsize=fontsize)\n",
    "    plt.tight_layout()\n",
Felix Reichel's avatar
Felix Reichel committed
    "    savename = \"fig3_growth_rate_lithium\"\n",
    "    savepath = os.path.join(savefolder,savename)\n",
    "    plt.savefig(savepath+\".pdf\", dpi=900, format='pdf')"
   ],
   "metadata": {
    "collapsed": false,
    "pycharm": {
     "name": "#%%\n"
    }
   }
   "execution_count": 7,
   "outputs": [],
   "source": [],
   "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
}