HBonds.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 509,
   "id": "71942024",
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import pandas as pd\n",
    "import seaborn as sns\n",
    "import MDAnalysis as mda\n",
    "from MDAnalysis.analysis.hydrogenbonds.hbond_analysis import HydrogenBondAnalysis as HBA\n",
    "from MDAnalysis.lib.nsgrid import FastNS\n",
    "import tqdm.auto as tqdm\n",
    "from MDAnalysis.lib.distances import calc_angles\n",
    "from MDAnalysis.lib.distances import calc_bonds\n",
    "from joblib import Parallel, delayed\n",
    "from numba import njit\n",
    "import matplotlib.pyplot as plt"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 224,
   "id": "9f49c106",
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "959b5a08",
   "metadata": {},
   "outputs": [],
   "source": [
    "u = mda.Universe(\"/home/dibyendu/Cage_Calculation/Systems/273k1bar/run1/npt1.gro\",\n",
    "             \"/home/dibyendu/Cage_Calculation/Systems/273k1bar/run1/npt1.xtc\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 225,
   "id": "d9dd52ce",
   "metadata": {},
   "outputs": [],
   "source": [
    "liq = mda.Universe(\"/home/dibyendu/DM_Projects/Sample_Structures/Ice/liquid.gro\", reset_offset = True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 537,
   "id": "10a758f6",
   "metadata": {},
   "outputs": [],
   "source": [
    "ih = mda.Universe(\"/home/dibyendu/DM_Projects/Sample_Structures/Ice/ice1h.gro\", reset_offset = True)\n",
    "ic = mda.Universe(\"/home/dibyendu/DM_Projects/Sample_Structures/Ice/cs2.gro\", reset_offset = True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 227,
   "id": "77e4f6cc",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "application/vnd.jupyter.widget-view+json": {
       "model_id": "588cbfa752164605aa8b137a215782ea",
       "version_major": 2,
       "version_minor": 0
      },
      "text/plain": [
       "  0%|          | 0/1 [00:00<?, ?it/s]"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "/home/dibyendu/Soft/anaconda3/envs/Torch/lib/python3.10/site-packages/MDAnalysis/analysis/base.py:438: UserWarning: Reader has no dt information, set to 1.0 ps\n",
      "  self.times[i] = ts.time\n"
     ]
    },
    {
     "data": {
      "text/plain": [
       "<MDAnalysis.analysis.hydrogenbonds.hbond_analysis.HydrogenBondAnalysis at 0x2b7ebac39570>"
      ]
     },
     "execution_count": 227,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "hbondsih = HBA(\n",
    "  universe=ih,\n",
    "  donors_sel='resname SOL and name OW',\n",
    "  hydrogens_sel='resname SOL and name HW1 HW2',\n",
    "  acceptors_sel='resname SOL and name OW',\n",
    "  d_h_a_angle_cutoff=150,\n",
    "  d_h_cutoff=1.2,\n",
    "  d_a_cutoff=3.5,\n",
    ")\n",
    "hbondsih.run(verbose = True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 538,
   "id": "aae36442",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "application/vnd.jupyter.widget-view+json": {
       "model_id": "489be9e74e124eb6b3232822e08d1ae1",
       "version_major": 2,
       "version_minor": 0
      },
      "text/plain": [
       "  0%|          | 0/1 [00:00<?, ?it/s]"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "/home/dibyendu/Soft/anaconda3/envs/Torch/lib/python3.10/site-packages/MDAnalysis/analysis/base.py:438: UserWarning: Reader has no dt information, set to 1.0 ps\n",
      "  self.times[i] = ts.time\n"
     ]
    },
    {
     "data": {
      "text/plain": [
       "<MDAnalysis.analysis.hydrogenbonds.hbond_analysis.HydrogenBondAnalysis at 0x2b7ed8cfd1e0>"
      ]
     },
     "execution_count": 538,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "hbondsic = HBA(\n",
    "  universe=ic,\n",
    "  donors_sel='resname ICE and name OW',\n",
    "  hydrogens_sel='resname ICE and name HW1 HW2',\n",
    "  acceptors_sel='resname ICE and name OW',\n",
    "  d_h_a_angle_cutoff=150,\n",
    "  d_h_cutoff=1.2,\n",
    "  d_a_cutoff=3.5,\n",
    ")\n",
    "hbondsic.run(verbose = True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 228,
   "id": "cb11169f",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "application/vnd.jupyter.widget-view+json": {
       "model_id": "77bc5cd022f44ffebe45b4a0e5f0e5c2",
       "version_major": 2,
       "version_minor": 0
      },
      "text/plain": [
       "  0%|          | 0/1 [00:00<?, ?it/s]"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "text/plain": [
       "<MDAnalysis.analysis.hydrogenbonds.hbond_analysis.HydrogenBondAnalysis at 0x2b7ebac3bc40>"
      ]
     },
     "execution_count": 228,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "hbondsliq = HBA(\n",
    "  universe=liq,\n",
    "  donors_sel='resname SOL and name OW',\n",
    "  hydrogens_sel='resname SOL and name HW1 HW2',\n",
    "  acceptors_sel='resname SOL and name OW',\n",
    "  d_h_a_angle_cutoff=150,\n",
    "  d_h_cutoff=1.2,\n",
    "  d_a_cutoff=3.5,\n",
    ")\n",
    "hbondsliq.run(verbose = True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 230,
   "id": "643eec05",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<AxesSubplot:ylabel='Density'>"
      ]
     },
     "execution_count": 230,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAjcAAAGdCAYAAADuR1K7AAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjQuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/MnkTPAAAACXBIWXMAAA9hAAAPYQGoP6dpAABhkElEQVR4nO3deXxU5dn/8c9M9oQkEEI2kkBYZN/3TRQUxaVabdVaRa1acW2lPFZsa5+2T4v91Vq1KmpFqAsuFVBbl4KyqewQ9n1LQkgIIWQhZJ/z++MkgUASsszMmZl836/XvM7JmXtmrjNkJhf3ue77thmGYSAiIiLiI+xWByAiIiLiTEpuRERExKcouRERERGfouRGREREfIqSGxEREfEpSm5ERETEpyi5EREREZ+i5EZERER8ir/VAbibw+Hg2LFjhIeHY7PZrA5HREREmsAwDIqKikhISMBub7xvps0lN8eOHSMpKcnqMERERKQFMjIySExMbLRNm0tuwsPDAfPNiYiIsDgaERERaYrCwkKSkpJq/443xtLkZs6cOcyZM4cjR44A0K9fP55++mmmTp3a4GNWrlzJjBkz2LlzJwkJCTzxxBNMnz69ya9ZcykqIiJCyY2IiIiXaUpJiaUFxYmJiTzzzDNs3LiRjRs3MmnSJG644QZ27txZb/vDhw9zzTXXMGHCBFJTU3nqqad47LHHWLhwoZsjFxEREU9l87RVwaOiovjLX/7Cvffee8F9v/zlL/n000/ZvXt37bHp06ezdetW1qxZ06TnLywsJDIykoKCAvXciIiIeInm/P32mKHgVVVVvP/++xQXFzNmzJh626xZs4YpU6bUOXbVVVexceNGKioq6n1MWVkZhYWFdW4iIiLiuyxPbrZv3067du0ICgpi+vTpLF68mL59+9bbNjs7m9jY2DrHYmNjqaysJDc3t97HzJ49m8jIyNqbRkqJiIj4NsuTm169erFlyxbWrl3Lgw8+yF133cWuXbsabH9+IVHNVbWGCoxmzZpFQUFB7S0jI8N5wYuIiIjHsXwoeGBgID169ABg+PDhbNiwgRdeeIHXXnvtgrZxcXFkZ2fXOZaTk4O/vz8dO3as9/mDgoIICgpyfuAiIiLikSzvuTmfYRiUlZXVe9+YMWNYunRpnWNLlixh+PDhBAQEuCM8ERER8XCWJjdPPfUU33zzDUeOHGH79u386le/YsWKFfz4xz8GzEtK06ZNq20/ffp00tLSmDFjBrt37+bNN99k7ty5zJw506pTEBEREQ9j6WWp48ePc+edd5KVlUVkZCQDBw7kyy+/5MorrwQgKyuL9PT02vYpKSl8/vnnPP7447z88sskJCTw4osvcvPNN1t1CiIiIuJhPG6eG1fTPDciIiLexyvnuRERERFxBiU3IiIi4lOU3IiIiIhPsXyeGxERSx34Cjb9E0pOQcfuMH4GdOhidVQi0gpKbkSkbaqqgE8egW3vnz125BtIfReufx6G3GFZaCLSOrosJSJtj2HA5/9jJjY2O4x8AL7/GqRMBEcFfPoo7P3S6ihFpIWU3IhI27NpPmyaB9jg1nfhmv8Hg26DaZ/A4DvAcMBH98CpNKsjFZEWUHIjIm1L8UlY+ltz/4rfQu9rzt5ns5mXpJLHQsUZWPIrS0IUkdZRciMibcvKZ6CsAOIGwNjHLrzfLwCu/SvY/GD3v+HgMvfHKCKtouRGRNqOvEOwYa65P+WPYPerv11sXxj5U3N/+Z/cE5uIOI2SGxFpO9b/A4wq6D4Juk1svO2EGeAXCEc3QMZ698QnIk6h5EZE2oay05D6jrk/+qGLt28XAwNuMffXvOy6uETE6ZTciEjbsO19KCuEqO7QfXLTHjOmOgna/SkUHHVdbCLiVEpuRKRt2PRPczvyfrA38asvth90GW8ODd/2oetiExGnUnIjIr4vdz9kbwO7/9lLTU016FZzu+0Dc/I/EfF4Sm5ExPftWGRuu10OYR2b99g+3wO/IDixx0yQRMTjKbkREd9mGLDjI3O//83Nf3xIe+h1tbm/9QOnhSUirqPkRkR8W84uyN1n9r6cOxtxcwz4obnd/W9dmhLxAkpuRMS37f3C3Ha/HIIjW/Yc3SeBfzAUpMPxHc6LTURcQsmNiPi2A1+Z255Xtvw5AsPMeh04myyJiMdSciMivqsk/+zswj1akdzA2Utaez5r3fOIiMspuRER33VoubncQvQl0KFL657rkqsBG2RtgYJMZ0QnIi6i5EZEfNf+6ktSre21AXM5hsTh5v7Br1v/fCLiMkpuRMQ3GQYcWmHu92jicgsXU7Nsw8Flznk+EXEJJTci4pvy06HwqDkrcfJo5zxn90nm9tAKcFQ55zlFxOmU3IiIb0r7ztwmDDFHOzlD52EQFAElpyBrq3OeU0ScTsmNiPimmuSmyzjnPaefP6Rcau7r0pSIx1JyIyK+KW21uXVmcgPmZIAAB5c793lFxGmU3IiI7ynMgrxDYLND8ijnPnfKZeb26AaoLHPuc4uIUyi5ERHfk17daxPbv+VLLjSkY3cI6wRVZZC52bnPLSJOoeRGRHzP0U3mNsnJvTYANht0GWvu19T1iIhHUXIjIr4nc6O5rZl0z9lq6njS17jm+UWkVZTciIhvqao4O0y78zDXvEZNz036OqiqdM1riEiLKbkREd9yfCdUlpq1NlHdXfMaMX3N5y8vguPbXfMaItJiSm5ExLdkVtfbJAwFu4u+4ux+kFQ963HNkHMR8RhKbkTEt9QkN66qt6lRW1Ss5EbE0yi5ERHfUpPcuKrepkZNUXHaanORThHxGEpuRMR3lJ+B3H3mfsIQ175W/CAICIWSPDix17WvJSLNouRGRHxHzi4wHBAWA+Fxrn0t/8Czl740342IR1FyIyK+I3ubuY0b4J7XO/fSlIh4DCU3IuI7sqqTm/iB7nk9FRWLeCQlNyLiO7Kr55xxV89N52Fg84OiY5Cf4Z7XFJGLUnIjIr7BUWVO4AcQ56aem8Cws71EGevc85oiclFKbkTEN5w8AJUl5gimqG7ue92axTkz1rvvNUWkUUpuRMQ31FySiu1vziDsLrXJzVr3vaaINErJjYj4hppLUrH93Pu6NclN9g4oO+3e1xaReim5ERHfkLPb3Lo7uYnsDBGJYFSdnR1ZRCyl5EZEfEPOLnPbqbf7XztZdTcinkTJjYh4v7LTkJ9m7sf0cf/r19bdaMSUiCdQciMi3q9mbaewGAiLdv/rJ400t0fXg8Ph/tcXkTosTW5mz57NiBEjCA8PJyYmhhtvvJG9extfgG7FihXYbLYLbnv27HFT1CLicU5U19tY0WsDEDvAHIJeWgC5WkRTxGqWJjcrV67k4YcfZu3atSxdupTKykqmTJlCcXHxRR+7d+9esrKyam89e/Z0Q8Qi4pFqiolj+lrz+n7+5mzFAOkaEi5iNX8rX/zLL7+s8/O8efOIiYlh06ZNXHrppY0+NiYmhvbt27swOhHxGjXFxDEWFBPXSB4NR74xi4qH32NdHCLiWTU3BQUFAERFRV207ZAhQ4iPj2fy5MksX768wXZlZWUUFhbWuYmIj7G65wZUVCziQTwmuTEMgxkzZjB+/Hj69+/fYLv4+Hhef/11Fi5cyKJFi+jVqxeTJ09m1apV9bafPXs2kZGRtbekpCRXnYKIWKG0EIqyzP3oS6yLI3G4uc07CKdPWBeHiGAzDMOwOgiAhx9+mM8++4xvv/2WxMTEZj32+uuvx2az8emnn15wX1lZGWVlZbU/FxYWkpSUREFBAREREa2OW0QslrkJ/jHJHCn1P/utjeXl0WZx820LoPe11sYi4mMKCwuJjIxs0t9vj+i5efTRR/n0009Zvnx5sxMbgNGjR7N/f/1fakFBQURERNS5iYgPyT1gbqM9YFBBzZBwXZoSsZSlyY1hGDzyyCMsWrSIZcuWkZKS0qLnSU1NJT4+3snRiYhXOFn9H5uOPayNA8yiYoB0JTciVrJ0tNTDDz/MggUL+OSTTwgPDyc7OxuAyMhIQkJCAJg1axaZmZm89dZbADz//PN07dqVfv36UV5ezjvvvMPChQtZuHChZechIhbKrU5uPKLnprqo+FgqVJaBf5C18Yi0UZYmN3PmzAHgsssuq3N83rx53H333QBkZWWRnp5ee195eTkzZ84kMzOTkJAQ+vXrx2effcY111zjrrBFxJOcrL4s1dEDkpuobhAaDWdyIWvr2ctUIuJWHlNQ7C7NKUgSEQ/ncMCfEqCyBB7dDB27Wx0RvHc77P0MpvwfjH3U6mhEfIbXFRSLiLRI4VEzsbEHQPsuVkdjUlGxiOWU3IiI96qpt4lKMZdA8ATnFhW3rY5xEY+h5EZEvJcn1dvUiB9s9iQV58CpI1ZHI9ImKbkREe9VO1LKA4aB1wgIhoTB5r4uTYlYQsmNiHiv2jluPKjnBrTOlIjFlNyIiPfypNmJz1Wb3Ky3Ng6RNkrJjYh4p/Jic7QUeG7PzfGdUFpgbSwibZCSGxHxTicPmtuQDhDW0dpYzhceCx26AgYc3Wh1NCJtjpIbEfFOnlpvU0N1NyKWUXIjIt7JU+ttaii5EbGMkhsR8U6etBp4fWqSm6MbwVFlbSwibYySGxHxTic9vOcmpg8ERUD5abOwWETcRsmNiHgfwzh7WcpTa27sfpA43NzXpSkRt1JyIyLe53QOlBcBNnNdKU+VVLPO1Bpr4xBpY5TciIj3OXXY3EYmgn+QtbE0pssYc5u2RotoiriRkhsR8T551clNh66WhnFRnYebi2gWHdMimiJupORGRLxPTc+NJ1+SAggMhc5Dzf201dbGItKGKLkREe9T23Pj4ckNQJex5lbJjYjbKLkREe/jLT03AF3Gmdt0JTci7qLkRkS8jzf13CSNBJsd8g5BYZbV0Yi0CUpuRMS7lBXBmVxz3xt6boIjIW6Aua/eGxG3UHIjIt6lptcmJMpMHLxBsupuRNxJyY2IeBdvqrepoaJiEbdSciMi3sWb6m1q1CQ3ObvgTJ61sYi0AUpuRMS7eGPPTVg0RPcy97UUg4jLKbkREe/ijT03oEtTIm6k5EZEvIs39tzA2flulNyIuJySGxHxHpXlUHDU3Pe6npvqRTSztprD2UXEZZTciIj3KMgAwwH+IRAeZ3U0zROZCO27gFEFGeutjkbEpym5ERHvce5q4DabpaG0iOpuRNxCyY2IeA9vrbepoeRGxC2U3IiI9/DWkVI1aoqKMzdCRam1sYj4MCU3IuI9vL3nJqobtIuFqnLI3GR1NCI+S8mNiHgPb++5sdl0aUrEDZTciIh3MAw4dcTc99aeGzh7aUorhIu4jJIbEfEOp49DZQnY7BCZZHU0LVfTc5O+DqoqrY1FxEcpuRER73AqzdxGJIJ/oLWxtEanPhDcHiqKIXur1dGI+CQlNyLiHfKrk5v2ydbG0Vp2OyRXz1asuhsRl1ByIyLeoSa56dDF2jicQUXFIi6l5EZEvMMpH+m5gbqLaDoc1sYi4oOU3IiId8hPN7ftfaDnJn4gBIRBaT6c2G11NCI+R8mNiHgHX7os5RcASSPNfV2aEnE6JTci4vkcVVBw1Nz3hctScM6lqe+sjUPEBym5ERHPV3gMHJVgD4DweKujcY7aouI15gSFIuI0Sm5ExPPVXJKKTAS7n7WxOEvnYeAXCKezIe+Q1dGI+BQlNyLi+WqKiX2h3qZGQLCZ4IDqbkScTMmNiHi+2mHgPpTcgOa7EXERJTci4vlqh4H7SDFxjdrkRkXFIs6k5EZEPF/tMPCulobhdEmjzIVA89POjgYTkVazNLmZPXs2I0aMIDw8nJiYGG688Ub27t170cetXLmSYcOGERwcTLdu3Xj11VfdEK2IWMaXZic+V1A4xA8y99PWWBuLiA+xNLlZuXIlDz/8MGvXrmXp0qVUVlYyZcoUiouLG3zM4cOHueaaa5gwYQKpqak89dRTPPbYYyxcuNCNkYuI21SWQ9Exc9/Xam5A892IuIC/lS/+5Zdf1vl53rx5xMTEsGnTJi699NJ6H/Pqq6+SnJzM888/D0CfPn3YuHEjzz77LDfffLOrQxYRdys8CoYD/IOhXYzV0Thfl7Gw5iVIV8+NiLN4VM1NQUEBAFFRUQ22WbNmDVOmTKlz7KqrrmLjxo1UVFRc0L6srIzCwsI6NxHxIucWE9ts1sbiCsljzO2JPVCca20sIj7CY5IbwzCYMWMG48ePp3///g22y87OJjY2ts6x2NhYKisryc298Ith9uzZREZG1t6SkpKcHruIuJCvDgOvERoFnfqY++q9EXEKj0luHnnkEbZt28Z777130ba28/73ZlRPXX7+cYBZs2ZRUFBQe8vIyHBOwCLiHr46DPxcmu9GxKksrbmp8eijj/Lpp5+yatUqEhMTG20bFxdHdnZ2nWM5OTn4+/vTsWPHC9oHBQURFBTk1HhFxI18aTXwhnQZCxvnqqhYxEks7bkxDINHHnmERYsWsWzZMlJSUi76mDFjxrB06dI6x5YsWcLw4cMJCAhwVagiYhVfvywFZ3tusrdDaYG1sYj4AEuTm4cffph33nmHBQsWEB4eTnZ2NtnZ2ZSUlNS2mTVrFtOmTav9efr06aSlpTFjxgx2797Nm2++ydy5c5k5c6YVpyAirtYWLktFJECHFHNUWMZ6q6MR8XqWJjdz5syhoKCAyy67jPj4+NrbBx98UNsmKyuL9PT02p9TUlL4/PPPWbFiBYMHD+YPf/gDL774ooaBi/iiihJz1WzwvdmJz5c82twe3WBtHCI+wNKam5pC4MbMnz//gmMTJ05k8+bNLohIRDxKzZIEge0gpIO1sbha4nDY+p6SGxEn8JjRUiIiF/D1OW7OlTjC3B7dBA6HtbGIeDklNyLiuQqqp26IbAPzU8X0A/8QKCuAk/utjkbEqym5ERHPlV+T3DQ+RYRP8POHzkPNfV2aEmkVJTci4rlqam7at4GeGzDrbkDJjUgrKbkREc/Vli5LQd26GxFpMSU3IuK58ttYctO5uucmZyeUnbY2FhEvpuRGRDyTowoKM839tnJZKiIeIhLNyfyOpVodjYjXUnIjIp6pKAuMKrD7Q7tYq6NxH9XdiLSakhsR8Uw1xcQRncHuZ20s7lRbd7PR2jhEvJiSGxHxTDX1Nr68plR9apObDdCEWdxF5EJKbkTEMxVUz07cFua4OVf8QLAHQHHO2RmaRaRZLF1bSkSkQTWXpSweKWUYBscKSskuKKGwpJLyKgfhQf5EhASQ3DGUiOAA575gQAjEDYBjm83emw5dnPv8Im2AkhsR8Uy1l6Xcn9yUVVbx1a4cPtmSyeb0fHJPlzXYNrpdEP07RzCiaxRju3dkcFJ7bK1dBytxRHVysxEG/KB1zyXSBim5ERHPVOD+pRfOlFfyxjeHmffdYU6dqag97m+3Ed8+mMiQAAL87JwureTUmXJyT5eTe7qMFXtPsGLvCQASO4Rww+AE7hmXQnS7oJYFkjgC1r+mEVMiLaTkRkQ8j2Gcc1nK9QXFhmHw8ZZMnvliD8cLzV6auIhgbh7WmUm9Y+mXEEFwwIUjtopKKzh4opjU9FOsP5zHyn0nOHqqhJeXH+TNb49wz7iuPHhZd8Kbe+mqZjh49jaoLAP/FiZJIm2UkhsR8Twlp6C8eobeyM4ufalTxeU8tXg7X+zIBsyel/+5qhfXDUzAz9745aXw4AAGJ7VncFJ77hmXQkl5Fcv25PD6N4fYmpHPKysO8vn2LF758TD6JkQ0PagOXSE0Gs7kQtY2SBrRijMUaXs0WkpEPE9Nr01YJ7PA1kW2ZORz9Qur+GJHNv52G7+48hK+/sVEbhjc+aKJTX1CAv24dmA8Hz80ln9MG07n9iEcOXmGG1/5jsWpR5v+RDZb3SHhItIsSm5ExPO4YcHMxalHueW1NRwvLKNbpzAWPzSORyf3JMi/9RMG2mw2ruwby38eHc/lvTpRXulgxodb+ffWY01/kppLU5mazE+kuZTciIjnyXddMbFhGLzw1X4e/2Ar5ZUOrugTy6ePjGdAYqTTX6tDWCBz7xrBj0clYxjw+AdbWL43p2kPrum5yVDPjUhzKbkREc9T4JrZiR0Og9/9exd/+2ofAA9d1p3X7xxGuyDXlR/a7TZ+f0N/rh+UQKXD4MF3NnEgpwkrfnceCtjMyQyLjrssPhFfpORGRDyPCy5LORwGsxZtZ/7qIwD87/V9eeLq3thbUFvTXH52G8/dMoix3TtSWuHg5x+kUl7paPxBQeHQqbe5r0tTIs2i5EZEPI+TL0sZhsHTn+7gg40Z2G3w3C2DuHtcilOeu6kC/Oz87dbBtA8NYEdmIc8t3XfxB9WuEK7kRqQ5lNyIiOepGS3lpNmJZ3+xh3fWpmOzwV9vGcRNQ61Zryo2IphnbhoIwGurDrLxSF7jD1BRsUiLKLkREc9SUWouGglOuSz15reHeX3VIQBmf38A3x9i7UKcV/eP4wfDEjEM+N9/78ThaGTl7841yc1mcFS5J0ARH6DkRkQ8S2GmuQ0Ig5AOrXqqJTuz+cNnuwD45dW9uW2k62c7boonp/YmPMifHZmFfLS5kflvYvqY70P5aTix130Bing5JTci4lny081t+yRzMrsW2nmsgMfeT8Uw4PZRyUyf2M1JAbZedLsgHp3cA4C//Hcvp8sq629o96seNYUuTYk0g5IbEfEsThgplX+mnAfe3kRphYNLL+nE77/Xr/UrdTvZ3WNT6NoxlBNFZby64mDDDTsPM7eaqVikyVqU3Bw+fNjZcYiImGoXzGxZbUyVw+Cx97dw9FQJyVGh/P22Ifj7ed7/4wL97Tw5tQ+AuQp5cXn9DWuXYdjkpshEvF+LPvE9evTg8ssv55133qG0tNTZMYlIW1YzDLyFI6X+vmw/q/adIDjAzmt3DiMytJkrcrvRVf1i6RMfQXF5FW9+18B/GmtGTJ3YDWVF7gtOxIu1KLnZunUrQ4YM4Re/+AVxcXE88MADrF+/3tmxiUhb1IrLUhuO5PHi1/sBmH3TAPrEN2MlbgvYbDZ+Vl17M/+7IxScqbiwUXic+V4YDjiW6uYIRbxTi5Kb/v3789xzz5GZmcm8efPIzs5m/Pjx9OvXj+eee44TJ044O04RaStamNwUnKng5+9vwWHATUM7Wz7ku6mm9I2jd1w4RWWVzG2o96a27kZFxSJN0aoL0f7+/nz/+9/nww8/5M9//jMHDx5k5syZJCYmMm3aNLKyspwVp4i0BQ4HFFQPBW/mZalff7KDzPwSunQM5fc39HdBcK5ht9t4bHJPAOZ/d5ji+kZO1U7mp7obkaZoVXKzceNGHnroIeLj43nuueeYOXMmBw8eZNmyZWRmZnLDDTc4K04RaQtOHwdHBdj8oF1ckx/2+fYs/r31GH52Gy/cNsSlC2G6wtX94ujaMZTC0ko+2lTPvDe1RcUbwGhk0j8RAVqY3Dz33HMMGDCAsWPHcuzYMd566y3S0tL4v//7P1JSUhg3bhyvvfYamzdvdna8IuLLai5JRXQGv6YlKCdPl/Gbj3cA5irfg5Pauyg417Hbbdw73lzr6s3vDlN1/qzF8YPA7m8mfwWNTPonIkALk5s5c+Zw++23k56ezscff8x1112H3V73qZKTk5k7d65TghSRNuLcCfya6OlPdnKyuJzeceE8OqmniwJzvZuHJRIZEkDayTN8tft43TsDQiC2n7mv+W5ELqpFyc3SpUv55S9/SVxc3W5jwzBITze/nAIDA7nrrrtaH6GItB3NnOPmq13H+Wx7Fn52G8/+cBCB/p43n01ThQb68+NR5vIQc7+pp7C45tKU6m5ELqpF3wTdu3cnNzf3guN5eXmkpKS0OigRaaOaMVKquKyS3366E4D7J3Sjf+dIV0bmFneN7UqAn431R/LYmpFf986aRTQ1YkrkolqU3BgNFLSdPn2a4ODgVgUkIm1YMybwe/6rfWTml5DYIYSfTfbey1Hnio0I5vpBCQDM/fa83puaEVNZW6CqnvlwRKRWs4YUzJgxAzAnnnr66acJDQ2tva+qqop169YxePBgpwYoIm1IEy9L7cku5M3vjgDwhxv6ExLo5+LA3Ofe8Sks2pzJZ9uzeHJqbxLah5h3RHWH4PZQmg/Hd0DCECvDFPFozUpuUlPN2TENw2D79u0EBgbW3hcYGMigQYOYOXOmcyMUkbajNrlpuOfGMAx+/+9dVDkMpvaP4/LeMW4Kzj36JUQytntHVh88yT9XH2HWNeb6U9jt5mR+B782L00puRFpULOSm+XLlwNwzz338MILLxAR4dlTm4uIFyktgLICcz+ic4PNlu46zuqDJwn0t/NUzR9+H3PfhBRWHzzJgvXpPDq559l5exKHn01uRt5vbZAiHqxFNTfz5s1TYiMizlUzM3FIBwhqV2+Tssoq/vj5bgDun5BCUlRove283WWXxNCtUxhFpZV8uCHj7B21I6ZUVCzSmCb33Nx0003Mnz+fiIgIbrrppkbbLlq0qNWBiUgbU1id3DRSb/P2mjTSTp6hU3gQD17Ww02BuZ/dbuMn41L49cc7mL/6CHeN7Yqf3XZ2jamTB+BMHoRGWRuoiIdqcs9NZGQkNputdr+xm4hIs11kGHhRaQUvLz8AwMwpl3jdEgvNdfPQRNqHBpCed86kfqFRENXN3M/UDPAiDWnyt8O8efPq3RcRcYqLjJSa++1hTp2poFunMG4e6h0rfrdGSKAft49M5pUVB3nz28Nc1a960tTEEZB3yLw01fMKa4MU8VAtqrkpKSnhzJkztT+npaXx/PPPs2TJEqcFJiJtTCPJzanict6onrV3xpWX4O/nvTMRN8e0MV3xt9tYdziPHZnVxdaazE/kolr0DXHDDTfw1ltvAZCfn8/IkSP561//yg033MCcOXOcGqCItBGNJDevrjrI6bJK+sRHcE3/eDcHZp24yGCuHWie75s1k/olVtfdZG7UCuEiDWhRcrN582YmTJgAwEcffURcXBxpaWm89dZbvPjii04NUETaiAZqbk4Vl/P2mjTArLWx223ujsxSPxlnLmnz723HyCkshdgB4BcEJafMy1MicoEWJTdnzpwhPDwcgCVLlnDTTTdht9sZPXo0aWlpTg1QRNoARxUUHjP3z+u5mb/6CGfKq+iXEMEkH5uwrykGJbVneJcOVFQZvL02DfwDIX6QeadWCBepV4uSmx49evDxxx+TkZHBf//7X6ZMmQJATk6O5r8RkeY7fRwclWD3h3axZw+XVTJ/9REAHryse+2IzbbmJ+PN3pt316VTWlF1dp0p1d2I1KtFyc3TTz/NzJkz6dq1K6NGjWLMmDGA2YszZEjTpwRftWoV119/PQkJCdhsNj7++ONG269YsQKbzXbBbc+ePS05DRHxFDX1NuEJYD+7TtR769IpKKkgJTqMqW2o1uZ8U/rG0rl9CHnF5Xycmnk2udFkfiL1alFy84Mf/ID09HQ2btzIl19+WXt88uTJ/O1vf2vy8xQXFzNo0CBeeumlZr3+3r17ycrKqr317OkbKwKLtFm19TZnL0mVVVbxxrdmTckDl3YzJ7Fro/z97NwzritgDok3aibzy94OFSXWBSbioVo8C1ZcXBxxcXF1jo0cObJZzzF16lSmTp3a7NeOiYmhffv2zX6ciHioekZKfbYti+OFZcSEB/H9oQ2vNdVW3DIiib8t3cf+nNOsPB7CZWGdoPgEZG2D5FFWhyfiUVrUc1NcXMxvfvMbxo4dS48ePejWrVudm6sNGTKE+Ph4Jk+eXLuYZ0PKysooLCyscxMRD3NecmMYBm9+Zw59vmtsV4L8/Rp6ZJsRERzAbSOTAXjj2yNaZ0qkES3qubnvvvtYuXIld955J/Hx8W4r8ouPj+f1119n2LBhlJWV8fbbbzN58mRWrFjBpZdeWu9jZs+eze9+9zu3xCciLVRQd12pjWmn2JFZSJC/nR9V/0EXuGdcV+avPsK3B3I5PqE/sXwOGethzMNWhybiUVqU3HzxxRd89tlnjBs3ztnxNKpXr1706tWr9ucxY8aQkZHBs88+22ByM2vWLGbMmFH7c2FhIUlJ9a9dIyIWOW+Om5oJ624a2pmosECrovI4iR1CuWZAPP/eeox/HU/gEdBwcJF6tOiyVIcOHYiK8ozVaEePHs3+/fsbvD8oKIiIiIg6NxHxMOdcljp66gz/3ZkNwN1jUywMyjPdP8F8T147EIlh8zNXU695/0QEaGFy84c//IGnn366zvpSVklNTSU+vu0OERXxeuXFUJJn7kcm8u66dBwGjOvRkV5x4dbG5oEGJrZndLcoihxBZIdUjxTNWGdtUCIepkWXpf76179y8OBBYmNj6dq1KwEBAXXu37x5c5Oe5/Tp0xw4cKD258OHD7NlyxaioqJITk5m1qxZZGZm1q5j9fzzz9O1a1f69etHeXk577zzDgsXLmThwoUtOQ0R8QQ19TZBkZT7t+PDDeYlqjtHd7UuJg93/4RurD2Ux/Lirtxu2wMZG6D/zVaHJeIxWpTc3HjjjU558Y0bN3L55ZfX/lxTG3PXXXcxf/58srKySE9Pr72/vLycmTNnkpmZSUhICP369eOzzz7jmmuucUo8ImKBc+a4+XJnNieLy4mNCOKKPm1vqYWmurxXDN07hbH2ZA9uD0Q9NyLnsRlG21pWtrCwkMjISAoKClR/I+IJNv0T/v0Y9LyKW08/zrrDefxsck8ev/ISqyPzaO+vT+fvi5bxXfDPMOz+2J7MgMBQq8MScZnm/P1uUc0NQH5+Pm+88QazZs0iL8+8Xr5582YyMzNb+pQi0hZVF8MWBMaw7nAedhvcNlIjGi/mxiGdKQtLINvogM1RCcdSrQ5JxGO0KLnZtm0bl1xyCX/+85959tlnyc/PB2Dx4sXMmjXLmfGJiK+rTm42nAoDYHKfWOIjQ6yMyCsEB/gxbWwKmxxmUbGhS1MitVqU3MyYMYO7776b/fv3ExwcXHt86tSprFq1ymnBiUgbUF1z89Uxc2DCHaO7WBmNV7ljdBe22cy5v07t/dbiaEQ8R4uSmw0bNvDAAw9ccLxz585kZ2e3OigRaUOqe24OlnUgOSqUCT2iLQ7Ie0SFBdKh13gAAo5thLZVQinSoBYlN8HBwfWu0bR37146derU6qBEpI1wOMxJ6IBjRkduH5WMvQ2v/t0SV11xFWVGAOGOAtL2b7M6HBGP0KLk5oYbbuD3v/89FRUVANhsNtLT03nyySe5+WbNtSAiTXQmF6rKqTJs5PlF8cNhiRd/jNSREtuBjGDz0tS6VV9aHI2IZ2hRcvPss89y4sQJYmJiKCkpYeLEifTo0YPw8HD++Mc/OjtGEfFV1fU2x+nAlf2T6NguyOKAvFO7nmMBcKSv40RRmcXRiFivRZP4RURE8O2337J8+XI2bdqEw+Fg6NChXHHFFc6OT0R8WPGJNMKAY0Y0Px6l1b9bKrbfpbDjdQazj7fXHGHGlF4Xf5CID2t2cuNwOJg/fz6LFi3iyJEj2Gw2UlJSiIuLwzAMbDZdLxeRptmxcwejgKKgWC5L8YzFeL2RLWkUAJfYjnL/ml08eFkPQgL9LI5KxDrNuixlGAbf+973uO+++8jMzGTAgAH069ePtLQ07r77br7//e+7Kk4R8TGGYZBxZB8AnRJ76D9GrdEuBqNDCnabQdeyPXy0WauES9vWrJ6b+fPns2rVKr7++us6a0IBLFu2jBtvvJG33nqLadOmOTVIEfE93x04SVhpNvhBz569rQ7H69mSRsGpwwyz72PuN4f48UiNPJO2q1k9N++99x5PPfXUBYkNwKRJk3jyySd59913nRaciPiut9YcIcF2EoCgKNXbtFrSSABG+h/gyMkzrNx3wuKARKzTrORm27ZtXH311Q3eP3XqVLZu3drqoETEt2Xml/DV7uMk2HLNA5EaAt5q1XU3w+wHsOPgn2uOWBuPiIWaldzk5eURGxvb4P2xsbGcOnWq1UGJiG97b106AUY5nWzVk4EquWm9mD4QGE6Q4wy97EdZsfcEh3OLrY5KxBLNSm6qqqrw92+4TMfPz4/KyspWByUivqussor3N6QTX31JioAwCOlgbVC+wO4HicMA+FFcFgBvr0mzMiIRyzSroNgwDO6++26CguqfaKusTJNHiUjjvtyRTe7pcka1K4RKzF4bjZRyjqRRcGgFUyLSefrYKP61KYNfTLmEsKAWTWkm4rWa9Rt/1113XbSNRkqJSGNqehNuTDFgP7ok5UzVRcWxhVtJib6bw7nFLE7N1Err0uY0K7mZN2+eq+IQkTZg17FCNqadwt9uY0x0qZIbZ+s8HLBhyzvE/ZeG89SSYt5ac4Qfj0rWPELSprRobSkRkZZ4fdVBAK7qH0e7UrMuhMgkCyPyMSHtzcJi4Iboo4QG+rHv+GnWHsqzNi4RN1NyIyJukXaymE+3HgNg+qXdoaB6Fl313DhX9aWpsOOb+P6QzgD8c/URCwMScT8lNyLiFq+uPIjDgImXdGJAYqSSG1epnu+G9HVMG9MVgCW7ssnML7EuJhE3U3IjIi53LL+EjzaZycyjk3qAYSi5cZWa5OZYKr2iAxndLQqHAQvWaVi4tB1KbkTE5V5fdYiKKoNRKVEM7xoFZ/KgsronISLB2uB8TVQ3COsEVWVwbEtt782HG49SUeWwNjYRN1FyIyIulXu6jPc3pAPwyKQe5sGCDHPbLhb86583S1rIZoPk0eZ++hqu6BNLdLtAThSV8fXuHGtjE3ETJTci4lJzvz1MaYWDQYmRjO8RbR7UJSnXSh5jbtPXEuhv5wfDzBFp761PtzAoEfdRciMiLlNwpqJ20r5HJvU8O9eKkhvXSqruuclYCw4Ht40wk5tV+0+QkXfGwsBE3EPJjYi4zPzVRzhdVknvuHAm9445e0fNZSnNceMa8QPBPwRKTsHJ/XSNDmN8j2gMAz7YkGF1dCIup+RGRFzidFkl81YfBuChy3tgt58zQ656blzLLwASh5v76WsA+NHIZAA+3JihwmLxeUpuRMQl3l2bRv6ZClKiw7h2QHzdOwszza2SG9c5p+4G4Mq+sXQMCySnqIxle1RYLL5NyY2IOF1pRRX/+MbstXnwsu742c9b10g9N66XXDOZn9lzE+hv5wfDzfdbhcXi65TciIjTfbgxg9zTZXRuH1K7BECtyjIoyjb3VXPjOokjwWaHU0dq3+8fjTAvTa3cp8Ji8W1KbkTEqcorHby28hAAD0zsRoDfeV8zBUcBAwJCIbSj+wNsK4IjILafuV99aaprdBjjenTEMMwEVMRXKbkREaf6ODWTzPwSOoUHccvwenpm8qsvibRPNiecE9c5r+4GzhYWf7Ahg0oVFouPUnIjIk5T5TCYs/IgAPdPSCE4wO/CRucmN+JatYtorq49NKVvnAqLxecpuRERp/lsexaHc4tpHxrAj0d1qb+Rkhv36TLO3GZvh9ICoLqweJgKi8W3KbkREadwOAxeXnYAgJ+MSyEsyL/+hjXJjYqJXS8iHqK6g+GAtDW1h2+tnrF4xb4TZOaXWBWdiMsouRERp/hq93H2Hi+iXZA/d1WvRF2vmtmJ1XPjHl3Hm9sj39Qe6tapHWO6ddSMxeKzlNyISKsZhsHLy81emzvHdCEyNKDhxrWXpRq4bCXO1XWCuT3ybZ3DPxpVPWOxCovFBym5EZFW+/ZALluPFhAcYOfe8SkNN6wsh8Jj5r56btyja03dzTYoya89fFW/WDqEBpBdWMqKvSesiU3ERZTciEirvVRda/OjkclEtwtquGFh9Rw3/iEQFu2e4Nq6iISzdTfnDAkP8vdTYbH4LCU3ItIqG47kse5wHgF+Nn56abfGG2uOG2vUU3cDcFv1nDfL9+ZwTIXF4kOU3IhIq9TU2vxgWCLxkSGNN65NbjRSyq1SLjW359XddO/UjlEpUTg0Y7H4GCU3ItJie7OLWLH3BDYbPHBp94s/QHPcWKNL/XU3ALefU1hc5TDcHJiIayi5EZEWe+Mbcw2pq/vF0TU67OIPyNcwcEtExEPHHtV1N2vq3HVVvzjahwZwrKCUVftUWCy+QcmNiLRITmEpH2/JBLh4rU0N9dxYp7bupu6lqeAAP24eahYWL1BhsfgIJTci0iLzVx+hospgRNcODEnu0LQHaY4b69TMd3N41QV3/WikWQO1bE8O2QWl7oxKxCWU3IhIsxWXVfLO2jQA7p/QxF6bynIo0hw3lqnpucneDsUn69zVIyackV2jqHIY/EuFxeIDlNyISLN9sCGDwtJKUqLDuKJPbNMeVJhp1nz4B0NYJ9cGKBcKj4OYfoABh5ZfcPePRpm9N++rsFh8gJIbEWmWyioHc789DMB9E1Kw25s4X825C2ZqjhtrdL/c3B68MLmZ2j+eyJAAMvNL+Ga/CovFu1ma3KxatYrrr7+ehIQEbDYbH3/88UUfs3LlSoYNG0ZwcDDdunXj1VdfdX2gIlLrix3ZZOaX0DEssLYQtUlUTGy9HpPN7cGvwajbOxMc4MdNQzsDsGCdCovFu1ma3BQXFzNo0CBeeumlJrU/fPgw11xzDRMmTCA1NZWnnnqKxx57jIULF7o4UhEBc4HM11eZw7+njelKcIBf0x+s1cCtlzzWvCxYlAUn9lxw9+3VMxZ/tfs4mZqxWLyYv5UvPnXqVKZOndrk9q+++irJyck8//zzAPTp04eNGzfy7LPPcvPNN7soShGpsf5wHtszCwjyt3PnmGaOeFLPjfUCgs0J/Q5+DQe+gpg+de7uGRvO2O4dWX3wJO+sTeOXV/e2KFCR1vGqmps1a9YwZcqUOseuuuoqNm7cSEVFRb2PKSsro7CwsM5NRFrmrTXmCKmbhiYSFRbYvAcrufEMPau/Q/f9t9677x7bFTAX0yytqHJTUCLO5VXJTXZ2NrGxdUdmxMbGUllZSW5ubr2PmT17NpGRkbW3pCStaSPSElkFJXy5MxuAu8a2YJ4aJTee4ZLq5CZ9zQVLMQBM7hNLYocQ8s9U8En1JI0i3sarkhsA23mjLIzqorjzj9eYNWsWBQUFtbeMDM3hINISC9alU+UwGJUSRe+4iOY9uKrCHAoO5mgpsU5UN4juBY5K8/LUefzsNqZVX3Kc992R2u9YEW/iVclNXFwc2dnZdY7l5OTg7+9Px44d631MUFAQERERdW4i0jxllVW8Vz01/13Vly2apSDj7Bw37Zo4L464Tq+rze3eL+u9+9bhyYQE+LEnu4jvDpyst42IJ/Oq5GbMmDEsXbq0zrElS5YwfPhwAgICLIpKxPd9sT2b3NPlxEUEc2XfFiQneea8OHToCnav+trxTZdUJzf7l0BV5QV3R4YGcOsIs4ft1ZUH3RmZiFNY+i1z+vRptmzZwpYtWwBzqPeWLVtITzf/hzhr1iymTZtW23769OmkpaUxY8YMdu/ezZtvvsncuXOZOXOmFeGLtBn/XHMEgB+PSibArwVfG6fMx9Ohq7NCktZIHAmhHaE0H9K+rbfJveNT8LPb+PZALtuPFrg3PpFWsjS52bhxI0OGDGHIkCEAzJgxgyFDhvD0008DkJWVVZvoAKSkpPD555+zYsUKBg8ezB/+8AdefPFFDQMXcaFtR/NJTc8nwM/GbSNbWAx8qqbnJsV5gUnL+flD72vN/V2f1NskKSqU7w1KANR7I97H0nluLrvsskaL1ebPn3/BsYkTJ7J582YXRiUi56oZ/n3tgHg6hQe17Elqem6ilNx4jL43wOa3YPd/4JpnwX7hhIwPTOzG4tRMvtiRxZHcYrpGh1kQqEjz6eK3iDQor7icT7eaK3lPa0khce0THTG3uizlObpeCsGRUJwD6WvrbdI7LoLLe3XCYcDr3xxyc4AiLafkRkQa9P6GdMorHQzoHMmQpPYtexLDOKfmRj03HsM/EHpVX5rauajBZg9e1gOAjzYdJaeo1B2RibSakhsRqVdllYN315o1b3eO6dLgXFIXdeYklBcBNk3g52kGVNcr7lgIlWX1NhnRtQNDk9tTXulg3ndH3BebSCsouRGRetUsnhgVFlhbWNoiNb02EQnm2kbiObpdDuHxUHLKHBZeD5vNVtt7886aNApL61/qRsSTKLkRkXrNX30EgB+NTGre6t/nO3eOG/Esdj8YeIu5v+W9BptN7h1Dz5h2FJVVsmBdeoPtRDyFkhsRucDurELWHsrDz27jjtEtWEfqXBoG7tkG/cjc7v8vFNe/Rp/dbuOBid0BeOObw1pQUzyekhsRucA/q3ttru4XR3xkSOueTBP4ebaYPhA/2FxrasfCBpvdMDiBzu1DyD1dxgcbtEafeDYlNyJSx6nicj6uXg26RetIna/mspTmuPFcNb03WxY02CTAz86Dl5m9N6+uPEh5pcMdkYm0iJIbEanjg40ZlFY46BsfwYiuHVr/hBoG7vkG/ADs/pC1BXJ2N9jsB8MSiY0IIquglEWbj7ovPpFmUnIjIrUqqxy8XT0j8d1ju7Z8+HeNihIoMicB1GUpDxYWDT2nmPtb3m2wWXCAH/dP6AbAnJUHqaxS7414JiU3IlLrq905ZOaX0CE0gO8NbsXw7xr51SNrgiIgNKr1zyeuM+QOc5v6LlQ0PFnf7aOSiQoLJO3kGf6zLctNwYk0j5IbEalVU0h828jk1g3/rlE7DLwLtLYXSFyr51UQkQglebBzcYPNQgP9uXe8eYnxpeUHcDgaXh9QxCpKbkQEgD3Zhaw5dNI5w79rqN7Ge/j5w/B7zP0NbzTadNqYLkQE+3Mg5zT/3ZnthuBEmkfJjYgA8I9VZi/LVf1i6dy+lcO/a5zSBH5eZeg0sAdA5kY4ltpgs/DgAO4eZyasf192AMNQ7414FiU3IkJ2QSmfbjWHf9cUjDpFTc+NhoF7h3Yx0O9Gc/8ivTf3jO1KWKAfu7IKWb43x/WxiTSDkhsRYd53h6moMhiZEsWQZCcM/66hpRe8z4j7zO32j+BMXoPNOoQF1l6+fPFr9d6IZ1FyI9LGFZZW8G71ekHTJzqx18bhUM2NN0oaBbEDoLK00WHhAPdOSCHI386WjHxWHzzppgBFLk7JjUgbt2BdOqfLKukZ047LLolx3hMXZUFVGdj8IDLRec8rrmWzwcjq3pv1r0NVZYNNY8KD+dHIZABeWnbAHdGJNImSG5E2rKS8ije+MS8d/fTSbtjtThyunXfQ3HboCn4Bzntecb0Bt0BIlDlP0Z7/NNr0p5d2I8DPxppDJ9l4pOHLWCLupORGpA17b306uafLSOwQwo1DOjv3yU9W/0++Yw/nPq+4XmDo2dqbNS832jShfQg3DzV75uasOOjqyESaRMmNSBtVWlHFa6vMP0YPXtadAD8nfx2crP5Dp+TGO424D/wC4eh6yFjfaNOfXtoNmw2+3pPDvuNFbgpQpGFKbkTaqH9tOsrxwjLiI4P5wTAX1MTU9tx0d/5zi+uFx8LAW8z9NS812rRbp3Zc1TcOgNdXHXJ1ZCIXpeRGpA0qq6xiznIz+Zg+sTtB/k5YauF8uizl/UY/bG53//vsyLcGPFA90u6TLZlkFZS4ODCRxim5EWmD3l+fwbGCUmIjgrh1RJLzX6Cq4uwfQ/XceK/YvtB9EhgOWPtqo02HJHdgVEoUFVUGb3572E0BitRPyY1IG1NSXsVL1b02j07q6ZwFMs+Xnw6OSvAPgXAnrC4u1hnziLlNfRtK8httOv0yM5FdsC6dgjMVLg5MpGFKbkTamLfWHOFEkTlC6pbhLui1gbr1NnZ9zXi17pMgpi+Un4bN/2y06WWXdKJXbDjF5VW8sy7NTQGKXEjfOiJtSFFpBa+uNEcx/WxyTwL9XfQVoGJi32GzwZjq2pt1r5mXHBtsaqutvZn33RFKK6rcEaHIBZTciLQhr608xKkzFXTrFMb3nT2vzbk0DNy3DPghhMVAYSbsWNho0+sHJZAQGUzu6TIWbc50U4AidSm5EWkjsgtKeeNbc5juL6/ujb+z57U5l0ZK+Rb/IBg93dz/9m/mumENCPCzc2/1yvL/+OYQVQ4tqCnup+RGpI14/qt9lFY4GN6lA1P6xrr2xdRz43tG3AdBkXBiD+z9rNGmt41IIjIkgMO5xSzZme2mAEXOUnIj0gbsO17EhxszAJh1TW9sNieuIXW+8jNQeNTcV3LjO4IjYeT95v6qZ8FouEcmLMifaWO6APDqyoMYjbQVcQUlNyJtwJ+/2IPDgKv7xTGsS5RrXyyveobakA4Q6uLXEvca/aA5vD9rCxxc1mjTu8Z2JcjfztajBaw9pAU1xb2U3Ij4uLWHTvL1nhz87DaeuLqX619Q9Ta+Kywaht1t7n/z10abRrcL4ofDzWU9atYwE3EXJTciPswwDGZ/vhuAH41Molundq5/USU3vm3so2APgLTvIG1No03vn9ANuw1W7D3B7qxCNwUoouRGxKd9tj2LrUcLCA3042eTL3HPi9YWE2uOG58U2RkG327ur/hTo027dAxj6oB4QAtqinspuRHxUeWVDv7y370A/PTSbnQKD3LPC6vnxvddOtPsvTm8Cg5/02jTBy41h4V/uvUYR0+dcUd0IkpuRHzVgnVppJ08Q3S7IO6vnnfELZTc+L72yTB0mrm//E+NjpwamNiesd07UuUwmKsFNcVNlNyI+KDC0gpeXGYmGY9f2ZOwIH/3vPCZPCipHhkT5caEStxvwi/ALwjSV8OhFY02nT7RvET5/voMThWXuyE4aeuU3Ij4oNdWHiSvuJxuncK41VWLY9Ynd7+5jegMgWHue11xv8jOMPwec/8ivTcTekbTNz6Ckooq3l6rBTXF9ZTciPiY7ILS2u5/ly+zcL4Te8xtJzcMORfrjZ9hzntzdD0c+KrBZucuqDl/tRbUFNdTciPiY/621I3LLJyvNrnp497XFWuEx8KIe8395X9stPfm2gHxdG4fQl5xOf/adNRNAUpbpeRGxIfsO17Evza5aZmF+uSYc+qo56YNGf84BITBsVTY92WDzfz97Nw/IQWAf6w6RGVVw4tvirSWkhsRH+LWZRbqc8Icek6Mem7ajLBoGPWAub/8j42uGH7LiCQ6hAaQnneGL7WgpriQkhsRH7HmoJuXWThfST4UHTP31XPTtox9FIIiIHs77P6kwWahgf5MG9MVgFeWa0FNcR0lNyI+wDAMnvnCzcssnK+m1yY8wVxBWtqO0CgY84i5v+yPUFXZYNO7x3YlLNCPXVmF/Fe9N+IiSm5EfIAlyyycr6aYOKa3Na8v1hrzEIR2hJP7Ydv7DTbrEBbIT8abtTd/W7ofh0O9N+J8Sm5EvFx5pYP/96XZa/LApd3dt8zC+TRSqm0LCjeHhgOseAYqyxpset/4boQH+7P3eBH/2Z7lpgClLVFyI+Ll3l2XRnqeuczCfdWjUSyhkVIy4l7zsmRBBmya32CzyNAA7htvznvz/Ff7NHJKnE7JjYgXKyyt4O9WLLNQH42UkoAQmPiEub/qL1Be3GDTn4zvSvvQAA6dKObTrcfcFKC0FUpuRLyYZcssnE8jpaTGkDugQwoUn4B1rzbYLDw4gJ9Wrxj+wtf7qVDvjTiR5cnNK6+8QkpKCsHBwQwbNoxvvvmmwbYrVqzAZrNdcNuzZ48bIxbxDJYus3A+jZSSGn4BcPlT5v53L5iJbwPuGtOVjmGBpJ08w6LNmrVYnMfS5OaDDz7g5z//Ob/61a9ITU1lwoQJTJ06lfT09EYft3fvXrKysmpvPXv2dFPEIp7D0mUWzqeRUnKu/jdDTF8oLYDVf2+wWViQPw9eZq4Y/uLXByivVO+NOIelyc1zzz3Hvffey3333UefPn14/vnnSUpKYs6cOY0+LiYmhri4uNqbn5+fmyIW8Qx7s89dZqGP+5dZOJ9GSsm57H4w6dfm/to5cDqnwaY/HtWFTuFBZOaXsGCdVgwX57AsuSkvL2fTpk1MmTKlzvEpU6awevXqRh87ZMgQ4uPjmTx5MsuXL2+0bVlZGYWFhXVuIt7uz1+eu8xCB6vD0UgpuVCva6DzMKgohm+ea7BZSKAfP5ts9r6/8PV+Ckoq3BWh+DDLkpvc3FyqqqqIja3bnR4bG0t2dv2zVsbHx/P666+zcOFCFi1aRK9evZg8eTKrVq1q8HVmz55NZGRk7S0pycKiSxEnWHPwJMusXGahPhopJeez2WDy0+b+xrmQn9Fg09tGJNEjph2nzlTwyvIDbgpQfJnlBcXnd6cbhtFgF3uvXr24//77GTp0KGPGjOGVV17h2muv5dlnn23w+WfNmkVBQUHtLSOj4Q+YiKdzOAxmVy+zcPvIZGuWWTifRkpJQ7pdBimXQlU5rPxzg838/ew8dY1ZrzXvuyNk5J1xU4DiqyxLbqKjo/Hz87uglyYnJ+eC3pzGjB49mv379zd4f1BQEBEREXVuIt7qs+1ZbDtaQFigH49N9pBCeo2UksZMqu692bIAchv+rr68Vwxju3ekvMpRm8CLtJRlyU1gYCDDhg1j6dKldY4vXbqUsWPHNvl5UlNTiY+Pd3Z4Ih6nvNLBX/5rJhI/tXKZhfMd325uY/tZG4d4pqQRcMlUMKpg+Z8abGaz2fj1tX2x2+Dz7dl8dyDXjUGKr7H0stSMGTN44403ePPNN9m9ezePP/446enpTJ8+HTAvKU2bNq22/fPPP8/HH3/M/v372blzJ7NmzWLhwoU88sgjVp2CiNvULLPQKdziZRbOl12d3MQNsDYO8Vw1I6d2LoKsbQ0265sQwZ2juwDwv5/u1MR+0mIWztUOt956KydPnuT3v/89WVlZ9O/fn88//5wuXcxf7qysrDpz3pSXlzNz5kwyMzMJCQmhX79+fPbZZ1xzzTVWnYKIWxSWVvDi12aX/s+vsHiZhfPVJDfxA62NQzxXXH/o/wPY8REs+z/48YcNNp1xZS/+vS2L/Tmn+efqI9w3oZsbAxVfYTMMo02tN19YWEhkZCQFBQWqvxGv8afPd/P6qkN06xTGkp9fau1sxOeqqoTZnaGyFB7dDB27Wx2ReKqTB+GlEeblqbv+AykTGmz6wYZ0frlwO2GBfiyZMZHO7UPcGKh4qub8/faQb0gRacjh3GLmfWcus/Cba/t6TmIDkHfQTGwCwsz1hEQa0rE7DL/H3P/8f6Cq4flsfjgsieFdOlBcXsWvF2+njf0fXJzAg74lRaQ+f/xsFxVVBhMv6cTlvWOsDqeu2nqb/mDX14lcxOW/gtCOcGJ3o4tq2u02nrl5IIF+dpbvPaFVw6XZ9G0k4sFW7TvBV7tz8Lfb+M11HjhBXtZWc6tiYmmK0Ci44nfm/vI/Qd6hBpv2iGnHo5N6AGZx8cnTZe6IUHyEkhsRD1VZ5eAP/9kFwLQxXekRE25xRPXQSClprsE/hq4ToOIMfPIoOBoeEfXAxO70jgvn1JkKfl/9WRBpCiU3Ih7q3XXp7M85TYfQgNq1dzyKYUDWFnM/TiOlpInsdvje3yEgFNK+hdUvNtg00N/On28eiN0Gn2w5xrI9x90YqHgzJTciHuhUcTnPLd0HwC+m9CIyNMDiiOpx6jCUnAK/QIjtb3U04k2iUuCq6gn9vv49pK9rsOmgpPb8ZJxZrP7rxTsoKtXCmnJxSm5EPNDfvtpHQUkFvePC+dHIZKvDqV/mZnMbNwD8A62NRbzPsLvNuW+MKvhwGhQcbbDpjCmXkBQVwrGCUv7vP1qaQS5OyY2Ih0lNP8Xba9MAePr6vvjZ619I1nLHUs1twlBr4xDvZLPB9c9DTF84nQ3v3gKlBfU2DQ305y8/GITNBh9szOC/O7PrbSdSQ8mNiAcpr3Tw5MLtGAbcNLQzY7tHWx1Sw2p6bjoruZEWCgqH2z+AsBjI2Qnv/ABKC+ttOrpbR35aPVvxrEXbySkqdWek4mWU3Ih4kNdXHWTv8SKiwgL59bV9rQ6nYVWVZ4uJ1XMjrdE+Ge74CILbw9H18PaNUFR/4fCMKZfQJz6CvOJynvhomyb3kwYpuRHxEHuyC3nx6wMAPH1dX6LCPLiOJXevOZQ3sB1Ee+BILvEu8YNg2idmgpO5Cf4xCbJ3XNAsyN+PF24bTKC/nRV7T/DOuvQLn0sEJTciHqGssoqfv7+F8ioHk3rHcMPgBKtDatzRjeY2fjDY/SwNRXxEwmC472vo2AMKj8LcKbDn8wuaXRIbzpNX9wbM2bsPnjjt5kDFGyi5EfEAf12yjz3ZRXQMC+TPNw/EZvPQIuIaGdVDd5NHWRuH+JboHnDfV5ByKVQUw/s/gq//AI6qOs3uHtuV8T2iKa1w8LP3UymrrGrgCaWtUnIjYrFv9p/gH9+Y09A/c/NAOoUHWRxRE6SvMbfJY6yNQ3xPSAe4YxGMuN/8+ZtnzTqc0ydqm9jtNp794SDahwawI7OQ2Z/vsSZW8VhKbkQslJlfwmPvpWIY8KORyVzZN9bqkC7udE71mkA2SBxhdTTii/wC4Npn4ea55orzh1fBaxMgbU1tk7jIYJ67ZRAA81cf4YvtWVZFKx5IyY2IRcoqq3jo3c2cOlNB/84R/PZ6Dx4dda70teY2pi+EtLc0FPFxA34AP10O0b2gKAvmXwvfvWgu/QFM6h3LAxPN4eFPfLSNtJPFVkYrHkTJjYgFDMPgt5/sZGtGPpEhAcz58TCCA7ykMFf1NuJOnXrB/ctgwA/N2YyX/gY+uAPKigCYOaUXw7t0oKiskocXbFb9jQBKbkQs8fqqQ7y/IQO7DZ6/bTBJUaFWh9R0qrcRdwtqBzf9A679q7mW2Z7/wJtXQ346AX52/n77EDpU19/88TMtzyBKbkTc7ssdWTzzpVkA+Zvr+nJ5rxiLI2qG0gI4tsXcV3Ij7mSzwYj74J4vzBmNj+8w58PJWE98ZAjP3ToYgLfWpPHJlkxrYxXLKbkRcaM1B0/y2PtbMAyYNqYLd4/tanVIzXPkO/PSQFR3aJ9kdTTSFiUON+twYgdA8QmzDmfrB1zeK4aHLusOmPU324/Wv06VtA1KbkTcZNvRfO5/ayPllQ6u6BPD09f19fz5bM53aIW57XaZlVFIWxeZCD/5EnpfB1XlsPin8NXv+MWVPbm8VyfKKh3c/9ZGrT/Vhim5EXGDAzlF3PXmek6XVTK6WxQv3T4Ufz8v/PgpuRFPEdQObnkbxs8wf/72Ofz+NY0Xbr6EHjHtyC4s5b5/bqS4rNLaOMUSXvjtKuJdjp46w51z13PqTAUDEyN5464R3jMy6lyFWeaaUtggZYLV0YiA3Q5X/Ba+/1ptoXHEguuY9/14osIC2Xa0gAff3Ux5pcPqSMXNlNyIuFBOUSl3zl1PVkEpPWLaMf+ekbQL8rc6rJap6bVJGGzOIiviKQbdBnf9G0KjIXs7SR9O4aOJuYQE+LFq3wme+GgrVQ6tIN6WKLkRcZGcolJu/8c6DucW07l9CG/fO9KzV/q+mH1fmNseV1gbh0h9kkeb8+HED4KSU3RbNp1lPT4g0l7Cx1uOMfNfSnDaEiU3Ii6QU1TKj15fy4Gc08RHBrPg/lHER4ZYHVbLVZbBga/N/V7XWBuLSEM6dIF7v6quw7ERf3gRayJ/w+V+W1mcmsnjH2yhokqXqNoCJTciTpZTaCY2B08UEx8ZzPs/HU2XjmFWh9U6h7+B8tMQHg/xg62ORqRh/oFmHc49n0P7ZEJLjjEv4M88H/AK32zdwz3zNlBYWmF1lOJiSm5EnCinsJQf/cNMbBJ8JbEB2Pu5ub3karOIU8TTdRkLD66B0Q8BNm70+5avg/6HqEOf8INXviMj74zVEYoL6VtKxEmOnjrDba+fTWze85XExuGAvdX1NrokJd4kqB1cPRvu+wpi+hJlK+LFwJd58tRvuffFxfx3Z7bVEYqLKLkRcYK92UXcPGc1h6qLh30msQFIXw1FxyAoElIutToakeZLHA4/XQmX/wrDL5BJfltYZMzguwV/4unF2zituXB8jpIbkVbaeCSPH766muOFZVwS246PHhzjO4kNwLYPzG3f70FAsLWxiLSUfyBMfALb9G9xJI2ina2U3wf8k++l3stP//ouy/fmWB2hOJGSG5FWWLbnOHfMXUdhaSVDk9vz4QNjvHtU1Pkqy2DXJ+b+wFutjUXEGTr1wn7Pl3DNs1T6hzHcvo95ZTNY+9bTzHh/E3nF5VZHKE6g5EakBQzD4J21adz/1iZKKxxc3qsT7943mvahXjyPTX32/ddcCTyiM3QZZ3U0Is5ht8PI+/F/dD1V3a8kyFbJrID3uGPXT/nJX9/j3XVpVGrIuFdTciPSTGWVVcxatJ1ff7yDKofBTUM78/q04YQEeuGSChezaZ65HfBDjZIS3xOZiN8d/4IbXqYqoB1D7Qd4v+oXZH/6e258YRnf7D9hdYTSQjbDMNrUlI2FhYVERkZSUFBARESE1eGIlzleWMr0dzaRmp6PzQZPXNWb6RO7ed/q3k2Rux9eGg7Y4GdboENXiwMScaH8DByfPor90HIAjhrRvFR5I3ndbuSJ6wfTIybc4gClOX+/ldyINNGmtFNMf2cTJ4rKiAj25++3D2XiJZ2sDst1vvglrHsVLpkKt79vdTQirmcYsGMhjv/+GvvpLADyjTC+cIyktMskrrp8MgkpvcHug720XkDJTSOU3EhzVTkMXl15kL8t3Uelw+CS2Ha8fudwukb70Iio85Xkw/MDoKwQ7lgEPSZbHZGI+1SUwMZ5VKx+mYCio3XvsgVCVDcCortB+2SITDK30ZeYN12+dZnm/P320uWJRdzjSG4xT3y0jfVH8gC4bmA8z9w80HtX9m6qta+YiU2nPtDtcqujEXGvgBAY8xABox6AwyvJ3fQxxfu/IbY8g2DK4eQe83a+oAjoPAySRpoLeSaONCcSFLfz8W9okZapqHLwj28O8cJX+ymrdBAW6MfvbujPzUM7+2Z9zblKTsHaOeb+ZU/qf6LSdtn9oPskortPIhpYdyCHt778huJje+hsyyXF/yTjoovpGZiHf+4e8z8Eh5abNwCbH8QPhOSx0HkoxPaDjj3BT396XU3vsMg5DMNg6a7jPPPFHg7lFgMwoWc0f/r+AJKiQi2Ozk2+e8H8ko7pB32+Z3U0Ih5jVI8YRj58Eyv3neDZJXt5N7MQMiEs0I+7xiTyQK8yIk+mQsZ6SFsDBelwLNW81fALgk69IG6AmezE9DVv7WLA1//j5EaquRHBTGpW7DvBS8sOsCntFAAdwwKZdU2fttFbUyN3P7wyBhwVcNt70FtrSYnUxzAMluw6zgtf7WdXViEAoYF+3DmmCz+d0I2O7YKg4KiZ5GSsheztcHwnlJ+u/wlDOpiXgWN6192GRSvpqaaC4kYouZFzFZZW8ElqJgvWZ7C7+gsq0N/OfeNTePCy7oQHB1gcoRsZBrz9fbNLvecUuP1DfamKXIRhGHy1O4cXvt7HjkzzOyQkwExy7p/QjU7hQWcbOxyQnwbHd5iJTvZ2yNkNpw6D0cCkgX5BEB4L7eIgvPrWLtbcRiSYCVB4XJv4rCq5aYSSGzEMg83p+by3Pp3/bDtGaYX5pRIS4Mcdo5O5f0I3YiLa4BpK6/8Bn880v0wfXgtR3ayOSMRrGIbBsj05vPD1frYdLQAgOMDOj0d14aeXdiO2se+UihKz1/TEHsjZBTl74MRuOJUGNOFPdHB78zJX90lwydUQ08cnkx0lN41QctM2GYbB3uNF/GdrFv/edoy0k2dq7+sZ047bRiZz89DOvrd8QlNl74B/TIKqMrj6GRj9oNURiXglwzBYsfcEz3+9n60Z+QAE+Nm4bmAC945PoX/nyKY/WUUpnD5u3oqyoCjbvNUcO5UGeQcv7PWJTIYhd8Cwu81eHx+h5KYRSm7aloMnTvOfrVn8Z9sx9uecvdYdHGDn2gEJ3D4qiaHJHdpOTU19CrNg7pVQkAE9r4LbP/DJ//WJuJNhGKzan8tLy/az4cip2uMjU6L4ybiuTOodS6C/E0YiVpTCyf2QsQ72LYHDK6Gy1LzP7g99b4BxPzdHbXk5JTeNUHLj28orHWxKO8U3+0+wfO+J2joagEA/OxN7deK6gfFc0SeWMF+fq6YpzuTBW98zr/1HdYd7l0JYR6ujEvEp247mM/fbw3y2LYtKh/knt0NoANcNTODGIZ0Zmtzeef/BKj8Dez6DDf8wE54afa6Hy2aZI7S8lJKbRii58R1VDoOcolL2HT/NrmOFbDiSx9pDJzlTXlXbxt9uY3zPaK4bmMCVfWOJDGlDBcIXU5RtFhDn7ILQaLhvqepsRFwou6CUf645wkebjnKiqKz2eJeOoUztH8/lvToxtEsHAvycNLdU1lZzaocdi6it3el7ozl/VUwf57yGGym5aYSSG89S5TA4ll9C2skzHMsvoaCkgsLSCnNbUsGZ8irKKh2UVVZvK8z9kvIqcorKav8XdK7odoFM6NmJ8T2imdQ7hg5hbbSOpjHpa+Ffd5vX8dvFwZ2LIbav1VGJtAlVDoPvDuTycWomX+7MrvMfsvBgfy7t2YmJl3RiVLcokqNCW9+rk7MbVv4Zdi6uPmCD/jebPTnRPVr33G6k5KYRSm6s43AYHD5ZzPajBWw7WsD2zHx2ZBZSUlF18Qc3wM9uo0tUKH0SIhjYOZLxPaPpExeB3a6akXpVlJhfct+9CEYVRPcya2yiUqyOTKRNOlNeyVe7c1i+J4eV+06QV1xe5/7YiCBGdI1iVEoUQ5I7cElseMtrdY7vghWzYfen5s82Owy8DSY+4RXfAV6V3Lzyyiv85S9/ISsri379+vH8888zYcKEBtuvXLmSGTNmsHPnThISEnjiiSeYPn16k19PyY17VDkM0vPOsPNYAduPFrD1qJnInC6rvKBtoJ+d5I6hdG4fQvvQACJDzFtEcAChQX4E+fsR5G8nyN9OcIBf7TYmIoiY8GD8lMhcXEUpbHsfVv4/KMw0j/X/AVz/gta+EfEQVQ6DbUfzWbH3BN8eyGXb0Xwqqur+iQ7ws9ErLpz+CZH0S4igW6d2pESHERcR3PT/1GVtg+V/gn1fmD/b/c2enJEPQOIwJ5+V83hNcvPBBx9w55138sorrzBu3Dhee+013njjDXbt2kVycvIF7Q8fPkz//v25//77eeCBB/juu+946KGHeO+997j55pub9JpKbpzDMAwKSys5XljK8cJSsgtKycwv4UDOaQ7knOZQbjHllRdOShXkb6dfQgQDE9szMDGSgYmRpES3U4LiClWV5rTvuz6Gre/DmVzzeERnmPr/oM91loYnIo0rrahiS0Y+6w/nseFIHtuOFlBQUlFv2+AAO107hpHYIZTYiCDiIoKJjQgmNjK4ej+IyJCAupe4jm6C5X+Eg1+fPZYwFAbfDr2mQmSii8+webwmuRk1ahRDhw5lzpw5tcf69OnDjTfeyOzZsy9o/8tf/pJPP/2U3bt31x6bPn06W7duZc2aNU16TW9LbgzDoMphUFFlUF7loLLKQUWVQUWVg/IqBxVVDiqr76uorHtf5XntyirM/fLqGhZz67hge+F9Vec87uy2qp56l3MF+dvpFRduJjGd2zMgMZKeMe3wd1axnJgJTHmRudhlwVHIT4e8Q+baNpmboaL4bNuIRBjzMAz/CQS0wUkKRbycYRgcPVXCjswCdhwrYHdWEUdyi0nPO1Nv/eH5AvxsRIUFEhUWRHS7QDpW7/dx7Gd4zkd0yfoCu+Ns8uSI6Y+t6zhssX3NteY6dIGQKMsW/mzO32/LxsKWl5ezadMmnnzyyTrHp0yZwurVq+t9zJo1a5gyZUqdY1dddRVz586loqKCgIALR8KUlZVRVna2Kr2gwJw5srCw8IK2rXGquJzH3kvFYRgYgMMwfxEdhlG7bxhU/2zum+2q2zjMNpUOMyGpdBiUVxlUVjnw5KqoiGB/YiOC6RQRRGx4ECnRYXSPaUe36HYktA+5oEfmTHED66q0dR/cAWWnwVFlTshlVG8djuqfq284zDYVpVBWBFWljT9vYAR0vxz6fg+6Tza/lErKzZuIeJ1IfxjXJYxxXcJqj1VUOWoHZmQVlnCioIycojKOF5VyoqiMnMJS8ksqKQOyzhSTdaK+Z/4eUUzke36rudxvC4NtB/HL2A4Z2y9oWUQYJYRSYfOj0uZPFQE4bH44MP/j+lrEY5RFJPPS7UOdeu41f7eb0idjWXKTm5tLVVUVsbF1Z0+MjY0lOzu73sdkZ2fX276yspLc3Fzi4+MveMzs2bP53e9+d8HxpKSkVkQv59ppdQDSiCLgneqbiEjDMoCtwB8u2rLoIvffBcDbLprovKioiMjIxmd6tnwWs/OHuBmG0eiwt/ra13e8xqxZs5gxY0btzw6Hg7y8PDp27Fj7mMLCQpKSksjIyPCKS1XOpvPX+ev82+75g94Dnb93nL9hGBQVFZGQkHDRtpYlN9HR0fj5+V3QS5OTk3NB70yNuLi4etv7+/vTsWP9s6oGBQURFBRU51j79u3rbRsREeHR/7CupvPX+ev82+75g94Dnb/nn//FemxqWFbZGRgYyLBhw1i6dGmd40uXLmXs2LH1PmbMmDEXtF+yZAnDhw+vt95GRERE2h5Lh63MmDGDN954gzfffJPdu3fz+OOPk56eXjtvzaxZs5g2bVpt++nTp5OWlsaMGTPYvXs3b775JnPnzmXmzJlWnYKIiIh4GEtrbm699VZOnjzJ73//e7Kysujfvz+ff/45Xbp0ASArK4v09PTa9ikpKXz++ec8/vjjvPzyyyQkJPDiiy82eY6bhgQFBfHb3/72gstXbYXOX+ev82+75w96D3T+vnf+ls9QLCIiIuJMmk1NREREfIqSGxEREfEpSm5ERETEpyi5EREREZ/ic8nNqlWruP7660lISMBms/Hxxx832v7uu+/GZrNdcOvXr19tm/nz59fbprT0Iuv6uNns2bMZMWIE4eHhxMTEcOONN7J3796LPm7lypUMGzaM4OBgunXrxquvvnpBm4ULF9K3b1+CgoLo27cvixcvdsUptFpL3oNFixZx5ZVX0qlTJyIiIhgzZgz//e9/67Tx5d+BFStW1Htue/bsqdPOG34HWnL+vvQdMGfOHAYOHFg7GduYMWP44osvGn2ML33+m3v+vvTZr9Hc98CXPv/n8rnkpri4mEGDBvHSSy81qf0LL7xAVlZW7S0jI4OoqCh++MMf1mkXERFRp11WVhbBwZ61svLKlSt5+OGHWbt2LUuXLqWyspIpU6ZQXFzc4GMOHz7MNddcw4QJE0hNTeWpp57iscceY+HChbVt1qxZw6233sqdd97J1q1bufPOO7nllltYt26dO06rWVryHqxatYorr7ySzz//nE2bNnH55Zdz/fXXk5qaWqedr/4O1Ni7d2+dc+vZs2ftfd7yO9CS8/el74DExESeeeYZNm7cyMaNG5k0aRI33HADO3fWvwKcr33+m3v+vvTZr9Hc96CGL3z+6zB8GGAsXry4WY9ZvHixYbPZjCNHjtQemzdvnhEZGenc4NwgJyfHAIyVK1c22OaJJ54wevfuXefYAw88YIwePbr251tuucW4+uqr67S56qqrjNtuu825AbtAU96D+vTt29f43e9+V/uzL/8OLF++3ACMU6dONdjGW38HWvLv70vfAYZhGB06dDDeeOONeu/z9c+/YTR+/vXxlc/+uRp7D3z18+9zPTetNXfuXK644oraiQRrnD59mi5dupCYmMh11113QWbviQoKCgCIiopqsM2aNWuYMmVKnWNXXXUVGzdupKKiotE2q1evdnLEzteU9+B8DoeDoqKiCx7jq78DNYYMGUJ8fDyTJ09m+fLlde7z1t+Blvz7+8p3QFVVFe+//z7FxcWMGTOm3ja+/Plvyvmfz5c++9C898DXPv9Kbs6RlZXFF198wX333VfneO/evZk/fz6ffvop7733HsHBwYwbN479+/dbFOnFGYbBjBkzGD9+PP3792+wXXZ29gULlcbGxlJZWUlubm6jbc5fxNTTNPU9ON9f//pXiouLueWWW2qP+fLvQHx8PK+//joLFy5k0aJF9OrVi8mTJ7Nq1araNt74O9CSf39f+A7Yvn077dq1IygoiOnTp7N48WL69u1bb1tf/Pw35/zP5yuf/ea8B776+ddlqXP86U9/Mjp27GiUlZU12q6qqsoYNGiQ8eijj7YyQtd56KGHjC5duhgZGRmNtuvZs6fxpz/9qc6xb7/91gCMrKwswzAMIyAgwFiwYEGdNu+8844RFBTk3KCdrKnvwbkWLFhghIaGGkuXLm20nS/9DtTnuuuuM66//vran73xd6Al5+8L3wFlZWXG/v37jQ0bNhhPPvmkER0dbezcubPetr74+W/O+Z/Llz77LX0PavjC5189N9UMw+DNN9/kzjvvJDAwsNG2drudESNGeGzm/uijj/Lpp5+yfPlyEhMTG20bFxd3Qfadk5ODv78/HTt2bLTN+Zm8J2nOe1Djgw8+4N577+XDDz/kiiuuaLStL/0O1Gf06NF1zs3bfgdacv6+8h0QGBhIjx49GD58OLNnz2bQoEG88MIL9bb1xc9/c86/hi999qFl78G5vP3zD7osVWvlypUcOHCAe++996JtDcNgy5YtxMfHuyGypjMMg0ceeYRFixaxbNkyUlJSLvqYMWPGsHTp0jrHlixZwvDhwwkICGi0zdixY50XvJO05D0AeO+997j77rtZsGAB1157bZNex1d+B+qTmppa59y85XegNefvC98B9TEMg7Kysnrv87XPf30aO3/wnc9+Yy72HpzPWz//dVjSX+RCRUVFRmpqqpGammoAxnPPPWekpqYaaWlphmEYxpNPPmnceeedFzzujjvuMEaNGlXvc/7v//6v8eWXXxoHDx40UlNTjXvuucfw9/c31q1b59Jzaa4HH3zQiIyMNFasWGFkZWXV3s6cOVPb5vzzP3TokBEaGmo8/vjjxq5du4y5c+caAQEBxkcffVTb5rvvvjP8/PyMZ555xti9e7fxzDPPGP7+/sbatWvden5N0ZL3YMGCBYa/v7/x8ssv13lMfn5+bRtf/h3429/+ZixevNjYt2+fsWPHDuPJJ580AGPhwoW1bbzld6Al51/DF74DZs2aZaxatco4fPiwsW3bNuOpp54y7Ha7sWTJEsMwfP/z39zz96XPfo3mvge+9Pk/l88lNzXD2s6/3XXXXYZhGMZdd91lTJw4sc5j8vPzjZCQEOP111+v9zl//vOfG8nJyUZgYKDRqVMnY8qUKcbq1atdfCbNV995A8a8efNq29R3/itWrDCGDBliBAYGGl27djXmzJlzwXP/61//Mnr16mUEBAQYvXv3rvOL70la8h5MnDix0d8Zw/Dt34E///nPRvfu3Y3g4GCjQ4cOxvjx443PPvvsguf2ht+Bln4GfOU74Cc/+YnRpUuX2jgnT55c+0fNMHz/89/c8/elz36N5r4HvvT5P5fNMAzDVb1CIiIiIu6mmhsRERHxKUpuRERExKcouRERERGfouRGREREfIqSGxEREfEpSm5ERETEpyi5EREREZ+i5EZERER8ipIbERER8SlKbkRERMSnKLkRERERn6LkRkRERHzK/wcBqM4zxVavSgAAAABJRU5ErkJggg==",
      "text/plain": [
       "<Figure size 640x480 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "sns.kdeplot(hbondsih.results['hbonds'][:,4])\n",
    "sns.kdeplot(hbondsliq.results['hbonds'][:,4])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 510,
   "id": "dd9217d8",
   "metadata": {},
   "outputs": [],
   "source": [
    "@njit\n",
    "def localHBs(tid, data):\n",
    "    nd = data[(data[:,1] == tid) | (data[:,3] == tid)]\n",
    "    data1 = nd[:,1]\n",
    "    data2 = nd[:,3]\n",
    "    return data1, data2"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 489,
   "id": "04f4fa01",
   "metadata": {},
   "outputs": [],
   "source": [
    "@njit\n",
    "def neighbour(tid, arr):\n",
    "    return np.array(list(arr[arr[:,0] == tid][:,1]) + list(arr[arr[:,1] == tid][:,0]))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 519,
   "id": "86ac4661",
   "metadata": {},
   "outputs": [],
   "source": [
    "@njit\n",
    "def localHBsInfo(tid, arr,hbonds, posi):\n",
    "    arr_neighbour = neighbour(tid, arr)\n",
    "    #arr_neighbour = neighbours(tid, arr)\n",
    "    h_neighbour = 4*arr_neighbour\n",
    "\n",
    "    my_list = []\n",
    "    dt = np.column_stack(localHBs(tid, hbonds)).astype(np.int32)\n",
    "    for data in dt:\n",
    "        my_list.append(data)\n",
    "    for nns in h_neighbour:\n",
    "        list_pair = np.column_stack(localHBs(nns,hbonds)).astype(np.int32)\n",
    "        for pair in list_pair:\n",
    "            if arrayChacker(my_list, pair):\n",
    "                d1 = np.linalg.norm(posi[tid] - posi[pair[0]])\n",
    "                d2 = np.linalg.norm(posi[tid] - posi[pair[1]])\n",
    "                if d1 < 5.0 and d2 < 5.0:\n",
    "\n",
    "                    my_list.append(pair)\n",
    "            \n",
    "    return my_list"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 520,
   "id": "3ca44f8d",
   "metadata": {},
   "outputs": [],
   "source": [
    "@njit\n",
    "def arrayChacker(datalist, pair):\n",
    "    for data in datalist:\n",
    "        if (data == pair).all():\n",
    "            return True\n",
    "    return False"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 521,
   "id": "1246afcd",
   "metadata": {},
   "outputs": [],
   "source": [
    "def localHBDistribution(universe, hbond_obj) :\n",
    "    pos = universe.select_atoms(\"name OW\").positions\n",
    "    posi = universe.select_atoms(\"all\").positions\n",
    "    cutoff = 3.6 \n",
    "    \n",
    "    gridsearch = FastNS(cutoff, pos, universe.dimensions, pbc=True)\n",
    "\n",
    "    results = gridsearch.self_search()\n",
    "    arr = results.get_pairs()\n",
    "    hbih_nums = Parallel(n_jobs=-1, backend = \"threading\")(delayed(localHBsInfo)(i*4,arr,hbond_obj, posi) for i in tqdm.trange(len(pos)))\n",
    "    number_list = list(map(lambda x : len(set(tuple(i) for i in x)), hbih_nums))\n",
    "    return np.array(number_list)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "11b47cfd",
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": 522,
   "id": "bd32de61",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "application/vnd.jupyter.widget-view+json": {
       "model_id": "66bde622c80947748557060e67dd397f",
       "version_major": 2,
       "version_minor": 0
      },
      "text/plain": [
       "  0%|          | 0/11664 [00:00<?, ?it/s]"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "hbih_nums = localHBDistribution(ih, hbondsih.results.hbonds)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 523,
   "id": "dbe05e5b",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "application/vnd.jupyter.widget-view+json": {
       "model_id": "efa807622d7648108eee149ffef456dd",
       "version_major": 2,
       "version_minor": 0
      },
      "text/plain": [
       "  0%|          | 0/11222 [00:00<?, ?it/s]"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "hbliq_nums = localHBDistribution(liq, hbondsliq.results.hbonds)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 539,
   "id": "c5a4e14d",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "application/vnd.jupyter.widget-view+json": {
       "model_id": "67b54c12c79b40f68a28fbe3d86f0eae",
       "version_major": 2,
       "version_minor": 0
      },
      "text/plain": [
       "  0%|          | 0/13600 [00:00<?, ?it/s]"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "hbic_nums = localHBDistribution(ic, hbondsic.results.hbonds)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 546,
   "id": "279a6146",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<AxesSubplot:ylabel='Density'>"
      ]
     },
     "execution_count": 546,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAkkAAAGdCAYAAAAGx+eQAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjQuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/MnkTPAAAACXBIWXMAAA9hAAAPYQGoP6dpAAA6rUlEQVR4nO3de3TU1b3//9ckk4sCSbk1IcoliC1XLyQakhqpFYNBW1E8Rs85Ea36bX4962DIYVUurVKsBhUtUm4HBTVakbYB6zpGSbQSsYxUMCC1qYdKJIDJSYM1USm5zf79kcwnGWYSkhAyn3Gej7VGzJ49n8/eDp+Vt+99cxhjjAAAAOAlLNANAAAAsCOCJAAAAD8IkgAAAPwgSAIAAPCDIAkAAMAPgiQAAAA/CJIAAAD8IEgCAADwwxnoBgQrt9utTz/9VIMGDZLD4Qh0cwAAQDcYY/TFF18oISFBYWFd54oIknrp008/1ciRIwPdDAAA0AtHjhzR+eef32UdgqReGjRokKTW/8gxMTEBbg0AAOiO+vp6jRw50vo93hWCpF7yDLHFxMQQJAEAEGS6M1WGidsAAAB+ECQBAAD4QZAEAADgB0ESAACAHwEPktauXavExERFR0crKSlJO3fu7LJ+aWmpkpKSFB0drbFjx2r9+vVe72/dulXJycn6xje+oQEDBuiSSy7R888/71Vn6dKlcjgcXq/4+Pg+7xsAAAheAQ2StmzZotzcXC1ZskRlZWVKT09XZmamKisr/davqKjQrFmzlJ6errKyMi1evFjz5s1TYWGhVWfIkCFasmSJXC6XPvjgA91555268847tX37dq9rTZo0SVVVVdbrwIEDZ7WvAAAguDiMMSZQN09JSdHUqVO1bt06q2zChAmaPXu28vPzferfd999euWVV1ReXm6V5eTkaP/+/XK5XJ3eZ+rUqbruuuv04IMPSmrNJL388svat29fr9teX1+v2NhY1dXVsQUAAABBoie/vwOWSWpsbNTevXuVkZHhVZ6RkaFdu3b5/YzL5fKpP3PmTO3Zs0dNTU0+9Y0xevPNN/XRRx/pyiuv9Hrv4MGDSkhIUGJiom699VYdOnSoy/Y2NDSovr7e6wUAAL6+AhYk1dbWqqWlRXFxcV7lcXFxqq6u9vuZ6upqv/Wbm5tVW1trldXV1WngwIGKjIzUddddp1/96le65pprrPdTUlJUUFCg7du366mnnlJ1dbXS0tJ0/PjxTtubn5+v2NhY68WRJAAAfL0FfOL2qTteGmO63AXTX/1TywcNGqR9+/bpvffe00MPPaS8vDzt2LHDej8zM1Nz5szRlClTNGPGDL366quSpOeee67T+y5atEh1dXXW68iRI93uIwAACD4BO5Zk2LBhCg8P98ka1dTU+GSLPOLj4/3WdzqdGjp0qFUWFhamcePGSZIuueQSlZeXKz8/X9/97nf9XnfAgAGaMmWKDh482Gl7o6KiFBUV1Z2uAQCAr4GAZZIiIyOVlJSkkpISr/KSkhKlpaX5/UxqaqpP/eLiYiUnJysiIqLTexlj1NDQ0On7DQ0NKi8v14gRI3rQAwAA8HUW0ANu8/LylJ2dreTkZKWmpmrDhg2qrKxUTk6OpNYhrmPHjqmgoEBS60q21atXKy8vT/fcc49cLpc2btyozZs3W9fMz89XcnKyLrjgAjU2NqqoqEgFBQVeK+gWLFig73//+xo1apRqamr0i1/8QvX19Zo7d27//gcAAAC2FdAgKSsrS8ePH9eyZctUVVWlyZMnq6ioSKNHj5YkVVVVee2ZlJiYqKKiIs2fP19r1qxRQkKCVq1apTlz5lh1vvrqK/34xz/W0aNHdc4552j8+PF64YUXlJWVZdU5evSobrvtNtXW1mr48OGaNm2a3n33Xeu+AAAgcF79oErP/LFCV1w4TLkzvhWwdgR0n6Rgxj5JAACcHRvfqdCD//MX/eDiBK267dI+vXZQ7JMEAADgjyd/E9b5Yvd+QZAEAABsxTPGFdbFlkD9gSAJAADYitsTJZFJAgAAaOeZLE0mCQAAoAM3c5IAAAB8tY+2kUkCAACwWKvbAhylECQBAABbcXsyScxJAgAAaOeZkxTgKUkESQAAwF7YJwkAAMAPdtwGAADwgzlJAAAAfpi27SQDHCMRJAEAAHtxs08SAACAL3bcBgAA8Mezui3AURJBEgAAsBX2SQIAAPDDsLoNAADAl9vaTDKw7SBIAgAAtmINtxEkAQAA+OJYEgAAgA7aM0kESQAAABZWtwEAAPhhrInbZJIAAAAsrG4DAADww7C6DQAAwBebSQIAAPjRfsAtQRIAAIDFbWWSAtsOgiQAAGArRp5MUmDbQZAEAABsxZqTFOCdkgiSAACArXB2GwAAgB9sJgkAAOAHmSQAAAA/BjZ9pgmOwxrYUBPQdhAkAQAAW0mtL9JrUYt0ycfrA9oOgiQAAGAvxt36J3OSAAAA2jmsICmwYQpBEgAAsCeCJAAAgA7aMkkm1Ifb1q5dq8TEREVHRyspKUk7d+7ssn5paamSkpIUHR2tsWPHav1670ldW7duVXJysr7xjW9owIABuuSSS/T888+f8X0BAED/8Ay3OUI5k7Rlyxbl5uZqyZIlKisrU3p6ujIzM1VZWem3fkVFhWbNmqX09HSVlZVp8eLFmjdvngoLC606Q4YM0ZIlS+RyufTBBx/ozjvv1J133qnt27f3+r4AAKAfWZmkwAZJDmM8+1r2v5SUFE2dOlXr1q2zyiZMmKDZs2crPz/fp/59992nV155ReXl5VZZTk6O9u/fL5fL1el9pk6dquuuu04PPvhgr+7rT319vWJjY1VXV6eYmJhufQYAAJze/zzxI11f/5L+d+zt+tbtv+rTa/fk93fAQrTGxkbt3btXGRkZXuUZGRnatWuX38+4XC6f+jNnztSePXvU1NTkU98YozfffFMfffSRrrzyyl7fV5IaGhpUX1/v9QIAAGeBZ7gtVA+4ra2tVUtLi+Li4rzK4+LiVF1d7fcz1dXVfus3NzertrbWKqurq9PAgQMVGRmp6667Tr/61a90zTXX9Pq+kpSfn6/Y2FjrNXLkyB71FwAAdI9DnsPbQnhOkiQ5Tpm5bozxKTtd/VPLBw0apH379um9997TQw89pLy8PO3YseOM7rto0SLV1dVZryNHjnTZLwAA0Es22SfJGagbDxs2TOHh4T7Zm5qaGp8sj0d8fLzf+k6nU0OHDrXKwsLCNG7cOEnSJZdcovLycuXn5+u73/1ur+4rSVFRUYqKiupRHwEAQM85PNOlQ3V1W2RkpJKSklRSUuJVXlJSorS0NL+fSU1N9alfXFys5ORkRUREdHovY4waGhp6fV8AANCf2jJJAR7wClgmSZLy8vKUnZ2t5ORkpaamasOGDaqsrFROTo6k1iGuY8eOqaCgQFLrSrbVq1crLy9P99xzj1wulzZu3KjNmzdb18zPz1dycrIuuOACNTY2qqioSAUFBV4r2U53XwAAEDh2OZYkoEFSVlaWjh8/rmXLlqmqqkqTJ09WUVGRRo8eLUmqqqry2rsoMTFRRUVFmj9/vtasWaOEhAStWrVKc+bMsep89dVX+vGPf6yjR4/qnHPO0fjx4/XCCy8oKyur2/cFAACB45m43dVc4X5pRyD3SQpm7JMEAMDZ8foj/6Zr//k/+tuE/9C4rIf79NpBsU8SAACAf23DbaG+BQAAAEBH7avbwgPaDoIkAABgKw5rdVuI7rgNAADgj2d1m4PhNgAAgHbWsSQBXt1GkAQAAOwl1HfcBgAA8M8em0kSJAEAAFsJsw6vZ3UbAACAxVrdFsacJAAAAEv7xG2G2wAAACzWFgAESQAAAB15jpVluA0AAMDiOZbEEcbEbQAAAIs1cZvNJAEAANpZE7fJJAEAALRj4jYAAIAfnkySg+E2AACAdmSSAAAA/LAySWEESQAAABZPkGTIJAEAALTzbAHAcBsAAEAHDLcBAAD4YU3cDnCYQpAEAABsJczKJLEFAAAAQAdtO2472HEbAADA0n4sCcNtAAAAlrC2OUlhBEkAAADt2o8lIUgCAACwePZJYk4SAABAB9bqtsAubiNIAgAAdsPEbQAAAB+eTFIYw20AAADtrLPbyCQBAAC0CzOezSQJkgAAACwccAsAAOBHmGe4jUwSAABAK+MZapPkCPAeAARJAADANtymQyYpjNVtAAAAklozSdYWAARJAAAArVozSaxukyStXbtWiYmJio6OVlJSknbu3Nll/dLSUiUlJSk6Olpjx47V+vXrvd5/6qmnlJ6ersGDB2vw4MGaMWOG/vSnP3nVWbp0qRwOh9crPj6+z/sGAAB6xm2MHA7P6rYQnpO0ZcsW5ebmasmSJSorK1N6eroyMzNVWVnpt35FRYVmzZql9PR0lZWVafHixZo3b54KCwutOjt27NBtt92mt956Sy6XS6NGjVJGRoaOHTvmda1JkyapqqrKeh04cOCs9hUAAHSPZ05SWIAzSQ7TcRp5P0tJSdHUqVO1bt06q2zChAmaPXu28vPzferfd999euWVV1ReXm6V5eTkaP/+/XK5XH7v0dLSosGDB2v16tW6/fbbJbVmkl5++WXt27ev122vr69XbGys6urqFBMT0+vrAACAdicam1X30IUa4fhM//zhWzpn1NQ+vX5Pfn8HLERrbGzU3r17lZGR4VWekZGhXbt2+f2My+XyqT9z5kzt2bNHTU1Nfj9z4sQJNTU1aciQIV7lBw8eVEJCghITE3Xrrbfq0KFDZ9AbAADQFzqubgvZidu1tbVqaWlRXFycV3lcXJyqq6v9fqa6utpv/ebmZtXW1vr9zMKFC3XeeedpxowZVllKSooKCgq0fft2PfXUU6qurlZaWpqOHz/eaXsbGhpUX1/v9QIAAH2r4+q2kN9x+9SNoowxXW4e5a++v3JJevTRR7V582Zt3bpV0dHRVnlmZqbmzJmjKVOmaMaMGXr11VclSc8991yn983Pz1dsbKz1Gjly5Ok7BwAAesRtOhxLEqqr24YNG6bw8HCfrFFNTY1PtsgjPj7eb32n06mhQ4d6la9YsUIPP/ywiouLddFFF3XZlgEDBmjKlCk6ePBgp3UWLVqkuro663XkyJEurwkAAHquNZMU4seSREZGKikpSSUlJV7lJSUlSktL8/uZ1NRUn/rFxcVKTk5WRESEVfbYY4/pwQcf1Ouvv67k5OTTtqWhoUHl5eUaMWJEp3WioqIUExPj9QIAAH3LdNgnKaSH2/Ly8vT0009r06ZNKi8v1/z581VZWamcnBxJrdkbz4o0qXUl2+HDh5WXl6fy8nJt2rRJGzdu1IIFC6w6jz76qH76059q06ZNGjNmjKqrq1VdXa0vv/zSqrNgwQKVlpaqoqJCu3fv1s0336z6+nrNnTu3/zoPAAB8uL123A5skOQM5M2zsrJ0/PhxLVu2TFVVVZo8ebKKioo0evRoSVJVVZXXnkmJiYkqKirS/PnztWbNGiUkJGjVqlWaM2eOVWft2rVqbGzUzTff7HWvBx54QEuXLpUkHT16VLfddptqa2s1fPhwTZs2Te+++651XwAAEBimwz8DPdwW0H2Sghn7JAEA0Pdqvjipc1eM1kDHSWnePmlIYp9ePyj2SQIAADiV4ew2AAAAX6bDZpIESQAAAG3cxlj7JKmLfRP7A0ESAACwjY6r28gkAQAAtDEddtwmSAIAAGjDxG0AAAA/3G63whwESQAAAF6MOmzfSJAEAADQyu1uCXQTLARJAADANkzHIIlMEgAAQBvDcBsAAIAPt9vd/gNBEgAAQCuCJAAAAH8IkgAAAHy5DRO3AQAAfBivTBIH3AIAALRidRsAAIAvt2FOEgAAgC+vzSQZbgMAAJDUvgVAiw1ClMC3AAAAwKNtuM0osFkkiSAJAADYiCFIAgAA8OXZAoAgCQAAoAO3u3ULADdBEgAAQAdtO24bG4QogW8BAACAR9ucJDJJAAAAHbitOUmBD1EC3wIAAAAPz+q2wCeSCJIAAIB9GDdzkgAAAHyYtgNu2QIAAACgI8OcJAAAAB+ezSTdAT7cViJIAgAANmLIJAEAAPjyBEl2QJAEAADso211m9sGIUrgWwAAANCmbXEbq9sAAAC8WJtJBj5ECXwLAAAA2rRvJkkmCQAAoAPPZpKBD1EC3wIAAIA2ViaJfZIAAADate+TFKRBUkVFRZ81YO3atUpMTFR0dLSSkpK0c+fOLuuXlpYqKSlJ0dHRGjt2rNavX+/1/lNPPaX09HQNHjxYgwcP1owZM/SnP/3pjO8LAAD6gefstmCduD1u3DhdddVVeuGFF3Ty5Mle33zLli3Kzc3VkiVLVFZWpvT0dGVmZqqystJv/YqKCs2aNUvp6ekqKyvT4sWLNW/ePBUWFlp1duzYodtuu01vvfWWXC6XRo0apYyMDB07dqzX9wUAAP3DTpkkh/Ect9sDf/7zn7Vp0yb9+te/VkNDg7KysnTXXXfp8ssv79F1UlJSNHXqVK1bt84qmzBhgmbPnq38/Hyf+vfdd59eeeUVlZeXW2U5OTnav3+/XC6X33u0tLRo8ODBWr16tW6//fZe3def+vp6xcbGqq6uTjExMd36DAAA6NqO136j7+6+R5URYzVqSVmfX78nv797lUmaPHmynnjiCR07dkzPPPOMqqurdcUVV2jSpEl64okn9Pe///2012hsbNTevXuVkZHhVZ6RkaFdu3b5/YzL5fKpP3PmTO3Zs0dNTU1+P3PixAk1NTVpyJAhvb6vJDU0NKi+vt7rBQAA+pbDcyxJsE/cdjqduvHGG/Wb3/xGjzzyiD7++GMtWLBA559/vm6//XZVVVV1+tna2lq1tLQoLi7OqzwuLk7V1dV+P1NdXe23fnNzs2pra/1+ZuHChTrvvPM0Y8aMXt9XkvLz8xUbG2u9Ro4c2WldAADQO3YabjujIGnPnj368Y9/rBEjRuiJJ57QggUL9PHHH+sPf/iDjh07phtuuOG013CcEikaY3zKTlffX7kkPfroo9q8ebO2bt2q6OjoM7rvokWLVFdXZ72OHDnSaV0AANBLbk+QFPiJ287efOiJJ57QM888o48++kizZs1SQUGBZs2apbCw1g4lJibqv//7vzV+/PhOrzFs2DCFh4f7ZG9qamp8sjwe8fHxfus7nU4NHTrUq3zFihV6+OGH9cYbb+iiiy46o/tKUlRUlKKiojp9HwAAnDnj2UwyWIfb1q1bp3/9139VZWWlXn75ZV1//fVWgOQxatQobdy4sdNrREZGKikpSSUlJV7lJSUlSktL8/uZ1NRUn/rFxcVKTk5WRESEVfbYY4/pwQcf1Ouvv67k5OQzvi8AAOgnJsgzSSUlJRo1apRPYGSM0ZEjRzRq1ChFRkZq7ty5XV4nLy9P2dnZSk5OVmpqqjZs2KDKykrl5ORIah3iOnbsmAoKCiS1rmRbvXq18vLydM8998jlcmnjxo3avHmzdc1HH31UP/vZz/Tiiy9qzJgxVsZo4MCBGjhwYLfuCwAAAsRGE7d7FSRdcMEFqqqq0je/+U2v8s8++0yJiYlqaWnp1nWysrJ0/PhxLVu2TFVVVZo8ebKKioo0evRoSVJVVZXX3kWJiYkqKirS/PnztWbNGiUkJGjVqlWaM2eOVWft2rVqbGzUzTff7HWvBx54QEuXLu3WfQEAQIDYaE5Sr/ZJCgsLU3V1tU+QdPjwYU2cOFFfffVVnzXQrtgnCQCAvveHrU/pex8s0N/OmaJx973T59fvye/vHmWS8vLyJLWuDLv//vt17rnnWu+1tLRo9+7duuSSS3reYgAAALVvAWCH42V7FCSVlbXufGmM0YEDBxQZGWm9FxkZqYsvvlgLFizo2xYCAICQ4dlM0g6r23oUJL311luSpDvvvFNPPvkkw0wAAKBPeWYB2WFOUq8mbj/zzDN93Q4AAID2LQCCKZN000036dlnn1VMTIxuuummLutu3br1jBsGAABCUDDOSYqNjbWO7YiNjT1rDQIAACHMyiQFUZDUcYiN4TYAAHBWWJmkwA+39SpM++c//6kTJ05YPx8+fFgrV65UcXFxnzUMAACEIBvNSepVkHTDDTdYR4V8/vnnuvzyy/X444/rhhtu0Lp16/q0gQAAIHTYaXVbr1rw/vvvKz09XZL0u9/9TvHx8Tp8+LAKCgq0atWqPm0gAAAIIZ6DQGwwJ6lXLThx4oQGDRokSSouLtZNN92ksLAwTZs2TYcPH+7TBgIAgBAS7MNt48aN08svv6wjR45o+/btysjIkCTV1NSwwSQAAOg909L6Z7Bmku6//34tWLBAY8aMUUpKilJTUyW1ZpUuvfTSPm0gAAAIHQ7PcJsNVrf1asftm2++WVdccYWqqqp08cUXW+VXX321brzxxj5rHAAACC2eA27tMHG7V0GSJMXHxys+Pt6r7PLLLz/jBgEAgBBmo4nbvQqSvvrqKy1fvlxvvvmmampq5Ha7vd4/dOhQnzQOAACEGBtN3O5VkHT33XertLRU2dnZGjFihHVcCQAAwBnx7LgdrJmk1157Ta+++qq+853v9HV7AABAKLPRAbe9asHgwYM1ZMiQvm4LAAAIcZ7VbcYGg1S9CpIefPBB3X///V7ntwEAAJy5IB9ue/zxx/Xxxx8rLi5OY8aMUUREhNf777//fp80DgAAhBYT7KvbZs+e3cfNAAAAkBzBvk/SAw880NftAAAAUPtwW+AnJfU6TPv888/19NNPa9GiRfrss88ktQ6zHTt2rM8aBwAAQkywbwHwwQcfaMaMGYqNjdUnn3yie+65R0OGDNG2bdt0+PBhFRQU9HU7AQBACLDObgvWTFJeXp7uuOMOHTx4UNHR0VZ5Zmam3n777T5rHAAACDHBvk/Se++9px/96Ec+5eedd56qq6vPuFEAACBE2Wh1W69aEB0drfr6ep/yjz76SMOHDz/jRgEAgBBlnd0WpEHSDTfcoGXLlqmpqUmS5HA4VFlZqYULF2rOnDl92kAAABBKgnx124oVK/T3v/9d3/zmN/XPf/5T06dP17hx4zRo0CA99NBDfd1GAAAQIhw2Gm7r1eq2mJgYvfPOO3rrrbe0d+9eud1uTZ06VTNmzOjr9gEAgFBiTdwOfCapx0GS2+3Ws88+q61bt+qTTz6Rw+FQYmKi4uPjZYyRwwbpMQAAEJwcNjq7rUctMMboBz/4ge6++24dO3ZMU6ZM0aRJk3T48GHdcccduvHGG89WOwEAQAhwBOtmks8++6zefvttvfnmm7rqqqu83vvDH/6g2bNnq6CgQLfffnufNhIAAIQGz5QkOwRJPWrB5s2btXjxYp8ASZK+973vaeHChfr1r3/dZ40DAAChJWiH2z744ANde+21nb6fmZmp/fv3n3GjAABAiDJBugXAZ599pri4uE7fj4uL0z/+8Y8zbhQAAAhNQZtJamlpkdPZ+TSm8PBwNTc3n3GjAABAiArWfZKMMbrjjjsUFRXl9/2GhoY+aRQAAAhNDhsNt/UoSJo7d+5p67CyDQAA9J4nkxQe2Gaoh0HSM8880+cNWLt2rR577DFVVVVp0qRJWrlypdLT0zutX1paqry8PH344YdKSEjQT37yE+Xk5Fjvf/jhh7r//vu1d+9eHT58WL/85S+Vm5vrdY2lS5fq5z//uVdZXFycqqur+7RvAACgZ+yUSQrogN+WLVuUm5urJUuWqKysTOnp6crMzFRlZaXf+hUVFZo1a5bS09NVVlamxYsXa968eSosLLTqnDhxQmPHjtXy5csVHx/f6b0nTZqkqqoq63XgwIE+7x8AAOiZoD+7ra888cQTuuuuu3T33XdLklauXKnt27dr3bp1ys/P96m/fv16jRo1SitXrpQkTZgwQXv27NGKFSs0Z84cSdJll12myy67TJK0cOHCTu/tdDq7DKIAAEAgBOnqtr7U2NiovXv3KiMjw6s8IyNDu3bt8vsZl8vlU3/mzJnas2ePmpqaenT/gwcPKiEhQYmJibr11lt16NChLus3NDSovr7e6wUAAPqWJ5Nkh7NgAxYk1dbWqqWlxWffpa7mBlVXV/ut39zcrNra2m7fOyUlRQUFBdq+fbueeuopVVdXKy0tTcePH+/0M/n5+YqNjbVeI0eO7Pb9AABA9wTtPklnw6mRojGmy+jRX31/5V3JzMzUnDlzNGXKFM2YMUOvvvqqJOm5557r9DOLFi1SXV2d9Tpy5Ei37wcAALqLOUkaNmyYwsPDfbJGNTU1ne7qHR8f77e+0+nU0KFDe92WAQMGaMqUKTp48GCndaKiojrdHwoAAPSN9tVtgQ+SAtaCyMhIJSUlqaSkxKu8pKREaWlpfj+TmprqU7+4uFjJycmKiIjodVsaGhpUXl6uESNG9PoaAADgzNlpdVtAW5CXl6enn35amzZtUnl5uebPn6/Kykpr36NFixZ5bU6Zk5Ojw4cPKy8vT+Xl5dq0aZM2btyoBQsWWHUaGxu1b98+7du3T42NjTp27Jj27dunv/3tb1adBQsWqLS0VBUVFdq9e7duvvlm1dfXd2uzTAAAcDa1ZpIcYYEPkgK6BUBWVpaOHz+uZcuWqaqqSpMnT1ZRUZFGjx4tSaqqqvLaMykxMVFFRUWaP3++1qxZo4SEBK1atcpa/i9Jn376qS699FLr5xUrVmjFihWaPn26duzYIUk6evSobrvtNtXW1mr48OGaNm2a3n33Xeu+AAAgMKxMkgK/us1hjNUa9EB9fb1iY2NVV1enmJiYQDcHAICvhffyZ+qyhnd1YOoyTfnBvX1+/Z78/g58LgsAAMCDOUkAAAC+HDbaAiDwLQAAAGjj2UzSQZAEAADQji0AAAAA/GjPJAV+dRtBEgAAsA3ObgMAAPDHcyarDTaTDHwLAAAA2rC6DQAAwA/PAbdkkgAAADoIszJJ4YFtiAiSAACAjVir2wLcDokgCQAA2Ig1JymMTBIAAIDFMyeJidsAAAAdeIbZmLgNAADQAWe3AQAA+OGwNpMM/NRtgiQAAGAbHEsCAADgh2d1m4N9kgAAANqFWXOSGG4DAADwQSYJAACgA2ufJLYAAAAAaOc5u43VbQAAAB2wTxIAAIAfYaxuAwAA8OVguA0AAMBX+z5JgQ9RAt8CAACANtaO22EMtwEAAFjCyCQBAAD4aj/gNvAhSuBbAAAA0IZjSQAAAPzwTNwOY04SAABAO0+QJDJJAAAA7dr3SSKTBAAAYGk/uy3wIUrgWwAAANDGM3E7jGNJAAAA2rXvuM2cJAAAAItnuE0MtwEAALQyxrRvAcCO2wAAAK3cpuPEbeYkAQAASPJkktp23Ga4TVq7dq0SExMVHR2tpKQk7dy5s8v6paWlSkpKUnR0tMaOHav169d7vf/hhx9qzpw5GjNmjBwOh1auXNkn9wUAAGeXdyYp4CFKYIOkLVu2KDc3V0uWLFFZWZnS09OVmZmpyspKv/UrKio0a9Yspaenq6ysTIsXL9a8efNUWFho1Tlx4oTGjh2r5cuXKz4+vk/uCwAAzj4j0x4k2WB1m8OYtuN2AyAlJUVTp07VunXrrLIJEyZo9uzZys/P96l/33336ZVXXlF5eblVlpOTo/3798vlcvnUHzNmjHJzc5Wbm3tG9/Wnvr5esbGxqqurU0xMTLc+AwAAOneysVnRDw+VJH01768aMGREn9+jJ7+/A5ZJamxs1N69e5WRkeFVnpGRoV27dvn9jMvl8qk/c+ZM7dmzR01NTWftvpLU0NCg+vp6rxcAAOg7buO2/j2kh9tqa2vV0tKiuLg4r/K4uDhVV1f7/Ux1dbXf+s3NzaqtrT1r95Wk/Px8xcbGWq+RI0d2634AAKB7jLt9cMvBjtu+Y47GmC7HIf3V91fe1/ddtGiR6urqrNeRI0d6dD8AANA1t7vF+nc7ZJKcgbrxsGHDFB4e7pO9qamp8cnyeMTHx/ut73Q6NXTo0LN2X0mKiopSVFRUt+4BAAB6zu3uONwW+InbAQvTIiMjlZSUpJKSEq/ykpISpaWl+f1MamqqT/3i4mIlJycrIiLirN0XAAD0gw5Bkh123A5YJkmS8vLylJ2dreTkZKWmpmrDhg2qrKxUTk6OpNYhrmPHjqmgoEBS60q21atXKy8vT/fcc49cLpc2btyozZs3W9dsbGzUX/7yF+vfjx07pn379mngwIEaN25ct+4LAAD6X8fhtjAb7Lgd0CApKytLx48f17Jly1RVVaXJkyerqKhIo0ePliRVVVV57V2UmJiooqIizZ8/X2vWrFFCQoJWrVqlOXPmWHU+/fRTXXrppdbPK1as0IoVKzR9+nTt2LGjW/cFAAD9z2t1mw0ySQHdJymYsU8SAAB96/jxWg391QWSJLPk/+SIiO7zewTFPkkAAAAdea9uC/HhNvhRUy7tf0mKOU9K+X+Bbg0AAP3GdJi4LRsMtwW+BfD2WYX0x5XSB1sC3RIAAPqV6TAnSQrhLQDQifDI1j9bGgPbDgAA+pnpMNwmGxxwS5BkN+Ft+z0RJAEAQoxnM0m3cRAkwQ8ySQCAEOVZcO+2wVCbRJBkP05PkNQU2HYAANDf2jJJhiAJfpFJAgCEKM8WAG6bhCf2aAXaeYKk5obAtgMAgP7WNtxml12uCZJsZu/RLyVJjY0ESQCA0OI2rZkkY4M9kiSCJNs53tA6DhvmZk4SACC0eDaTZOI2/HI6o1r/VLOVdgQAICQYgiR0ISIyqv0HJm8DAEKI2+2Zk2SP8MQerYDFSZAEAAhVxrO6jUwS/HBGRrf/wF5JAIAQ4rYOuCVIgh8RERFqMW1/OcgkAQBCCBO30aXI8DA1ydn6A3slAQBCifHsuG2P8MQerYAl0hmmRnkOuWW4DQAQOgyr29CVKGeYGj2ZJIbbAAAhxLjJJKELkc4Ow20tDLcBAEKHJ5NkHGSS4EdkeJiaTLgkqaWJTBIAIHQwcRtdiugwJ6mZ89sAAKGkbZ8ktgCAXx1XtzU3ESQBAEKHaTuOy22T8MQerYAlItxhTdxubjwZ4NYAANB/jLs1k2TIJMEfh8OhZkdbkMQ+SQCAEOLJJBEkoVPNbXOSWpiTBAAIJdbqNnuEJ/ZoBby0ONqCJFa3AQBCSPs+SWSS0ImWtuG2libmJAEAQoe1TxJBEjrTEtaaSXI3k0kCAIQQhttwOi2OSEmSmy0AAAAhxLO6jc0k0Sm3J5PE2W0AgBDE2W3olCdIMmSSAAAhxNonibPb0Bk3c5IAACGofXWbPcITe7QCXtxhrXOSxHAbACCkmEA3wAtBkg2Z8LbhNjJJAIAQ0n4siT3CE3u0At48c5JamJMEAAghbAGA0zHhnuG2psA2BACA/sTZbTitcOYkAQBCjyGThNPxZJIcZJIAACGEs9twWg4rSGJOEgAglLjb/iRIkiStXbtWiYmJio6OVlJSknbu3Nll/dLSUiUlJSk6Olpjx47V+vXrfeoUFhZq4sSJioqK0sSJE7Vt2zav95cuXSqHw+H1io+P79N+nQmHsy1IcpNJAgCEEE8mic0kpS1btig3N1dLlixRWVmZ0tPTlZmZqcrKSr/1KyoqNGvWLKWnp6usrEyLFy/WvHnzVFhYaNVxuVzKyspSdna29u/fr+zsbN1yyy3avXu317UmTZqkqqoq63XgwIGz2tceacskhREkAQBCiWEzScsTTzyhu+66S3fffbcmTJiglStXauTIkVq3bp3f+uvXr9eoUaO0cuVKTZgwQXfffbd++MMfasWKFVadlStX6pprrtGiRYs0fvx4LVq0SFdffbVWrlzpdS2n06n4+HjrNXz48LPZ1R4JI5MEAAhBxrO6LdQnbjc2Nmrv3r3KyMjwKs/IyNCuXbv8fsblcvnUnzlzpvbs2aOmpqYu65x6zYMHDyohIUGJiYm69dZbdejQoS7b29DQoPr6eq/X2eKIiJIkhbG6DQAQSgwTtyVJtbW1amlpUVxcnFd5XFycqqur/X6murrab/3m5mbV1tZ2WafjNVNSUlRQUKDt27frqaeeUnV1tdLS0nT8+PFO25ufn6/Y2FjrNXLkyB71tyc8E7fDDJkkAEDoYAuAUzhOmZxljPEpO139U8tPd83MzEzNmTNHU6ZM0YwZM/Tqq69Kkp577rlO77to0SLV1dVZryNHjpymZ70X1pZJCme4DQAQQhxue61ucwbqxsOGDVN4eLhP1qimpsYnE+QRHx/vt77T6dTQoUO7rNPZNSVpwIABmjJlig4ePNhpnaioKEVFRXXZp74S3jYnKdw098v9AACwA8NwW6vIyEglJSWppKTEq7ykpERpaWl+P5OamupTv7i4WMnJyYqIiOiyTmfXlFrnG5WXl2vEiBG96UqfC/dkkgxzkgAAIYThtnZ5eXl6+umntWnTJpWXl2v+/PmqrKxUTk6OpNYhrttvv92qn5OTo8OHDysvL0/l5eXatGmTNm7cqAULFlh17r33XhUXF+uRRx7RX//6Vz3yyCN64403lJuba9VZsGCBSktLVVFRod27d+vmm29WfX295s6d229974onSHKSSQIAhBTT9qc9gqSADbdJUlZWlo4fP65ly5apqqpKkydPVlFRkUaPHi1Jqqqq8tozKTExUUVFRZo/f77WrFmjhIQErVq1SnPmzLHqpKWl6aWXXtJPf/pT/exnP9MFF1ygLVu2KCUlxapz9OhR3XbbbaqtrdXw4cM1bdo0vfvuu9Z9A609k8ScJABA6GifuG2P4TaH8cx8Ro/U19crNjZWdXV1iomJ6dNr735/r1Je+Z7+qWids/T/+vTaAADY1Z9eytflf12uPQO/q+QFvz8r9+jJ72975LPgxelszSRFiEwSACB0eDJJdlndRpBkQ87ItjlJarHOsQEA4GvPZsNtBEk2FBEV3f4DeyUBAEKFsdfEbXu0Al4iIjvsx8TRJACAUMEWADgdZ0SHTFIzQRIAIERYmSSG29CJyMgINZu2r4ZMEgAgZJBJwmlEhoepqW0LK9PSEODWAADQTzyr2wiS0JlIZ3uQ1NRIJgkAECI4uw2nE+UMU2NbkNTceDLArQEAoJ9YmSSCJHQiIrw9SGpqZLgNABAiGG7D6YSHOazhNjJJAICQ0ba6jYnb6FJzW5DU0kQmCQAQIthMEt3R7IiQxHAbACCEcCwJusMTJJFJAgCEDA64RXdYw23NBEkAgNDgEBO30Q3uME8miX2SAAAhgtVt6I6WtuE2dxOr2wAAocGwug3d4ZmT5OaAWwBAiHB4VrcRJKErnuE2N8NtAIBQwcRtdIcVJHHALQAgRLRP3CZIQhdMW5Bk2AIAABAqmLiN7nCHRUqSDHOSAAChwtpM0h7hiT1aAR+eTJJamgLbEAAA+g0Tt9ENJrwtk9RCJgkAEBpY3YZuac8kMScJABAiDBO30R3OqNY/GW4DAIQKawsAe4Qn9mgFfIV7MkkMtwEAQoODTBK6pW1OkoNMEgAgZDAnCd3RFiSFuckkAQBCg4MgCd3haJuT5GC4DQAQKtwMt6EbHM624TY3w20AgFDh2UwyPMDtaEWQZFNhTs9wG0ESACA0ePZJcpBJQlfagySG2wAAoYKz29ANnjlJ4aY5wC0BAKB/sAUAuiUsoi1IYrgNABAi2ofb7BGe2KMV8BHeNtwWbgiSAAChwrT90x7hiT1aAR9OTyaJIAkAECo8w21h9ghP7NEK+AiPbA2SnARJAIAQ4dlM0iF7zElyBroB8C+8LZPkFBO3gf70+YlGHT5+QpI0+NxIjRp6boBbBIQOh80ySQRJNmUFSaxuA/rNFyeblPHLt1XzRYNV9vxdlyv9wuEBbBUQOhxsAeBt7dq1SkxMVHR0tJKSkrRz584u65eWliopKUnR0dEaO3as1q9f71OnsLBQEydOVFRUlCZOnKht27ad8X37m2dOUoQYbgP6y6Z3PlHNFw06NzJcQwa0Lp54bPtHMm0rbgCcXZ7VbQRJkrZs2aLc3FwtWbJEZWVlSk9PV2ZmpiorK/3Wr6io0KxZs5Senq6ysjItXrxY8+bNU2FhoVXH5XIpKytL2dnZ2r9/v7Kzs3XLLbdo9+7dvb5vIDgjo1v/ZLgN6Bd1J5r09DuHJEmP3nyRSuZfqXMjw/XB0Tq9UV4T4NYBocJeQZLDBPB/kVJSUjR16lStW7fOKpswYYJmz56t/Px8n/r33XefXnnlFZWXl1tlOTk52r9/v1wulyQpKytL9fX1eu2116w61157rQYPHqzNmzf36r7+1NfXKzY2VnV1dYqJielZx7vh2NFKnff0lNYf7v+HbcZnga+rx4s/0jN/+ED/MuSQfnZlrMIcDhV9VK+H/jJcMSMu0Kv/eYXCwuwxmRT4uvrL8qs08eT72nXxcqXd+P+dlXv05Pd3wOYkNTY2au/evVq4cKFXeUZGhnbt2uX3My6XSxkZGV5lM2fO1MaNG9XU1KSIiAi5XC7Nnz/fp87KlSt7fV9JamhoUEND+zyFuro6Sa3/sc+Gk41Nqm9ojV8PPnC53DaZ6Q98XaU2NWquKhX5jxZ9+fvWsiskveaQPv40Xnvuj7XLJsDA11ZC82HVy+irf548a79fPdftTo4oYEFSbW2tWlpaFBcX51UeFxen6upqv5+prq72W7+5uVm1tbUaMWJEp3U81+zNfSUpPz9fP//5z33KR44c2Xkn+8zefrgHgM59EegGACEmr+119nzxxReKjY3tsk7AV7edetKvMabL03/91T+1vDvX7Ol9Fy1apLy89i/M7Xbrs88+09ChQ/v8tOL6+nqNHDlSR44cOStDeXbwde8j/Qtu9C/4fd37SP96zxijL774QgkJCaetG7AgadiwYQoPD/fJ3tTU1PhkeTzi4+P91nc6nRo6dGiXdTzX7M19JSkqKkpRUVFeZd/4xjc672AfiImJ+Vr+5e/o695H+hfc6F/w+7r3kf71zukySB4Bmw0cGRmppKQklZSUeJWXlJQoLS3N72dSU1N96hcXFys5OVkRERFd1vFcszf3BQAAoSegw215eXnKzs5WcnKyUlNTtWHDBlVWVionJ0dS6xDXsWPHVFBQIKl1Jdvq1auVl5ene+65Ry6XSxs3brRWrUnSvffeqyuvvFKPPPKIbrjhBv3+97/XG2+8oXfeeafb9wUAAJAJsDVr1pjRo0ebyMhIM3XqVFNaWmq9N3fuXDN9+nSv+jt27DCXXnqpiYyMNGPGjDHr1q3zueZvf/tb8+1vf9tERESY8ePHm8LCwh7dN9BOnjxpHnjgAXPy5MlAN+Ws+br3kf4FN/oX/L7ufaR//SOg+yQBAADYFTsUAgAA+EGQBAAA4AdBEgAAgB8ESQAAAH4QJJ1l69at00UXXWRtiJWamup1+K4/paWlSkpKUnR0tMaOHav169f71CksLNTEiRMVFRWliRMnatu2bWerC13qaf+2bt2qa665RsOHD7fqb9++3avOs88+K4fD4fM6efLk2e6Oj572b8eOHX7b/te//tWrnl2+P6nnfbzjjjv89nHSpElWHTt9hx3l5+fL4XAoNze3y3rB9Ayeqjt9DLbnsKPu9C8Yn0OP7vQv2J7BpUuX+rQjPj6+y8/Y5RkkSDrLzj//fC1fvlx79uzRnj179L3vfU833HCDPvzwQ7/1KyoqNGvWLKWnp6usrEyLFy/WvHnzVFhYaNVxuVzKyspSdna29u/fr+zsbN1yyy3avXt3f3XL0tP+vf3227rmmmtUVFSkvXv36qqrrtL3v/99lZWVedWLiYlRVVWV1ys6Oro/uuSlp/3z+Oijj7zafuGFF1rv2en7k3rexyeffNKrb0eOHNGQIUP0L//yL1717PIderz33nvasGGDLrrooi7rBdsz2FF3+xhsz6FHd/vnEUzPodT9/gXjMzhp0iSvdhw4cKDTurZ6BgO6AUGIGjx4sHn66af9vveTn/zEjB8/3qvsRz/6kZk2bZr18y233GKuvfZarzozZ840t956a983the66p8/EydOND//+c+tn5955hkTGxt7FlrWN7rq31tvvWUkmX/84x+dft7u358xPfsOt23bZhwOh/nkk0+sMrt9h1988YW58MILTUlJiZk+fbq59957O60brM9gT/roj92fw570LxifwzP5/uz+DD7wwAPm4osv7nZ9Oz2DZJL6UUtLi1566SV99dVXSk1N9VvH5XIpIyPDq2zmzJnas2ePmpqauqyza9eus9PwbupO/07ldrv1xRdfaMiQIV7lX375pUaPHq3zzz9f119/vc//4QZCT/p36aWXasSIEbr66qv11ltveb1n1+9P6t13uHHjRs2YMUOjR4/2KrfTd/gf//Efuu666zRjxozT1g3WZ7AnfTxVMDyHvelfMD2HZ/L9BcMzePDgQSUkJCgxMVG33nqrDh061GldOz2DAT2WJFQcOHBAqampOnnypAYOHKht27Zp4sSJfutWV1f7HLQbFxen5uZm1dbWasSIEZ3WOfXQ3v7Sk/6d6vHHH9dXX32lW265xSobP368nn32WU2ZMkX19fV68skn9Z3vfEf79+/3Spf3l570b8SIEdqwYYOSkpLU0NCg559/XldffbV27NihK6+8UlLn33Ggvj+p999hVVWVXnvtNb344ote5Xb6Dl966SW9//77eu+997pVPxifwZ728VR2fw572r9gew7P5PsLhmcwJSVFBQUF+ta3vqX/+7//0y9+8QulpaXpww8/tA6n78hWz2Cf5qXgV0NDgzl48KB57733zMKFC82wYcPMhx9+6LfuhRdeaB5++GGvsnfeecdIMlVVVcYYYyIiIsyLL77oVeeFF14wUVFRZ6cDp9GT/nX04osvmnPPPdeUlJR0Wa+lpcVcfPHF5j//8z/7qsk90tv+eVx//fXm+9//vvWz3b4/Y3rfx4cfftgMHTrUNDQ0dFkvUN9hZWWl+eY3v2n27dtnlZ1uKCPYnsHe9LEjuz+HZ9o/D7s+h2faP7s/g/58+eWXJi4uzjz++ON+37fTM8hwWz+IjIzUuHHjlJycrPz8fF188cV68skn/daNj4/3iYRramrkdDqtiLuzOqdG1f2lJ/3z2LJli+666y795je/OW16OSwsTJdddpkOHjzYl83utt70r6Np06Z5td1u35/Uuz4aY7Rp0yZlZ2crMjKyy7qB+g737t2rmpoaJSUlyel0yul0qrS0VKtWrZLT6VRLS4vPZ4LtGexNHz2C4Tk8k/51ZNfn8Ez6FwzPoD8DBgzQlClTOm2LnZ5BgqQAMMaooaHB73upqakqKSnxKisuLlZycrIiIiK6rJOWlnZ2GtxDXfVPkjZv3qw77rhDL774oq677rpuXW/fvn0aMWJEXzaz107Xv1OVlZV5td3u35/UvT6Wlpbqb3/7m+66665uXS8Q3+HVV1+tAwcOaN++fdYrOTlZ//Zv/6Z9+/YpPDzc5zPB9gz2po9S8DyHve3fqez6HJ5J/4LhGfSnoaFB5eXlnbbFVs9gn+al4GPRokXm7bffNhUVFeaDDz4wixcvNmFhYaa4uNgYY8zChQtNdna2Vf/QoUPm3HPPNfPnzzd/+ctfzMaNG01ERIT53e9+Z9X54x//aMLDw83y5ctNeXm5Wb58uXE6nebdd9+1ff9efPFF43Q6zZo1a0xVVZX1+vzzz606S5cuNa+//rr5+OOPTVlZmbnzzjuN0+k0u3fvtn3/fvnLX5pt27aZ//3f/zV//vOfzcKFC40kU1hYaNWx0/dnTM/76PHv//7vJiUlxe817fQdnurUoYxgfwb9OV0fg+05PNXp+heMz2FHp+ufR7A8g//1X/9lduzYYQ4dOmTeffddc/3115tBgwZZq/Hs/AwSJJ1lP/zhD83o0aNNZGSkGT58uLn66qutXz7GGDN37lwzffp0r8/s2LHDXHrppSYyMtKMGTPGrFu3zue6v/3tb823v/1tExERYcaPH+/18PennvZv+vTpRpLPa+7cuVad3NxcM2rUKOuaGRkZZteuXf3Yq3Y97d8jjzxiLrjgAhMdHW0GDx5srrjiCvPqq6/6XNcu358xvfs7+vnnn5tzzjnHbNiwwe817fQdnurUX0DB/gz6c7o+BttzeKrT9S8Yn8OOuvN3NJiewaysLDNixAgTERFhEhISzE033eQ159HOz6DDGGP6NjcFAAAQ/JiTBAAA4AdBEgAAgB8ESQAAAH4QJAEAAPhBkAQAAOAHQRIAAIAfBEkAAAB+ECQBAAD4QZAEAADgB0ESAACAHwRJAAAAfhAkAQAA+PH/A91tgjts5b7jAAAAAElFTkSuQmCC",
      "text/plain": [
       "<Figure size 640x480 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "#sns.kdeplot(hbliq_nums, bw_method=0.5)\n",
    "sns.kdeplot(hbih_nums, bw_method = 0.01)\n",
    "sns.kdeplot(hbic_nums, bw_method = 0.01)\n",
    "#plt.xlim([4.9, 5.1])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "9caa0d72",
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3 (ipykernel)",
   "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.10.9"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}