{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Chapter 5: Point-Neuron Network Models (with PointNet)\n",
    "\n",
    "In this chapter we will create a heterogeneous network of point-model neurons and use the PointNet simulator which will run the network using the NEST simulator. As with the previous BioNet examples will create both a internal recurrently-connected network of different node types, and an external network of \"virtual\" neurons that will drive the firing of the internal neurons. And we'll show how to drive network activity by using a current clamp.\n",
    "\n",
    "PointNet, like BioNet and the other simulators, use the SONATA data format for representing networks, setting up simulations and saving results. Thus the tools used to build and display biophysically detailed networks in the previous chapters will work just the same. \n",
    "\n",
    "Requirements:\n",
    "* bmtk\n",
    "* NEST 2.11+"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 1. Building the network\n",
    "\n",
    "There are two ways of generating a network of point-neurons. Either we can take the existing biophysical network created in the previous chapters and make some minor adjustments to the neuron models being used. Or we can build a new network from scratch using the BMTK Builder.\n",
    "\n",
    "### Converting networks\n",
    "We want to take the BioNet V1 network and change parameters so that the individual neurons are using point models. Luckily there parameters are stored in the node and edge \"types\" csv files, thus we can easily change them with a simple text editor (emacs, vi, sublime-text, etc). Here is an example of the old *V1_node_types.csv*: "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>node_type_id</th>\n",
       "      <th>ei</th>\n",
       "      <th>morphology_file</th>\n",
       "      <th>model_processing</th>\n",
       "      <th>pop_name</th>\n",
       "      <th>location</th>\n",
       "      <th>model_template</th>\n",
       "      <th>model_type</th>\n",
       "      <th>dynamics_params</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>100</td>\n",
       "      <td>e</td>\n",
       "      <td>Scnn1a.swc</td>\n",
       "      <td>aibs_perisomatic</td>\n",
       "      <td>Scnn1a</td>\n",
       "      <td>L4</td>\n",
       "      <td>ctdb:Biophys1.hoc</td>\n",
       "      <td>biophysical</td>\n",
       "      <td>472363762_fit.json</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>101</td>\n",
       "      <td>i</td>\n",
       "      <td>Pvalb.swc</td>\n",
       "      <td>aibs_perisomatic</td>\n",
       "      <td>PV</td>\n",
       "      <td>L4</td>\n",
       "      <td>ctdb:Biophys1.hoc</td>\n",
       "      <td>biophysical</td>\n",
       "      <td>472912177_fit.json</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>102</td>\n",
       "      <td>e</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>LIF_exc</td>\n",
       "      <td>VisL4</td>\n",
       "      <td>nrn:IntFire1</td>\n",
       "      <td>point_process</td>\n",
       "      <td>IntFire1_exc_1.json</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>103</td>\n",
       "      <td>i</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>LIF_inh</td>\n",
       "      <td>VisL4</td>\n",
       "      <td>nrn:IntFire1</td>\n",
       "      <td>point_process</td>\n",
       "      <td>IntFire1_inh_1.json</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "   node_type_id ei morphology_file  model_processing pop_name location  \\\n",
       "0           100  e      Scnn1a.swc  aibs_perisomatic   Scnn1a       L4   \n",
       "1           101  i       Pvalb.swc  aibs_perisomatic       PV       L4   \n",
       "2           102  e             NaN               NaN  LIF_exc    VisL4   \n",
       "3           103  i             NaN               NaN  LIF_inh    VisL4   \n",
       "\n",
       "      model_template     model_type      dynamics_params  \n",
       "0  ctdb:Biophys1.hoc    biophysical   472363762_fit.json  \n",
       "1  ctdb:Biophys1.hoc    biophysical   472912177_fit.json  \n",
       "2       nrn:IntFire1  point_process  IntFire1_exc_1.json  \n",
       "3       nrn:IntFire1  point_process  IntFire1_inh_1.json  "
      ]
     },
     "execution_count": 1,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "import pandas as pd\n",
    "\n",
    "pd.read_csv('sources/chapter05/converted_network/V1_node_types_bionet.csv', sep=' ')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "and here is the *V1_node_types.csv* used for PointNet:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>node_type_id</th>\n",
       "      <th>ei</th>\n",
       "      <th>morphology_file</th>\n",
       "      <th>model_processing</th>\n",
       "      <th>pop_name</th>\n",
       "      <th>location</th>\n",
       "      <th>model_template</th>\n",
       "      <th>model_type</th>\n",
       "      <th>dynamics_params</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>100</td>\n",
       "      <td>e</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>Scnn1a</td>\n",
       "      <td>L4</td>\n",
       "      <td>nest:iaf_psc_alpha</td>\n",
       "      <td>point_process</td>\n",
       "      <td>472363762_point.json</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>101</td>\n",
       "      <td>i</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>PV</td>\n",
       "      <td>L4</td>\n",
       "      <td>nest:iaf_psc_alpha</td>\n",
       "      <td>point_process</td>\n",
       "      <td>472912177_point.json</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>102</td>\n",
       "      <td>e</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>LIF_exc</td>\n",
       "      <td>VisL4</td>\n",
       "      <td>nest:iaf_psc_alpha</td>\n",
       "      <td>point_process</td>\n",
       "      <td>IntFire1_exc_point.json</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>103</td>\n",
       "      <td>i</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>LIF_inh</td>\n",
       "      <td>VisL4</td>\n",
       "      <td>nest:iaf_psc_alpha</td>\n",
       "      <td>point_process</td>\n",
       "      <td>IntFire1_inh_point.json</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "   node_type_id ei  morphology_file  model_processing pop_name location  \\\n",
       "0           100  e              NaN               NaN   Scnn1a       L4   \n",
       "1           101  i              NaN               NaN       PV       L4   \n",
       "2           102  e              NaN               NaN  LIF_exc    VisL4   \n",
       "3           103  i              NaN               NaN  LIF_inh    VisL4   \n",
       "\n",
       "       model_template     model_type          dynamics_params  \n",
       "0  nest:iaf_psc_alpha  point_process     472363762_point.json  \n",
       "1  nest:iaf_psc_alpha  point_process     472912177_point.json  \n",
       "2  nest:iaf_psc_alpha  point_process  IntFire1_exc_point.json  \n",
       "3  nest:iaf_psc_alpha  point_process  IntFire1_inh_point.json  "
      ]
     },
     "execution_count": 2,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "pd.read_csv('sources/chapter05/converted_network/V1_node_types.csv', sep=' ')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Changes:\n",
    "* **model_type** - PointNet will not support the \"biophysical\" model_type and only support \"point_process\" neuron models.\n",
    "* **model_template** - nrn:IntFire1 and ctdb:Biophys1.hoc are special directives for running NEURON based models. Instead we replaced them with the \"nest:\\<nest-model\\>\" directive (note we can replace iaf_psc_alpha with any valid NEST model).\n",
    "* **dynamics_params** - We have new json parameters files for the new NEST based models.\n",
    "* **model_processing** - \"aibs_perisomatic\" is a special command for adjusting the morphology of biophysical models, and since our NEST-based models do not have a morphology we set it to none which tells the simulator to use the models as-is (note: you can implement custom model_processing functions for PointNet that will be explained later)."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We must also adjust the *edges_types.csv* files:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {
    "scrolled": false
   },
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>edge_type_id</th>\n",
       "      <th>target_query</th>\n",
       "      <th>source_query</th>\n",
       "      <th>syn_weight</th>\n",
       "      <th>dynamics_params</th>\n",
       "      <th>distance_range</th>\n",
       "      <th>delay</th>\n",
       "      <th>target_sections</th>\n",
       "      <th>weight_function</th>\n",
       "      <th>model_template</th>\n",
       "      <th>weight_sigma</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>100</td>\n",
       "      <td>pop_name=='Scnn1a'</td>\n",
       "      <td>ei=='e'</td>\n",
       "      <td>50.00000</td>\n",
       "      <td>ExcToExc.json</td>\n",
       "      <td>NaN</td>\n",
       "      <td>2.0</td>\n",
       "      <td>NaN</td>\n",
       "      <td>gaussianLL</td>\n",
       "      <td>static_synapse</td>\n",
       "      <td>50.0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>101</td>\n",
       "      <td>pop_name=='LIF_exc'</td>\n",
       "      <td>ei=='e'</td>\n",
       "      <td>50.00000</td>\n",
       "      <td>instanteneousExc.json</td>\n",
       "      <td>NaN</td>\n",
       "      <td>2.0</td>\n",
       "      <td>NaN</td>\n",
       "      <td>gaussianLL</td>\n",
       "      <td>static_synapse</td>\n",
       "      <td>50.0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>102</td>\n",
       "      <td>model_type=='biophysical'&amp;ei=='i'</td>\n",
       "      <td>ei=='i'</td>\n",
       "      <td>50.00000</td>\n",
       "      <td>InhToInh.json</td>\n",
       "      <td>NaN</td>\n",
       "      <td>2.0</td>\n",
       "      <td>NaN</td>\n",
       "      <td>wmax</td>\n",
       "      <td>static_synapse</td>\n",
       "      <td>NaN</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>103</td>\n",
       "      <td>model_type=='point_process'&amp;ei=='i'</td>\n",
       "      <td>ei=='i'</td>\n",
       "      <td>50.00000</td>\n",
       "      <td>instanteneousInh.json</td>\n",
       "      <td>NaN</td>\n",
       "      <td>2.0</td>\n",
       "      <td>NaN</td>\n",
       "      <td>wmax</td>\n",
       "      <td>static_synapse</td>\n",
       "      <td>NaN</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <td>104</td>\n",
       "      <td>model_type=='biophysical'&amp;ei=='e'</td>\n",
       "      <td>ei=='i'</td>\n",
       "      <td>50.00000</td>\n",
       "      <td>InhToExc.json</td>\n",
       "      <td>NaN</td>\n",
       "      <td>2.0</td>\n",
       "      <td>NaN</td>\n",
       "      <td>wmax</td>\n",
       "      <td>static_synapse</td>\n",
       "      <td>NaN</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>5</th>\n",
       "      <td>105</td>\n",
       "      <td>model_type=='point_process'&amp;ei=='e'</td>\n",
       "      <td>ei=='i'</td>\n",
       "      <td>30.00000</td>\n",
       "      <td>instanteneousInh.json</td>\n",
       "      <td>NaN</td>\n",
       "      <td>2.0</td>\n",
       "      <td>NaN</td>\n",
       "      <td>wmax</td>\n",
       "      <td>static_synapse</td>\n",
       "      <td>NaN</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>6</th>\n",
       "      <td>106</td>\n",
       "      <td>pop_name=='PV'</td>\n",
       "      <td>ei=='e'</td>\n",
       "      <td>0.00035</td>\n",
       "      <td>ExcToInh.json</td>\n",
       "      <td>NaN</td>\n",
       "      <td>2.0</td>\n",
       "      <td>NaN</td>\n",
       "      <td>wmax</td>\n",
       "      <td>static_synapse</td>\n",
       "      <td>NaN</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>7</th>\n",
       "      <td>107</td>\n",
       "      <td>pop_name=='LIF_inh'</td>\n",
       "      <td>ei=='e'</td>\n",
       "      <td>50.00000</td>\n",
       "      <td>instanteneousExc.json</td>\n",
       "      <td>NaN</td>\n",
       "      <td>2.0</td>\n",
       "      <td>NaN</td>\n",
       "      <td>wmax</td>\n",
       "      <td>static_synapse</td>\n",
       "      <td>NaN</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "   edge_type_id                         target_query source_query  syn_weight  \\\n",
       "0           100                   pop_name=='Scnn1a'      ei=='e'    50.00000   \n",
       "1           101                  pop_name=='LIF_exc'      ei=='e'    50.00000   \n",
       "2           102    model_type=='biophysical'&ei=='i'      ei=='i'    50.00000   \n",
       "3           103  model_type=='point_process'&ei=='i'      ei=='i'    50.00000   \n",
       "4           104    model_type=='biophysical'&ei=='e'      ei=='i'    50.00000   \n",
       "5           105  model_type=='point_process'&ei=='e'      ei=='i'    30.00000   \n",
       "6           106                       pop_name=='PV'      ei=='e'     0.00035   \n",
       "7           107                  pop_name=='LIF_inh'      ei=='e'    50.00000   \n",
       "\n",
       "         dynamics_params  distance_range  delay  target_sections  \\\n",
       "0          ExcToExc.json             NaN    2.0              NaN   \n",
       "1  instanteneousExc.json             NaN    2.0              NaN   \n",
       "2          InhToInh.json             NaN    2.0              NaN   \n",
       "3  instanteneousInh.json             NaN    2.0              NaN   \n",
       "4          InhToExc.json             NaN    2.0              NaN   \n",
       "5  instanteneousInh.json             NaN    2.0              NaN   \n",
       "6          ExcToInh.json             NaN    2.0              NaN   \n",
       "7  instanteneousExc.json             NaN    2.0              NaN   \n",
       "\n",
       "  weight_function  model_template  weight_sigma  \n",
       "0      gaussianLL  static_synapse          50.0  \n",
       "1      gaussianLL  static_synapse          50.0  \n",
       "2            wmax  static_synapse           NaN  \n",
       "3            wmax  static_synapse           NaN  \n",
       "4            wmax  static_synapse           NaN  \n",
       "5            wmax  static_synapse           NaN  \n",
       "6            wmax  static_synapse           NaN  \n",
       "7            wmax  static_synapse           NaN  "
      ]
     },
     "execution_count": 3,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "pd.read_csv('sources/chapter05/converted_network/V1_V1_edge_types.csv', sep=' ')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "* **model_template** has been changed to use a NEST based model type (static_synapse)\n",
    "* Use different **dynamics_parameter** files \n",
    "* It's important to readjust **syn_weight** as values appropiate for NEURON based models are oftern wrong for NEST based models. \n",
    "\n",
    "Notice we don't have to change any of the hdf5 files. The network topology remains the same making it a powerful tool for comparing networks of different levels of resolution."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Building a model from scratch.\n",
    "\n",
    "We can use the BMTK Network Builder to create new network files just for point-based modeling\n",
    "\n",
    "#### V1 Network\n",
    "\n",
    "First lets build a \"V1\" network of 300 cells, split into 4 different populations "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "from bmtk.builder.networks import NetworkBuilder\n",
    "from bmtk.builder.auxi.node_params import positions_columinar\n",
    "\n",
    "net = NetworkBuilder(\"V1\")\n",
    "net.add_nodes(N=80,  # Create a population of 80 neurons\n",
    "              positions=positions_columinar(N=80, center=[0, 50.0, 0], max_radius=30.0, height=100.0),\n",
    "              pop_name='Scnn1a', location='VisL4', ei='e',  # optional parameters\n",
    "              model_type='point_process',  # Tells the simulator to use point-based neurons\n",
    "              model_template='nest:iaf_psc_alpha',  # tells the simulator to use NEST iaf_psc_alpha models\n",
    "              dynamics_params='472363762_point.json'  # File containing iaf_psc_alpha mdoel parameters\n",
    "             )\n",
    "\n",
    "net.add_nodes(N=20, pop_name='PV', location='VisL4', ei='i',\n",
    "              positions=positions_columinar(N=20, center=[0, 50.0, 0], max_radius=30.0, height=100.0),\n",
    "              model_type='point_process',\n",
    "              model_template='nest:iaf_psc_alpha',\n",
    "              dynamics_params='472912177_point.json')\n",
    "\n",
    "net.add_nodes(N=200, pop_name='LIF_exc', location='L4', ei='e',\n",
    "              positions=positions_columinar(N=200, center=[0, 50.0, 0], min_radius=30.0, max_radius=60.0, height=100.0),\n",
    "              model_type='point_process',\n",
    "              model_template='nest:iaf_psc_alpha',\n",
    "              dynamics_params='IntFire1_exc_point.json')\n",
    "\n",
    "net.add_nodes(N=100, pop_name='LIF_inh', location='L4', ei='i',\n",
    "              positions=positions_columinar(N=100, center=[0, 50.0, 0], min_radius=30.0, max_radius=60.0, height=100.0),\n",
    "              model_type='point_process',\n",
    "              model_template='nest:iaf_psc_alpha',\n",
    "              dynamics_params='IntFire1_inh_point.json')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We can now go ahead and created synaptic connections then build and save our network."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<bmtk.builder.connection_map.ConnectionMap at 0x7f1c796b4130>"
      ]
     },
     "execution_count": 5,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "from bmtk.builder.auxi.edge_connectors import distance_connector\n",
    "\n",
    "## E-to-E connections\n",
    "net.add_edges(source={'ei': 'e'}, target={'pop_name': 'Scnn1a'},\n",
    "              connection_rule=distance_connector,\n",
    "              connection_params={'d_weight_min': 0.0, 'd_weight_max': 0.34, 'd_max': 300.0, 'nsyn_min': 3, 'nsyn_max': 7},\n",
    "              syn_weight=3.0,\n",
    "              delay=2.0,\n",
    "              dynamics_params='ExcToExc.json',\n",
    "              model_template='static_synapse')\n",
    "\n",
    "net.add_edges(source={'ei': 'e'}, target={'pop_name': 'LIF_exc'},\n",
    "              connection_rule=distance_connector,\n",
    "              connection_params={'d_weight_min': 0.0, 'd_weight_max': 0.34, 'd_max': 300.0, 'nsyn_min': 3, 'nsyn_max': 7},\n",
    "              syn_weight=3.0,\n",
    "              delay=2.0,\n",
    "              dynamics_params='instanteneousExc.json',\n",
    "              model_template='static_synapse')\n",
    "\n",
    "\n",
    "### Generating I-to-I connections\n",
    "net.add_edges(source={'ei': 'i'}, target={'pop_name': 'PV'},\n",
    "              connection_rule=distance_connector,\n",
    "              connection_params={'d_weight_min': 0.0, 'd_weight_max': 1.0, 'd_max': 160.0, 'nsyn_min': 3, 'nsyn_max': 7},\n",
    "              syn_weight=-3.0,\n",
    "              delay=2.0,\n",
    "              dynamics_params='InhToInh.json',\n",
    "              model_template='static_synapse')\n",
    "\n",
    "net.add_edges(source={'ei': 'i'}, target={'ei': 'i', 'pop_name': 'LIF_inh'},\n",
    "              connection_rule=distance_connector,\n",
    "              connection_params={'d_weight_min': 0.0, 'd_weight_max': 1.0, 'd_max': 160.0, 'nsyn_min': 3, 'nsyn_max': 7},\n",
    "              syn_weight=-3.0,\n",
    "              delay=2.0,\n",
    "              dynamics_params='instanteneousInh.json',\n",
    "              model_template='static_synapse')\n",
    "\n",
    "### Generating I-to-E connections\n",
    "net.add_edges(source={'ei': 'i'}, target={'ei': 'e', 'pop_name': 'Scnn1a'},\n",
    "              connection_rule=distance_connector,\n",
    "              connection_params={'d_weight_min': 0.0, 'd_weight_max': 1.0, 'd_max': 160.0, 'nsyn_min': 3, 'nsyn_max': 7},\n",
    "              syn_weight=-3.0,\n",
    "              delay=2.0,\n",
    "              dynamics_params='InhToExc.json',\n",
    "              model_template='static_synapse')\n",
    "\n",
    "net.add_edges(source={'ei': 'i'}, target={'ei': 'e', 'pop_name': 'LIF_exc'},\n",
    "              connection_rule=distance_connector,\n",
    "              connection_params={'d_weight_min': 0.0, 'd_weight_max': 1.0, 'd_max': 160.0, 'nsyn_min': 3, 'nsyn_max': 7},\n",
    "              syn_weight=-3.0,\n",
    "              delay=2.0,\n",
    "              dynamics_params='instanteneousInh.json',\n",
    "              model_template='static_synapse')\n",
    "\n",
    "### Generating E-to-I connections\n",
    "net.add_edges(source={'ei': 'e'}, target={'pop_name': 'PV'},\n",
    "              connection_rule=distance_connector,\n",
    "              connection_params={'d_weight_min': 0.0, 'd_weight_max': 0.26, 'd_max': 300.0, 'nsyn_min': 3, 'nsyn_max': 7},\n",
    "              syn_weight=3.0,\n",
    "              delay=2.0,\n",
    "              dynamics_params='ExcToInh.json',\n",
    "              model_template='static_synapse')\n",
    "\n",
    "\n",
    "net.add_edges(source={'ei': 'e'}, target={'pop_name': 'LIF_inh'},\n",
    "              connection_rule=distance_connector,\n",
    "              connection_params={'d_weight_min': 0.0, 'd_weight_max': 0.26, 'd_max': 300.0, 'nsyn_min': 3, 'nsyn_max': 7},\n",
    "              syn_weight=3.0,\n",
    "              delay=2.0,\n",
    "              dynamics_params='instanteneousExc.json',\n",
    "              model_template='static_synapse')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [],
   "source": [
    "net.build()\n",
    "net.save_nodes(output_dir='sim_ch05/network')\n",
    "net.save_edges(output_dir='sim_ch05/network')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "### Building external network\n",
    "\n",
    "Next we want to create an external network of \"virtual cells\" with spike-trains that will synapse onto our V1 cells and drive activity. We will call this external network \"LGN\" and contains 500 excitatory cells."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [],
   "source": [
    "lgn = NetworkBuilder('LGN')\n",
    "lgn.add_nodes(N=500, pop_name='tON', potential='exc', model_type='virtual')\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We will use a special function for setting the number of synapses between the LGN --> V1 cells. The select_source_cells function will be called during the build process."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<bmtk.builder.connection_map.ConnectionMap at 0x7f1c442da0a0>"
      ]
     },
     "execution_count": 8,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "import numpy as np\n",
    "\n",
    "def select_source_cells(sources, target, nsources_min=10, nsources_max=30, nsyns_min=3, nsyns_max=12):\n",
    "    total_sources = len(sources)\n",
    "    nsources = np.random.randint(nsources_min, nsources_max)\n",
    "    selected_sources = np.random.choice(total_sources, nsources, replace=False)\n",
    "    syns = np.zeros(total_sources)\n",
    "    syns[selected_sources] = np.random.randint(nsyns_min, nsyns_max, size=nsources)\n",
    "    return syns\n",
    "\n",
    "lgn.add_edges(source=lgn.nodes(), target=net.nodes(pop_name='Scnn1a'),\n",
    "              iterator='all_to_one',\n",
    "              connection_rule=select_source_cells,\n",
    "              connection_params={'nsources_min': 10, 'nsources_max': 25},\n",
    "              syn_weight=15.0,\n",
    "              delay=2.0,\n",
    "              dynamics_params='ExcToExc.json',\n",
    "              model_template='static_synapse')\n",
    "\n",
    "lgn.add_edges(source=lgn.nodes(), target=net.nodes(pop_name='PV1'),\n",
    "              connection_rule=select_source_cells,\n",
    "              connection_params={'nsources_min': 15, 'nsources_max': 30},\n",
    "              iterator='all_to_one',\n",
    "              syn_weight=5.0,\n",
    "              delay=2.0,\n",
    "              dynamics_params='ExcToInh.json',\n",
    "              model_template='static_synapse')\n",
    "\n",
    "lgn.add_edges(source=lgn.nodes(),  target=net.nodes(pop_name='LIF_exc'),\n",
    "              connection_rule=select_source_cells,\n",
    "              connection_params={'nsources_min': 10, 'nsources_max': 25},\n",
    "              iterator='all_to_one',\n",
    "              syn_weight= 10.0,\n",
    "              delay=2.0,\n",
    "              dynamics_params='instanteneousExc.json',\n",
    "              model_template='static_synapse')\n",
    "\n",
    "lgn.add_edges(source=lgn.nodes(),  target=net.nodes(pop_name='LIF_inh'),\n",
    "              connection_rule=select_source_cells,\n",
    "              connection_params={'nsources_min': 15, 'nsources_max': 30},\n",
    "              iterator='all_to_one',\n",
    "              syn_weight=10.0,\n",
    "              delay=2.0,\n",
    "              dynamics_params='instanteneousExc.json',\n",
    "              model_template='static_synapse')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Finally we build and save our lgn network."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [],
   "source": [
    "lgn.build()\n",
    "lgn.save_nodes(output_dir='sim_ch05/network')\n",
    "lgn.save_edges(output_dir='sim_ch05/network')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 2. Setting up PointNet Environment\n",
    "\n",
    "#### Directory Structure\n",
    "\n",
    "Before running a simulation, we will need to create the runtime environment, including parameter files, run-script and configuration files. If using the tutorial these files will already be in place. Otherwise we can use a command-line:\n",
    "```bash\n",
    "$ python -m bmtk.utils.sim_setup        \\\n",
    "    --network sim_ch05/network/         \\\n",
    "    --include-examples                  \\\n",
    "    --tstop 3000.0                      \\\n",
    "    pointnet sim_ch05/\n",
    "\n",
    "```\n",
    "\n",
    "or"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [],
   "source": [
    "from bmtk.utils.sim_setup import build_env_pointnet\n",
    "\n",
    "build_env_pointnet(base_dir='sim_ch05',      \n",
    "                   network_dir='sim_ch05/network',\n",
    "                   tstop=3000.0, \n",
    "                   dt=0.01,\n",
    "                   include_examples=True,         # Copies components files\n",
    "                  )"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The network files are written to **circuit_config.json** and the simulation parameters are set in **simulation_config**. The simulation time is set to run for 3000.0 ms (tstop). We also specify a membrane-report to record V_m property of 4 cells (gids 0, 80, 100, 300 - one from each cell-type). In general, all the parameters needed to setup and start a simulation are found in the config files, and adjusting network/simulation conditions can be done by editing these json files in a text editor."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### lgn input\n",
    "\n",
    "We need to provide our LGN external network cells with spike-trains so they can activate our recurrent network. Previously we showed how to do this by generating csv files. We can also use NWB files, which are a common format for saving electrophysiological data in neuroscience.\n",
    "\n",
    "We can use any NWB file generated experimentally or computationally, but for this example we will use a preexsting one. First download the file:\n",
    "```bash\n",
    " $ cd sim_ch05\n",
    " $ wget https://github.com/AllenInstitute/bmtk/raw/develop/docs/examples/spikes_inputs/lgn_spikes.h5\n",
    "```\n",
    "Then we must edit the simulation_config.json file to tell the simulator to find the nwb file and which network to associate it with.\n",
    "\n",
    "```json\n",
    "\n",
    "\"inputs\": {\n",
    "    \"LGN_spikes\": {\n",
    "        \"input_type\": \"spikes\",\n",
    "        \"module\": \"h5\",\n",
    "        \"input_file\": \"$BASE_DIR/lgn_spikes.h5\",\n",
    "        \"node_set\": \"LGN\"\n",
    "    }\n",
    "},\n",
    "```\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 3. Running the simulation\n",
    "\n",
    "The call to sim_setup created a file run_pointnet.py which we can run directly in a command line:\n",
    "```bash\n",
    "$ python run_pointnet.py config.json\n",
    "```\n",
    "or if you have mpi setup:\n",
    "\n",
    "```bash\n",
    "$ mpirun -np $NCORES python run_pointnet.py config.json\n",
    "```\n",
    "\n",
    "Or we can run it directly"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "2020-04-27 14:24:33,004 [INFO] Created log file\n"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "INFO:NestIOUtils:Created log file\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "2020-04-27 14:24:33,026 [INFO] Batch processing nodes for LGN/0.\n"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "INFO:NestIOUtils:Batch processing nodes for LGN/0.\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "2020-04-27 14:24:33,048 [INFO] Batch processing nodes for V1/0.\n"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "INFO:NestIOUtils:Batch processing nodes for V1/0.\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "2020-04-27 14:24:33,096 [INFO] Setting up output directory\n"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "INFO:NestIOUtils:Setting up output directory\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "2020-04-27 14:24:33,099 [INFO] Building cells.\n"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "INFO:NestIOUtils:Building cells.\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "2020-04-27 14:24:33,116 [INFO] Building recurrent connections\n"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "INFO:NestIOUtils:Building recurrent connections\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "2020-04-27 14:24:33,167 [INFO] Build virtual cell stimulations for LGN_spikes\n"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "INFO:NestIOUtils:Build virtual cell stimulations for LGN_spikes\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "2020-04-27 14:24:59,771 [INFO] Network created.\n"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "INFO:NestIOUtils:Network created.\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "2020-04-27 14:24:59,790 [INFO] Starting Simulation\n"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "INFO:NestIOUtils:Starting Simulation\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "2020-04-27 14:25:01,560 [INFO] Simulation finished, finalizing results.\n"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "INFO:NestIOUtils:Simulation finished, finalizing results.\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "2020-04-27 14:25:01,581 [INFO] Done.\n"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "INFO:NestIOUtils:Done.\n"
     ]
    }
   ],
   "source": [
    "from bmtk.simulator import pointnet\n",
    "\n",
    "configure = pointnet.Config.from_json('sim_ch05/simulation_config.json')\n",
    "configure.build_env()\n",
    "network = pointnet.PointNetwork.from_config(configure)\n",
    "sim = pointnet.PointSimulator.from_config(configure, network)\n",
    "sim.run()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 4. Analyzing results\n",
    "\n",
    "Results of the simulation, as specified in the config, are saved into the output directory. Using the analyzer functions, we can do things like plot the raster plot"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYUAAAEGCAYAAACKB4k+AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4yLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+j8jraAAAgAElEQVR4nO2de3ReVZn/vw9J2vRHbVosVQqWhlK1lh+1pUsU8pvFUAQG8LrU4ec4MLZDB1QqM0ukODqEWS6o46DS8TaFVEEHGW9rQC5yiRcIP6m0lAAlQtqGFlskLZBgNEnT9Pn98Z7z9uTk3M/eZ+/zvs9nray873nP2fvZz372efZ9EzNDEARBEADgCNMCCIIgCPYgTkEQBEGoIk5BEARBqCJOQRAEQagiTkEQBEGo0mhagLzMnj2b58+fb1oMQRCEUrFly5b9zHy0/3rpncL8+fOxefNm02IIgiCUCiLaFXRduo8EQRCEKuIUBEEQhCqFOAUiaiCirUR0l/O9lYg2EdF2IvpvIpriXJ/qfN/u/D6/CPkEQRCECkW1FD4NoMfz/UsAvsrMJwJ4FcAq5/oqAK8617/q3CcIgiAUhHanQETHATgfwM3OdwJwJoAfO7fcAuD9zuf3Od/h/L7CuV8QBEEogCJaCl8D8FkAh5zvrwcwwMwHne+/B3Cs8/lYAC8AgPP7oHP/BIhoNRFtJqLN+/bt0ym7YCkdXX04fd0v0NHVZ1oUQagptDoFIroAQD8zb1EZLjNvYOblzLz86KMnTbMV6oCNXX3YMzCMjeIUaoqyOntb5c4il+6WwukA3ktEzwO4HZVuoxsBzCQid43EcQD2OJ/3AHgTADi/twB4WbOMQglZ2daKY2dOw8q2VtOiCAopq7O3Ve4scml1Csx8NTMfx8zzAVwI4BfM/DcAfgngQ85tFwO4w/l8p/Mdzu+/YDnwQQhgVVsrHll7JlaJU6gpyursbZU7i1xU1DuXiM4A8BlmvoCITkCl5XAUgK0APsbMo0TUDOB7AJYCeAXAhcy8Myrc5cuXs6xoFgRBSAcRbWHm5f7rhS1eY+ZfMfMFzuedzPwOZj6RmT/MzKPO9RHn+4nO75EOQRDKhq19z1koY1pUy5w0vKD7bNWfrGhWhJvBq2/dbGVGF4GtRm6TXG4f73V39xQqjw4d2NqPHoVqmZOGF3SfrfoTp6AIN4M7e/oLzWibXnjrO3uxZ2AY6zt7TYsyAZsK38q2VjQQYZy5UHl06GBlWyuaG4/A3sFhrL610oXb0dWHJe33Y8m191tZK1bd9580vMVzZ6CBCIvnztAmiyrEKSjCzeAVi+YUmtE2vfDAvv+WYFPhW9XWis+dv6hweXToYFVbK8bGGcxAZ08/gIo9Do6MYXB4zMpaseoJCknD27b3NYwzY9ve17TJoorSb51tC6vaWo1k7sq2Vmzs6rPihbfmrIXWyOLFVN6EYUIeXXGuWDQHnT39WLFoDoCKPa5/sBcgTLADm+zUBK5ehkYPoqOrzyp79FPY7CNdyOyjStPcLXA2G1sYXvkBlDotQm3S0dVX6RblSuVnVVtr6nJ3+rpfYM/AMI6dOQ2PrD1zQthR4egq38ZnHwn6sKVpnhWv/GVPi1CbbOzqw+DwGAZHDneLpbXVsC68uHCKLhPiFGoAm/rMs+CVv+xpEWqTlW2taJnWhJbmpqptprXVsDGEsHDcwfnFc2cUWiak+yiEsnfJqEb0UT5U5Fkt5XvZ0hLW3aQK6T5KiXRjTET0UT5U5Fkt5XvZ0mKq1SxOwYN3LnXQvOKiZbBlbjfgNJ+bm6qzJwQ7iLIRFS+VuDCCFm3aZLdeVL9k/elUnW5TU1al+8iDt7kGQGvTzWYZwtDdnK03VHRnmM4TN353QZ6NdqsLv+5N50VapPsoATYMeNogQxi2yVN2VHRnmM6ToEWbpmUqCn86ayXd0lIoCWUbJBPikTwVTNqAtBRKTtkGyYR4bN3mQCgOG8u1OIWSUCtNU0GoV4IGom0s1+IULCFu5oLUKvWic8ZMkrB1xx+0c2ncM0Xrw9ZZSy55ZxvZ2CoIQpyCJbgGs/7BXqsLRpmJKsT+cw5UvqCSvAx0vjDCdi6Ne0anPEFhp4lTVf6kCccvn4ptLrLqWefUdXEKluAaDAilqE2UkagC6D/nQOVLMUkXgc5uBHeNScu0psTh65YnKOw0carKnzTh5J1tFNTaz6pnnfuFyewjy5AZKfpIsxslILu12oyqclLW8qbCVsNmH4lTsJiOrr7q3vRrViwEIC8qQQ95X45lfbnWMzIlNQW2DHj5+4JtPe5SsI+iB0FtH0S1pUxHYYuM4hQCsMXAJ/UFW3rcpWAfqvb6T4qNUyu96CjTqgd7s8iow5HIcZwB2HJ0YNARijbIJdhPWhvOe1ynbUee+tFRpv0vcfdzVj1kkdErgyr918yYgvRpCkkQO1Gvg3rVadbBXpX6yhNWzY8p2NLlkwQTfYemF1D540m7mEoVXjsJSq+teaOSoLKSRwYVZc+W/vQ0cninmKZZXJpHX375dCxqrRmnYHufphdVDizPwhudciWRJe1iqijS6MFrJ0HpVa2DJLKljTNvX7bKRVRh4XllTVIBsMWx6CoDXtnyvKtUOpQwGtrb21MHbhMbNmxoX716NZbNm4VVba1YNm9W9beOrj5cfttWAMCyebMmfc+LG96mnS/jyh9144YHnsO9T72Ib/96Z2wcO/qHsNInb1ouv20r9gwMY0f/UKKaQpI4VciVhO7dA2huasClZyzIHVcaPfjtJCi9WXQQZluubI/07se2vYNYd+/vJt2TNk5verfuHkhlA1HdDVnzPqjseWXdNzSK0YOHsKN/qHrNfU5F/N64vLrw5sm3f7UDV9z+BLbtHcR7lsyNDKd79wAOjB/ClMYjlJUDr2xL583Cw7378dBz+3DTw32p48mqJ68MAHDvXXfM/sKVV3zRf1/NjCkEEXYIRktzE6Y3N+bu0/MfMOKliIM2yji3XEecNugh7ICVjq4+XHd3D8aZJxxEk8c28ixcKvogGP9aG7emqyN+fz560/qHwZFqHuy4/rzIcLzPeQd/VS2Sc3XgkkQXKmzUL8PmG/8Boy/2kv++0juF6ce9hW+8/eeRK1QXz52BbXtfq/4fGj2IweExZYVz8dwZeOi5fRgdP4Q3z5mOodHxUgy6mTgpysbTqVTIFFRow+zPpG2YHhQuMn5vXJt2vozOnn6sWDQHGy6aNLYa+pwOJ9bRVVlzNDo2jqlNDVizYuEkXUQ5OBVydHT14bIPnnngQH/f1Ek/MnOp/6a88UQ+7fpOjuK06zv5+Kvuqt5388M7+bTrO/nmh3dGPlcWsqYnrx6yPG9K91Hx6pLJb3f1yiW3PMYnrL2bL7nlMdOiMHO6/FZhG1nCKOKdBWAzB7xTa7ql4FKrU/BcOVS1fNJiY60/DBOy2mInpllw9T2Ju25009HVhy/e9QwYQEtzE7rbz9YeZxbbK8J2jExJJaJmIvotEXUT0TYiuta53kpEm4hoOxH9NxFNca5Pdb5vd36fHxfHW9/4ulilqZ62Zcv012rfJMP4edJpKXoKYlJZVcolZ2BUWLFoDhqIsGLRHGVhZs2njV19hzcEmNSbrkeeLOXEpO1obSkQEQE4kpmHiKgJQBeATwP4JwA/ZebbiejbALqZ+VtE9AkAJzPzpUR0IYAPMPNfR8VhYkO81bduTtw/qZMy10RtbWXYKpcwkaz55Pbng4E1Z03uyy9aHpMYaSk4XVdDztcm548BnAngx871WwC83/n8Puc7nN9XOI5FKXlrg9v2voZxZmzb+5piySokla/MNVFb15XYKpcJim7NZV1vkoZVba3ovuZsdLefrbTc1JLdaB9TIKIGAFsAnAjgGwC+DOBRZj7R+f1NAO5l5pOI6GkA5zLz753fdgA4lZn3+8JcDWA1AMybN++UXbt2pZIpr1fXXUMvY61DqD2KtkOx+2Ixts0FM48z89sBHAfgHQDeqiDMDcy8nJmXH3300amfz+rV3ZoMAK01dBtrHbZsQyAUR5Ad5rUD1X3vNlErZaTQ2UdE9C8AhgFcBeCNzHyQiN4FoJ2ZzyGi+5zPvyGiRgB/AHA0RwhZ5JiCrbMIikBqcfkQO1DzvM2ULW2mZh8dTUQznc/TALwbQA+AXwL4kHPbxQDucD7f6XyH8/svohxC0WSpydgyUykvZa/FmUbsQM3zNlMradM9++hkVAaOG1BxQD9k5n8lohMA3A7gKABbAXyMmUeJqBnA9wAsBfAKgAuZeWdUHLYfx1krNcQgajltqhFd2YGt+WBCLjmjuQBsNThd2NZcrjf9ly29Nshrm826mJCr5s9TsIFa6SJIim3N5XrTf9nSa4O8ttmsi01yyXGcCrHlGM+i8B/BaLomWG/6L1t6bZDX1mNDbZJLuo8EZdjaNBcEYTLSfSRox6YmsCAI2RCnoBmTC1qKjrsM227UU34Uha4FbUXoy6Y8sUUWcQqaMTm4VkTcriGvvnXzJIM2YeRBcXqv1Xp+AJN1oDsf8qYr7HnV5zYH2aruPInTvS226UWcgmZMdqkUEbdryJ09/ZMM2oSRB8XpvVbr+QFM1oHufEibLv+LMux5Ffrypj3IVnXnSZzu/bbZ0tyEodGDRlsLMtAspMI/w8g9g3d0fBxTGxombEds8gxo79GXQLpzjFXKYWImVlAemV4f4KXICQlB51kXeSxqnO51H7sZhSxeE5QQZLQ2zjoyLZPp+G3GNidlE0XqRmYfFUCavltbBpXSEtTcNtklE6bHLDKpzBOZiXUYv17LMCHBFDboRloKCklTO5SapBpU6lHyRA+iVzuRlkIBpKkdSk1SDSr1KHmiB9FruZCWgiK8A5ybdr4CELBmhbozYNPKYesgYz1hg+6DBt5N2oENA/CiiwrSUtCMd7rb4MgYBofHsNEzN1rH2EFQ2EVPR8wiY5770t6bFB1hptG9LjuJmjKcFhUy2rBO5P5nXpokg4kxPlee6+7usWpsUZxCBoIMyG0ir1g0By3NTWiZ1lTdAExXIfCG7cq0eO6MCU110033pOlPoycdOtURZhrd67ITr12qnPOfVx5T60TI+UzOd5f1nb3YMzCM9Z29hcrTQIRxZuML1rzILqkZ8BYOt9kXtcuhrp0hvbtOujIBmDCYZ3r3xaQ7Y6bZQVPHbps6wkyje107iKrM/ywy+rtITNrjqrZWrO/sxeDwGGZMa5ooB/v+FyQPoO/94CVNV5WMKWTAdF9gEDbKpIJaTVe9YNvMozB7qnU7C8qHsDGF2JYCES2L+p2ZH88saUlJWtsp0tBMtwh0EdQqU01R+RQVT62+lGw4Q8FLrZaTONLkQ5IxhRucv28A2ARgA4CbnM/fyCFnzWN6kLdIdA3UFbl/k+58ioqnVm0laDGWjQs3TYwpFKmHNIviYp0CM/8lM/8lgBcBLGPm5cx8CoClAPbklraG0flCs61g6Xqpxb1UVOihqMHPqHhMTwgoElW2orQM5BxTyCKLrRWBxGMKRLSNmRfHXSsaW9YpFI3NfbWA3g3ovGkHYJUehHhUdZWpLAN5Zcoii+kuQxXrFJ4kopuJ6Azn7yYAT6oT0X5sqp3bVrNc1dZa7bd0m+JFbNXs14Pth+gUKZ9N9upF1f4+KstAXpmyyGLDPkdBpGkpNAO4DMBfOJceAvAtZh7RJFsiimwp2FY7tw1XPy3NTZje3GikBmQyj5LEXaR8Yq9CFLlbCsw8wsxfZeYPOH9fNe0Qisa22nkYpmqIrn7WnLVQeQ0o6oS3IBlsPUSnSPmKjMvWVomQntiWAhH9kJk/QkRPIWAYhplP1iVcEup1TCGKWqwhumlyV4DWUtpqgVq0uVonT0vh087/CwC8J+BPsIyytGjSoHK7hnqjiP23/NurCNHY3LJStqKZiH7DzO9SElgKbGopmJ5NIAhB6KzFSwshGzborYhdUpsVhlVKbJ13bAM214xqHZ0tR9Mn3JUVm1vzKlsKjzNz5JYYOrClpdDhTMUEY8Lh9UIFG2pGgjm8rWh3U7qWaU3ovuZs06LVLXKegmY2dvVhcHgM05sbxSEEYHPNSNDPhFa0gR1JheSodAo06QLRm4jol0T0DBFtI6JPO9ePIqIHiKjX+T/LuU5EtJ6IthPRk3Gb8dmE96Xnbx6baC7b1kQPW6hjm5wmsFEHqmXylo81Zy2sTl0WklOUnaTqPiKi4wEsZOYHiWgagEZm/qPz20nM/LTv/mMAHMPMjxPR6wBsAfB+AH8H4BVmXkdEawHMYuariOg8AJcDOA/AqQBuZOZTo2SypfvIi7+rxETXSVm6a8oip05s1EF1IeK0JkyfamYhYlkoaoKJajvJ3X1ERJcA+DGA/3QuHQfgf9zf/Q7Bufaiu7W24zx6ABwL4H0AbnFuuwUVRwHn+q1c4VEAMx3HooSiPK2/q8RE10lZumvKIqdObNSBKxMYMnkihqImmBRlJ2m2uXgCwDsAbGLmpc61p5j5fyd8fj4qW2OcBGA3M890rhOAV5l5JhHdBWAdM3c5v3UCuIqZN/vCWg1gNQDMmzfvlF27diVKg401MkGwGZlmHU9ZdZT5kB0Po8x8oPIOB4ioEQmHiohoOoCfALiCmV9zwwAAZmYiSjXkxMwbUDnXAcuXL0/8rG0HfgjqKGvBtJ16PZQmDbWmozQDzb8mos8BmEZE7wbwIwA/i3uIiJpQcQj/xcw/dS6/5HYLOf/7net7ALzJ8/hxUHhmQ1G7EqropkobRtj9urvM4s42UB1/WHiyRqQ86CofWXapzXo2h42TA1SRximsBbAPwFMA/gHAPQA+H/WA0zXUAaCHmb/i+elOABc7ny8GcIfn+kXOLKR3Ahhk5hdTyGgFKl5QacMIu1/3y9IbflBcquMPC8/Gfvl6Ie0LUlf5SBKu/544+00aTi05iTS7pB5i5puY+cPM/CHnc1zXzekA/hbAmUT0hPN3HoB1AN5NRL0AznK+AxVHsxPAdlSO/PxE2gTZgIoXVNowwu7X/bKMOttAR/xh4dm6N309kPYlr6t8ZNmlNs5+k4aT1dHZ6EyS7JIauDuqSy3sklrv/dH1nn4bKVOeJJW1TGlKS9q0ufcPjR7E4PCYkckvYQPNSZzC8c7HTzr/v+f8/xgq48RrlUmZARVOod5nJdV7+m2kFvOkFtOUFRsOpMq8ToGZdzHzLgDvZubPMvNTzt9VAGpi4xLT/dGmm5Bu+hfPnWFdU1YFpvWbBdM2qYOsaSoq/5Ie5KQCnQdS5SXtOoVPMvMjzvfTAHyTmd+uUb5YbFzRnBZbalC2yKGaWk1XvVBU/tXbQU4qNsRbBeCbRPQ8Ee0C8E0AK1UJWHby1GZsqRXqkENHLS/tdERb9CtkQ1X+xdli0EFOquy3TK3V1FtnE1ELADDzoBaJUmJLS0Fqo8Ho0EtQmKJ/IY4sNqLKrmy0TxV7H7UQ0VcAdALoJKIbXAdRT4R5/KJqo7bVOJLWvlTppaOrDy//aRQEYPHcGdVrQyMH0TKtSVoDJcCUDWexRVX2W6bWapoxhZ8AeBqHN7L7WwBLmPmDmmRLRNEtBdMe33T8puVx4wNgdBdaITuSX3agYkxhATNfw8w7nb9rAZygTsRysLKtFS3NTRgaPah1O4eo+MNqHCZqYEHy6JTD1b+3VVDrs6fKTJAt2FhrNtV6sa3lD6RrKfwGwJWeHUxPB/DvzPwujfLFYmJMwdY+bRtkMCmHLemvV4IWcNmYJ7bI2dHVh+vu7jE200lFS+EyAN9wZh89D+DrqOyBVHcUsZ2DKrnqSQ5b0l+vBG31YGOe2CLnxq4+jDOjgcgq/aRpKUwF8CEACwDMBDCIyormf9UnXjy2zD4CansZfx5EL/WBN58BTMrzIu0gKi5b7NG0HJm3ufAE8HMAAwAeBzDuXmfmG1QJmQWbnILOJmhHVx/Wd/YCDKw5a2GpXq42diEIejHdxSo2F4+K7qPjmPlCZv43Zr7B/VMoY6npcDa3amnWMy1yY1cfBofHMDgypm0bbF2DXmmb5kUPvtV6fCqIkjnpYHKRXTQ2b6mh6kwHXaRpKWwA8B/M/JRekdJhoqVgYqCqiJaCLbUrU9NcazU+FUTJXMb0hFFEWvxxeL8DKEyXKloKbQC2ENGzRPQkET1FRE+qE7E8rO/sxZ6B4cpL2kF3LWhVWyu6rzkb3e1na+s6sqVGX/SgX1HbKKiOTzfe9ITJXGsLB4toYag600EXaVoKxwddd3ZQNYaJlsKS9vsxODKGluYmdLfXxEaxmailGqIKak0fSdJTa2nOShn1kLul4G6h7f9TK2Y5WHPWwuq2t7ajs7/ShlqNTeTRh+lxjayLzHRsY1LkNtlSFiaTekM827Bp9pGNmOqvFCYSN/3Q9LiGLTXdorfJNp1ek6gYUxBKiG39lfVK3Bm+psdRbLGNouSwJb02Ii0FQSgA0wuVBMFPWEuh0YQwgr3Iy0sPq0SfQkmQ7qOM2LDIRAdx3Rx+alUPglCviFPISNqXZ1lI29daZj2IQxNsxaRtilPISK0OVK1qa8Uja89M3NVRZj2U2aHVGqYdtBv/6ls3W1FRMGmb4hQSEGSwaV+etqGqEJZVD7W2Ehcw/2INIqlMph20G/8Dz7w0abcCE5isbIlTSIBpg9VBLaYpDRu7+jA4MobpUxtL59DCsDFPk8pkusXpxj+1yXklGp6UabKyJU4hAaYNVge1mKY01GL6bUxTUplMtzjd+K88562l2a1AF7JOQRAEoQ6RFc2CIAhCLOIUBEEQhCriFCzGxtkkRVGGtJdBRqH+yLszslanQEQbiaifiJ72XDuKiB4gol7n/yznOhHReiLa7hzis0yHTGUqyDbMJjGlLxvSHkcZZBT0Y9s7xWuX6x90DgR7MPkUW90the8CONd3bS2ATmZeCKDT+Q4AfwVgofO3GsC3dAiUpCDbksk2zCYx9eIzkfa0+W5D/gjFEWYftlUOJtglORcp8pEJaJ99RETzAdzFzCc5358FcAYzv0hExwD4FTO/hYj+0/n8A/99UeGnnX2UZMO3ethr3dXD4rkzsG3va6H6qKcN8tx8byDC585fpDy99aTLWiTsvWBzvkbJFjb7yIRTGGDmmc5nAvAqM88korsArGPmLue3TgBXMfOkNz4RrUalNYF58+adsmuX2gPgbM5kVXhfgOPMNe0Ak9LR1Yfr7u7Rpo96qGzUMrX2XrBy62xmZiJK7ZWYeQOADUClpaBarnrY5nhlW+uklkK94+a5W/BV4+pcdF1O6uG9ANRh95EgCIJg1+K1OwFc7Hy+GMAdnusXObOQ3glgMM4hCIIgCGrR2n1ERD8AcAaA2UT0ewDXAFgH4IdEtArALgAfcW6/B8B5ALYD+DOAj+uUTRAEQZiMVqfAzP835KcVAfcygE/qlEcQBEGIRlY0C4IgCFXq2ikkXaxky2I2XZQhfUXmlY36yCOTTemxSRZVmEqTrnjr2ikkXYnovS/vviIuWZ/VYQhRegiKz0QhSJtX193dEylfVBpsW6EKpM8j72/X3d2jLT1pbcFNx/rO3ppxDqbsRVe8de0Ukm5T4L3PmxF5MiXo2bAC5r2uwxCi9BDkENd3Ht5PpaiC7cq4eO6MyDhXtrVWF+RF6ShKj1m3r9DpLJPmUdBv48xoINKyPiKpE3Zx0wGGdY43K968KbLCpGubFTlkJyXeVY2bdr6Mzp5+rFg0BxsumjTdN3E47oKYsBWv3uveBVBFLKTxyum+AFqamzC9uRFDowcxODwWu0JX5UrQJKuCk8Tnv0eFjKZWLEfJrnsVbtZV4LW2OtilTKvWjW1zoRuTi9dUG0BYQSmyAEXtiZT1RapSTzpe5l4Z8+x7VIsvuiwOthZJmsYy6UKcggaSbipXJtyXI6FydnnLtCZ0X3N27HOmaqt+h5M1rjz7HnV0VbYoBgFrViw0YgO6dKzCoZt4Uaoum2n0EJReXZWZrHR09eGyD5554EB/31T/b3U9ppAX97DvbXtfM9Y/qroP0+2nnNrkmEbCOkNUv7bOQ9n9/apZx1xWtbXic+cvytRHu7GrD4MjYxgcHjPWR65r0FFFv7WJgVg3zs6e/sBJImlJo4eg9PqvmZ7MsLGrD9TQOCXot4b29vaCxVHLhg0b2levXp3qmY6uPlx+21YAwLJ5s5TIsaN/CCvbWrF190Bs2Crjv/y2rdgzMIwd/UOxL92geP3Xls2bhVVtrZg2pRE7+odw6RkLEsvo6iDufpXpd+X1hpNUjiRhJaV79wCamxoC9ZU2vVn1kzXdUXKoqs2mKR+q2NE/hFNPOAojY4cmjIk90rsf05sbU8Wf1jaC8sJ/TXV+peWu279z4AtXXvFF//W67D7K2xRMGnbYoHCeJnmeZmhQvGkHbgFkelEEDVgn7b833dTOS9ouLp2DlXntJS8mB2LDugjLbF95ZLdpQzzj5G0KugQ1RxfPnYEGIiyeOyP02TxNcn+YSbtmOrr6MDRyEC3TmibEm0QWFdNwvc/FTRv16zVvU9sfnsoutyRhpe3i0nmiWxpdJpEjrS5NnlbndhG2NDdhaPTgJPtKOq3WJnR0Q9Vl91HepmBHVx8+vvExPNDzEl4bOTih62bdvb/D4MgYnt//5wlNV29cebopguRJwuW3bcW+oVHMnj4VX//o4eOvk8rixrl03qxU8btdIYvnzsDI2CEsnjsD9297KVQ3rqz+LrEsaXbjfrh3P/YNjaJ79wC+/+iu6vckXW5xeGV1vwOYlN8AqgU3Tod57cNLUFdUUl0mkSNN92XSMHWybN4sfP/RXZPy/5He/RhnztS1lJYk3bhpyFo27r3rjtlB3Ud16RTCCMqYICN2X7AA0ECENWctnPC7a2AjY4fwyNozMxlYmJHkKVRBxhMWj/e6e7iId8whafzuS8PVxbp7fzfhe1g4XlmzptmNu7mxAbOnT8WB8UPYNzRa/R6ni6QF1ZXVrbUFvSC9L8+vf3RZ5jxM+/Lwv7TT6jLOPlxnn7aSkFR+FcQ5xmXzZmF6c2O13EY5eBUEOdK0ztXFm59pdHv5bVvx4m/vafjCZ9Zc6/9NnIIHN2O6X6jUKIFw5boDi5855y0TMtE1sDWe2koAABLKSURBVLyDSFmNJIywl8HHNz5WrUVfdsaC0PhV1mTidBNl6N7vSQYt3cHyr390GaY0HjHhe1wLJUke+PUalbasNuFNc5TjCSOPLYbpwL2+/48HMKXxCCx1Xq5Zw8tKErtM4hj95TZIz3FxpSkjUQPRWQfj0+o2bKC5rp1CWA3iwMFDkd0LW3cPYNve13DpGQsCX5iqmshJCnPemte3f70DowcPobmpYYJT8MefxODCWloAJoyBZGlpBL2ot+4eiCy4QCWv3BdWklpVnhkiUfmexya8aV7Z1qpMpqSExbejf6ja+krzku9+YQAHDh7ClMYjcpeRpC9CNw3f/tUOXHH7E9i2dxDvWTI3tEXsfcb9HlaBSitLUJ5s3T1QtdUsjt+fziRdg//y2X98or29fYP/t7qcfeQSNhMi7ewQkzMq8satcqVmkm06VM62AibPhPLGBSA03jJtR2Dz7Jgssulc4R7Hgqvvqe4FteP681LJsuTa+zE4PBa6oDNPPpnYyiZs9hGYudR/p5xyCmfl5od38mnXd/LND+/M9VzWcFRgMu6kshQpozeuqHht0lu9YVL3l9zyGJ+w9m6+5JbHUsuiU24TOgGwmQPeqcZf6nn/8jgF3dTCi0dXGvKE6z57yS2P8cnX3Mcnt9+XOJyi80SHY7LlRZaEJPGrkjEunLy/Z407abhFl7Uwp1CT6xSK3L42Kl4Vc4jzpEXFs+422WFp8Mbhfl596+bYeFVsO97Z0596ewl/vLr029HVhyXt9+OLdz0Ten5A0LbkSeQI011QGDrmsSeV1V0sFrdFvCoZ49YUeWVRrauoZ4Nsbsm192NJ+/0Tyo2uLenXP3g43CTUpFMwta+IP17T+8ZkfdZbgMCITEPQwjbvfjNh+HXjdSjeAhP17IpFc9DS3DRpQV4UefdK8r5M4l4EgyNjYFSmLQedH+CVRcWisqAwdCwWSyrrxq7gsxx0lJO4cPyy5NWV36lEPRtkc4PDYxgcGZtQbtyyBlJ81gT5/sdQs7OP4kbhdc2XVjG/3itf2nngXjbtfLm6iO49S+aGxrFp58tYd+/vAFR08fHvPIbhsXEQAVee+9bA6Zte/AvbohamuXhnW0yY3bTvTxg9eAijBw+Fzr5w9fqeJXNx2RkLcJmz31CSPA3KE2+exYWRZiaQd9ryqSe8ftK9flm6dw/gwHj8rJwou4qKQ6XN7+gfwuK5MybYTdh9a86avHNsWDlJurAr6ZqiKFny6CrN+o+g37pfGAABOMTA0nkzMTJ2aNLU6SC7ypKH3qnY3meuvfbaF2X2kQfbZ5+okC8uDO8ZAt79YJa034/BkTG0NDehuz1+22wVsnn3pSECZkxtCnyZpAlThVx+dM4E0m2TqsMvam+krHt2qZbDiwo7yJIGlemWvY98mNyDJQkq5AsKw9vs9XbDeO9bc9ZCHDtzGtactTB3Ovy48S+eO2OSbNOnNqJlWhM+f/7b0N1+durCpktnXnRuA+7dN0sHqm1eRxkKCjPpNZXElR0VdpAlDUW8t+qypRDl5aN+W33rZnT29OPEOUdiaHR8Qh/h4rkzsKnvFYCRuoYbJ6M3jm17Xws9OCTPWgJ/GHFxxcUX9lvStQwdXdEH17jhT5/agO39f5qQJ0l3j42Lo2hUHTUadj+QbIdbvw3E5W9QuG5Z8R5Va/N6iyRkqaUH6RIA1nf2Vt8VQPqdh7Pq0vvc3/+fE6Sl4JJmpoCXzp5+jDPj2ZeGAgdXvYNHKmX0D+CGDeQmGQCMq2kkjSsuvrDfwuIPHIyLmFnkhv/sS0OT8iSpLuLiKJoktcC0A+Pu/dfd3VOdhZJkgDhN3gfJ5JaVzp7+zLLbRpZaepAuwwaa0+glqy6TPFeTA81JBmO8A2XegdaoHSy37R3E8/v/jDe/YTqamxqqg6vduwfAxGhsIBw5pRFL582MHIBLu9GafwA3aiA3z75CwOHB6RPnHInmpobYQeOwQdowPYYNyPm3w3D16j24Jij8uS3NGPjz2IQ82bp7AA8/tz/00Bs3Dx5+bj+IgCOnxudZXrIOgvufD5t4EBW+u9EbEcCM0IkHXlybPnJqY6AO/ZMY/HntlpUVi+ZMiCvPbrdBactCVHhRv4XlT5x8/skXS+fNQvcLA2hurNhn2p2HXbJuF+LGded3bqyfgeakzbywgdakuE2xodGDGBweS7ztRdDvJprWrhwt05owfWrjhKl6WQay8h5eFLaNQJie88jhvyduC4O8RMmkotsvanJAVX8jBzE4Eq+/JPElzeuk3XhR9xQ5QJ6lbKqWL+m7IG+8dTXQnLSZFzbQmhT//GL3+bj4g3430bR25fDOoc8zkJXm2cD0su+/7964NRNp5Jh0T0jcqoibQ5+32y9qLro7KOpOIEiSP1lsOIik3XhR9xQ5QJ6lbKqWL+m7QNugc9Ay5zL9FbnNRdSeR3mXqBe5HYCOcNMu8096TZV8eeTXjW79m8S/11AQtsruolo+ndttpAEh21zUZPdRWtI211qamzC9uXHC/UFdMVm7geJm7iQ92zgrQfHnaUKn6WoIm5WR5LPpGS06uwDDwjbR7eiNc9POlyfNMvLe56490W2zUbLqmGGWR+95uyrjbCHprLHSzD4ionOJ6Fki2k5Ea4uIM21zLWgZelBXjGp5VrZFn22siqD48zSh03Q1uLMy1j/YGzrDJW7miyl0yhIWtuo4k+xr5I0zaJaR975xp9Kp22bD0DXDLJfeU3RVBuVH2F5GaWeNhWHV7CMiagDwcwDnALgewPprr732ofb29n1hz6g6eS3J6L87+yBoGXp19k7fK9VZBaoPNvGfDqVzC4+ok6GSLOX3yuAeXpLkRKmHeytZ7c4cCjobOuyzjllDaUlymEvWk7WiDrvJmn6/raQ9tGZK4xGBs4xcuncPVGd4Bc0iKyLP/LPYggjb8iWKrHr3bjvh2kJYvEH5EXcwVtIZg6WYfURE7wLQzsznON+vBgBmvj7sGRXdR2nxN9/c7y//aRQjY4e0zWDxxw0Ed53o2gIgrsns/T1oFpOOxVneZ/wLhHR2McXJ6R7mAgBfuOBt+PJ9v8PI2CE0Nx2B1x85VVn+BHWPePWxaecrkV0nQYsGk+oqSTdj1HYmpraYCZI770zEoPD93ThR8bokOfArTxnxPlOW2UfHAnjB8/33zrUJENFqItpMRJv37QttRGjD3wRzv4+OHarcoNHPJuk60TUrIa7p6f09SIYkcmXZPiBsgZDOLqa4MFcsmjPhXtc2Rp0anKr8Ceoe8eojruvEL0sa/SfpZvSH73Yp+XdOLZIgufPORAwK39+NExZvA1WmjAXpJCg/8pSRJGWgMXGoFsHMGwBsACothaLjd2vC3imoQTXVIuL2fnZxu2t0xx31e5AMuuVK0lJQHWdYmBsuWh46IKtSDyvbWqsthSB7dFsKYXLmkSVIB/5r/vD9NmKCILlV50mQPYbFC+ifLBFnr16k+0gQBKEOKUv30WMAFhJRKxFNAXAhgDsNyyQIglA3WNV9xMwHiehTAO4D0ABgIzNvMyyWIAhC3WCVUwAAZr4HwD2m5RAEQahHbOs+EgRBEAwiTkEQBEGoYtXsoywQ0T4Au0zLIQiCUDKOZ+aj/RdL7xQEQRAEdUj3kSAIglBFnIIgCIJQRZyCIAiCUEWcgiAIglBFnIIgCIJQRZyCIAiCUEWcgiAIglBFnIIgCIJQRZyCIAiCUEWcgiAIglBFnIIgCIJQRZyCIAiCUEWcgiAIglBFnIIgCIJQRZyCIAiCUEWcgiAIglBFnIIgCIJQRZyCIAiCUKXRtAB5mT17Ns+fP9+0GIIgCKViy5Yt+4POaC69U5g/fz42b95sVoa1d0/4/vy68w1JIgiCkAwi2hV0XbqPBEEQhCqlbykIgiAkxduqlxZ9MNJSEARBEKqIUxAEQRCqSPeRIAhCgdg+MUWcgiAImbH9BZcGGW+oIN1HgiAIQhVpKQiCoAWpeZcTaSkIgiAIVcQpCIIgCFWk+0jIjXQTJEP0JARhm12IU0iBf6aFUDxRBShv4VI9k8YWe7HhpVNLs5RqHXEKQmmx5aUrCLWEOIU6R2pwQhBiF8mwoRWmGnEKJSOqdqzTKHW/JGqxcCWlHl/A9ZbfZWrVilPQgAqDN1Voksar08iLfEmWqbAK6VGRv0ntUYXd2uAs69IpRGWeDZniJ4th21L7DJM9jTx5C7a8+NUi+jSPzvdUXToFIRkqnJHKsG3FhopEmfVpg/5qlSy6FaegmaxNTxVx5b2vzJhKY5Hx2tiCLCr9qrtqinjONEnlrhunUNaMVEW9p7+s6M4322vpNoytmQxDZ3hh1I1TiCJM2ToyQV7O5UfyUB/1VnuPImvLP6/zJGbOFYBpph6zkI+5+GuBv4UNIAtCWsSWzCN5oJZdX7pgCzMv91+v6ZaCGI6gCrEloV4ofUuBiPYB2GVaDguYDWC/aSEsQXQxEdHHYUQXhzmemY/2Xyy9UxAqENHmoKZgPSK6mIjo4zCii3jkPAVBEAShijgFQRAEoYo4hdphg2kBLEJ0MRHRx2FEFzHImIIgCIJQRVoKgiAIQhVxCoIgCEIVcQoWQ0QbiaifiJ72XDuKiB4gol7n/yznOhHReiLaTkRPEtEyzzMXO/f3EtHFJtKSlxBdtBPRHiJ6wvk7z/Pb1Y4uniWiczzXz3WubSeitUWnQwVE9CYi+iURPUNE24jo0871urONCF3UpW0ogZnlz9I/AH8BYBmApz3X/g3AWufzWgBfcj6fB+BeAATgnQA2OdePArDT+T/L+TzLdNoU6aIdwGcC7n0bgG4AUwG0AtgBoMH52wHgBABTnHveZjptGXRxDIBlzufXAXjOSXPd2UaELurSNlT8SUvBYpj5IQCv+C6/D8AtzudbALzfc/1WrvAogJlEdAyAcwA8wMyvMPOrAB4AcK5+6dUSoosw3gfgdmYeZeY+ANsBvMP5287MO5n5AIDbnXtLBTO/yMyPO5//CKAHwLGoQ9uI0EUYNW0bKhCnUD7ewMwvOp//AOANzudjAbzgue/3zrWw67XCp5wukY1udwnqSBdENB/AUgCbUOe24dMFUOe2kRVxCiWGK+3hep5T/C0ACwC8HcCLAG4wK06xENF0AD8BcAUzv+b9rd5sI0AXdW0beRCnUD5ecpr+cP73O9f3AHiT577jnGth10sPM7/EzOPMfAjATah0AQB1oAsiakLlJfhfzPxT53Jd2kaQLurZNvIiTqF83AnAnSVyMYA7PNcvcmaavBPAoNOVcB+As4loltOEPtu5VnrcF6DDBwC4M5PuBHAhEU0lolYACwH8FsBjABYSUSsRTQFwoXNvqSAiAtABoIeZv+L5qe5sI0wX9WobSjA90i1/4X8AfoBK03cMlT7OVQBeD6ATQC+ABwEc5dxLAL6BygyKpwAs94SzEpUBte0APm46XQp18T0nrU+iUoCP8dz/z44ungXwV57r56EyQ2UHgH82na6MumhDpWvoSQBPOH/n1aNtROiiLm1DxZ9scyEIgiBUke4jQRAEoYo4BUEQBKGKOAVBEAShijgFQRAEoYo4BUEQBKGKOAWhbiGimUT0Cc/3uUT0Y01xvZ+I/kVBOP9ORGeqkEkQgpApqULd4uyVcxczn1RAXP8PwHuZeX/OcI4HcBMzn61GMkGYiLQUhHpmHYAFzn77Xyai+e55DUT0d0T0P865BM8T0aeI6J+IaCsRPUpERzn3LSCinxPRFiJ6mIje6o+EiN4MYNR1CET0XSL6lhPOTiI6w9m0rYeIvuvc0+Dc9zQRPUVE/wgAzLwLwOuJ6I3FqEioNxpNCyAIBlkL4CRmfjtQbTl4OQmVXTebUVnxexUzLyWirwK4CMDXUDkI/lJm7iWiUwF8E4C/e+d0AI/7rs0C8C4A70Vlxe3pAP4ewGNE9HZU9vc/1m3FENFMz7OPO/f/JFuyBSEccQqCEM4vubJH/x+JaBDAz5zrTwE42dmZ8zQAP6pswQOgcniLn2MA7PNd+xkzMxE9BeAlZn4KAIhoG4D5AH4N4AQi+g8AdwO43/NsP4C5eRMnCEGIUxCEcEY9nw95vh9CpewcAWDAbWlEMAygJSRsb7jVsJn5VSJagspBOJcC+Agq+xQBlZbLcIp0CEJiZExBqGf+iMoRjpngyr79fUT0YaB6FvKSgFt7AJyYJmwimg3gCGb+CYDPo3IUqcubcXjXT0FQijgFoW5h5pcBPOIM5n45YzB/A2AVEXUD2IbgIxwfArCUPH1MCTgWwK+I6AkA3wdwNVA9O+BEAJszyisIkciUVEEoACK6EZVxhAdzhvMBVA6q/4IayQRhItJSEIRiuA7A/1IQTiPkaElBI9JSEARBEKpIS0EQBEGoIk5BEARBqCJOQRAEQagiTkEQBEGoIk5BEARBqPL/AW7VxzKcRudqAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x288 with 2 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "from bmtk.analyzer.spike_trains import plot_raster, plot_rates\n",
    "\n",
    "plot_raster(config_file='sim_ch05/simulation_config.json')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Or we can plot the rates of the different populations"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "No handles with labels found to put in legend.\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYIAAAEHCAYAAACjh0HiAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4yLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+j8jraAAAgAElEQVR4nO2deXxcdbn/38/smcm+dUla2tCW0pa9lEVAFlkElIuCgiji9Sfo1auoeF2uV9Er3ntdcecCIihXkUUFBdl3sEBbupdC6Zp0yb5nkpnJ9/fHOWcySSbJZJkkk3ner1denTlz5swzp8n5nGf5Po8YY1AURVGyF9dUG6AoiqJMLSoEiqIoWY4KgaIoSpajQqAoipLlqBAoiqJkOZ6pNmC0lJaWmgULFky1GYqiKBnF2rVr640xZcleyzghWLBgAWvWrJlqMxRFUTIKEdkz1GsaGlIURclyVAgURVGyHBUCRVGULCfjcgSKoijZTiQSobq6mnA4POi1QCBAZWUlXq835eOpECiKomQY1dXV5OXlsWDBAkQkvt0YQ0NDA9XV1SxcuDDl42loSFEUJcMIh8OUlJT0EwEAEaGkpCSppzAcKgSKoigZyEARGGn7cGSNEOyu7+Cul3fT0hmZalMURVGmFVkjBFsPtPLNh7ZwoLVrqk1RFEWZVmSNEOT43AB09sSm2BJFUZTxM9RQsbEMG8seIfBaQtClQqAoSoYTCARoaGgYdNF3qoYCgcCojpe28lERCQDPA377c+43xnxzwD5+4LfACUAD8EFjzO502BNUj0BRlBlCZWUl1dXV1NXVDXrNWUcwGtK5jqAbONsY0y4iXuBFEfm7MWZ1wj4fB5qMMYtE5Argf4APpsOYPiGIpuPwiqIok4bX6x3VOoGRSFtoyFi020+99s/A4NUlwF324/uBc2QstU8pkOOzNC8cUY9AURQlkbTmCETELSLrgVrgCWPMKwN2qQD2ARhjokALUJLkONeKyBoRWZPMFUqFoFdDQ4qiKMlIqxAYY2LGmGOBSmCViKwY43FuNcasNMasLCtLOldhRLRqSFEUJTmTUjVkjGkGngEuGPBSDTAPQEQ8QAFW0njC8XtcuESrhhRFUQaSNiEQkTIRKbQf5wDnAm8M2O0h4KP248uAp81YimBTs4egz6MegaIoygDSWTU0B7hLRNxYgnOvMeZvIvJtYI0x5iHg18DvRGQH0AhckUZ7yPG56Ypo1ZCiKEoiaRMCY8xG4Lgk27+R8DgMXJ4uGwaS43WrR6AoijKArFlZDNZaAs0RKIqi9CerhMAKDakQKIqiJJJVQhD0aWhIURRlIFklBDlerRpSFEUZSFYJgZUj0KohRVGURLJOCNQjUBRF6U9WCUHI76G9Wz0CRVGURLJKCEpyfXT2xLSEVFEUJYGsEoLSXD8A9e3dU2yJoijK9CGrhKDMFoI6FQJFUZQ4WSUEcY+gTYVAURTFIbuEIM8HQH17zxRboiiKMn3IKiEoCWmOQFEUZSBZJQQ+j4uCHK8KgaIoSgJZJQQApbk+FQJFUZQEsk4I5hbmsK+xa6rNUBRFmTZknRAsLA2xu76DNE3EVBRFyTiyUgjauqM0dGjlkKIoCmShECwoDQGwu75jii1RFEWZHmSdECwssYRg1wQKweaaFmK9GmpSFCUzyTohmFuYA8CBlvCEHG/tniYu/tmL/PrFnRNyPEVRlMkm64TA53ER8rlp6YpMyPG27m8B4O1aDTUpipKZZJ0QABQGfTR3jiwE+5u7+MIf19M5zFQzJ+lcFPJNmH2KoiiTSdqEQETmicgzIrJVRLaIyOeS7HOmiLSIyHr75xvpsieRghwvLV0jVw397Om3+NPrNfxt44Eh96lustYkRGO9E2afoijKZOJJ47GjwBeNMetEJA9YKyJPGGO2DtjvBWPMxWm0YxCWEIzsEfjclk4eaB46n+AknScq1KQoijLZpM0jMMYcMMassx+3AduAinR93mgoDHpTCg21hq2Q0NYDLUPus6ehE1AhUBQlc5mUHIGILACOA15J8vIpIrJBRP4uIsuHeP+1IrJGRNbU1dWN257CoJfmFC7cNc1W2GdzTWvS140x8RCTCoGiKJlK2oVARHKBB4DrjTEDr6jrgMOMMccAPwP+kuwYxphbjTErjTEry8rKxm1Tfo6Xls7IiG0mauz4/8HWML1J1gmEI71EYtZ2FQJFUTKVtAqBiHixROD/jDF/Gvi6MabVGNNuP34E8IpIaTptAijM8dET6yUc6Uvw1raF2dfYGX8ejfVysDVMrt9DrNckvdC3hfu2taoQKIqSoaSzakiAXwPbjDE/GmKf2fZ+iMgq256GdNnkUBj0AtCcUDm06qanOP17z8Sf728OE+s1HDOvAICGjsGtq50cQlmeXz0CRVEylnR6BO8APgKcnVAeeqGIfFJEPmnvcxmwWUQ2AD8FrjCT0Ba0IMcWgmESxptqrATx6YutUFSy8ZaOR1BZlENHT4yIlpAqipKBpK181BjzIiAj7PNz4OfpsmEoCm0hGO4uflNNC163cOrhJQA0JulW6ngElUVBXt/bzDcf2sJN/7QC28lRFEXJCLJyZXFBsL9HkCwRvLmmhSWz8phTYPUmakgy1czxCC46ajaHl4X4/St7+dL9G2nSFteKomQQWSkEhUGrHYRT+tnU2XfhdkRhT2MHi8pzKQp6ERkqNGR5BMfMK+TR68/gyDn53L+2mr9vPpjur6AoijJhZKUQFAwIDdW29d3th6MxALp6egn6PHjcLoqCvqTJYscjyAt48bpd3POJk633RmJptV9RFGUiyUohCPnceFwSDw3VJQhBZ491Ee+OxAh4rdNTHPLRkMQjaO2K4hLreGB1NgXoiWrSWFGUzCErhUBE+q0uTvQIumwhCEdjBLzWBT7kc8cFIpG2cIS8gDeeHFYhUBQlE8lKIYC+1cUA9QmJ4I6eKLFeQyRmCHj67vSTXdzbwlHyAn2FV26X4HYJPTENDSmKkjlkrRAUJnQgTUwWd/bECNsxfic05PO46EmyRqC6qYuSXH+/bT63K952QlEUJRPIXiEI+uIrixPLPbv6CYHtEbgHewQtnRHW7m3itEUl/bYP5T0oiqJMV7JWCApy+lpRNyWsMO7siRG2L+T9PIKEi3tzZw/n3fwcsV7D2Utn9Tuuz+OiW4VAUZQMIquFwAkNNXf2UJ5nhXg6e6KDPQKPu19oaGN1C4dauzmmsoBj5xX2O24y70FRFGU6k7VCUBT00RaOEon10tQZoaLIWkGcmCPwe5KHhg61WhPLfnrlcbhd/dtJ+IfIJyiKokxXslYIinOt1cVNHT00d/YwtzBRCJKEhmKJLautKqPyvMCg41phJK0aUhQlc8haISgNWUJQ195Nc2eECkcIuqN0DwgN+QfkCOrauskLeMixF5Il4tXQkKIoGUbWCoFT9rmnoZNor6Es14/XLXRGYvE2E305gsGhISenMJChSk0VRVGmK1ksBJZH8HZtO2ANq8nxuu3y0QGhIffg0FCysFB8X/UIFEXJILJWCEpD1h391gPWGOXikI+gz9O/aihhZXGs18R/3q5rZ1b+MB6BCoGiKBlE1gpBfo4Hj0v4++aDuASOn19E0O+mo59HMLiZ3BfuXU9zZ4TZ9pyCgeg6AkVRMo2sFQIRochOGK9aWExRyEfQ5x6wsrgvNASWEOyq7wDg46ctTHpczREoipJppDSqUkSKgLlAF7DbGDMjrnRO++n3H18JQNBrh4aSJIsBumMx2sNRLjp6DmVDJIv9miNQFCXDGFIIRKQA+DRwJeAD6oAAMEtEVgO/NMY8MylWpolrTl3A23XtXHaCLQR+N40dPfHQkN/Tt44ALI+gvTtKnn9o/dQcgaIomcZwHsH9wG+B040xzYkviMgJwEdEpMoY8+t0GphObnzv8n7Pgz431U0xuiMx/B5XfM6Af4AQ5I4kBBoaUhQlgxjyimaMOXeY19YCa9Ni0RSS4/XEcwTOxR+sRWJgjaDs7ImRGxhGCDQ0pChKhjFislhEnhKRCwdsuzV9Jk0dIb/bLh/tjecHoC9Z3NRhNakb0SNQIVAUJYNIpWpoIfBlEflmwraVI71JROaJyDMislVEtojI55LsIyLyUxHZISIbReT4Udg+4eT4rPLRgeEfJ0fgDLDPG8Yj8LpdRHsNvb06nEZRlMwgFSFoBs7BShL/1U4ip0IU+KIxZhlwMvBpEVk2YJ93A4vtn2uBX6V47LQQ9HroifZS395NYdAb3+4IQaM9wCbX7036/sR9NU+gKEqmkIoQiDEmaoz5F+AB4EWgfKQ3GWMOGGPW2Y/bgG1AxYDdLgF+ayxWA4UiMmdU32ACCdpN5Pa3dFEU9MW3Oxd3Z5LZcDkCvwqBoigZRipCcIvzwBhzJ3AN8PhoPkREFgDHAa8MeKkC2JfwvJrBYoGIXCsia0RkTV1d3Wg+elQE/bYQNIcpTBQCO0fQ2Ol4BIO7jsb3TagwGor27iif+O2a+OI0RVGUqWRIIRCRYhEpBu5zHtvPdwE3pPoBIpKL5Ulcb4xpHYuRxphbjTErjTEry8rKxnKIlHA8glivoSghNOQfTWjIPbIQPLn1EE9sPcQPH98+bpsVRVHGy3DrCNYCBhBgDrDffoy9vWqkg4uIF0sE/s8Y86cku9QA8xKeV9rbpoQcb9/pcNpPQEKyuH3k0NBQHsHBljCFQS8Br5t9jZ0AeAZMN1MURZkKhvQIjDELjTFVxpiFwDbnsbN9pAOLtRrr1/Z7fzTEbg8BV9vVQycDLcaYA2P5IhNBKCHkkyxH0OcRjCwEiY3njDGc/F9P8cm7raUXm2paAGsojqIoylSTUq8hLA9gtLwD+AiwSUTW29u+BswHMMbcAjwCXAjsADqBj43hcyaMoC9RCBKqhuxwz177Tj6UZDKZg9O6OtEj6Oixehc9u93Kb2zZb0XIdtVpjkBRlKknVSEYNcaYF+kLJQ21j8HqZzQtSAwNFSbxCLqjvZy7bBYe99A5dr/dsTScMLe4pSsSf9wT7WV/Sxcugf0tYbp6YklHXiqKokwWwzWd+0LC0/IBzxkm3JOxzCvumzFQFOrzCEI+D+ctm0VrOMI33zNwKUR/nBXJ3ZE+j6Cls08I9jd3YQycuKCYV3c3Ut/ezbzi4ER9BUVRlFEznEeQl/D4tgHPZyR5AS/Xv2sxNz/5FnMSBs+4XMKtV4+4mBroqzByZhpAf49gzZ4mABbNyuXV3Y109sRQFEWZSoZrOvetyTRkunD9u5bwL2cuioeDRovfzhEkJosTheDlHfUALC7PBaCzJzpWUxVFUSaE4dYRfN0eSDPU62eLyMXpMWtqGasIQN9Us0SPoDXcJwRvHGwD4PAySwi61CNQFGWKGS40tAn4m4iEgXX0DaZZDBwLPAl8N+0WZhjJPILWBI/g7bp2yvL88V5GQ4WGVu9sAODkqpJ0maooigIMHxp6EHhQRBZjlYLOAVqBu4FrjTFdk2NiZpHMI0gMDXVHeynP85NjJ5U7I4OFwBjDFbeuBmD3f1+UTnMVRVFGLh81xrwFvDUJtswIhsoRuAScztQhvydeMhpO4hHsqG2PP+7siRL0pa3KV1EUJaWmc8ooGKpqaHZ+IP485HPHL+4Dk8XbDrRy7o+f7/dcURQlnagQTDAul+BzuwZ5BKV5ftx2b6GgzxNfxTwwNPTXDfsB+OQ7Dwfgn+9cQ12btqJQFCV9qBCkAb/X1c8j2NvQyez8QDwvEPS58XtciPSvGtpc08Ivn32bk6uK+fIFR3D64lJauiI8te3QpH8HRVGyh1RmFi+x5xZvtp8fLSJfT79pmYvf4457BG3hCDvrOziqoiCeFwj5PYgIQa+7X9XQZ+95HYCLj56LiHDXx1aR5/fEm9QpiqKkg1Q8gtuArwIRAGPMRuCKdBqV6QS8Lrptj8BpMLeisiDuETiCkOPzxIWgNRxhZ10HHz9tIR8++TDACjMtr8hnswqBoihpJBUhCBpjXh2wTZfDDoPfY+UI9jV2xstAj6ooiM8fcLqXBn1uuuxk8VZbME5bXNrvWEdVFLDtYBsRHX2pKEqaSEUI6kXkcOxW1CJyGTBlMwMygYDXTTgS49k3rbbTV66aT2muP97L26kYCvr6QkPOXf9RFQX9jrWiooCeaC/bD7Zxx4u7+q1JUBRFmQhSKVD/NHArsFREarBGVV6VVqsyHMcj2FzdQmHQy3cvXQFAr7GkwKkYCnjddNkhpF31HRQFvZTm+vsdyxGGGx/awpo9TbR0Rfj8uUsm66soipIFpCIExhjzLhEJAS5jTJuILEy3YZmM4xFs3t/CURUFWMPaEoTAP9gjaO6K9BuP6bCgJESu3xPvWpo3zJhMRVGUsZBKaOgBAGNMhzGmzd52f/pMynz8Hhft3VHePNTG8rl9oR7jrCxOyBE4QtDSGaEwxzvoWC6XsHxufvx5JDaWYXGKoihDM9xgmqXAcqBARN6X8FI+VvM5ZQgCXjc1zV1EYqbfsBtHCJwcQY7PE19v0NzVQ3le8tP67xcdyVPbavnJU2/R0a15ekVRJpbh4gxHABcDhcB7Era3AZ9Ip1GZjt/joi1sXbBLQn0xfzMgR2CtI7D2a+6MsKQ8+eyfoysLObqykN+8tIt2FQJFUSaYVLqPnmKM+cck2pTxOOMqAUpz++L+fU3n+tYTJIaGCoKDQ0OJ5Po9KgSKokw4qWQeXxeRT2OFieKxC2PMP6fNqgwnUQhKEqqAnGRxTkL5aFdPjEisl7buKIU5g5PFiYT8Hg0NKYoy4aSSLP4dMBs4H3gOqMQKDylDUFnUlxcoSfAI8u1ksDfefM5NtNfQ0N4DEB9WMxQh9QgURUkDqXgEi4wxl4vIJcaYu0Tk98AL6TYsk1lYGoo/zvP3neLfXHMij2w6QFme5SU4nsH+FmvGz0hCkBdQIVAUZeJJxSNwlrI2i8gKoAAoH+lNInKHiNQ6zeqSvH6miLSIyHr75xupmz29qbLnEQPxNQQA84qDXPfOw+PbnN5DB5rDABQkKR9NJOQbHBqqbQ3T2NEzIXYripKdpOIR3GoPsf868BCQC/xHCu+7E/g58Nth9nnBGHNxCsfKKBJDQ8PhVA8diHsEqeQI+s8vWPXdpxCBXf+lIy0VRRkbqYyqvN1++DxQBSAi81N43/MismA8xmUqXndqYx6cLqT7bY8g2YKyRHL9btrCg3sNGV1jpijKOBhWCETkFKACeN4YUysiRwNfAU4H5k3A558iIhuA/cANxpgtQ9hxLXAtwPz5I2rQtOBP/3JqPPQzFIM9ghGEIOChoyeGMaZfyElRFGU8DHnrKiLfB+4A3g88LCLfAR4HXgEWT8BnrwMOM8YcA/wM+MtQOxpjbjXGrDTGrCwrK5uAj04/x88v4sg5+cPu0ycEYUQgLzBy1VCs18SH3oQjgwffK4qijJbhPIKLgOOMMWE7R7APWGGM2T0RH2yMaU14/IiI/FJESo0x9RNx/Ewgx2ud/gMtXeQHvPGZxkORa1cgtXdHCXjdNHf2hYnUS1AUZawMF8wOG2PCAMaYJuCtiRIBABGZLfaVS0RW2bY0TNTxMwHHIzjU2j1iWAisqiGAdrt9RVNnX7WQlpUqijJWhvMIqkTkoYTnCxOfG2PeO9yBReQPwJlAqYhUA98EvPZ7bwEuAz4lIlGgC7jCmOxKezpCACMnisHKEUDfRT9RCJo7IyOGlhRFUZIxnBBcMuD5D0dzYGPMlSO8/nOs8tKsJSdBCApGKB2FvtBQR3eU5s4ebnp4W/y1lq7IhGTvFUXJPoZrOvfcZBqSjSRWFaXiEYQScgR/eHUfW/bH0yz98gWKoiijIbWCdyUteNwufPaag5FWFYO1jgAsIXj6jUP9XtNZxkq20NtriMZ6p9qMGYUKwRQT8Fr/Bakki3P91j7VTV2stUdXOjR3aZsJJTu4+Gcvsvjrf2f1zqyqLUkrOgB3inFmFKTiEThzDB7eeIBeA7+86nj8Hhcfv2uNhoaUGU1rOMIX791AT7SXrQeskOjaPU2cXFUyxZbNDEYUAhH5KzCwmqcFWAP8r1NiqoyNT5xexau7GzhtcemI+zrlo1sPtFKW5+eC5bNxuYSA16WhIWVGs3Z3E09s7R8O3VnXwbq9TRQHfSxI6PirjJ5UPIKdQBnwB/v5B7HmESwBbgM+kh7TsoPPvWsxqS7UdrkkPvD+HYeX4LIXoBXm+Gju1NCQMnPZWd/R7/mKinx21rfzvl++rE0XJ4BUhOBUY8yJCc//KiKvGWNOFJGkvYGU9OH3uOjsiTG/pO8OqDDo1dCQMqPZVd9OQY6Xzp4okZjh6MpCfv/KXkCbLk4EqQhBrojMN8bshXjnUafhvt6GTjLOYrLKwr5W1wU5Xpo1NKTMYHbVd7CwNMRtV68kHInx2JaDU23SjCKVqqEvAi+KyDMi8izWdLIbRCQE3JVO45TBRGLW7c/cBCEoDHppUY9AmcHsquugqixEWZ6fecVBzjlyVr/Xs6wpwYQzohAYYx7BCmJfD3wOOMIY87AxpsMYc3O6DVSSU5Ew/KYwx6flo8qMJRrr5UBrmMqiYHzbwtIQ/37hkfHnTrHElv0tdPZo363Rkuo6ghOA5cAxwAdE5Or0maSkwpyCQPyx5giUmUxTZwRjoCy3fxuWT5xRxU+vPA6A2rZuNuxr5qKfvsh/PfLGVJiZ0aRSPvo74HBgPeA0wDcMP4JSSRMiVnIs4E3sU+SlO9pLOBLrt11RZgINHd0AlOT6B71WnmdtO9Qa5tbndwLw9Bu1/OfkmTcjSCVZvBJYlm2dQacrz91wFnXt/ZduFOZYd0otXREVAmXG0dBuhT1LQoMbM87KtzzjA81h1uy2VtvXNHdR3dTZL5SkDE8qoaHNwOx0G6KkxvySICccVtxvm9OeQsNDykykvt3xCAYLQUVhDiGfm/vW7qMrEuOaUxcAsLmmddC+ytCkIgSlwFYReUxEHnJ+0m2Ykjpltnu8v7lrii1RlImnzyMYHBryeVycvriM12xv4L3HzgWsclMldVIJDd2YbiOU8bF0dh4Am2taOGtp+RRboygTS2NHD26XDNmP69xls3h0y0H8HhfHVBZSludnV337JFuZ2YwoBDqXYPqTF/BSVRpiU03LVJuiKBNOQ0c3xSFfvKXKQC45di7FuT7mFARwu4SFpSF21qlHMBqGDA2JyIv2v20i0prw0yYiGoCbZqyoKODxrYd4ZNOBqTZFUSaU+vaepIliB4/bxVlHlLN0dj4Ah5eFNDQ0SoYUAmPMafa/ecaY/ISfPGNM/uSZqKTCxUfPAeDGh7boKktlRrGnoaPfSvqRmFOQQ0NHjw6vGQXDJotFxC0iujojAzhv+Wx+ePkx1LZ1DxpaoyiZSmdPlB217ayYm/q9pzMCNhxVIUiVYYXAGBMDttuN5pRpzplHlAFw2S3/YLPmC5QZwLYDrfQaK/SZKs7Uv66e2Ah7Kg6plI8WAVtE5CktH53elOT6ueG8JQDsa+ycYmsUZfxsqrZuaI6qHI0Q2B5BRIUgVVIpH/2PtFuhTBgXHz2XHzz+Jp16N6TMAF7c0UBFYQ6z8wMj72yjQjB60lY+KiJ3ABcDtcaYFUleF+AnwIVAJ3CNMWbdWD5L6SPos/4IOvWPQMlwwpEYL+2o5/KVlViXi9SI5wgimiNIlXSWj94JXDDM6+/Gam+9GLgW+FXqZitDkWMLQVg9AiXD+dvGA3RFYoNmD4yE4xF06c1QygznEVwFVvnoWA5sjHleRBYMs8slwG/tZnarRaRQROYYY7QQfhwE7QH3GhpSMpneXsP3Hn2DY+YVcvqi0lG9N8dn3d9qaCh1hksW/9l5ICIPpOGzK4B9Cc+r7W2DEJFrRWSNiKypq6tLgykzB7dL8HlcdEZ0OIeSuXT0RKlt6+bio+YMuaJ4KPwe9QhGy3BCkHj2q9JtyHAYY241xqw0xqwsKyubSlMygqDPraVzSkbjXMQDvtG3VY+HR1UIUmY4ITBDPJ4oaoB5Cc8r7W3KOAl63RoaUjKabjvRmzOG+RpaNTR6hhOCY5zkMHB0GnoNPQRcLRYnAy2aH5gYctQjUDKcuEfgTXWabh9aNTR6hkwWG2PGNepKRP4AnAmUikg18E3Aax/7FuARrNLRHVjlox8bz+cpfQR9Hh3grWQ0zo3M2DwCe2WxegQpk8qCsjFhjLlyhNcN8Ol0fX42k+Nz6x+BktGE4x7BGITASRarV5wyo/e7lGlPjldDQ9lKJNZLbVt45B2nOV3jEAKXXTkXjurfQKqoEMxAgj5NFmcrv3hmB6tueoq9DZnda8qJ748lRwDWzZAuqkwdFYIZSI4KQdayo9Ya0fiDx7en/J7711ZzxveemVZzLJzQ0FhyBGAJiCaLUydtOQJl6ghqjiBryQtYf9Lr9qY+k+KG+zYA0NETI9c/PS4J48kRgB0e1b+BlFGPYAaiVUPZS3u3dfGrae5KuY7e67bWjjZ19KTNrtHSNW6PwK3rCEaBCsEMJMfrJhzppa6tm9ZwZKrNUSaRjm7rBsAY2JviTAqPy7oMNHdOn9+V8SSLnfepR5A6KgQzkPwcLwAn3vQkR9/4OC/vqJ9ii5TJor07SshusbCzLrUB7h7HI+icPh6BE9/3e8Z2icr1e9jf3EVv7/TJe0xnVAhmIJcdX8lNl67gaxcuBeCNg21TbJEyWXT2RFluj3X85N1rqW4a2Svwuq3LwPQSghgBr2vUDecc3nd8BW/XdfDn17VrTSqoEMxACoJerjrpMD5xehVet1Db1j3VJimTREd3jFn5AS49zmrk+/KOhhHf43FNvxyBJQRjb25w6XEVVJWG+Mt6FYJUUCGYwYgI5XkBalszf4GRkhrt3VFy/W6+f9nReN3CzvqRw0NxIZhOOYKe2JgTxWD97r9r2SxW72ygvVsLJ0ZChWCGU5bnV48gi+jsjhL0efC4XcwvDrKrvn3E9zhJ1ebpFBqK9o7LIwA4c0kZkZhhze7GCbJq5qJCMMOZle/nUBZ4BHsbOjnxpifZUdvGwZYwJ333Sc7+wbPjbrVxzW9e5QePpb44ayrp7TV09MQI2WsBqspy2TWER3Dva/v4wC3/IByJxe+Y72Sn110AACAASURBVPrHHh6cJqGUrp7xhYYAKouCANTpjdCIqBDMcMrzAlnhETyzvZa6tm5+8tQOth1s5VBrNzvrO/jHzrFXTBljeHZ7HT9/ZscEWpo+Ou07+1y/dQGtKg2xq76DT/9+HQ+sreb6e17nzpd2Adb5enV3I798ZgeRWF9lzb1r9g0+8BTgJIvHQ3GuDxh7EjwS6+Xrf9nENx7cTDQ2s1cpT49lhEraKM/z09IV4dVdjaxaWDzV5kw4De3d7Gvqile+/HXDfuYUBOKvP7mtlrOXjm74uUNrV2bFlp01BM7c6nOOnMVzb9bxxJZDPLzRGvXxl/X7WTI7L15a+tOnLZGbV5zDvsYuKgpzpsDywYQj48sRAIR8bnxuF40dg3MfL++opzUc4fzlsxFJXpm0s66Du1fvBeCqkw7jiNljGt+eEahHMMNZYv/yfuw3r86Yu5rGjh5qW8P09hqu+91arv71K/3u+m59ficApy8u5dVdY48PVzf3lV62ZcDCPGdBmNMmYtXCYh69/gyOnV8IWPOs5xQE+NWzb7OroYOLjp4Tf+/15yyhqixER/fgUFok1jvpfYhauiLxENdYERGKQz4aO/o8YmMMD6yt5kO3v8In717HhuqWId/f3t33f17TnNlN/EZChWCGc/7y2Xz7kuV09MRYuyf1/jPTlce3HOT4/3yCVd99igt+8jxr9jTRGo5S3dRJjtfNecusu//CoJcls/Koaeoa80Vsf3NfbmV3/fS+EIQjMc6/+XmAQRfQSvsuv6o0xEVHzeGFt+rpifZy2qJSikNW+CQ34CHX76FjQGuScCTG8m88xrf+unUSvoVFb69hT2MnC0qC4z5WUcjXzyP4ygOb+OJ9G3DblVKbqpuHfG9buO9c1DTP7DybCkEW8L7jK/G6hUe3HJxqU8bN1gPWlNTLTqjkzUPt+OyQ0PaDbRQFvZxyeAkAbhEqCnPoisTGXBa5v7kr/vjy/32ZlmlUXjkQ57zA4NW4c20hmFuYw9lHlse3LywNcfriUgCiMUPQ5+aFt+q5+GcvcMHNz3PBzc9z3o+fpyfWy50v707/l7Cpae6iJ9rLwtLccR+rOOSNe4tr9zTxxzX7OOuIMl768tkUh3xsqhnOI0gQgqauIfebCagQZAG5fg8XHTWH37y0mx8+vp1YBi+7r2nqojzPz02XruCT7zyc/3jPMsASgsKgj4WlIQA6eqLxC2DiBd3hH283sH6fdTe4dk9T0hLDmuYuvG7hnKXlhCO9PLHtULq+1rjZbF/Q3nXkLE44rKjfaxVFOfF/Vy0o5qqT5nPZCZUcO6+Qb713OdecuoCzlpaR6/cQ6zVsrmllXnGQ+cVBltqhxcqiycsdOJVOVWWhcR+rOOSn0V4o97rdkfVHHziW2QUBVlQUcO+a6vjvwUCcnIvXLUl/h2YSKgRZwlfefSQAP3t6B28cbB1h7+nL/pYu5hbm4Pe4+cq7l3JKlZUA7+iJURTyUmXfRYYjvfGL197GTvY1dlLf3h2/y7vyttX80y9eAuD9v3qZy275x6DPqm0NM7sgwO0fXcmsfD8Prq9h9c4GNte0TKve/Z09Ue55dR/FIR+3XX3CoNCQI4gVhTl43C5uuvQofnD5MQS8bgqDPm5873KCPk/8fXl+D7ddvZJb7Z9/OfNwDraEJ+0GIi4EpRMgBEFvXAgOtYbxe1wUBq1eXE4Y8arbVvPCW3WD3uuEhpbMyqNmAoRgZ137tP3bUyHIEmYXBPjjtScD0JSkiiJT2N8cjt/hApTn91UIFQZ98dcuWD47fgG86eFtnP69Z1j5nSf58O2v9Dve23VDL7hqC0cpyPEiIpy3bDYvvFXPFbeu5uKfvcjmmunzB/3ff3+DrQdaOX5+YdIKmMPLQrgEFpUPH2pxhKAw5O23vaIoh2ivmbQRmLsbOgj63JTl+cd9rKKQj5auCNFYL7Vt3ZTn++Pn6MMnH8Zj159BR0+Mj/z6Vf7xdv92HM5NwxGz89gzzolvnT1RLrj5BS64+QX+Mg37H2n5aBZRZCcGG6fRCtLR0NtrqGnu4txlfeWgeX5PfBpVUdCL2yWs/uo5FAa9+D0ucrzufndz6/c19+tT/6PH34w/jvWaeBIRoDUcIc9vXRS/euFSLjxqDtVNnXzp/o3sbezkqMqC+L4fveNVjqoowOt28avnduBxufjfj5zAOxaVpuVcOBhjeHzLIRaV5/KDy49Juk9lUZBnbjiT+cXDJ1+drqVFQV+/7U5JaU1TFzleN1ff8Srzi4M880YtvQZuunQF7zu+clzfI9Zr+Oc7X+OkqmJqW7uZnR8YsqxzNJSEnLUEEWpbuynPC/R7/YjZeTz++TM478fP8/jWg/EcE0B7OErQ52b53AL+tK6G2tZwvxuPVPnlszv41bNv02NX7X3+3vV8+YGNBLxu/njdySydnT+ObzgxqEeQRTh/4NOplcBo+M7D2+iJ9vardReR+J2j8/1mFwQIeN2ICIcNqDzxuqXfStOHNx2IP24c0HStLRyNT/wK+jyccngJ7zrSEqHEu+NIrJfn3rQWnv34yTc5dl4h7d1Rtu5Pv9ewuaaVg61hrjujisIBF/BEDisJjXhhdTyC/MAAj8A+3995eBtX3f4KG6tb+NvGAyyalUevMWwcpgQzVf7w6l6ee7OOl3bUc6g1THn++L0B6Lv5aers4VBbmFlJjrtkVh5nLy3nyW2H+rWttvo2eTjK7uY6XGJ5KPY1dvK9R7fHw0y326G2K1fNp6Urwmu7p0clnwpBFuHERjMxNFTX1s0d9qrYUxPu2gDed1wly+fmc1qSu++FA+LMMdurALjxPctYtbAYxwkY2IrDEoL+F8XCoBef28Wh1m5e3lHPA2ur2ZcwAOasI8q45cMn4BImZSjQj57YTq7fwzlHjm3RXCLO+gPfgKqjw0pCnLGkjHAkRqzXcN0ZVZy/fBY/v/I4SkK+cTd1a+roic9YrmnqskI4eaO/805GsS2OjR091CXxCBwuOXYu+xq7uH9ddXxbW3eU3ICH5XPzERleCDp7oty3Zh/bB7R8HxhueseiUr50/lK+cfEycrxudqU4MyLdpDU0JCIXAD8B3MDtxpj/HvD6NcD3ASdo9nNjzO3ptCmb8bpd5Pk906rvfKrct9ZqffDwZ09j8az+Kzw/f+4SPn/ukqTvc4TgoqPncMbiUr78wCY22rXjJ1WVcM07FrJubxPv++XL1LV1U9sWpiDHi9/jprUrEvcIHBwPZEdtG7c89zYAH1hphUXu/vhJnGaXY+b6Pf3q0CcSp0Xz3oZOntlex5fOPyK+HmA8OB6BU5Lr4PO4+O0/r0r6ntyAh/Zxfs9HtxykuTPCaYvsBYBC0jv3seC0mahu6qKtOzpk3uE9R8/lthd28n+r9/CBlfMAKzSU57eS6IvLc1kzxN17R3eUf/3D6zz9Ri15fg/3feoU8gNeikO+QUnmHDv85nIJC0pDSZsCrt3TyBfv3cADnzqVktyJOQ8jkTYhEBE38AvgXKAaeE1EHjLGDFyZ8kdjzGfSZYfSn8KQN+NCQ/et2cf3Ht3OrHw/y+aMLp7qhIaMMfHk8et7LSEoty8Kzr//+fBWdtZ1cHJVMb//fyfT3hONT3tLpDzfz5PbauPP711j3UXOLey728wLeNPiEby+t4lLf/kyt1+9knDUynW8c0nZhBzbWX/gHcVUsFy/Z9wewYGWMCJw9tJyXrSn6U20R+CE6WYNEeN3uYQT5hfxwLoajDGICB3d0bg4nnlEOb95aVc8XOTQE+3lnd9/hvr2Hs5ZWs7zb9Vxwc0vALB0dh4rKgqYnR/gYJLGj1WlIbbsH+xl/OSpHexu6GRnfUfmCwGwCthhjNkJICL3AJcAk7dEURlEUdA3rfrOp8IG+w7+jmtOHHUC0ckbBLzueKz7qTdq8bgk/ppzl+j039m6v5X2nijGQH5g8J+IIxwlIR/nHFmeIAR9uYv8HO+EegSv7W7kGw9uocaeOHbPa3tZVJ6Hz+1iyayJ6YETtePjzjD7VAhNgOdT1xamJOTvl8+ZqByBkzd5+W1LYJYO0y9oYWmI9u4ode1WCKm9O8r8kGXT2UvLufX5nbz4Vj0XrJgdf89ruxupb+/hylXz+M9LVvDGwTa27G/h5bcbeHD9fkSEuYUB/u8TJw3ytKrKQjy65SA90d54OG77wTaef9MqZZ3MQUHpFIIKILGVYTVwUpL93i8iZwBvAp83xgxqfygi1wLXAsyfPz8NpmYPhUFfxnkEu+o7OGZeIcvnFoy88wDOWlrOdWdU8f9Or6Ik5OOaUxdwqDXMioqC+BhEv8fN585ZzIGWLkJ+D795aTcHW6w7uIGhIegrK7x85TzOXTYrLgSJbZPzAh5au8YvuKt3NhCJ9bJ6ZwNvHmrj/OWz2FTTwiu7GtlR284Rs/MGxfTHSk909HOC8wIeDrSMr6z0UGs35Xn+fkI6UR6Bz+MiL+DhjYNteN0yrGhWlVnltV++fyOLynOpae5i2VzLAz12ntWv6a1Dbf2E4Klttfg8Lv7j4mV43C5WVBSwoqKAyqIgD67fz7YDrVx89BwOLxtcurt4Vh6xXsObh9pYYSekf/3iTkTAmMkdHTrV5aN/Bf5gjOkWkeuAu4CzB+5kjLkVuBVg5cqV02clTwZSFPSyO4WpVdOJXXUdnFxVMvKOSfC6XXz1wiPjz2987/Kk+zk5BqcfvzPneWCyGODKVfOpb+vhM2cvItfv4cMnz6c70r+hX37A069X0Vi54tbVAJxcVcycggC/vOoEXtpRz6fuXkttWzeX2/HsieCspeW4XcJVJx2W8nty/ePPEdTa1TwLSkIcVhIkHImxZNb420s4FId8tIWjI4qmk096Znsdr+xqRIDj5lurtANeN6W5/kEx/7V7GjlhflG846vDioSbloohVmU71Uhb9rewoqKA7QfbeGBdDZcdX8l9a6sn1XNPpxDUAIm/pZX0JYUBMMYkptRvB76XRnsUrFDJwDLJ6UxXT4z9LeFB1T/pwgkfvWH37hlYSglw8dFzufjoufHn3/mnowbtkxfw0tbdNmj7aEhcvbx6Z18b8XcsKmXjjeeP69jJqCjM4e3vXjiq9+T6vePOERxq7Wb5nAJyfG6e+9JZ4zpWMqL2vIVjKguH3c/xSFZU5PO3fz190OsVhYF+QmCMYWd9R3w+dCIFQS+HlQTZ09DJ4UP0TDqsOEie38OmmhY+sNLwjQc3kxfw8LULj+TBDftnTGjoNWCxiCzEEoArgA8l7iAic4wxTiH3e4FtabRHwYqHt3dH6eyJDrqLmY7sbrC8lwWTJQT23VufRzC2c5QfGH/svL69/4WgcprMCkgkN2Ali3t7TTzUNhqisV4a2rsnrEooGd12yOsTp1cNu5/bJbz8lbMpHSJBW1GUE/+9AGjo6KEtHB3yJuXe606hprmLoyuShzRdLmF5RT7Pbq/jhO88SWNHDzdduoKikI/ioG9SQ0NpW0dgjIkCnwEew7rA32uM2SIi3xaR99q7fVZEtojIBuCzwDXpskexcKomalszY2pZtd31caRVsRNFeV4Aj0vicwzGKgR5AStZPJ6eRE7Pnbn2oJ05hRMTN59I8uwKmoHtq1PBGMNNj2yj10DZGFbspsqdHzuROz92Yko3E3MLc4YMH80tyGF/c19bc+f/ZyghmJUf4Pj5RXjcQ19mP7ByHtVNXTR19vD1i47kihOtHGhh0Jt0oE66SOstoTHmEeCRAdu+kfD4q8BX02mD0h+n4qW2rXvS7rLHg9P1cag460TjdglXrJrH63ubKc31j/lz8wJWJ8/OhBnCieysa+cfOxs4b9nspLXtkVgvNz9ptb84qaqEP79eE6/qmU7k2kLZ3j148d1IPLmtlt+8tJs5BYF488B0sGKIO/LRUlGUQzjSy69f3IXf62aD3bW0ahztsi89roL1+5pZtbC4X7ixOORj9c4GdtV3TEpYdPrHBpQJxSnLm6wGYuOlprkLv8cV7xkzGSSL+Y8W56LYFrZq0Y2xVjQbYy0q+soDm3h1dyP3vLqPv3z6Hf16HAH8ffNBXn67gZDPzT+/YyF/fr2Gs44oT/ZRU4pTU9+RQp6gob2bgNcdF8bfv7KHecU5PP3FM+OjRqczzhqW7zzcF8Eez80CWAsUv33JikHbfR4X7d1R/ukXL7Hhm+eN+fipokKQZcyyy/IOZUhoqKbZmqM7EQ3IJpMCeyFaY0eP1cr6hV3c9Ej/FFhprp9NNS388bV9fOik/mXRT249RNDn5rWvv4ugz8PO7144phh8unGEYKR8yBNbD/HJu9dSkOPlkc+ezuyCADvrOzimsjAjRAAsz2zDN8+Ll9mC5fkNFPGJwKkYapmAEuRUyIz/AWXCcHrlTHeP4NHNB7jud2uoburqV1+eKTgtn7cfsqqPXtnVQEVhTr8Oobd/dCUnLSzme4+9QWdCjL2rJ8az22u58Kg58YT+dBQB6AsNNbT38JUHNvLwxgOD9on1Gm58aAtgCePp33uav27Yz77GzgmZOTCZFOR4Kcvzx38S145MJP/z/qPI8boJ+tyTMvtChSDLcHrl1E1jj6ClK8In717HY1sOsWFfc79uo5nC4WUhAl4Xm6otIdhU08KJC4q47IRKHv/8GXzmrEUcXVHAR09dQHNnJL6q+a1DbVx522paw1E+eOLErRFIF0tn51ES8vE/j77BPa/t49O/Xzdon/3NXdQ0d/Gfl6zglg8fTyRmuOlhK0m8cAKmkM1Els7O5wvnLqGzJ0ZrV3p6ViWioaEsZHZBYNjVoJFYL6/tbiQSG/lOpDzPz5Gj7P+TjM01LURivRw3v4iX7H4zDisqpr5f+2jxuF0sm5PP5poWatvCHGrtjictl8zK44bzjwD6RkDub+5i8axcrv3dWnbVd3DtGVWcuCB9CdSJIi/g5csXLOXfHtgIJK+g2WlX1ywqz2XVwmI++c7D4w37JmIu8UzFqRI70NpFQXB0ifjRokKQhcwtzIl34EzGg+v3c8N9G1I6lsclrP2Pc+Mx8bHQ0N7NxT97EYDHrj8jXpY3vzjI/uYuPnhiZrYVOXZeEXe/sodfv2i1zz5u/uAFTfGhL81dPL7lELvqO7jjmpWcvXT8baUni8tOqOSe1/aybm9z0iZ9u+wpcI5InL98VlwIJmIu8UxlToH1u3GgOZz24TUqBFnI3MIAj20OD7kIaFd9O26XcO91JwNDx6bX7Wnipke2cdy3H8cAJSE/T33xnSmLwsMbD/DdR7bxX+/rq9I5/+bnAZidH+Cv/3oaMLg/fqbw8dMX8vtX9/C/z+3klKoSjp9fNGif4pCPgNfF/uYu9jd34fO4OH3xxHQTnSxcLuG3Hz+JT929lhfeqmfJ1/9OxJ7GVZrrZ3F5Lnl+D6V2S+jj5hfx2PVnEPC6kq7cVizm2OtHdtS2c9bS9FaMqRBkIZWFOfTEeqlv7046em9/c5jZ+QFOOGz40ERBjje+IGjZnHy2Hmjlc/e8zkdPWRD/xT3UGuamh7clbcn87Hary+I9r+0d9Nq84pxxeRnTgYrCHO762Cpe2dXIFavmJa18srpT5lDT3EVjRw9HzsnPmCqaRHL9HqpKQ7zwVj090V4+c9YiXAJ3/WMPL7/dwNGVBf2+/xHDdAFVLGblB1hUnst//X0bpy4qGVPTxVTJvN84Zdw4VTjVAxpoOdQ0daVUG5242vdz71oMWBf3j935Wnz7t/66hUe3HKSpo2fQz0l275xHNh0cdOyGDOqHNBwnVZXw2XMWD9tNs6Iwh801rWyuaeWoDMyHODgtnwuDXm44/wi+cN4R/ODyYzh2XmF82IuSOm6X8MdrT0ZEeGTT4GqsiUQ9gizEucjXNHUlDVfUNHfFG5wNR2LI5pTDB3cH3dfYySObDvLZcxbzhSEmiH3pvg3ct9Zq4/yJ0xdy2wtWPP1jpy4Y8fNnCkdVFPDCW1aC/NTD0zvsPp0U2QnNUEIPq3OXzeLcZZmT75hulOT6OXFBEU9tq+VL5y9N2+eoEGQh1gIta1LXe46Z2++1aKyXg63hftO2hsPrFiIxQ37Ay/ziIHsbOwna4/iefsOa4pWsO6PDMfMK40Lw7xct42t2y+hMW0A2Hr50/hF84vQq3G7J6Ji5Myg+5E9PbX228o7DS/nhE2/S1ROLj7qcaDQ0lIXkBbxcelwFv1u9mwMt/cND+5q6iPWalBdxvfBvZ/PsDWcC8MCnTuXyEyrp7Ilx8nef4vuPbaeqNDRsr5SjBvSBEZGsEgGwvnNRyJfRIgB9oaGcDOhqm0k4f4vJxl1OFCoEWcoHV84jEjPxhUwOP33qLXxuF2ekWLkyuyAQb15XlufnnUdY7zvYGuZdR5b3GwqTDE0azjxCabprzVZm29VDB8c5CW44VLqzFKcpWuJQke5ojL9t3M+HVs1n3hjbPieuAr75iuNG3D/gdXPeslmsXDA4V6FkFsdWFjKnIMAXzztiqk2ZUTit4w+l0SNQIchSnGZhiWMG3zzYTiRmWLVwbGMhoa+S6FNnHp7ye269euWYP0+ZPhQEvfzjq+dMtRkzDmdoj4aGlAknsY+8w6aaFmB8LR1Kcv28+rVz+Lfz9a5QUSaCvICXkM/Na7saicZ6R37DGFAhyFKcyg5HCDbXtPC1P28iL+AZ9zSw8vxA1iV8FSWddPTEeOqN2kGtzCcKFYIsxe9x4/O44n3kH1jn1PJX6UVcUaYZZ9lFGNekaX2N5giymDy/h/Zuq/XDlppWjp9fyGfPWTzFVimKMpCffeh4eo1JW4mxegRZTG7AQ3s4Sm+vYcv+lkE1/YqiTA9y/Z60rjNRIchicv0e2sJR3verl+noibFchUBRshINDWUxIb+HfU2dvHmonaKgl/O0J4yiZCXqEWQxeX4Pbx6yhob86APHxlsEKIqSXaRVCETkAhHZLiI7ROQrSV73i8gf7ddfEZEF6bRH6Y+zlgCSjxhUFCU7SJsQiIgb+AXwbmAZcKWILBuw28eBJmPMIuDHwP+kyx5lMM7qYq9b4rNzFUXJPtLpEawCdhhjdhpjeoB7gEsG7HMJcJf9+H7gHNEi9knD8QjmFQfxZOBULEVRJoZ0/vVXAPsSnlfb25LuY4yJAi3AoEY3InKtiKwRkTV1dXVpMjf7OGfpLE5aWMzVJx821aYoijKFZETVkDHmVuBWgJUrV5opNmfGsGphMX+87pSpNkNRlCkmnR5BDZA4qLTS3pZ0HxHxAAVAQxptUhRFUQaQTiF4DVgsIgtFxAdcATw0YJ+HgI/ajy8DnjbG6B2/oijKJJK20JAxJioinwEeA9zAHcaYLSLybWCNMeYh4NfA70RkB9CIJRaKoijKJJLWHIEx5hHgkQHbvpHwOAxcnk4bFEVRlOHRmkFFUZQsR4VAURQly1EhUBRFyXJUCBRFUbIcybRqTRGpA/aM8e2lQP0EmjORTFfb1K7RMV3tgulrm9o1OsZq12HGmLJkL2ScEIwHEVljjFk51XYkY7rapnaNjulqF0xf29Su0ZEOuzQ0pCiKkuWoECiKomQ52SYEt061AcMwXW1Tu0bHdLULpq9tatfomHC7sipHoCiKogwm2zwCRVEUZQAqBIqiKFlO1giBiFwgIttFZIeIfGWKbdktIptEZL2IrLG3FYvIEyLylv1v0STZcoeI1IrI5oRtSW0Ri5/a53CjiBw/yXbdKCI19nlbLyIXJrz2Vduu7SJyfhrtmiciz4jIVhHZIiKfs7dP6Tkbxq4pPWciEhCRV0Vkg23Xt+ztC0XkFfvz/2i3qkdE/PbzHfbrCybZrjtFZFfC+TrW3j5pv/v257lF5HUR+Zv9PL3nyxgz43+w2mC/DVQBPmADsGwK7dkNlA7Y9j3gK/bjrwD/M0m2nAEcD2weyRbgQuDvgAAnA69Msl03Ajck2XeZ/X/qBxba/9fuNNk1BzjefpwHvGl//pSes2HsmtJzZn/vXPuxF3jFPg/3AlfY228BPmU//hfgFvvxFcAf03S+hrLrTuCyJPtP2u++/XlfAH4P/M1+ntbzlS0ewSpghzFmpzGmB7gHuGSKbRrIJcBd9uO7gH+ajA81xjyPNQsiFVsuAX5rLFYDhSIyZxLtGopLgHuMMd3GmF3ADqz/83TYdcAYs85+3AZsw5q9PaXnbBi7hmJSzpn9vdvtp177xwBnA/fb2weeL+c83g+cIyIyiXYNxaT97otIJXARcLv9XEjz+coWIagA9iU8r2b4P5J0Y4DHRWStiFxrb5tljDlgPz4IzJoa04a1ZTqcx8/YrvkdCeGzKbHLdsOPw7qbnDbnbIBdMMXnzA5zrAdqgSewvI9mY0w0yWfH7bJfbwFKJsMuY4xzvm6yz9ePRcQ/0K4kNk80NwP/BvTaz0tI8/nKFiGYbpxmjDkeeDfwaRE5I/FFY/l506KudzrZAvwKOBw4FjgA/HCqDBGRXOAB4HpjTGvia1N5zpLYNeXnzBgTM8YcizW3fBWwdLJtSMZAu0RkBfBVLPtOBIqBL0+mTSJyMVBrjFk7mZ+bLUJQA8xLeF5pb5sSjDE19r+1wJ+x/jgOOa6m/W/tVNk3jC1Teh6NMYfsP95e4Db6QhmTapeIeLEutv9njPmTvXnKz1kyu6bLObNtaQaeAU7BCq04ExITPztul/16AdAwSXZdYIfYjDGmG/gNk3++3gG8V0R2Y4WwzwZ+QprPV7YIwWvAYjvz7sNKqjw0FYaISEhE8pzHwHnAZtuej9q7fRR4cCrssxnKloeAq+0KipOBloRwSNoZEJO9FOu8OXZdYVdQLAQWA6+myQbBmrW9zRjzo4SXpvScDWXXVJ8zESkTkUL7cQ5wLlb+4hngMnu3gefLOY+XAU/bHtZk2PVGgpgLVhw+8Xyl/f/RGPNVY0ylMWYB1nXqaWPMVaT7fE1kpns6/2Bl/d/Eik/++xTaUYVVrbEBeKP+6QAAAoRJREFU2OLYghXXewp4C3gSKJ4ke/6AFTKIYMUePz6ULVgVE7+wz+EmYOUk2/U7+3M32n8AcxL2/3fbru3Au9No12lYYZ+NwHr758KpPmfD2DWl5ww4Gnjd/vzNwDcS/g5exUpS3wf47e0B+/kO+/WqSbbraft8bQbupq+yaNJ+9xNsPJO+qqG0ni9tMaEoipLlZEtoSFEURRkCFQJFUZQsR4VAURQly1EhUBRFyXJUCBRFUbIcFQJFUZQsR4VAUcaBiCyQhFbZo3jfXBG5f4jXnhWRleO3TlFSwzPyLoqiTDTGmP30rRRVlClFPQIl67Hv6reJyG32kJLHRSRHRI4VkdV2J8o/S9+wmRPEGmiyAfh0wnHcIvJ9EXnNfs91I3zmZvtxjojcY9vwZyAn3d9ZURJRIVAUi8XAL4wxy4Fm4P3Ab4EvG2OOxmor8E17398A/2qMOWbAMT6O1YPmRKzulZ+w+/iMxKeATmPMkfZnnDDub6Moo0CFQFEsdhlj1tuP12K1bi40xjxnb7sLOMNuVFZorME5YPXycTgPqzHZeqxZACVYAjMSZ2D1tcEYsxGr/42iTBqaI1AUi+6ExzGgcAzHECxP4bGJMUlRJgf1CBQlOS1Ak4icbj//CPCcsXrXN4vIafb2qxLe8xjwKXsuACKyxG41PhLPAx+y37MCqzOmokwa6hEoytB8FLhFRILATuBj9vaPAXeIiAEeT9j/dmABsM7uZ19HarOnfwX8RkS2YfXqn9TpVIqibagVRVGyHA0NKYqiZDkaGlKUNCIiR9G/sgig2xhz0lTYoyjJ0NCQoihKlqOhIUVRlCxHhUBRFCXLUSFQFEXJclQIFEVRspz/Dx6qPUHzhBOSAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "plot_rates(config_file='sim_ch05/simulation_config.json')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "In our simulation_config.json in the reports section, we can see we also record the V_m (i.e membrane potential) of a select sample of cells. By default these files are written to an hdf5 file with the same name as the report (membrane_potential.h5), and we can use the analyzer to show the time course of some of these cells."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "## 5. Additional Information"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Recording membrane potential  (and other cell properties).\n",
    "\n",
    "By default bmtk will automatically record and save the spike-trains of our (non-virtual) cells. Sometimes it can be useful to record other properties besides spikes, like the membrane voltage traces, of all or a subsection of the cells. To do so we only need to open the _simulation_config.json_ with our prefered text editor and insert the following section:\n",
    "\n",
    "```json\n",
    "  \"reports\": {\n",
    "    \"membrane_potential\": {\n",
    "      \"cells\": {\n",
    "          \"population\": \"V1\", \n",
    "          \"node_id\": [0, 20, 60, 80, 100]\n",
    "      },\n",
    "      \"variable_name\": \"V_m\",\n",
    "      \"module\": \"multimeter_report\",\n",
    "      \"sections\": \"soma\"\n",
    "    }\n",
    "  },\n",
    "```\n",
    "\n",
    "This will tell bmtk to record the \"V_m\" variable for selected cells with ids 0, 20, 60, 80, 100 - and save the results in a file _output/membrane_potential.h5_. After re-running the simulation the file will be generated, [according to the SONATA format](https://github.com/AllenInstitute/sonata/blob/master/docs/SONATA_DEVELOPER_GUIDE.md#frame-oriented-node-element-recordings), which we can read using hdf5 or bmtk like below:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXkAAAEWCAYAAACDoeeyAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADt0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4yLjByYzEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy/xvVyzAAAgAElEQVR4nOydd5jc1Ln/v2c0fbbvuhfWHeMKNsY0ExLAGAgthJBLCXCBC5f09gshBZLrwCXc3OQSSiCQ2AkJISEJIQkl9GIbbIOxjbFxW9u7Llu8faqk8/tDOpojjTR9tsyez/Pss7szGulII73nPW8llFIIBAKBoDxxDfYABAKBQFA6hJAXCASCMkYIeYFAIChjhJAXCASCMkYIeYFAIChjhJAXCASCMkYIeYFAIChjhJAXCASCMkYIecGIhBDyHCHkBzavX0QIOUwIcQ/GuASCYiOEvGCksgrAVYQQYnn9agCPU0rlQRiTQFB0hJAXjFT+CqAewOnsBUJILYALAKx2+hAh5NeEkAcIIc8SQvoIIW8RQsYSQn5KCOkkhGwnhBxf+uELBNkhhLxgREIpjQB4EsA13MuXA9hOKX0/w8cvB/AdAA0AYgDWAnhX//9PAH5S9AELBHkihLxgJLMKwGWEEL/+/zX6a5n4C6V0I6U0CuAvAKKU0tWUUgXAHwAITV4wZBBCXjBioZS+CaAdwMWEkGkAlgD4XRYfPcL9HbH5v6JogxQICkREEAhGOquhafCzADxPKT2SYXuBYFghNHnBSGc1gLMA3IjsTDUCwbBCCHnBiIZS2gRgDYAQgL8N7mgEguJDRGcogUAgKF+EJi8QCARljBDyAoEFQsgHeqKT9efKwR6bQJArwlwjEAgEZcyQCqFsaGigjY2Ngz0MgUAgGFZs3LixnVI6yu69ISXkGxsbsWHDhsEehkAgEAwrCCH7nN4TNnmBQCAoY4SQFwgEgjJGCHmBQCAoY4SQFwgEgjJGCHmBQCAoY4SQFwgEgjJGCHmBQCAoY4SQH+as3d2B1t5oVttSSvGX95oRl1UAQEJR8eyWQ+Cznne19mHjvs68x6OqFE9vaoGiavt8Y2cb9neEAQCHu6P44GB3XvvdcbgXW1u6jfN4elMLEop2Hrvb+vDO3qMpn3lx2xF0hxP46Egvth/uyel4G5qOYldrHzY3d+Ff245gza52471oQsHTm1pM121rSzc+PKQd41B3BG/sbMv6WM2dYazd3QEA2NvejzW72kEpxVMbmxGTlZzGnQ2tvVG8ubMdqkrxyvZWWLPeX/rwCFp7sruneF7Z3oqWrggicQX/2HwIW1u60dTej7isplyv3mgCz209VPC5OPHaR21o6YoUZV8dfTH8Y3PqWFu6InhlR2vO+1uzux3NneFiDC0rhlQylCB3PvvIOkyoCeCtb30847bPbT2Mr/zhfexp68fXzpmFn724Ez9/ZRceu3YxPn7sGADAWT95DQDQdPf5eY3nDxsO4LY/b8HR/jiuO3UKrn70HXgkgp0rz8PSu17Ka9/bD/fg3J++YXz2mc2H8KUnNqG5M4Jbz5yOT/xP6piP9ERxw+oNOH1GA97Y2Z7zcS97aG3Ka+zz9zy3A4+9tRe1QS+WzdSSDC+4701jm/N+9gY6w4msj3fGj1+FolI03X0+zrz3VQDAA1eegK/98X00dWjfVTG55P41aOmK4PufPA53PrMN9332eHxywXgA2gT676s2oLE+iFe/cWbW+6SU4rpfr8e4aj+WzRiFP2w4YLz3nx+bhgde3Y2g142zj9Pus689+T5e2HYEr3z9Y5jSECrq+QHA5x57B1V+NzbfsbzgfV10/1to7ozgE7PPhd8jGa8v/9/X0ReTc76f/+2Rt+F1u/DRf60oeGzZIDT5MiBbjaWtLwYA6AzHAQAHdG2iO5Io2lhae/Rj9MeN1xJKan2kvpiMrnDc9Nqf323GgaOpGs5lD5oFbluvdox2/Xzs6I/JAIDmzuJoczxMCwvH7bXsznBu15Oteng69HM72h9Pea9Q2P2yT19hsesJALI+lgM5Xrde/Xof6Ylib3u/6T22r3BcNl5j27DVWCnoicqZN8oCdg9Zy3z1xfLfP1tNDwRCyI8g2I3llTRthMkWFyFFO4aiasdwuUiKGYDniofX4qyfvG783x+T8dUn38dVj76dsq31YWLHcLucx80EZ5pN8oYJQo9Ugp3rsInRI5XuEZXZdeTOg103Kcd7oqtfm9hqg94UE1NCv+/4c2HHSfcd5otqM2kWA4rC9zsYBSGFuWYYk+sNwzQSn0d72FRDEBbvQZO5hzeWRlvZ2mK2kR/RbcBMu8zmGJLLWQAqlI2j+EKSaZ9SKWYQHSaASzmRMCWaP498z+2oviqrDnrQb1nhsH3yAj1hTNTF/37iJVwdFIpcogkoHUKTH8bYLfPT8es1Tab/VcqEZTE1+aQAjiYyOw0TigpKKR57a6/xWnNnGAt/8AI2NKU6UwFAUZITCe8Q5ZGV4p+bcXy19Fo20+TdJTwGWxHxWnu+GjYzAdYGvQhbVl5M6HrdyXORjfMr/vdTKhNQMZTwdIpPqRBCfhiTr1IwVXd0KQOkyaeTF/0xGdsO9eC36/Ybr727vwtd4QR+9VZT2mNILoI9FvsvwxBWJRAislI6UwODCSpPSVcLqROh8VqO1y2ia+8hnzvFvGaYCTkhzyaxYt571n2X8vvpiebnx2K+qlKOzYoQ8sMYNUfVYmJtAABQ4XObPl8aTZ4YmrzPLTlu3xeTTY4/APDqAoY5CKeN0ial+pA35RhOV0AuwQTGMEwNpbSXD4gmnzoRGiugHK8bPylZHdLsPR+vyevXsBh2bqexFHuC50c6/44X8toHU0oaSxBR5IQQ8sOYXIV8jyWKphTOSd6ZxzR55gOwoy8mmyINvG6XISSYRsjMImyS4jXQqEOEi0pLp83JhlO0hJq8jVO02Cg2vg32/eU68cc556rV7pzUrFPNNaUgbuPoHSr06RE/NQHPgB1z6F0FQdawhzRbpYs5xNjjZUTXlFiT91s0eT76oT8mp9gpjXHqgpppZuxTfHQNO8bYKr9pH0yIFPPcGKVyvPLXxZhISuCYZNhF0uRrk0+ksbHb2ciN77QEst7wARRZyLP7MZ0vLFMwhPFMeJxXt8VGCPlhTK4hkNabk2m7xQ2hTAqJaMJek49yIXY9UbOQj8uqEePOnherZmjS5PV9WQVuISF6mR5Up/j4QuGvi8zMHyXU5Itpk0+kEazbD/c6HrsUGKajEmnyTgEFz209jJPvejlt/gb7jnnTVakRQn4Yo+ZpbmEyzLDJl8DxKrlcRry09Ybuj3FCPpIwtrvu1EYASecUC4NkcdZs3LwAj8TtoxWY2cFq788Gu+QtHpbUU2wtlE2KABAfSJu8K9Umn2toYzZ2cP5ysWNnewk/ONiNc/73NdtkuZSxyKVxurOxMqWkwuc2TcKPvbUXh3ui2GEzqTHYdyw0eUFWMCFI8hTShk3e4S6glGYVBmm3T5MmbzHXRDhNuLM/jpi+XbVup+zQhfy+jjB+8Mw2Q+AxktqmCwcdsn3ZOHa29uU0fsCsUdtRKk1e5swaA6HJKzZae942+QK159beqOn8rbz0YSs+OtKHax57J+uxFMNcY5dYlTS5mPfPlKVImvsjGYwgNHlBFiTNLbl9jkU0qCwZxmGS+PIfNuGSB9bktG/eBOCkyYcTyRC7znDC0IyYM2pfRzIs8rG39nI2ed0myoXIPffBYe09i1pttfNPH12R9TlkmtisQr5YWYx8Ek8yFLU4jyilFGt3d5iEVjqbfM7RNXLm3AG760QpxZeeeA9LVr6E/33xI8fPBr2aomAtmWAHc7wWQ5PnvxM2fD5qjD8l9jymK3cQ0z/rFUJekA1MSOdrU2crASfn5NObDuLDQz05JZcwAUwI8PnfvQcg9YbmzTXdkYTxUNZV+AAA65vMVTCtGl4iC23TKqhZ2Gg2xBLZny9QPPsyH3HCJshiCYM3d7Xjs4+sMyWd2Wai5plExvaVz6349KaDAIA/v9viuA0TtqMrfRn3l865qaoU7+3vzHpitkteYq9ZNXm2y940MfRRkQwlyIVcHKe809Vqk8/0+VwKmDGBt7ctqXFZNXl+OdsT1WzyHokY5horTPCwcTMhzE8+1kfW+nBmK4af23oIP3tpZ1bbGo7hIoUD8ufDrlGx3CV79O9jN/e9GPcEd4x8k8iyiZaxe4vfPt1EzJzx2ZQs+MN6rQJmyJu6v79vOYRLHliDl7dnVyLYVEjMosn7PZLpnNjEmK4wGvvsQJawEbVrhjG5hFDaVb1T03yenxR6ozIaKjJrUNrntOPs4xxk1iU8v5ztDidQG/TC55YQ8to7o6wPNntQuvRqj3babixHXwLj5t++m/Nn2MqiUHiHL/NnFEsYvL1Xq1fPa592K5BCbfLpyHQu6VYtbPXXG5VBKU3rh9p0oAsAMK7Gn/Ieq9ufbXVSk7lGF+lJTV4yrQiyMdfk6uMqBkKTH8aw+ysbTZ6/uaxx8nb0c2Vh0y0/rTDBEeGO57Msm9n+GuuDaO+LISYr8LldqNUzWp0wbKK6KYOFqo2q8JkEyGNv7sUdz2yz/3ARMR76HM07TvCafCbnb678c4vmu7CzyfPqaKEZr/lQqWvw/els2bpgVVSasf5LnX4f2WnybFWZbd2nhM2xnByvzBST1lzDJm/Luua+l3bm3VAnE0LID2OUHByvdg+GEcZmc7/zJpXeHOpys33ymajWKIffv6PVqZlQG8D7zd0IxzQhP21UBX540ZyMx2APygsfHAGQfKgZP/j7tpTPlHJ1/Mgbe/L+7H6u6qbMrQjYyqtY4z5zltbcpI/zh9gJZj6ZLReY4zW9rTv1PZVS9OkKBT+2lP1zY81Ux51N/rwgXbenAx8e6jGEfLoJhcfO8Wpo8m6zuYbd831ZmGt4KKX42Us78dibTVmNKVeEkB/G5FJgjK/xzR5Ea1kE/gENm4R87po8r4laAy6O6I1FJtcFAWhOQabtf2rRRACa4F50TK15fPpv9qCwGPipo0IlqYGSCXa51uzWKmHOHJN9BA/jjV3JNoFxmXJ/F9dBx74XvnGHXay6XKBNPlfCccW4jukEL389MgnoiI0gveLhdVjxszeMc862wJjd9+Dk2GXHTacUReVUM1x/XIGs0rzun2wQQn4Yk0uBsaiNSSEp5JnQT77HP0i5dNgxNHn9eEGvhA8P9eLJ9cl2cMeNr8KxYyvx76dNAQC09sYMbT/odePuS+fhiZuW4tYzpzmcS9IpObE2gEAWiSWldHQxf0WlP/d6JKZEJE6TTzoyi+vU5Wu9263kjPLDOYZuxpXMKw92nF9yKx8mEMdU+RBJKI6TWy6afDonMBtntqtTfgVsTYbyeVymY2Ql5G0mIKZE5XP/ZIMQ8sOYXGqmd9ikWlvtkqbojkR+5hpDk9c/Xx3wYEtLN7751GZjm+5IAlUBD8ZWB4zX+NIHVyyZjJljKlFluemZwGMTSFxRDWfdIDTcMR56wzyQxyB4ocZH6RRbk2f742u92zle8y3Tm20tmv6YjP/6x4fG/30xTcAtPqYOAHDPc9vT7l/bR3p/RbqMZXYdCtHkYzaavKyoxsoy3b7ZZ/kRMvNOpb80cTBCyA9jcgmhtNPG2QNpLRcA5G+uUTkh73W7YDey3qiMKr8bQY9kRPbYZQBWOYRUMtMTpZq9PxsfYbHMOXaCvKMv/z6svZzQ5e2/1izfQmGCz1aT565N3jZ5ff9KGilPkewAxmAKxMXHTwAAHOy2j3qJZWmuUVWatmRCXM5Nk+ePxb573ibPeE+P6AFSq73y2K2o2bNZIYS8wIpdgSknIgn+ZtV+Wx9I2RTCV6gmr8LvdtmGusVkBT6PBJeLGLHRdjXnnTQb/kExNPmsR1gYpuU7paCUGkLeOoZsNHsmRLxul0WTL250DRNuTHMG7MMejXsqz+gaJcPkZG1yzu6t6oAHS6bUOU6YCUU1HOy9aYR8pnBW9v2lE8Q8rNm98f/RMN7eq3Us46Nr2nUtvrE+mD5OnlNQGOweqMwhYS8XhJAfxig5xDTb1Vuxajz8A2KOk89ek2djiiYUxyJMsYQKn26DrzSEvI0m72Cj5CcgLQY/8/kXy5xjDZfsicoF9RTlIzHMNnnnyKd8YEK4mxOyMZvEHJZdnGsVSjZBpROyq9Y0pdjTmZCv8LlR4XObQnfN46eoCWr3QzpNnjfVJFeqyddYYl+2igs/KVEAl/9iLf617QhcJJn/QSk1etyePK0efTHZsQ6PnU3+//Tku4BDnkihiGSoYQy7j7JZWfMhkcnaNWYJwgt2ptERktkGysM+F5NV1IY8sHvm44pq2OAr/G6g276xSJC76T0SsS09zAp4DZRNno8aokhWzLQbA6XpE9UopVi1dp/xPx/tVOxm1IYGywk3u7Bavn1jLhiafJr4879vPpRiWvxv3QZfHfSgwufGnjYnIa+iLujFHvSnF/JpQoUB4Kj+fWVrk7cK5UPdmrlJchHTd9sT0cY0sVaLGOuNyrZ5H3Zx8hv2aWU8sgkgyAehyQ9jcslOtBPURi6M/pt3bjGNvDrgSbs8tqIYQl5xrAIYl1XjPRZRYGeu4U09TvvyuiXTw+akQZVKky9EGJvsszQ5aQP29vJCsBtnzCYWP2+bvP65TCUemPnjP5ZNNb0+rsqv94e1Vyjismpo8umia8wrCX11wY2JNRzPNmLMpBxxp2ZdMbT2RuFzuzBGb17jtFKwThq8Y9eaNFgshJAfxuQSJx/O0SbP/q4JeLJOHOE/l1BoigBmxOVkVExFGnMNj8ftshV4XsOsoL2Xi//A9jgZzBQmTZ6aH9IUm3yGY1njuXNt55gL2caxG12vcrTJs8lVzmATZ7Zwq0BzuQgq/W6Tz4BHW/1JCHikrM01/GcZrBRGXFazKjGQbhvCmQnbemMYXxMw6i/ZrRQo5bJ1bXIDrN3NioUQ8sOYXIpJ2dW4ZrZK9ps3g7B9Vwc8aTP4nMYEONci4UMfWURBJiHvpMl7JHMEj1MxtWzFZ6Zsd+tDn1IILQdBzYR8pd8NCmpbu7xYcj8uqylp+MljpH7v+SZDZdLku3WzBv99X6YnwIW8bkQTqu1qLKFofpx02j5gNtfYrVD5iTUb5SWSSDVzWqHQJo/aoAdV+v1s59i1M4+xVcmPL5tf9HaSDCHkhzG51P4O2yw7rfZTXgtjAr866M2YfGLeR3KfPofwRkWl8EqaJleh1xfJtFT1urXEE6sgTDq/tP+bOjLXG09Hppom5smSGgLE53YBlFqW8dnti2l/6cIPCyWhqKgJ2NcG+sEz29CiN1/JJWKLhwn3TGWXmfDjnfJ3XqiVsmDRVHbfYVxW4ZFcqPS70wYCpOsna6U3KmPu95/HI687l6UwTeqWU+Pv7aP9cdQGvUbYr52yYVc/ij1bpYqRB4SQH9YkOztlp8lbtWFrdA2vhZk0+Zxs8qnhjXYwRyuLKMioyevvW7Uhr9s8kfCldHmy0bDZBHLtKY2OKeZWEwsz17Dx5ZLezx76gEeyncC0cWe9O0dUlZqiU6zsae/Ht/RktWRnrxzb/6nZmWvi/KSow64dE3Rn/eT11P0rFB43wbhqf9oKknbRNdaEJubQ//36/eiLyVj5zw/hBP99//yVXbbbUErRH5dR4XcbQt7OXGMXI8+erVCJwieBAoU8IeQOQkgLIWST/nMe995thJBdhJAdhJDlhQ9VYCUXTb4/LiPo025u9hhYZQqvhbG/a3Qhn60ZQraYa4hDeCObcJiwdzLHsAfSK7lAkWouMTR5/X+mQa26fklW4+VhmnR9yGsb6WBth8jb5H1ubXzpbPRWDCGvn2O2lRFzJcE50Z1gSUrMVJJrHXumIGTbGpHX5FkkT7q0/oSswitJmFQbNFYdtttlocmzePvdemvI8dXOtnB+5baai4RaOrXOdGdHEwoCHokz16QqRuZ7R7tezBQ6ZIW8zv9SShfqP/8EAELIcQCuADAHwLkAHiCEDFzn2hGCnIMmH44rCFoEl2rY5PX9OUTXZFPeNfm5pKBysgEDSe2NPeBO4m3W2EoAyYfAqiF5JWKaSCJxGUGvhGUzGrIaL4+RQewijlLO2hicCRA2SeWiyUe49HgKe39AMcQ+m3hqg86lnJmWmW+XK3beXeHsQhN5BzeLoqpKY7KIKyo8boKaoCdtIpO5kYx2LnzhN0CbxIHkd5muNn3EoYz0Tz9zPHccbTLweyRU+NxwEXtN3m4CZJp8unMvlFKZay4C8ASlNEYp3QtgF4DcVStBWowCZVmoXZG4gqAuKJkWkVq7xkaT15f42Uat8EKiPuRz1AiZkGdjd4osuf282bhq6WSMqfKBUorfcNoUvx92Tv1xBUGvO6/m5pnaKVIKNHGlgSm4ptG6z4DvQZpp8cO0RLZacboGW1u6sXHfUZxy10uOjcvTwb7X2pCzpsy0zHTlp52glKJdz1TNpger5CK211ixcQCz/ccVLey2KuBBLE1kjF2oqJMmz4R8S1cEW1vsa7nbNZ+55PgJGFvtN93b0YQKv0cCIQRVAfuJyJR1rv9mQr7CV5riZEBxhPznCSGbCSGPEUJYbdgJAA5w2zTrrwmKSC4FysJxOaXzkpERCC09n7ensvR0ZmPMNoyStytb67zzMJssMYS8/XaLG+vwXxfPM7azNtPwWJy74ZhsSqJiZCO0mJCRXM45tLv0JT4jZrHJO9lt7YiYbPLU1lyz/VAPLrjvTXzqwbU42B3FGr2zUS4wTT5ddy8mNNMV93Li9Z3tOW3vZJZaOKnG+JuP6FJUCkphOF4B51h52cYmnyrktevQxhXtu+C+N233Z7eCtU4kiqpNQszEV+X32Mbh720Pp7yWNNeUztCRUcgTQl4khGy1+bkIwIMApgFYCOAQgP/JdQCEkJsIIRsIIRva2toyf0BgkIvjNaxruEBSi+A1p1t/9y7+8/Fk6zuZc7wCmcu7Wj8HAPUVXkdhyYQ80+iysflTpMYSW23ymiZvicHOUqnn8w7sxkOhVZz0uZMRPQlekwc1auRr26c/J5ZYZdjkbY650zKp7G3vw3/8ZkOK2SgdTMiPTWN7ZuYaJY9WhuEcHPMMu8VSpd+Duy6dBwBY8IMXjNfZxON1u4y8Cielw85cZnW81ld4TftN99loQklRom5bcax+DsTYBkiaJ6sCbtvomj9tTOq9RtNv5ni16WJVLDIKeUrpWZTSuTY/T1NKj1BKFUqpCuARJE0yLQAmcbuZqL9mt/+HKaWLKaWLR40aVej5jCiYUMgmBf1QdzSlyp2aVOXxzy2HTSYZRaWQXMSoLWMn5GOygic3HHAsi+tUewZIar51ugkhUwgZO8N+i13TWukyEldSnFg+t5RV5ihbhUgu4lgF8mh/3LDpAlx0jZRbWznAosnDPrrGWrXx/ld24/kPjuDJDQdStnWCaZ4VPjdObKxNu00+NvlMn3ByqtvBN/Nm14NdY48eJw84mw/tmrtbNe/qgMd29dtkY2qKyaopEuiqpZON0gWMiMWBXuW3N9fMHV8NAJjSEDJe64vKmh2/RDHyQOHRNeO4fy8BsFX/+28AriCE+AghUwDMAPBOIccSpCJzmmc6uvRU7lam/THZnubplHUhzyYGu4Sov79/CN/802b8aWOz7T78HvsqlACMOPnPnDgZ3zl/Nq49ZUrac2DjtmqNXkucfH881VyTLpSTJ2muIaZOWoyv//F9HOqOGjVJKKVcdI0WBmnXLs4JXsgD5rIGDKdwwVxqy/ARQH+8+RTcft5s2+1Ov+dlw3GaS1KX3bXisXPAO0Vd8UL+sD7BGX4PiThq8s9uOYT7X9lla26yvuZzuzBKN101VPiwWO9AdtgyoQLaufGRQOy+5WG+Fb/JXGMTQikrhgmTjejd/Z05hSjnQ6E2+XsIIVsIIZsBnAngKwBAKf0AwJMAtgF4DsCtlNKBb1Ne5igsBT3Dt8gExRkOESd2ZgJFVeHmSgHb3YisDsn+o6m2RsC+Hg2DCV6v24UbTp+aURCzycKa7eiRiGkiCcdSzTU+tysrm7zKTZp2q5O/vKctRvlVB2+uAXJr9sELB0rtHa9OAiAXfZuN0ZPhRjlwNIJ393fmsGcNu/hvnlyqK1YFkteWOXH5a+x0P97y+Lv48fM7zJq8g03e53ZhVKUm5CfVBXD3p+YDSBYvY7D8AruYfh4+SgrQVgr2IZSa3Z6f3jZxdehLRUGGIErp1WneWwlgZSH7F6SHKSiZHK/MfjtjjBaOaDVd2Nk3DU3e4aGKxBU88OpuAM71PdIlOGWrXfNQaA5kr9vFJSHpsf/6Ex1OyCn2zWyPpXLXM13IaIjzbSRNCSSllk06InEFP9NLzLKlei61a3KpJxS3TETpFn5MYOUyicQ4m7SdwLcrOc3GYF2RBDzJ7+6Arjzwk1QojdIBJB2l/G6tQt7rdhmrsYBHMsxvViHPnPx8zgR/L7FziFpWZJV+t60mH0koWl4ITd6vx45N7YBWbETG6zCGmVAymWtYfK5T5xm7xseKSjVN3iGagbfFW+3GDJ9HMh4w6zyUKcPVCvt4JKGYooSsBcXCMSVFc2SJSplgKxoXSS+sg5xJIWbU4dHGkW1Vyp2tvanHz8EenssSPyFnp8kD2ce587Dm1E79A3gh+X+fPR4PXXWC8b+1Rs6ssZW4/tRk71/AbJNnq6h1e47a+jBYlUmP5OLi5G2EvB4aHPBIqAp44CLmstFA0jHOl9ywu2/DcbPjNeRzIxxXUsYXjWsJU7wqr1KaNgqtGAghP4z53xc/AuD8cDFYaKTHZd8PNRI3PwSSixiafMAjwUW0TNIdh5OCiRcG1u45DB/X7SiYp3ZtxdqMhO8M9e7+TnT0x1Mcr940ZiMe3lyTzs7MC4OETI0GKED2mnwHtw/2zDvVrmEZmefOGYumu8/HuGp/Tpp8zKLJZ0MucfKGMHTYP/99XbhgPM6dO844Z+vEI7kIvvfJ41DhcxsCO26jyf/+nf2Y+u1/pkx27L7kz9U68Qa9biMxzO+VILm02PYui7OUafK8T8GUxKWfBTO7scmMhUNalaeInhXLE5NV214KxWW0RPwAACAASURBVEQI+WEKn26dKRmKOZ6cyuhazS1eyQVFoXpjBIKQz40HX92N5T99Hfv04lG8BubkHPR7XMYEc6yeucofI1copYgkVNODwguWSx9YA8DsvANYolJmqcXXUk9nrjn7uDH6eIC4osDDTTR2NmE72m1CIO00UwCYOkqro8OuuVaJMXdNPp9rng0xWYFHIo73YbpmGE6ri9qQxxDY7P71uV0pWdv7O8wKBpuAvZLLsDmx82emoQqf27DJsxFX+NwpwQVs8uJ7udpdw1++sVfbTh8bU2isXa6iieQqk3Kv+bNUQvJFCPlhSmuvvYnEDsOm6aBpWTUOj6Rp8qxIFS809+gFwJgmOWN0BbrCCVsbpM8tGQ/oRQvHm4p+5Vqrg8mPmEWT99nUrJ85xjyh+CQXdrf1p615Apija5wE9HNfPt20/1hCC7FjY3DS5Nt6Y6aJhk+WYZ91svSM0zV5ZrMPeaWcunWl2uSLG66n9fOVHPdr53hlmzopHrVBr2Ej5801LhfBTz+z0Niuvc88WTJtnJ88EpakwRAn5JNVIFMTmOw0eX5VyM7hnSbW89WsyYdjVk1ev1e414QmL3CE1+QyxYDLilmTs25tTdyQXASKqhoPBS/kD3ZrgrI/JqPK78bXzpkJwD7G2OdxGcf2SC6cNKUeAPDJBePzskNSmLUhAMbDCgosmVKHCp8b584da/ocK6h48f1vmV7vjSawak2TMQmqWYSk8k5dCmqKozY1hUDye2ntjeLElS8ajuqeaAI//Ps2AMCzXzrd2N7J8dqgnyPLUwp4Jdv+AE5YY/nt+I8zzJ2aMt1TkbiCJStfxHNbD+uN2Z33nU6Td6p2WR1IhiEmHa/a9zKDUxasygXrYeuWSEqcPBP8FT7JuG/4JtrWhiWGJu/geLXC7kumyVtXW6yBPQDjIUzXC7lYCCE/TMkl/Zxt69QI4mGbetqy7ngFzFr37X/Zij9tbEZ/XEGFz40JNVpiyKHuaIpJxOeWjLZwHilpujmeS1/PFt7xyguN6aMrDNtoTFax6JjUZB8mEK1Zok9tbMb3//YB/r75IACzJu9Epd9taHDhmIK/vX8QTR1hY3x2jtddR7Ss1ee2HgYAbGxKhinOHldljN+pTG/IUtsm5HWbOn1lwrqSszu721bYx8470dwZRmtvDPc8tx3RhJo2XNZ+AtBG4aTJ8wlF1pXIGC7ruTuSMK2euiJxeCRzbRz2Prv/g1436nSbPHOaanXqU80rgLOQt47cr7/HFAFrQbLkqi+Z5a2tgoQmL7AhW9svv621wQbDWoKWQhMoTNhZs1F//vJOhOMygj53suxBVE6JDvG5Xcax3RIxnLCZWuw5oVKKnohsWlmwB5BCM+XkErXDzDesOXOmAmWAecLr6E+1q5tKDeuXgz3sLHnt7b1HbfftZK6ZpJdKmDZa02ADXslkCqCUpi1zkI0mDwDzJlSnjN0Jdk4qpYVp8g5jqgq4DfOJNTqIr6bZE5FN2nxnOGGsDpjSEVdUeCSCU6bVG5+fOiqEhgof/uOMaQC0yDOrkI8ZUUOcuSbNNWQaPCvpbbXJa6s+dr/SZG19ockL7OCF/B83Nqd9yJlwcW6hRzChJmD+jJIU8ta4cwotKSnkc5tCLK0p8XwSkldycV2Hcr/tCCFo7YkhklDQyKWFa+9pv2Oyarv0dZJXxgrHEqfupMhLLmJqN8hWCCdNqTPGYFs7RX+NxZM/9Jpmtll32ydM43eqG/OxWaPx6OcW46tna6axoFcyhOy7+zvxl/dacOLKF43Emnuf34EfPLONO775+3eaw75zfqo2f7Q/butnYIJVMbRRydHOzoQ8HxPvFCfPqEyjyUsugj/ctNQYx+1/2WJ8rjuc0BPkkvvafzQMFyG499ML8K+vLEPAK6HS78GG75yFCxeM14+X6sxmmjy/SvHZxMlb32NKiNUmH5MVk00+miEqqVgIIT9MsZprfrtun8OWMGnTgFZJcX1TUpvsDCdQG/Lg82dOB6BpcYpKje2t8fWSiyAc06paVhi1RBIpQp53xHk4IZ9r/1AG06wm1KQW2mINPexS6J2coWx/zInJR9fYYV3RsPjwG0+fqo/BvmkIe40JjRmjK+D3uFIKhvHzQyW3Ygh4JHxi9hhDkw163QjHZUQTCi59YA2++uT7AIAtzZqQ//kru/DYW3tTzj9TCGU11zmKjf2EH/4L33t6a8q2LKOTUqRo8tbrx2zVdtfVaUxVfjdisoqYrKTU7AeAk6bWoz7kRU8kgec/OGK8zvcPZj6cf2w+hJisIuh1GwmBVip8npS2gklNns/LcL6GLKmNZVzbavKe5GqahekKTV5gS8IiuA53O0fbyAqFiyRDLR97ay8+/dBa431FpQh53fj68lm49pRG7TMqNTRua0iiRAj6YjKCXje8bhd8bhd6Y7JtA2ZG0CsZmmo+5hr+EzWW5hdJzcjeieWUoMRq+jANzkiGchDyzDTFJi++LAGzq9tNKOxhZkIj6HMbTmh+/Lwmv2RKnfG3VRAGvBIiCSWlCJZTU3Gr49Lp6lvNdmzcT6xPLYbGNHlKNVuzSZO3rNT8dpq8/ttpQmUlrp/behgJvemHNTpMc86m+ib41Va654In5NUiwfjvz1pdErDa5O3Hzla+fC4Dq3PER4MlQzSFJi+wwWoWsIaSWbfNlO3I25pZbXO3TXQNoNms++MyKnzJNO6+aKq5hifglQybfK79Q63Y9SrVtDbVdunLXys+yalDb3TB4qONKpQO9gxr+vkWvdEEEwLMzsqEqWET1gVHBxcSaKfB/nXTQePvdGWBgx5NIHVYMjStkxmbvOKyChdxtn8zeCFPafpGMSwiy84mb12pGULe5vjWCYHB7rkvPbHJSOayKgeVAY8xUZv2yfmeDmUp5Nlqg49astPk7b63mqAHv7vxJON/ZpMP2+yLL1NtmIOEJi+ww/pAp0veSShUu/HTKNBW26jMh1DqZoqpo0KoD3nR0R9DW2/M0Kgr9OScdGn5Qa/bmATycrxyH7G2sSMkufy11eS5a8MnvLTrjlNDk89grmHFs9gc8PL2VgC6Jk+02jcJhab4Pvjv5hev7UZcTm2qbiVdHgETSFYttS8qmyY0VhYgoaiWmivm82PfvdVBmq7NHhOGcVk1Jlem2VqFOZsEzTZ5tq1Tgl7yPP62SSsMZ71mVX433tAblvzw4rkYrYdFuvWidS9vb8Xv39kPAPjvT81zPBcgeb35qCVDk+eunbltofZ7wcQanDItWfzPK7ngdhGTJs8LeT4azLr/UiCE/DBFttjk7erPGNuqakY7OK+FU8CkybMHYEJNAJccPwHtfXFEEypO16tassiEdP1Ng14paZMvUJO3K+gkqypUal/igfdf8NopS7bhnYiAc3SNUyNsJsQSlmgJdlReyD+x/oA2EaRJu9f26azdse/DqqV2RxIm7bG1J2YcP91Kjk1qhBDco1dkBKitKYTBzrU7ktAdiqmNuRnsPbt70Gl1MXVU0rm+Xg85tV4T/j6YWBMwrgubDPpiMv72vrY6umD+eMdzAZJ2dDvt21xqOHW8Vl8NIcTkHNf2ZdbaKahtiGYpEEJ+mGIVqHb13vltNTuls6BngorJN1a7Bkg6ARWVmgQdy/wMed3oz6DJa+YaswM4F/ixW5fMhBBT6rsVXpNnQr4vJhtp8+w1FgnkqMnrQsX6LrOz7tNT7JkguOqXbyMcl01C3u3SyhjzwsIusSld2CETSN/mokoALYuW3xcr1RtX1LRRIbxQvvzESYZGzGvy1lo+7P6TVYqj/XFT1q/1+2XHtpvcnaJrlk6tNwqVWffD4BO2JtYGjExT3ibPyJRhza43HxFjbe1o/Zth1/Am5HMbmnxcVo17rMInmaLB7M6r2AghP0yxCnm7sgLJbSk8GcoRnzCZSyKi2krBqskrKjVFYLCeoSGfG/1xOW2CFq/JZ1MNMV/s7JtmIa9dp5c+TEZkGJq8kfFqv+8qR03eUvVSnzA3N3fj9Y/aTMcfVxPQJl138iB2q7BAHnHnmiafnOxZnaFElpo8D6Xme8pqn+e/656obBJ+Vjs7E2L8cZiTOV3zk+Mnm5PmrGamuJwcw4TagOHwzEeJMMw13PVjeRf8YfnryO5nq8+K7a8/LmNLczdmfudZPK37W9gzI2zygoxY29PZ9ZRkyIoKt6XhtZXP6VE1TGNO6J8BkrU4VGrW5JltOORzIxxT0mryfjenyRfQ6uzjx45OOQ/+Xzv7Jm8OuX7VelNp5CWNdYbGmim6pspvtskbx/S4TM1M+PPb2x5GXFZR4XPj1On16I/JuiaffLBtNfk0jTasFT0ZPRZzDcudiFtt8pbPWU0mneE4Xth2xCTYrULeamLyuZPNMBwdr9x1SRbNcxZB/KR6+eKJKe+z7+GLH5+OoNdtCEtFpTm3MQw4mGv4yCnAfG59hnaeOvmz+kKsrs3v3tZ8A3wzddbKMpTmuy4GQsgPU6zhiuG44mgTT3Ax73ZMaQiZHrbemIzth3sNJxPTxGSV2mqzIa+Evlh6m7zLRYxJIFOUhx1rd2sONrtOOnyUSSatKJpQ8eCruw1n68S6ALrCCVBKjWgYp+gaZ5u8ZDS4AMyCYndbn5EEU+Fzozea4GrQa9hp8unstE4TQE8kYdoXuy7W6Cprn1zrpJtQNBPMPq7CozWG3BrC63U7J7sZ5hqJ1+Qz50xUcWaQeRNTS2EwM8l8/b3XP2oDoK2g0ik9dtjZ5KO6Js8rL7wfIFncLHXSZbkM7Bli0W/1FV5EEwpe2HYEz245BMC5z0OxKO3eBSXDTqD2RmXbwl8J3Qbs9DjxNsGDXKVGZkNlWuf46oCtoGNNEviHwVZOsoiKPDR5tnJZefHclPPY09Zn/J2NfVNWVcP2OqEmAFml+vi1950cr07mGo/kMpUOZqGZgJZ4NmN0BbxuFyr9HsNB7eWEm11sfTY2eSs90YSpFgtfxZH3AXToAufTiybijxubHb+PXa3J65pqrtFCRXlfCCs5bY3KYRMvL/zZ/ZuuThB/ve1WaN9cfix8bgmn2bS15IX8504+xvEYjGS9GXNEjN8jGXXtr156jGnyZTZ3OyEd8kl48cMOw2nMqAt60aRPns/qtYxyrciaK0KTHya8sbPN1J6MCb1TpyeTapxC3mRdkw/53Bhf7cc1lpv+qqXJ/0/lHhimZc2dUIW7Lp2Huz41Dw0hH6yEvBL647KpwNZtK44FAFwwfxxO0hN7rjxpMiQXMUrn5sK/nzYFdSEvzjx2NCr9HtSHvEaG7mdOnGRsxz+E8yZU4/LFE/HVs2ealsk+t2RoYSwevTuSSNrkXcBnl0zCpLoAvvBx7RgukvRbjKsOpAjh605tNP5eOq0eXsmFExtr0dTRjx1HeuFza12NtFWD2YF30fETTPtyEa1D0vhqP75y1syUa8ELeX4+6o4kENGF1MTaYFLIK9SUSHTajFHwSi5ctmgiKn1uXH+avYOTN2ulaPIqNcXy8+fj87gQ8GiVHkdX+jC+2g+v24V/W5L8nuZPrIHP7cLFC83nzsNrzXYrm7HVftx16TzjvXs/vSBlm4euOgHf/+Qcx2Mw7Mw1TJM/4ZhaeCSCi483R+gsnzNW/57rYMXJpOaWXFgwsdr0mrVsSLERmvwwoK03hqsffQenz2jAb/5dS7pgmtTjNyzFv7YdwY2rNzg6XxOKCrfLBb9HwprbPoFoQsHqtVoZhOe/vAyzuIYeVy89Bq09Udz38i7DFkkIwWeXTAaQ1DD5eOGQz6076jQB8/gNJ+HU6dpk8fN/S7Z6u3zxJFy+OPmg58ItH5uGWz42zfh/43fPNv7+zImT8f+e2mIaHwA884XTjL+/+IkZaPzWPwBotlQX0bZl1Qi7IwnDXOMiBHddOt/47NfOmWUay/iaAD784bnG/gDgq+fMgs8j4cfP78Alx4/H6uuX4K5nP8T6pk50hbvRUOFDpc9tmFN4oXjhgvGoD3lx5S/fBgDsuet8AMAavbaNFd5c84urFuHkafX41VtN+Mm/PjJMNBNqAobjNS4rpu5VZ8wchY9WrgAAbLlzecr+/3rrqVjxszdMPQusjakTsoqQ1w2Xnh/gc7vwp5tPxmUPrUWV34OXf/gxAFpCGCEEO354rslxOmtsJXb81wrb82PwZpBswgwvWzRRq21UH8Ta3R144NXdmNJQ4ehj4XEKofR5XDhhci12rjwv5TNnHTfGuI7pxm7l6c+fhssfWot3mo4i4JEy9mguFCHkhzjd4QROXPkiAODNXe3G6zJnZ2UmFLsO8YCeDMM95LxJw+5mZI2Nw/HU/XkkF+6+dB4WcOWCWc9TtpIoxLGaLywhiyUs2XHSlDq8vfco+mIyJEmLZWbXrjuSyKrUcDpuPH0qRlX4cNECTTvlTVvtfTFUcpqp1eGYTihYYdEcly2aiHPmaLXz6yu074wlSI2v8SMcVxCXVSQUalvTxwk2lrbeGHxuF2KymqJAsASrkE/LkfB7JIzTi9zxJRmYYM+nUQkv2LMd/9X6qnTp1HqsmDvOpMCkPZYeBsvf84V0baq0yeX4+b8db/zNotRKbY8HhJAf8rRx5Qr48q9GFis4Ie+gycuW5Bv+gbMTLuwGdIqWuULX6hksOoDZQfMtQFYInzphAlat3ZeSDcvzxE1LcfH9b6EnmoDP7ULAKxnn2hWOg51upsbojJ9+ZqHJ9ut1u3A5ZzqyhhLyE5A13tpOKDgR9Lrx5v87E6Mrk+YSdg8c7GJCXhO4rN56VQ7CpFKPFlGp1pSluTNiG0LpkVwIeiX0RmUEvRIm1ATw4lfPwBRLldBikGvCkEdyYZ7FLJIOl94N7L6Xd+GrZ88EIUS3yedn0bZ7rvgwZdZIvC7N/VoshJAf4lj7rzLinCbPhIeTTT6hUgQdIlrs7IF+o+Z1djDHEevKk08p4UL52vJZOHfuOFNDCSuEaA2be6NaTXq/R8K4ak0YHuyKGk7rbBX5i493ticDqQ61qiJp8oBmc7fb96HuCAgBxlYlfQ1OtXKc4LXLCp8bIV2Q87AaPcm+ptrv6aMrUApK3QeVpyciozroQTShoMbB2Z4Ju0mVv66sJIhdHaZiIxyvQ5x+h4bNWoKLJo3YA+4UNiYrqmMylJ29kpkrstWe2ETBVhKDYa6p8ntw8rT6jNuxuuFRvcNUbdADr+TCkd5oxto1ufLJBePwqROS8d18tIg1CihXIW+FafKHuqMIed3GCqU7ksiqQB2P5CLG6sznkfSooFRzjUev0QKUPkIkX406F+68UHPQHurRooQK6b8asFGe+PLR7PtyipQqJkLID3HCDpq8rCbNNUGvZk+0Nj1g5PqQnzFrFK5eegy+e0F2LeFYstRgmmuypVKvG87aCBJCEFdU/OK1PUaCT7bmmozH8ntw76eTDtx0mjyr72Jt3pItbALZ296PvpjM+WkSiOWoyQNJrZNFBdmFUHoll1FIbHye486WUtd3AZL1cphvqxCbvFVP2PS9s01mUqbBF+teS4cw1wxx7LvLSKaStoQQw/Foh6ykT4ay4nNL+OHFc7Penmlxnbq5JlOFxcGECaxoQk3Rnt/dp8U0ZxONkS3swa70u002ebt4/qduOQWN9cGU17PBmr/AO5StjvdsqPJ7cKQnlhTyMWsylKZk/PjT8/HLN/bimLr8xp0t+WrUucBPjEBhmvyJjXVorA+iqSOMm5ZNTemBwJznpY6sAYSQH/Iwb/9x46qw7VAPOvsTGFstpdQjqfC5HYuUJdRUTf47589OW544F5iQZ7W9S1mbplCqAh6E4wp6owmjENcTNy3FFQ+vw6t6xqRTxmu+/OOLp2FUpc/kiLVrfG3XhDxbeCE/tSFk2JK7I4mUsgbZ0FDhw87WPngkLYmr01K3Xau/48Ip0xpMZXZLhVPceTGxmj21OPn8NPlJdUG8+o0zoarUNjGQBVGk8yEVi6H7NAoAJON2r1yqRbRc9tAaHO2Pp5hgKvSCSHZoWpf5Trvh9Km4VU8mKhRmv2WafKmr6hUCE4ZHemJGvPnSqfVYMLHaqPVSbOVqzvhqjK70m1YOuQrdTPAa4b+dNNkw33SFExkLlNnBEtYSimo4q3niafw8xeThqxfhlGn1Ja/vAqRGqRWiyTNcLmIbPnrOnDG49pRGfN2Sg1EKhu7TKACQFPKT9GiK5s4IXvzwiG6TNycksQexN5rAg6/uNiJztHrypfuqmZbFSvcOZU2ePch9MdkUbjl7XJXxdzHNNTz8d1BsIc/TUOHTNXA3OsPxvDT58+aNAwDsPNKnm7jsHa+l5pw5Y/G7G5fmFWefK2wS7o4koKpau75SRfUEvW7cceEcU1XXUjF0n0YBACASl0GIltzC6NHD4viHjO82v3rtPvz3c9vx/AdabYxsSg0XgtftgldyJc01w0CTB8wdppxS9EtFKVc7zElZF/Lqqz6as0BmTsi2vhgq/e6UBiIJhZrKJZcDbkkrItcTkQ0n/ED4AkrN8D+DMqc/riDokVDP1YzpYc40t8VcE2OavPabVUYcCK0r5JMKa+83QFSZSiUnrwkv8AfCcVyKieSey7RInsYGbdVXE/QapQlynVQm1GrRMopKUeX3IC6rpsYh+ZiAhgNVfje6Iwmu9d/AxeeXivL7lsqMcFxB0Oc2p8jr2hkfjx7iHK/WOhxadE1pv2o+0WMoR9fw15GvE14bGlghXwpN/vLFk7DhO2fh2LGa6akuqEXIALmfk88t4dYzp+HxG04yzBhX6bV1AL1G/RD+nvOlKuBBTzSR7NokNHmBlb3t/XhabzxcDMJxLWWctxN39MVsHa/WhtQshlmLrimtdl2nrzQ8kr2jaajAC3nWKAVIppkDpbPJ85TKJMRX26z0e4yywvl8/99YfixOnd5gCPn1TZ3GvTVQNvmBpjrgEZq8ID3/9fdt+NITm2yLe+VDOK4Yjs3bz9OSkzr64kYIG4MJeUopovqy+r0DnUgoKigtvHl2JkbpBbKG+oPPC3le0KareVNMWAasUwOSYsLb0r0FCCu+89Gv1zRBUSlUOvS/63yoCniMBDJAaPICG17a3goARgOFQonEFcP8cuOyqbhg/ji098WQUKhpuVzhd0OlWpehmK7Bv7e/C4++uRcASu4kYz6DgXBaFgIb38JJ5k5DAyF0AeDuT83D29/+RN7x17nAm9AKWcmxCpeA1jeWNfwoN8croN0HPUKTFzjBFxNj4YSF0q+baxhjq/w43BPVa8SnVpPsjcqmnqFv7+kAkFoRsdgwQTAYdWty5Z1vfwKP33CS6TVmky+1ndkjuQYkAQYwl1EoZPI9flIN7r50HibVBdDSGTGEfFna5P2auaacNHmR8VpEvv3nLcbfThUhcyUSV4zMTAAYVxNANKGivS9mMtewUrW90YRhrgEA1u+71PVkmC04bNOUeqgx2kbIhrwSbj5jGs7X48PLgcoiOcMJIbhiyWQ8s/kgOsPxrJpwD1eqAx70xxXDvzUQNXNKTfl9S3nw1q52x1rsufDn95IO12LsD9A0eb4c8CQ9tM1qrmEPdI9Fkz/crZmNSv1AMk1+OAh5Owgh+NaKY3OqQT7UqeCqHhajSmRN0Ktl0DJzTRkKeVZfiGU/D+Xs7WwZ/mdQIL3RBK785du4afWGou63WJp8OKaY2r2dNqPBEO68aaSKM9dEZRULJ9VgxdyxOHCUCfnSavJVOTS9EAwMfCOSYgj52qBWw4Y1Hh/K+RD5wnwzTMgLTR4AIeQLhJDthJAPCCH3cK/fRgjZRQjZQQhJbSQ5RGDFiNbtOVrwvvweF64/VWuKbM0QzJf+uGx6QINet9F5x85cw5xGfo8LE2sDRk/RUjv6BqKNmSA3eHNNviWMear8HvREk9mgQ93Jng9MWWntyS+JbChS0JNJCDkTwEUAFlBKY4SQ0frrxwG4AsAcAOMBvEgImUkpHXJreafyvPmQUCiCXglBr1QUTV5RKaIJNaWxwNhqP3Yc6U0pawDomnxCQV3Ii2Pqk23YSn2zFtr0QlB8mFZa4XObyjbkS6XfA0WlXC/f4S8ArVRxBewAc9P04Uqh39ItAO6mlMYAgFLaqr9+EYAnKKUxSuleALsALCnwWCXBqTxvrigqhaI38tA0nsKFPNPCrUKedcnhX+cdr+G4gpDXjVGcw7bUy86KEncGEuQOa0zxsVmjirI/Zq9u79NqFA1Et6aBhikrh3VN3q495nCj0G9pJoDTCSFvE0JeI4ScqL8+AcABbrtm/bUUCCE3EUI2EEI2tLW1FTic3LGWUM0XPna4KpDaSScfwvoqw9pKzNpXE9CiQ1xEO59wTEbAKxk9S4GB0OSFTX6oMa46gN/dcBLu/fSCouyPfcftehZtOZprmLLCzDWBMrDJZ5ymCCEvAhhr89bt+ufrACwFcCKAJwkhU3MZAKX0YQAPA8DixYuz7R1dNHo5c42q0rxT2vnY4coiafIsVtdveZhYggbfOox1h+qNJhBOaAlUJiFfYq1LaPJDk1OmF6+hB9Ny243Ik+EvAK0wm/zhnqimOA2DvI9MZHwyKaVnOb1HCLkFwJ8ppRTAO4QQFUADgBYAk7hNJ+qvDTl4c01nOI56rvZHLrCIA6/bhSq/21jSFoKTg6vCCJc0TyRaw2UZ4ZhWCqEuyGvypX0gc21jlkgk0NzcjGg0WqIRjWz8fj8mTpwIj6d4KywWwdXWVz7hhVbYs6XSoV0yOxcKVb/+CuBMAK8QQmYC8AJoB/A3AL8jhPwEmuN1BoB3CjxWSejjele29+Uv5PkEkaqAB3va+wsem1Nm4aUnTMCjb+7FshlmWyvfJCLolUxldQfigfzO+bMxb0J2cebNzc2orKxEY2PjkC5oNhyhlKKjowPNzc2YMmVK0fbLzDUdugJTDtmgViQXQdArIRxXipa1PtgUKuQfA/AYIWQrgDiAz+la/QeEJrAYgwAAIABJREFUkCcBbAMgA7h1KEbWAGZNvq03hlljK/PaD58gwlKjCyUh22cWzhlfjaa7z7f9zCs7NL9G0CuZtOuBiPe94fTsLXXRaFQI+BJBCEF9fT2K7eNipoykJl9+5hpAyykIx5WCeu4OJQoS8pTSOICrHN5bCWBlIfsfCPh49ra+/E0HcSWZIFKt98SklKYIMUq1Cn7ZmDeMfWaphW8/3Gv8bW18PBSX1kLAl45SXFvDJl/G5hog2eP33Dl2rsjhR3l+S1nS2hNFW1/MqA1zqDt/Ic+bVmqCWjzxw6/vSSk5fOcz27Dsnldy2me2mYWnTKs3/raGXZar1jVYNDU1Ye7cuTl95qGHHsLq1atz2teqVaswY8YMzJgxA6tWrcprrMUi6JXgkQhae8o3ugYAqB7+wUJGhzvlcRZ50NwZxmn/rQnbExtrkVBUowRAPvCmFWYLv+vZ7ZBcxGTG+PWaJgBaXH0mbd5w5mZZI2T19Usw/fZnAaQmcQR9QsgPNjfffHNO2x89ehR33nknNmzYAEIIFi1ahAsvvBC1tYNjRiCEoDboRWtvmQt5/Xe5lOooz28pC3Zwpg1WQfH37+w3OunkCm9aqeEcnpubu223b8/iOLkWguJb/LEkjmtOPgYhr1SWxaQKoampCbNnz8aNN96IOXPm4JxzzkEkok3ymzZtwtKlSzF//nxccskl6OzsBABs3LgRCxYswIIFC3D//fcb+1IUBd/4xjdw4oknYv78+fjFL35he8w77rgD9957b9p98Tz//PM4++yzUVdXh9raWpx99tl47rnninkZcoaF5RJSHolCdlBdlS+X3I/y/JaygMWgA5pTslP3pP9r2xFcsWRyzvvjTSsBTzJ0kbVLs9IZjmesK57Io0aIi2jhX0yT/8FFc3HHJ+dk/fnB4M5nPsC2gz1F3edx46vw/QznvXPnTvz+97/HI488gssvvxxPPfUUrrrqKlxzzTW47777cMYZZ+B73/se7rzzTvz0pz/Fddddh5///OdYtmwZvvGNbxj7efTRR1FdXY3169cjFovh1FNPxTnnnJM2ssVpXzwtLS2YNCkZiTxx4kS0tAxuJDIT8tUBT85hs8MFZq4pl1IdI1a94zvPR+IKvnLWTABAR39qfHtnfxyHulNNOf0xGd/56xa09kRNphW+yxDT8K0ctTmOlXgedbtZr1XeJl8OCR2lYMqUKVi4cCEAYNGiRWhqakJ3dze6urpwxhlnAAA+97nP4fXXX0dXVxe6urqwbNkyAMDVV19t7OeFF17A6tWrsXDhQpx00kno6OjAzp07HY+bbl9DnbG6YlIu4YXpqBqgbmGlpjymqjyIc5q85CL40lkz8Ms39hglRnk++8g6bD/cmxK2+MbONvx23X4oKsVZs8cA0ARyDdcUutNBmHf2Z35IEjna5AHgrNmj8cT6A8NKC8mkcZcKny+ZEyFJkmGuyRVKKe677z4sX24utnr77bfjH//4BwDNBJQrEyZMwKuvvmr839zcjI997GN5jbFYTKwtvJrlUKcm6EFHfxyhMvFjjVhNngn506Y34Nvnaw2ya0NedIZThTILTeS1fyCp9fdEZJP9vJ4rJ+CksR/tz8Emn0MvzW+fPxu/uvZETKwNZv0ZQZLq6mrU1tbijTfeAAD85je/wRlnnIGamhrU1NTgzTffBAA8/vjjxmeWL1+OBx98EImENnF/9NFH6O/vx8qVK7Fp06YUAZ9uXzzLly/HCy+8gM7OTnR2duKFF15ImUgGmlP1Mgk/uGhomwAL4YErF+HzZ07HqDwTI4caw0fdKzLMJn//lScY5pW6kDetGaW1J4ZJdUnhua8jDEAzhzDTitftMjlArZOGV3IhrqiGDyAd8Tw68FT5PTjz2NFZby9IZdWqVbj55psRDocxdepU/OpXvwIA/OpXv8L1118PQgjOOeccY/sbbrgBTU1NOOGEE0ApxahRo/DXv/417TGc9sVTV1eH7373uzjxRK3u3/e+9z3U1dUV6Szz46Sp9Xj+y8swc0zFoI6jlMwaW4lZY2cN9jCKxogV8kyA8gkd9SFv2lj5o/1xk5B/+PU9ALTGI9Zwx+e/vAw/+dcOvLDtiBEuSSmFrGrb2a0YUsYo5y7kBdnR2NiIrVu3Gv9//etfN/5euHAh1q1bl/KZRYsW4f333zf+v+cerUeOy+XCj370I/zoRz9Ke8w77rgj476sXH/99bj++uvTn8wAk29WuGBwGLHSI5ZItXfbmWtYOBWQ3vTSpxcLY3a8WWMrccq0BlCaFOhxRQULtnGy1fOwejiFNGEWCAQjmxErPeKKCreLmCJP6kNedPTHTYKdb0zNC3mtZIH2d2d/wugwxbfBq9Vt80ygR+NJZ2825pp8QigFAoGAZ8RKj7isptTeqA15EZdVk2Dn2wPyQj6hUCOetqM/ht6oDK/bZSofwBywzEHLOj0BQFcW5pqEosJFci/jKxAIBIwRLeStGjKrv84Lc5MmzwlmFmnTUOFFNKGitTdm1Ntm1Fr2xwv5o9nY5BVV2OMFAkFBjFgJYivkQ3ZCnmsq0s8Lec2UwhokN3X0p6RB11dYhLw+YUyoCaArqzh5KuzxAoGgIEasBIkrqUKe2dB5LTuasLfJMyE/rlpLDtnXEU5JQHLS5MdV+9Ebk00JWfZjVMqmO41AIBgcRqwEiclKipZcb3GUAsC7+7oAaJ3p+cibmC6wx+ua/NH+eIqQ97pdqPS5DSHPJozxNdrE0BVJb7JJyDTrMsOCgWWgSg1/85vfxJw5czB79mx88YtfNIICNm7ciHnz5mH69Omm1wUCKyNWyGuOV3Pacq2NuYZVizxj5ihTRIyhydck07wrfam1LuoqvCnmmnE12sSQqbRBwma1IRi+3Hzzzbjmmmuy3n7NmjV46623sHnzZmzduhXr16/Ha6+9BgC45ZZb8Mgjj2Dnzp3YuXPnoFenFAxdRqwEidnY5Kv8brhdJMUsU+l3o6HCZ9Lw44a5JllJcoJNXY/aoDfFXDNeN/FkSogSjtfSMRxKDRNCEI1GEY/HEYvFkEgkMGbMGBw6dAg9PT1YunQpCCG45pprMmbYCkYuIzbj1U7IE0JQq5c2UFWKfUfD6IkmUOX3oE5PlFJVCpeLGJp8Q4UPFT43+mIyZo1JzQSsD3nR0qUJj4jFXJMpISqhqPC4RoCQf/ZbwOEtxd3n2HnAirvTbjLUSw2ffPLJOPPMMzFu3DhQSvH5z38es2fPxoYNGzBx4kRju6FQglgwdBkBEsSe3qicEvIIaGGUR/vjePyd/Tjz3lexdncHqgMe1Aa9UCnQo2e2shBKn9uF5XovyMn1qUXBJtUFsf9oGJRSwybPtP9MCVExWYXPM2K/opIz1EsN79q1Cx9++CGam5vR0tKCl19+2SicJhBky4jV5HsiCcy2qcGhUIoXth2B36PZ6w91R3FMfdAUXlkT9BplEXxuCd+9YDYWTKrGksbU4lHTRlcgHFdwuCdq2OQNTT6DuSaWUOEfCb1ZM2jcpWKolxr+y1/+gqVLl6KiQisGtmLFCqxduxZXX301mpubje2am5sxYcKEvMYuKH9GrJrYE03YNgVgwpwXwBU+T7JEgf46M9d43S7UBL245uRG2+YcLPrmcHfUMNfUBDwIeKSM5pqorAhNfoAZSqWGJ0+ejNdeew2yLCORSOC1117D7NmzMW7cOFRVVWHdunWglGL16tW46KKLin4tBOXBiNTkFZVq5hobIf+Fj0/H1Y++g60tyd6sFT6Jy4Y1m2syRb+MrtSEfFtvDJGEAp/bBZeL6Db+DOaahIr6kBDyA81QKTV82WWX4eWXX8a8efNACMG5556LT37ykwCABx54ANdeey0ikQhWrFiBFStWFOnsBeXGiBTyfVEti9XOJr/4GM3kwgvgkM+N2pA2ITDtm2ny/gya9ugqzSRw0282wiu5ENSrVNYEPY7mml2tvaj0exCTlZQwT0FxGA6lhiVJcozUWbx4sWn8AoETI1LIM+dptY0mH/BKRjVKht8jJW3yzFyTYI7X9EK4PuQFIVpz4Liiol5v8l3n0IUKAM76yesIeiVUBzxGQ26BQCDIh7KzBVBKce2v3sFf37MPKfvnlkP487vae06NekdX+U3/bz/cg4BHgs/tStHkrZUsrbgll2mbgIdp8l5bmzyLwAnHFRzqjmJXa1/a/QsEAkE6yk6T74nKeHVHG17d0YaLj0+NOPjPx981/rbT5AFgTJUPHx4C5k+sxubmbpw2fRQIIab2gNkKeQAgSDpkFT39vC7osbXJswxbBl87RyAQCHKl7IT8oW7nMDhrfY8qv72QH6tr8lMaQlh13RJjMuCzV+OyCslFTP1cnSBc0A37fE3Qi+5IArKimvbBlzYGgMeuPTHj/ocrWuMVUZunFIhaNgJG2Zlrmo8mhbz1Ro9Zqj5WBeznOGauCXgk1Ia8RmhkQ6UP7X0xPL2pBT9/ZReyFU/8dr2605fZ+LsjZm0+YhHyfNmEcsLv96Ojo0MIoxJAKUVHRwf8/vK8dwS5UXaa/A2rNxh/d4YThjAFkn1dGenMNYC5yQcAjK3yYcfhHtz1z+0AAFnNTkAFfW70W4R3Ldc1qr4imZRjPWa5aroTJ05Ec3Mz2traBnsoZYnf7zeVPhCMXMpOyPMc7o6ahbxsFqAhr/3ps5h4q+lkTJUfbb0xo058tqy6bgmeercZj76513htgl6JsqUrgplczRurkC9XPB5P2touAoGgOJSduSbgkdCgd2Q63BPBGzvbcNH9b6E/JiNq0eTtMlSBpLlmQo25quSYKj9U6vw5J44bX4XvXnAc6kNeXHdqIwBgYq1W5+b+l3eZto3GR4aQFwgEA0NZafKyoiKSUHDhgvH4w4YDeG1HG1at3QcA2NrSbZhIMnHC5Bo8eOUJ+Pjs0abXmUO2rTdm97GMbPzu2cbfo3QTzYHOsGkbpslfMH8czj5uTF7HEQgEAkZZafJduhNz1thKuAgMAQ8AO1v7sg5HJIRgxbxxKYlOYzknaH3Ii7/eemreY3W5CL7w8elo642ZxsWE/HcvOA4XLRRFpwQCQWGUlZBn4YmjKn0YVekzvbenrd8w15w+owG3nzc75/2P4ZKkls8di4WTagoYrTYZqRSmhCc2RlYFUyAQCAqhrIT8nrZ+AFp44thqsz39YFfE0Ji/9IkZuHHZ1Jz3Xx/ywq3b40NFKDcwbZRWQrapo994jY0xIIS8QCAoAmVlk7/5txsBaElL46r8eB/AadMbQAhwsDsp5PPVkl0uglGVPhzqjiLgEJmTC6xd4MGuZGx/JK7ARSAaeAsEgqJQVpo8Y0JtwLCf14W8mFgb0DT5HEoROMFi64NF0OSr/B5U+txo6eSEfEJBwCOVbXy8QCAYWAoS8oSQPxBCNuk/TYSQTdx7txFCdhFCdhBClqfbT7Go8Llx/alTUB3wGEJepRTjqwNo74sb2aWF2LtZKYRimGsAbUJq6Yoa/0cSiqg8KRAIikZBNgdK6WfY34SQ/wHQrf99HIArAMwBMB7Ai4SQmZTSkgWBJxQVfTEZNUFNCLNyAL1R2Wi319Su2b4L6bZUqdegL4a5BtBaAbZw5ppoXBFOV4FAUDSKYq4hmm3hcgC/11+6CMATlNIYpXQvgF0AlhTjWHYkYmF8uH0b3EgKeRYJ0xNNGLbvPW1aFEshQrTSptFIIUyoCZht8gkh5AUCQfEolk3+dABHKKWsRf0EAAe495v110rCe/96HPP/eCqOIUdQ4dOEMCtn4JVcRubqbj36ppDm2CF9//0xuZAhG4yvCaA7kkCfvr9wXCmKvV8gEAiALMw1hJAXAYy1eet2SunT+t+fRVKLzwlCyE0AbgK0xsX5UFdTCwAIIoagbkaZMboC31g+C5eeMAH1IS1mfv9RLbu0kMiV/1g2DZubu3H+/HF574OHj7CZOabScLwKBAJBMcgo5CmlZ6V7nxDiBnApgEXcyy0AJnH/T9Rfs9v/wwAeBoDFixfnVXe2oU7ryxoiUYT0HqqEENx65nSnMedzGADA5PognvnCaXl/3opRqKxTF/Jxxai9IxAIBIVSDHPNWQC2U0qbudf+BuAKQoiPEDIFwAwA7xThWLZUV1UDAAKcJj9cmFCjFSpr1mvYhOPysDsHgUAwdCmGkL8CFlMNpfQDAE8C2AbgOQC3ljKyhnhDAIAQzKWFeR648gQAwOxxVaUaRl6MrvTB63YZpqRoQhWOV4FAUDQKVhkppdc6vL4SwMpC958VupD//GnjMKUhZLvJefPGoenu8wdkOLngchFMrgtiXwevyQshLxAIikN5ZLzqQv7YuuEpHBvrQ5yQF9E1AoGgeJSVkEeiP/12Q5TG+iD2He2HolLEZGGuEQgExaM8hLzkBYgExMOZtx2CNDaEEE2oRslhFiEkEAgEhVIeQp4QTZuPD09NfuoobSXy0vYjAIDJdcHBHI5AICgjykPIA5qQH6bmmumjtbryr2xvBQAcU2/vPBYIBIJcKR8h7wkOW01+VIUP1QEP1jd1AgAaKnwZPiEQCATZUT5C3hsatjZ5Qghmja00/q/Vi6wJBAJBoZSPkPdVArHewR5F3hynJ2lV+d1wS+XztQgEgsGlfKSJrwqIdQ/2KPLmWF2T74+XLDFYIBCMQMpHyPurgejwFfKs3IKi5lWjTSAQCGwpMyHfM9ijyBveJi8QCATFonzKHfqrgFgPQKkWNz/M8HskXHnSZMwcI4S9QCAoHmUk5KsBqgLxPs0JOwxZecm8wR6CQCAoM8rHXOPTSwgPY7u8QCAQFJvyEfJ+rXGIEPICgUCQpIyEPNPkh6/zVSAQCIpNGQl5ockLBAKBlTIS8jXa75jQ5AUCgYBRPkJeOF4FAoEghfIR8oZNvmtwxyEQCARDiPIR8m4f4PYLx6tAIBBwlI+QB/TSBkKTFwgEAkZ5CflgA9DfMdijEAgEgiFDeQn5UAPQ3zbYoxAIBIIhQ5kJ+VFCyAsEAgFHGQr59sEehUAgEAwZykvIV4wC4r1AIjLYIxEIBIIhQXkJ+dAo7bfQ5gUCgQBA2Qr51sEdh0AgEAwRykvIV4zRfvceHtxxCAQCwRChvIR8zWTtd9eBwR2HQCAQDBHKS8gH6wF3AOjaP9gjEQgEgiFBeQl5QoCaSUC3EPICgUAAlJuQB4DqScJcIxAIBDrlJ+RrJgHdQsg7oshAzyHza6oKRDoHZzyCwWH9o8Cz3xL9F0YA5SfkqycB4Q4g3j/YIxmavP0g8JNjgY7dydfe/TXw343A4S2DNSrBQPOPr2r3wlv/N9gjEZSY8hPyNcdov4XJxp6tT2m/1z2QfG33y9rv5vUDPx7BwBPrS/797qrBG4dgQCg/IV8/Tfvd/tHgjmMooqrA0b3a33vfSL5O9NugdfvAj0kw8LBkwTFztYJ+4aODOx5BSSk/Id8wU/vdJgRWCm3btaYq9TOA9h3J2vt9euXOva8P3tgEA0efLuQnLdF+d+4dvLEISk5BQp4QspAQso4QsokQsoEQskR/nRBC/o8QsosQspkQckJxhpsFvgrNZHPkgwE75LDhwDrt92lfSf5/6H1g/xrt/7YPgaY3B2dsgoGDCfnx+mPZc3DwxiIoOYVq8vcAuJNSuhDA9/T/AWAFgBn6z00AHizwOLkxbgFw8N0BPeSw4PAWwFcNzL0UkLzA/rXAL5Zp7y26DgAB9q0Z1CEKBgBmrhl/vPa7u2XwxiIoOYUKeQqgSv+7GgBTCS4CsJpqrANQQwgZV+CxsueYU7Ss1+7mATvkkKe7RXO6jp4NeALAhEXAOm7unbIMqJuiTQTb/ga899vBG6vADKXA2ge0VVcxYGa6hpmArwro2Fmc/Y40VAXo3DfYo8hIoUL+ywB+TAg5AOBeALfpr08AwIe3NOuvpUAIuUk39WxoaytSV6djTtF+71tbnP2VA0/doMVE1zZq/886D1Bl7e9b1mra/Zi5QPMG4MmrgadvFWGoQ4VD7wPP3wY88vHi7C9yVBPubq826R/ZVpz9jjSe+RLws/lA+9CeJDMKeULIi4SQrTY/FwG4BcBXKKWTAHwFwKO5DoBS+jCldDGldPGoUaNyPwM7xszVbuJ9bxVnf+UAs7vPOFv7fdyF2m9PCBhznPb32HlAL2efPfDOwI0vEz0HgQ//Djx2LvDoOUC0Z7BHNHCwEFdVBrb/s/D9hY8CgVrt79GzgdYPtNXCSGL/OmDt/UCsN7/Pqwrw3m+0vw++l98+lER+n8uRjEKeUnoWpXSuzc/TAD4H4M/6pn8EoLvr0QJgErebifprA4NLAhpPA3a9OLJu3q4DZicapdpqhkXPnPIFYN5l2t+1jcAFPwWu44TGuAXm/f3mYuDvXy3pkLPmhe8Cf7hS8yMceFt7SAeSw1uAjb8e+Pupu0U7ricIeCuBP14L9B4pbJ99h5O9Fxpmaiu8gcx4bt4APHS6NlkPRj5Lz0HgseXA898GXrsn+Xp3MyDHs9tHK7f66diV+xi2PgWsHDsg1oZCzTUHAZyh//1xAGzd8jcA1+hRNksB/P/2zjxMiuLu459id4HlWnZF7nNBTuUKLHIpigJiBFGiBK9XRUwEb00UL3wTjSeJryYaFYkajBgB8eAQFVAJCijXrrCwXIpBkNuLu94/vt12z+7sBbvs7FCf55lnqrtneqq6qn/9u6pmt7V2c7QTlBptBml5g2MtDMqSlwbBuDaQPUPb6+bChAHwWAtttxwQ+fkuV0L9jsF2/ShJUIvHK7++MPb8FxY9n/9NsmE+rHpH5f88KY28qALT2iC9s/tovUe7sdZ/WHqzdl+/SuZ51pTCP1uSvPEb2LVRGWMjZsOhffDPC45uOYJvs+HEVir7kwd3bjjqqhaZ6bfDN8v1sJ7zwLH7XVA7X7ko2PaVou1r4YkOkccO7st/jIblypGsevv5y7LMjkHa8tEK+WuAx40xy4AHUSYNwHRgHZADPAdcd5S/U3zanCfNZ9Fzx/ynjzk7N+jG2bFO2/8aBh/8MdJdldpU1k1BVAu5y/qFbr7vi/AnLO/eDe/cCivfjF6/fwyEV4fLf/nu3dLIw0srFMSmxcoIOXcc9PsjVEgMMkQOH4I17ykA9uJ58Ewv3biHDx29xuuzZ3Mwue6LaSVzzqKw/4dACLQdJNdKj+thSybkvH9k5/xxB3y/BU5sre2f/4OhFAKIhw9Hzq4F9cuWLGj9SzjlIljxetHHQUFkToHnz5b7JVqQ+tNnYWyKBPk3K7QsedPewTpXX34iobtuTqBdP9sHJl2q8v4f5F7xJ45t/A9Urw+NexQv+JrzPjzaQr8DpXPdc3FUQt5a+7G19hfW2g7W2m7W2s+8/dZaO8pa29xae4q1dnHJVLcYVKoGXa9S55dUVkKsMmUkLHxW5f5/0vuHj+qVUAlqtYLTf1+0c7UfJo2+x2i4bKr2bV+rAf7CAA363L7EsKade6bx7Ht1Y/m8dH5Q3uT5/AszkZe/qv8JaHu+lpNOTg3cC4vGw8QLFQDzGdcG/jcNHm+p4xDUP2tq4dcgN76WWylFwqA4fxSfOQXeux/mPAh/PgVm3wdr5xTNivEDokOehdPvULnvWMVRcqe6HtwHy1/Tv6IVZHn5/v2GXfSe6mvynrBZ8FcJw5l3Ft9n/N8lQQ7++o/gzevhTw10nh+2K2tr7QeyRloOgL73wOED8GRnLZx3NHzytMbTQ40VpJ42OjhmLXw8LtiuVkdjO7VJ0Le+Nl65Jrx3nx4WW7+AVW9LMXm8NfyhFjySLkt5zbvQpLvOURxN/vOXNMvYx//9UkxySCy1M8cCvW6Gpa/AlGvh6llQOaWsa1RybFqsG7teB5m9PhnXSEvJnAxrZqk8uhgB1Av+HpTTvCUitq+BlW9J+wYJgIGPSuCCtBF/4H4X8srNfQjmP6Hyr1+VhbFnE9RoIA3vq4USmO/cAhf/U9ZXNDb+RzdU1RO0nZwWaFRhy6Fpb2h8qh5uPu/cAhWrQmoz1X/bamg3pOjXw28fQLvztdbLA3Xh2g/zxjDC/LgDJg6Frz+L3D//L3r1GQN9ojx4Fz4nP3n66cFM1HodoIKnjyUkQuNueYX8vEfgo8eC7UunQIu+ec+/eIIW8WvUTduVUyTYtudIMM8ao/2f/E2Crv8DMOP3SrGtUguan6k6hNmxXpZkzmyoXg+GvQIv/jJ0/b6EeQ/D8kneb9aE5mdASkNleWVPlyLW8BfBdw78BDPvgJwPIK0pnHU/NChgTuWBH/VuvQfcN8tlMdZuq3H83WbFoNoNUZuNkXX7/RbY/6P6uEYDOOMumHadfttnzgOwzw/0W41jgB43qO57vpaiklgx//qBHmTr5kr7b9JDlveXC+DVS3SeUYugVouCz3EExLeQT06FC5+Hf14IEwbC+U9DvfaFfy+W2bFeAVFfA6hYLTiW0hgSkqDDxUqJfPtmaNz9yH8rpaHWtXnbmyGbVFXnXvQcHNoP5z2hmyWcQvbTLr1vXg5zPati0FPQ6hzoMFw3/MBHdANuWhS0I+f96EL+xx3SqE6+INhXJU0PlY//DBs+glOvg27XKpiYUAnS0qXVnnodvPIrWD0L1szOe72Kiq/ldrlKWtxPO+HlC+DW7LwCzyd7RiDgOwzXA6r9xbIsNs6XQO5xPVSsEnznwF6YfpvKY/7rrTNkgrRXnyY95I7zfevLJkUKeFCwtkVfCcu9e6B6HWnZGz/WA6ZCQvDZRhnKFFnysgT5FW/B0omw4CkJxGX/0gtg8N+g0yWRv7XkZQl4kDCd93Dk8R3rYPXMYLvHaI0tkODNnq5+PLEV2EP6zZVvqw2gPwGaOBRuWQkH9ypzzlcwQJr6zo2yQttfpDpMG6XYT5hmp0FyzWC7pnddMyerjjWbQMfhsmbC80SypsqSbNRVE8i+/lwKRf2Ocp1h5fbx183Kj81LtaxIxgg4+UIpQVlRWcjSAAAPdUlEQVRTZC0AZL8DtW4s+BxHQHwLeYD0PjD8NZhyDfy9NzQ6FZr1lrZUrY73VK8A2JAJ7b1bGyrn3k8RvpO7XNzvkPc7q96ODJLt/15C5Kz7ILFSsD8hCQY/VdjVKZgKCTKrs6dDx0vgnEeU5VEhUe6hriPk35w2Sp8/oUXgRlnr+Yxvywl8/UNCk68aZsC8h4Lt/AKmfoCrcY9gX3Kabgjfsmg3JFIQdhyuF0DzvpHB0h+26VqGhURh7NooDbV+R7h9DXz2Irx1g6wS/3etlVZYva62t62GCklwx0ZZEj7n/UUxhIkXKq21xVnBsS2ZQXnVdAmeGvUhqXJkfdpdoCWCJwyEm5bD+/drPaKrZuq3/vHLYIxMvVZxhKETAjdTq3Miz3fWWLkfAHrfopTariMk5Bc8JWWp/TAtTbwjiv/cjwVVrycBu3omVK0NVWvpAb1urgLF547TA+XENsF3q9fRvbjwOQm9gz9Bw656wFWvrzH843bdv//+Hz08z3lYD3WfH7fD/u/UP7710vpcuYa++0aunORUPfzD1O8EJgHe9Fw77YdpXNRpp7RSkwAt+2v8974VTr89b9vrtNP76pkac35qajQ2eIsCNvVmmdcLJT20OjcIgpcw8S/kQR1//Wew+AUN+I8eD8y68kjDDLjiTbkN8jP7S4pBT8HmJZHCqPsoCfmZd0oz9DnhpGCW8bYcqFY3Mpgb0YauQbl+JwXjDh+K1DBBN1hSVc3Q9anhT5420u5qFDCZuk674IHTZwzMfVCCqEb9ApsNKCNo3VxptQkhU7zWSV4b10jIHzoIr18p11GlFLhhiY6lpUcKeJ8mPXS+1bMir2vY7bboeS0oV+fkvN8/oTn0vVdaf+ZkuQuGTpBQBbliFo+XC2GtF+DLng4YWTt1T8l7jUbOlcXmu7JSm8qts/sraHY6nPOQtPVvs3V85dtaOqTnTRLIzc+ESyYrNrL7K+h8GZx5DzzeSg8KkDDP/dsAGSNhxu/kWswYofWTEivBoCc9a2QvNBrvtQG5+cJC3o8DhYVkcqq0ZZBFd/hQ3gd7rRZwc6ZcVFlTgwe0379Nesj1tG21FJho1Ouods0ao9eF44M0ZdC4NhVkaSyeoIC3f0+07A+9bpFidN5fop+/BDg+hDyo03vfqteBvYEf+addBCqzCQ2EXGXwtqOVc3/maM5VhO+nNdPSBPdszysUS5qqJ0QKIghuJl/AN+kpoWBMEOTevqZg87XxqXpYZYxUVsMbv5Ff2E/tA7kZvpgGbQdHarMNukgIpjQqWMCDbvQvF8AZY/j5GmdNhW6/KfjaWauMIJ/ktKDsZ6ZsydLksuWTgtjAvt3w4SOyNFqH/NJhKlaRsAuvlLpqOrw3Vtf25AvkigJZndHwNchFz0vj9Ce5gTTag3thxWuBL3nFv2VZtB0U3Yqp30mC3scYXbv5T0AHzwdd9xRYN0+ur0mXBNdi15f6foUK0ta3rZYQNgZ63wYzPA04rVn0tmRcIy04MTmIPYRJqqyY2oGf4OUheRdU8/3nqflowsbk71arUV/WfdZUSPTGmK/xNztN3w2PyWjnvmqWrIa3bpRbK3u6rN7KKTC+v6wMkBXc/8HI7551X/7nLiGOHyEfJqmyOq6gzisP5DdwSxtjNFhnjVGaqj+h6t17ZDpbK4Gdn5ADZT+N8Py4/oNh6xfSNL9ZrkDvskm6QTJGRH639UC5b7qPKryu9TvCiPdU9mc3zhoj98EZYyI/uyVLwnHLF5EPqJTGQaYRKCaQli4tv+NwuZ1qtYLrPoFH0+HTZ/S5/LQ/kPa3bq6yTiqnwBu/VYrkkGehdmt99/OXpIVGwxfym5dJc69UPTjWpKfe53gxkb73yaVTp62smaJy1ljoc2fwgO06QsJwYkhT3bZGyyT4qZgt++nl022kYhC7v5I7Jz+iWTy5SUrWdfH9+9kz1Kf++Kndtqgti6TDxYqfdL5c220G6WHV4ddF+36VNGnvC/6qCZggKyCtucZvo25w2u2yiAoLzpYCx6eQdxw9HS9ReuDpvwv2VTlB6XFffSph76/tXxi+MFwzWybt+nnBsdRmka4akFC8akbx61ypOgx8TG6O3JPkdm+Cp3vk/U56HwUbU3ItvVT3FFkZj50kDfTyadJCm/QMAmkFPYRqt5HW92i63BJ7d0HPcRLwAJ0u1Ss/KqfIKt26Cs6+P/JYSgNZBLs2SvD1uln9Va128WIRxkRaUE17SVCtnwen/U7Whh9L8YV8NC4qwX+fqt1Ggd4Ff1OqpM+lk4vXtjDJqXBhaD5NYkXoenXxz9PlSliSpPH/zq1SgJLTpCAU5SFWSsTfn4Y4jg3JNeGa9+VX9PF9wlOvlcD3g5+FUbGq/PdLJ0YK+PqdpQGVJBnXKEtm/bzIJXbfvTvvZ0/qL+GdW8CDzuEz5GmlNQL0vFGuoN9vkFDNj7aD5W4CmOX9dtN8XDP50fde+PUrgQ85TPuLoUZDGPqChF/1OkcuBMMMehIuewPOvEsPmqII+ZKk02V6Dwt4iAzMlxWdL4er34ULntf2/u/k8ipDAQ9OyDtKkgbeBJudGyToqqQV+PEI+oTykjsMh1tXw8g5edP1SgI/IOeb/d+ulhvitNth7G7lSkNk2mZu0vvA3VvhtjWRefeNMpT9UVCWBSiGcXOmBPG+3bIMCnooFJcz74JbsqT5liSpTZTjDurfPV6g/VgJ+co1ArfPKaElCMKpqGVN+195/89A6YzfYuLcNY6S48RWcq/sXC9Nsjh0uVLuiVVvy5efkFQ6dQQF2kwFpdeBJjhVSIIML2Oj+2hNvGkeZTJRmMRKRy+Ye96grC9/pnJ5omptBVkTKqp8rKhUXRlSXUfIYmw1oPDvHGsGPqqYRjgvv4xwQt5Rchgj7XvnhvwzHQoiIan4s1GPhIREBXi//0ZB4szJcFK/ILWtYpW8GUWlRbdrI9MByxP+9UppFD0rprQ4+w9y7TX4ReAmizUSkmJCwIMT8o6SJjm1cFdFLFCtjjT5zcukFbYpIBPIER1fe88vNbK0aDUgNrX3GMX55B3HJ9XrSsj7E6WOleYeT/gxl7RCpvM7yhQn5B3HJ6nNlOOdPRNqtyvZoOfxwsF9evdnijpiEifkHccnLftpnZRNCzWz0VF8Ol+mBfD8GbGOmMT55B3HJ+HVOcO5/o6ik5auRdEcMY0T8o7jk4pVNcV/zyblvDsccYoT8o7jl9JcvdPhiBGcT97hcDjiGCfkHQ6HI45xQt7hcDjiGCfkHQ6HI45xQt7hcDjiGCfkHQ6HI45xQt7hcDjiGCfkHQ6HI44x1tqyrsPPGGO+BTYe4ddrAdtKsDpliWtLbBIvbYmXdoBri08Ta+2J0Q7ElJA/Gowxi621Xcq6HiWBa0tsEi9tiZd2gGtLUXDuGofD4YhjnJB3OByOOCaehPyzZV2BEsS1JTaJl7bESzvAtaVQ4sYn73A4HI68xJMm73A4HI5cOCHvcDgccUxcCHljzABjTLYxJscYc0dZ16coGGM2GGNWGGOWGmMWe/vSjDGzjTFrvPdUb78xxvyf177lxpjOZVjvF4wxW40xmaF9xa63MeYK7/NrjDFXxFBbxhpjvvb6ZakxZmDo2J1eW7KNMf1D+8t8/BljGhlj5hhjvjDGZBljbvT2l6u+KaAd5a5fjDGVjTELjTHLvLbc7+1vZoz51KvXJGNMRW9/JW87xzvetLA2Fglrbbl+AQnAWiAdqAgsA9qWdb2KUO8NQK1c+x4B7vDKdwAPe+WBwAzAAKcCn5ZhvU8DOgOZR1pvIA1Y572neuXUGGnLWOC2KJ9t642tSkAzb8wlxMr4A+oBnb1ydWC1V+dy1TcFtKPc9Yt3bat55STgU+9avwYM8/Y/A/zWK18HPOOVhwGTCmpjUesRD5p8BpBjrV1nrd0PvAoMLuM6HSmDgRe98ovA+aH9L1nxCVDTGFOvLCporf0Q2JFrd3Hr3R+Yba3dYa3dCcwGBpR+7SPJpy35MRh41Vq7z1q7HshBYy8mxp+1drO19nOv/B2wEmhAOeubAtqRHzHbL961/d7bTPJeFjgTeN3bn7tP/L56HehrjDHk38YiEQ9CvgHwVWh7EwUPiljBAu8aYz4zxoz09tWx1m72yt8AdbxyrLexuPWO9faM9lwYL/juDcpRWzwzvxPSHMtt3+RqB5TDfjHGJBhjlgJb0QNzLbDLWnswSr1+rrN3fDdwAkfZlngQ8uWVXtbazsA5wChjzGnhg1Z2WrnLby2v9Q7xNNAc6AhsBh4v2+oUD2NMNWAycJO1dk/4WHnqmyjtKJf9Yq09ZK3tCDRE2nfrY12HeBDyXwONQtsNvX0xjbX2a+99KzAVDYAtvhvGe9/qfTzW21jcesdse6y1W7wb8zDwHIFZHPNtMcYkIcE40Vo7xdtd7vomWjvKc78AWGt3AXOA7sg1lhilXj/X2TueAmznKNsSD0J+EXCSF7GuiAIWb5ZxnQrEGFPVGFPdLwP9gExUbz+b4Qpgmld+E7jcy4g4FdgdMsFjgeLWexbQzxiT6pnd/bx9ZU6uWMcQ1C+gtgzzMiCaAScBC4mR8ef5bscDK62140KHylXf5NeO8tgvxpgTjTE1vXIycDaKMcwBhnofy90nfl8NBT7wrK/82lg0jmW0ubReKFNgNfJ33VXW9SlCfdNRtHwZkOXXGfnf3gfWAO8BaTaI0v/Va98KoEsZ1v1fyFw+gHyDVx9JvYGrUAApB7gyhtryslfX5d7NVS/0+bu8tmQD58TS+AN6IVfMcmCp9xpY3vqmgHaUu34B2gNLvDpnAvd6+9ORkM4B/g1U8vZX9rZzvOPphbWxKC+3rIHD4XDEMfHgrnE4HA5HPjgh73A4HHGME/IOh8MRxzgh73A4HHGME/IOh8MRxzgh73A4HHGME/IOh8MRxzgh73CEMMY0NcasMsb8wxiz2hgz0RhzljFmvre+eoYx5vTQuuZL/NnLDkcs4iZDORwhvJUPc9Dqh1loevwyNBt2EHAlWqv8IWvtfG8hrb02WFXQ4YgpnCbvcORlvbV2hdViWFnA+1ba0AqgKTAfGGeMuQGo6QS8I5ZxQt7hyMu+UPlwaPswkGitfQgYASQD840xx3z5WIejqCQW/hGHwxHGGNPcWrsCWGGM6YrWCF9VxtVyOKLiNHmHo/jcZIzJNMYsRytYzijrCjkc+eECrw6HwxHHOE3e4XA44hgn5B0OhyOOcULe4XA44hgn5B0OhyOOcULe4XA44hgn5B0OhyOOcULe4XA44pj/B8nd5ICFwkKhAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "from bmtk.analyzer.cell_vars import plot_report\n",
    "\n",
    "plot_report(config_file='sim_ch05/simulation_config.json', node_ids=[0, 80])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Perturbation simulations\n",
    "\n",
    "A common use case is to mimick a network where the firing rate of select cells are being depressed or stimulated using current clamps or optogenetic methods. The best way to do this is by applying a current clamp in a similar manner used in [bionet](http://localhost:32000/notebooks/Workspace/bmtk/docs/tutorial/01_single_cell_clamped.ipynb).\n",
    "\n",
    "The point_120cells network in the bmtk examples consists of 120 recurrently connected excitatory and inhibitory cells with randomized inputs. Normally the results without any type of perturbations will look as follows:\n",
    "<img src=\"images/raster_120cells_orig.png\">\n",
    "\n",
    "To artifically excite some of the cells we can open the config.json file and add the following section to the \"inputs\":\n",
    "```\n",
    "\"exc_perturbation\": {\n",
    "   \"input_type\": \"current_clamp\",\n",
    "   \"module\": \"IClamp\",\n",
    "   \"node_set\": {\n",
    "     \"population\": \"cortex\",\n",
    "     \"node_ids\": [20, 21, 22, 23, 24, 25, 26, 27, 28, 29, 30, 31, 32, 33, 34, 35, 36, 37, 38, 39]\n",
    "   },\n",
    "   \"amp\": 230.0,\n",
    "   \"delay\": 1.0,\n",
    "   \"duration\": 3000.0\n",
    "},\n",
    "```\n",
    "\n",
    "This will apply a large step current for the duration of the simulation for cells 20 through 39. In this case we knew the exact node_ids of the cells we want to excite. However it is also possible to get a subset using a specific property, for example if we wanted to only select inhibitory L2/3 cells:\n",
    "```json\n",
    "   \"node_set\": {\n",
    "     \"population\": \"cortex\",\n",
    "     \"ei\": \"i\",\n",
    "     \"location\": \"L2/3\"\n",
    "   }\n",
    "```\n",
    "\n",
    "To simulate depression of cells 40 to 49 we will add another current_clamp but with a large negative current:\n",
    "```json\n",
    "\"inh_perturbation\": {\n",
    "  \"input_type\": \"current_clamp\",\n",
    "  \"module\": \"IClamp\",\n",
    "  \"node_set\": {\n",
    "    \"population\": \"cortex\",\n",
    "    \"node_ids\": [40, 41, 42, 43, 44, 45, 46, 47, 48, 49]\n",
    "  },\n",
    "  \"amp\": -230.0,\n",
    "  \"delay\": 1.0,\n",
    "  \"duration\": 3000.0\n",
    "}\n",
    "```\n",
    "\n",
    "After we rerun the simulation:\n",
    "<img src=\"images/raster_120cells_perturbed.png\">\n",
    "\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.8.2"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}
