{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Registration across Macaque MRI templates"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Authors:** Nikoloz Sirmpilatze (nsirmpilatze@dpz.eu) & Chris Klink (c.klink@nin.knaw.nl)                \n",
    "\n",
    "**Last updated:** May 11, 2020     \n",
    "\n",
    "**Requirements:**    \n",
    "* _python_ >= 3.7\n",
    "* _nipype_ >= 1.2.0\n",
    "* _nilearn_ >= 0.5.2\n",
    "    * Used only for visualisation\n",
    "* _nibabel_ >= 2.3.3\n",
    "* _joblib_ >= 0.14.1\n",
    "    * Used only for parallel processing, not necessary when registrations are done serially\n",
    "* _ANTs_ >= 2.4.0\n",
    "    * _antsRegistration_,  _antsApplyTransforms_ and _antsAverageImages_ need to be in your path as executables    \n",
    "\n",
    "**Citation**: Sirmpilatze, Nikoloz and Klink, P. Christiaan (2020). RheMAP: Non-linear warps between common rhesus macaque brain templates (Version 1.2)[Data set]. Zenodo. https://doi.org/10.5281/zenodo.3786357     "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The following templates are used: \n",
    "1. [NMT v1.2](https://afni.nimh.nih.gov/pub/dist/doc/htmldoc/nonhuman/macaque_tempatl/template_nmtv1.html)\n",
    "2. [NMT v1.3](https://afni.nimh.nih.gov/pub/dist/doc/htmldoc/nonhuman/macaque_tempatl/template_nmtv1.html)\n",
    "3. [NMT v2.0](https://afni.nimh.nih.gov/pub/dist/doc/htmldoc/nonhuman/macaque_tempatl/template_nmtv2.html)\n",
    "4. [D99](https://afni.nimh.nih.gov/Macaque)\n",
    "5. [INIA19](https://www.nitrc.org/projects/inia19/https://www.nitrc.org/projects/inia19/)\n",
    "6. [MNI macaque](http://www.bic.mni.mcgill.ca/ServicesAtlases/Macaque)\n",
    "7. [Yerkes19](https://github.com/Washington-University/NHPPipelines)\n",
    "8. [ONPRC18](https://www.nitrc.org/projects/onprc18_atlas)     \n",
    "9. [F99](https://fsl.fmrib.ox.ac.uk/fsl/fslwiki/XTRACT)\n",
    "\n",
    "Within this notebook, they are abbreviated as *NMTv12*, *NMTv13*, *NMTv20_sym*, *NMTv20_asym*, *NMTv20_05mm_sym_brain*, and *NMTv20_05mm_asym*, *D99*, *INIA*, *MNI*, *YRK*, *ONPRC18*, and *F99*.     \n",
    "\n",
    "**NB!** We do not provide copies of the actual templates (licenses often forbids redistribution), but instead suggest you follow the links above and get them at the source. We do offer the warp files and warped templates that will be produced by this workflow. They can be downloaded from [Zenodo](https://doi.org/10.5281/zenodo.3786357) or [G-NODE GIN](https://gin.g-node.org/ChrisKlink/RheMAP). "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "If you want to emulate this code, you can consider setting up your templates in the following folder structure:    \n",
    "\n",
    "|--- RheMAP   \n",
    "  |--- notebooks     \n",
    "  |--- templates     \n",
    "    |--- D99     \n",
    "      |--- D99_atlas_1.2a.nii.gz     \n",
    "      |--- D99_atlas_1.2a_in_MNI.nii.gz     \n",
    "      |--- D99_template.nii.gz      \n",
    "    |--- F99     \n",
    "      |--- struct.nii.gz     \n",
    "      |--- struct_brain.nii.gz \n",
    "    |--- INIA   \n",
    "      |--- inia19-t1-brain_truncated.nii.gz      \n",
    "    |--- MNI     \n",
    "      |--- macaque_25_model-MNI_brain.nii.gz      \n",
    "    |--- NMT     \n",
    "      |--- NMT_v1.2     \n",
    "        |--- NMT_SS.nii.gz      \n",
    "      |--- NMT_v1.3     \n",
    "        |--- NMT_SS.nii.gz     \n",
    "      |--- NMT_v2.0     \n",
    "        |--- NMT_v2.0_asym_05mm_SS.nii.gz  \n",
    "        |--- NMT_v2.0_asym_SS.nii.gz  \n",
    "        |--- NMT_v2.0_sym_05mm_SS.nii.gz  \n",
    "        |--- NMT_v2.0_sym_SS.nii.gz  \n",
    "    |--- ONPRC18     \n",
    "      |--- ONPRC18_T1W.nii.gz  \n",
    "    |--- YRK      \n",
    "      |--- MacaqueYerkes19_T1w_0.5mm_brain.nii.gz       \n",
    "\n",
    "After downloading the warp files and warped templates from [Zenodo](https://zenodo.org/record/3776856#.XqqfI3UzZjE), we suggest you include them like this:      \n",
    "\n",
    "|--- RheMAP   \n",
    "  |--- notebooks     \n",
    "  |--- templates     \n",
    "  |--- warps       \n",
    "    |--- final     \n",
    "      |--- D99_to_INIA_CompositeWarp.nii.gz     \n",
    "      |--- D99_to_MNI_CompositeWarp.nii.gz      \n",
    "      |--- etc   \n",
    "    |--- linear     \n",
    "      |--- D99_to_INIA_affine_0GenericAffine.mat     \n",
    "      |--- D99_to_MNI_affine_0GenericAffine.mat     \n",
    "      |--- etc   \n",
    "    |--- nonlinear     \n",
    "      |--- D99_to_INIA_1InverseWarp.nii.gz     \n",
    "      |--- D99_to_INIA_1Warp.nii.gz     \n",
    "      |--- D99_to_MNI_1InverseWarp.nii.gz     \n",
    "      |--- D99_to_MNI_1Warp.nii.gz  \n",
    "      |--- etc   \n",
    "  |--- warped_templates       \n",
    "    |--- final     \n",
    "      |--- D99_in_INIA_composite.nii.gz     \n",
    "      |--- D99_in_MNI_composite.nii.gz       \n",
    "      |--- etc   \n",
    "    |--- linear     \n",
    "      |--- D99_in_INIA_linear.nii.gz    \n",
    "      |--- D99_in_MNI_linear.nii.gz          \n",
    "      |--- etc   \n",
    "    |--- nonlinear     \n",
    "      |--- D99_in_INIA_linear+SyN.nii.gz     \n",
    "      |--- D99_in_MNI_linear+SyN.nii.gz     \n",
    "      |--- etc  "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Based on the warp files and warped templates you could of course reconstruct the original templates with something like the following.     \n",
    "\n",
    "On the command line:     \n",
    "```bash\n",
    "antsApplyTransforms -i <TEMPLATE1_in_TEMPLATE2_composite.nii.gz> \\\n",
    "                    -r <TEMPLATE2_in_TEMPLATE1.nii.gz> \\        \n",
    "                    -o <RECONSTRUCTED_ORIGINAL_TEMPLATE1.nii.gz> \\\n",
    "                    -t [<TEMPLATE1_to_TEMPLATE2_CompositeWarp>,1] \\\n",
    "                    -n Linear \\\n",
    "                    -d 3\n",
    "```    \n",
    "\n",
    "In NiPype:     \n",
    "```python\n",
    "import nipype.interfaces.ants as ants    \n",
    "ants.ApplyTransforms(\n",
    "            input_image=<TEMPLATE1_in_TEMPLATE2_composite.nii.gz>,\n",
    "            reference_image=<TEMPLATE2_in_TEMPLATE1.nii.gz>,        \n",
    "            output_image=<RECONSTRUCTED_ORIGINAL_TEMPLATE1.nii.gz>,\n",
    "            transforms=<TEMPLATE1_to_TEMPLATE2_CompositeWarp>,\n",
    "            invert_transform_flags=True,\n",
    "            interpolation='Linear',\n",
    "            dimension=3)\n",
    "```"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The end goal is to generate warps between each unique pair of the templates (forward and backward).\n",
    "* forwards (A to B, e.g. *NMTv1.2_to_D99*)\n",
    "* backwards (B to A, e.g. *D99_to_NMTv1.2*)\n",
    "\n",
    "![warps](RegisterTemplates.png)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Step 0: Preparations"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 0a. Import required libraries"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "import os\n",
    "import glob\n",
    "import time\n",
    "import shutil as sh\n",
    "import nibabel as nb\n",
    "\n",
    "from itertools import combinations\n",
    "from matplotlib import pyplot as plt\n",
    "\n",
    "import nipype.interfaces.fsl as fsl # nipype interface for FSL\n",
    "import nipype.interfaces.ants as ants # nipype interface for ANTs\n",
    "from nilearn import plotting # Plotting function from nilearn"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "# test if ANTs is detected\n",
    "!which antsRegistration"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 0b. Define relative paths to template files"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The skull-stripped isotropic volumetric images are used for registration.\n",
    "* For **NMT**, **D99**, and **YRK** the provided skull-stripped brains were used\n",
    "* In **INIA** the brain stem extends further down the spinal cord compared to the other temlplates, so the braisn stem was truncated at a level similar to the others\n",
    "* No skull-stripped image was provided with **MNI**, so the brain was segmented semi-manually using ITK-SNAP\n",
    "* For **YRK** we used the version provided together with [NHPPipelines](https://github.com/Washington-University/NHPPipelines)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "/media/DOCS/Documents/PRIME-RE/RheMAP/GIN/RheMAP/\n"
     ]
    }
   ],
   "source": [
    "# ============================================\n",
    "# NB! THIS COULD BE DIFFERENT FOR EVERY USER\n",
    "# ============================================\n",
    "BASE_path     = os.path.dirname(os.getcwd()) + '/'   # repo base folder\n",
    "print(BASE_path)\n",
    "\n",
    "# these follow the directory structure as outline above\n",
    "TEMPLATE_path = BASE_path + 'templates/'                   # templates base folder\n",
    "NMTv12_path   = TEMPLATE_path + 'NMT/NMT_v1.2/'\n",
    "NMTv13_path   = TEMPLATE_path + 'NMT/NMT_v1.3/'\n",
    "NMTv20_path   = TEMPLATE_path + 'NMT/NMT_v2.0/'\n",
    "D99_path      = TEMPLATE_path + 'D99/'\n",
    "F99_path      = TEMPLATE_path + 'F99/'\n",
    "INIA_path     = TEMPLATE_path + 'INIA/'\n",
    "MNI_path      = TEMPLATE_path + 'MNI/'\n",
    "YRK_path      = TEMPLATE_path + 'YRK/'\n",
    "ONPRC18_path  = TEMPLATE_path + 'ONPRC18/'"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [],
   "source": [
    "NMTv12_brain = NMTv12_path + 'NMT_SS.nii.gz'\n",
    "NMTv13_brain = NMTv13_path + 'NMT_SS.nii.gz'\n",
    "NMTv20_sym_brain = NMTv20_path + 'NMT_v2.0_sym_SS.nii.gz'\n",
    "NMTv20_05mm_sym_brain = NMTv20_path + 'NMT_v2.0_sym_05mm_SS.nii.gz'\n",
    "NMTv20_asym_brain = NMTv20_path + 'NMT_v2.0_asym_SS.nii.gz'\n",
    "NMTv20_05mm_asym_brain = NMTv20_path + 'NMT_v2.0_asym_05mm_SS.nii.gz'\n",
    "D99_brain    = D99_path + 'D99_template.nii.gz'\n",
    "F99_brain    = F99_path + 'struct_brain.nii.gz'\n",
    "INIA_brain   = INIA_path + 'inia19-t1-brain_truncated.nii.gz'\n",
    "MNI_brain    = MNI_path + 'macaque_25_model-MNI_brain.nii.gz'\n",
    "YRK_brain    = YRK_path + 'MacaqueYerkes19_T1w_0.5mm_brain.nii.gz'\n",
    "ONPRC18_brain = ONPRC18_path + 'ONPRC18_T1W.nii.gz'\n",
    "\n",
    "temp_names   = ['NMTv1.2', 'NMTv1.3',\n",
    "                'NMTv2.0-sym','NMTv2.0-0.5mm-sym',\n",
    "                'NMTv2.0-asym','NMTv2.0-0.5mm-asym',\n",
    "                'D99', 'F99', 'INIA', 'MNI', 'YRK', 'ONPRC18']\n",
    "temp_brains  = [NMTv12_brain, NMTv13_brain, \n",
    "                NMTv20_sym_brain, NMTv20_05mm_sym_brain,\n",
    "                NMTv20_asym_brain,NMTv20_05mm_asym_brain,\n",
    "                D99_brain, F99_brain, INIA_brain, MNI_brain, YRK_brain, ONPRC18_brain]\n",
    "\n",
    "# strength of contrast in the plot_anat function (\"dim\")\n",
    "contrasts = [0]*len(temp_brains) # create a list of zeros because this should exist\n",
    "contrasts = [-1, 0, 0, 0, 0, 0, 0, -1, -1, -0.5, 0, 0] # specify to improve display"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "In the following step you can create a subset of the templates to include in creating the warps. This can be useful when you add a new template and want to warp it to templates that have pre-existing warps (to each other). This way you can generate the new unique temlates without redoing the prexisting ones."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [],
   "source": [
    "# include_pairs_with   = ['NMTv2.0-sym',\n",
    "#                         'NMTv2.0-0.5mm-sym',\n",
    "#                         'NMTv2.0-asym',\n",
    "#                         'NMTv2.0-0.5mm-asym']\n",
    "include_pairs_with = temp_names # all included\n",
    "include_pairs_with = ['F99']"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The headers of NMT > v1.3 has both qform_code and sform_code set to 5. This is a relatively new definition ([see NITRC](https://www.nitrc.org/forum/forum.php?thread_id=10029&forum_id=1942)). Older version of nibabel (< 2.4.0) only recognized values 0-4 and reset the 5 to 0 which causes weird plotting artefacts (registrations will still be ok). As a workaround, you can use the code below to change the header value from 5 to something like 2 ('aligned'). A better option would of course be to update nibabel.\n",
    "```\n",
    "if NMTv13_brain in locals():\n",
    "    sh.copyfile(NMTv13_brain,NMTv13_brain + '.bak'); # back up the original file\n",
    "    img = nb.load(NMTv13_brain)\n",
    "    hdr = img.header\n",
    "    hdr.set_qform(hdr.get_qform(),code=2)\n",
    "    hdr.set_sform(hdr.get_sform(),code=2)\n",
    "    nb.save(img, NMTv13_brain) \n",
    "```"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 0c. Plot template brains"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {
    "tags": []
   },
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAqgAAAEYCAYAAACQmgS0AAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/NK7nSAAAACXBIWXMAAA9hAAAPYQGoP6dpAADTRElEQVR4nOx9aZCkWXXdyT0rM6sqKyu7urq7epmenp5Vw8DADCHEKrQRSMYEQkZGYZAUQpasxQpZhBABEQIRDi2WZEuWbSGzCslYEIMCa8GKAbFIMAszwDDMPr13dXVWdVZWZlZmVi7+kT6vzrv5ZdMzdHdV97wT0dFZmd/yvve95b5zz70vNhgMBggICAgICAgICAjYJohvdQECAgICAgICAgICFMFADQgICAgICAgI2FYIBmpAQEBAQEBAQMC2QjBQAwICAgICAgICthWSF3rgsWPHUKlULmVZLhvK5TL27du31cUICAgICAgICAiIwAUZqMeOHcP111+PVqt1qctzWZDNZvHoo48GIzUgICAgICAgYBviglz8lUrlqjFOAaDVal01bHBAQEBAQEBAwNWGoEENCAgICAgICAjYVggGakBAQEBAQEBAwLbCBQdJBQQEBAQEBFwY4vHA/wQERKHf71/QcaEHBQQEBAQEBAQEbCsEAzUgICAgICAgIGBb4aIaqIPB4Lz/PvvZz46cc+edd+Kuu+7C2bNnsb6+jkcffRTvec97MDExMfY+e/bswZ/8yZ/g6NGjaLVaOHnyJD7wgQ9g//79F/NxAgICAgICAgICtgCXRIP6wQ9+MPL7Rx55xPv7x3/8x/GhD30IyWQS9913H44dO4YXvvCFeOc734nXvva1eOlLX4p6ve6dc/PNN+Puu+/G3NwcnnrqKXz605/GoUOH8Ja3vAWve93r8NKXvhQPPfTQpXisgICAgICAgICAy4DYYDAYfLuDvvrVr+L222//thfjpWKx2Lc9ds+ePXj88ccxMTGBt771rc6oTafT+MhHPoI3vvGN+JM/+RP83M/9nHfegw8+iOc973n4sz/7M7ztbW9Dr9cDAPzyL/8yfv/3fx8PPfQQbr31Vny7x7r//vvxghe84NuWMyAgICAg4JkiBEkFBERj2wdJveUtb8HExAQ+85nPeIxrp9PBz//8z6PRaOCnfuqnUCqV3G8veclL8LznPQ8rKyv45V/+ZWecAsAf/MEf4J577sEtt9yC1772tZfzUQICAgICAgICAi4itsxAJSP7uc99buS3SqWChx9+GOl0Gq95zWtGzrnvvvtGXP8A8I//+I8AgH/xL/7FJShxQEBAQEBAQEDA5cCWGaj5fB4AcO7cucjfV1ZWAADPe97zvqNzAgICAgICAgICrixsmYF69uxZABgbec/vDxw48B2dExAQEBAQEBAQcGXhkhio49JMTU9Pu2Pojn/Tm96EVCrlnX/nnXfihhtuAABMTk667z//+c8DAF70ohfhxhtv9M7J5/N4wxveMHJOQEBAQEBAQEDAlYXLmmaq0+m4z3/+53+O3/iN38D+/fvxqU99Cr/6q7+KY8eO4SUveQn+9E//FBsbG0ilUl6012OPPYa/+qu/whve8AZ86lOfws/8zM/g3nvvxaFDh/CHf/iHzgC+0AixgICAgICAgICA7YdLYqC+9a1v/bbHNJtNvPa1r8WnP/1p/NAP/RB+6Id+yP129OhR/Kf/9J/w9re/fURv+tM//dOYnZ3FK1/5Si/xf6PRwDve8Q78zu/8zliNakBAQEBAwJWIvXv3olwun/eYSqWC48ePX6YSBVypuJC2pNiqdnVJDNQLxUMPPYQbbrgBP/qjP4oXvvCFSCaT+NrXvoaPfexjeOc73wkA+OY3v+mds7q6ile96lX4gR/4AbzqVa/C9PQ0jhw5go997GO46aabIs8JCAgICAi4UrF371488sgjmJiYwAc/+EG89a1vxb333osXvvCF3nHr6+u44YYbgpEaMBZRbYlIJBLYsWMHXv7yl+M973kPrrvuOgBb16621EAFgFarhY985CP4yEc+4n3/6le/GkB0GioA+Pu//3v8/d//vffdL/7iL573nICAgICAgCsN5XL5vNt/ExMTEyiXy8FADRiLqLb0gQ98ADfccANarRa+9KUv4bd+67fw2c9+Fo888ghmZma2rF1tuYEahZe97GW4/fbb8dBDD+Gf/umfLuicmZkZ/Jt/82/QbrfxoQ996BKXMCAgICAgICDgysctt9zi2PhXvOIV6PV6ePe734277rrrgiSblwpbuhfb8573PCQSCe+75z//+fjYxz6Gfr+PX/iFXxg557rrrhuJ0t+xYwc++clPolwu433vex9Onjx5ScsdEBAQEBAQEHA1gsbqmTNntrQcW8qg/sEf/AFuuukmPPjgg6hUKjhw4ADuvPNO9Pt9vO1tb4t01f/4j/84fu3Xfg333XcfTp48iVKphJe+9KXI5XL4wAc+gPe85z2X/0ECAgICAgICAq4CPP300wCAw4cPb2k5ttRA/ehHP4o3v/nNuO2221AsFnH27Fn85V/+JX7nd34HX/va1yLPufvuu3Hbbbfh9ttvx5133ol6vY4vfOEL+G//7b/hrrvuurwPEBAQEBAQEBBwBaPX66Hb7ToN6nvf+1687GUvw4/8yI9sabkuqoEai8We0fF/9md/hj/7sz97Rud84QtfwBe+8IVndE5AQEBAQEBAQMAoXvziF3t/33jjjfjUpz6FZHJrw5S2VIMaEBAQEBAQEBCwdfjwhz+Me++9F3fffTfe9ra34Vvf+hbe9KY3bXWxtmcUf0BAQEBAQEBAwKXHjTfe6AKjXvnKV6LX6+H973+/27lzqxAY1ICAgICAgICAAADAb//2b2NmZgbvete7tnTr+MCgBgQEBAQEXGG4++67ceTIkZHvt9KgCLg6MDMzg1//9V/Hr/3ar+FjH/sY3vzmN29JOYKBGhAQEBAQcIXh7W9/e+T3N99882UuScDViF/4hV/AH/3RH+E3f/M3t0yPGlz8AQEBAQEBVwje8pa3YDAYjP2XyWS2uogBVwjYlqg/VWSzWRw9ehSPPfbYyIZKlwvBQA0ICAgICAgICNhWCAZqQEBAQEBAQEDAtkIwUAMCAgICArYxKpUK1tfXv+1x6+vrqFQql6FEAVcqLrQtKbaqXV1QkFS5XEY2m0Wr1brU5bksyGazKJfLW12MgICAgICAb4vjx4/jhhtu+LbzVqVSwfHjxy9TqQKuRFxoW1JsVbuKDQaDwYUceOzYsatmZVYul7Fv376tLkZAQEBAwFWKeDw4KAMConChqdAu2EANCAgICAgIuDAEAzUgIBoXaqCGHhQQEBAQEBAQELCtEAzUgICAgICAgICAbYVgoAYEBAQEBAQEBGwrBAM1ICAgICAgICBgWyEYqAEBAQEBAQEBAdsKwUANCAgICAgICAjYVggGakBAQEBAQEBAwLZCMFADAgICAgICAgK2FYKBGhAQEBAQEBAQsK0QDNSAgICAgICAgIBthWCgBgQEBAQEBAQEbCsEAzUgICAgICAgIGBbIRioAQEBAQEBAQEB2wrBQA0ICAgICAgICNhWCAZqQEBAQEBAQEDAtkIwUAMCAgICAgICArYVkhd6YDwebNmAgICAgICAgIBLj2B1BgQEBAQEBAQEbCtcMIMaEBAQcLUgmdwc+vr9PtLpNIChp6jb7bpjOp2OOy6dTru/eXyz2bxcRQ4ICAh4TiEYqFuAdDqNbrfrZBMqn4jH494kWSgUsLKysiXlDAi4GsD+RaM0l8t5hme323W/dbtd9Pt9AECn00E2m/X6KX/L5XKo1+vIZrMAhkYurxGPx1Gv1y/DkwUEBARcvQgu/oCAgICAgICAgG2FwKBeQsTjcecKVJYGGLKo/K1WqzmWJp1Oe4xqrVZzx/F8vWaz2XR/FwoFx9x0u13H1EaVq9VqXdRnDQjYDlBXvQ3sZJvXPkFWVBlVsqTpdNrrJ/F43OvDrVYLuVzOO4fHTU1Nub9VBtDv90fKFeQCAQEBAaOIDQaDwYUcGKL4x6NQKACAN5nRjagT17jPOunxN37X7/e9CVXPowHK75LJpHesvS6hxjAAZLNZVKtV93uUURsQsNUoFArI5XKuvXe7Xc/o7Pf7zuXebDZdO6fGVI1JLuR4vF08Eslk0hv72Dco0eH9+/2+K5f22WQy6f3Gv/X6PFa1sCwHDedkMhlkAwEBAVcFdAw8H7bMQN27dy/K5fJ5j6lUKjh+/PhFve93ikKh4JiObDY7omVjPZHN5N+5XM595nnjDEjAZ4J4PzKyOkHrpGbLomg2m24y7nQ6IywqGZ9Wq+U9o15Hr7/dcCHtidiO7SogGlNTU55RCIwu0oBhu1WGk4tGvYb2I9WZAnDGazweRyqVAgBkMhkkk0n0ej13LV5/MBig1+uBw2e323XntVotJBIJd33tQ9ao1gWmXp/fq36W16HhGuWdsUFcem09R48NCAgIuJy4UAN1S1z8e/fuxSOPPIKJiQl88IMfxFvf+lbce++9eOELX+gdt76+jhtuuCEYEwHnBdvT//pf/wtvfetb3feJRAI7duzAy1/+crznPe/BddddByC0q4CAgICAgO2OLTFQy+UyJiYmvu1xExMTKJfLW25IJJNJx1C0Wi1PM6ZaN2U3k8kk5ubmPAYkFosBGDKog8HAMS5kY/hZWZvZ2VnHuDQaDQwGA0xOTrpjV1dXXdkSiYQ7j99tbGwA8JlQMj/8bmVlxbEp3W4XtVrNMVDKvEbp9bYDbHv6wAc+gBtuuAGtVgtf+tKX8Fu/9Vv47Gc/i0ceeQQzMzPbpl091zE3N+c+M1NFMplEOp1GrVYDMOxv6i4HfI0omc+pqSmP/cxms8hkMgCAdruNVCrltXfNkjE9Pe3aOFlTotfreX/zc6/XQyqVcv1rYmLCfc7n8yPHE61Wy/2eTCbRaDTcb6pFTyaTnodD2VT2P/1d60hZWtW0K3PMY1WWEBAQELCdEIKkxoATFuC76XK5nPuczWbRarWce7xQKLjJQCdAez1OopzQ7CTW7XbdxHXkyBFnfGWzWZw7dw7tdhvA0EAtFAru/Gw26+5LY3h6ehrAcELl/Wjo0uW5c+dON5lzcqtUKgDgGaR8Vk21Q0NiO+GWW25xbPwrXvEK9Ho9vPvd78Zdd93lMawBlxdc2FE7Td2zup6pD9WUTaVSyV0jl8u5/gPA9QW2TS7edLFWLBadyx3w++Lk5CQSiYTXb3heIpHw+qZeU69HqHtev4vFYu67yclJJwvo9XrI5/PuWvl83ltQplIpd6zVuPb7/chAR6ao4/1UIsSFti4wVfbDa/PYcflgAwICAi4HQuRTwFUPGqtnzpzZ4pIEBAQEBAQEXAgCgxqBbDbrBVGk02kvOpifi8Wi5+JnoAYwZHqKxaJjMpWRoBtfXfNkgsjakHldXV3F8vKyu0YymXTH7t69G9ls1rE63W7XlY0MDFkWdV3yHDJRysrSTXrgwAEAQ0aJ11hdXfUCvzTNjg0G2U54+umnAQCHDx/e4pI8t1AsFt3nfr/v2HYbAKTHFQoFpNNp16YLhQLOnTvnzmu1Wh7bx2vMzc1hZmbGtUFlOOkuZ5/qdruOmdzY2PDkN/QkAEOGU1nUVqvljrPHx2Ixd812u+3uz/7Dfsm+z+sPBgN3LQ124vfKtvIzpQX0rGicqw2AtKnubCor1iWZV/UAWUmBSgW2a18PCAi4ehAM1P8PdfupHovpojggX3vttc4N12w2vXQ22WzWRZJrJDAwnGB0Mut2u95vev/BYOAmmOnpaXe/jY0NxGIxp91LpVJot9vuOqqD43X4mxrEnEBpGGiaKRqh6mLlBFcqlbBz504nP6hWqy71TTweHzHstwq9Xs8ZM1/60pfw3ve+Fy972cvwIz/yI1tWpucC5ufnnZY0m8167QjYdPFPTU0hnU57elE9zhpq2jey2ay36NPfNjY2PPe/alCTyaQXZa8a7SiNKe/FfMLjwIWdXp/PwPInk0lPQqNZAbSciUTCK7OVHii0P2u52+22k/7wfmrQz87OuvvbrB+6+G61Wq4vs29rNgX7vgICAgIuNp7TBqoOxsoKdLtdxwwWCgXk83k3EZHNBODyLnIgn5yc9CbWQqHgJgNlZuLxOPL5/IimjJ8BPxhjZmYGAHDy5Emk02mPOcnn896ExwCMRqOBVCrljOd2u+0+p1IpJJNJL10WtaqDwQD1et0rD5kvGuOsm927dztDulareZPawYMHnbHSarUu68YAL37xi72/b7zxRnzqU586b1qvgGcPtivLrgOb7Xlubs61m1gs5hmTamzZ1E5qbAH+Yq5er3s6yVgs5q6pRiW10nosWceJiQmP1bSMZjabdX1fr7mxsYFEIuEZjhr0yDK2223PiNbyd7tdz3uSyWS8IEdNXaUsLxlg/q1BlryGPp9lWnm/fD7vGdJ2ganbuMbjcS8nq03PpXVrmdeAgICAZ4OgQQ246vDhD38Y9957L+6++2687W1vw7e+9S286U1v2upiBQQEBAQEBFwgnrOUkm48oO4rYJi2SHVxqntT9oBpW8is5PN5z80HwLEX6ipkFLMm41fNmOrlGo2Gp1VNpVIjbKxqS5XN1c0ANMUM2WKNOlZ3a6fTceXJZrOOJcnlciNsFp+vWCyiXq+76P9jx445BoaR/5Zdu1S48cYbXWDUK1/5SvR6Pbz//e/HX/3VX+ENb3jDJb33cwlky2zaNf62sLDgZawgc1coFFAsFj3WjZpNsqfaNvlbo9FAPp937adQKHgazUwm4+kkrUxF2Un2IfVe8Dc9r9vtuswA2k/J9Gq0vCbSt9dUllRBLwifVfu6Hk/ZCqH6clt/Wi6VCbAcek0rpeDYNjU15foyJRgcE3XsouZVtej67jTNVdjKNSAg4JngOWWg6uRTKpW8XIvlctmli5qYmPC0mzoxtFotN9Gtra0hFou5gZuuc2A44OvkzQGc1wTgBS/ws01VMzEx4Sa76elpxGIxLx+kuh3t/uPW3aYGos17qIEWuVzOTXiJRMIZoYVCwXNr5vN5N6FtbGwgm826uqhWq1hcXAQwzDup5dJyXg792m//9m/jE5/4BN71rnfh9a9/fdi29zuALoKYgkzTrlGDPTMzg/X1dayurgIYGlp79uwBMHTxJxIJ18bU+Op2u4jFYq7NqYE1PT3tBQLaNFAAvH6kxnIymXSGpqZOYy5Qu70pj8vlct6Cyu7epO3YGrY8Xo1c2+9s+1cjWoOkrAa01Wo5WY7dqUrHCJX/2BzJmv6K8giVA+juclHbt/L57JavNm8t61oN4Gw2GwzWgICA82LbGKh33303jhw5MvJ9SCAd8J1iZmYGv/7rv45f+7Vfw8c+9jG8+c1v3uoiBQQEBAQEBJwH28ZAffvb3x75/c0333xRrs+IXxq8i4uLjiHYs2cPpqen3Yp/YmLCc5l1u13nopuYmPBc6vV63XPHWRejBhooo2BTtfA37qSjAQbKBCkjQvZCpQMavDAuSIEMJq9rWRZNk5PNZjE7O+vur4yVMk/JZNJjVHu9nmON6vU6jh8/7n4rFovuc7PZvCyLkF/4hV/AH/3RH+E3f/M3gx71WSAqlRrbHlOSAZvR5mfOnEEmk3FJ9jUdGiPntb+xDdBVz/7Gc4HNoMSofeZZFpZNJSU8T2U09GhEufgt268eDw1itDszEdov6PKOYu05JigDbX/XrBw2kEzrh/VF1tWymARlQgBGMokAm5uH6Jig4yHP0+srE6pML+uWf2taOsqYWH+BTQ0ICLCIDdS3ex5cTLfo85//fNx///0XdOztt9+OBx544Fnfi4N+Op1GvV53f8/Pz3sTQrFY9CYiqyfTici61OyuUfysk5saqBy01U1G8F48tlqtunKePn0aGxsb3haR/X7fudX1flp+wHfr894asaufNfJ4YmJi5FpqTOukr393u123lWgmk/F2DqpWqyMTJ+/xbCJ/n0l7Ir7TdvVcAd+rGiSaLmpmZsa5oWu1mrewyefzrh2p3pF5PHmsji38XndsUoNUjVD9zLydusuZ6kN5Pv/X1Et20aXQZ9cFIN3Veg/9TWUB9vraD61RalM5qUE8rv+yfvmbps6yKaEUugjgfWmgqtHL9xRl9DIjA9FsNke2YNXFjW7BajXylzPTR0BAwNbhQkmpbcOgXkxosnwOgM1mE8ViEfPz8wCGgzMH3HK5jHQ67Q2WZFetwaTMDw0ynWA18CmbzUYmvraMjg74zWYT6XTaYxR0krK6O2VilZHI5XJjmQ07gabTafe3Tuz6HFH1m0wmPSZKA8iSyaTTIz799NMuMAYYutyZfJ37orPerNYu4PIjyqBR7Ta1j5lMxtvXfnp62v1GI5SMqmpOmWR+fX0dgN8XeT7bmBo/UX1NA/+0D6vhafXOUQwnEbUYjTL2+FnZwShviQ2Y0jady+W8slijVo1SPY7eD/a1TqfjWFFqxrVe1LAFNo1Z7deJRMLLD9vtdl3wFs+xAaBaZr7nqakpN3bS88R3qPXC+tKFacitGhAQoAjRIgEBAQEBAQEBAdsKW8KgVioVrK+vu4jRcVhfX3epTi4UyiJ2Oh1Pkzk/P++xFGT4pqamPBZPWZR0Ou2xIFYbms1mI3ehOp80QHee0u+BzTRW6k7Xe2rEs9XdWehvnU7HS6rN1E+Av7MO7606Q4Ll1e+0XrQOlW3asWOH23EKGLJru3btAjB8X6dOnfLKHSUhOB8utD0Rz6ZdPRdA5kt3gMrlcl56sWw26yLzV1dXkcvlsHPnTgCjid9zuZxjz5TdXF9fR7/fd9KAbreLtbU1AMDZs2cxMTGBvXv3uvuxDfR6vZF0Zbwmy6cMvm2LbKs2rRwAj3HU7/RYqw3nfVhP1gXP8li9KK/BcmnaOu13VmOrDLCWX8tI5tW60vmbSi2AzewdPE936SLW19eRSqU8V79uDKDaeEoMiEKh4P7mpgjA8F32ej3n5dGx2tZtYFMDAp6b2BINKgDs3bvXGYjjUKlUnI7xQqBGFzAsM136qVQKa2trnraSEyuAEaPM6tzUNWUNWZ3cbJBG1CAbFRCghqTNM0jQLc6gpX6/7z2zTf+i17WuRGu8qmtUNXpqWPBvnWDtdomqQeUznzhxwm3ByN9oTGazWWxsbDgDdnFxcUTCcCFG6oW0J+KZtqurFSrPAHyNIX87cOAAGo2GC6apVqueNKVcLrtgOABeWiQaeFGw2kR+PnXqlNd2ZmZm3E5qdpc1DRCkLMe60HmejgvqaqYxadNO6XH6t+3f2vfV0FRpkWra7ZigOlRbJ1Y/aqUX+t50gWnTX9lFso5zlFkwMI3vT7WlNGrVuIz6zDKq7CGVSnm/685bGxsbnpyH9ccdwqLSWlnjPCAg4MrDttegHj9+PBgJARcNoT0FBAQEBARcPbjig6Ti8bgLvqHLj2zCzMyMcz3VajXE43HHXlp3sI3mtUm+lWG0qaSURdFgI7sHvU2/EpVKyqae0WPW19fRbrddIEm/3/ei89PptGNjyKQoO6MuTr22BnPRlagBHwpld5VVZ92rfIF1zRRDLDczHPCZNGBmZmYGZ86cATB8ZzbgJKSjefYoFArOjR+Px7GysgJgs03TTT47O+v61GAwQKvVcu8nlUo5FqxQKHh9RoNolCkEfJZP+4z9e2FhAclk0skIzp0759j1YrGIfD7vsYVENptFrVaLDLYj02qzb2hZx8ltlHHsdruujpTR4zX0PGVktV+oq94GQll3vPVa2b6o947KQMK/tc/Ysug7m56e9rIb8DeOMyrxUS+KZVE1cEo39uh2N3cFW11dRSaTcb9ls1l3PwZtstwqOeEzBBY1IODqx5a5+L9TcLC0bq+pqSnP5Ugkk0msr687gzWfz2P37t3uWuoC13ypHNDVLWcNVDU8FeqO0nyfdnDVa3KCUQ2ZusoBuHKzbJqmivXCnINRRrc1xq2GTY3QZrM5Eg2tz6HuPP2t1Wo5ScLa2hrK5bL3zjj5MGuAXpOT3+LiIo4cOeIZ1lFa4IBvD2qzVZ6iWR1mZmbcAqLX641INfhbo9Fw75HGItuHbn1rJSq6CDpfpLweT9CQZso13iOfz3vXUa24pnyjPtQakby3/jZOm2p/i4rO1z5jZUC6wFVNq0qGrDve3lMNWT1Pn4fGpM3lqp+17+i99Teb1souKoj19XUkEgmXkooZGgC4HNHWYAVGtztW17+N8G80Gm6RUq/XPSlVMFQvLq6//nocP37c9T+Oz9wVMCDgYuBCXfxXnIFqg5Z08J+bm0OxWHQDZCaTcQNwrVbzmJlut4sdO3YAGKbOGTfgMzF4VGAQP48LbOh0Os4Q0+85gekEYCdPNfhOnjwJAC5tDyffvXv3epONav5YL1EBJfz/fPlTbT3wszVKef16ve5Nvt1u1wUiTUxMeKlvuLc3y2l1vNbIffLJJ119KtMV2NRvDzKmADwPQjab9QLVWq2Wx4Lx3dHwUmac7UaT5gN+X6RBo4sS24fUyLABhFGGUqfTQb1ed/l0U6mU02PToFEmkUgmkx5rz+OBzYWjeky0DesCULcL5jm2//E8whqy9ndt+zbFmtXK6jNpHSlLyv6jBrLW6fm05/Z/20f1fdkFhZabxibHHI5lmuy/0Wh4f1sWlt8RDKJjWbhNtV3sBFw49u/fD8DfvOHMmTNee6In8NChQwCG75c7PkYFuAZcfLAvFwoFL9iYMQDse1fSu7jQPrs9rM6AgICAgICAgICA/48rhkFV/ZMyEOVy2TE6k5OT3lak1mWlLrTV1VWXiFp3keJ56mLUvy1jOk7XRpDtATaj9xktGxW9G4/HXRQrAO8zWQaWdX5+HqVSaYSxIaxb02Y40BQ9lik9XxRwlD5WGVpguOsVMTs7OzaFDpkamw2Ax+nONEeOHHH1yXQ9YfeZ8dD3yGwMzFzBPgMMUzsx4hoYslWWuaOLv1gsOgkNGcUoto5pnpQdtP3GRrPzGjaCW3cr4/HAkO3hb4VCwduuWF38ZHJ1Vzn9TftxVB/WPqRtMyobgF5Xv9cxSRlhKxGy0fhR7ZssqPXWsL6iGGmeZzXzysIqQ2xT3VmXP5kcq8NXUDsaNUYyNRUzBeiuVjxOt2GlBtpmDqHnRusg4PxgdouojSVKpZKXeaFQKKDRaLi/G42Gexfs38yU89WvfvWyPsfViuc973neRjZW565jXrfb9bYP5rjebDa3dSrFq8rFn06nvZcAbOYwLZVKrnNNTU15g6VO0AxK4qCrL487SWlDsAaqnej1t6jJgGliaJS2Wi3npqI2Kyp4IJlMol6vewO2bvvY6/XcfVKpFMrlshfwZfOn6mRsNaj6jKxfqzfV56O+TaHGar/fd/Var9edIbNjx44RI1z1dFHaO17b6njpXlpaWvLKr66P5yqscabb5M7NzbkAJPsbJ6Tl5WUAm+nMgOE7KJfLbhKyqZ3UHah9KJfLjRw7LsjHDrL2Nz2PGkRek++9Uqmg2+26cloNKsun9aSfo4xL9hOVReizWvmB1plFVDCVdZ2fz11tZUDjdK1R97Waeb2mPqv+HhXgpmm07HMRnU7HG68YRAUMA+50ZzE1ULvdrhvzNKUV4G+Tm0wmUa1W3bH6fHT7B+nPKO688063HXQymUSpVHLGzO7du53RybmEdXrmzBnMzs66PkQ9OAA89thjmJqacn12ZmYG3/zmNy/L81xtuPbaawHAjcEcqyYmJtxCLZlM4lWvepVz5R85cgSFQgFHjx4FMCTdOA7wfXEcfOqppy7Tk1wYgos/ICAgICAgICDgisS2ZlBp/asrKh6PY+/evS44IhaLuRUGGZuoNDT8XdnAs2fPAhgGU83NzXlMhkb0RwUY8Xo2i4AN/CDDwwT8wHCF32g0vITnvE4ikRiJjOZxdnVFhpRMWKFQGMv8skwAXKoqfUaN2mQdsGzKSNigEGWXTp8+7RjUnTt3eu9P2SYtI12LUW5Au8pqtVru/o899piXfkbr+rkKZcqAoUeBgVDFYtHLYjEzM+MCqMjAsb7PnDnj+kY8Hsf+/fs9ZkXd+Ar1WNj0SvY4wHeXWxbRBhUCm5KCqCA9YJj1gW3gwIEDjiEi+6gR/hpJz/uyLFGppvhZswyoh4IeBL0fy2jrRb0O49hI2wdsX1A2V+vUsrJ8JpZF+54m5rdZF6zXSoOrms3miDeG5Wu3204qEo/HR6Q/ZD65MYDdDACA+45s68bGhmu3PIaMn3p0WBbWHVn1gGH6NjLMDJDkbm2lUsm9d/Ypnbempqa8TV60PVUqFRfhr/1quzF22xEMUltbW3NzJT2bHLu0Tufm5nD48GEcPHgQAPDJT37S81o1m03HrrZaLW9TlTvuuAP/8A//cNme7dvhqnDx86XpILN///6RHYOoJbVaKE3RpAMvMKwguiva7TZmZmZGDCoeZydaq82yOi5eQ138UdrNqKhfPoMaurzf2bNnsby8POKW0+toFLO6/jXnII1TNUo5Gdn3rMaC1czV63VXtpWVFaysrDi3/uzsrGeQWlemTpRTU1OeFMJGcNs6BYYu6kceecR1wGw267mfniuI0v8xDdns7Kw36U9MTLj3ows5G5HebDZdSrNms4nDhw97EhBNT2V1rgTbirr4bbuP0mDzWDXQtL3rIkSPo4HIHLrLy8uuHubn50dkJda1b6Pz9RnU8LcLUDVeVZepaeL0+vZ96bPrImFjY8MtvIGh9s+6X4nV1VVXl4lEAvF43Bl5qVTKW1RqOifNP8y/tVy6SO33+94WpvzMMXCcG193IQM2DdRutzuStYCwu1pphH8sFvPuMRgMnOwnmUx6W85qm3sujg0A8KpXvQoAcM8997j57ZZbbvH6JhC9m5z+pn1Ax+NWq+Xe6fHjx911MpkMMpkMHn744Uv5eFcsbr/9dk9/fc011wAYjolnz551kj/duv306dP4ru/6LvduVlZWUKvVUCqVAAzn0aWlJQDDeXR1dRWPP/44gOEYeNNNNwEA7rrrrsvzkOfBFW+gFotFb6W2sLAAYDj52q0Drc5S95zXgU/TqmiATaVSQTKZxNzcHAA/oIMTclQKKk5eWjc2jyMnLbvfPcsHjAYz8Tt7zbW1NVQqFXdt7oWuabXIXvDanBxsjkgtdzab9Yxz+zx2EtHyM2ipUqlgZmbGpe5SsO50cWB1gXo/vhfeW+tN0+mcPXvWdUB9RxqYdrVCDXxtc/Pz8271rZs6MD0bkU6nR/avZ/tYWVlxqc0SiQR2797tMeuqU7a5drXd2mA7y9JFaZJptFyIJlmNVRqovOba2ppjdjKZDHbs2OH1BS2jGsE6tuRyubEsp53QrSHdaDS83zU3sxqtVtOt7LIaj+fOnfN0tTTUgGF753awmUzGC67q9Xqehg3YNBJtmifqQnkdgsyrnqcGaa/XcySBGpNMOaVGjxrOg8HA052yfAwEVY8PwTGdZe/1eq6ul5aWvHadyWScnprs4XMJr3vd63DPPfcAGD7/HXfcAWCTMIhil61HwXrQ7JbgulFMLpdzzOni4iLS6bRrR9RJBgCHDx9GvV53bXJ2dhbf+73fC2A4j545c8aN29Zjq14QLorV3rHeKS7Yjx8/jgMHDgAADh06hI997GOX+CnPj6BBDQgICAgICAgIuCKxLbc6pZuGq4D5+XkXnWu1bcr4xeNxL32OshpkMJVt4mqx1+uhUqm41Eh79uzxNDfqurdbeCqDZF3XyhrpapWuSavh0/OUodEUW4PBwF3r4MGDmJqaimS34/G4twsUy6fl5Hk2JY+VDRBW09vv9x1beeDAgZEdqVRuoOyW6tusG1+lB9zuUNlyntfpdDA7O+s0r+fOnfMY4qtZj6paY11BHz58GNls1rk/i8Wi05mSYdR3wnZcKBQ8djqdTjtGjJozG2UPjLr4AXhs9zhXPXdCsjpGXl8jgwGMsLDKQJKFsK70dDrtpcPSKHHteyyHlsVqXKOi5a1WVRkB1ZerhoywkqCoTBvAkAEkO7iwsOB5MjSJfblcHmFONJKez93r9UaeXd3x+j+PJfR+duxSFlN/y2Qy3t/caAQYtlvdsGMwGDjmld/pNqg22b8ey/pbWFhAr9dzrGk8HnesUb1ef87pIk+dOuWkDTfddNPIrmGKKHc/P2uf03PZb7Xt7du3D8DwvVQqFbcxzvXXX49HH330UjzmFYE77rjD8zZ2Oh1XVwcPHvS8hhx7AX/c5FirXkT1NKj8gtIASgdarZaTwujmLdsd28rFr3kKO52OMzb379/vUd46OcTjm3uI2y1LAXiTDeC7F9W12el03GTXaDRcvkgaoFEaNXZOHax1n3urN+P9mNt03DWti1MnirNnz7p7zM7OekajDWDRCVXLws+sXwqzWUatm263631WnDhxwp1HXbDKH6zGkcfqxBQ1WFqjRgX7+s7q9bozUKvVqtNN0vjm8/H7qwGckNmei8Wim4QnJiaQSqW8oDlt79o3Wq2Wt6+8LgQ6nY7n4l9YWIjMA1woFEZcftpWrFwkarFmz7OyGW23utjkdbTM2g9VKsLnjkqTRIzrb/F43EvBQxc1jS8afxMTEyOBXnzmqOfjcVpmq/vUcU7Ps4FHWn4uwtWY1WAjrT8aiSx/Mpl0z0qpEOBLAQBfXmBzM/M7nmfftWpHWV6eo2OH3t/mW1YDXXWxyWTSkxRsbGx4i4FisficydX5kpe8BN/85jddMPGNN97o5RO286Qu8KIkNlESGLr3o/KprqysYHl52bmYc7mcG48rlYrTSl7NeMELXuAWbsDQzQ4Mx80dO3Y4Q3T37t0jhIrV6et3KrHSOBcly2jI6nj/9a9/HcCwr9955534+Mc/fnEf+BkguPgDAgICAgICAgKuSGwbF7+6mVutFkqlkotOU3aCqz+NOrcMZhRrA4y6J2w0PH9LpVKOjmegiQrB7fVsJDAwygTpPvKUMNj9s/WzRv/rCjWbzbqVF1es1iWjZRsnI2Ddso4sg6Tsk617Ms29Xs8xp5YNsnURFTnNe9vUVfZvPV7rMx6POyalUCi4qO1qteplBrgaoOxSt7uZ2H5hYcGxTZlMxmM0tb3rd4Df3/i7strqBlbY92yzWNj3rp4OZQCiMldoOZQZ1awZ2h+U3Wc5lGEfJ7ex9Wojk20UvwaT0T1GlkmDx3Q8USmKBnIq+23d8kxdpfWgXiWbGcCOGfY7YDjW8PoTExNevUxOTnqRxBodr6mdoqLqyaq2223vnfB4QrMG6Dvf2NjwJEGxWMxz32s54/G4Y6JisZjHrip7G4vFvPaqWQmy2ax33tWKH/zBHwQA3Hfffchmsy4BfFSmDLY1BgPyeyvbATZZbA2KsuO53oP9hH3n6NGjLgVSLpfD/Py8ayflctm96z179nieg1wuh3Q6jfe///0XoXYuH2666SbPc1er1bBnzx7vGLbHarXqMZ/6v3qG+M50/AP8sVwlVjpO53I5HDp0CADw8MMP45FHHsEb3/hGANhSJvXbYdsYqGoIMe2RRqsRUROtDmTjouE5GKr7WCcI7VzZbNYZYbVaDYlEwhl0OtFxUoqaeO1kODU1NUJr24mXiHJp8bh8Pj8216nVs1kDgZ+npqacHILPZF2NaqDo5FetVp3xvrCw4Ll+1P2jk6saVywLQRmClVuwHFov6trk+9OIYRrutVoNtVrNHVssFq/oyH59V3yPzKGXSqU8zRIQvVOSzUZh34H2P94T2DQCtY9ZF1OUy5wLkqiydDodFAoFL7+pamO1/dkFii27zS6gfVrrzOqq7UJVJwYti+Y31Umdizq9n074dtcqfQatS/1bjTL2H3XV81lzuZzXLxV2sa2fu90uNjY2PEONn2OxGCYnJ90z6C5PlI3QYNUxqVgsunRSfBaet76+jlar5Rnk6mbevXu3Z1xqrlN99n6/74xNyhd0nFEwjzSfgQayjhVXK17/+te73aJarRZe8IIXRBIkUYSKxgTYdHGdTsf9biP6rSzFtnPW/7XXXut2Umw0GlhdXXXyg+PHj7t8zSsrKzh06JD37rvdLn7u537OlfN3f/d3v8OaunRgtiGbxWP37t1OpsitxXUMIuy7scfZHS/17yj7Q69NMmn37t1YXFx0Zd3O2HINKldayeRmHsG5uTmUSiW3MlA2jIZe1ERrJ0XLIFpmRrVsVnyserdz5865zj0zMxM5ofD6WhYtn17T5gXViY/3VAaV6Ha7Xt4z3lPZ3XGMsZbPGhaWaavVapGTORMzc7DSAC12JN6fW12yLqzBoBPvOJ0r60jrQo3lTqfjVuXWqNatUJX1vZKgbIRODIcOHXIDf6FQcEwEA/80oE8DzjQATRcT/f5wC1v2v2Qy6fRS3W4Xu3bt8vSqBPuufc/8bpx3IYqxt6yMPoOW017HBghYHSownFTHbUXKa0RNFrbcyozaMcNOFAodJ8b1SXsvjl16TTV21XDWZ6URqP1E9ffKeGt9Mn2NBlkpS6qGn6Z9IutqU1YBm3lQ2a6Wl5c9Jnl2dtYtrngPYHMc03dix0EbW0DoIlb1qNxamu/hK1/5Cq423HHHHXjwwQcBANdddx3279/vsfPaXtTwUaOU/VZ/U88EsNlH2PajmHz7t36u1WpoNBou9ZS20dtuuw3FYtGNN3bjGO0jH/zgB59lTT17vPnNb8axY8cADFM4tVotp7NVEooaf46R3E4d2KwLSzQBmwttze+t56mXhe80Kk84x8MoD1a328W3vvUtd//nPe95l51FDRrUgICAgICAgICAKxJb6uLXFTt1p8BwpyOmpwB8ZpBMGS3wqETaSnOrlkVhV+h0O/IeauHPzs46F/HS0pIrp5af19HrK9NrpQH8DhhNaG5ZIJsCxK6KohgZdVOyfvX+1o2jjOY4mcTi4iKmpqY85kzvoytxfUc2klcjD/k86vrT+0a5f1nm06dPOzZnMBi4lSRXsowUTSaTkS6q7Q6b2B7wNafAkBnSOtftX6172rpGtQ1oVH+z2XTvg+ysspjjXIXadtmOLKPK8/Q3vYbd6ELbN9kFbVd6/DhJgU3+3+/7O1IpuxSl6Y6KUubzaP9WljQq+taWy/ZfvTelB/o86oFRqPuWdaD9SPu2pp7T9kHtPs/TtFKamN+CdavR+bz+YDDwJEnZbNZj6ZPJpOfWJzY2NkaYZY050K1W7UYAWvZkMulYWabUuhrTz73sZS8DADz44INubtq1a9fIWK5sqsYERHkBtR3aPmc3moiSinFc1zFM21oul3Mu/pMnT7oNV9rtNvbt2+fmW45DypTzmnfccYfbiOBi4tZbb3V19fDDD+O1r32tNzbq/N/v9117S6fTrh/s378fg8HA0+Kql0PHCGvDsH54TetZ0b5gxxrVClv7SsedYrGIp59+GsCQ3X3zm98MAPjoRz/6nVbfRcWWGqjq2i6VSi6107lz55BKpby8hXxJzWbTcx9b97TdLYovr9lsevS4TY+jEgPrJu33+06/Ua/XXW45ps1QFydBY8C69fV+41wj+kx2oLYGqcW4ABYOJvp8Opioq49lIJjKiQEj+rzW5WpdRbyWlkUHKxrH1v3EstjzaGByn226CDX1TavVQj6fd4NDrVYbmdS3O2xaJObMKxaLWF5e9tJ62cmE562srDhjdufOnd41O52OCyrjhK+GDHWEzKcX5W6lkRu1ILR9wuo1x7n/rYtfwTashqBeUw0zu9CxBjEnGbbTKFcY+6nKJAgaelF5lW2ZrSGuRq3Wi60jrQ8F5UlW58oyWncuYdOCqT6U99MczzpWaptQw4+/85qNRmPEUOSx09PT3u5XxWLRc8ETTE9lF9G8ty6uWq2W06eyXLxmu912ZaAe0i6WrwY89thjAIbvkH3akiEqTwF8Q8cuyvgd4Ne7PZYLfzVQdSywRqXKYfTYXbt24ezZs+4a1WrVu2+pVPIMLxqvmsLpYqJQKLhr/6t/9a+QTqedREzHlVe96lU4ceIEvvzlLwMYGrY0ZNlGdbzVcUvrzaaHUkJFCSqrL+fiWdNMafl0nrUpxfbv34/l5WUAwNNPPz0SwLVdEFz8AQEBAQEBAQEB2wpbxqDS2ufKYXJy0jGYXOlyFZNMJp07wDIMeg1lFQBffG3ZPssKWWZDWQ3dyUQZlUqlgkKh4O1epTvuaEJ6y5Ly2rymHqdJ03UVpq5wHmtdBcp6KKL+VnZXXezqcuh2u96q0Arqx6X7oouSn9UdawX3WjdWwtBqtVz05/r6uvttfn7eq1/rzqjX6y5q8pFHHvF22Njubj4+C+u5WCw6BnVlZQWFQsG9r5mZGa/9A5tBLMvLy15wnzIoNhpUxfb1et2tynn9KLYT8KUCeg0bCKWI6gv6G88HRndVi5IK8Debkoq/RQUlWklJVN8nolKu8Zhx7d96dZQRUlYjSv7CMltPg0oNothalkfPs8FV1lOhKWq0r+sYQDeulYvoeZZpAzZZUbo9daxKpVLodrtOQhLlKeL4PzExMeIJi5I7MNpc+zq9P7z31cagXn/99S6gcceOHSNePeseJjQQx8rfAD/ThPWuaP/T69r2ej4Jkd6j1Wq5iPKHHnoIe/fudR6ipaUlb2e5YrHopFvZbPaiu/nvvPNOAHCeOc5xWp+si4ceeghPP/00br31Vlc2ZUlVHqTncd6K8oBYmZLO/7bPq+wKgGc3WLZVx1y+G2aC+eY3v4knnngCAPCWt7xlS4LPxmHLDFQaYnv37gUw1Lqpjkpz2mmkXDab9SLpdSIiba6Tjg7IasDaSVJpdRsFqy/aDoDWoIrS1dn7Wdee7vhgdUM6MejxLItO4HbytpOIDiTWYFdoCiCm/OI1tANZHaNKNmyUaK1W8wxt1aCqwbixseH+7na7nqZydnbWufPsAkffZTabRaFQcK6gvXv3uohR/m6jf7cbtH3s3LnTuVe5m5HuAkawjfHYRCLhjuv3+57+UBdv3NKT7YWuH3t9nmf7kL5zdU3xGMDXrtH4idLH8ZyoyFS70LG6bqsB12vbCXNcOdX1Zt1juj2rHWv0+Ch3qUYja71Tj8b7aH3xGP5mZQhRcggaCVH9gvWlUfxa51auo89vx0ud7HTRru2h0+m4namA4aTP56G7X/OkanYAlRGw7eozWiIC2DTqtU2wXI1GYyTbwJWOV7ziFfjqV7/qnpGLdmLcooxpythG6vX6WJlQ1LuPavP2N9sXtB8rGQMM3xsN0htvvHEk44rOZbqAzGQyaLVaeM1rXgNg2FfuvvvuC6y9aHBrYLYTllv7GCV+U1NT2Llzp7NN9Bltirpx8h9gVIOqYwLgL4I1fSLnRr2OXtuOqTq29Pt9txX01NSUy0yw3bZBvewGKiuAQVGcQKenp70gpW6364yReDzuiabPnj3rVjjUVfE4zdkG+Iyrshc6kFmDUVfonJSi8oqxEUY1Pq7W9Zlsqhk1pO02dIRqdZQpBXztCsum6Sl4bV7XTnZqQPI6mgZGyxB1TXuMZYx1QFRDWrWk586dw/r6upuolCXld3zXc3NzXiowqyvVd611D8AZvWfPnvWYo+3EpqoBl81mvaBBHaBLpVLkBM1BifU+OTnpBiHbRrVNq+EFDINb+BuDsKIWHpYlsGWxxq0agzaQx7KkOrhaxl0NNZuXUVk/O+jrs9vv9X76HnK53NhUUspCab+1fdIyqNqftF/oooljkPYnqz1j21WjmkaaTlrWI2GNSz6/3k/rxG4QoHVl9a7WkNG/G42Gm/hZfr0Wf9O2TugkbdkofrZptNTDk0wm0Wg0rjgt+vnw0EMPodPpuByiuqFNlLfBLubswg/YDJbTtmj7qvXiRcUW6Duxx1niSK9fKBS8HOhWb62L0uXlZczOzro95guFgtuY4Mknnxypr5tuusmV/6mnnnLf33nnna7tMb2azqM6xrVaLS9IuNvtYmZmxh2rY5X2dRucq8eqXUJSRkkhO/6pbtwuVK0nSBfzOt/3+303B+7YscPV1xNPPIHv//7vx2c+85mR+tsKBA1qQEBAQEBAQEDAtsJlZVAtE7N//36nE7S0trq0+v2+Y1qz2Syq1aqj1VdXV91vXEFytWJZE12pKJUdpWvlZ7stqR5v9WjK4ERFRqprzep+tAwa8Wjd2FYDajVVuhJTdsFKDKyuUFfXukq1rkUtj7pb1F3DZ+U1qI3R+tDnn52d9VbljLxNpVIjGkeu+sjWqkuFjBJXp1puZohoNBoj7Np2gbrFWq2WW5mfOXPGk7So607ZMkZ0ahJ/Ze01tYxtq9Vq1f09MTHh3K18/5a9529abqul0uOjmFd9P3qeZWus/k1dhVHaMJbNsi7W/aj1bplEfrasBMHjlEGNYm/JVilLOi5Nl/UKWPmOPne/3x/ZdpXH2X6hdatlU20f68v2dcDXCup7ATZZS9XTWRck21IyubmbUywWc3IVQt3v3W7Xi8DncZlMxmtnls2zO+/pu8xkMpcs8vtygprHJ554wu0KCPi6UttmtU2wXtQTp8drH2efs14B/T/K/W513xY65+i1qEnWNEtWmqOpyvT91+v180o4Hn744cjv1eu0sbGBcrns7lGv170E/MqmUramnig7j0cx2HYutPViPQRqN+iYE2Vj2DFWv9fztH9PTk568sZHHnlkbB1eblx2A5UVv7CwgOnp6chJje4AHRh1gsxmsy4twvLyspcGSXfP0Yl1YmJiZKK1ncdqsIDNyStK7Az47rVut+tpzcYFRLDzRrlY7f86kafTaS9VCgcAnVR0gNA61c/sONYdwefhs2g9EFoeNRDotlFjxbp09fnZWcrlMrrdrptsrdtRy6CdigOcagltvkqWJZfLuUCrUqmExcVFb5LeDrtM2aC1ubk5b09zTtD5fH6sy53vhufl8/kRt5o1voDhIKzpSnSgpVFr3yXvp31GF0HsX1FSEZZZXec6sLZaLS+tjDUs1QVvNdCqj7ZtP8rtzXJGGZp8Pp0QiSjN6zgtaTwe93bMUyNC+xPLw3Jo/9QxkL9FGfi8n/ZjnRhV9qILZfYX60rneUoS9Pt9bytoa9SoFEjHJ31O5ktlndpcq4PBwDPWNZWUuivtYl7rU7WErHO7LfCVCLqnKVXiQhYYnTfGpZLSxb5N6aZjtRqdPI9gm4+ShdgAIDtOaLvQvlQqlbwgTWpTrVEMDNvMTTfd5NriyZMn3ftl6q3z4SUveQkAYG1tzUkIJycnkU6n3UKm2Wxifn7e28lSxw6Vtelz8Pl07tLn17qxhqQSSzq+USusxiTLofdU6KJZbR8t69TUlPcOo1L8bRW2NA9qQEBAQEBAwIWhWCx6hoWNHCe4SNEFqSUyovSJXNxbfbgaXpYA0oWmQkkYJR54rtVNE5qlxHq7NGB6eXkZTz31FE6cOOHqhmz97bffjvvvv39sPd5xxx2uPkqlkvucSCSwsrLiDN3l5WVXz7y/xi0wvkKfmXWhxrvNMGID0XTxqSyt1hP/p0FuF+F20W8JMuvp0vbA51tdXUU2m8Xi4uLYurucuCwGqrpzWRG7d+/2mMko0bRlHIHNilYXANkxroToRqIcgJ81HZS6OJTdiIJ1ZSu7alkiwq5m7EpTGRLrTrNuK2V+EomEi4yl6yHKbcvzbL1qeaLKxsEoyv3CY6MYLOt+tey0sm32s3Y0LbMOsLb8Vm6g14zabIBtYnJyEq1Wy6UqiUqxshVQdj8ej2N+ft4NtolEwpXftlEV1DebTTQaDSdnUCaN9cM2Z0X6ul+9Mlds3zq4WckJoV4JDo5RrizLhOpxlNNEMSY8Lsolzu/1eS0zycFbB33LLGkdKRur99Bj9fmVnVSGT8+x0eg6cathwGCmqDGJ9x0nWbAGh30+7aN2/NU6sy5BIplMemOm7jymYwLdsTa9EDBsY8pw9no9z2Ogf+t4pOMhrxklzwA2I7JZzlar5WUKCAgI2N645AaqDtCdTsfLL2bd19YdoIO1DkK6bVg6ncaOHTsADKOdNbVJPp93BmoymfSiHAF/IrYTn7q3otyUWr5xRpl10+n3ej91kXEw5vNqCipGNKp+dpyLtVareRHITInFutBjo7S66o61Zdd3Zo1Xnei0Pq2xaTWvUfWqk7wti40ejtLiat1Ql1mtVjExMeHpY7fDNqhqKOmuZfxNNdiqty0UCu7ZV1dXkUgkPCNRI70BRC48bE5Q1QXSjRX1DvjetN/qYklhFxrar9Q4phtL3X88Tg14Xmec9EYXoKrb0vP12ChXGJmMqAweNNTH9X29l23bOh7qe19ZWRkxNtkmbSRvNpt144A1ZFVuYMs4ziXL79UQVImVXkfri/WgUht9BtW56nvlYlv7uF7f/mYZpagFLa9vd7DiedpXrlRo3bNeot61lVcokWK3wVTphx33bXu044fVYGpfOV/GBH0O5rblNTTNIdlCXpeSEmA437daLbz4xS8GMBzbyaZOTExg9+7dOHXqlLsn834mEgk0Gg03Nubzefe8u3fvdluAAsOsL5OTk569oXlcreZZ7RT7t74LbcM6l3EMtQtvwsqY+B3vQbAfKsmnsTk6DrRaLTe/LC4ubqtMF5fcQNUBpFwuOwO13W67FDiEVqZdtbPSovRgmhJCByHVfnGSspMcrxXFggKjRphlRGzuMTtYRnV0siPsfNaFEo/7mjidvNvttjuW2kF7bV6zWq1G6nGZSmQcy2InYjU6tT6sG8jq8PiZuWOt7k/rU3VkWp+WKbLGpJ6nLit9L7VabWQwYrtYWlpy7+j666/Ho48+iq0Cy7Fjxw60221Pn0cNnnUdqfHVarWwa9euSGOJE8Y4Q0oXF8pecZCNGmitMaJGTafT8bbFjWLp+U6sntJeV8uoxor2W57PBakdnNXQ1UGfbVG1lzpGWGNO+5Omj9Lxw+pmO52O96x8X81mE+122y0ItG5paKmRRcOr1+t5i1jN7xmLxbwJjpplYDNoThcRls1Wo14nwHFph7gQVQ+QjtF6HZtjVscBHYOY4zhqTFaXKe9nFzO8jjK0zKs9zkt2pUC3BWX710WFJQiiCAOr41VdtPV8RJFG9jerAQZGAwD1vUTNFQRd33aeVq8Gx6aJiQlcc801ePWrXw0AOHbsGP7+7/8ewHAcmJycxOHDhwHA0zfzOoxj0Xaxe/du9Pt9FyiUSqW89m1JEiWQrF2gbnw7/lgDXetZx+moxW1UABzBhQbHPy2T9XYqQcH3XygUsLGxgYMHDwKAl45rK7D9QpkDAgICAgICAgKe07gsy8kobdvExATW1tYci2oj3ixryb+5ulDXl+oIxyXOB3zxdRSNrau9cSkZrFvZuuz0OOve0xWPskaW0bH3sGW2KYD0nrpKs2ynXenpStgyE1H6QK68LaOqZbVpjfTeyuQQvL6yOrrq1P/1vVv3tK5WuSJkPZ0+fRqpVArAkI1cXl72ysZrnD59GluFZDLpEvNnMhmsra259pHP50dS/ej7Wl1ddecpA6jvw4rm+Tsw6iLqdrtO88q2GRVdbr0EWj62fZUWWFYxiom1EaZWD6s7m+k4wchfKw8B4L5X7wLbIBODE9qHo9hcsoPVatXzguj1u92u08VvbGyg0Wg4ljQWi7lnGgwGXrJ6K7XJZDIeA8ixstVqeemW7A5JjUbD0+Ure5RIJNx5qVTKebSy2SxyuZxrg9qf+P4IlU+wvWlbUnfxuPFYxyqta362zLVl/qPYcO6QFOVRYju1TNqVBp1PmH7LuvX12Ci9OOcN67UDRlNQ6W/EOPe/zUZj5w7br6ICgaz7msdx96aVlRX322AwwAte8AK87GUvAwCcOnXKBfd87nOfw8zMjPt7enrajZPFYhEbGxvuOqVSCa961asADBnUL37xi84Ts7q6isXFRbcVa7Va9Tagsd5HHZtYfr4rK/fT+udnep3U+8Mxh8dZKVOUjWHlF9Y20P6tNka5XMaZM2fc+L/VuKQGqrp2gKFRykbCXW50cOZLqlarnptKBzJ2QtWM6Yu2E5o1cgidvPVetuzsTOM6up3AogwAXk87ve6bayd5Pde68PP5/Ij7NeqeOuHzOnp/7VhaF+wc2gn12ez9rLEbpfGyLg2ts6iyWL1dlHu52+16mlEN1IjH42i32y79WCqV8napyWQyTpOqi5Gt1KfRBQ8M23Sv1/MMECsdYZ9QLVWxWByZTPhM1G1HSQM4kGkQka0LXbjwGjSudAHF85iOJSqfb1Sf4fVteip9nqmpKW9f7Fqt5rnLo/SxwObCivVrjVLtm6rjs8++uLjogutWV1e9RYTtv/x7YmICxWLRvaNEIuHaXiaTQSqVcmNiv993Uo58Pu8ZkFoX2if4rDp+2e0KtZ+oO6/RaODcuXPue2uwzs3NAYDbuUzd+Non9XmjAuV04a86OFtnVmaghte4McGSALoVtcpBGo0Gut3ueSO7rwRwZzdgOK5tbGxEuoM5nmgdss+ptEyPBzYJIDuPRrn4+b2+X51/dO7SMV7lIvzNPoMaU51Ox7mcH374YWcsAsN0lfv27QPgx2tQG8oFCcd9vQf73Z49e5yOlQGa7EPMg0pQIwuMuvGz2eyIZlbbtPZf67onzpeeT+sDGA1Q7Pf7nlwomUw6Ykb7bzI5jE3gInYwGLh64mKX5x08eHBL3fxXtiAnICAgICDgOYKJiQkvuCWTyUTmE+ZCS1MS0dA5e/asx5DpYqFUKo0w4kpKWNJH94bXBTANRCU6dHFSr9c9I4wLRsZVaNyFPtett97qGWgnTpzwsj/8/M//PIChoZXJZJxHYN++fe7+7XYbZ8+e9YiJV77ylQCGW23/5//8n129zc/Pe2moOp2Ot0EG4Adv6sJNDW/9bBcAaoAy0EzP0/gbJSeazSZWV1fd4ncwGLj8re1226uXbrc7olNnG+j1ei6wnL/zN82xuxW4pAYqrX19obTwld3ib7rytStm/awdShkVy5JGudg11ZReQ92DgL97inWVKzupLIBlhsalSer3+yiVSl6gj7IGutrUlde5c+c81zVXWnoPywozIlxdLFYIbyNBeW2WTevTpsWxYngdoPgbV9ZR8grLbukqmK49ZU11taj/awRprVbD2tqae9fT09Mj8g12ZL0mMBzIjh07hsuNcrnsBgN6EFiX+Xzea5u5XM6t/jc2NrwUWrataqJuyzorlFlbXV31GEa6u3hNZTuVRdy5c6f3HnUg1jru9/sjrJ6yc5ax0c+ag1Dd4+pRAYbvmc9A9zWZDzIPfAaFuhtrtRpOnTrl8i6ScWBdLywseBIGXr9er3vtSz00g8HAlZmTOFnLbrfr9Vd1x1t3eFR+SJ6nHhmtdzKfnLTX19fdxLS6uuoZNMrMr62tYe/evR4bpOPquDGBZY3y1JDhtgEfhLJM1o2qQT3KpvJ7jYDmNRjxHRAQcOXgkrv4O52OpzOlQcoJdnZ2duQ860rWXI8cdKwuEtg0CO0qkudZd7UO8DqZ6gBs3VT2+mpQ2UnEulv1b3U/6e4U3D7NakmAoXtHd31gkmYtT9SEpnVqy2HLZifsKJdOlHaWGhdriACbbhNCJ1cuTPQeqjONMkr5rFYjRw3pxsaGk5DYumDZaQxmMhksLy+7Z1dj7HJCt3uly8WuxoHNd0CjtFKpOEPcave0PdAAtbITntftdt2ikVo+YKj5KpfLkQuPZDKJ9fV1d+zk5OSIzpSwi0dei7/ZBaEaHepi63a77h3p4g0YvmfWhbrxo1zH/EyXsGZC4ALlyJEjADbTfl1zzTXO7V0oFJBMJp02TiU7pVLJ2wFHc7BqvfCeuvMX3yv1xLqoVSNc3YB2nLMaNi5o2+02UqmU+20wGLi6zWazaDQarm57vZ47bnFxEd1uF9dcc42rP7sw0HdAcEzQY3UhYsfZcRlRrARDsyfoeJFMJkfczFHu7ysZukBjnm/Nl8z5dX19HbFYzC1AtB0MBgOvveZyOSctAeAWMDxPUxvqYqRer3uLRtVh27nR5msuFosj6ZoAuBRTusjodDquL1q74L//9/+OT3ziE+45fvVXfxXAsM+urq66SH1lJbmbIOvqoYcewu///u8DAH75l38ZX//6110dkKHUetV+rmkKbaYLS/woYaPjqC7AOOdF6Uq5AFN2WV31Ok+k02lsbGy496oL42Qy6S22gU2isF6vI5PJuHRdHMe2CpecQZ2amnLUshoqnNwoYi6Xy972fCrwtayNsqbKDOqgxetYIzQq2IP3sg0s6jzL5trUKdoIxunSeC/V2fDZba4+nWzIjqnhYicnHSC0U1jDU+vJsjE6IdigJa1Tq53RZ1SGzLo7lOGx71oNUj57lMaWZeZ9db/62dnZkYABXXzouy4Wi5iennbXiMfjuP766wHgsqeconehVqthbm7O06Bqu6rX626Bp8EGNo2Utg0uAvSdqxbx9OnT3jvj+1haWkK73XYDvdY5NaCaQomw2mJrTFrW3/ZnXfQpg84FnEUul0OxWPQmUts+1cjXiVP/f/zxx52ReMMNN6BcLjtmNB6PO4bxqaee8tJH9Xo9x4Qmk0kvWEKZbNXNZrNZpFIpt5CyBhvPJ/gbU/RpvYwba4DNyceym7FYzF3fpqDShdpgMEC1WsXx48cBDPNJni9Xoi2z9dwAm5ph9YDYRYyOJTp22AUuz1tfX/dYZ8A3mO++++6xZb5ScPr0aezatQsAPGOFf+uCRrX2uuClIcv2o4tMLoy0fnW+tdA507ZD/axGKA3ZKN2yXWB1u12cOXPGi7tQQzuZTOLd7343gGEfpq3Rbre99lyv111dFYtFFAoFZ7ynUin87u/+LoChS19z96ZSKZc7Gxh6MXVXKR0Prb1hDe0oLwHg6/pZP/pZF6LW6N3Y2PDGf80BPDEx4elulQRQfa5KA3gPkhVb7XUIaaYCAgICAgICAgK2FS4JgzpuNZTJZLyVdyKRcL9VKhUvFdDMzMxI9BqwGQmqjBs/c1/cce54y6JYd7i6onQFZ5k4XbWOi/DnuVGMKsus97MR9uNWVMDmqob3jhJUW/2musxs3VjXcNSKlte0ejDr0rPuCH6vOltdoau7hufp9W396nnq3pqZmRnZDSUqVRbvodeZn58HMGRQu92uYye3CqobBXyBe7/fR7vddvVl24qypupyUzc2MKxvei/ogtb7Kpu6trbmXGzKvKjOEvB3IWEda1nGpUDT9wH4rL0y4ZZlU51puVweSUenwQTW82Bdy7zH85//fE+vXK1WXT2phrdQKCCdTru2ozpTu7uMuvVVPsFE8rym3SHOyiQI7rKkY6LWiXpnlDlh2fh8ZKC0TsjmZjIZj0VRd3GtVnMsEhl8K2/ib9ZdqWOuZbkVOpbYXY+srpVMUCaTQbvd9lJJ6TWuBmSzWS/YRXd9KxaLXl1oQFUmk3F1wWN4rqZB4zxl+6POowSPs3MssDnGWxYV2PS0RGX/4JjB9tTr9TAYDBwT3G63vXFMx8JcLufu12g0nAyH99QI95mZGSftymQyrk6WlpZcdg3Wkcof6Gng86p9oEwzn0E9IXZO1XlN5yw7P9m+xfNKpRLOnj3rWNJsNjviTbP6bGDYTtRrqePV8vIyBoMB9u7d635j5pKtwCV18cfjw52cWGnFYtFz2+vApZVUr9dx+vRppxPUVCs0WOwL5TWs4afX1xdv3RhqxGgHiTLe9PlsQ9PfrIGtg6VeR8vFe1vNlf5mJz+9rh0s7EDDskQZ73qPKHc8n1WNfKsTi6pfqynWSYtlUZ2gSia0Y2mdcltP7jsPbA54U1NTI5pbq//RgYsDUKFQcLqey4U777zTlV0nEP4DhgaJ7iJUqVSchsi69HXwVuOemrOTJ08CgNuNBhjWieqgO53OyP7p+r7sokvdtjRqmApG27G+fy27GlT2/fB3nqd9UzWnXOTZ96xliNIjclHHa1YqFfcMNKp10Kc2jX1Y+75q8Wygn8qXWDeFQgEnTpxwk2+z2fSMVLtLkk0vo+5Du4iLkivxOKs9Jrrdrru/jjMnTpxwcgS9Hsuofd2m/dO+Zxf++r9dmNr7EHwG3k/7KtOzqQHHfvOZz3xm5FpXEthmlITo9XqYmZlx86M1iHQc1/dgdY6pVMqNJ6lUCun0Zs5qusl1fLZGjyU+gNE+rn2A5Yp6v+wvPG9jY8OLKs9kMm6hRJ0054DBYOCleRoMBq481rC0QdFqEFupgs0JalPu2XgP3sNqpdWmUFKq0+m4cm9sbKBer7v3rWnTOOZofWvw39TU1Mj8a20SQlOA6edkMolKpeKIn1arhZtuugnAMMXX5UZIMxUQEBAQELCNQQNiYmLCGX8TExOYmZnxDEZdGKnRbz0r+psuRqanp0cICsAnZjQvcKVScYsEa/RqPEW9Xndafy4sNT2TNWxV262ZGdRwZA5YMqPcAhgYLk643SkwJLno8dnY2HCLLv5NMGe2kgXJ5GZKpna7PULeqMHKoMlGo4FMJjOiH1Wop0aPO3v2rLu/ZgrRd0nkcjm3CIvyYqiHhOACRduKPsPc3JxjTR999NGRvMqXE5ecQS0UCi7yVVdGgL9q1sY8NTWFc+fOebu2sHEXi0VvhdPv971Vkw2gUuh5+sKiVhq66tfjlXnk3+rKto03qlFxZaUUu7Kb+pt1TUelYVKm0j6HNjwreFb2QjvLONkCy6V1My4oS1fIrVZrJJLfBs9oIJatXw2k4HETExOOKSW0IylTpM/Kd6vviQNHqVRCrVYbYXIuJVgu7ioEbOY65MCpjGKr1cL6+rq3+lU2SoMRrDu80Wi4YKdisegGXQ7KXDVrWiuyFXoPHawbjYY30ZB9nJiYwPT0tMfoRGWI4P9W8qETo7qMlcEtlUqerEOZeO3rVm7A7/j/kSNH8OSTTwIYstUcr8hIaF/QdmsDOfXalp2IGpMYCEcJgTLgli2xDIyVLOi9tb1EeUe0XehErB4YZerIUvOayriwLFFpn+ji13Jal76tN56nXibty1Y2Zo0yzeeYSCTw+c9/HgEBAVcmLrqBqgM601dEDYgcdKI0at1uFzt27HADVLVadXqRdruNnTt3egmAdQBUI8y6wwF4A2vUBMlj1B1i3ZoW1l2t99ZB17oU9FhrTEXVSzzuJ/jlRBg1AUQZWdaQVdcpo3ZpTKrBSsOPCYutqyLqfjp5c4Wo6WCiDFlbn1zZcdLhrjzAZlSouqL0Oa2OR1OcaLvTyc/uxLSwsODSbVwqqBFKqYEap/yNbmBG3HLBtr6+7t5dsVgcawzF43F3Do+17ZbvWY0auvisS5ZYX1/3IuH5HriNcZQ8RqNr9f76Wdt0lAst6jw7htgcvlHpxu677z4cOXIEBw4cADDMg6vPSsOX99C2r/223+97W/jqgsvqtLVMfGe8n7JR2je0XUZlFFEjXtPpWUMd8CNzta2MWyhOTk5iY2PDi6TW88ZJmTgG6t9Raej0eVkm1URrSiSbyUTvB/h6QmsMX8kgMzc3N+eNgXYhrvOiQtuFlbx0OpupIK2On65x3QyA4w3bqiUGeJzdBVFT1+mWpbt37x6Zf7X/q4vd6ovj8birG5URsP1SylStVj05iGrpdZHfbrfR7Xa9rDqq9S0UCt4zqRxqaWnJ2SkzMzOeJtSSYFr//X7fSYfYnvlMtVrNyyJSLBYdCcByauyALuqUGbW2g46ddmFKYpG/jRt/LwcuCYOqrOj09LSj1jU1CuAbJzogWhZvbm7OTd6VSgWLi4suP5dOIGT0VC+j2ijbuPSl2OsQZCT0OkSn08Hi4qI3qahOa2ZmxjUeqw/h8YQapDqYWKOX99VniNJ4EVaHx3t1u12cOnUKAHD8+HFvVwkK7IGhy4fHVSoVHDx4cKQxaz1pvdmyROnU+v1hAJWmteEgVK/XkUwmXVqVUqnksalR27ACo7o4rQcLLWepVEKxWHSM/OXQoqr435ZJtwnlcXx3arSzvNwzmvWgBiPddjopRenT+Pf5jCod9HS3EgBu68xkMomnn37aBRFNTk567KPC6lOj/lboxKIeCm3/1oWlTHmr1cJDDz0EAHjiiSdwyy23eEE/rBeb39DWybitJO0CzOrirZdl9+7dAIZaT23TusizsKnZeD+76Yj2AU7oWjZrtGi74nH5fB5TU1OeG1LrQc/X78jk6ru0zL/1bBCqw1fjn88WdV/OJ//wD/8QWWdXMnSxwHZsN/DQxRA9bWqwcDxjW41aZOj74bHal5rNphtvisXiyJbB2l6j+jmvXa1WvfmHnh3A77vUoKonR59Jxypl3dlGOIYmEglv0a/9MGrRpJ4tO3dYuQGvUyqVnPfFGoh2MW2vo3aQPofmuWaMhL5z3aHLel71GW3d2LRf1uPCRceePXu8ervcuHqWmAEBAQEBAQEBAVcFLomLnysE0tbqGrKR5FH6Lus+a7VaLlJxYmICS0tLbvtBTbROliRKIwn4CZ/VHadlJ3SFrit/pd9PnDgxwnwxgrRareL06dPOdahuTa5Ko+5nVzdRmkzrNtW/1TWirJHWr7q0gc0tHYHhqrxcLrsV5OLiosum0Gw2ceTIEbdK1G1jWQ6tN10xW1cfV+H9ft+rw2q16qXk2blzp3tXek3LhqubgvfT9qOMkta9ukV4v8u1m9Sdd97pmIGNjQ1v21+tE7KnwPDZ9u7d61xy/X7f9YVKpYJCoeBFYhN8Zuve42dlH7UPUwpij9fjWO6ZmRnX1xuNBs6dO+dcbIPBwGNb9Jq27/M7/j9Ot1itVj0WQtlPy9xpnzp16hSefvppAJtBCLZv8DxlJZTp4JigbCivUavVvET9lJXwGuop0jFqfn7eBSfYFHL6PDw3KtNHVCYT1YnbZ1SWHtgcq5n03h6vx7KcVlqhx1nNvNbtOAbVamtVHsT2qO2FfWViYuKqZE/JsAObm5AAmx6EKEbbSiEsS2o9XOqK1nbONqhtW7cM1uvYdqftVdndZDLpeRR1nmK7t6yl3RqZ1+n3+87bp8FNtj50nNKtwnltm9pPd1bSuSSbzTommmmsNM4jKhuNloP3VQ+ZxlKw3mg7VSoVZ2/o1s28pj6v1Zfb+Vjvr89qbQP14OrzbgUuiYHKh+PWh1GaRTs4W5esVq6NKpybm/PcmrpNmnVz285s84Ty3tqBrO6s399M8dLtdp3OZHJy0jN0mQKIv+3Zswf79u1z99DOy2fU+/N+tj7V0NKJlga//m4DC3hP7eCsW8okisWiMzqPHTvmpZlot9tuAMjn855RAMCbiHVgs8+k51SrVU/CoK7ScrnsDC5OttpmCD6nDpbqvrb6S72GDip6LPVEVq98qUADgOWncdDr9TzDSI3XvXv3YmZmxhvA6R47efIkFhcXR9w1hOq11XAAfFd6VMQu669erzsJxtramvfOWUZg873ozi40tDOZDBqNhlt06kDLdqRjhr7rWq3mjlVDkH1fB2SdEHXnrCeffNLVdS6X89oYA3QIXVhaQ9pqvsZJJmy7tSl7tK1yYqrVat7uN3ay0XGBvwO+fp1lJuhy160xdVxtt9seEWBlQdpnrWGr99N+b2UR48ZAHY+5MIh6jn6/73JjAv7i7Wo0TgFfSqPzDce1qMUjMDom6jE2XsFuC2wXRlbPzfN0fte2xven71SPs9pGNZD0vFar5aWLspHwmp5P+4QN6E0mk05qyO/VvlCNv24LyzIqeaD1axep/N0u3PQZaSOpTaFjhxqI5XJ5pD9pPeri18qhovqffrYLT9anykF27dqFo0ePjhx3uRDSTAUEBAQEBGxTJBIJ561SltKSGpZlt8yf9cop88dr0qi1OXY1wHVcFhmbE9VqsbW81kCziy9dgKdSKc8bojrzqPLo86vxqAusVCrlFjkavGeN8VKp5AxbQo1JfX7rodLMRNZbrNexRqXdzlTriTpU/qYeY128R2m89T42wFP1+WrYKqlx55134itf+QouJy6JgcqHZcNSVk1fiu0I40TVVsQcj8c99zLdYuvr65iYmIhcUdpOqR2ELySKDudqgowtsMkKdbtdZDIZb0eKa6+9FsDQLaOrXbuatW5mQoXVLL8+v6aQ4bFRDA9/U9cN76dJw4FhYz106BAA4MCBA6hUKs7lv7y87JI45/N5b3WlnUXfY9QzMQWYfUbrQs7lctixYweA4b7T2kbUFaTPxWfQKHHLYtndP/S9aDJ6PiPfE98n0xBdTFSrVW+w1bosFApeSq39+/cDwAjTAfii+/X1dZeMX1fbZE2UNSA4KbCNM+E5MGRvW62Wx46yPSSTSS8aVoPYmICe8pter+fSKa2vr6NcLnssGPsU+2VUOW3AUDy+mdbKZhtgnbKcALznI2Ni2Z+oerESG15T2RNlB62hYMdALZ/2b2WVLLNq62RlZcXdTwOY9Lq8nxoGnU7H259b2dxMJhPpyWHqPn4XFYQVlRmAdRI1/vN6URsKNJtNtNtt1+Y6nY4XJAjAazuXe9IMCAi4PLhoBqoOvHbQ10GMLixgdBC22iadeG0UrE5EHHDr9TpOnTrlNDtR+pwotxhdgFEGY71e91YUOrnNzs6i0+k4F+vCwsKIC0bvbVczdiWsdaKGFzUo/Ntu7WbzI/Kaen9daXHFqO4SddnNzc05978aq9VqFQsLC+7+dktFfWZ1m2od8DeCRqe2EU1boq5ZvQZdOFEuQj67Gu42fY7eT98D06No+S4Fnve85408G+8bi8U8Y1IXY/qMBD9rDkgAeOqpp1wGBEb8RrmD1tbWsLGx4Z61Xq8742B6etpL7aI602w2i3K57BkgqqdUPWyz2fS2Vu31eq4fnTt3zn2OxWJeBLm2Y74bZTv4Xp966ins3r3bS5VFsK9Rd9poNNw4NDU15RmUyjrYCHfbLqwxqyyEXSRrux3nzuQzApsaa9Ve6sJJs02sra2NSBai0kXxfjouW5ewGtlqvNptWKNcxiy/Mj42hY3WF//xGnbXMv4di8Xc8w0GA3S7XddernbjdN++fVhbW8M111wDwJfmkPmy8xh/U9LFMqv6nnSs5txjFxV2buZv42RQdkGm+nBguFiM+o335pjDxatmEdDsGhqfYtu21kc6nXbEx8mTJ9FsNt1cvbi46K7JNF4cqyYmJtBqtdxCuFAouN9KpZJHEFlNr8pUdGy0WYqUeLHzjSXxtE55vL5HTR1p+6W+G8oMot6N1f8ePHjQ1dPlxkUzULVBkjXRFAj8WwdZO3haZs/qQAl9ialUyqVBmJiYGGEtdXJRI1EnEU7W2pg0dZSWpdPpuKChZDLppZbgdYHR3KbWOLZuBW2sarjoXuCcvHlssVj09sktlUpu4rWGlq1POymrllQb7Pz8vEsV9Nhjj6FSqbhj6/W6t+CwRqJlUfmbsoP2f12MFAoFz4C0W8RG6ZoIq0G0LqSoSdoaJvV6fWSbu4uFer2OQqHgLa70uQaDgdP+Tk5OeoPyON3x5OSkE/UDw/qifqhUKmFyctJjz3Q7u/X1dY815WC+c+dOJJNJt0jSvJ1nzpwZ0a7yOdrtNubm5ry65UQyGAywf/9+937q9borV7Vaxfr6ujOItB1Rb6WLRX3/jz32mKunffv2uTJz62Qyv7lczpsQdIKwri/A9+TwfnSLWq8PMGr42YBBQidpnqeaU37Hsuj4lE5vbvPKMYPPqi5bhdWk6QRuDRFlN3O5XOQCkM9mDc9x/V7LyfLZxPoAvO1KeawypqlUCl/60pdGnu9qBN3E6gGyi4pxOzLZOdMaNjYwCth0v+vihOn+7D1YFjuv8Dz9Xz1vUeXhZ6bGY1/tdrs4ffq08/Dt2LHD/VYsFnHq1CkXs8AYAmBT0sBgUgDOs8S6YdnK5bIzwhYWFlAul3HkyBEAw5iMbrfrxsPp6Wmsrq4CGHpu1Wukc3q1WnVufn2XUfWmeWb5vfYZNd7r9bonVdD+pt7VWq3mBTtHzffWmNf3FWUwr6ys4PDhw3jsscdwuRD/9ocEBAQEBAQEBAQEXD5cNAZVrXOyIWQprV4JGHXvR1Hb49zV8XjcrbwZSQwMXe6avkJXm2RerCZPj1NWwJZFaXQil8thamrKWxmq69gyC9bNx79XVla83XjUzUidH+vo2muvdb+dPn3aY8zOnDnjmIdYLIZsNuvehWpAc7mc2yGEf2v9a91oXRw6dAjFYtExU81mc4QhV02lrrSV1VHWrd8f7mCkyYat20hF4VrGfr/v/WZX9prtwGrftL0qrN4yiiG7GEgkEs5dCQxZI3oCJicnvd2idEcS9gXLAALD9r+2tua1R7L7rVbLc6Vb7NixwzGcljm27h9lO7rdLnbu3AlgyBCwbSQSCc9jks1mXVmYKYHX2bVrl2MD6/U6HnvsMffuDh486LmrVadp20o2m8VTTz0FYJhKShlh7unN83R3FqttjdJRA34aFzKV1v1m68vKXGxWBT1GWVIyWcpw8rdisegxUgyeIFQjqvdknauLVutSYd28WkfK2pOVVcZH72dZU7L2g8EAvV7PeQn6/b57X2zTHMt4LADcf//9eC6A2V+YRSUqoMXKwfQ3MuDqTdHdoFQqoCwcx+ooVpSwXkn1Kuh8wO+AzaAhIplMOomKdVOToeRxu3btcu1kz549I7tkke1UVrTVaqFYLLo2RKkey6LzU6FQ8FIZdjod3HTTTQCG46eOpzp3nDhxYkS2OC6Jv3qIKd1h3VSrVa/+1XthZURRzKfWOY9fWVlxDC/r2Maq6Jxo5Y067tCbpXPq5cJF16Bms1k3mFhthuqROAhpZY8z4Ky7otPpeOlrOLlxK07eVwdOdVuwbLbj2XsQlUrFez7d2UEHgWaz6Rpzs9lEo9HwyqnGo27Jt76+7jortWxsdPl83l2flD3LmUgksLa25hmTWr/6/JOTk84A4W5RmutMc2fqQKLa1W63i2Kx6Drh0tKSN9lZ2YQ1Cm198/rqjldQJ8u6qFQqrsNxkFGNjdW4WlctYV1RtsPTULRau4uJXq/nGcdqrA4GA891F6VX1mdVF/iOHTucEaB60WRymF6Ihme323WD/t69e72gwXg8juPHjwMYtvFisehkHp1OxwvqUX2YLgCnpqa8iYcyAj47sGmg6fNRX8l71Go1b5cnbVfWRc3vAHi5bOPxOBKJhOf+s1pV7d9WUqCwiyD9ns+qkw/fo9XQE2po2vFJr6OaOS5YWDZNR2Vd7MDm+Mt2ob9r+240Gu63TCbjpejT966TJCe6qLzNtVrNC7jLZDKuDXARqC5+m69V55F77rkHzyVwIXn06FHcfPPNngGn/TRqHgNGU5ipkcUFjE2nxmsC0ekOea72DztG6lhq9cc6T3a7XdcHVTpTqVSwsbGBhYUFAMO2WyqVXPkWFhbcNSuVCg4fPuzmhNOnT7t7MLhT+wzHBNYbpXp79+511y+VSojHN/XQzAPK8uh8kcvlcOrUKadPLRaL7jnK5TKazaY3/lGOV6/XvXnUSnt0jNP+TANVySum1+M9WO5isYh6vT4i0+Ez2QW1vie1vfT+uVzOi4e5HLjoUfw6+avGgn+rgN4Gy/B41VK4gpprasQzJ1q+QDWMbNSzThR2kNVBnmWmToz3185CllRXgjpRWkOPZWC2AWppJicn3SRMQ1obJZOdcx92DS5JJpNuMwDLZDQaDddAlRHLZrNIpVLeu2DD4wqaZS+VSh4LpgNLOp12AVTNZtNLOBw1wNlOwO+V/bQMj2qlzpw546W/sAYCn5UDp+5vrFpFXb1bfbNGDAPjt0h9tmBdnj17Fjt27BjR2vGebNOAz57Zv5Wp5zMyaC8Wi7n33mg0vMUb95xmWaanp917LRaLrm5PnjyJfr/v2E+7io7FYt7qn8fROFVmnIb/3NzcCIPB9tftdrGwsODawNra2ohhpIO5aq7q9bqbdA4dOuQFAdBg5nl2gI7SJLO9Wy8MoW1OzysUCl5aHjXarH6T3+n/hJZT+ySZaR3ndD9yvbb2ES4motp0p9PxMqDoOG4XnxrcZHW4DIDj/TRDg+ZUBjAyBml2FACuvTz88MMj5b3awWdXjyTga+s5DthFFqELeEW/3/fmJusVVCOYCxNtmxrMpnEByqCyb+jiVdsoPQT2t06ng3379nlzjNVW6v3K5bIb5zWA99ixY1heXnYLoV6v53SkGxsbmJ+f9zacsRp07Ztzc3Nef+Zv9AQxi0+9Xvc2WNGNAzgf8/tGo+Ft76p6fLv4tTaT3VRAxye7MY8SG1HzL2G9yTpfsq3k83msrKy4+r4cxmrQoAYEBAQEBAQEBGwrXHQNqnWR2ihOqzvV8zSSTV0JSrmTqSALUCwWHWvClRbvZyNylQnl78CmllJXmzbCnqsGdUfG43Evl6WuPLm6oluh293c2iwWi3m6kcnJSbdiZNQmr1kqlRy7urS0hHq97jQ6xWIR+/fv9yKerRaM16nX6y6KkbIEZWT0szJf6ibhrhbKYNH1sbS05LkhlfmMetf2+6jIbNYpn2l6etqxyTt37vRcylGpSqJ0PFxVKjOlbhjFOFnCdwJtOzYhM+suFouhWCx6WS2UxbZuPqvD5LtTLenGxoYn89B+0ul0sLy87K65b98+9+7YZmzdAkPdar1ed66rZrPppABs72yb2qaZP1f1w7x3o9FwKap4TTIPNhOHfceacSKKCbVbJGrdqvcmyt3Pe0RJLHgdHXf0GjbbAMH7apvQbA1WwmIlMwrt96pF1OfheKxlsRKZKEkUtx9WVsdKeXSsVmmXjqnKkvZ6PVSrVS8Di543GAyek8wpoZlirFzKegittpCwTJiOv3YuVk27jpd0U0dp9qlFttpWXu/hhx92LOH09DQeeeQRAEMPSiqVwunTp0fuf80113jsLl3Mqhm3/Uw9MSrdoawJgMdCT09P48CBA+6ahULBe2att3q97nmsbPuen593WQTOnTvnyfZSqZSXoo+ZCGq1mpdKz74HOybY96YeC2VU9Txg2IetdJDlVjvGyut0blAJUj6fRzabdUz05WBQL0mQlAb16ORjqWutfJvr0hoGWmGVSsXl6dQ92a1AWQ0k6+IHRt1gOqirgVoul0caMDDUwHS7/g4VmqKICad5TRoLqVQK+XzeGZ7qUuUzaTlZrgMHDniGrt3eUcvGwUv1MrxOpVJBq9Xygiw0mEblCOpyP3XqlKez0Y4yPz/vJQ/P5XLu+gzqiNLVWE2eGl8sm16TxvmJEycwOTnpnimTyXiuZjvQqitW68suoHQbRetav5jgJEzoYmZhYQHT09Pu2dQ44EIuqm1zQaLBQZQrJJNJ5PN5b+tde74GEamERl26mpP1fCm7GAigmlB9Bs2fqu+YejDeb2pqyslIarWa907Pnj3rXMa6kGNZrIZdJzZ1E+rkrUYoy6ruvqhnJVS3pZo9e6x1iasWTScmO1ZZIyVKo8/n1AWwXYxFtSteT/VmlhTQ96XbpeozaTAa4KeN00VrJpPx9pVvtVpu8fnNb34Tz3WwnigHI1Rmom0T8NMc2QBCbU82BaGVrgCjG1KoEapBUjo/KDgvMs1dPp9317jllltwyy234NOf/jQA4MEHH8QNN9zgPZNNgaXtW2UD2q9VClAul9HpdJw+VCVNpVLJM4Jte7Xac+0Xduy12l4SAK9+9avx6KOPug1e0um0W8jzfUYt+DSNlNYHy6H3p32l9gvLw/6swV92zB+3aNb5WRcyvV4PxWLRG2cvNS66BjUgICAgICDg2YOGFZk4wmrRdQMF63W0zKv1JPFYGwuiC1IaRXodvb8at2o8tVotzM/POwav0Wi4477ne74Hr371q/HRj34UwNCYJOHE5PeqpdTNOeyCZ9xi0gYLqtbfGmEaR2KzEPA++twsC72LLE+1WnX5VF/zmteg0+m4zUE6nY4jIGKxGCYnJz0jUJ+DuWD5m2KcR4d/a+AZ78t6UyPfkoVap1qP6oHJ5/M4fvy4Y/cvBy56FH88HnduGxuMoJViXe5aWXb7zNXVVcc2ra6ueml3GBjF66tbzrJluguFrihVhMyyKZOgUcXqKqSEQF2HdnXJhl+v193qKpFIYGpqyq1uLLujqxal4llunsfGM062oHWaTCbdyieZTHoJ/vV8BjrpM/KatVoNKysrXmSkBq0VCgXXCXWwYtSisrJWiqAMt75PTRLd7Xady+jkyZNeijGNJreDsa0jTTliV5a6dad1o15MkKWnu+Ts2bMjbnz+TTdT1LMAm+22Xq+j0Wi4yOhkMukFH8ViMc8FxjaQTCYxNTXlflO2ndGvfK/a9yiN4eSikyXfq7KmvAbTs/H4arXqGHwyGBrFz8+NRsPri/l83pvY1JWuLIeyLiybnVT1eA30sBIQZTeUebETt3Wj23EPGB131I1vWalut+vqj14OrSPLylpDhd/bMYHHkX1i39bMAOyvdmwDhsyUujHVoBoMBl6Wkx07dnjjdrVaxYkTJwAMJUKXO4VNQEDA9sZFZ1A7nY5z/1g31flcZDZaFvBXKkqlT05OusHSpk8Yp4fSVRn/1slAVzG6MqWuRo9V7Zw1GHQiWltbc4M3c1oCw8lU87Wqy4VuMJ1U7EpL7797927P2CSs8az0P7Wkx44dA+C75Thp2Z2lgE05RdQOWISu4FTaYbWDWl82RZQaiOru1fe3Z8+ekclW24+tU30+1eZYN5mmvbGLg4sJGhzc6YSaZQAjados86H9qtvtelkeNJJ9MBg45iAej2N5edkt9G655Ra32qeOjKxNqVRyhh/bt6ZIUXmN5vBVLVer1XLbiLKc3IKYkcF8J6VSyVvIAZspqI4fP+5+SyQSOHfunHv2a665xi26aExatzuh71hdgboQ0Prls9v6V4ZG24bVgmmbtUyH1Q9rux4nRbGR1TY3o5ZV27TVGsbjce8eVoOq7UqzJ2h9RBmSPC8Wi40sIHSxz3Z05MgRVKvVy5625kqBasU3NjZcHU5NTY20HxvVz8/WPaxt3kbR2wWwzQqhC3pC5xfAZzeB4YJSF3Rc+P3lX/4lfu/3fs+NRTt27HD9OB7301NZvayNs1CGzxIK2q9szl3to7wPr6kLNytj0nmGC1bO7Tt37nS66Wq1it27d3sSOz7j3NycZwtZhtbKLawMRzMK2EWmnmvHCO3rWm98TkLHR0swWmnapcZFN1AtFa5shWrNrGGjhiaZR9Wh6F7Teh/VWdjUQ/w96m/b0K0BpQyIdkLVwTEPqj6vTgzVahXnzp0DMBRI28nHdjT+ZjudNh573tLSUqQRxWdSfaK6X3K5nNMnqpaPz646VzUYtY7JmvI8ZTvVIO50Op7mx+oWx70XGqsaOKL1wPsS1ujQcms90EDSegI2k+XrJK2d/GKjVqt55ef7YKJ+1VrpgKH59RjEBAwNuB07drjyHzhwwL3XfD6PTqfjAhOWlpZcvTLPqb4LNeiZGxDY1CEDmwYpJzPd8ld1csBoyjftY2pwsQ/xfjMzM441P3PmjEvPBcBb5HFiGaev5LX5t9UhE2owchIZ597TyUC9RTqJ8V5RhgI/2/yf/E3/1+tbg9d6JLT+LRM7bkHNc8jS9no9b9zudDbTr6nem/fTfLt6vePHj7sJrd1uu/aXzWaDcXoekIVWrwPgGzNsV9qe1OgZF9gGbMYv8Hsdj/Wa/KzBT4Smm7LXj8fjnhv73Llzjrw4cuQIut3NPMzjSCz+rW3Wzn/8ztYN5ylq7jWnMIOCoxZq9EzYxWnUHMtjOR6urq5idnYWAPC+970P2WzWBVfv2LHD23hFxwEa5YRdaFvZgcbaWAJHy2mNaevB1EWAdf+rHcO2VygUkMvlXNDq5UA03RAQEBAQEBAQEBCwRbjoDGo6nXYWtnXfWq2hrhx0tZdOp92uUICfgJb6qyhKPEo2oL/p6p4rTD1PVyLK1ul2ovr9ysrKyKqRz95qtbCxsYFrrrkGwJCleuKJJwAM0/jYVZOudOxOD0rpWz2qrq70+S1LpM/HuledK59raWkJsVjMMZ6rq6ueG1pX11pPvJ4yTLrZQaVScQybsjpcvUax14S6CDWa1NZ9VFom3k8TxmsQgDKog8EAq6urkVurXipopgONiFbmWt9ltzvM4sBsBolEwpVzcnISu3bt8lbDmq5mbm7O6cOXlpbcBg/1eh2nTp3yZDNa57t37/ZSmSnbsnv3bicNUJTLZVSrVU9KoTKZcrns/a3b6WmQAN1hfFZNB3c+nbGylmTilY1U9sWOJTomWfei1otmuxjnFmN/iZIesO1HSQjYvpWp1Pqy8hpbRs3QoeOqfX96nkbZa/1xrNJj2c97vR4SiYS7XywWcywpPQSaPYSfo2RJAZug563f7zstL3BhEffAZvu0TDwwGsWv791GhvOaUQwmj7dZQIBNORjH5EQi4bw9yWTSY1ety9hKvlTeYtulzchhvRq6y6T10mkMiOrMtb74m84t2ke0rnbu3OkYb/UsAsNtqHVu1jpU5pNjfZRHkeXQd6N1bmWSqte39aJeOVse1ivvqRKCXC6HXbt2AQAeffRRXGpctBlYH5aD1fr6OhKJhJsUU6mUl7pEXduaLgrwK9S6O9W4s1S14nyuPKs9rFarXoCAuhF0IFetICcK7SCqqykUCnjssccADAdn3S0qKmCB17f6Njv52ChGXkcNd/u/ui5p5HAQLBaLrmMdOHAAlUrFPYfqD7vdrmdYRC0ICG3YuVwO1WrVvcOoNBVRbnx2fl3U2GAXLYt10/C9NJvNEQ0Vj9UtaFdXV7GxseHV56VEt9t1W4rOzs46F2oikXCGHQBvF6Rms+l2LmH5qfk+cOCAC2oitF6TyaTTtq2trbk6oQG6uLjozlODWBcbuu+0HUyBzTqrVqtYXFx0v01NTbl3x93JOE5oLuFWq4W5uTmnc221Wk4rvby8jOnpaedG04nDanOtUW9dalo/4xZ2NogynU57k5oaanpNvdfU1JQXnKkGNicK3k/bKQ0FfR7ty2pw28DQcW7fqHRA1iXM89rtthc9bqVMbHP1et0LklISQrXMALC4uHhJJTNXE3Sc0tRDumDieBhlzAC+oWEXO3qslbHZa2gb0jmW+Yyt7jUKut99q9XyAlwZE8FrKNnA/sJ76JbMNJ51/rcLRr2uzg92LrGLUoWSIZYYseMK+wVlWnTr67ug7aEyQh1/bOqwcRIIljVKS8vFSVS7YV1Ytz4wXFTaeVXlcFaWdqlx+e4UEBAQEBAQ8G2hTJw1CNRAsdkqVIcP+ESFGpK6GFLNJYki9TBokLA1lqI8mVoW9XaQeTt9+rS3JagaoLy3llW9E/qZz6/PatNRKSOqpJayy1Y/HqXb1pgMXbCqB/P06dPuGWu1GtbW1rB///6Rd5PNZiMNbR7XarW8haOyyQoa72pAjtOyKvOqhAKfzxrlqi+2Xm8SBpcDF81AjYrsPnv2LGZmZrzViXV/KUuh4t94fDOS3roq7Mu2K8io1Z6FrkSU3bS/AdHyAH6vrEMymXSMAYOEyIrl83kvolqjG+0qya6KtfNq+pdCoTAiB7BMkK03lrNQKLio6lqt5tz4nU4HCwsL3jMyAIL1q6l8tHNYZltdMcrSsoPa52PZdLWsnZXn6ntQNkoHUd2piPfR/3nsYDBwx5FRZudUJvJSgeWvVCreLkj0OgBDKQnrnC5VtrlMJuPaGAMLxyWh7vc3AwrJlAPAU089hX379uHWW29157CelpaWsLi46N4XM1Do9aKCfAB4qdQ0D+DCwoI3EcTjcddvVlZWHOvIeuH7ufnmmz0GV9sF3V3KEigzqQFVwCazzHJEtUe6D6M8BpaV5Xe2HmgIRE1GZDV0ArLsbZSBwTKrLEjZmGQyGVkuOz7SkwIMXfOJRMJJR7Qf6vfAsM1df/31AIYM6okTJzwpB8eVe+6557JOZgEBAVcXLjqDanWI586dc8ZPPp8f2YGEg6y6vmjt62rIrs4IpaPtgG9d6GrQAJs6KN5P86/aqMUoY5d0d1Q+U056O3fuBDB0H3Pgtq5xfZ6oSUTdn4uLi65eJiYmcM0110Tqxji5qjRC66Lf77vdOx5++GHPHW/dFtSOav5QXofPbhMdJ5ObqapYT2owWPeKTsgqYbDuBqsb1lWhbs+qOlO+D2BzRxFto2tra67MeuzlAO917NgxZ6TFYjEMBgNneHa7Xc/9n81mvR3KdFEAbEoooiQXPGZubs4Ze0tLS3jiiSdGovN53NLSkpfflu2YBiifQVPgUEeq/UGNKC0n86kCm+5qRnhns1kXiUuNmdUXE1bLqW1fF2jjXJ8sm7YDq/PUdmqThuuxdryyC1Ceo+Wy7JRCmRT2CXVZaj1zgWav1ekMUwBqfmKmVWMuVR5bKBSce7Lb7WJ9fd1t1Wh3viuXy87N/+ijj7qdvwKePdg3mIOYcgu7RabVfRI2lsG6rpXoOJ8m27bNKNlJlDyG4722NSUXdu7c6eRGOqf1+31voWZlOfq8HA+0XytZpAs1qyW3mTDs/VU2YHWnBJ9R5ViamSKVSo1oe7VsCitxiHoH9n9mSdKx074v/q27SvIZlRzU8zR3sUoIOHdeaumb4qIbqJqgfWZmBhsbG16CaWU+G42G61QqlOaAri9XK1pdEraDaKe1Ogu+GGCY6kLzMOqAbhukulnUYKQ2UH/TVDrApuG2d+9eZ4B0u92R9E3aeWyaD2WgU6nUiI5LO7N2aOv+Gdfp5ufnnbFuDVqrb7H6NvubGiHa6ZTlsdcYp8floKoMU5Q+ir/x+nRJqSGgW362Wi3HCGndT0xMIJ/Pb0kAR7e7mYjd7tSRTCYdi12tVtFut93gPj097dWlsr6aDJ/tlIOZGppzc3P46le/6s7TtFLNZhMHDhxw5506dcoxr2TCNUG/tsWpqSnXxzQQioMlNa/axqhV5aKoWq16fV/f+fkWdsrQRmnO7MQVtajlIG4NSvtu+JtCB3/LjGofsedZxjTKG2QDOVRPzvGJfTuVSjk2noYrdc86MdGY0By6LNvq6ipKpZIXAKrM+Oc//3kEXFwomVGr1ZyBms/nR7SU6v62m6hEedCazaZnVOm8SL20kgvaTqM8AOp9UKiXUOfp6elpxGIx5zE6evSoK49dfNk5B8AFLS55HZ3zonTZtqz0HkUZxLyOEmcaE5BIJFy5Y7EYpqenvTHI3pPjYZSOV+8/bmFq9af1et3T4Ko3zRr2fBZgMy81P2u9qa1FIudyzo8hzVRAQEBAQEBAQMC2wiUJkqI1TpaKFjhX+fwul8t5Ef4EV0LWnQZssmOWriYsG6crCnXjk8Xkb5a6Vjper6MuZ7IOURFvTKDLVdL8/LzHTmo5LROkrkObNH95edmVKZ/PeytI6xpQ2HqyWk9dlelqT10BLJeW12r3lImyq26WVZk1W2Z7Da17m6ZM2ScrrNc2o6vHfD6PfD7vXGhnz551rGMikXA7gVxudLtd5xqdmZlBLBbzmC+6UHu9nuf+X1hY8Fb/i4uLLp2LppNhxK26gvUd3HrrrR77qkyEMjP79u1z5Txx4oR3f2BTErOwsICnnnrKfT8/P++lktINJgqFgvtN2waP1X5ivSLKBlgNqkofbCSulteyq9oXtP3bYAoti02YrTrxb5dM3XoTeP10Oh2Z8sxqXVUetbi4iEQi4Vjuqakpx743Gg2POUun044xLRaLyOfzrs6OHj3qeYoA4Ctf+QoCLg+o3WWfpRyn3++7cYFeKZWm2awSUeM451Bl+QmbTcXOpyrPivI+6Hxpx3LraeTfmkaL2XS0LyWTSc+jaRlAvYeex3vxf5ueyrrV9ZraD1UeocdywxTKFXXcarfbXsabZrPpbTVuteda3zpXaoAa61Dfm7rjrTRP3zHvw7oYN+daF7/WYTI53K5ax/tLjUsaxd9sNjE5Oelt3cjdI8rlsrd/uhofgO+StpoUG6yglanGlXUB64spFouuY/PFRhmMlnKPmsh0IrTGFe9h3RF8Fv6mBqntENoRJycnvb3WtUx2otVj7ECibhutaw5AWm/j3B9a94xK1HcRNYDwPLtbT5TGqdVqedGmOlBRi6juLbsY0bKpW3ScW3WrocFak5OTrt7PnTvn7eJTKBSwd+9eAJt1BAyNoYMHD7r9zRcXF52rnBOCbuGrW9/q4L60tOTJSlZWVpw7bmpqysufqgFGWo98Byq30fauaeXsrnEaUKXX4Pvn/bS9sV3aPsl7W92Utk1dZFpXmB3kCbZbbX/67DrRad1a6YENFLV91KbNAYaTuG4je8cdd7gcy9/4xjdGdvvh/SYnJ5HJZHDdddcBGBpBTIN39OjRkWcM2DpwgZFOp700cxbjFu1RBpcdA+0ORfxN+zTbq47ddpHG36zhouOuyouKxSKOHj3qpCbFYnHEWLbkAs+1u0XquKLn2fHASgG1bqzL35IbuVxuxFZgmWyGAcbbxGIxpFKpEQkCy6zjn/ZzjgFqN1gXvyW0ohbGjMHgeK82DPOJa7uxskF9Jo13sLrnS42QZiogICAgIGAbQgNiAH/Bxu2LNRm7Eh1qoEXlheZvvV7P86zogo3Mo7KPtmxaRhqplrCwBuD111+PI0eOABgaqNQ3Ly8vo1QqOcY4ijzQRag1Qu3zW6PUloVQg0wXkbymanJpvJ44cWJkkUoDlcHR+t7UyO12u94zWsZSg52tQc1j1RsMwAtIJyHFY9Swp8dFy8M6iSKZlN21dXmpcUkN1H5/uAuGNggyMUyVYylpYDMCPGr1wUpXhlLpfA3YOJ/LW68ZxfQp86gdVt0PgO+C107PgUN/s0FF1j3PZ7CdSX9Tt611oyj7aBtzlEg8qlFypRsV7GRTGGmdkYnSTmPZIG3g+lnrRQdcrlCjZANcgWod6krSpqey9U5WP5vd3BPcMtyXG3zPS0tLmJycdOXXAMJer4fdu3d7qb4YbNRqtbCwsICFhQV3HZ2UcrmcczOdOHFihH3k/efm5lxdnDp1CvF43LmGtb8x+bZOgup+K5fLjgli++BxKysrrmylUmnE/aTMqwbkaHtUr4Ay/jyPz8BUWFo2O3FpJKxOeOqC1/ZnZUg6RljWX8eLKBe/emfsO+EkpnkTuUMM+9pDDz3kyvj6178e/X7ften77rvP81rt2LHD3YPsaUBAQMB2xCV38Xe7XW+7MZ1s1SCwbvRxrnoOypomSN0IOqHZNBtqGNlUQvF43HMdWvrd6tJ4nE1npJPPiRMn3DPqDgzWnWgNKTXudAVl5Qa5XM6LXtd64gQaZaRZHWk8vrlFK+tPy2oNdzXedTJXvYxeX78ndIWoxkS9Xve0kGo8a4qhZrPp5eC07mWbOkMZAusisnKMrQLrpNUabpPLvlIoFLz0bMqGlEol58av1WqoVqtuETg/P+/0oouLiyM6aP42NTWFcrnsfqvX6y7jxMLCgreTTTwe97JytFotb4c0gv2E74tlBEYlH+l02kUtVyoVLC0tedpmPg9X/TYZN3/TvxXMLKLyF33WdDodGe1qjV89jwakppjTcUjHGdsvVdpgF6C6yNK+n81mvcW3unPtYrfVark6K5VKTvJx9OhRZDKZEHV/BYGLw4MHDwIYvndqildXV51WGRidt3TRqQQBr0vZ0MbGhlvEWJe+3ZUPgMfY2vnA7lZodZ/6G6VC6XTa06uvrKx4ecmVMLISNisd0gWeJV50MakLT51To3S1/N5eBxgm56cGtVgsuucvl8tYWlry3o1NR6ljjZW8EZrFh+VSRlPLsrGxMcJ+Ktuqz2/nOjuP6mLebmlqbadLiUtqoFLLoWyP1Yiwkeh2fRyY1cDQxqTMJABnXCWTm+l4CDXe9G9Nn2BZGytgtgaoGsc8n8fyxS8tLaHX63lMpnZWa5xbZsZOcDxOBwROyGokjtPx6jPQcFfmi52s2+1icXHRY1u17tSNo0Y9n4kDS61W8zRHeh2WQT9bkTjvvbq66t6p3VvYanFsOiDeX98hn4FBee122wWKFAqFbZPDcWlpybmJksmkM1a5fSQZzaefftpt/Tk/P49KpeICLBYWFlwO0Wq16thMYFiXhw8fdveamppybN2JEydc3VH/pnpV3YZ0ZWXF8y5omrVarTZ2sw1dWDGAi+XM5XLOYE2n044NZMoplq1Wq7n7sQ2x3drALasxZzviIkfbnxqhLAPfgz6DLm40sGBckBdh2622ae0X8Xjc05Ap9N7WJahs9fz8vKuLXq/ntKoBVwY493GB2mg03FidTCaRSCS8+c966QglUmjY0EDt9XruvGaz6c1HzHmrhpZ6n5LJzVzYNremXSzq9qVq9KphWyqVcPbsWbeg5QYbOh/Zcd4a5cBmv9Nya/+3i1nrirdxHVpWjqNnzpxxx7Ae+Xlpacnr2zYPqXpZ1EaiDWHtAX2/WtaNjQ0vaI5jszV01QtmtcdRaR55rnrdnn766UjC61Jie0SHBAQEBAQEBAQEBPx/XFIGlTpATRod5bZ3hZGIQ9V6nk/g3O1u7gqUyWQ8zasKkVmOKE0oV2j6t+o8oyQGvIYyJM1m09Ng6spXn4PnaTnVhdDv992qpV6vj2hkNcJa62acKJx/qwvZMppchc/OzmJ5ednpGqlnZNkAeNsaaqS8TaSvbnvL7mo5owT8ALC+vo5z585hx44d7jdtI1anrMy2rjyVgae7V7ddJTt5//33Yzug0+mgXq+7tqueB7qxyDCk02n3LLVaDcVi0dsmlOxGuVxGq9Vyx+rKuFKpoFarudV3qVTy3Oi7d+92LrjFxUXXJpikX9uxglpywGcCgE39Ku/B+83Pz3t9Q6Uxqhfm3xqgoG1M3e9kQdTToX3BBinYXa1YFia/5vX5HZ9vnBvQyoJ0n2uWh1DmYpwWXdkWXmec9ly1t5lMBq1Wy7UrrcuA7Yt+v+8yt+Tzec9LqO1seXnZeV3YpnXO0f6nHqj19XWvr/B3YDOrS1TEOzNJqIfBun91TuVv7Cu8V7PZdKysjdpn1g71Mug99P7qFeWzWG8FYfuyllPnMcoRomJZ6Hnk3KmeINaplQRpuQk7bxcKhZF5HdjcyVGfP5VKubppt9seYxultdd60XFEGVsdr5rNppt/n3zyyZE6uNS45FEh6vau1+veLjnaSDUXGD+r61phdVtEvV7HxMSEa+zqDrQ7UFmXn2pEaGxp+W3n0nuT8tfgGz4LtT38nefpddTdyYFE3Z/W/aITlbo/tKGzI0e56lkWXrdQKDg3MQBcd911Lmft8vKy+03dovouWDbtkGqAzM/Pj6SBUqNGDUYAbjBeX193O3IAowEsCn12Ggh6rLprl5eX3f127dqFs2fPYruh0+m4PMJ03QObbUe36aWxSvkHDZBOp+MtJubn593C48SJE04Dxu/5vkql0oh29aabbgIw3IGNEoJyuewFJuqg2+12vfQsOllSGqLSHxvhymMrlYo7rlgsejIMO8jq/yopYLtRF6MatlH6d5bZllsXOjyGZVFDluOF3WVHxzU9l9dXSY66/tSw1efg8+kCTw31XC7n3le73R6RxgRsf6isK5VKuTkGGAZQUq5E+RswGjuh79zmCNV+S7c/x2DGI9g5h9exuk/r4tf5meedO3cO1WoVu3btAgA8/vjjHpkxGAw8eZh1z/NZWK4oTeT5ZDXAqMtdYXWt+r+OOdlsFmtra67O9+7d68amer3u6pJQKYRKLjgH6jNpPeq8qdIJboHLBYvaN3bRrc/J9hRldFuJYbfbdePHOEP9UiKkmQoICAgICNiGUKMP2NSdAkPNI3NaAv5ChQtZXfDogl2NvlQq5QxSbgRCYqVer7vth3mu9VrpgotlIQGiBIkm4x8MBiMZYXg/u/mNGpBqzNEI0wWeGnYW44wrrSdrtFqPhNWEqjG7vr7uPsdisZFAJI2/sV4WDY60+lEb+MU6jsVi6Ha77l0NBgNvYZHNZr2NOWxbYP3bRYYubJrN5sg7vZy45AaqCo43Njbci+j1esjn816j1ZWJrhqV1rfiZ8Xy8rLXSHQlRsbGuieAUfZPVxd0T0cFJFgmaH193TXK3bt3o1AoOHZOV2yW7dGOQ2mDNgQtl7rAo2QS1m0RJVvQa/Ce/K1SqTh3MDDcf5c7E+XzeSwsLDj2ThNI28hAK1BXNk0HVftMzWbTib7ZRnRgtRGbUQMR3Tc2fxuwOVByL/uvfvWrXifeDiB7RhH+zp07HWuyuLiIM2fO4PrrrwfgBxH1+31P9lAul71oWAAuMKpSqbj6IUuhO9KQpWPAHN36c3NzLtCGgyxZg1Kp5DEByq7a9qzvUlOLkTHleblczrHd9GTwWJtSTic2df9bF791k+mkp2Vlv7cBmjxGJ0d9Pr12s9nEsWPHXN0q46EBZSyHXp/vhcdYlsUGVGldknE/duyY84asra15wYsBAQEB2xmXhUHl4L24uOgSyU5MTHhaB3WlafQ14NPcHJg54KvxyhQ8XFGoi4HX0Khz3s+u5KyBWKvVvBUWr8GoYtXyUK8Rj8exsrISOSlS02NdFSyzdZXaiEIb1aeTE8Hn03LrxKSTq35PVy9dwbwPz7HbouozqNtCJ1Qahboq43lczWmUpF4/kUiM7DI07r3os1n3Ep+hUqng5MmTIyvj7QbNN7q4uIhDhw4BGKaWqdfrzqjWiHHNFwoM64fGK3VNbEua9olabV00EPv27UOlUvEWekxhRDlB1EKLC0xtH9blHpVmhdpLdfnrFqnK/GiaMf6ufXlcxgddLNq0c8om2EWwslBcRKim3WaKAIbRvJVKxUuPp5rQfD7v6dR16+eVlRV861vfGqmjPXv2oFQqjWjTWa6VlRVvEby2tgZg2PYXFhZC/tMrCOxHbE+DwcC1l2w2i42NDedKbrfbrk1kMpmRjCjq2td5c3193cvM0W63vW2Wx2WlsJpTO2/qWK6Zd/gdxzDmRAfg6dJ5Py27Elecj/QZNROPrUcbta/jlo5NVsqj85OWZXJyEpVKxS04H3/8cTdOA0ONN9nuRCLhiJ1SqeSNlZp9pFqtenXMetXPPHZmZsY7t9fruXFmenrak80x1gPY3Nqcv2taSTK0ytLqe7vqGFTAF/nS8JmdnUUsFnMv19LHVtCsncKmiGJDm5yc9Ab7Xq/nKG/bOawGw95bNWN6Xr/fdx2r2+3i5MmTLshm3759Lpik3W67bSkBv2GpMc57aIe0E70NjrDGrE5GltFR41b1cGrAqSi7XC6jWq26QW/Pnj1ux49ut4unn37aaUKZAkSfSVlLXUTYiV8NEKs3Zl0wTZcuVKxRHaWxVcNY3yEwTA1Sq9U8Efh2hPYHHWh27tyJnTt3eppUDoo6iPNvnlupVLzBXQ09thvVrupCbmpqyl2nWCy6yWxpaQmFQsEZkNquWX59r/oetd+qJjmdTnu6ax2A6bGw9cTzlCW1ek67sLJ5eKPYSAvtQ7YOLUPMhcHx48dx7tw5b1Gpk91gMHDsuOoJdewAhm2A7fj06dPehKsbaDQaDSSTSechOHbsmHu2lZWVEBh1hYGLMraver3uvCIcw/muNc0R51Udn9l+1tfX0e12HbOeSqXceRMTE97Cb2Njwxs/td9yzI0iV2jMcHzlwprIZrPOMxeLxTyvp44pUeOIdYcTajzqfMvfbL1q/9F5UxeY7GdRwU65XA6ZTMbdZ2JiwgvUzeVyXkwK7xeLxbyYj6g5PsrQtpphjrdK7tDeefLJJ73gWh37ad/oxik6bq+srLjALyUDozy2lxpBgxoQEBAQELANwUUQjZBGo+EMi6mpKY/VV+9hq9VCr9fzcqYSExMTqNVqznjSZPAMwuJCKZFIoN1ujwSgAtFJ7Akazwz0JIsPDMmbRCLhzhsMBt7icXJy0jO61GuoRlJUIFBURgxeV3/TxaY+PzEuyEoXmDyG9VgoFJxhBwzrWZ+b3uP19XUsLi46QsDeV3WuSiRZg5TBsjRmGeUPDBnUjY0Nb2HMctpsO1rGRCKBWCzm2o0GiG8FLouBqi/09OnTADZX/WTquPIH4CqH0MZjO6ytbN2pKpPJuJXk2tqa1xCVbaGbwArM9d7aSMigrq6ueuxLrVZzDG4+n/c6ukbqs/xaFoVljXRAsCJtLZuyVFGrK3Vn6qpNMxj0ej3v/mfOnPGiEROJhFtRHTt2zHUyugmiOhZ1yOoetcxTFGvF1Rs1jqlUytUndaxRTDp1w1pPR48edeW3GsTtCGUxW62WN/Ax9RPgu/EZVU/WQpnxXC6HXC7n6ZDZjqempjx3dTKZ9Fx+VnZB7N6922srqklmG1NGQ9+1tj+d2GybVqG/7S+2n2gbst4EG+FqJUNEVJQwz1M3pdWmx+NxV2eLi4vONc9xTfualitKhsLrr6ysOKalUCg4FrtcLiOZTLp+ODk56c5lxgx6PUqlkrdLT0BAQMCVgsvKoOqAfPr0aczMzLgBmKtBImo1CGzqVdQw4QS2vLyMXq/nVio2JYeu6CYmJtygTi2jneC0zCz38vJy5GoQGBpzGmEJjGph9PpqpKkhbHdZ0hWT6pGs61zdmHZFqJobPo9O7kQikfD2V6cOlHW2trbmrj0xMeGMo1qt5ukR1Vi1xqc1unXSjloRUlKQyWRcm6BGyrotgM0Vob6LcbkDtzPU7a15aVdXV7Fnzx4AwKFDh5wBsrKygrm5OS+1E41c9hM1IHVBpoE1GnxEN5HKM/TdqcvbGnvab+Lx+Mh+91avps8c5X5jIOU4NkWvY7Xbeoxtj+djYVQOAPg7zug4pDlaFxcXXdujFp26XatV1WA+XQhydx/dTYznraysoNVquV3C1G3faDTwla98BQFXF3TRqTtH5XI51y4mJibcONftDneK4jiZTqed+7darbrvgWGb5JjJeYH9aHp62uuDUTKeKDKF4wvn4unpafdbpVLx5h9gU5IQi8WQz+fduEIGVYkWgmWxAYJaNpZb60nHEEKDE3W8i5q/dY4jYcLy6PiXz+ddPU9OTrpnzGQyLgAWgBfIqsGdvKbKC/v9viPBqtUqqtWqG++V1FheXkY8HsfevXsBbDLhwHB8P3r0qCPzstmsa1MbGxvemLTVkqDg4g8ICAgICNimUENrZWXFW4jPzMw471cmk3FGSKvVwsbGhjOgqDvlb+122xE7mgi/Wq26xRSPVaMpSlduA16BzQWieio0a8jy8rJbxMViMbc4SyQSnofNBkbyO/0tSj9uM9XU6/XI/MEsq5JFajyz7tWLoveo1+uu/geDgef91aC1VqvlyLh4PI5rrrnGaYAbjYarN271HBU8Ho/H0Wg0nEG6urqKxcVFT6qg9W8DRtkWOp0Ozp0758VfqPGsC+Wt9jRumYFKbQOZxqWlJee2pGtQE9CqlkbF4IVCwbGi/D9q9cPcblFMGhugsk0a8Vav171r2xUOnyGTyYyI1NkodWcY62JWxkYZUy2bltfmpNPgk6ggGF7Hri51oKF78tSpU2g0Gp5Ugsjn80gkEs6tXy6XXVmr1arbYxnwI6xZh2TmNChseXkZ3W7XreBrtdrIblz8W13GLIt2Qr6j5eVl5HI591uxWHTPSlf/lQCtSz5bsVhEq9Vyz5FMJr3dvgC496MSiGKxiGq16phY3Umq2+16A3in03GyiqmpKedSJuymEspa2sAEdfFrkJ6yrXYQ1LbK3bGAYXtrNpveblhW7qJtXrOBKAujTAf7TFSAZlQKKRuwx+tUq1XHQAPAz/zMzwAADh48iJ/92Z91xsDc3Jwrvw3Q04AXTjB8Xy984Qs9dqrVarmo4JWVFS/rRkBAQMDVgMtuoNqUKESpVHI0c7fbdRMxMDQ86Oali50TUr1e96Kc0+nNrR/V8Go0Gmi3224CUBo7lUp5om1NXbGysoJqtertsqOa2mPHjnmrJLqgM5kM5ubmnNGtK026G9QwtUaouh3V4LaRdFFuUl7Prjr1flER/wBw0003ebrGWq3m6qnRaGB9fd2blFV/uLS05AwpRbVaRbvd9twIRL1e9wxb1fkxolLftV6TvwMYWUX3+333LtTg1ojUKwX6ftbW1lybBka1nRrRrW2culI1sNg+CoWCp0/VFGG8Pq+pDAKvr1IEXfSoHMW2R17LfubiRbMI0Fimq1FdbLpYU91x1EJSXZM2D6pCyxnF5vBzt9t1RqLuHlWr1fDxj38cAPDmN78Z5XLZLWLX19c9XfDCwoLXn1iW2dlZb2w7duyYV/5CoeDYsrm5uWCgXuWI0loXCgWsra15qck4vzFAJmrcm52dxalTp9w4q8wbt+vUIKlYLObamiV4lN1dX193MoJGo4Hp6WnvOhz/19fXvfk2mUw6giKVSiGZTI6wnFZyx3LrsdrHbdaexcVFR5BwAW/nVV7f2ig6PvT7fTc3cutPlj0Wi2F5ednVo5ZfxyObB3ppacnNbZlMBv1+38vSohKAVCrl3ikX+lFkFxfPx48fBzDUppMZpTSAY+zKyoorL+V9UTEAW4EtdfH3+303AGcyGdfRms0mVldXXQWeOnXKC0DQoA01JnVCBPz0PEyloS9ed2DIZrOuEyr7t76+jpmZGW9S1E4Tj8ednqRYLLrGv3v3bk/LoknEeS6hGh8+m3YYNV6tgarsj9Xq6LHqbqEe1er+WJZSqeT0bd3u5lZnJ06cwMzMjJdgXPM9JpNJ9w6feOIJN5jl83mkUin3m2qKz5075wJOeD+rlY1yE+VyOU9TaXfM0FRSu3btwv33348rFdaN1m63vcA8PufU1JTHmmazWZfvlC4v9iNtVzRstc2pgViv173BXPVoNu0TwckrKviNfUhTktlJQAPsaHy1223s379/JI0Yr6kaOJsORfuiThbWWFUoQ8rztA+trKzgxIkTAIaTnmplOfb84R/+oef202jaRCKB2267zTNYOVGdOHHCG8uUgbaeBc1pG3B1Qr1nzLVtdaa6e1G9Xsf8/LznjeIYmU6nUavVvNSGRLPZRDqddtrRdruNWCw2EnwMbI5FvAfzkAPD/rG+vu6uo/lb2+2217dsnEi1WnW6es6jUWQK50qN0SB0sQgMY150jtD5z+o+1b6gwap9kX2eOWj5XMlk0ot+V93txsaGZzzqwldd/MxoQMO+1+s5+0KJMH1XOudrfEq1WnV1TCII2CRwVHerzz4YDLZNYGXQoAYEBAQEBFwBaDabnjGlzL0u9PP5vMsyAwyNIA0+1VzUUcYZjcaNjQ1vM4Ber+cM4nq97nkl6/W6u9bMzAyy2az7bXl52ds+VfOCApspmKiBVK+ZXaCrjM16DXUR2ul03EJRM13U63WUSqWRRTH/t4t39dpopp6ZmRkkk0lvm1h6eqvVqpcSrNvtOjZz586dI8GXmq+0Wq16TCyPm5yc9LyPzOyhxIIawcCmMXvq1KkRQkr/pvGaSCTw1a9+FdsFW26g6opP84JVq1XnPovH44751EbI85TRsZ2X55EVUkZRX2KpVHKd8sknn3Qv/cYbb8TGxobbslRTMqXTaTSbTde4Dh486NwIXNloZ1JWyjI8VkMX5dKxO3VEMVTWrRmlSeWK1Kb94TWV4s/lcrjhhhvctfv9vhPRqz6QgwHrRlerXKHqgEh3T61WG9EA2hWxuvb5nW6WYMFj6Brm/1cqqtWqa1fU8up7VncWsynwe7qjLGup7Yq7EmkQBc9jAIW2FZWqnE9LyiAC/maj6vkO7Q5NqtFmQAcwHNiLxaKXB3Kca07bFNufltvez7L2vIb+r2MG+zfZlHvuucebZKgL5oTOCb7b7eKaa64BMNzYQ8uSTCZdf+ImC2SPV1dXPTZXJ/ilpSXs378fwJWlsQ4ICAg4H7bcQFVoSgN1gVvX/fmgBsv5rtFsNr0JVVeCpVLJGceM0mO6Bp0Iqc88ePAgAD/vpDWWVedioYYmJ/Lzpbax7o6oXS54T3sMy6aGhk7KfA47gQNDV+JTTz3lbbPJ+3EHIE6otVrNTdLdbhfr6+uRqbe+XcqfgCFYZ9x67/rrr3ff60JLg8N0JyndDQTwg6T6/b5r77ymzTOsumB9fxoBDGz24Xq9jq997Wuetoz9S7Xf9vnYb8hKxGIxx0rMz897UhV9HhqZqhG1+YU1Elfb2TgdF8ulrjWVm5RKJdfGn3jiCZd7NJ1OO+0570P26MCBAzhw4ACAoYFtWQ+t51wu52nYWcaTJ0+i1Wq5FFQA8PDDDyPguQH2qYmJCWxsbHgLf/at2dlZrK2teW5dnWMUlOPwGrlczkvRp7sz6v0BP0WRMp/5fN4LYE0kEiNjActdKBScbIGLZd3pz2rNo6RpLLvOW/F43PWfqPlOF6CWldUFudZdv7+Z+o2yKCUBeB3mP2Yd0OXPZ1S5kGr84/E4ZmZmXBpBDb48efLkSOwAz4l6Ritfss+k4zYXv9uJPQW2mYEaEBAQEBAQMB6PPfYYAOD6669Hv7+ZJ7rX6zljxbrRo7wbaqCpd02PXV9f93TtmUzGueNVc8rrqB5ar5NMJkcCprhoS6VSnpGloJcyKr82j1dvy/kMS9WuaxqqZrPp5V2111e9Zjwe97KfpNNp98w2KFljYFSHrhp7PiNB7yYX6EoWMbaCoPdWN1+x7znKE8sFiQbMbmUg1PmwbQ3UKAbxO7lGFFQkrKln+v2+W92dO3cOiUTCsT+ao4zBI1zJaHCT3VPXBlvYFFC68gH8VZE2Mj3WBovoiolsalT0I1eauhLVjqz1rYFQ/FtT5vBzsVhEs9l0q+1ms+lWfnTh0k0d3JDPHJwUWq0WarWatwuUykLK5bK345BGrGqEv7Y/plVj++COVADcO+R1dEMNti1lTTWIjRotQjfeWF1ddX1Kgwv4mcfu2rXLTQBkLDS4g2XRACJgdJtDhcqE7E50WjdWaqIMDWUDnKySySR27doFYLiZAl3/yWQSmUwGz3/+8wEMN1bg+6nVap7GTFkVsmEahMa9y9fX13Hw4EF87nOfQ0BAQMDVim1roF4OWDeffuZEziTCuk2optwplUqegaCTm13B2YlunO5t3Hn8W92hiijjVV2eep4ayXo/RjtbyQGP0/RfwKZhwON4D51YaQBdiSmetgto1DN7ATMr3HDDDc5ALBaLmJqachphZUloWPJvzTfLXKe8zokTJ7y0I/V63Vts0FgtFosj+lFN8bawsODaUbVadW3x1KlT3s4mqVTKGXrFYtGTKWjbpKsvSjtNuYmyBLrYs/ptNZzHpZKiq09T8BB097PcsVjMJd7u9XquzzBX72233eaeT/OzWpZJdcDApnGcy+VcJHWxWAzGaYBjwjTjBdsWJTzsj1G7QUVJvOj6VfJC4yd0MTkxMeGdm8/nHQnBVIq8TqVS8RLlt9ttL0e4kjxR0eiE5t62LGm323VjGhlTmxaOn60cgP1V85Xzmul02kmgNM6C5daxhHXK8ZQLVSW2NHCJ5WF2lXq9jmw269hmZX5TqZQ3xvL+quVXaP0rIcW6uBLm4+e0garQpNnqUtjY2MDk5KQXjceGXa1WUSwWR7RuwKjm1BqTasxFGZqEdkALspvqqtFoR2WbNMLSGsE6ebNR6wARFe3IsmkAU7fbddqZ1dVVL+VIwMUB3ToPPfQQgOG7PHz4MIBNLakNIgRGt0FVtpOTAAc6DX6bmprC1NSUO88GF9qNFWgc8zjdrYb9ZN++fR5Tb/MF28lUGVO7RbBtw1FuS153nB7P6rq1/nQhqeMCy8JJ4frrr3caO44LwNBwLZVK3lanmiJMJ06r/dUIbdUIf/3rX0dAwJNPPultS6xaxk6n48V02JgIG+zIvkEDST2BuhCbnp72NKl2S22212Kx6PVPJT2SySS+67u+yws2Vm28LmY532r/ZB/kedp/VH+v44Ht89lsNlIDqotfXlMN9Gw26+Ul5TEAPI9Hp9PxPFErKysjKSB5j7W1NeclYhwLy7OysuKur9ua8h7dbtcZt1FuevUUbfWuUM8GwUANCAgICAi4AqH5TPk3EY/HPdZSM8xowB+/AzYNMvX2lctlJ7kB4OWzTiQSLvi1Wq26jVrOnTuHfD7vAqiATZlSsVj0dqdTY7nf73ub4QAY+dsGAvE6movY5ni1wZ3qbbTeFQWvr+xr1M6NvKcatsBmjmJlt23w0trammNFKdlSXS29V5r3m9fJZDKYmZmJLLviSjROgWCgemADKhaLrqNXq1WcPn3arWjy+bxrhEwiTNgUTzpAAKNRh8qWaAewjKWNVLbMkN5TYaUD2rEsg2vTTEXtOMTz2GEqlYo3EMTjcSfgD9H4lw5TU1NuQHvooYfc5LF///4R+YUyCppmymo0dcIql8uu/dfrdRSLRccaaFokamF5zZWVFdceisUiCoWCNwmwLzDdGyPgi8XiCPOpac6U1bBtVSdSfS5NVWW9EJbN0a13bVS9RtZr7kROqsoQab3TDToYDDwXoTLJLJdKfTSaWF2PS0tLnqszICAg4GpHMFAjYPNOnj592q1Sut2u29lhz549KJVKnhvT5kyMmiStFkXdmDbIQydMTorWsB0X7ASMau0I68ZUqEtFJ3MGhXGVzB24gGEi5rNnz3r5aoNr/9JgaWnJGYzVahXf+MY3AGy6ivibBrTR0GIbyGazzl1E44npj5LJpJcCRvOuqsFG41Q12brLmrrx1XXe6XQwNzfnjrW7WNmARRuZq33EtvdxTIciKjUXr6PRvTxOo2RVq83zgc10dPysgYA7duzwXH1aR5orVmUJ7D9//dd/7Y690vP5Blwa2MURobp/7dM2CFAlNouLi1hcXPSOaTabXgS+zls2VRHHBs4NymhyDJmfnx/ZvUnd3/Z/JXAsg6hsoxIr7I86Lug9VMZm2Vqduzjn6hxMkOSxi1Z+1gVmuVz2dPzcaQvY1PIDm2kEn3jiCQBD3SkXppQY8X4LCwses301IhioAQEBAQEBVyiiNMxKcACjbnSFelrm5+fRarW8RVMymXQBevV63W06wdRR9Bbs27fPETtcrGrQFO9LXasuBK2uXIN87SKV4LFRGTvi8bi32NTfovSpGpxoNfA2UFnLqvIALQsX3Voe6tBXVlaQy+W8DCs07KvVKs6cOePKlkgkvAWy6mXL5bLTDF+tCAbqGJB5KhQKaLVaLgAin8+7DlmtVlEul73OxIZVr9dH3JPK2tjAJ5soXNMIacdSl3vUIGRTSUUNSuw8GjQyLhBK6yIeH+Zh42C1urrq2NQzZ86MRH8GXDpwVa1pio4cOYJGo4EXvehFAOBNNM1mE3Nzc44ZWFpaciv4HTt2uG3zgCGLQOF9s9n0mENNR0V9lLZduu07nQ4qlYprD81m0zGABw4cQC6X8wZeHmc9Ccp0kEGNCi60rnZlUG0ksmYeYJ9SVlbbcDKZ9KLutXxkl4Bhf9cIYkWlUsE999wDALjlllu8fqfP/j/+x/9w59x00004cuSIF7wZEBAQ8FxCMFC/DRj9zoliZWXFpcdoNBp48skncd111wHw3Ziq1bPg5BnlgqTL3uZ15Gfr4rQuButuiDI8ub+xrrata0KfX7MWnD592hnoiUTC7RzFPYGvVDH2lQrNzgAMc2VyAXHw4EHPZazG5Llz51zGBW4fqtHqNFZbrRbm5+e99sgUVzRW1Zikayoej3s5WlXiUiwWR1Lb2EwUKo2xOu5x7VblKPoZGI3iV92nXlPd/1o+/RsYLgwqlYqLqK9Wq87dZvP8Pvnkk3jyyScBDHd8euELX+jq7+Mf/ziiEHaGCnimsPONZSbHgSQMMMw6oempVDMNDF3O3BiAeZV5X2b8APy8ycC336bbEjnKrkYtRPWz6ra1rPY61k2vUjsbta/3GBeTQXbZPldUWTS9nN1lTmUSi4uLSCQSXlQ/GVObiJ+7NzJI7WpEMFAvADqZqj6nXC5jfX3dMUOHDh3yGEttvGq82UlYdaTUmKgmjhM0DU7LLqleSFkihQ3oOF/kYjabdSlKpqamnNFx5MgRtNttZwAtLi56k34wTrcGakjF43GnKZ2YmHCLp1ar5XTCgK/7bLVaWFlZcRqpYrHoJhgyjBo0pQsk1ZKtrKy4qN1rr73WC8qqVquOXWWuXUJdamThzxfAZ9PlKFQjavubprUh22k30lD9Wbfb9XTlWmf1eh3Hjh3DU089BcCfuM6Hhx9+OBifAZcFdpGmhpxCF7nqbgaG80+pVHK69kKh4AwvLkA55ywtLXn6ajWQtT9aQ1EXjCyfGtY2fRyh2lALq13VY2mE2gUof1NPoA2w1LGLxrTWnabA0jnXBjPbQE2t83K57IxYzTt98OBBz6vC61/IuHOlIhioAQEBAQEBVyGUQRxHIExMTDhdKTBcaNIIu+aaa7Bnzx5nBKnkjNfVz9abp3+r8Wo1sZYZ1eCuqA0t9B5ElHZU60EDRG0daVCmLZuWT8tBljQq+NgSNloX9DypR3N5edmVrVgsugW07v63b98+9Ho95405efIkbrzxRlzNCAbqswAZ04mJCRSLRafDrNfrTrvHTqgrQXWTRkXO8zjVh0ZJBXTlyXMIdY3afZI1EtImJ7eRy+wUp06dci5LuhLItDHiMGB7wOblO3r0qMtZODk5ibW1NfcOG40Gjh8/DmA4CK6vrzsXtbr1VlZWnPYYgJcoOh6P4+zZs07LWq/XcfDgQQCbbnzVtar7S9lIAJ5eVNvjt0ujpn1g3KYS/E6zDygrqqxsvV73WA8rheFxtVoNDzzwgOtj09PTePTRR0feSUBAQEDAs0MwUL8DnD592iXLBeBc4cAwGtKmwOBEGTWRcjLlFnFRxqyuEIHNXYX0Ow0oUaNUgzHoKo1KZcVrcIs2TW3T7XZx9uxZb+eggO0Fq8Pke9y7dy9arZaX1JntttvtYmpqygsAYltuNBrYsWOHax/VatVjQvL5vAsg3LFjh/t87NgxL6hIAwZ5D9WBWne8/kZYViYqxYtKY7RN67HsY7y+LtY0PzD7L/uiLhy/+MUv4ujRoyNBXQEB2xHna5+JRMKRELVaDZ1Ox9t6WNlPBv8Cm6mj2HdsOjqNgVDW0rr6lcFkPmadm1TPOU4KMO45bW5y/q/3UPD+1u1v5XHAaF5my+jaoGXej1pSnsetn4Fhyi99Rr1GqVRCNpt1pMPp06dx3XXXOVLgakQwUAMCAgICAp6j0GwbNJzUCLXxEnY7UTXe6GkhAaKGJ+/BxaoaeppL3Bp53P4T8AORaUjaHKq8P8/nb1a/qvcgou6vx9rgZjU8bV5X1dbqLlCM/yDR0+v1cO2117py6yKZx/P/qakp3HbbbQCAz3zmM/jMZz7jAqpuv/123H///biaEAzU7wDNZhNnzpzB/v37AfjM0+LiIsrlsnP5W9YG2FxVaiAU85pF7fgTtd+4diZdafb7fS+A6nx78SqjW6lUcOzYMaeJ2djYcKwRg2kCtjeiBuzl5WUvzVS5XHardmBzpymez+NmZmZQKpW8veB1siCLDwyZGLajs2fPolAoRKaSIruvui4bVT8uildhpTIqqdHzbBCC3cpRodG+ZHI12vfzn/88AOCf/umfMD8/j3w+DwBOFxYQEBAQcHEQDNTvELVazblHdu7c6fSoq6ur3lagCwsLLjKPEymNAF1dWteDuldUp8r/dYWnhmaUG0XP5yoYGBrTdAUfP37cc1tsbGw4HaHqDwO2P9R4q9VqLhoUGBqQdPfv3LnTC6LodDqureZyuZHcp4zUz2QyXhS/7iQ1MzPjufW0XTOBt2o7bbStXczpM2lfGfe8fA49Lmp3KpUaANGBFzTOH3zwQZfPlAvP4NoPuFKh7J6V2HDhxUVgu90GAPc9YXXvnFdWVlZGpGzK0qo0TdPM2evZNE46H3KsiMpUQCInKkUcoeORjj9Ri2fVz9uUeITKq9LptLd9tI6vNnUVx18tvx2T+Jtq5K+//no88MADbny6GseiYKB+h+j3+17OOKbSWVxcxMbGhgsoosEKDNNFUOPD88blXtTJ2251xok8KkWIzSdnXRUnTpxw+SyPHDnibeWaTCadzoXHBFzZ6HQ6OHbsGPbt2wcAaLfbboGk2wMCQ92Zbtuprip1lZVKJayvr3sRt2y33J1G2UcNINSB1qZjs1vk2lRQquVS2BzBqptT/ag+D4O1bF0QlUoFDz30EADgsccec3kdu92uywEcEHAlQg0m9Sgkk8mRrXypV5+amvKMNRv/oO5//V1d+nSjqybc3l/nQ5vPVedI9eBoXIXVsOuGH7yHptlTA1H/1jLwnoRNgWUj9y1hZOUGVndL2PehUE/snj170Gq18PTTTwMYjlUcnyiZuNIRDNSAgICAgIDnEKweUw3JUqnkGXaDwcBJ0JSx1EUg4AdC0VC1i0tgVOKmQZHWMLMbcNhFqA2U1AWxllV1tTRyx20UoBlGKHfSwFB9Dq07/Z3yJ61XLacNGrXPPS740gY079ixAydPnnTPHBX4dSUjGKgXEadOnXJR74cPHwaw2fDOnj2LM2fOAACeeuop7NmzB9dffz0AeMnAdWDg+WR+6CZQd4iuTG00siYjBnzG6YEHHnAR16lUyqUKGgwGGAwGgTm9CtFqtVwbyGQyrj1yACZrXigU3AqcAzDZ1Wq16iad5eVlTE9Pu8koFou5SaBQKHiTh81RqG1VGVObeNq2d14L8N1mNrig1WpFTo6EppXSsqn+ttPp4Bvf+AZOnz4NYJhWji7OXq8X9NgBAQEBlxDBQL1EeOyxx7y/5+bmvImXkx4w1JJwVUT3R5TRSYOUE2ixWEStVovcScLqVdVNQeNjcnISwDCVECdx7owTcHWCWmNgc8u9tbU1LC8ve8YdPzNvKjVo/X7fyUH27t2Lubk5r63ZFb+mh1E9mE0XpdDUTlE5em0uVv7fbDYj9Wi8l+4WZV2BPG9lZcUtMiuVCur1upPjtNttp78NQVEBVyro3tb+qNtpptNpL4AylUq5323UuurFdcFIxlIXujY9nF5Hy6L9M2qnJKtRjUqUH+WS14WwTWVnd46yZVU5hD02iolVzT7gb7VMUkmfWcchfT6yrcBw3tasCbyWsrRELpcbkUBdiQgGakBAQEBAwFUO3QZbDSt1ac/NzaFQKDjygl4VaxTy/3FBunptwDcI6X7Xa2rWmiiPB2Fd82pc6n25kKUmU41YJXz4HGrcWQ2s3TLV7mpldbgshy709RzWmWpiNa4k6pn1swZ/ataRQqEwot+/0hEM1MuEpaUll1CXbA6ZoKWlJW/vc+10NtBJA6GslkfTStH9yGOnpqac2/7MmTPodrtulaysWsBzB7oJg+rO6vX6CBPKdnX48GEcOnQIwGZUKaGTBwdOuzmEXkthmQ+VtRCc5LTdc0Cm3ouTMOBPCoVCwWMwdMLRgIJKpeKkD4lEAsVi0d2j0WiEvhJwxUKzZigTp4E+c3NzyGazSKVSAIaSr16vFzkfWSZWEaWFZP9jVHzU1qOWEbUGngYJ2Wh4vW+32/WO5Xe8v8IauDYan99HPWOURpSwWQy0HMyOwmcadw2bSk/lVzyXsiNlfq+WiP5goF5GqPv82muvdZGRx48fdxpQNl4VaWvHUtfIww8/jHK5jPn5eQCjbpROp4Pv+Z7vAQDcddddLkfr5OQkut3uSMqQgOcurJ6SbYoswJ49ewAAhw4dGpEC6ALJMgU244R+b/ONErpAU1cY2zevyW1YeY4O+jZ4odls4qabbgIAHDhwAB//+McBDHeESiaTXgouPs/s7CxOnDjh5A3BOA0ICAi4fAgG6hbhySefdJPizMyMS0F1+vRpzM/PO9dEoVBwKyyyRMy7ury8jF6v5ybsYrHoGK3V1VXs2LHD5W188MEHsWvXLgCbKUO+/vWvX45HDbgCoSmU9u3b53Y6sUyoMiE25Yp1xxE2SEpTxTDNS1TUr93OVPMT8jPvv7i46KJbE4kEDh06hDvuuMOVUxGPx70tYBk8eM8994TcvwFXBbQfcrGo2kXqrLlw5G+xWAwAIlnNbrc71o1td1nSMnAM0QWqasnV5W09Kzai315fPTh6PuBHxtvUVVoX1mup11MXuu4oxevbHKk2fVVUWVWfz1RVNqsBf7PSiFwu5707vc7VgGCgBgQEBAQEXMWwwTe6Q5omv0+lUiNZLTY2NkaChngNe009T41QNYhpIKrxqPpUTRdFQ3dcSig1AtVYprEWtXsiz4lKuE9j2UqLiCjjOcrQth4cPZ71a5Pz63l6Dz5DoVBAtVr1ZFDJZBLT09MAhtI9vsdarXbe3SOvFAQDdQtBt+rk5KTTg/Z6PS9xPkXqwKbmj6taADh58qTTzGUyGXf8xMQEjh8/7higAwcOuHMefvjhS/dQAVcdjh07hu/+7u92f1shvjITGhRhAxh0QLc6UE0XZQdsK3fR4AYbxc9dVZ566ilvAlxaWsLf/d3fuevoM3z3d3+36ycnT5503oyAgICAgK1DbDAYDC7kwKuFMt6uYN5UpvmgEdpoNEYmZM3hODU15f6uVqvOsB0MBp7BYNNeBQQ8G7zhDW9wn6mJVkOQbZEBCsoiqCFr3XPqJrOpomwErg0a5HGPP/64S/80PT3ttKOtVisYnQHPadj0UN1u17mGi8Wi65vz8/MjAYWNRsNF9ZdKJS/NobKo2h87nQ46nY5LT8W/geHYoPOW3d6b1+L/dgyxKRTHsZ31et3TpBOa0pHn27RUyqiqPt4uzm26Lv1eg5Yph+D9dGy00gQd75rNpjfG1et1Jw2sVCpear3HH3/cyfcefvjhEVnBdsKFbihwwQzq3r173VaFVzoqlQqOHz++1cUIuIrw/Oc/f6uLACC07YCAgFGoPtHm0ux2u5iYmAAwqrns9/tIJBIjubiBTRezLi6tdpNQYoULW/WMqLtbAyPttaKyi9j0VSyblQroOdlsdiTtFp/XBnSqkcctXFlWNUJ1EWCzG9j6V8PXyga0PjSXKXOeq3GnnijdnCRKYnAl4oIN1EceecQ14suNEydO4Hd/93fxwAMP4Gtf+xpWV1fxgQ98AG95y1tGju10Onjve9+Lj3zkIzh58iR27dqFH//xH8e73vUuV/719XXccMMN22oiV4Zzfn7ea7xcFbXbbeRyOdx6660Ahq7XVCqF7/3e7wUA/J//83+cuz+Xy4XdoC4j7r///kt6/fe///349Kc/jQcffBBnzpzBnj178H3f931417ve5YLfgMvTtqMCF/g3B+Eo5oKDLZkWTcZvJygiKk2MThD87vTp03j88ccdE9RoNNz9rpZ9qQMCvhOogaj5M2OxmJOYsQ/bvezVKLUGkk0TxXsoa6pGnzVAlcHksZbBVCbSGoHj0tfZrU41mEtlRXZjAGs8j0vJxQwmUUY5z1Hj3epBo/S51WrVq1/duY+789EbtLKy4m1qks1msby8HFkXVyou2EDdKuMUAJ544gn8+Z//OW677Ta85jWvwV/8xV+MPfZNb3oT/uZv/gbvete78KIXvQj//M//jPe+97345je/ib/+678GMHyWcrm8rQxUhUZQR8Hu9vSVr3zlUhYnYBvg3e9+N175ylfife97H/bs2YNHH30U73nPe/CpT30KDzzwAHbu3Ang0rdtnSC497ONQOVnu0OMdbNpXkROMrVabWRnKZ2A4vG4t6vNN77xDQDDrBg6CQaXfkBAQMCVjSsiSOplL3sZzp49CwC47777xhqoX/7yl/HJT34Sv/d7v4df+ZVfAQC8+tWvRjKZxDve8Q783//7f/F93/d9l63cAQEXCw888ADm5ubc3y9/+cvxghe8AC960Yvwp3/6p3jnO995Wcpx11134Ud+5EcAjKZ90s+WeVHjkcaquhEJMhhqBGvgVb/f9wKhGNxEFiEYpgEB0YjSYQLD4FqrI9d+q67zqAWp9nkuHgE/Ab91aev9LaPKa2k5ozwpdhEMwNvExmYOsCwofxtXTj2W0EwBrIsottKmhNJ4EF5Hx0OWm6yz/m0zHFSrVXeepqG89tprXZB0rVa7KpL1PyMe+Atf+AJisVikgfjhD38YsVgM995770UrHHGhdPWXvvQlAMBrXvMa7/vXvva1AIBPfOITF7dgAQH/H61WC89//vNx6NAhrK6uuu8XFxcxPz+PV7ziFW5f+2cDNU6J22+/HYlEYtt6AgICArYPGNBoPRKxWMx9f/r0abRaLTSbTefOZ+AiDUKmqKKBxnOZZUOvz2NpSOVyOS/nsc2Vqptx6CJX76/H0phrNpveFq56Lv/xOM1RGvWP141KccVFuaaDUr0nr8EyqMyBx9HzRM1otVp1//gbwU17+K/RaCAWiyEWi2H37t0olUqujtPpNObn551E8Gpw8z8jBvWlL30pnv/85+OP//iP8aY3vcn77Y/+6I/wohe9CC960Yvc9mgXVIDkxSNx2RCo2SRCYvqAS41sNouPf/zjuP322/GTP/mT+MQnPoF+v49//a//NQaDAf7iL/4CiUTiovaNf/zHf0Sv18PNN998MR7hgkGpzA/+4A+OBBTYxPnKElim1Q7c/KzX1Kh+pofiJhb79+93ZQoyl4CAgICrC8/YOvzFX/xFvPWtb8WDDz6I2267DQBw77334t5778WHPvQhAMCHPvQhvPWtb72g611glqsLArcy/NKXvoRrrrnGff/FL34RAJyAOCDgUuC6667D+9//fvzYj/0Y/vAP/xArKyv43Oc+h7/7u79zgUwXq2+sra3h537u57B371785E/+5EUp/zPF3/3d3+G1r33tyA4xwGiaFbuaV6NUo1bj8binQ63X686Nn06n8epXvxqve93rAAwDCLiVb0BAwIVhamoKBw4ccP1K2cper4d6ve4tGAuFghdlrsGNmreYrCAwXFzncjmvX49z8duFbafT8XZHUre7utyjAo+iMgrwey23XfhqUJgtmzKydoe7ccn47XUBeG57DfDU+iALzEW41inLommydMzVgKpsNntV7IL3jA3UN73pTXj729+OP/7jP8af/umfAgD+y3/5L9ixYwd+7Md+DADwwz/8w8/Y1R/FLD1TdvWHfuiHcOjQIbz97W/Hzp078aIXvQhf/vKX8Y53vAOJROKqoLwDtjfe+MY34nOf+xz+w3/4D+j1enjHO97h6Z6fTd+waLVaeP3rX4+jR4/i7rvv9nRfAQEBAeeDBhkCQ6OKebeZ21QNLTVK1RNi85dqBo6oFEeqJQei00fRGLXufT2Of9NdrxpYNWzV8I4yLKPsC7r49Z61Ws17fith0L/VkFc3e7/fd9fgOWqwE3Z3LjVymVeW92g2m2i1Wmg0GgCG+nu+26jduq5EPGMDNZPJ4G1vext+7/d+D7/zO7+DjY0NfPzjH8ev/MqvOFd6qVRy229dKKKYpWfKrqbTafzt3/4tfuInfgLf//3fDwDI5/N43/veh/e85z3Ys2fPM7peQMCzwU/+5E/iT/7kT5BOp/GLv/iL3m/Ppm8o2u02/uW//Jf44he/iE9/+tO48847v9Pifkf49Kc/7TTfUftyRyXxV+aU32sKqn6/7zJZNBoNzM7OAhhuZlEoFHDLLbdc6scKCLhqkUwm0W633UYWAFyaqUwmg1wuN8Lu0dBS5pHf27Rz9l72N2uY2TyoysraICSb0srqP/WzDcy0EiOblF8/2+1MlRlWw14ZzHq97rG7drc8Gqi2HpQxpvSJ5SBLCmym+LLp81KplCsbg6SeswwqAPzbf/tv8R//43/E//yf/xOtVgvdbhc/+7M/635/Nm7Mi8EsAcChQ4fwz//8zzh58iRWVlZw7bXXYnV1Fb/0S7+El73sZd/x9QMCzodGo4Gf+ImfwOHDh3HmzBn89E//ND71qU+5378TF3+73cbrXvc6fPazn8WnPvUpl/92q/E3f/M37vOrX/1qAKNR/DZwwUbY6pa+586dc4PrxMSEi1r94Ac/eEmfIyAgICBg++BZGai7du3Cj/7oj+K//tf/ik6ngx/+4R/Gvn373O/PxticnZ11TMnFwJ49exxj+s53vhP5fB4/9VM/ddGuHxAQhZ/92Z/FsWPHcM899+CRRx7BG97wBvz+7/8+/v2///cAnv1CjMzp3XffjU9+8pP4gR/4gYtd9IuCf/iHf4j8/hWveIXnblO2YWVlxRmh1Wo1JNYPCLiEeOyxx3D48GHn1lcvB5lC3QZU00DZpPrKPuqOS3S92xRyPEfZRXVFM6WVSgWy2ezIxgG8npUiaLmUXdUtQa3EQGUCZFc5BilzynN0Byo9pl6vu+twu1hlf5VBVWmATWulwaVWZ6t63GKxiHg87hbziUTiinfpWzzrEPpf+qVfcu7FD3zgA95vF9vYBIC/+qu/ArCZpP6+++5zegvdH/y3f/u3MT8/j3379uHMmTP4+Mc/jrvuugsf+chHgos/4JLi/e9/Pz760Y/iAx/4AG6++WbcfPPN+Hf/7t/h7W9/O17ykpfgjjvueNZ94w1veAP+9m//Fr/xG7+B2dlZfPnLX3a/TU1NuQDBgICAgG+HdrvtpEatVsu5iavVKpLJpNvvnTpTzTeschzVv1s3vrrVmXKJnzXYR3enYr5Q1XxqXlBr2GrqKZUKMOWVbhWqsDpTlQboeeN21uLz63Nns1l33srKimewl0olz+XOXaEI3SWvWq1699cNTfR+1LHyPvpOuS3qlY7Y4DsIo7/mmmswMTGBhx9++GKWKRJc7UVBH+E3f/M38eEPfxgnTpzAxMQEXvziF+Md73gHXvrSl3rn3H777XjggQcuWXkDnlv42te+hjvvvBNvfOMbPVd0u93GS17yEiwvL+OBBx5wWqNnivO1/5e//OX43Oc+5/4ObTsgIOB82L9/vzNmVldXUS6XAQyzg5TL5REDlVB9Kg3U8wVEqVFlk+ZHbfVpDVQGCak8aCsMVDWYGWFvt15m/lhezxqoWk9qoCoT/GwMVG5OcubMGfdOT506haWlpZEsB9sFUW0mCs+aQf3617+OI0eO4I//+I+f7SWeES7Ujn7Xu96Fd73rXZe4NAEBPr7ru74rcjDMZDK47777vuPrX8x0bAEBAQE0oKanp5HP5wEAGxsbOHHihAugmpmZ8Vz1tVrNywJQr9fd32rMarJ6nqcspLqxrUu/VCp56el01yUanoRGvavxSuY1ygjldaJ2eeLxGrSk5yaTSef+p3FO49VKDprNprtOLpdz5/HZNX0Uy8JAU/7NHTR5XrPZdGXtdrvo9Xo4d+4cgOFcw2tmMpmrgkF9xgbqk08+iaNHj+Id73gHdu3ahbe85S2XoFgBAQEBAQEBlwJHjx7F9ddfD2Do5dEo/WQy6QwjRvWrsTMuLZJN3cQ0UPxN8xnzd8CPsFejEdg03tQoU2guUs0nStYzykC1kfLdbtfLH6oscTqddrtU8RpWr6qGr3rIlP3sdDruGSqVive3ppxiZgSeNzEx4RYAdmtpXosZGDQFFuBrgq9UPGMD9T3veQ8+8pGP4MYbb8T//t//eyR1QkBAQEBAQMD2BtNCdrtdrK+vu+/VyDxz5gwKhYJzHWseTuof1fW8srIycg37N4OHaKDRKOY19VgyqGpcaqqpcUFalBCoYanuf5sCSoOpVC/LZ1RZAa+j26by/lo3apTr/dLptBf8pGwyy0kZgd0kYWpqyj1/s9lEo9FwGVDq9TomJiYADAPhrgY8YwP1gx/84FWR7uX+++/f6iIEBFwShLYdELD1uBpcrAEBW4lnrUENCAgICAgIuDLx9a9/HQBw8803OwZ1cnISqVTKBWVms1kUCgXnZtYtROmOJvtXq9UcY1gsFj1mUpP91+t1dDodF4gVFfikQUP1et372zKomhJK3eHVajVyEwFlR/k32ddz585h7969XrltKiuym7x/VDQ+2VX+pts383yVIzDNHplXZYZZp61WC9ls1pW1Wq16+aOLxeJV4dZXBAM1ICAgICDgOYp6ve5c+PV6HZlMxhk+3ERDdZeas1Rd3GqscacmDa6iYUZZgKag4mcb2AQMA4xopFnjtd/vO+Nubm7Ou7/qNfWaNmVUoVDwjqNxTYzLFGB30pqamnIGMlNnaRouPkOlUsHU1JT3/FbCQKkEsJkRIZlMolqtYm1tzV0vm826+8disasuD2rwQQQEBAQEBAQEBGwrBAY1ICAgICDgOYqjR4/i8OHDADbd5Ezc3+12sbq66lz+jOoHhgxivV730ixpvlDAj9QnE7h7925vRyQbba9ufAZJaRomRS6Xc8xsPB53n+1uVRqUZFNH8Z68HtNn8W/LxGpg+Pl20tKUWGSUgaGMQMuqz8vzbLl4jY2NDccYF4tFrK+ve5kPrpbgKCIYqAEBAQEBAQHIZrOIxWKe0ZPJZFzaqcnJyfPuHkXDrtvtehHouVzOGVtMhs9z6/W6M/R4HI09zbHK36MkALyONZBVYqBQOYJG7duk+HTVMxm+ygEINcbVWAY2DUx1vc/MzKBer7t7FQoFdz/qc5muamVlxTPO2+029u/fD2C4wQIXEqzXqw1X3xMFBAQEBAQEXDDIvN16662eMUVmjwZTt9t1DF6pVPLYxmw26z43Gg0UCgVnsE1NTXmJ6TXXKLDJPna7XeRyOfc3A4PUwFTNazqddnpNZVqz2axnSKoelrpRZVR57PLysgvUAjaDlpTR1fsrE6zbjmazWaTTafeblrvT6WBqagonTpwA4OdvpeGu13n00UcBAAsLC5icnHTGsoUuHK4WBA1qQEBAQEBAQEDAtkJgUAMCAgICAgLw9a9/HTfddJNjQgEgkUi4BPDJZBKNRgPAJkvKqPJsNuu2SNVUTcDQxU7GkluEErp9p90ulCyluu5Ztk6ng06ng+XlZQBD+YFqUHWHJgVZVnXx856Tk5MolUoeK6pZBpRNrVarnuQgHo87ppcMMZlQ1ZKSdVW3PtHtdlEsFh1LXa1WcfPNN7t6On36tNs5inXOsrMeriYEAzUgICAgICDAgS799fV1b6epVCrlDLTBYIDFxUVn6DUaDaePpJtetxPVvKNqoKoBCPiGHvWqauDRWKVxNjs7CwBYXFx05aabXs/j/Wu1mpfftNPpOJnA+vr6/2vvfkKjuqI4jv+eeWMmmdiETFJDJaIu1IUEqlDEglBKoZZ2pVkKKnShuHCh4MY/CzdqN+66UVNQ3AhFSqFQF3VRFCxYSRcVxH+J2qaTOJP/Y+ZluhjO7X2TaBOdmEny/YAkkzfz3huDeObce86Jtbry21vZNW3vqD9Fy85v1xsYGIgtufvTogqFggs27bE91x+PKkktLS1uK0AikYgt/9uHhWKxKEnq6+vTYkOACgAAJEm9vb1qbW2VVMp2BkEQ64tqQVsqlVIqlXKB3tjYmP766y93Hn8fqF8IZEGnBaF+b1ELbP0qej/Q9INHq7i3jK6f+QzDcErlfvkYUr+4ye4liiJls1m1tbXFruePN/ULqvz3aAVV0n89Sv3sqt+jtaenxwWWYRi6IrTW1lZNTk4ql8tJKn0I8IukJiYmtHLlSklSPp9XJpPR48ePZ/JrXZDYgwoAAICqQgYVAABIKmUCLaMXhuGUtlLl04rseFNTk8t89vf3q7293S3H+xOnRkdHY/s+m5qaYkv8YRjGxonaz+yrLYFnMhnXG1QqtW/yJz6V9xD1e6IuW7Ys1rrKn86UzWbd9ZuammI9XP3epv5YVftq41vHx8dj07P8rQCFQkHFYlFRFLnzWMZaKrWPsuxqEARuSb9YLKqhocFljHt6emITpxajoGh/EwAAoCLK+2UuJBaQJZNJpdNppVIpSaUA0QKrmpoaRVHkArZUKuVeVygUFEWRG6Fa3j90YGDABXN+oY8VEJXv87Tz+sv/z549UyKRcEVDDQ0NWrNmjaSpraT8Jvq259XuZ3x8XH///bckaWJiwm0PkErjU/2m+plMxgWm1sTf7tVv4F8oFNTX1zdt8ZdtL7DrDw0NxbZQJBIJd79jY2OxJv65XM59ePD3qi40/h7k11m4/4IAAACwKLHEDwAAHMv2DQ4OanJy0lXxJ5NJl0EsFAoKgsBlMLPZrJtsVFtbqyAIYi2h/CyotY+yY+UTl/xtBFbUJJUyb3b9FStWKJvNTjmv8SdZNTQ0xJb//Uzs5OSky2BacZO9LpFIKJ1OuyIqv7jKMrR+YZTfOmpkZMQVO1kmWSoVNw0NDbmWXDU1NbHMcz6fd0v8dXV17j5zuZxrP7VUEKACAACnvBepBUnZbDa259Tfk1pXVxdbSi8Wiy7QGhkZiQWStbW1LujzA1C7ti1f29K4BX52Djt/Y2PjtJOsLHj0+5naNcIwVH19feycdr0oimIBY7FYjI1F9QNS20f74sULd0/2vEQiodra2lgLLtvnmkgkFASBamtr3ets24S18bKepplMZtHvM30d9qACAFBhC3kP6utY8FZfX69kMuke+8VUURTJDy38gLeurs6N9CxnAa9f0JRMJt25/Kzo8PBw7DptbW3uHkZHR2OtnZYvX+5+H/7PpFL/0D///FNSqSgplUq5wQSNjY2x9lF+AB6GofL5vMsg+/txi8Wi8vm8y8xaZtmOhWEYyzZbNnVoaEi9vb2xAq/FiD2oAAAAWJDIoAIAUGGLNYPqs2b59r1V5luW1FhrJf+xn0G17y2D6U+L8jOYExMT7jy259UmSX3wwQduqX5yclLLli1zGU+/aX95q6znz5+7PbZNTU3q7+93XQvstX6bKLuXfD4fG2cqyd2bvb8gCNwxW9JPJpMaGRlx18xms8pkMrFzL3YzfZ/sQQUAALPmtzoqFAou0Kqvr1dra6srDgrD0C2F+z1Npf8CPakU4BUKBddKaXh4WPl83gWBw8PD7nV+H1WpNBbUlvj9yVOSYkVZw8PDGhwcdMvq/taEJ0+eKJPJuNfZSFS/56m/pSEIglgbKAtQbS+r3Wt/f7/bqzo6Ouq2IJilEpjO1uL/iAcAAIAFhSV+AAAqbCks8ZfzG9P7LZn8ZvtW6GRFQn4GdXx8PFYkZROf7PHLly9jM+39xy0tLdq6dau7D38pf3x83P0+Xr58qcHBwViW1ubZ+4VV9ri8yb4t95cXMtm2AmPdAuz6lm2ebhrXUjPTjDEBKgAAFbYUA9S3FYZhLCC0x/4eUD/Q9IPCQqGg9vZ2SdKGDRtiFf/lQW4mk9GjR48klXq9+uNK/f2x1nfUAk0LvKX/tg2U3w/+HwEqAADzhAB1bqxevVqSpjStD8PQZTXT6bTa2tq0YsUKSaWMqgVFuVxOPT096uvrc6/1C51m6lWZVPw/2kwBAABgQSKDCgBAhZFBnV/lQwTs9/HkyZP5vC2IJX4AAOYNASowPZb4AQAAsCARoAIAAKCqEKACAACgqhCgAgAAoKqE//8UAAAwG+3t7WppaZnv23hrmUxGPT09830bWIKo4gcAoMLGxsZUV1f3zq/b29urb775Rnfu3NHdu3eVy+V08eJF7dmzZ9rnX79+XceOHdPdu3dVX1+vL7/8UmfOnNH7778vqfQ+Nm7cSJCKiqGKHwCAeTIfwakk3b9/X5cvX9by5cv1xRdfvPa5N27c0I4dO7Ry5Updu3ZN586d0/Xr1/Xpp58qn89LKr2PxZAJxsLDEj8AAIvE9u3b9c8//0iSfvvtN125cuWVzz1y5IjWr1+vq1evuvnza9eu1ccff6wLFy5o//797+SegemQQQUAYI4EQfDKP48ePar49WY6IODp06e6ffu2du/e7YJTSdq2bZvWr1+v77//vuL3BswGGVQAAObIzZs3Y4/Hxsa0e/duRVGk5uZmFYtFRVE0o3P5geTb+uOPPyRJHR0dU451dHTo119/rdi1gDdBgAoAwBzZunWr+z6KIu3cuVO5XE43btzQe++9p66uLu3du3dG56pkTXN/f78kqbm5ecqx5uZmdxyYLwSoAAC8AwcPHtSPP/6oH374QZs3b5YkffXVV7p9+/aszjNd1vVNs6tBEMzq58C7QoAKAMAcO3XqlL799ludP39en3/+uft5c3OzGhsbZ3Wu7777bkrWdbbZ1XQ6LUnTZkoHBgamzawC7xIBKgAAc6irq0vHjh3TyZMntW/fvtix6YLNV7Eg9E2yruU2bdokSeru7p7Sjqq7u9sdB+YLASoAAHPkp59+0tdff619+/bpxIkTU46/SbCZTqddBvRNrVq1Sh999JEuXbqkw4cPq6amRpJ069Yt3bt3T4cOHXqr8wNviwAVAIA58PDhQ3V2dmrdunXau3evbt26FTv+4YcfViTYLHf16lVJ0oMHDySV+qE2NDRIknbt2uWed/r0aX322Wfq7OzUgQMH1NfXp6NHj2rTpk0zzuoCc4VRpwAAzIFffvlFn3zyySuPP3z4UGvWrKn4dV9X4FT+X/7PP/+s48eP6/fff3ejTs+ePetGnUrSli1bdOfOnYrfJ5ammY46JUAFAACvRICKSpppgMokKQAAAFQVAlQAAABUFZb4AQAAUFXIoAIAAKCqEKACAACgqhCgAgAAoKoQoAIAAKCqEKACAACgqhCgAgAAoKoQoAIAAKCqEKACAACgqhCgAgAAoKr8C+PFQcjac98BAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 660x260 with 4 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "for name, brain,con in zip(temp_names, temp_brains, contrasts):\n",
    "    if name in include_pairs_with:\n",
    "        display = plotting.plot_anat(brain, display_mode='ortho', title=name, \n",
    "                  draw_cross=False, annotate=True, dim=con);\n",
    "        plt.savefig(name + '_template.png'); # save png\n",
    "        plt.draw()\n",
    "        #display.close()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 0d. Parallel or serial    \n",
    "Registration can be a time consuming process. If you have the resources, you may opt to parallelize some of the steps. You can configure that here."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [],
   "source": [
    "do_parallel = True  # Set to True for parallel processing\n",
    "par_numjobs = 5     # How many parallel processes? This is only used when do_parallel is True\n",
    "if do_parallel is True:\n",
    "    from joblib import Parallel, delayed"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Step 1: Linear registration\n",
    "* This is performed using [ANTs](http://stnava.github.io/ANTs/) **Affine** registration (12 degrees of freedom: 3 translations + 3 rotations + 3 scalings + 3 shears)\n",
    "* The _antsRegistration_ function is called through its [nipype wrapper](https://nipype.readthedocs.io/en/latest/interfaces/generated/interfaces.ants/registration.html#registration).\n",
    "* [Anatomy of an antsRegistration call](https://github.com/ANTsX/ANTs/wiki/Anatomy-of-an-antsRegistration-call) does a good job of explaining what the various parameters mean."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 1a. Define function"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [],
   "source": [
    "def register_linear(moving, target,\n",
    "                    moving_prefix, target_prefix,\n",
    "                    include_pairs_with):\n",
    "    \"\"\"Linearly registers moving image to target image\n",
    "\n",
    "    Parameters\n",
    "    ----------\n",
    "    moving: 3D NIFTI image file\n",
    "    target: 3D NIFTI image file\n",
    "    moving_prefix: string\n",
    "        one of the template names \n",
    "        is used in the naming of output files\n",
    "    target_prefix: string\n",
    "        one of the template names \n",
    "        is used in the naming of output files\n",
    "    include_pairs_with: list\n",
    "        list of templates to include in execution\n",
    "        will only run registration if 'moving' or 'target' is on this list\n",
    "        \n",
    "    \n",
    "    Return\n",
    "    ------\n",
    "    Linear transformation matrix\n",
    "        file named as '{moving_prefix}_to_{target_prefix}_0GenericAffine.mat'\n",
    "    Transformed image\n",
    "        moving image transformed into target image space\n",
    "        file named as: {moving_prefix}_in_{target_prefix}_linear.nii.gz\n",
    "    \"\"\"\n",
    "    if moving_prefix in include_pairs_with or target_prefix in include_pairs_with: \n",
    "        rgstr = ants.registration.Registration(\n",
    "               dimension=3,\n",
    "               float=False,\n",
    "               fixed_image=target,\n",
    "               moving_image=moving,\n",
    "               output_transform_prefix='{0}_to_{1}_affine_'.format(moving_prefix, target_prefix),\n",
    "               output_warped_image='{0}_in_{1}_linear.nii.gz'.format(moving_prefix, target_prefix),\n",
    "               initial_moving_transform_com=0,\n",
    "               winsorize_lower_quantile=0.05,\n",
    "               winsorize_upper_quantile=0.95,\n",
    "               interpolation='Linear',\n",
    "               use_histogram_matching=[True],\n",
    "               transforms=['Affine'],\n",
    "               transform_parameters=[(0.1,)],\n",
    "               metric=['MI'],\n",
    "               metric_weight=[1],\n",
    "               radius_or_number_of_bins=[32],\n",
    "               sampling_strategy=['Regular'],\n",
    "               sampling_percentage=[0.2],\n",
    "               number_of_iterations=[[1000, 500, 250, 100]],\n",
    "               convergence_threshold=[1e-6],\n",
    "               convergence_window_size=[10],\n",
    "               shrink_factors=[[8, 4, 2, 1]],\n",
    "               smoothing_sigmas=[[3, 2, 1, 0]],\n",
    "               sigma_units=['vox'])\n",
    "    \n",
    "        rgstr.run()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 1b. Apply function to register all unique template pairs"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\n",
      "Performing linear registration in parallel...\n"
     ]
    },
    {
     "ename": "OSError",
     "evalue": "No command \"antsRegistration\" found on host chris-mh16. Please check that the corresponding package is installed.",
     "output_type": "error",
     "traceback": [
      "\u001b[0;31m---------------------------------------------------------------------------\u001b[0m",
      "\u001b[0;31mRemoteTraceback\u001b[0m                           Traceback (most recent call last)",
      "\u001b[0;31mRemoteTraceback\u001b[0m: \n\"\"\"\nTraceback (most recent call last):\n  File \"/home/chris/anaconda3/envs/py3.7/lib/python3.7/multiprocessing/pool.py\", line 121, in worker\n    result = (True, func(*args, **kwds))\n  File \"/home/chris/anaconda3/envs/py3.7/lib/python3.7/site-packages/joblib/_parallel_backends.py\", line 620, in __call__\n    return self.func(*args, **kwargs)\n  File \"/home/chris/anaconda3/envs/py3.7/lib/python3.7/site-packages/joblib/parallel.py\", line 289, in __call__\n    for func, args, kwargs in self.items]\n  File \"/home/chris/anaconda3/envs/py3.7/lib/python3.7/site-packages/joblib/parallel.py\", line 289, in <listcomp>\n    for func, args, kwargs in self.items]\n  File \"/tmp/ipykernel_13807/957643870.py\", line 56, in register_linear\n    rgstr.run()\n  File \"/home/chris/anaconda3/envs/py3.7/lib/python3.7/site-packages/nipype/interfaces/base/core.py\", line 398, in run\n    runtime = self._run_interface(runtime)\n  File \"/home/chris/anaconda3/envs/py3.7/lib/python3.7/site-packages/nipype/interfaces/ants/registration.py\", line 1020, in _run_interface\n    runtime = super(Registration, self)._run_interface(runtime)\n  File \"/home/chris/anaconda3/envs/py3.7/lib/python3.7/site-packages/nipype/interfaces/base/core.py\", line 756, in _run_interface\n    % (executable_name, runtime.hostname)\nOSError: No command \"antsRegistration\" found on host chris-mh16. Please check that the corresponding package is installed.\n\"\"\"",
      "\nThe above exception was the direct cause of the following exception:\n",
      "\u001b[0;31mOSError\u001b[0m                                   Traceback (most recent call last)",
      "\u001b[0;32m/tmp/ipykernel_13807/232977590.py\u001b[0m in \u001b[0;36m<module>\u001b[0;34m\u001b[0m\n\u001b[1;32m      5\u001b[0m         \u001b[0mdelayed\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mregister_linear\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m      6\u001b[0m         \u001b[0;34m(\u001b[0m\u001b[0mbrains\u001b[0m\u001b[0;34m[\u001b[0m\u001b[0;36m0\u001b[0m\u001b[0;34m]\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0mbrains\u001b[0m\u001b[0;34m[\u001b[0m\u001b[0;36m1\u001b[0m\u001b[0;34m]\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0mnames\u001b[0m\u001b[0;34m[\u001b[0m\u001b[0;36m0\u001b[0m\u001b[0;34m]\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0mnames\u001b[0m\u001b[0;34m[\u001b[0m\u001b[0;36m1\u001b[0m\u001b[0;34m]\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0minclude_pairs_with\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m----> 7\u001b[0;31m         for names, brains in zip(combinations(temp_names, 2), combinations(temp_brains, 2)))\n\u001b[0m\u001b[1;32m      8\u001b[0m     \u001b[0mend_time\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mtime\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mtime\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m      9\u001b[0m     \u001b[0mprint\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m'Linear registration completed in  {0:.2f} minutes'\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mformat\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mend_time\u001b[0m \u001b[0;34m-\u001b[0m \u001b[0mstart_time\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m/\u001b[0m\u001b[0;36m60\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n",
      "\u001b[0;32m~/anaconda3/envs/py3.7/lib/python3.7/site-packages/joblib/parallel.py\u001b[0m in \u001b[0;36m__call__\u001b[0;34m(self, iterable)\u001b[0m\n\u001b[1;32m   1096\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m   1097\u001b[0m             \u001b[0;32mwith\u001b[0m \u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0m_backend\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mretrieval_context\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m-> 1098\u001b[0;31m                 \u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mretrieve\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m   1099\u001b[0m             \u001b[0;31m# Make sure that we get a last message telling us we are done\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m   1100\u001b[0m             \u001b[0melapsed_time\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mtime\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mtime\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m)\u001b[0m \u001b[0;34m-\u001b[0m \u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0m_start_time\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n",
      "\u001b[0;32m~/anaconda3/envs/py3.7/lib/python3.7/site-packages/joblib/parallel.py\u001b[0m in \u001b[0;36mretrieve\u001b[0;34m(self)\u001b[0m\n\u001b[1;32m    973\u001b[0m             \u001b[0;32mtry\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m    974\u001b[0m                 \u001b[0;32mif\u001b[0m \u001b[0mgetattr\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0m_backend\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0;34m'supports_timeout'\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0;32mFalse\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m--> 975\u001b[0;31m                     \u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0m_output\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mextend\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mjob\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mget\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mtimeout\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mtimeout\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m    976\u001b[0m                 \u001b[0;32melse\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m    977\u001b[0m                     \u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0m_output\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mextend\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mjob\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mget\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n",
      "\u001b[0;32m~/anaconda3/envs/py3.7/lib/python3.7/multiprocessing/pool.py\u001b[0m in \u001b[0;36mget\u001b[0;34m(self, timeout)\u001b[0m\n\u001b[1;32m    655\u001b[0m             \u001b[0;32mreturn\u001b[0m \u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0m_value\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m    656\u001b[0m         \u001b[0;32melse\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m--> 657\u001b[0;31m             \u001b[0;32mraise\u001b[0m \u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0m_value\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m    658\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m    659\u001b[0m     \u001b[0;32mdef\u001b[0m \u001b[0m_set\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mself\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mi\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mobj\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n",
      "\u001b[0;31mOSError\u001b[0m: No command \"antsRegistration\" found on host chris-mh16. Please check that the corresponding package is installed."
     ]
    }
   ],
   "source": [
    "if do_parallel: # parallel\n",
    "    print('\\nPerforming linear registration in parallel...')\n",
    "    start_time = time.time()\n",
    "    Parallel(n_jobs=par_numjobs, backend='multiprocessing')(\n",
    "        delayed(register_linear)\n",
    "        (brains[0],brains[1],names[0],names[1],include_pairs_with) \n",
    "        for names, brains in zip(combinations(temp_names, 2), combinations(temp_brains, 2)))\n",
    "    end_time = time.time()\n",
    "    print('Linear registration completed in  {0:.2f} minutes'.format((end_time - start_time)/60))\n",
    "        \n",
    "else: # serial\n",
    "    for names, brains in zip(combinations(temp_names, 2), combinations(temp_brains, 2)):\n",
    "        \n",
    "        # explicit here for clarity\n",
    "        mov_pref = names[0]\n",
    "        targ_pref = names[1]\n",
    "        mov_brain = brains[0]\n",
    "        targ_brain = brains[1]\n",
    "    \n",
    "        print('\\nLinearly registering {0} to {1} template...'.format(mov_pref, targ_pref))\n",
    "        start_time = time.time()\n",
    "        register_linear(mov_brain, targ_brain, mov_pref, targ_pref, include_pairs_with)\n",
    "        end_time = time.time()\n",
    "        print('Linear registration completed in  {0:.2f} minutes'.format((end_time - start_time)/60))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 1c. Visualize linear registration results"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "for names, brains in zip(combinations(temp_names, 2), combinations(temp_brains, 2)):\n",
    "    if names[0] in include_pairs_with or names[1] in include_pairs_with:\n",
    "        mov_pref = names[0]\n",
    "        targ_pref = names[1]\n",
    "        mov_brain = brains[0]\n",
    "        targ_brain = brains[1]\n",
    "        targ_index = temp_names.index(targ_pref)\n",
    "\n",
    "        display = plotting.plot_anat(targ_brain, display_mode='ortho',\n",
    "                            title='{0} edges on {1} volume'.format(mov_pref, targ_pref),\n",
    "                            draw_cross=False, annotate=False, dim=contrasts[targ_index])\n",
    "        brain = '{0}_in_{1}_linear.nii.gz'.format(mov_pref, targ_pref)\n",
    "        #display.add_edges(brain)\n",
    "        display.add_contours(brain)\n",
    "        plt.savefig('Linear_{0}_on_{1}.png'.format(mov_pref, targ_pref))\n",
    "        plt.draw()\n",
    "        display.close()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Step 2: Non-linear registration\n",
    "* Symmetric Diffeomorphic ([SyN](https://www.ncbi.nlm.nih.gov/pmc/articles/PMC2276735/)) registration is performed using [ANTs](http://stnava.github.io/ANTs/).\n",
    "* The starting point for this step is the transformed image from the previous (linear) step\n",
    "* The _antsRegistration_ function is called through its [nipype wrapper](https://nipype.readthedocs.io/en/latest/interfaces/generated/interfaces.ants/registration.html#registration).\n",
    "* [Anatomy of an antsRegistration call](https://github.com/ANTsX/ANTs/wiki/Anatomy-of-an-antsRegistration-call) does a good job of explaining what the various parameters mean."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 2a. Define function"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def register_SyN(moving, target,\n",
    "                 moving_prefix, target_prefix,\n",
    "                 include_pairs_with):\n",
    "    \"\"\"Non-linearly registers moving image to target image using ANTs SyN\n",
    "\n",
    "    Parameters\n",
    "    ----------\n",
    "    Same as register_linear function\n",
    "\n",
    "    Return\n",
    "    ------\n",
    "    Non-linear forward warp\n",
    "        file named as '{moving_prefix}_to_{target_prefix}_1Warp.nii.gz'\n",
    "    Non-linear backward warp (inverse of the above)\n",
    "        file named as '{moving_prefix}_to_{target_prefix}_1InverseWarp.nii.gz'\n",
    "    Transformed image\n",
    "        moving image transformed into target image space\n",
    "        file named as: {moving_prefix}_in_{target_prefix}_linear+SyN.nii.gz\n",
    "    \"\"\"\n",
    "    if moving_prefix in include_pairs_with or target_prefix in include_pairs_with: \n",
    "        rgstr = ants.registration.Registration(\n",
    "               dimension=3,\n",
    "               float=False,\n",
    "               fixed_image=target,\n",
    "               moving_image=moving,\n",
    "               output_transform_prefix='{0}_to_{1}_'.format(moving_prefix, target_prefix),\n",
    "               output_warped_image='{0}_in_{1}_linear+SyN.nii.gz'.format(moving_prefix, target_prefix),\n",
    "               initial_moving_transform_com=0,\n",
    "               winsorize_lower_quantile=0.05,\n",
    "               winsorize_upper_quantile=0.95,\n",
    "               interpolation='Linear',\n",
    "               use_histogram_matching=[True],\n",
    "               transforms=['SyN'],\n",
    "               transform_parameters=[(0.1,)],\n",
    "               metric=['CC'],\n",
    "               metric_weight=[1],\n",
    "               radius_or_number_of_bins=[4],\n",
    "               sampling_strategy=['Regular'],\n",
    "               sampling_percentage=[0.2],\n",
    "               number_of_iterations=[[500, 200, 50]],\n",
    "               convergence_threshold=[1e-6],\n",
    "               convergence_window_size=[10],\n",
    "               shrink_factors=[[8, 4, 2]],\n",
    "               smoothing_sigmas=[[3, 2, 1]],\n",
    "               sigma_units=['vox'])\n",
    "\n",
    "        rgstr.run()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 2b. Apply function to register all unique template pairs    "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "if do_parallel: # parallel\n",
    "    print('\\nPerforming non-linear registration in parallel...')\n",
    "    start_time = time.time()\n",
    "    Parallel(n_jobs=par_numjobs, backend='multiprocessing')(\n",
    "        delayed(register_SyN)\n",
    "        ('{0}_in_{1}_linear.nii.gz'.format(names[0], names[1]), brains[1], names[0], names[1], include_pairs_with) \n",
    "        for names, brains in zip(combinations(temp_names, 2), combinations(temp_brains, 2)))\n",
    "    end_time = time.time()\n",
    "    print('Non-linear registration completed in  {0:.2f} minutes'.format((end_time - start_time)/60))\n",
    "    \n",
    "else: # serial\n",
    "    for names, brains in zip(combinations(temp_names, 2), combinations(temp_brains, 2)):\n",
    "\n",
    "        # explicit here for clarity\n",
    "        mov_pref = names[0]\n",
    "        targ_pref = names[1]\n",
    "        mov_brain = '{0}_in_{1}_linear.nii.gz'.format(mov_pref, targ_pref)\n",
    "        targ_brain = brains[1]\n",
    "    \n",
    "        print('\\nNon-linearly registering {0} to {1} template...'.format(mov_pref, targ_pref))\n",
    "        start_time = time.time()\n",
    "        register_SyN(mov_brain, targ_brain, mov_pref, targ_pref, include_pairs_with)\n",
    "        end_time = time.time()\n",
    "        print('Non-linear registration completed in  {0:.2f} minutes'.format((end_time - start_time)/60))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 2c: Visualize non-linear registration results"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "for names, brains in zip(combinations(temp_names, 2), combinations(temp_brains, 2)):\n",
    "    if names[0] in include_pairs_with or names[1] in include_pairs_with:\n",
    "        mov_pref = names[0]\n",
    "        targ_pref = names[1]\n",
    "        targ_brain = brains[1]\n",
    "        targ_index = temp_names.index(targ_pref)\n",
    "\n",
    "        display = plotting.plot_anat(targ_brain, display_mode='ortho',\n",
    "                            title='{0} edges on {1} volume'.format(mov_pref, targ_pref),\n",
    "                            draw_cross=False, annotate=False, dim=contrasts[targ_index])\n",
    "        brain = '{0}_in_{1}_linear+SyN.nii.gz'.format(mov_pref, targ_pref)\n",
    "        #display.add_edges(brain)\n",
    "        display.add_contours(brain)\n",
    "        plt.savefig('Nonlinear_{0}_on_{1}.png'.format(mov_pref, targ_pref))\n",
    "        plt.draw()\n",
    "        display.close()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "toc-hr-collapsed": false
   },
   "source": [
    "## Step 3: Derive composite warps by combining linear and non-linear steps\n",
    "* This is done with the _antsApplyTransforms_ function, using its [nipype wrapper](https://nipype.readthedocs.io/en/latest/interfaces/generated/interfaces.ants/resampling.html#applytransforms).\n",
    "* **Caution:** the order with which transforms are passed is akin to linear algebra notation, so first transform is applied last\n",
    "* Both forward (templateA-to-templateB) and backward (templateB-to-templateA) warps are computed"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 3a. Define function"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def combine_transforms(moving, target,\n",
    "                       moving_prefix, target_prefix,\n",
    "                       include_pairs_with):\n",
    "\n",
    "    \"\"\"Combines linear transform with non-linear warp to generate a composite warp\n",
    "\n",
    "    Parameters\n",
    "    ----------\n",
    "    Same as register_linear function\n",
    "\n",
    "    Return\n",
    "    ------\n",
    "    Composite forward warp\n",
    "        file named as '{moving_prefix}_to_{target_prefix}_CompositeWarp.nii.gz'\n",
    "    Composite backward warp (inverse of the above)\n",
    "        file named as '{target_prefix}_to_{moving_prefix}_CompositeWarp.nii.gz'\n",
    "    \"\"\"\n",
    "    if moving_prefix in include_pairs_with or target_prefix in include_pairs_with: \n",
    "        forward = ants.ApplyTransforms(\n",
    "                input_image=moving,\n",
    "                output_image='{0}_to_{1}_CompositeWarp.nii.gz'.format(moving_prefix, target_prefix),\n",
    "                reference_image=target,\n",
    "                transforms=['{0}_to_{1}_1Warp.nii.gz'.format(moving_prefix, target_prefix),\n",
    "                            '{0}_to_{1}_affine_0GenericAffine.mat'.format(moving_prefix, target_prefix)],\n",
    "                dimension=3,\n",
    "                print_out_composite_warp_file=True)\n",
    "        forward.run()\n",
    "\n",
    "        backward = ants.ApplyTransforms(\n",
    "                input_image=target,\n",
    "                output_image='{1}_to_{0}_CompositeWarp.nii.gz'.format(moving_prefix, target_prefix),\n",
    "                reference_image=moving,\n",
    "                transforms=['{0}_to_{1}_affine_0GenericAffine.mat'.format(moving_prefix, target_prefix),\n",
    "                            '{0}_to_{1}_1InverseWarp.nii.gz'.format(moving_prefix, target_prefix)],\n",
    "                invert_transform_flags=[True, False], # invert the linear transform\n",
    "                dimension=3,\n",
    "                print_out_composite_warp_file=True)\n",
    "\n",
    "        backward.run()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 3b. Apply function to generate all composite warp pairs"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "if do_parallel: # parallel\n",
    "    print('\\nCreating composite warps in parallel...')\n",
    "    start_time = time.time()\n",
    "    Parallel(n_jobs=par_numjobs, backend='multiprocessing')(\n",
    "        delayed(combine_transforms)\n",
    "        (brains[0], brains[1], names[0], names[1], include_pairs_with) \n",
    "        for names, brains in zip(combinations(temp_names, 2), combinations(temp_brains, 2)))\n",
    "    end_time = time.time()\n",
    "    print('Composite warps computed in {0:.2f} minutes'.format((end_time - start_time)/60))\n",
    "    \n",
    "else: # serial\n",
    "    for names, brains in zip(combinations(temp_names, 2), combinations(temp_brains, 2)):\n",
    "\n",
    "        # make explicit here\n",
    "        mov_pref = names[0]\n",
    "        targ_pref = names[1]\n",
    "        mov_brain = brains[0]\n",
    "        targ_brain = brains[1]\n",
    "\n",
    "        print('\\nCreating composite warps between {0} and {1} templates...'.format(mov_pref, targ_pref))\n",
    "        start_time = time.time()\n",
    "        combine_transforms(mov_brain, targ_brain, mov_pref, targ_pref, include_pairs_with)\n",
    "        end_time = time.time()\n",
    "        print('Composite warps computed in {0:.2f} minutes'.format((end_time - start_time)/60))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 3c. Visually assess the derived composite warp pairs\n",
    "* First use the derived composite warps to transform brains in both directions"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "for names, brains in zip(combinations(temp_names, 2), combinations(temp_brains, 2)):\n",
    "    if names[0] in include_pairs_with or names[1] in include_pairs_with:\n",
    "        mov_pref = names[0]\n",
    "        targ_pref = names[1]\n",
    "        mov_brain = brains[0]\n",
    "        targ_brain = brains[1]\n",
    "\n",
    "        forward = ants.ApplyTransforms(\n",
    "                input_image=mov_brain,\n",
    "                output_image='{0}_in_{1}_composite.nii.gz'.format(mov_pref, targ_pref),\n",
    "                reference_image=targ_brain,\n",
    "                transforms=['{0}_to_{1}_CompositeWarp.nii.gz'.format(mov_pref, targ_pref)],\n",
    "                dimension=3)\n",
    "        backward = ants.ApplyTransforms(\n",
    "                input_image=targ_brain,\n",
    "                output_image='{1}_in_{0}_composite.nii.gz'.format(mov_pref, targ_pref),\n",
    "                reference_image=mov_brain,\n",
    "                transforms=['{1}_to_{0}_CompositeWarp.nii.gz'.format(mov_pref, targ_pref)],\n",
    "                dimension=3)\n",
    "\n",
    "        print('\\nTransforming {0} into {1} space...'.format(mov_pref, targ_pref))\n",
    "        start_time = time.time()\n",
    "        forward.run()\n",
    "        end_time = time.time()\n",
    "        print('Transformation complete in {0:.2f} minutes'.format((end_time - start_time)/60))\n",
    "        print('Transforming {1} into {0} space...'.format(mov_pref, targ_pref))\n",
    "        start_time = time.time()\n",
    "        backward.run()\n",
    "        end_time = time.time()\n",
    "        print('Transformation complete in {0:.2f} minutes'.format((end_time - start_time)/60))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "* Then, overlay the transformed (by composite warps) brains on target brains for visualization\n",
    "* two plots (forward/backward) will be generated for each template pair"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "for names, brains in zip(combinations(temp_names, 2), combinations(temp_brains, 2)):\n",
    "    if names[0] in include_pairs_with or names[1] in include_pairs_with:\n",
    "        mov_pref = names[0]\n",
    "        targ_pref = names[1]\n",
    "        mov_brain = brains[0]\n",
    "        targ_brain = brains[1]\n",
    "        mov_index = temp_names.index(mov_pref)\n",
    "        targ_index = temp_names.index(targ_pref)\n",
    "\n",
    "        display = plotting.plot_anat(targ_brain, display_mode='ortho',\n",
    "                            title='{0} edges on {1} volume'.format(mov_pref, targ_pref),\n",
    "                            draw_cross=False, annotate=False, dim=contrasts[targ_index])\n",
    "        brain = '{0}_in_{1}_composite.nii.gz'.format(mov_pref, targ_pref)\n",
    "        #display.add_edges(brain)\n",
    "        display.add_contours(brain)\n",
    "        plt.savefig('Composite_{0}_on_{1}.png'.format(mov_pref, targ_pref))\n",
    "        plt.draw()\n",
    "        display.close()\n",
    "\n",
    "        display = plotting.plot_anat(mov_brain, display_mode='ortho',\n",
    "                            title='{1} edges on {0} volume'.format(mov_pref, targ_pref),\n",
    "                            draw_cross=False, annotate=False, dim=contrasts[mov_index])\n",
    "        brain = '{1}_in_{0}_composite.nii.gz'.format(mov_pref, targ_pref)\n",
    "        #display.add_edges(brain)\n",
    "        display.add_contours(brain)\n",
    "        plt.savefig('Composite_{1}_on_{0}.png'.format(mov_pref, targ_pref))\n",
    "        plt.draw()\n",
    "        plt.show()\n",
    "        display.close()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 4. Clean-up\n",
    "* Create subfolders for linear and non-linear warps, and for the warped volumes\n",
    "* Put all files in sensible subfolders"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "os.makedirs(BASE_path + \"warps/linear\" , exist_ok=True)\n",
    "os.makedirs(BASE_path + \"warps/nonlinear\" , exist_ok=True)\n",
    "os.makedirs(BASE_path + \"warps/final\" , exist_ok=True)\n",
    "\n",
    "[os.rename(os.getcwd() + '/' + f, BASE_path + \"warps/linear/\" + f) for f in glob.glob('*_affine_*Affine.mat')];\n",
    "[os.rename(os.getcwd() + '/' + f, BASE_path + \"warps/nonlinear/\" + f) for f in glob.glob('*_to_*_1Warp.nii.gz')];\n",
    "[os.rename(os.getcwd() + '/' + f, BASE_path + \"warps/nonlinear/\" + f) for f in glob.glob('*_to_*_1InverseWarp.nii.gz')];\n",
    "[os.rename(os.getcwd() + '/' + f, BASE_path + \"warps/final/\" + f) for f in glob.glob('*_to_*_CompositeWarp.nii.gz')];\n",
    "[os.remove(f) for f in glob.glob('*Affine.mat')]; # spurious linear warps (generated by nonlinear)\n",
    "\n",
    "os.makedirs(BASE_path + \"warped_templates/linear\" , exist_ok=True)\n",
    "os.makedirs(BASE_path + \"warped_templates/nonlinear\" , exist_ok=True)\n",
    "os.makedirs(BASE_path + \"warped_templates/final\" , exist_ok=True)\n",
    "\n",
    "[os.rename(os.getcwd() + '/' +  f, BASE_path + \"warped_templates/linear/\" + f) for f in glob.glob('*_in_*linear.nii.gz')];\n",
    "[os.rename(os.getcwd() + '/' +  f, BASE_path + \"warped_templates/nonlinear/\" + f) for f in glob.glob('*_in_*linear+SyN.nii.gz')];\n",
    "[os.rename(os.getcwd() + '/' +  f, BASE_path + \"warped_templates/final/\" + f) for f in glob.glob('*_in_*composite.nii.gz')];\n",
    "\n",
    "os.makedirs(BASE_path + \"images/templates\" , exist_ok=True)\n",
    "os.makedirs(BASE_path + \"images/linear_reg\" , exist_ok=True)\n",
    "os.makedirs(BASE_path + \"images/nonlinear_reg\" , exist_ok=True)\n",
    "os.makedirs(BASE_path + \"images/all_warp_pairs\" , exist_ok=True)\n",
    "\n",
    "[os.rename(os.getcwd() + '/' + f, BASE_path + \"images/templates/\" + f) for f in glob.glob('*_template.png')];\n",
    "[os.rename(os.getcwd() + '/' + f, BASE_path + \"images/linear_reg/\" + f) for f in glob.glob('Linear*.png')];\n",
    "[os.rename(os.getcwd() + '/' + f, BASE_path + \"images/nonlinear_reg/\" + f) for f in glob.glob('Nonlinear*.png')];\n",
    "[os.rename(os.getcwd() + '/' + f, BASE_path + \"images/all_warp_pairs/\" + f) for f in glob.glob('Composite*.png')];"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 5. We will now create the average 'meta'-template in NMTv2.0 symmetric space (just because we can)\n",
    "* Get all the templates (in NMTv2.0 sym space) & average them into a 'MetaTemplate'\n",
    "* Plot the result"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "wtp = BASE_path + \"warped_templates/final/\"\n",
    "image_list = [NMTv20_sym_brain]\n",
    "for file in os.listdir(wtp):\n",
    "    if file.split('_')[2] == 'NMTv2.0-sym':\n",
    "        image_list.append(wtp + file)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "avg = ants.AverageImages(\n",
    "    dimension=3,\n",
    "    output_average_image=\"MetaTemplate_in_NMTv2.0-sym.nii.gz\",\n",
    "    normalize=True,\n",
    "    images=image_list\n",
    ");\n",
    "avg.run();\n",
    "\n",
    "plotting.plot_anat(\"MetaTemplate_in_NMTv2.0-sym.nii.gz\", display_mode='ortho', \n",
    "          title=\"MetaTemplate in NMTv2.0-sym\", draw_cross=False, annotate=True, dim=con);\n",
    "plt.savefig('MetaTemplate_in_NMTv2.0-sym.png')\n",
    "plt.draw()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "os.makedirs(BASE_path + \"warped_templates/metatemplate\" , exist_ok=True)\n",
    "[os.rename(os.getcwd() + '/' +  f, BASE_path + \"warped_templates/metatemplate/\" + f) for f in glob.glob('MetaTemplate*.nii.gz')];\n",
    "os.makedirs(BASE_path + \"images/metatemplate\" , exist_ok=True)\n",
    "[os.rename(os.getcwd() + '/' + f, BASE_path + \"images/metatemplate/\" + f) for f in glob.glob('MetaTemplate*.png')];"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3 (ipykernel)",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.7.12"
  },
  "toc-showcode": false,
  "toc-showmarkdowntxt": false
 },
 "nbformat": 4,
 "nbformat_minor": 4
}