{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "3e316db481ea72fe",
   "metadata": {},
   "source": [
    "# Construction and Electrophysiological Simulation of Thalamic Neuron Models\n",
    "\n",
    "This example demonstrates how to implement thalamic neuron models using the `braincell` framework, including four typical thalamic neuron types: `HTC`, `RTC`, `IN`, and `TRN`.\n",
    "By constructing these neurons with complex combinations of ion channels, you will learn how to simulate the electrophysiological properties of thalamic neurons and understand the functional significance of differences in channel composition among neuron types.\n",
    "\n",
    "## Preparation\n",
    "First, ensure that the necessary libraries (`braincell`, `brainstate`, `brainunit`, `matplotlib`) are installed, and import the required modules:\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "id": "3127da63817f773a",
   "metadata": {
    "ExecuteTime": {
     "end_time": "2025-10-13T09:25:02.108097Z",
     "start_time": "2025-10-13T09:25:02.105610Z"
    },
    "execution": {
     "iopub.execute_input": "2026-05-25T10:28:13.691337Z",
     "iopub.status.busy": "2026-05-25T10:28:13.691045Z",
     "iopub.status.idle": "2026-05-25T10:28:17.192852Z",
     "shell.execute_reply": "2026-05-25T10:28:17.191545Z"
    }
   },
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "An NVIDIA GPU may be present on this machine, but a CUDA-enabled jaxlib is not installed. Falling back to cpu.\n"
     ]
    }
   ],
   "source": [
    "import time  # for timing the simulation\n",
    "\n",
    "import brainstate\n",
    "import braintools\n",
    "import brainunit as u\n",
    "import matplotlib.pyplot as plt\n",
    "import braincell"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "de625e637a78e3fa",
   "metadata": {},
   "source": [
    "## Code Explanation\n",
    "\n",
    "### Base Class for Thalamic Neurons\n",
    "\n",
    "Define a base class for all thalamic neurons, encapsulating common functionalities:\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "id": "442d98ad44653814",
   "metadata": {
    "ExecuteTime": {
     "end_time": "2025-10-13T09:25:02.121103Z",
     "start_time": "2025-10-13T09:25:02.117631Z"
    },
    "execution": {
     "iopub.execute_input": "2026-05-25T10:28:17.196477Z",
     "iopub.status.busy": "2026-05-25T10:28:17.195950Z",
     "iopub.status.idle": "2026-05-25T10:28:17.201866Z",
     "shell.execute_reply": "2026-05-25T10:28:17.200798Z"
    }
   },
   "outputs": [],
   "source": [
    "class ThalamusNeuron(braincell.SingleCompartment):\n",
    "    def compute_derivative(self, I_ext=0. * u.nA):\n",
    "        I_ext = self.sum_current_inputs(I_ext, self.V.value) * self.area_factor\n",
    "        for key, ch in self.nodes(braincell.IonChannel, allowed_hierarchy=(1, 1)).items():\n",
    "            I_ext = I_ext + ch.current(self.V.value)\n",
    "        self.V.derivative = I_ext / self.C\n",
    "        for key, node in self.nodes(braincell.IonChannel, allowed_hierarchy=(1, 1)).items():\n",
    "            node.compute_derivative(self.V.value)\n",
    "\n",
    "    def step_run(self, t, inp):\n",
    "        # Define update rule at each time step\n",
    "        with brainstate.environ.context(t=t):\n",
    "            self.update(inp)\n",
    "            return self.V.value\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "a8dfd148e1bd6ac0",
   "metadata": {},
   "source": [
    "### Implementation of Four Thalamic Neuron Types\n",
    "\n",
    "In the thalamus, different neuron types achieve distinct functions through specific combinations of ion channels. The following describes their channel composition:\n",
    "\n",
    "#### HTC\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "id": "d856d4866382d6dc",
   "metadata": {
    "ExecuteTime": {
     "end_time": "2025-10-13T09:25:02.135968Z",
     "start_time": "2025-10-13T09:25:02.130762Z"
    },
    "execution": {
     "iopub.execute_input": "2026-05-25T10:28:17.204608Z",
     "iopub.status.busy": "2026-05-25T10:28:17.204386Z",
     "iopub.status.idle": "2026-05-25T10:28:17.210058Z",
     "shell.execute_reply": "2026-05-25T10:28:17.209299Z"
    }
   },
   "outputs": [],
   "source": [
    "class HTC(ThalamusNeuron):\n",
    "    def __init__(\n",
    "        self,\n",
    "        size,\n",
    "        gKL=0.01 * (u.mS / u.cm **2),  # Potassium leak channel conductance\n",
    "        V_initializer=braintools.init.Constant(-65. * u.mV),  # Initial membrane potential\n",
    "        solver: str = 'ind_exp_euler'  # Integration method\n",
    "    ):\n",
    "        super().__init__(size, V_initializer=V_initializer, V_th=20. * u.mV, solver=solver)\n",
    "\n",
    "        # Membrane area parameter\n",
    "        self.area_factor = 1e-3 / (2.9e-4 * u.cm** 2)\n",
    "\n",
    "        # Sodium channel\n",
    "        self.na = braincell.ion.SodiumFixed(size, E=50. * u.mV)  # Sodium reversal potential 50 mV\n",
    "        self.na.add(INa=braincell.channel.Na_Ba2002(size, V_sh=-30 * u.mV))\n",
    "\n",
    "        # Potassium channel\n",
    "        self.k = braincell.ion.PotassiumFixed(size, E=-90. * u.mV)  # Potassium reversal potential -90 mV\n",
    "        self.k.add(IKL=braincell.channel.K_Leak(size, g_max=gKL))  # Potassium leak current\n",
    "        self.k.add(IDR=braincell.channel.KDR_Ba2002(size, V_sh=-30. * u.mV, q10=2.0, temp=u.celsius2kelvin(16.)))  # Delayed rectifier potassium current\n",
    "\n",
    "        # Calcium channel\n",
    "        self.ca = braincell.ion.CalciumDetailed(\n",
    "            size,\n",
    "            C_rest=5e-5 * u.mM,  # Resting calcium concentration\n",
    "            tau=10. * u.ms,  # Calcium decay time constant\n",
    "            d=0.5 * u.um  # Calcium diffusion distance\n",
    "        )\n",
    "        self.ca.add(ICaL=braincell.channel.CaL_IS2008(size, g_max=0.5 * (u.mS / u.cm **2)))  # L-type calcium channel\n",
    "        self.ca.add(ICaN=braincell.channel.CaN_IS2008(size, g_max=0.5 * (u.mS / u.cm** 2)))  # N-type calcium channel\n",
    "        self.ca.add(ICaT=braincell.channel.CaT_HM1992(size, g_max=2.1 * (u.mS / u.cm **2)))  # T-type calcium channel (low-threshold)\n",
    "        self.ca.add(ICaHT=braincell.channel.CaHT_HM1992(size, g_max=3.0 * (u.mS / u.cm** 2)))  # High-threshold calcium channel\n",
    "\n",
    "        # Calcium-activated potassium channel (IAHP)\n",
    "        self.kca = braincell.MixIons(self.k, self.ca)  # Mix potassium and calcium ions\n",
    "        self.kca.add(IAHP=braincell.channel.AHP_De1994(size, g_max=0.3 * (u.mS / u.cm **2)))\n",
    "\n",
    "        # Hyperpolarization-activated cation current (Ih) and leak current (IL)\n",
    "        self.Ih = braincell.channel.HCN_HM1992(size, g_max=0.01 * (u.mS / u.cm** 2), E=-43 * u.mV)  # Regulates resting potential and rhythm\n",
    "        self.IL = braincell.channel.IL(size, g_max=0.0075 * (u.mS / u.cm **2), E=-70 * u.mV)  # Background leak current\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "c519f384afa6ff61",
   "metadata": {},
   "source": [
    "#### RTC"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "id": "d1145a1a4b5adb37",
   "metadata": {
    "ExecuteTime": {
     "end_time": "2025-10-13T09:25:02.167383Z",
     "start_time": "2025-10-13T09:25:02.162904Z"
    },
    "execution": {
     "iopub.execute_input": "2026-05-25T10:28:17.212129Z",
     "iopub.status.busy": "2026-05-25T10:28:17.211723Z",
     "iopub.status.idle": "2026-05-25T10:28:17.217538Z",
     "shell.execute_reply": "2026-05-25T10:28:17.216756Z"
    }
   },
   "outputs": [],
   "source": [
    "class RTC(ThalamusNeuron):\n",
    "    def __init__(\n",
    "        self,\n",
    "        size,\n",
    "        gKL=0.01 * (u.mS / u.cm** 2),\n",
    "        V_initializer=braintools.init.Constant(-65. * u.mV),\n",
    "        solver: str = 'ind_exp_euler'\n",
    "    ):\n",
    "        super().__init__(size, V_initializer=V_initializer, V_th=20 * u.mV, solver=solver)\n",
    "\n",
    "        self.area_factor = 1e-3 / (2.9e-4 * u.cm **2)  # Membrane area parameter\n",
    "\n",
    "        # Sodium channel\n",
    "        self.na = braincell.ion.SodiumFixed(size)\n",
    "        self.na.add(INa=braincell.channel.Na_Ba2002(size, V_sh=-40 * u.mV))\n",
    "\n",
    "        # Potassium channel\n",
    "        self.k = braincell.ion.PotassiumFixed(size, E=-90. * u.mV)\n",
    "        self.k.add(IDR=braincell.channel.KDR_Ba2002(size, V_sh=-40 * u.mV, q10=2.0, temp=u.celsius2kelvin(16.)))  # Gating shift\n",
    "        self.k.add(IKL=braincell.channel.K_Leak(size, g_max=gKL))\n",
    "\n",
    "        # Calcium channel\n",
    "        self.ca = braincell.ion.CalciumDetailed(size, C_rest=5e-5 * u.mM, tau=10. * u.ms, d=0.5 * u.um)\n",
    "        self.ca.add(ICaL=braincell.channel.CaL_IS2008(size, g_max=0.3 * (u.mS / u.cm** 2)))\n",
    "        self.ca.add(ICaN=braincell.channel.CaN_IS2008(size, g_max=0.6 * (u.mS / u.cm **2)))\n",
    "        self.ca.add(ICaT=braincell.channel.CaT_HM1992(size, g_max=2.1 * (u.mS / u.cm** 2)))\n",
    "        self.ca.add(ICaHT=braincell.channel.CaHT_HM1992(size, g_max=0.6 * (u.mS / u.cm **2)))\n",
    "\n",
    "        # Calcium-activated potassium channel (IAHP)\n",
    "        self.kca = braincell.MixIons(self.k, self.ca)\n",
    "        self.kca.add(IAHP=braincell.channel.AHP_De1994(size, g_max=0.1 * (u.mS / u.cm** 2)))\n",
    "\n",
    "        # Ih and IL currents\n",
    "        self.Ih = braincell.channel.HCN_HM1992(size, g_max=0.01 * (u.mS / u.cm **2), E=-43 * u.mV)\n",
    "        self.IL = braincell.channel.IL(size, g_max=0.0075 * (u.mS / u.cm** 2), E=-70 * u.mV)\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "17c464d4af82fcbb",
   "metadata": {},
   "source": [
    "#### IN"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "id": "b8848c7317a31960",
   "metadata": {
    "ExecuteTime": {
     "end_time": "2025-10-13T09:25:02.181844Z",
     "start_time": "2025-10-13T09:25:02.177974Z"
    },
    "execution": {
     "iopub.execute_input": "2026-05-25T10:28:17.219624Z",
     "iopub.status.busy": "2026-05-25T10:28:17.219392Z",
     "iopub.status.idle": "2026-05-25T10:28:17.223901Z",
     "shell.execute_reply": "2026-05-25T10:28:17.223395Z"
    }
   },
   "outputs": [],
   "source": [
    "class IN(ThalamusNeuron):\n",
    "    def __init__(\n",
    "        self,\n",
    "        size,\n",
    "        gKL=0.01 * (u.mS / u.cm **2),\n",
    "        V_initializer=braintools.init.Constant(-70. * u.mV),  # More hyperpolarized initial membrane potential\n",
    "        solver: str = 'ind_exp_euler'\n",
    "    ):\n",
    "        super().__init__(size, V_initializer=V_initializer, V_th=20. * u.mV, solver=solver)\n",
    "\n",
    "        self.area_factor = 1e-3 / (1.7e-4 * u.cm** 2)\n",
    "\n",
    "        # Sodium channel\n",
    "        self.na = braincell.ion.SodiumFixed(size)\n",
    "        self.na.add(INa=braincell.channel.Na_Ba2002(size, V_sh=-30 * u.mV))\n",
    "\n",
    "        # Potassium channel\n",
    "        self.k = braincell.ion.PotassiumFixed(size, E=-90. * u.mV)\n",
    "        self.k.add(IDR=braincell.channel.KDR_Ba2002(size, V_sh=-30 * u.mV, q10=2.0, temp=u.celsius2kelvin(16.)))\n",
    "        self.k.add(IKL=braincell.channel.K_Leak(size, g_max=gKL))\n",
    "\n",
    "        # Calcium channel\n",
    "        self.ca = braincell.ion.CalciumDetailed(size, C_rest=5e-5 * u.mM, tau=10. * u.ms, d=0.5 * u.um)\n",
    "        self.ca.add(ICaN=braincell.channel.CaN_IS2008(size, g_max=0.1 * (u.mS / u.cm **2)))\n",
    "        self.ca.add(ICaHT=braincell.channel.CaHT_HM1992(size, g_max=2.5 * (u.mS / u.cm** 2)))\n",
    "\n",
    "        # Calcium-activated potassium channel (IAHP)\n",
    "        self.kca = braincell.MixIons(self.k, self.ca)\n",
    "        self.kca.add(IAHP=braincell.channel.AHP_De1994(size, g_max=0.2 * (u.mS / u.cm **2)))\n",
    "\n",
    "        # Leak (IL) and hyperpolarization-activated current (Ih)\n",
    "        self.IL = braincell.channel.IL(size, g_max=0.0075 * (u.mS / u.cm** 2), E=-60 * u.mV)\n",
    "        self.Ih = braincell.channel.HCN_HM1992(size, g_max=0.05 * (u.mS / u.cm **2), E=-43 * u.mV)\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "8de99a6ba6eced59",
   "metadata": {},
   "source": [
    "#### TRN"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "id": "4d459cf6562f41f3",
   "metadata": {
    "ExecuteTime": {
     "end_time": "2025-10-13T09:25:02.197006Z",
     "start_time": "2025-10-13T09:25:02.192835Z"
    },
    "execution": {
     "iopub.execute_input": "2026-05-25T10:28:17.225950Z",
     "iopub.status.busy": "2026-05-25T10:28:17.225806Z",
     "iopub.status.idle": "2026-05-25T10:28:17.231318Z",
     "shell.execute_reply": "2026-05-25T10:28:17.230335Z"
    }
   },
   "outputs": [],
   "source": [
    "class TRN(ThalamusNeuron):\n",
    "    def __init__(\n",
    "        self,\n",
    "        size,\n",
    "        gKL=0.01 * (u.mS / u.cm** 2),\n",
    "        V_initializer=braintools.init.Constant(-70. * u.mV),\n",
    "        gl=0.0075,  # Leak conductance coefficient\n",
    "        solver: str = 'ind_exp_euler'\n",
    "    ):\n",
    "        super().__init__(size, V_initializer=V_initializer, V_th=20. * u.mV, solver=solver)\n",
    "\n",
    "        self.area_factor = 1e-3 / (1.43e-4 * u.cm **2)\n",
    "\n",
    "        # Sodium channel\n",
    "        self.na = braincell.ion.SodiumFixed(size)\n",
    "        self.na.add(INa=braincell.channel.Na_Ba2002(size, V_sh=-40 * u.mV))\n",
    "\n",
    "        # Potassium channel\n",
    "        self.k = braincell.ion.PotassiumFixed(size, E=-90. * u.mV)\n",
    "        self.k.add(IDR=braincell.channel.KDR_Ba2002(size, V_sh=-40 * u.mV, q10=2.0, temp=u.celsius2kelvin(16.)))\n",
    "        self.k.add(IKL=braincell.channel.K_Leak(size, g_max=gKL))\n",
    "\n",
    "        # Calcium channel\n",
    "        self.ca = braincell.ion.CalciumDetailed(size, C_rest=5e-5 * u.mM, tau=100. * u.ms, d=0.5 * u.um)\n",
    "        self.ca.add(ICaN=braincell.channel.CaN_IS2008(size, g_max=0.2 * (u.mS / u.cm** 2)))\n",
    "        self.ca.add(ICaT=braincell.channel.CaT_HP1992(size, g_max=1.3 * (u.mS / u.cm **2)))\n",
    "\n",
    "        # Calcium-activated potassium channel (IAHP)\n",
    "        self.kca = braincell.MixIons(self.k, self.ca)\n",
    "        self.kca.add(IAHP=braincell.channel.AHP_De1994(size, g_max=0.2 * (u.mS / u.cm** 2)))\n",
    "\n",
    "        # Leak current (IL)\n",
    "        self.IL = braincell.channel.IL(size, g_max=gl * (u.mS / u.cm **2), E=-60 * u.mV)\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "b5405d5120eedd72",
   "metadata": {},
   "source": [
    "### Simulating Thalamic Neuron Electrical Activity\n",
    "\n",
    "Define the `try_trn_neuron` function to set input currents, simulation parameters, run the simulation, and visualize the results:\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "id": "5b1dff1ed3d6567c",
   "metadata": {
    "ExecuteTime": {
     "end_time": "2025-10-13T09:25:03.733730Z",
     "start_time": "2025-10-13T09:25:02.225284Z"
    },
    "execution": {
     "iopub.execute_input": "2026-05-25T10:28:17.234004Z",
     "iopub.status.busy": "2026-05-25T10:28:17.233654Z",
     "iopub.status.idle": "2026-05-25T10:28:22.058046Z",
     "shell.execute_reply": "2026-05-25T10:28:22.056997Z"
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Simulation elapsed time: 1.6359 seconds\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAj4AAAHHCAYAAAC/R1LgAAAAOnRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjEwLjgsIGh0dHBzOi8vbWF0cGxvdGxpYi5vcmcvwVt1zgAAAAlwSFlzAAAPYQAAD2EBqD+naQAAeIRJREFUeJzt3Xd4k1X/BvA7u7ulG2jpAGRvBBEBBaQ4UAS37KFiQRn6e0VfwI0T9XXxOsDxojjALQqyFEFENsimUCgdtKW7TZrk/P5IkyZt0iZp0iTN/bmuXm2f52mek3Rwc873nCMRQggQERER+QGppxtARERE1FwYfIiIiMhvMPgQERGR32DwISIiIr/B4ENERER+g8GHiIiI/AaDDxEREfkNBh8iIiLyGww+RERE5DcYfIiIiMhvMPiQV/nwww8hkUjw999/Wz1/9dVXo3v37gCAJ554AhKJpNG3q6++2vT1W7Zswbhx4xAfHw+lUonY2FiMGTMGa9eubY6nR35s+/bteOKJJ1BUVNTs93777bfx4YcfNvt9Wzq+rr5J7ukGEDlr3Lhx6NChg+nzsrIyzJo1C7fccgvGjRtnOh4XFwcAWLJkCZ566il07NgR9913H5KSklBQUICffvoJ48ePx6pVq3D33Xc3+/Mg/7B9+3Y8+eSTmDJlCiIiIpr13m+//Taio6MxZcqUZr1vS8fX1Tcx+JDP6tmzJ3r27Gn6PD8/H7NmzULPnj0xYcIEi2u/+uorPPXUU7j11lvx6aefQqFQmM498sgj+OWXX1BdXd1sbXeFiooKBAUFeez+5eXlCA4O9tj9iczZ+n3QarXQ6/VQKpUeaBV5Iw51kV9YtGgRIiMjsWLFCovQY5SWloYbb7yxwceQSCSYPXs2vvnmG3Tv3h0qlQrdunXDzz//XO/arKwsTJs2DXFxcabrVqxYYXGNcVjvzJkzFse3bNkCiUSCLVu2mI4Zh/h2796NoUOHIigoCI899hgAIC8vD9OnT0dcXBwCAgLQq1cvfPTRRxaPeebMGUgkErz88st499130b59e6hUKlx++eXYtWtXg8/bvK1bt27FAw88gNjYWCQkJJjOr1u3DkOGDEFwcDBCQ0Nxww034PDhwxaPkZOTg6lTpyIhIQEqlQqtW7fGzTffbPH8k5OTceONN2L9+vXo3bs3AgIC0LVrV6tDkadPn8Ztt92GyMhIBAUF4YorrsCPP/5o9bX84osv8OyzzyIhIQEBAQEYMWIETp48aXHtiRMnMH78eMTHxyMgIAAJCQm48847UVxcbHHd//73P/Tr1w+BgYGIjIzEnXfeiXPnzjX4+j3xxBN45JFHAAApKSmmYVjjc9dqtXj66adN35fk5GQ89thjUKvVDT6uPa9rcnIyDh8+jK1bt1od/i0qKsLcuXORmJgIlUqFDh064IUXXoBerzddY/7z8+qrryIpKQmBgYEYNmwYDh061GgbjfeZN28ekpOToVKpkJCQgEmTJiE/Px+Aa34fzNv52muvmV7Pf/75BwBw9OhR3HrrrYiMjERAQAD69++P7777zuJ+xnb88ccfmD9/PmJiYhAcHIxbbrkFFy9eNF3X2OtK3os9PuSViouLTX8QzTnTK3PixAkcPXoU06ZNQ2hoaJPatW3bNqxduxYPPPAAQkND8Z///Afjx49HZmYmoqKiAAC5ubm44oorTEEpJiYG69atw/Tp01FSUoK5c+c6de+CggJcd911uPPOOzFhwgTExcWhsrISV199NU6ePInZs2cjJSUFX375JaZMmYKioiI89NBDFo/x6aeforS0FPfddx8kEglefPFFjBs3DqdPn7YaCOt64IEHEBMTg8WLF6O8vBwA8Mknn2Dy5MlIS0vDCy+8gIqKCrzzzju46qqrsHfvXiQnJwMAxo8fj8OHD2POnDlITk5GXl4eNmzYgMzMTNM1gOH7dccdd+D+++/H5MmTsXLlStx22234+eefce2115pe4yuvvBIVFRV48MEHERUVhY8++gg33XQTvvrqK9xyyy0W7X7++echlUrx8MMPo7i4GC+++CLuuece7Ny5EwCg0WiQlpYGtVqNOXPmID4+HllZWfjhhx9QVFSE8PBwAMCzzz6LRYsW4fbbb8eMGTNw8eJFvPHGGxg6dCj27t1rcwhr3LhxOH78OD777DO8+uqriI6OBgDExMQAAGbMmIGPPvoIt956KxYsWICdO3di6dKlOHLkCL7++usGvyeNva6vvfYa5syZg5CQEDz++OMAaod/KyoqMGzYMGRlZeG+++5Du3btsH37dixcuBDZ2dl47bXXLO718ccfo7S0FOnp6aiqqsLrr7+O4cOH4+DBg6bHtKasrAxDhgzBkSNHMG3aNPTt2xf5+fn47rvvcP78edPr4Qhrvw9GK1euRFVVFe69916oVCpERkbi8OHDGDx4MNq2bYtHH30UwcHB+OKLLzB27FisWbOm3s/MnDlz0KpVKyxZsgRnzpzBa6+9htmzZ+Pzzz8HgAZfV/JygsiLrFy5UgBo8K1bt25Wv/bixYsCgFiyZInF8W+//VYAEK+++mqT2gZAKJVKcfLkSdOx/fv3CwDijTfeMB2bPn26aN26tcjPz7f4+jvvvFOEh4eLiooKi+eakZFhcd3mzZsFALF582bTsWHDhgkAYvny5RbXvvbaawKA+N///mc6ptFoxKBBg0RISIgoKSkRQgiRkZEhAIioqChRWFhoutb42nz//fcNPndjW6+66iqh1WpNx0tLS0VERISYOXOmxfU5OTkiPDzcdPzSpUsCgHjppZcavE9SUpIAINasWWM6VlxcLFq3bi369OljOjZ37lwBQPz+++8WbUlJSRHJyclCp9MJIWpfyy5dugi1Wm269vXXXxcAxMGDB4UQQuzdu1cAEF9++aXNtp05c0bIZDLx7LPPWhw/ePCgkMvl9Y7X9dJLL1n9fu/bt08AEDNmzLA4/vDDDwsAYtOmTTYf097XtVu3bmLYsGH1jj/99NMiODhYHD9+3OL4o48+KmQymcjMzBRC1P78BAYGivPnz5uu27lzpwAg5s2b1+D9Fy9eLACItWvX1jun1+uFEK75fTC2MywsTOTl5VmcGzFihOjRo4eoqqqyuPeVV14pOnbsaDpmbMfIkSNNbRNCiHnz5gmZTCaKiopMx2y9ruTdONRFXumtt97Chg0b6r2Z1/TYq6SkBACa3NsDACNHjkT79u1Nn/fs2RNhYWE4ffo0AEAIgTVr1mDMmDEQQiA/P9/0lpaWhuLiYuzZs8epe6tUKkydOtXi2E8//YT4+HjcddddpmMKhQIPPvggysrKsHXrVovr77jjDrRq1cr0+ZAhQwDA1P7GzJw5EzKZzPT5hg0bUFRUhLvuusviucpkMgwcOBCbN28GAAQGBkKpVGLLli24dOlSg/do06aNxf++w8LCMGnSJOzduxc5OTmm5z1gwABcddVVputCQkJw77334syZM6ahDaOpU6da1HjUfd7GHp1ffvkFFRUVVtu1du1a6PV63H777RbPNT4+Hh07djQ9V0f99NNPAID58+dbHF+wYAEA1Bu+M+fI62rNl19+iSFDhqBVq1YWz2nkyJHQ6XT47bffLK4fO3Ys2rZta/p8wIABGDhwoOk52LJmzRr06tWrXq8KYBhCdoa13wej8ePHm3rTAKCwsBCbNm3C7bffjtLSUtPzLCgoQFpaGk6cOIGsrCyLx7j33nst2jZkyBDodDqcPXvWqfaS9+BQF3mlAQMGoH///vWOG/9AOyIsLAwAUFpa2uR2tWvXzmqbjP/oXLx4EUVFRXj33Xfx7rvvWn2MvLw8p+7dtm3begWaZ8+eRceOHSGVWv4fpkuXLqbzDbXfGILs/UczJSXF4vMTJ04AAIYPH271euNrr1Kp8MILL2DBggWIi4vDFVdcgRtvvBGTJk1CfHy8xdd06NCh3j+Gl112GQBDrUl8fDzOnj2LgQMH1ruf+fM2LnsANP68U1JSMH/+fCxbtgyrVq3CkCFDcNNNN2HChAmmUHTixAkIIdCxY0erz9WeoUJrzp49C6lUajFDEQDi4+MRERHR4D+0jryu1pw4cQIHDhywCAnm6v6sWnvul112Gb744osG73Pq1CmMHz++0fY4wtrvg1Hdn9OTJ09CCIFFixZh0aJFVr8mLy/PItQ19XeFvBeDD7V4nTt3BgAcPHiwyY9l3tthTggBAKaC0AkTJmDy5MlWrzX2Wtn6n65Op7N6PDAw0KG2WtNY+xtTtw3G5/vJJ59Y/YdWLq/9EzN37lyMGTMG33zzDX755RcsWrQIS5cuxaZNm9CnTx97n4JT7Hner7zyCqZMmYJvv/0W69evx4MPPoilS5fizz//REJCAvR6PSQSCdatW2f18UJCQprURmd7Ppryuur1elx77bX4v//7P6vnjYGzObjy98HWz+nDDz+MtLQ0q19TN3g29XeFvBeDD7V4l112GTp16oRvv/0Wr7/+epP/gWpITEwMQkNDodPpMHLkyAavNf4Psu6Cdo50pSclJeHAgQPQ6/UWvT5Hjx41nXcn47BfbGxso8/XeP2CBQuwYMECnDhxAr1798Yrr7yC//3vf6ZrjP87N/+H8Pjx4wBgKoJOSkrCsWPH6j1+U593jx490KNHD/z73//G9u3bMXjwYCxfvhzPPPMM2rdvDyEEUlJSnAoEtv5hT0pKgl6vx4kTJ0w9VoChgLuoqMiu59LY62rr3u3bt0dZWZld3zugtofP3PHjxy2K023dp7HZX674fbAlNTUVgKFXzt7nag9nwyp5Fmt8yC88+eSTKCgowIwZM6DVauudX79+PX744Ycm30cmk2H8+PFYs2aN1T/05tNhjaHBvI5Cp9PZHCKz5vrrr0dOTo5ppglgmBr9xhtvICQkBMOGDXPmadgtLS0NYWFheO6556zOuDM+34qKClRVVVmca9++PUJDQ+tN2b5w4YLFTKaSkhJ8/PHH6N27t6lX6frrr8dff/2FHTt2mK4rLy/Hu+++i+TkZHTt2tWh51FSUlLv56JHjx6QSqWm9o0bNw4ymQxPPvlkvf/1CyFQUFDQ4D2Max7V/Yf9+uuvB4B6M6iWLVsGALjhhhtsPqa9r2twcLDVFaNvv/127NixA7/88ku9c0VFRfVek2+++caiFuavv/7Czp07cd1119lsI2Coudm/f7/VGWrG19IVvw+2xMbG4uqrr8Z///tfZGdn1ztv/nvpCFuvK3k39viQX7jjjjtw8OBBPPvss9i7dy/uuusu08rNP//8MzZu3IhPP/3UJfd6/vnnsXnzZgwcOBAzZ85E165dUVhYiD179uDXX39FYWEhAKBbt2644oorsHDhQhQWFiIyMhKrV6+2Gsxsuffee/Hf//4XU6ZMwe7du5GcnIyvvvoKf/zxB1577TWXFHQ3JCwsDO+88w4mTpyIvn374s4770RMTAwyMzPx448/YvDgwXjzzTdx/PhxjBgxArfffju6du0KuVyOr7/+Grm5ubjzzjstHvOyyy7D9OnTsWvXLsTFxWHFihXIzc3FypUrTdc8+uij+Oyzz3DdddfhwQcfRGRkJD766CNkZGRgzZo19WqeGrNp0ybMnj0bt912Gy677DJotVp88sknpiALGP5hfuaZZ7Bw4UKcOXMGY8eORWhoKDIyMvD111/j3nvvxcMPP2zzHv369QMAPP7447jzzjuhUCgwZswY9OrVC5MnT8a7776LoqIiDBs2DH/99Rc++ugjjB07Ftdcc43Nx7T3de3Xrx/eeecdPPPMM+jQoQNiY2MxfPhwPPLII/juu+9w4403YsqUKejXrx/Ky8tx8OBBfPXVVzhz5ozFVPMOHTrgqquuwqxZs6BWq/Haa68hKirK5lCZ0SOPPIKvvvoKt912G6ZNm4Z+/fqhsLAQ3333HZYvX45evXq55PehIW+99Rauuuoq9OjRAzNnzkRqaipyc3OxY8cOnD9/Hvv373f4MW29ruTlPDKXjMgG41TSXbt2WT0/bNgwh6ezm9u4caO4+eabRWxsrJDL5SImJkaMGTNGfPvtt422DYBIT0+vdzwpKUlMnjzZ4lhubq5IT08XiYmJQqFQiPj4eDFixAjx7rvvWlx36tQpMXLkSKFSqURcXJx47LHHxIYNG6xO37X1vHNzc8XUqVNFdHS0UCqVokePHmLlypUW1xin+Vqb9tzYayZE49+XzZs3i7S0NBEeHi4CAgJE+/btxZQpU8Tff/8thBAiPz9fpKeni86dO4vg4GARHh4uBg4cKL744guLx0lKShI33HCD+OWXX0TPnj2FSqUSnTt3tjrN/NSpU+LWW28VERERIiAgQAwYMED88MMP9doFK9PUja+H8XU6ffq0mDZtmmjfvr0ICAgQkZGR4pprrhG//vprvfuuWbNGXHXVVSI4OFgEBweLzp07i/T0dHHs2LEGX0MhDNPH27ZtK6RSqcXU7erqavHkk0+KlJQUoVAoRGJioli4cKHF1Gtr7H1dc3JyxA033CBCQ0MFAIsp2KWlpWLhwoWiQ4cOQqlUiujoaHHllVeKl19+WWg0GovX66WXXhKvvPKKSExMFCqVSgwZMkTs37+/0ecthBAFBQVi9uzZom3btkKpVIqEhAQxefJki2Ufmvr70NDPufHxJ02aJOLj44VCoRBt27YVN954o/jqq69M19j6Wbc2rb6h15W8l0QIVmoRkXdITk5G9+7dXTLsSK5z5swZpKSk4KWXXmqwV4vIF7DGh4iIiPwGgw8RERH5DQYfIiIi8hus8SEiIiK/wR4fIiIi8hsMPkREROQ3uIBhHXq9HhcuXEBoaCiXIyciIvIRQgiUlpaiTZs2DS5iyuBTx4ULF5CYmOjpZhAREZETzp07h4SEBJvnGXzqMC7xf+7cOYSFhXm4NURERGSPkpISJCYmNrpVD4NPHcbhrbCwMAYfIiIiH9NYmQqLm4mIiMhvMPgQERGR32DwISIiIr/B4ENERER+g8GHiIiI/AaDDxEREfkNBh8iIiLyGww+RERE5DcYfIiIiMhvMPgQERGR32DwISIiIr/B4ENERER+g8GHyMOqqnXQ64Wnm0FE5BcYfIg8qKhCgyuWbsSED3Z6uilERH6BwYfIg/7KKERRRTW2nypAmVrr6eYQEbV4DD5EHpRVVGn6uLiy2oMtISLyDww+RB6UkV9u+risij0+RETuxuBD5EFZl2p7fDjURUTkfgw+RB5kPtRVrdN7sCWO0+sFLpaqPd0MIiKHMPgQeVCp2fCWRuu64COE+6fHP/XDP7j82V/x6c5Mt9+LiMhVGHyIPMg8oLiqx2fah7sw/JWtKKlyb7H0h9vPAACe/fEft96HiMiVGHyIvIQrenw0Wj02Hc1DRn45/jiR74JWNa5co2uW+xARuQKDD5EHmQ9IaVzQ46M360G6WObe+huZVOLWxycicgcGHyIv4coaHwAocfO6QKEBcrc+PhGROzD4EHmQeQ2yK3p8zB+vxM3rAgUra4NPcxRTExG5AoMPkZeodnGPT6mbg49SXvvnQ8tNVonIRzD4EHmQMKvy0bkgO5g/nrvXBVLIamt8XD1MR0TkLgw+RF5Cp7cdHkrtnJpuPuLk/uBT++eDwYeIfAWDD5EHmQcVWzllze7z6PHEeqz+y7GFAt0dRqQSsx4fH1t1moj8F4MPkQeZj27Z6vFZ8OV+AMCjaw869HjuDj7mw2rs8SEiX8HgQ+QlXNFpYj67yt29MOa9VWoGHyLyEQw+RB5kMdTl4inhzRlG2ONDRL6CwYfISzRU3GyvZh3qcvEaREREzYHBh8ijzKaz28gOEgd2hmjOWV321CcREXkbBh8iL2ErPChlzv2a6ty8qKB5PRE7fIjIVzD4EHmQPdPZ5Y5sBmr2ePpm3EZCyx4fIvIRDD5EXsJWj4/UgeBjsRJ0M24j0Zz3IiJqCgYfIg+yqJOx0UMjc6THx4y7s4hlbxWDDxH5BgYfIi9hKzzIHKhuNg8j7h5+8lTvEhFRUzD4EHmQZYGw9fDg2FBXreYsu+Hu7ETkKxh8iDzIPC7YCg+O9PiYc3dxs/nD6xl8iMhHMPgQeQlb4cGRGh/zHiR398LYE9qIiLwNgw+RB1nW5Nga6nLg8cw+bs5eGNb4EJGvYPAh8hK2hqacHepye4+PHfVJRETexmeDz/PPPw+JRIK5c+eajlVVVSE9PR1RUVEICQnB+PHjkZub67lGEjXCYmhK54Li5macYm65ZQWDDxH5Bp8MPrt27cJ///tf9OzZ0+L4vHnz8P333+PLL7/E1q1bceHCBYwbN85DrSRyjCt6fJp1irkdw3RERN7G54JPWVkZ7rnnHrz33nto1aqV6XhxcTE++OADLFu2DMOHD0e/fv2wcuVKbN++HX/++acHW0xkm12zupxcwLBZV25uxu0xiIiawueCT3p6Om644QaMHDnS4vju3btRXV1tcbxz585o164dduzYYfPx1Go1SkpKLN6IPMHmOj4Obc9e+6H7FzCspeMupUTkI+SeboAjVq9ejT179mDXrl31zuXk5ECpVCIiIsLieFxcHHJycmw+5tKlS/Hkk0+6uqlE9rFjU9GmbFmh1wuHaoQc0ZxT54mIXMVnenzOnTuHhx56CKtWrUJAQIDLHnfhwoUoLi42vZ07d85lj03UGHsKhJ1duRlw7xCUxdR5DnURkY/wmeCze/du5OXloW/fvpDL5ZDL5di6dSv+85//QC6XIy4uDhqNBkVFRRZfl5ubi/j4eJuPq1KpEBYWZvFG5Am2RqZkjox01ckfzVXnwx4fIvIVPjPUNWLECBw8eNDi2NSpU9G5c2f861//QmJiIhQKBTZu3Ijx48cDAI4dO4bMzEwMGjTIE00mapTFWjguHuoC3Bt8LKbO25iKT0TkbXwm+ISGhqJ79+4Wx4KDgxEVFWU6Pn36dMyfPx+RkZEICwvDnDlzMGjQIFxxxRWeaDKRQ1xR3CzqDHa5syfG/F7s8SEiX+Ezwccer776KqRSKcaPHw+1Wo20tDS8/fbbnm4WkU321MmYBx8hBCQNBKG6DyHcWeMjzD9m8CEi3+DTwWfLli0WnwcEBOCtt97CW2+95ZkGETWBrR4f85yjF47V/DRXjQ87fIjIV/hMcTNRS2TPFhPmPT6NBZm6Z90ZSMzbzlldROQrGHyIvISt8GDZ49NI8KlzvrmGoLhyMxH5CgYfIg+yZ28tiUWNj2OP31yBhLmHiHwFgw+RB1kOF1m/xrykp/EeH8vP3TvUVfvgehb5EJGPYPAh8hL2FTc7FjDcGUgsZ6S57TZERC7F4EPkQXZtWWGWfBzdd9SdRccsbiYiX8TgQ+QlbK/j0/g1Rs051GV5HwYfIvINDD5EnmTHdHbzKp9Gg0+dCe1u7fExuxeDDxH5CgYfIi9hOzyYBwwHH7O59upycAiOiMhTGHyIPMie6ewW13vRUJf5Q3PLCiLyFQw+RF7CVvCx6FlpdKjLUnMNQXGoi4h8BYMPkQfZs45PU6aNu3OvLg51EZEvYvAh8hL2hJTGanbqb1nRpCY1ovbBOdRFRL6CwYfIgyzW8bERHsxDRWP5ojmHuriODxH5IgYfIg+yZ9sHy6EuB1dubrZNSpvlNkRETcbgQ+QlbPf4NH6NtWsBd6/j0zz3ISJyJQYfIg+ynBJuz3T1xgJG3QUMnWuXPSyH4Bh8iMg3MPgQeRFrBc5NmdXVbJuUclYXEfkIBh8iD6rbUWJtKMuiDsjBoa7GhsZcpbnuQ0TUVAw+RF6ksZ4TR3tW3JlHzB+bQ11E5CsYfIi8SGM9J41vUurY9U1h2RPlttsQEbkUgw+RF7Fa4+PAejn1hrqaKZE0132IiJqKwYfIQ6wNDzVWjOxovnDrUJfZx5zOTkS+gsGHyEOsZQWrxc1woLi53nT25kk+zD1E5CsYfIi8iLUeH0eKiJtzqMtiuw0OdRGRj2DwIfIQa1HB+nR2s/MOzupqrjzCoS4i8hUMPkRepLGeE0eLm905zdyRzVOJiLwFgw+Rh1gvbrZyXRNqfNy5sCCLm4nIFzH4EHmRxoa6HM0X7t2rq/ZjrtxMRL6CwYfIQ6zW+PjQUJc51jYTka9g8CHyAnKpBID1YNOUTUrdO6uLu7MTke9h8CHyEPOsIKsJPlaDivnKzQ4GmWYb6mKXDxH5CAYfIi/QUPBxqLi5zmm37tVlcR+33YaIyKUYfIg8xDzQyBoY6jLXWMCot3JzMyUSDnURka9g8CHyEHuHuhzZpLQut+YeDnURkQ9i8CHyAvYWNzu6ZYV7h7rsH4IjIvIWDD5EXkAqMfb4NHxdY+frxg+3Bh9uUkpEPojBh8gLyBsc6nK+Z6W5anzY40NEvoLBh8hDzLOC1O51fLynxsdid3YGHyLyEQw+RF6g4R4f6x9bU7cGyL1DXWY9UQ7uGk9E5ClyR7/gyJEjWL16NX7//XecPXsWFRUViImJQZ8+fZCWlobx48dDpVK5o61ELYq16eyN9Zw43uPD6exERObs7vHZs2cPRo4ciT59+mDbtm0YOHAg5s6di6effhoTJkyAEAKPP/442rRpgxdeeAFqtdqd7SZqUUzr+FhdwLBWY9PG6xc3N7Fhdt6LQ11E5Cvs7vEZP348Hn74YXz11VeIiIiwed2OHTvw+uuv45VXXsFjjz3mijYStUiW6/gY/g9ifcsKYe1DuzTXrC4u40NEvsLu4HP8+HEoFIpGrxs0aBAGDRqE6urqJjWMqKUzzwqymr7XxldudnAdH67cTERkwe6hLoVCgTfffBNFRUV2X09E9qnt8al/ril7YjVXTwxXbiYiX+HQrC5jDc/dd9+NTZs2uatNVi1duhSXX345QkNDERsbi7Fjx+LYsWMW11RVVSE9PR1RUVEICQnB+PHjkZub26ztJLKXeS+JzFDiY7VWxrEtKyzPuyuQ1J895pbbEBG5nEPBJycnB8uXL0d2djauvfZapKSk4Omnn8a5c+fc1T6TrVu3Ij09HX/++Sc2bNiA6upqjBo1CuXl5aZr5s2bh++//x5ffvkltm7digsXLmDcuHFubxtRU8lrenwaG5pytGbHXUNQzbk1BhGRKzkUfAIDAzFp0iRs3rwZJ06cwMSJE/HBBx8gJSUFo0ePxpdffum22p6ff/4ZU6ZMQbdu3dCrVy98+OGHyMzMxO7duwEAxcXF+OCDD7Bs2TIMHz4c/fr1w8qVK7F9+3b8+eefbmkTUVNY1vg0sI6P+Z5Yjc3qqnPaXbOt6j4qcw8R+QqnFzBMTU3FU089hYyMDKxbtw5RUVGYMmUK2rZt68r22VRcXAwAiIyMBADs3r0b1dXVGDlypOmazp07o127dtixY4fNx1Gr1SgpKbF4I2puDa3j05TZU6zxISKy1OSVmyUSCeRyOSQSCYQQzTKbS6/XY+7cuRg8eDC6d+8OwDAMp1Qq6021j4uLQ05Ojs3HWrp0KcLDw01viYmJ7mw6kYnldPYG1vFxqMbHkrtmdTXnCtFERK7kdPA5d+4cnnrqKaSmpuLaa6/FhQsX8N577yE7O9uV7bMqPT0dhw4dwurVq5v8WAsXLkRxcbHprTnqlYjqctXKzc21OzuHuojIVzm0ZYVGo8HatWuxYsUKbNq0Ca1bt8bkyZMxbdo0pKamuquNFmbPno0ffvgBv/32GxISEkzH4+PjodFoUFRUZNHrk5ubi/j4eJuPp1KpuMUGeYa9PT5mHzvagWNterwrNFctERGRqzkUfOLj41FRUYEbb7wR33//PdLS0iCVNs8+p0IIzJkzB19//TW2bNmClJQUi/P9+vWDQqHAxo0bMX78eADAsWPHkJmZiUGDBjVLG4kcYbFXl8R2cbM5hxcwbKZAwqEuIvIVDgWff//735g4cSJiYmLc1R6b0tPT8emnn+Lbb79FaGioqW4nPDwcgYGBCA8Px/Tp0zF//nxERkYiLCwMc+bMwaBBg3DFFVc0e3uJHCGTGYe66p+z3AXdOzYpFXUGu4QwtFNSE+CIiLyVQ8Fn/vz59Y6VlZVBr7fsTw8LC2taq6x45513AABXX321xfGVK1diypQpAIBXX30VUqkU48ePh1qtRlpaGt5++22Xt4XIFSyKmyW2h7rMNZZ7mqvo2NrD6vQCchmDDxF5N4eCj1FGRgZmz56NLVu2oKqqynTc+D8+nU7nsgaaP3ZjAgIC8NZbb+Gtt95y+f2J3MlVxc11uavGxxrOaCciX+BU8JkwYQKEEFixYgXi4uLYvU3kBLsXMDSfzt7YAob1vrb50gjrfIjIFzgVfPbv34/du3ejU6dOrm4PkV9qaKjLYuVmh2d1Ne9QFxGRt3NqStbll1/O9W6Imshik1KZi9bxqTery7m2NaZucbPhXgw+ROT9nOrxef/993H//fcjKysL3bt3h0KhsDjfs2dPlzSOyF802ONjdsjR9XKaM4zom7GeiIjIWU4Fn4sXL+LUqVOYOnWq6Zhxywp3FTcTtTRWa3ys7dVl/rHDe3U141AXe3yIyAc4FXymTZuGPn364LPPPmNxM5ELSE0LGDZ8XePFzZbn3VbjY+UYh7qIyBc4FXzOnj2L7777Dh06dHB1e4j8huUmpYb31sKDxQKGDmYLd2URa7PF3LUhKhGRKzlV3Dx8+HDs37/f1W0h8ivG3hmJBJA2sFeXuUZ7VeruodWMYYRDXUTkC5zq8RkzZgzmzZuHgwcPokePHvWKm2+66SaXNI7IX5j26mqkxsfR4aTm2J1dIjH0LHE6OxH5AqeCz/333w8AeOqpp+qdY3EzkZ1qcoIEDe/Obp4yGp3OXufz5ihulkslqNYJtw2rERG5klPBp+7eXETkPIlEUlvc3Eh6cHQLCrd1wljUJxmCD3t8iMgXOFXjQ0RNZ33Lioava2wLirqnmyOMNDRMR0TkbewOPqtXr7b7Qc+dO4c//vjDqQYR+ZvGhrosZ3V5S41P7eMaC7Obc18wIiJn2R183nnnHXTp0gUvvvgijhw5Uu98cXExfvrpJ9x9993o27cvCgoKXNpQopbGPCfYO9TVWAdO3XV8mqPGp6HeKiIib2N3jc/WrVvx3Xff4Y033sDChQsRHByMuLg4BAQE4NKlS8jJyUF0dDSmTJmCQ4cOIS4uzp3tJmoxJBKzdXysblJay9G1cppjAUN5AzvLExF5G4eKm2+66SbcdNNNyM/Px7Zt23D27FlUVlYiOjoaffr0QZ8+fSCVsmyIyB4Ww0V29/g4OtTleLscZWw7V24mIl/g1Kyu6OhojB071sVNIfIvxpwggcRsuKjhTUobHeqquzu7u3p8zHeWlzL4EJHvYPcMkadJGg4P5j1D3rI7u/mj1u4zxuBDRN6PwYfIQ5wJD41OZ6/zua4ZskhtaHP/vYiImorBh8hDjCHGfDq71XV8zIe6HJw55a4p5tZmdXGoi4h8AYMPkYdJJLWLAFrfnb3240a3rKhz3n2zumo3WK1pOoe6iMgnMPgQeYjFOj52Tgn3mlld5vuMcVYXEfkQu2d1zZ8/3+4HXbZsmVONIfJHhlldho8bCw+NL2BY53o398JIJBKzVafdeisiIpewO/js3bvXruskxn5vIrKLRNJwcbN3bllRS8K9uojIh9gdfDZv3uzOdhD5HevbPjS2O7tjXT7uCiPCfKjLzt4qIiJvwBofIg8xFQijkeJm848dzBbuyiLmxc2mtrO4mYh8gFMrNwPA33//jS+++AKZmZnQaDQW59auXdvkhhH5C4lE4rbiZnfPtJJAUjvUxeBDRD7AqR6f1atX48orr8SRI0fw9ddfo7q6GocPH8amTZsQHh7u6jYStUgWQ12mOpmGr2t0gUMP7s7O3ENEvsCp4PPcc8/h1Vdfxffffw+lUonXX38dR48exe2334527dq5uo1ELZIxJ5gvYGh9d/baY/bmGOMcA7ft1WW6EaezE5FvcSr4nDp1CjfccAMAQKlUory8HBKJBPPmzcO7777r0gYStXgS163jYzxdG0aa3jzr96mtT+IChkTkS5wKPq1atUJpaSkAoG3btjh06BAAoKioCBUVFa5rHVELZrHDuYtWbjYyBSk398JIGtlglYjI2zhV3Dx06FBs2LABPXr0wG233YaHHnoImzZtwoYNGzBixAhXt5GoRZMAkNb8F8TqOj5mHze26Wi9Hh93DXVxry4i8lFOBZ8333wTVVVVAIDHH38cCoUC27dvx/jx4/Hvf//bpQ0kaqnMY4LMzkUA7d10tLnCiOWsLrfeiojIJZwKPpGRkaaPpVIpHn30UZc1iMhfmBYBtNj2oYmblNa8l7q57qa27YDMWEjNHh8i8gF2B5+SkhKEhYWZPm6I8ToiapxEYn9Njr29KsYg5e4s0tiMNCIib2N38GnVqhWys7MRGxuLiIgIq3tyCSEgkUig0+lc2kiilslKcbPVYGM+nd2xoS63bVlh1ibu1UVEvsTu4LNp0ybTEBf37SJqOsv9rhrapLT248ansxvOS928to7FMJ2bp84TEbmS3cFn2LBhpo9TUlKQmJhYr9dHCIFz5865rnVEfkAikdTuzt5IULE3XNQOPzWpaY3iUBcR+Rqn1vFJSUnBxYsX6x0vLCxESkpKkxtF5A8sZnU1uHJzrcbCRW1xs5t7fMw+5gKGRORLnAo+xlqeusrKyhAQENDkRhH5A8uhLsPH1np8zOt67A0ybq/xMWs81/EhIl/i0HT2+fPnAzB0zS9atAhBQUGmczqdDjt37kTv3r1d2kCilk4iqe2haXzLioYfy7SAodmsLlv/UWkKi33GuFcXEfkQh4LP3r17ARj+kB48eBBKpdJ0TqlUolevXnj44Ydd20KiFsp8ZpRSbujy0WjrF+VYrNxs53CS1Czn6EXtWjuuJpFwAUMi8i0OBR/jbK6pU6fi9ddf53o9RE1Q20EiQYBCBgBQa/X1emjMO1Ian85uOC8zSz56ISCDi3t8LLasqL0PEZG3c2rl5pUrV7q6HUR+SyKBKfgAhvBj/rk5e+uHpWbBSacXsPFwTVCzO7uEs7qIyLc4VdxcXl6ORYsW4corr0SHDh2Qmppq8eZpb731FpKTkxEQEICBAwfir7/+8nSTiOox7yAJkNf+KlZV6+pcZ39xc90an7r3cRXzwmwuYEhEvsSpHp8ZM2Zg69atmDhxIlq3bu3ywsmm+PzzzzF//nwsX74cAwcOxGuvvYa0tDQcO3YMsbGxnm4ekYmxxkcCQC6TQi6VQKsXqKzWIcLG1zg6qwtwTyAxFTebL2DIHh8i8gFOBZ9169bhxx9/xODBg13dniZbtmwZZs6cialTpwIAli9fjh9//BErVqzgZqrklYz/bwhQyFCm1qKq2rJK2GIdHyeGutxRe2Pe42PMWMw9ROQLnBrqatWqlcUO7d5Co9Fg9+7dGDlypOmYVCrFyJEjsWPHDqtfo1arUVJSYvFG1Bzq5pEAheHXse5Ql3nysXd3doviZjckEmFW42PvBqtERN7AqeDz9NNPY/HixaioqHB1e5okPz8fOp0OcXFxFsfj4uKQk5Nj9WuWLl2K8PBw01tiYmJzNJXIRFIz40olN1Qg1ws+Zuydzi6z6PFpQuNsMJ+RxqEuIvIlTg11vfLKKzh16hTi4uKQnJwMhUJhcX7Pnj0uaVxzWLhwoWlhRgAoKSlh+KFmVTvUZezxsT3U1Vinium8WdmdO7aSMN5HKml4g1UiIm/jVPAZO3asi5vhGtHR0ZDJZMjNzbU4npubi/j4eKtfo1KpoFKpmqN5RBbqD3XV9PhonZ/VZWTcPFSnF26p8TE+pvl0dg51EZEvcCr4LFmyxNXtcAmlUol+/fph48aNpnCm1+uxceNGzJ4927ONI6rDfFYXUBt81E0Y6jJfDdoYfNzZEyOBhD0+RORTnKrxAYCioiK8//77WLhwIQoLCwEYhriysrJc1jhnzJ8/H++99x4++ugjHDlyBLNmzUJ5eblplheRtzEuB+GKoa7ax6yt83HnUJcj+4wREXkDp3p8Dhw4gJEjRyI8PBxnzpzBzJkzERkZibVr1yIzMxMff/yxq9tptzvuuAMXL17E4sWLkZOTg969e+Pnn3+uV/BM5Gl1Q0ygwvDrWK7R2vwaexcwBNxbe2PeW8Xd2YnIlzjV4zN//nxMmTIFJ06cQEBAgOn49ddfj99++81ljXPW7NmzcfbsWajVauzcuRMDBw70dJOIGhUWUBN81JbBxzxP2DudXSqRmNbXccsChqYeHw51EZFvcSr47Nq1C/fdd1+9423btrU5bZyILNWNCSE1wae0qk7wgXlxc8MblQorRcfuWcenlvE+WgYfIvIBTgUflUpldaG/48ePIyYmpsmNIvIH5iEFAEJtBZ86eaKhDpzaFZUlbp1tZRGwuI4PEfkQp4LPTTfdhKeeegrV1dUADN3dmZmZ+Ne//oXx48e7tIFELZ0x+ISoDOthNRZ8GhruMp9m7s6i49q9usxXbnb5bYiIXM6p4PPKK6+grKwMsbGxqKysxLBhw9ChQweEhobi2WefdXUbiVqkujnB2ONTpq62OF436DTUg2NeeyN3Z3Gzec+Sca8u9vgQkQ9walZXeHg4NmzYgD/++AP79+9HWVkZ+vbta7FHFhE1zDw8ALaHuuoGH73lbHer10rNe2LcEkjMhrpkUjfeh4jItZwKPh9//DHuuOMODB482GKHdo1Gg9WrV2PSpEkuayBRS1e3xqdMXTf4WF7fYI+P8THh3mnmtVtW1O7VxeJmIvIFTg11TZ06FcXFxfWOl5aWcqFAIrtZBgVjjU+ZWY+PtRlcDfWsCFOPj8RsAcMmN7Qevam3Cqjp8OE6PkTkE5wKPkII02qz5s6fP4/w8PAmN4rIHwiz8ADU9viUmAUfaxmnoVoaixWV3VrjU9t4rtxMRL7EoaGuPn36QCKRQCKRYMSIEZDLa79cp9MhIyMDo0ePdnkjiVoy438iQlT1i5ut9aI0NKSkNytuNk0zd8dQV817rtxMRL7GoeBj3Phz3759SEtLQ0hIiOmcUqlEcnIyp7MT2aluTAgLMAx1VVXrUa3TQyGTWg0TDQUM860kpG5cWJArNxORr3Io+Bh3ZU9OTsYdd9xhsV0FETmm7lBXsEpmOleu1iIiSGl1scKGAobxlNRsOrt7Vm6uv1cXi5uJyBc4Natr8uTJAIDdu3fjyJEjAIBu3bqhT58+rmsZUQtXd+VmuUwKlVwKtVaPMieDD8wXMHRnT4xZLRFXbiYiX+JU8MnLy8Odd96JLVu2ICIiAgBQVFSEa665BqtXr+a2FUR2MO+dMQpRyaHWalCu1tVcUxsm5FIJtHphd4+PcWFBt2xZUfNeAonZys0MPkTk/Zya1TVnzhyUlpbi8OHDKCwsRGFhIQ4dOoSSkhI8+OCDrm4jUYtkPvXcKKhmuMu4lo9F8JE1HjDMZ1u5dZNSsx4fdw6pERG5mlM9Pj///DN+/fVXdOnSxXSsa9eueOuttzBq1CiXNY6oJdObhQejYKXhV7JCo7W4BgDkUikAfYMBw7zHxzTN3C09PrWP6c4iaiIiV3Oqx0ev10OhUNQ7rlAooG9oPX0iMtFb6fExTmkvr+nxMV/AUCFrPGBYm2bu1r26LBZKZPAhIu/nVPAZPnw4HnroIVy4cMF0LCsrC/PmzcOIESNc1jiilswUfMx+C4NMa/kYa3xqz8nt2BOrdvjMzcGn5r35fbiODxH5AqeCz5tvvomSkhIkJyejffv2aN++PVJSUlBSUoI33njD1W0kapGs9/gYanxqh7rMenzsCBjmPTHuXFFZbz57jD0+RORDnKrxSUxMxJ49e7Bx40bTdPYuXbpwd3YiBxhHhc23fwlSWm5Uah4wZLLGA4b59XJ39sSY7SxvLLpm7iEiX+Bw8Pn888/x3XffQaPRYMSIEZgzZ4472kXU4unNhqWM6tf4oOYa+2pprE4zd8smpbVtN/b4aFnfR0Q+wKHg88477yA9PR0dO3ZEYGAg1q5di1OnTuGll15yV/uIWixr6/gEKQ1DXXXX8ZHauSCh+fUyN87qMrZBKpWYTZt3+W2IiFzOoRqfN998E0uWLMGxY8ewb98+fPTRR3j77bfd1TaiFk1Y6fFRyQ3BR1PTTaM3G1KyJ8iYr6/jznV8jAFLxlldRORjHAo+p0+fNm1XAQB33303tFotsrOzXd4wopbOfCd1I5XC8CuprjYEH53OcJFcZl/PivmiiO5cX8c4fGbe48OVm4nIFzgUfNRqNYKDg2u/WCqFUqlEZWWlyxtG1NJZq/FRyWuCj9Yw1FVdk3JkZgGjoVoa80URjVtWuKPHR2fe48Pd2YnIhzhc3Lxo0SIEBQWZPtdoNHj22WcRHh5uOrZs2TLXtI6oBbM2nd001KU1hBttTY+PQia1a72c2lPu3UPLGKYMgcxwjMGHiHyBQ8Fn6NChOHbsmMWxK6+8EqdPnzZ9bt5tT0S2CSvFzUpTj09N8Knp3ZFLzdflsf2YxpAjl0qgc2NPjHlxs7Fd3KuLiHyBQ8Fny5YtbmoGkf+pXbnZvMfHcqjL2OMjl0pM6/LoGhjq0lkZGnPLUJexx8d8hWjW+BCRD3Bq5WYiarra6ey1x1T1enyMxc1Su9blMV0vde8mpaYaH6mUNT5E5FMYfIg8xGqNj6JujU/tUJc909mNs8BkMvcWHZt6fKTu3ROMiMjVGHyIPMTaOj5KmWWPj87U42O2NUQDAcNqj4871/GxM5AREXkLBh8iD9FZ2avLtI6PaTq7MchIzbaGaKDHR29lCMqdKzebrRckRG2YIyLyVgw+RB7S0Do+GlOPT81Ql8y+YmXzHp9mKW426/ExP05E5K2cDj6///47JkyYgEGDBiErKwsA8Mknn2Dbtm0uaxxRSyasruNjOdRVrXOsWNlyfR33b1Iqk0hMu8YD7lklmojIlZwKPmvWrEFaWhoCAwOxd+9eqNVqAEBxcTGee+45lzaQqKWytkmpcQFD05YVZrO65HYUEVv0+BjX13HLUJfhvbROj4877kVE5EpOBZ9nnnkGy5cvx3vvvQeFQmE6PnjwYOzZs8dljSNqyYwhwXzNT+MChsZNSqvNZ3XZEXzM1/GxZzd3ZxnbpTAbgnPXvYiIXMmp4HPs2DEMHTq03vHw8HAUFRU1tU1EfsFaj49xVpdOL6DV6U09PwEKmV1BxlqPjzuGn4w1SEqZ1KL9DW2gSkTkDZwKPvHx8Th58mS949u2bUNqamqTG0XkD/R6K8XNitpfSY1Oj6qa2V0BCmntpqMNreNjrPGRSU17aLmjuNnYI6WU1w7BAZzSTkTez6ngM3PmTDz00EPYuXMnJBIJLly4gFWrVuHhhx/GrFmzXN1GohbJfFVmI6XZxxqtHlXVNcFHLoNMKrX4ugYfU+reTUpNPT5yqcWWGw3tHE9E5A0c3p0dAB599FHo9XqMGDECFRUVGDp0KFQqFR5++GHMmTPH1W0kapHMV2U2ksukkEoMw2CG4GO4RqWQmep3Ghrqqray0rM7enyMs86MxdgyqQQ6veBQFxF5PaeCj0QiweOPP45HHnkEJ0+eRFlZGbp27YqQkBBXt4+oxart8ZFYHFfJZais1kFt3uOjkKKq2nB9Q0HGvCbInQsYmvf4AIZp7ToIDnURkddzKvgYKZVKdO3a1VVtIfIrpqnqUssRZ6VcahZ8aoOMsTenoXBhrAlSyd27eaipxqdmaE4qBaBzT+8SEZErORV8ysvL8fzzz2Pjxo3Iy8uDvk7/9unTp13SOKKWrFpff6gLMJvSrjUrbpbLUCHRAmg4yKjNhsbcGnxq2mVsqyG86bmAIRF5PaeCz4wZM7B161ZMnDgRrVu3tthriIjsY76TurnajUp1FkNd9kxnV5v1+Lhzk9K6Q13G7MZ1fIjI2zkVfNatW4cff/wRgwcPdnV7iPyG+Qwsc8Yp7Rqt5To+9uyCbiw6Nq/xccdqyqaia2ONjxvvRUTkSk5NZ2/VqhUiIyNd3RYiv6I1DXXVqfGR1a7ebN7jY+wZaqiOxjyQyNzY41OuMQy7hagM/3dy57AaEZErORV8nn76aSxevBgVFRWubo9VZ86cwfTp05GSkoLAwEC0b98eS5YsgUajsbjuwIEDGDJkCAICApCYmIgXX3yxWdpH5AytzkaPj3Gj0mo9StWGgBGsktu1ErPFUJdpVpdr2w0AZVW17QIYfIjIdzg11PXKK6/g1KlTiIuLQ3JyssV+XQBcvl/X0aNHodfr8d///hcdOnTAoUOHMHPmTJSXl+Pll18GAJSUlGDUqFEYOXIkli9fjoMHD2LatGmIiIjAvffe69L2ELmCtQUMgdq1cTQ6vSlghKjktcNJDQaf2uJmuR3XO6tMXafHx429S0REruRU8Bk7dqyLm9Gw0aNHY/To0abPU1NTcezYMbzzzjum4LNq1SpoNBqsWLECSqUS3bp1w759+7Bs2TIGH/JKOhs1PuazuowBIzRAUVusbKOORghhKjq26PFxx1BXneDjzlWiiYhcyangs2TJEle3w2HFxcUWdUY7duzA0KFDoVQqTcfS0tLwwgsv4NKlS2jVqpXVx1Gr1VCr1abPS0pK3NdoIjOmVZbrzuqS187qKq2qBgCEBshNAclWkDH29gB1iqFdHHz0eoFyjWFIre5QF9fxISJv51SNj6edPHkSb7zxBu677z7TsZycHMTFxVlcZ/w8JyfH5mMtXboU4eHhprfExET3NJqoDls9PiorPT4hKnmjPTjGQmjjYxhH0FzdC2MsbDa2C+BQFxH5DqeCj06nw8svv4wBAwYgPj4ekZGRFm/2evTRRyGRSBp8O3r0qMXXZGVlYfTo0bjtttswc+ZMZ5pvYeHChSguLja9nTt3rsmPSWSPap31Gh9jj09JldZ0TWiA3KyA2PrjGUOSUi6FQua+dXzK1YaAJZUYZpsBHOoiIt/h1FDXk08+iffffx8LFizAv//9bzz++OM4c+YMvvnmGyxevNjux1mwYAGmTJnS4DWpqammjy9cuIBrrrkGV155Jd59912L6+Lj45Gbm2txzPh5fHy8zcdXqVRQqVR2t5nIVYybjsrq1vjUBKGCMsOsRYkECFbWzuqytVaOMZDUnWLu6rV1LlUY2tUqSGlavLR2Q1SX3oqIyOWcCj6rVq3Ce++9hxtuuAFPPPEE7rrrLrRv3x49e/bEn3/+iQcffNCux4mJiUFMTIxd12ZlZeGaa65Bv379sHLlSkjrrH0yaNAgPP7446iurjbNMtuwYQM6depks76HyJOMs7oUNmp8CsoNtWchSsMwlzHI2JrOXmaa+m6YFeau4ubCckPwiQyuradjjw8R+QqnhrpycnLQo0cPAEBISAiKi4sBADfeeCN+/PFH17WuRlZWFq6++mq0a9cOL7/8Mi5evIicnByL2p27774bSqUS06dPx+HDh/H555/j9ddfx/z5813eHiJXMK7jI5Nan85uDBghAfYVEBtnWgUrDdc3VgztrIKadrUyCz7G0Tqu3ExE3s6pHp+EhARkZ2ejXbt2aN++PdavX4++ffti165dbhk22rBhA06ePImTJ08iISHB4pyo+UMbHh6O9evXIz09Hf369UN0dDQWL17MqezktYyBRGFjOnt+zVBXvSnjjQQfdxccX6oJPlHmwUfCWV1E5BucCj633HILNm7ciIEDB2LOnDmYMGECPvjgA2RmZmLevHmubiOmTJnSaC0QAPTs2RO///67y+9P5A7Vtmp8aoJPXkkVACAiyDB0K29kOMk0AyzAvWvrFFgZ6pJwVhcR+Qings/zzz9v+viOO+5AUlIStm/fjo4dO2LMmDEuaxxRS2bq8am3crOxxqe2iBhovAenXG19GwlX98IU1tQeWfT4cJNSIvIRDgef6upq3HfffVi0aBFSUlIAAFdccQWuuOIKlzeOqCUzrrJcdwFDY/AxigoxBIxGh7pqFhUMUdbZP8vVs7rKDYsqtrIy1GVrqj0RkbdwuLhZoVBgzZo17mgLkV/R1KQEpY11fIxMPT6NFBCX1t041BhGXLxL6cVSQ49PdEhtPZ/UTYslEhG5mlOzusaOHYtvvvnGxU0h8i/GHp+6QaduEDLW0hhnf2ltBJna4mZZzfXu6fHJLqkEALSJCDAd45YVROQrnKrx6dixI5566in88ccf6NevH4KDgy3O27uOD5E/sxV8VAobPT6NbFJqq8bHlcNPer1ATrGh6Do+PNB03F2rRBMRuZpTweeDDz5AREQEdu/ejd27d1uck0gkDD5EdjAOddWt6QlUyCw+r+3xMXxuq1elpKp2J3fD9cYw4rrkU1CuQbVOQCIBYkNrh7rc1btERORqTgWfjIwMV7eDyO+oq401PpZBJyJIafF5XJhhSEnaSI9PcaWm5usVlte7sBfG2NsTE6KymI3GdXyIyFc0eXd2IYRpEUEisp+puNlGMbORsZbGOPvLVpApqjDMtooItFz3x5VZJLvYUN/TOjzA4rjUDfciInIHp4PPBx98gO7duyMgIAABAQHo3r073n//fVe2jajF0umFKcDUCz7BCtPHAQopwgPt68EprjQEn/CgukNdLuzxKTHW91gGn8bqj4iIvIVTQ12LFy/GsmXLMGfOHAwaNAgAsGPHDsybNw+ZmZl46qmnXNpIopbGWNgMNNzjE6yU1+6A3kCQEUKgyBh8jEHJDcEn21jYHFYn+HBWFxH5CKeCzzvvvIP33nsPd911l+nYTTfdhJ49e2LOnDkMPkSNMA8+dYubFTIpYkJVuFiqRu/ECNPxhlZHrqrWmx7TWCPU2BYXzsgsrAAAJLQKsjhek804q4uIvJ5TQ13V1dXo379/veP9+vWDVqttcqOIWjq1zrDKskRSG1DMPX1zN1zVIRoLr+9iOmYcTtJaCRfGYS65VIJgpaFY2nxozFV1eJkFhuDTLsoy+HDLCiLyFU4Fn4kTJ+Kdd96pd/zdd9/FPffc0+RGEbV0pjV8ZFLTUJa50d1b438zBqJDbIjpWEPDSUVmM7rqDo0Bris6PltQDgBIqht8uI4PEfkIu4e65s+fb/pYIpHg/fffx/r16017dO3cuROZmZmYNGmS61tJ1MLYWrywIQ3ttm6c0RUWWFsYbR58dHpRbxd4RxVVaExrBbWLtAw+7toJnojI1ewOPnv37rX4vF+/fgCAU6dOAQCio6MRHR2Nw4cPu7B5RC2TrcULG2Kq2bGyZYVxqCvCRvBxxRDUmZphrthQFYKUln86uI4PEfkKu4PP5s2b3dkOIr9iPtRlr4YWMCw2ruETVH/HdMB6XZCjbA1zAeYzyJp8GyIit2ryAoZE5Dhnhroa2nvrUoWhxie8gaGupjIVNkcG1zvX2M7xRETewqnp7FVVVXjjjTewefNm5OXlQV9nL6A9e/a4pHFELVVltWFWV0Cdfbka0tDMqYJyQ/CJClbWux5wzRDU6XxDj09KdP0eH9NQF4MPEXk5p4LP9OnTsX79etx6660YMGCA1VkpRGRbpcYQfAKVjgcfrZUun/wyNQAg2mzjUPNaZlcMdR3PLQUAdIwLrXfOHYslEhG5g1PB54cffsBPP/2EwYMHu7o9RH7B2OMT5EjwMfWq1D+XX1a/x0cikUAmlUCnF03uidHpBU7mlQEALrMSfLhlBRH5CqdqfNq2bYvQ0Pp//IjIPhXGHh8nhrqs9aoUWOnxAVy3vk5mYQXUWj1Ucmm9qeyA2Sal7PEhIi/nVPB55ZVX8K9//Qtnz551dXuI/ELtUJf9na6yBtbKMQ11BVsGH2nNb3hTg49xmKt9TIjV9YBqV4lu0m2IiNzOqaGu/v37o6qqCqmpqQgKCoJCobA4X1hY6JLGEbVUpqEuF/T46PUCBTVDXdGhSotzcqkUgL7JwedETfDpFG+9p5ezuojIVzgVfO666y5kZWXhueeeQ1xcHIubiRxUoTGsgOxIcbPUxrBVSVW1qXg5MlhZ52sM75tae3Mkx1jYHGL1PLesICJf4VTw2b59O3bs2IFevXq5uj1EfqGiCbO6AEMvj7GuxjjMFRYgh0ous/o1Ta29OXi+GADQo2241fPcsoKIfIVTNT6dO3dGZWWlq9tC5DeqmjDUBVhOT79YWjPMFaKy8jXSetc76lK5BpmFhsULe7aNsN62mh4fV+0CT0TkLk4Fn+effx4LFizAli1bUFBQgJKSEos3ImpYk3t8zAJGQXlNYbPV4GN435QhqINZht6e5KgghAcprF7DdXyIyFc4NdQ1evRoAMCIESMsjgshIJFIoNPpmt4yohbMGHzqbvbZEPO9t8wDRm6JIfjEhFoJPi5YUfnA+SIAQM+ECNtt415dROQjnAo+3LCUqGlqp7M7sEmp2aXmtTQ5xYZh59bhAfW+RiarWe25CT0x+2vqe3omWK/vAbhlBRH5DqeCz7Bhw1zdDiK/Ulpl2E09VGV96MgauVny0elqA0Z2cRUAIN5a8JE0rbhZrxfYdcawPEXfpFY2r+NQFxH5Cqd3Z//9998xYcIEXHnllcjKygIAfPLJJ9i2bZvLGkfUUpVWGaazhwXaH3zM1w207PExBJ/W4YH1v6aJgeRYbimKKqoRpJTZnNEFADIXTZsnInI3p4LPmjVrkJaWhsDAQOzZswdqtaHGoLi4GM8995xLG0jUEpXU9PiEBdrf6SqRSEzhx7wHp6EeH3kTp5nvPF0AAOifHAmFzPafC25ZQUS+wqng88wzz2D58uV47733LFZtHjx4MPbs2eOyxhG1VCWVhh6f0AD7e3yA2uEuY5DR6wVyS4w9PvWDj61FD+3152nDMNcVqZENXtfU+xARNRengs+xY8cwdOjQesfDw8NRVFTU1DYRtWhV1TpoaqY/hQU4VmZnLPPR1tT45JerodULSCU2ZnU1YahLq9NjR02Pz8CUqAavNS2UyKEuIvJyTgWf+Ph4nDx5st7xbdu2ITU1tcmNImrJjMNcEgkQ7MB0dqD+7CljfU9MqMrqUJS8CYFk99lLKK6sRkSQAr0TI+xsl8O3ISJqVk4Fn5kzZ+Khhx7Czp07IZFIcOHCBaxatQoPP/wwZs2a5eo2ErUopmEuldxUG2OvusXKtfU99Qubza/X6hxPJJuO5gEArukUa3VH9obaRUTkrZyazv7oo49Cr9djxIgRqKiowNChQ6FSqfDwww9jzpw5rm4jUYtimsruYH0PYFasrLfs8WkdVr++B2ja+jq/HskFAIzoEtvotdydnYh8hVPBRyKR4PHHH8cjjzyCkydPoqysDF27dkVIiPWdm4moVnGlcUaX48FHVmeW1oWaxQutzegCAHnNPPNqB3t8Tl0sw6mL5ZBLJRh6WUyj17O4mYh8hVPBx0ipVKJr166uaguRXygoM24qqnT4a+sGjPOXDMEnoZX1oS5j3U+1g3tJfLvvAgDgqo7RCLOjZ6opRdRERM3JoeAzbdo0u65bsWKFU40h8gfGTUWjgh0PPqbZUzU55nzNrukJrYKsXq90IvgIIfDtPsOipGN7t7WvXdyygoh8hEPB58MPP0RSUhL69OkDwT9wRE7JN/X41J9+3hhj8NHWJB9jj09ipPUeH6XcEHw0Dgx17TtXhLMFFQhSyjCqW5xdX8PiZiLyFQ4Fn1mzZuGzzz5DRkYGpk6digkTJiAysuGFzYjIUn5ZTY9PE4KPXgiUq7UoKDeEKFs9PqahLq39PT5f/H0OADCqa5zdu8cbe3ycmDxGRNSsHJrO/tZbbyE7Oxv/93//h++//x6JiYm4/fbb8csvv7AHiMhOxhqfKCdqfEwBQ1/b2xMWIEe4jUJpY/DR2DnUlV+mxpo9hmGuuwcm2d0u48KK3LKCiLydw+v4qFQq3HXXXdiwYQP++ecfdOvWDQ888ACSk5NRVlbmjjYStSjGGh9nipvNi4jP1dT3JEZa7+0BAKW8ZlaXnT0+n+w4C41Wj16JEbg82fZu7HVJWeNDRD7C6d3ZAUAqlUIikUAIAZ1O56o2NUitVqN3796QSCTYt2+fxbkDBw5gyJAhCAgIQGJiIl588cVmaRORIy6WGoOP80NdOr3AuUs1wcfGMBfg2KyuqmodPvnzLABg5pAUSCT2L67IWV1E5CscDj5qtRqfffYZrr32Wlx22WU4ePAg3nzzTWRmZjbLOj7/93//hzZt2tQ7XlJSglGjRiEpKQm7d+/GSy+9hCeeeALvvvuu29tEZC+NVo+8muDT2sZqyw1R1RQrq7U6nCtsuLAZqJ3VZU9x89o9WSgs16BtRCBGd4t3qF2c1UVEvsKh4uYHHngAq1evRmJiIqZNm4bPPvsM0dHR7mpbPevWrcP69euxZs0arFu3zuLcqlWroNFosGLFCiiVSnTr1g379u3DsmXLcO+99zZbG4kakltSBSEMs62cmc4eqJQBACo0utoenwaGuhTGWV2NDHXp9QLv/34aADDtqhTIrez71RDO6iIiX+FQ8Fm+fDnatWuH1NRUbN26FVu3brV63dq1a13SOHO5ubmYOXMmvvnmGwQF1f9Dv2PHDgwdOhRKZe0/JmlpaXjhhRdw6dIltGplf70CkbsYC5LbRgQ6vE8XANMsqwqNtrbGxwVDXZuO5uF0fjlCA+S44/JEh9ulkBmn2TP4EJF3cyj4TJo0yaFxf1cRQmDKlCm4//770b9/f5w5c6beNTk5OUhJSbE4FhcXZzpnK/io1Wqo1WrT5yUlJa5rOFEdF4pqg48zgsx6fBpbwwcAlKYtKxoOPu/V9PbcPaAdQlSOL+gur5nW5cxmqEREzcnhBQxd6dFHH8ULL7zQ4DVHjhzB+vXrUVpaioULF7r0/gCwdOlSPPnkky5/XCJrsmqCT5sI63trNcYYfLIuVaJMbdjl3dYaPoD5Aoa2g8+B80XYmVEIuVSCKYOTnWqXo9PmiYg8pUl7dTXVggULMGXKlAavSU1NxaZNm7Bjxw6oVJazYPr374977rkHH330EeLj45Gbm2tx3vh5fLztQs2FCxdi/vz5ps9LSkqQmOh4Vz+RPYzDU20jbIeVhhiHug5dKK55nEAEKGQ2r68d6rLdE/Pe7xkAgDG92jhVcG24T81QF4MPEXk5jwafmJgYxMQ0vvPzf/7zHzzzzDOmzy9cuIC0tDR8/vnnGDhwIABg0KBBePzxx1FdXQ2FwrCY24YNG9CpU6cG63tUKlW9QEXkLqcuGta6ah8b7NTXG3t8DmcZhmRTYxp+HFNPjNb6chNZRZX46WA2AGDGkBSr19jDnoBFROQNPBp87NWuXTuLz43T5tu3b4+EhAQAwN13340nn3wS06dPx7/+9S8cOnQIr7/+Ol599dVmby+RNUIInMwzBJ8Osc4t/WAMPqU1w1wp0Q0HH2UjgeSj7Weg0wtc2T4K3dqEO9UmAJDbWUtERORpPhF87BEeHo7169cjPT0d/fr1Q3R0NBYvXsyp7OQ1LpapUVKlhVQCJEc52+Nj+Sub2kjwUchtB5IKjRar/8oEAEy/yvneHsC5XeCJiDzBJ4NPcnKy1b3Bevbsid9//90DLSJq3Km8cgCGdXcaqstpSEyo5bBsakzDPUdKmeE+1tbxWbsnCyVVWiRFBeGaTrFOtcfIuO6PXhjW8pE5MVWfiKg5NGnLCiKy39EcQ11Oh0bCSkPiwixng3VpHdbg9caVnquqLWt8hBD4aPsZAMDkQclOrSlkzljcDLDXh4i8G4MPUTM5cN4wE6tHgvO1NPFmwSc2VFWvB6iuIFXtuj/m9p4rwom8MgQqZLi1f4LT7TFSmK30zEUMicibMfgQNZP954sAAL0SIpx+jDYRAaYFBgekRDZ6fbBppWfL4PP1niwAQFq3OIQFKJxuj5HcrMfI3p3giYg8gcGHqBmUVFXj9EVDjU/PJvT4yGVSPHVzN4zsEof/S+vc6PXBNT0+5Rqt6ZhGq8f3By4AAMb1bXpvD2DYnd24qHu1nsGHiLyXTxY3E/mafZlFAAwLDkaFNG3dqHF9E+wOLKa9vdS1PT6/n7iIoopqxIaqMLiDazYZlkgkUEil0Oj0XMuHiLwae3yImsEfp/IBAFekRjXrfY1DXYZAYuiJ+fVIHgBgdPd4l86+4urNROQLGHyImsH2kwUAgMEdmjf4BCprp81XaHQQQmDLMUPwuaZz06aw1yXnWj5E5AM41EXkZkUVGtPeWq4aWrKXUi6FQiZBtU6gQqNF1qVKZBdXIUAhxSAX9z5x2woi8gXs8SFys01H8yAEcFlcSL11eJqDsc6nXK3FpqOGjXuv6hDt9CKKtii4bQUR+QAGHyI3++lgDgBgdPfWHrl/RJBhunpRRTU2HXXPMBfAHh8i8g0MPkRuVKbW4rcTFwEA13WP90gbomtmkR3LLcXec0UA0OQtKqzhRqVE5AsYfIjc6KeD2dBo9UiJDkbn+FCPtCGmJvh8tfs8hAA6x4eiTUSgy++jkBr+nGjZ40NEXozBh8iNPqvZ/fzWfgmQSDyzcadxW4u9NWsJjewS55b7NLQTPBGRt2DwIXKTI9kl2JtZBLlUgttcsB+Ws+r27gzv4vphLqC2xkfD4ENEXozBh8hN3v89AwAwqlscYkObfzaXUY+2tVtkxIcFNGmvsIYE1swSq7sTPBGRN2HwIXKDc4UV+GafYSPQe4e292hb+ie3QptwQ/C6f1iqS1drNsfgQ0S+gAsYErnB21tOQacXGNIxGr0TIzzalgCFDF+nD0ZOcRV6ubEtATWrRFdqGHyIyHsx+BC52NGcEny+y1DU/OCIjh5ujUFcWIDbF0809vhUVrPGh4i8F4e6iFxICIGnvv8HemFYt+fy5EhPN6nZ1AYf9vgQkfdi8CFyoe8PZGP7qQIo5VI8dn0XTzenWRk3RGWNDxF5MwYfIhfJK63C4m8PAQAeuLo9EiODPNyi5mXc+4s1PkTkzRh8iFxACIHH1h5EUUU1urUJQ/o1HTzdpGbHWV1E5AsYfIhc4L3fT+PXI3lQyqRYdntv02J+/iRQYXjOrPEhIm/mf3+diVxs+8l8PL/uKABg0Ziu6OShPbk8jTU+ROQLGHyImuBcYQVmf7YXegGM75uACQPbebpJHmOs8algjQ8ReTEGHyInFZSpMWnFXygs16BbmzA8e0t3j21E6g3CAhUAgOLKag+3hIjINi5gSOSECo0W0z76Gxn55WgbEYgVUy439Xj4q4ia4FNUYTv4VGp0OJFXinK1DmGBciREBCE8SNFcTSQiYvAhcpRGq0f6qj3Yf64IEUEKfDRtgNtXRfYFrYKUAICiCk29cxn55Xh5/TH8cigHWr2wOJccFYR+SZG4pnMMhl0Wg9AABiEich8GHyIHaLR6PLBqNzYfuwiVXIoPJvdHh9gQTzfLK0TU9NyUa3TQaPVQyg0j6ZuP5SF91R5T7U90iBLhgQoUV2qRX6bGmYIKnCmowJo956GQSTCofTRu7tUGad3jEaLinygici3+VSGykyH07MGvR/Kgkkvx/uT+6JfkP1tSNCYsQAGJBBACKKrUIDY0AAfOF2HW/3ajqlqPgSmRWDKmG7q0DjXVQhVXVGPvuUvYcaoAG/7Jxen8cvx2/CJ+O34Rj39zEGnd4jG2T1sM6RANuR8uEUBErsfgQ2QHjVaP9E/34NcjuVDKpXhvUn8M6Rjj6WZ5FalUgvBABYoqqlFYrkFUsAqPrjmIqmo9ru4Ug/cm9a+3vlF4kAJXd4rF1Z1isfD6LjiZV4afDmbj671ZyMgvx7f7LuDbfRcQHaLEzb3bYlzftujaOsyvi8iJqGkYfIgaodHqMfvTPdjwjyH0vD+pP4ZextBjTZvwQBRVVONcYSX+yijEP9klCAuQ4+Xbetm1qGOH2BA8OKIj5gzvgP3ni/HN3ix8v/8C8ss0+GBbBj7YloHO8aEY17ctxvZui1jWVhGRgyRCCNH4Zf6jpKQE4eHhKC4uRlhYmKebQx5WrTOEnl8O1/b0DGPosSl91R78eDAb9w1Lxac7M1FapcXTN3fDxEHJTj9mtU6P345fxNo9WdjwTy40Oj0AQCoBhnSMwbi+bTGqa7xpAUUi8k/2/vvNHh8iG4w9PetrenrendiPoacRKdHBAID/bj0NAOjWJgx3D0xq0mMqZFKM6BKHEV3iUFxRjR8OXsDaPVnYffYSth6/iK3HLyJEJcf1PeIxvm8CLk+OhFTKoTAiso49PnWwx4cAQK3VIb2mkFkpl+K/E/vhmk6xnm6W19t+Mh93v78TACCRAGtnXYk+7Vq55V4Z+eX4em8W1u45j/OXKk3HE1oFYlyftrilb4IpiBFRy2fvv98MPnUw+FBVtQ4PrNqDTUcNs7fe5fCW3XR6gXs//htbjl/EglGX4YGr3b9LvV4vsOtMIdbuycKPB7NRptaazvVtF4ExvdpgdPd4tA4PdHtbiMhzGHycxODj36qqdbjvk93YevwiAhRSvD/pclzVMdrTzfI5VdU6j6xkXanRYf0/OVi7Jwu/n7gI87US+7SLwPXdW2N093gkRgY1e9uIyL0YfJzE4OO/qqp1mPnx3/j9RD4CFTJ8MLk/ruzA0OOr8kqq8N3+C/j5UA52Z16C+V+6Hm3DMbp7PEZ0iUWnuFBOjydqARh8nMTg458qNTrM+HgX/jhZgCClDCumXI4rUqM83SxykdySKvxyOAfrDuZgZ0aBRU9Qm/AAXN05Ftd0isWV7aMQzNWiiXwSg4+TGHz8T4VGi2kf7sKfpwsRrJRh5dQBGJDCFZlbqvwyNTb8k4v1h3Ow/VQB1Fq96ZxSJsXA1Ehc1SEag9pHoWvrMK4YTeQjGHycxODjX8rVWkz9cBf+yihEiEqOD6dejv7JDD3+oqpahx2nC7DlaB42HcvDucJKi/OhKjkGpETiitQoXJEaha5twiDjVHkir8Tg4yQGH/9RptZiyoq/8PfZSwhVyfHR9AHo66ap1+T9hBA4dbEcW47l4c/TBdiZUYjSKq3FNcFKGXokhKN3Yiv0aReBPokRXD2ayEsw+DiJwcc/lFZVY8rKXdh99hJCA+T4ZPpA9E6M8HSzyIvo9AL/XCjBn6cLsON0Af7KKLSYKm/UJjwAPRMi0KV1GDrFh6JL61AktgriIopEzYzBx0kMPi1faVU1Jq/4C3syixAWIMf/ZgxEz4QITzeLvJxOL3AirxT7Mouw75zh7XhuqUWhtFGQUoZO8aHoHB+K1OgQJEcHIyU6CImRQVDJubUGkTsw+DiJwadlMw894YEK/G/6QPRICPd0s8hHlau1OHC+GIcvFONoTimO5pTgeG4ZNGYF0+akEqBNRCBSooORFBWE1uGBaBMRYHgfHoi4cBWDEZGTWmTw+fHHH/HUU0/hwIEDCAgIwLBhw/DNN9+YzmdmZmLWrFnYvHkzQkJCMHnyZCxduhRyuf3TUxl8Wq6SmtCztyb0rJoxEN3bMvSQa2l1epwpKMfRnFIcyylFRn45MvLLcSa/HOUaXaNfHx2iQuvwAESFKBEZrERUsBKRwSpEhRg/ViI8UIGQADlCVQoEKKRch4gILXCT0jVr1mDmzJl47rnnMHz4cGi1Whw6dMh0XqfT4YYbbkB8fDy2b9+O7OxsTJo0CQqFAs8995wHW07eoKSqGpM++Av7zjH0kHvJZVJ0iA1Fh9hQ3Niz9rgQAhfL1DiTX4Ez+eXILKzAheJKZBdVIbu4EheKq6DR6pFfpkZ+mdru+8mkEoSo5AhRyREaYHgfqJRBJZdCJZdBKZfWfCyt+dhwTCaVQCqRQCox7KsmlUggMX4OQCo1fC6EgE5v9iYE9HoBrd7wXifMPtYDeiGgr/kavTBsKaIXhuuEMNzr8uRIpHWLR2Sw0vXfAKJG+ESPj1arRXJyMp588klMnz7d6jXr1q3DjTfeiAsXLiAuLg4AsHz5cvzrX//CxYsXoVTa9wvGHp+Wxzz0RAQZhrcYesjbCCFQWK5BdnEVcoqrUFiuQUG5BoXl6pr3GhSUGd6XVFWjTK2F9//1tk0mlWBQahSu6xGPtG7xiA5RebpJ5ONaVI/Pnj17kJWVBalUij59+iAnJwe9e/fGSy+9hO7duwMAduzYgR49ephCDwCkpaVh1qxZOHz4MPr06WP1sdVqNdTq2v9dlZSUuPfJULOqG3pWzRiIbm0Yesj7SCQSRIWoEBWisiuY6/UCldU6lKm1KK3SokytRVmVFqVV1ais1kGj1UOt1de8r/1cXfO5Tm/ogdELQ+gSMPbWGN5DGAq6pVJDb5BcKoFUKoFMIoFcZugtkklr3mo+Np6XSg09R7Uf13wulaC0SouNR3NxKKsE207mY9vJfCz65hAGpkTh+p6tkdYtDrGhXCKA3Mcngs/p06cBAE888QSWLVuG5ORkvPLKK7j66qtx/PhxREZGIicnxyL0ADB9npOTY/Oxly5diieffNJ9jSePuVSuwZSVf2H/+WKGHmpxpFIJglVyBKvkiPOxzul5116GswXlWHcoBz8dzMaB88XYUbNswOJvD2FAciSu72HYUDaO6ySRi3l0LfZHH30UkppxZVtvR48ehV5vmCHx+OOPY/z48ejXrx9WrlwJiUSCL7/8skltWLhwIYqLi01v586dc8VTIw/LK6nCne/+if3ni9EqSIFPZ1zB0EPkRZKignH/sPb4bvZV+P3/rsFj13dG78QICAHszCjEku8O44qlG3Hb8u1YsS0D2cWVjT8okR082uOzYMECTJkypcFrUlNTkZ2dDQDo2rWr6bhKpUJqaioyMzMBAPHx8fjrr78svjY3N9d0zhaVSgWVimPLLcn5SxWY8P5OnCmoQFyYCv+bPhAd40I93SwisiExMgj3Dm2Pe4e2R1ZRJdYdzMZPB7OxJ7MIu85cwq4zl/DUD/+gV0I4ru0ah1Hd4tExNoSz2cgpHg0+MTExiImJafS6fv36QaVS4dixY7jqqqsAANXV1Thz5gySkpIAAIMGDcKzzz6LvLw8xMbGAgA2bNiAsLAwi8BELdvpi2WY8P5OXCiuQkKrQHw64wq0iwrydLOIyE5tIwIxY0gqZgxJRXZxJdYdNAyH7c68hP3ni7H/fDFeXn8cSVFBGNU1Dtd2jUe/pFbcQ43s5hOzugBg7ty5+Oqrr7BixQokJSXhpZdewvfff4+jR4+iVatW0Ol06N27N9q0aYMXX3wROTk5mDhxImbMmOHQdHbO6vJdh7KKMWXlX8gv06B9TDBWzbgC8eGsDyBqCS6WqrHxSC7W/5OLbSfzLRaJjApWYnjnWFzbNQ6DO0QjWOUT5avkYi1uAcPq6mosXLgQn3zyCSorKzFw4EC89tpr6Natm+mas2fPYtasWdiyZQuCg4MxefJkPP/881zA0A9sOZaH9FV7UK7RoVubMHw8bQCiOD2WqEUqV2vx2/GL2PBPLjYezUNxZbXpnEImQf+kSAzrFIOhHWPQpXUoh8T8RIsLPs2Fwcf3fL4rE499fQg6vcDgDlF4Z0I/hAUoPN0sImoG1To9dp0pxPrDudh0NA+ZhRUW52NDVRh6WQyGXhaDwe2j+B+iFozBx0kMPr5DCIFXfz2B/2w8AQAY17ctnh/XE0q5RycrEpEHnckvx9bjF7H1+EXsOFWAymrLbUI6xoZgYGokrkiNwsCUKMSEMgi1FAw+TmLw8Q0VGi0e+eoAfjxgmPE3Z3gHzL/2MnZpE5GJWqvD32cuYevxi/jt+EUczSmtd01qTDAGpkShf1Ir9G4XgZSoYEhZKO2TGHycxODj/c4VVmDmx3/jaE4p5FIJnhnbHXcOaOfpZhGRlyss1+CvjEL8eboAOzMKcTSnpN62H2EBcvRKjEBvszcOj/kGBh8nMfh4t+0n85H+6R5cqqhGdIgK70zoi8uTIz3dLCLyQUUVGuw6cwk7Txdg77kiHMoqhtpstphR6/AAdI4PRZfWYejcOgxd4kOREh0MuYzD6t6EwcdJDD7eSacXeHPTSby+8Tj0AuiZEI7/TuyH1uGBnm4aEbUQ1To9juWUYu+5IuzLLML+80U4mVdm9VqlXIrL4kLQISYEKdEhSI0JRkq04Y3T6T2DwcdJ7go+lRodAhRS1qA4Ibu4EnNX78POjEIAwK39EvDM2O4IUMg83DIiaulKq6pxLKcUR3JKcSS7BEezS3A0pxQVGp3Nr4kLUyElOhjtIoPQJiIQbSIC0bbmLT48gH+73ITBx0nuCj7TP9wFrV5g8ZiuaB8T4rLHbenWH87B/605gKKKagQrZXjmlu64pU+Cp5tFRH5Mrxc4d6kCR7JLkZFfjtMXy5CRX46M/HIUlGsa/froEBXiw1WIClYhOkSF6FAlooMN76OCVYgKUSI8UIHQAAVCVXKfKLZWa3UordKitEqLkspqw/uqapRWVaOkUmt4X3OsXK3F8gn9XN4RYO+/3+yPawYZ+eX4/UQ+NDo90l79DdOuSsGc4R0QyrVmbLpUrsET3x/Gt/suAAC6tw3DG3f1RUp0sIdbRkT+TiqVICkqGElR9f8eFVdU43S+IQhlXarEheJKnL9UiQtFlcgqqkRVtR75ZWrkl6ntvl+oSo7QADlCAxQIC5QjRCVHoFIGlVwGlVxqeFOYfSyXQSqVQCoBJDXtlUgkho8lEhjzhlanh1YvoNUJVOv10OpEzeeG49U6Paqq9ajUaFFZrUOl2ccVGh2qNDrTx9Zqoxqi1uo91vPFHp863NXjk5Ffjmd++Acbj+YBAGJCVXgkrRPG9WnLArk61h3MxqJvDyG/TAOpBJg5NBXzr70MKjm7h4nIdwkhUFRRjayiSuSVViG/VIP8cjXySzUoKDeEoYIyDfLLNCipqrbYlsNXhKrkCAtU1AQ1OcICDB8bjxk+V2B8v7Yu/5vOoS4nubu4efPRPDz5/WGcKTCsLpoUFYT0azrglj5tofDzAHSusAJP/fAPNvyTC8Cw0NhLt/VC78QIzzaMiMgDqqqNw0fVtcNINcNHaq0e6mo91FpDb4vhcx2qao7phCFoCQEICOj1gF4ICBiOA4BcKoVcJoFCJoVMKoFCJql3LFAhQ5BShoCa94EKGQKUMgQpZAhU1p4LDVAgRCX36GaxDD5Oao5ZXWqtDh/+cQbLt57CpQrDHjPxYQG4c0Ai7ry8nd9trFmp0eGdraewfOspaLR6yKQSzBrWHnNGdGAvDxER2YXBx0nNOZ29XK3Fqp1n8e5vp5FfZiiIk0klGNIxGtf3aI20rvEID2q5dUB6vcD3By7gxZ+PIauoEgAwuEMUnhjTDR3jQj3cOiIi8iUMPk7yxDo+aq0OPx/KwaqdmfirZso2AMilEgxqH4WrOkTjqo7R6BIf5hPV/Y0RQmDT0Ty89Msx0xLybcID8O8bu+K67vGc8k9ERA5j8HGSpxcwPHWxDD8eyMaPB7JxLNdyX5moYCX6J7cyLKeeEIEeCeE+NzNsx6kCvPTLUezJLAJgKIS7d2gqZgxJRaCSw1pEROQcBh8neTr4mDuZV4atxy9i24mL2JlRWG/BLIkESGwVhA6xIYa3mBC0jw1Gm4hAxIYGeLTITKPVI6e4CueLKnChqAoXiiqxM6MAf5wsAAAEKKSYOjgF9w1NRUSQ0mPtJCKiloHBx0neFHzMabR67D9vWEZ93/ki7D9XhPOXKm1eL5NKEBeqQpuIQMSFBSAiSIFWQUpEBCkQHmh4C1TKoJRZrv8gl0qhFwI6ISCEgF4YZgKoq/Uo12hRqdGhXKNDhVqLco0OxZXVKCxXo7Bcg4IyDQrLa94qNPU2/wMAhUyCuwe0Q/o1HRAb5l9F3ERE5D5cwLCFUcqluDw50mJDzoIyNY7nluHkxTKcyivDyTzDolm5JVXQ6gUuFFfhQnGVx9qskkvR1my59sTIQNzcuy0SI4M81iYiIvJvDD4+LCpEhUEhKgxqH2VxXKcXuFiqxoXiSmQXVSGvtApFFdUoqtCgqLLa8HFlNdTVOmiM6z9oa9Z+0AvIalb2NKz8aVj9UyU3rNkQrDS+N6wcGh6oQFSwEpHBSrQKNiy3HhmsRGyYClHBShYqExGRV2HwaYFkUgniwwMM6wG183RriIiIvId/LxVMREREfoXBh4iIiPwGgw8RERH5DQYfIiIi8hsMPkREROQ3GHyIiIjIbzD4EBERkd9g8CEiIiK/weBDREREfoPBh4iIiPwGgw8RERH5DQYfIiIi8hsMPkREROQ3GHyIiIjIb8g93QBvI4QAAJSUlHi4JURERGQv47/bxn/HbWHwqaO0tBQAkJiY6OGWEBERkaNKS0sRHh5u87xENBaN/Ixer8eFCxcQGhoKiUTissctKSlBYmIizp07h7CwMJc9rq/w9+cP8DXw9+cP8DXw9+cP8DVw5/MXQqC0tBRt2rSBVGq7koc9PnVIpVIkJCS47fHDwsL88ofdyN+fP8DXwN+fP8DXwN+fP8DXwF3Pv6GeHiMWNxMREZHfYPAhIiIiv8Hg00xUKhWWLFkClUrl6aZ4hL8/f4Cvgb8/f4Cvgb8/f4CvgTc8fxY3ExERkd9gjw8RERH5DQYfIiIi8hsMPkREROQ3GHyIiIjIbzD4NJO33noLycnJCAgIwMCBA/HXX395uklu8cQTT0AikVi8de7c2XS+qqoK6enpiIqKQkhICMaPH4/c3FwPtrhpfvvtN4wZMwZt2rSBRCLBN998Y3FeCIHFixejdevWCAwMxMiRI3HixAmLawoLC3HPPfcgLCwMERERmD59OsrKyprxWTRNY6/BlClT6v1MjB492uIaX34Nli5dissvvxyhoaGIjY3F2LFjcezYMYtr7Pm5z8zMxA033ICgoCDExsbikUcegVarbc6n4hR7nv/VV19d72fg/vvvt7jGV58/ALzzzjvo2bOnaVG+QYMGYd26dabzLfn7DzT+/L3t+8/g0ww+//xzzJ8/H0uWLMGePXvQq1cvpKWlIS8vz9NNc4tu3bohOzvb9LZt2zbTuXnz5uH777/Hl19+ia1bt+LChQsYN26cB1vbNOXl5ejVqxfeeustq+dffPFF/Oc//8Hy5cuxc+dOBAcHIy0tDVVVVaZr7rnnHhw+fBgbNmzADz/8gN9++w333ntvcz2FJmvsNQCA0aNHW/xMfPbZZxbnffk12Lp1K9LT0/Hnn39iw4YNqK6uxqhRo1BeXm66prGfe51OhxtuuAEajQbbt2/HRx99hA8//BCLFy/2xFNyiD3PHwBmzpxp8TPw4osvms758vMHgISEBDz//PPYvXs3/v77bwwfPhw333wzDh8+DKBlf/+Bxp8/4GXff0FuN2DAAJGenm76XKfTiTZt2oilS5d6sFXusWTJEtGrVy+r54qKioRCoRBffvml6diRI0cEALFjx45maqH7ABBff/216XO9Xi/i4+PFSy+9ZDpWVFQkVCqV+Oyzz4QQQvzzzz8CgNi1a5fpmnXr1gmJRCKysrKare2uUvc1EEKIyZMni5tvvtnm17S01yAvL08AEFu3bhVC2Pdz/9NPPwmpVCpycnJM17zzzjsiLCxMqNXq5n0CTVT3+QshxLBhw8RDDz1k82ta0vM3atWqlXj//ff97vtvZHz+Qnjf9589Pm6m0Wiwe/dujBw50nRMKpVi5MiR2LFjhwdb5j4nTpxAmzZtkJqainvuuQeZmZkAgN27d6O6utritejcuTPatWvXIl+LjIwM5OTkWDzf8PBwDBw40PR8d+zYgYiICPTv3990zciRIyGVSrFz585mb7O7bNmyBbGxsejUqRNmzZqFgoIC07mW9hoUFxcDACIjIwHY93O/Y8cO9OjRA3FxcaZr0tLSUFJSYvG/Zl9Q9/kbrVq1CtHR0ejevTsWLlyIiooK07mW9Px1Oh1Wr16N8vJyDBo0yO++/3Wfv5E3ff+5Samb5efnQ6fTWXxDASAuLg5Hjx71UKvcZ+DAgfjwww/RqVMnZGdn48knn8SQIUNw6NAh5OTkQKlUIiIiwuJr4uLikJOT45kGu5HxOVn73hvP5eTkIDY21uK8XC5HZGRki3lNRo8ejXHjxiElJQWnTp3CY489huuuuw47duyATCZrUa+BXq/H3LlzMXjwYHTv3h0A7Pq5z8nJsfpzYjznK6w9fwC4++67kZSUhDZt2uDAgQP417/+hWPHjmHt2rUAWsbzP3jwIAYNGoSqqiqEhITg66+/RteuXbFv3z6/+P7bev6A933/GXzIpa677jrTxz179sTAgQORlJSEL774AoGBgR5sGXnKnXfeafq4R48e6NmzJ9q3b48tW7ZgxIgRHmyZ66Wnp+PQoUMWdW3+xNbzN6/X6tGjB1q3bo0RI0bg1KlTaN++fXM30y06deqEffv2obi4GF999RUmT56MrVu3erpZzcbW8+/atavXff851OVm0dHRkMlk9Sr4c3NzER8f76FWNZ+IiAhcdtllOHnyJOLj46HRaFBUVGRxTUt9LYzPqaHvfXx8fL0id61Wi8LCwhb5mgBAamoqoqOjcfLkSQAt5zWYPXs2fvjhB2zevBkJCQmm4/b83MfHx1v9OTGe8wW2nr81AwcOBACLnwFff/5KpRIdOnRAv379sHTpUvTq1Quvv/6633z/bT1/azz9/WfwcTOlUol+/fph48aNpmN6vR4bN260GP9sqcrKynDq1Cm0bt0a/fr1g0KhsHgtjh07hszMzBb5WqSkpCA+Pt7i+ZaUlGDnzp2m5zto0CAUFRVh9+7dpms2bdoEvV5v+uPQ0pw/fx4FBQVo3bo1AN9/DYQQmD17Nr7++mts2rQJKSkpFuft+bkfNGgQDh48aBEAN2zYgLCwMNNwgbdq7Plbs2/fPgCw+Bnw1edvi16vh1qtbvHff1uMz98aj3//XV4uTfWsXr1aqFQq8eGHH4p//vlH3HvvvSIiIsKigr2lWLBggdiyZYvIyMgQf/zxhxg5cqSIjo4WeXl5Qggh7r//ftGuXTuxadMm8ffff4tBgwaJQYMGebjVzistLRV79+4Ve/fuFQDEsmXLxN69e8XZs2eFEEI8//zzIiIiQnz77bfiwIED4uabbxYpKSmisrLS9BijR48Wffr0ETt37hTbtm0THTt2FHfddZennpLDGnoNSktLxcMPPyx27NghMjIyxK+//ir69u0rOnbsKKqqqkyP4cuvwaxZs0R4eLjYsmWLyM7ONr1VVFSYrmns516r1Yru3buLUaNGiX379omff/5ZxMTEiIULF3riKTmksed/8uRJ8dRTT4m///5bZGRkiG+//VakpqaKoUOHmh7Dl5+/EEI8+uijYuvWrSIjI0McOHBAPProo0IikYj169cLIVr291+Ihp+/N37/GXyayRtvvCHatWsnlEqlGDBggPjzzz893SS3uOOOO0Tr1q2FUqkUbdu2FXfccYc4efKk6XxlZaV44IEHRKtWrURQUJC45ZZbRHZ2tgdb3DSbN28WAOq9TZ48WQhhmNK+aNEiERcXJ1QqlRgxYoQ4duyYxWMUFBSIu+66S4SEhIiwsDAxdepUUVpa6oFn45yGXoOKigoxatQoERMTIxQKhUhKShIzZ86sF/p9+TWw9twBiJUrV5qusefn/syZM+K6664TgYGBIjo6WixYsEBUV1c387NxXGPPPzMzUwwdOlRERkYKlUolOnToIB555BFRXFxs8Ti++vyFEGLatGkiKSlJKJVKERMTI0aMGGEKPUK07O+/EA0/f2/8/kuEEML1/UhERERE3oc1PkREROQ3GHyIiIjIbzD4EBERkd9g8CEiIiK/weBDREREfoPBh4iIiPwGgw8RERH5DQYfIiIi8hsMPkTkFaZMmYKxY8d67P4TJ07Ec88957bH/+eff5CQkIDy8nK33YOIGseVm4nI7SQSSYPnlyxZgnnz5kEIgYiIiOZplJn9+/dj+PDhOHv2LEJCQtx2n1tvvRW9evXCokWL3HYPImoYgw8RuV1OTo7p488//xyLFy/GsWPHTMdCQkLcGjgaM2PGDMjlcixfvtyt9/nxxx8xc+ZMZGZmQi6Xu/VeRGQdh7qIyO3i4+NNb+Hh4ZBIJBbHQkJC6g11XX311ZgzZw7mzp2LVq1aIS4uDu+99x7Ky8sxdepUhIaGokOHDli3bp3FvQ4dOoTrrrsOISEhiIuLw8SJE5Gfn2+zbTqdDl999RXGjBljcTw5ORnPPPMMJk2ahJCQECQlJeG7777DxYsXcfPNNyMkJAQ9e/bE33//bfqas2fPYsyYMWjVqhWCg4PRrVs3/PTTT6bz1157LQoLC7F169YmvqJE5CwGHyLyWh999BGio6Px119/Yc6cOZg1axZuu+02XHnlldizZw9GjRqFiRMnoqKiAgBQVFSE4cOHo0+fPvj777/x888/Izc3F7fffrvNexw4cADFxcXo379/vXOvvvoqBg8ejL179+KGG27AxIkTMWnSJEyYMAF79uxB+/btMWnSJBg7ztPT06FWq/Hbb7/h4MGDeOGFFyx6spRKJXr37o3ff//dxa8UEdnNLXu+ExHZsHLlShEeHl7v+OTJk8XNN99s+nzYsGHiqquuMn2u1WpFcHCwmDhxoulYdna2ACB27NghhBDi6aefFqNGjbJ43HPnzgkA4tixY1bb8/XXXwuZTCb0er3F8aSkJDFhwoR691q0aJHp2I4dOwQAkZ2dLYQQokePHuKJJ55o8PnfcsstYsqUKQ1eQ0Tuwx4fIvJaPXv2NH0sk8kQFRWFHj16mI7FxcUBAPLy8gAYipQ3b95sqhkKCQlB586dAQCnTp2yeo/KykqoVCqrBdjm9zfeq6H7P/jgg3jmmWcwePBgLFmyBAcOHKj3mIGBgaYeKiJqfgw+ROS1FAqFxecSicTimDGs6PV6AEBZWRnGjBmDffv2WbydOHECQ4cOtXqP6OhoVFRUQKPRNHh/470auv+MGTNw+vRpTJw4EQcPHkT//v3xxhtvWDxmYWEhYmJi7HsBiMjlGHyIqMXo27cvDh8+jOTkZHTo0MHiLTg42OrX9O7dG4BhnR1XSExMxP3334+1a9diwYIFeO+99yzOHzp0CH369HHJvYjIcQw+RNRipKeno7CwEHfddRd27dqFU6dO4ZdffsHUqVOh0+msfk1MTAz69u2Lbdu2Nfn+c+fOxS+//IKMjAzs2bMHmzdvRpcuXUznz5w5g6ysLIwcObLJ9yIi5zD4EFGL0aZNG/zxxx/Q6XQYNWoUevTogblz5yIiIgJSqe0/dzNmzMCqVauafH+dTof09HR06dIFo0ePxmWXXYa3337bdP6zzz7DqFGjkJSU1OR7EZFzuIAhEfm9yspKdOrUCZ9//jkGDRrklntoNBp07NgRn376KQYPHuyWexBR49jjQ0R+LzAwEB9//HGDCx02VWZmJh577DGGHiIPY48PERER+Q32+BAREZHfYPAhIiIiv8HgQ0RERH6DwYeIiIj8BoMPERER+Q0GHyIiIvIbDD5ERETkNxh8iIiIyG8w+BAREZHf+H/2p8TCUPwxeQAAAABJRU5ErkJggg==",
      "text/plain": [
       "<Figure size 640x480 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "def try_trn_neuron():\n",
    "    # Set the simulation time step\n",
    "    brainstate.environ.set(dt=0.02 * u.ms)\n",
    "\n",
    "    # Generate stepwise input current\n",
    "    I = braintools.input.section(\n",
    "        values=[0, 0.05, 0],  # Current amplitudes\n",
    "        durations=[50 * u.ms, 200 * u.ms, 100 * u.ms]  # Duration of each segment\n",
    "    ) * u.uA\n",
    "\n",
    "    # Generate simulation time series\n",
    "    times = u.math.arange(I.shape[0]) * brainstate.environ.get_dt()\n",
    "\n",
    "    # Select the neuron type to simulate (HTC here; can be replaced with RTC, IN, TRN)\n",
    "    neu = HTC(1, solver='ind_exp_euler')  # Use independent exponential Euler method\n",
    "    neu.init_state()  # Initialize neuron state\n",
    "\n",
    "    # Run the simulation and record elapsed time\n",
    "    t0 = time.time()\n",
    "    vs = brainstate.transform.for_loop(neu.step_run, times, I)  # Loop to update neuron state\n",
    "    t1 = time.time()\n",
    "    print(f\"Simulation elapsed time: {t1 - t0:.4f} seconds\")  # Output computation time\n",
    "\n",
    "    # Visualize membrane potential changes\n",
    "    plt.plot(\n",
    "        times.to_decimal(u.ms),  # x-axis: time\n",
    "        u.math.squeeze(vs.to_decimal(u.mV))  # y-axis: membrane potential\n",
    "    )\n",
    "    plt.xlabel('Time (ms)')\n",
    "    plt.ylabel('Membrane Potential (mV)')\n",
    "    plt.title('HTC neuron response to step current')\n",
    "    plt.show()\n",
    "\n",
    "\n",
    "# Execute the simulation\n",
    "if __name__ == '__main__':\n",
    "    try_trn_neuron()\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "9f34a308cd614c8c",
   "metadata": {},
   "source": [
    "## Results Interpretation\n",
    "\n",
    "After running the code, you will obtain the membrane potential trace of an HTC neuron under step current stimulation. Typical characteristics are as follows:\n",
    "\n",
    "- 0–50 ms (no current): The neuron remains at rest, with a stable membrane potential around -65 mV, maintained by the balance of leak current and Ih current.\n",
    "\n",
    "- 50–250 ms (0.05 uA current injection):\n",
    "   - The injected current depolarizes the membrane, reaching threshold and triggering action potentials.\n",
    "   - Due to the involvement of T-type calcium channels and IAHP currents, burst firing may occur, which is a hallmark of thalamic neurons participating in rhythmic oscillations.\n",
    "\n",
    "- 250–350 ms (current off):\n",
    "   - The membrane potential gradually returns to the resting level, but afterhyperpolarization may appear due to the sustained action of calcium-activated potassium channels.\n",
    "\n",
    "## Extended Exercises\n",
    "\n",
    "- Replace the neuron type, e.g., `neu = TRN(1)`, to compare how TRN and HTC neurons respond to the same current.\n",
    "- Modify the conductance of T-type calcium channels, e.g., `CaT_HM1992(g_max=1.0 * ...)`, to observe its effect on burst firing.\n",
    "- Extend the simulation time and analyze the spontaneous rhythmic activity of different thalamic neuron types.\n",
    "\n",
    "Through these exercises, you will gain a deeper understanding of thalamic neuron electrophysiological properties and their role in regulating brain rhythms.\n"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.13.11"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}