{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "bf6256ae",
   "metadata": {},
   "source": [
    "## 12/07/2022\n",
    "## Registration with QC metrics\n",
    "### starts after tif conversion"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "id": "c1337f9f",
   "metadata": {},
   "outputs": [],
   "source": [
    "import matplotlib.pyplot as plt\n",
    "import matplotlib.patches as patches\n",
    "from pathlib import Path\n",
    "import numpy as np\n",
    "import os\n",
    "import sys\n",
    "import glob\n",
    "import pandas as pd\n",
    "import xml.etree.ElementTree as et\n",
    "import datetime\n",
    "from skimage.io import imread, imsave\n",
    "from imageio import volread as imread\n",
    "\n",
    "import tifffile\n",
    "import pystackreg\n",
    "from pystackreg import StackReg\n",
    "from skimage.filters import threshold_otsu\n",
    "\n",
    "from pystackreg.util import to_uint16   # make sure version 0.2.5 (not anything below)\n",
    "#from ims_to_tiff import convert_to_tif   "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "id": "dc1d0c16",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "/mnt/disks/store/101222_D10_Coverslip1_Processed\n"
     ]
    }
   ],
   "source": [
    "print(os.getcwd())"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "id": "9adc6d6d",
   "metadata": {},
   "outputs": [],
   "source": [
    "# Define input params, directories\n",
    "allfiles = sorted(glob.glob('tif/*'))\n",
    "allCycles = []\n",
    "for x in allfiles:\n",
    "    if x.startswith('tif/Cycle_'):\n",
    "        allCycles.append(x.split('_')[-1])\n",
    "\n",
    "_allFOVs = sorted(glob.glob('tif/Cycle_0/*'))\n",
    "allFOVs = []\n",
    "for x in _allFOVs:\n",
    "    allFOVs.append(x.split('F')[-1][0:3])\n",
    "    \n",
    "    \n",
    "for i in allCycles[1:]:\n",
    "    if i == '0':\n",
    "        continue\n",
    "    try:\n",
    "        os.makedirs(f'tmat_Cyc_{i}')\n",
    "        os.makedirs(f'reg_bin_Cyc_{i}')\n",
    "        os.makedirs(f'reg_Cyc_{i}')   \n",
    "    except: \n",
    "        continue"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "id": "8d4fc6fe",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "100\n"
     ]
    }
   ],
   "source": [
    "for i in range(len(allFOVs)):\n",
    "    if allFOVs[i] == '100':\n",
    "        print(i)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "id": "35a1f0ef",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "['100',\n",
       " '101',\n",
       " '102',\n",
       " '103',\n",
       " '104',\n",
       " '105',\n",
       " '106',\n",
       " '107',\n",
       " '108',\n",
       " '109',\n",
       " '110',\n",
       " '111',\n",
       " '112',\n",
       " '113',\n",
       " '114',\n",
       " '115',\n",
       " '116',\n",
       " '117',\n",
       " '118',\n",
       " '119',\n",
       " '120',\n",
       " '121',\n",
       " '122',\n",
       " '123',\n",
       " '124',\n",
       " '125',\n",
       " '126',\n",
       " '127',\n",
       " '128',\n",
       " '129',\n",
       " '130',\n",
       " '131',\n",
       " '132',\n",
       " '133',\n",
       " '134',\n",
       " '135',\n",
       " '136',\n",
       " '137',\n",
       " '138',\n",
       " '139',\n",
       " '140',\n",
       " '141',\n",
       " '142',\n",
       " '143',\n",
       " '144',\n",
       " '145',\n",
       " '146',\n",
       " '147',\n",
       " '148',\n",
       " '149',\n",
       " '150',\n",
       " '151',\n",
       " '152',\n",
       " '153',\n",
       " '154',\n",
       " '155',\n",
       " '156',\n",
       " '157',\n",
       " '158',\n",
       " '159',\n",
       " '160',\n",
       " '161',\n",
       " '162',\n",
       " '163',\n",
       " '164',\n",
       " '165',\n",
       " '166',\n",
       " '167',\n",
       " '168',\n",
       " '169',\n",
       " '170',\n",
       " '171',\n",
       " '172',\n",
       " '173',\n",
       " '174',\n",
       " '175',\n",
       " '176',\n",
       " '177',\n",
       " '178',\n",
       " '179',\n",
       " '180',\n",
       " '181',\n",
       " '182',\n",
       " '183',\n",
       " '184',\n",
       " '185',\n",
       " '186',\n",
       " '187',\n",
       " '188',\n",
       " '189',\n",
       " '190',\n",
       " '191',\n",
       " '192',\n",
       " '193',\n",
       " '194',\n",
       " '195',\n",
       " '196',\n",
       " '197',\n",
       " '198',\n",
       " '199',\n",
       " '200',\n",
       " '201',\n",
       " '202',\n",
       " '203',\n",
       " '204',\n",
       " '205',\n",
       " '206',\n",
       " '207',\n",
       " '208',\n",
       " '209',\n",
       " '210',\n",
       " '211',\n",
       " '212',\n",
       " '213',\n",
       " '214',\n",
       " '215',\n",
       " '216',\n",
       " '217',\n",
       " '218',\n",
       " '219',\n",
       " '220',\n",
       " '221',\n",
       " '222',\n",
       " '223',\n",
       " '224']"
      ]
     },
     "execution_count": 11,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "allFOVs[100:]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "id": "c4069683",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "['0', '8', '9']"
      ]
     },
     "execution_count": 12,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "allCycles"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "f539d4a9",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "cycle 9   FOV 100 \n",
      "Saving Registered Images... reg_Cyc_9/Cycle_9_F100_reg.tif\n",
      "cycle 9   FOV 101 \n",
      "Saving Registered Images... reg_Cyc_9/Cycle_9_F101_reg.tif\n",
      "cycle 9   FOV 102 \n",
      "Saving Registered Images... reg_Cyc_9/Cycle_9_F102_reg.tif\n",
      "cycle 9   FOV 103 \n",
      "Saving Registered Images... reg_Cyc_9/Cycle_9_F103_reg.tif\n",
      "cycle 9   FOV 104 \n"
     ]
    }
   ],
   "source": [
    "# # loop through registration\n",
    "for c in allCycles[2:]:    \n",
    "    for fov in allFOVs[100:]: \n",
    "        # read images\n",
    "        ref_name = f'tif/Cycle_0/Cycle_F{fov}.tif' \n",
    "        mov_name = f'tif/Cycle_{c}/Cycle_{c}_F{fov}.tif'\n",
    "        ref = imread(ref_name).astype(np.uint16)\n",
    "        mov = imread(mov_name).astype(np.uint16)\n",
    "        \n",
    "        print(f'cycle {c}   FOV {fov} ')\n",
    "        \n",
    "        ref_max = ref.max(0)\n",
    "        ref_binary = ref_max[0] > threshold_otsu(ref_max[0]) # nuclei\n",
    "        mov_max = mov.max(0)\n",
    "        mov_binary = mov_max[0] > threshold_otsu(mov_max[0]) # nuclei\n",
    "        \n",
    "        sr = StackReg(StackReg.RIGID_BODY) # create instance of stackreg class\n",
    "        tmat = sr.register(ref_binary, mov_binary) # calculate tmat\n",
    "        out = sr.transform(mov_binary)\n",
    "        out = pystackreg.util.to_uint16(out) # make sure to convert later into bool. Save first then convert.\n",
    "        \n",
    "        # save binary images\n",
    "        fname_to_save = f'reg_bin_Cyc_{c}/Cycle_{c}_F{fov}_bin_reg.tif'\n",
    "        tifffile.imwrite(fname_to_save, out, imagej=True, photometric = 'minisblack', metadata={'axes':'YX'})\n",
    "        \n",
    "        # save tmat\n",
    "        np.save(f'tmat_Cyc_{c}/Cycle_{c}_F{fov}_tmat.npy', tmat)\n",
    "        \n",
    "        # REGISTRATION\n",
    "        reg = np.zeros(mov.shape, dtype=np.uint16)\n",
    "        for Z in range(mov.shape[0]):\n",
    "            for ch in range(mov.shape[1]): \n",
    "                reg[Z,ch,...] = sr.transform(mov[Z,ch,...], tmat=tmat)\n",
    "                reg[Z,ch,...] = to_uint16(reg[Z,ch,...])\n",
    "\n",
    "        fname_to_save = f'reg_Cyc_{c}' + f'/Cycle_{c}_F{fov}_reg.tif'\n",
    "        print('Saving Registered Images...', fname_to_save)\n",
    "        \n",
    "        # save registered images\n",
    "        tifffile.imwrite(fname_to_save, reg, imagej=True,\n",
    "                         photometric = 'minisblack',metadata={'axes':'ZCYX'})\n",
    "        "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "1d14a40d",
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "9e883b67",
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "id": "219a37ed",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAABJIAAAI4CAYAAAAxusOtAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/Il7ecAAAACXBIWXMAAAsTAAALEwEAmpwYAABjiklEQVR4nO3dd7wcVd3H8c9JDwkQAtISegcFpYmoiBSlNxXwoYkoNqQoKPqg2FCxgYoNUUBEqgiIjwoKioqiKE1AqkAInZBQ027O88eZm2xubpm7O7Oz5fN+va53d3b2zC93ifne35w5E2KMSJIkSZIkSUMZUXUBkiRJkiRJag82kiRJkiRJkpSLjSRJkiRJkiTlYiNJkiRJkiRJudhIkiRJkiRJUi42kiRJkiRJkpSLjSRJC4UQPhBCeCKE8EIIYfmq66lXVv/aFR37jSGEu6s4tiRJah8hhNeHEO7Ncss+VdczHFXnnRDCJ0MIZ1V1fKnbhRhj1TVIKkgI4UFgJaAHmAfcALw/xjgtx3tHA88B28QYby2zzux4fwC2AeYDs4HrgQ/FGB8r+9jZ8bcHfhpjnNrgOBFYL8Z4XxF1SZKkzpJlns2AlWOMc2q2/x64Msb4zex54ZkiG/MlIAKzgIuAE2KMPUUdY4jjfwZYN8Z4cANjbE8BmU1ScZyRJHWePWOME4FVgCeAb+d830rAOOCO4R4wJPX8/8lRWa3rAhOBr9UxRmlCCKOqrkGSJLWvEMKawBtJjZy9+ry8BnXkrgGOM1hm2SzLW28CDgDeXcQxi9BAhpRUIf/SSh0qxjgbuBTYuHdbCGFsCOFrIYSHs0vYvh9CGB9CWB/onZ48M4Rwbbb/tiGEf4QQZmXft60Z6w8hhFNCCH8hnelaO4SwYQjhmhDCjBDC3SGE/XPWOhO4HHh1zfgDjhVCWD6E8MsQwnNZXV8IIfy55vUYQlg3e7xbCOHOEMLzIYTpIYTjQwgTgF8Dq2bTyV8IIawaQvhMCOHSEMJPQwjPAe8KIWwdQvhrCGFmCOGxEMIZIYQx2djXZ4e8NRvjgBDC9iGER2pq2Sj7Wc0MIdwRQtir5rVzQgjfCSH8KqvvxhDCOnl+ZpIkqS0cCvwNOAc4rHdjCOF+YG3gl1mG+Gv20sJMke23RwjhlixH3BBC2LRmjAdDCB8PIdwGvDjUCbBsptNfWDxvDTb+5iGEm7OMckkI4aIQwhey1/rmnY9nOev5LLftGELYBfgkcED2Z7o127e/DHl4COGu7P0PhBDel+07WGb7ac3x98py1sxs/I36/JyODyHclmXai0II44b85CQNyEaS1KFCCEuRzjr9rWbzl4H1SQFiXWAK8OkY4z3AJtk+k2KMO4QQJgO/Ar4FLA98A/hVWHztpEOAI4GlgaeAa4CfASsCBwLfDSFszBCyMfcD7sueTxhirO8ALwIrk0LZYX3HrPEj4H0xxqWBVwLXxhhfBHYFHo0xTsy+Hs3235vUgJsEnE+6TPA4YAXgdcCOwAcBYozbZe/ZLBvjoj5/rtHAL4Grsz/Hh4HzQwgb1Ox2IPBZYLnsz3/KUD8vSZLUNg4l5YnzgbeGEFYCiDGuAzxMNpM8xvi6bP+FmSKE8Brgx8D7SFnsB8CVIYSxNeO/E9idlN/mD1ZICGFD0uyo3rw14PjZSbNfkBpgk4ELgH0HGHcD4ChgqyxvvRV4MMb4G+CLwEXZn2mzmrfVZsiHgCeBPYBlgMOB00IImw+R2XqPv35W37HAK4D/IzXoxtTstj+wC7AWsCnwrsF+VpIGZyNJ6jyXhxBmkq6D3xn4KqSpw6R/sI+LMc6IMT5P+sf9wAHG2R24N8Z4XoxxfozxAuA/wJ41+5wTY7wjCy67kELD2dn+NwM/B94xSK3fCiHMAp4mNWo+nG3fY6CxQggjgbcBJ8cYX4ox3gmcO8gx5gEbhxCWiTE+G2P81yD7Avw1xnh5jHFBjPHlGOM/Y4x/y+p4kBSy3jTEGL22IV2y9+UY49wY47XAVaTQ1+sXMca/Zz/D86k5SyhJktpXCOENpMvXLo4x/hO4H/ifYQxxJPCDGOONMcaeGOO5wBxSvuj1rRjjtBjjy4OM868QwovAXcAfgO/mGH8bYFQ2/rwY42XA3wcYvwcYS8pbo2OMD8YY7x/iz7YwQ2bj/yrGeH9M/kg6CffGIcbodQDwqxjjNTHGeaSlEsYD29bs860Y46Mxxhmkk3yvzjm2pH7YSJI6zz4xxkmk9Y6OAv4YQliZdIZmKeCf2bTfmcBvsu39WZV0hqjWQ6RZTL1qF/FeA3ht79jZ+AeRZg0N5OgY47KkM0PLAb2LKA421itIwab22IMtJv42YDfgoRDCH0MIrxtk3yXGCiGsH0K4KoTweEiXu32R1PTKY1VgWoxxQc22vj/Dx2sev0RqPEmSpPZ3GHB1jPHp7PnPGHwWdV9rAB/tk4dWI+WLXkPeUAXYnJQvDgBeC0zIMf6qwPS4+J2Z+j1WdsncscBngCdDCBeGEFbtb9+Bxgoh7BpC+FtISxrMJGW34eSthZk1y13TMG9JpbGRJHWo7MzSZaSzRG8gzfp5Gdgkxjgp+1o2W3yxP4+SAkat1YHptYepeTwN+GPN2JOy6ccfyFHr7cAXgO9kM6cGG+sp0p3eau/csdogY/8jxrg36dKyy4GL+6l9sbf0ef490kys9WKMy5Cu9Q9D/ZkyjwKrhcUXkez7M5QkSR0mhDCedDnVm7KTUY+TLpXfLISw2eDvXmgacEqfPLRUNku8V65bcGczfS4G/gp8Osf4jwFTslzWa7C89bMYY+8MrAicOkR9C7dnl+r9nDSTaKXshOj/sShvDfVnXCyzZjWvhnlLKo2NJKlDhWRv0kyfu7KzMz8kXXO+YrbPlBDCWwcY4v+A9UMI/xNCGBXSoo8bky7N6s9V2f6HhBBGZ19b1S52OIRzSXeO22uwsWK6Xe1lwGdCCEtl1/sfOsDPYEwI4aAQwrLZVOfngN7ZQU8Ay4cQlh2irqWz972QHatvY+wJ0mKZ/bmRdNbrY9mfYXvSpYEXDnFMSZLU3vYhnczbmHQZ1auBjYA/MUBuYclM8UPg/SGE12a5bkIIYfcQwtIN1PVl4L3ZbPXBxv9rVv9RWQ7cG9i6vwFDCBuEEHbIGkKzSScua/PWmmHwO7ONIV0a9xQwP4SwK/CWmteHymwXA7uHtMD3aOCjpEv0bhjqhyGpPjaSpM7zyxDCC6TmxynAYTHG3lvLfpy0wOLfssu0fgds0N8gMcZnSGsVfRR4BvgYsEfN9Oy++z9P+kf/QNKZocdJZ6PG9rd/P++fC3wT+FSOsY4Cls22n0daYHHOAEMfAjyY/XnfT7pEjhjjf7L3PZBN5x5oCvbxpPUMnicFrov6vP4Z4NxsjMXuUpf9mfYkLRL5NGlNgkOzY0uSpM51GHB2jPHhGOPjvV/AGcBBof87rH2GmkwRY7wJeG/2nmdJGe5djRSVzQK/HjhhsPGzDLMfcAQwEziYdKKvv7w1ltSgepqUzVYEPpG9dkn2/ZkQQr/rVGa572hSQ+hZUu66sub1QTNbjPHurL5vZzXsSVrEfG6OH4mkOoTFL3uVpPYTQjgVWDnGOJx1ByRJkpRTCOFG4PsxxrOrrkVStZyRJKnthBA2DCFsmk3D3pp0tuwXVdclSZLUKUIIbwohrJxd2nYY6eYov6m6LknV629KpSS1uqVJU5xXJV03/3XgikorkiRJ6iwbkC43mwA8ALw9xvhYtSVJagUtdWlbCGEX0hopI4GzYoxfrrgkSZKkjmcGkyRJebVMIymEMBK4B9gZeAT4B/DOGOOdlRYmSZLUwcxgkiRpOFppjaStgftijA9kK+xfCOxdcU2SJEmdzgwmSZJya6U1kqYA02qePwK8tu9OIYQjgSOzp1s0oS5JklShGGOouoYOZwaTJElLGCiDtVIjKZcY45nAmQAhhNa4Lk+SJKnDmcEkSRK01qVt04HVap5PzbZJkiSpPGYwSZKUWys1kv4BrBdCWCuEMAY4ELiy4pokSZI6nRlMkiTl1jKXtsUY54cQjgJ+S7r17I9jjHdUXJYkSVJHM4NJkqThCDG27yXuXp8vSVLnc7Ht1mMGkySp8w2UwVrp0jZJkiRJkiS1MBtJkiRJkiRJysVGkiRJkiRJknKxkSRJkiRJkqRcbCRJkiRJkiQpFxtJkiRJkiRJysVGkiRJkiRJUrcbAUzLHm86+G6SJEmSJEnqdlOAh4DRA+9iI0mSJEmSJKnbjAPuBEYB9wH3AgFYffC3jSq7LkmSpDxGAbcDC4BNKq5FkiSp480B9odR/4bb18mfwWwkSZKklhCADYEI3JZtO7C6ciRJkjpbBO6FsMHwMpiNJEmSVLkxwI3Z4wC8Kns8vppyJEmSOtdE4I/p4ZgRw89gNpIkSVKlxgJ/Al5dcR2SJEldYSSwef0ZzMW2JUlSZcYD1wJbVV2IJElSN1ga+E1jGcxGkiRJqsxc4HNVFyFJktQtRgPbNJbBbCRJkqTKjAaOr7oISZKkbvEcsF9jGcxGkiRJqswoYKeqi5AkSeoW44AjG8tgNpIkSVIllgLOr7oISZKkbjIGltqlsQxmI0mSJFViHnBZ1UVIkiR1i2WA02DeuxvLYDaSJElSJeYBP6+6CEmSpG4xB/gjzHtLYxnMRpIkSarMbOCDVRchSZLUDeYAlwMHNpbBbCRJkqTKzAfOqboISZKkLtNIBrORJEmSJEmSpFxsJEmSJEmSJCkXG0mSJKlS84DPAj3AyRXXIkmS1NFeAk5JD+vNYDaSJElSpeYDX6q6CEmSpG4wG/g8cGr9GcxGkiRJqtx84GvARIAAnFBpOZIkSZ1rDvCF9HCxDJaTjSRJklS5HuBE4BO9G1YBjgFWrKoiSZKkDjYXOKOfDJaDjSRJktQyIvC9CHwUOB1YrdJyJEmSOtNc4CPAD9LTCHwv51ttJEmSpJaxIMCHjoCzjqi6EkmSpA43jzQDHFgAfAg4K8fbRpVXkSRJ0jCMAA6BeBa8FxhTdT2SJEmdrgf4KXBwmpWUJ4M5I0mSJLWGkcA5i54eVlUdkiRJ3WI+8C7gokWbhspgNpIkSVL1RgD7LLl5ZpPLkCRJ6jo9wP8Av1i0aeYgu9tIkiRJ1RtJv6e/7m96IZIkSV1oAfA24Kr0dLAMZiNJkiRVrwf4btVFSJIkdbEI7JV9f27g3WwkSZKk6o0GflV1EZIkSeI64N6BX7aRJEmSqjUC2Ab4S9WFSJIkdbkI7Dj4LjaSJElStUYDpwAfq7oQSZIkDWVU1QVIkqQuN4d05mt21YVIkiRpKM5IkiRJkiRJUi42kiRJUvUicHfVRUiSJGkoNpIkSVL15pJuNytJkqSWZiNJkiRVbwxwTdVFSJIkaSg2kiRJUvXmAusAj1ZdiCRJkgZjI0mSJLWG+cDqVRchSZKkwdhIkiRJreWZqguQJEnSQGwkSZKk1tEDrAjMrLgOSZIk9ctGkiRJai0LgMlVFyFJkqT+2EiSJEmt6cWqC5AkSVJfNpIkSVLricDE7LskSZJaRoixfRNaCKF9i5ckSbnEGEPVNWhxZjBJkjrfQBnMGUmSJEmSJEnKxUaSJEmSJEmScrGRJEmSJEmSpFxsJEmSJEmSJCkXG0mSJEmSJEnKxUaSJEmSJEmScrGRJEmSJEmSJAjAAmCLgXexkSRJkiRJktStRpOaRwuAHlIzaRCjyq9IkiRJkiRJLWuI5lEtZyRJkiRJkiR1i2WAmcB4YD7w8vDebiNJkiS1hNGkHDMCmJd9bV5pRZIkSR3oOWCN7PtIGD1yeBnMRpIkSarMBNIJsbHAi9n32aRr70cxrFnWkiRJymMZmPAUzBxVXwZzjSRJklSZQJpdPYs0I4ma75IkSSpBgDC6/gzmjCRJklSJZYAnSc2ksRXXIkmS1C0azWA2kiRJUlOtADxMCi/jK65FkiSpW6ywAjz8PITpjWUwL22TJElNFYAppGaSJEmSmiMEmDKx8QzmjCRJktR0I0jTqiVJktQ8RWQwG0mSJEmSJEmdbCXgT8D6jQ9lI0mSJDXVM8BGVRchSZLUTZ6CZ/aFjf7U+FA2kiRJUlMtB/y+6iIkSZK6yQJY7m74/U6ND2UjSZIkNdWzwG5VFyFJktRlnp0Eu/2q8XFsJEmSpKaaBPy86iIkSZK6zKSR8PPVGx/HRpIkSWqqWcAhVRchSZLUZWbNgEO2bXwcG0mSJKmpeoC7qy5CkiSpy/T0wN0FhDAbSZIkqelmAVsC21VdiCRJUheZNQu23L6xDDaqqGIkSZLy6gH+CUysuhBJkqQu0tMD/7ylsQzmjCRJkiRJkqRu8QJpStJb6nu7M5IkSVJlXgZ2qnk+cgT89tfALsDVwAcrKUuSJKlz9cDLf4Kdxi/aNBL4bc6320iSJEmV6QF+nz0eCVwRYfdvwq9+ReowLVNZaZIkSR2rB/j9HGBvGHkFXAHsDvwqx3u9tE2SJLWEBcBZAd73PmDXqquRJEnqcAuAq2DB3nDWfvC+nG8LMcYyyypVCKF9i5ckSUsII+Ft8yEAFwHLbwkzboqh6rq0ODOYJEmdJQR42zsgXDR0BvPSNkmS1DIicGn2eDTwbIW1SJIkdYsY4dJLgQijDx48g9lIkiRJ1RsFnEmaipT5WVW1SJIkdaMFwOHws58Byw+8m5e2SZKk6o0FZvezfUuIXtrWcsxgkiR1vhj7z2Auti1JkiRJkqRcbCRJkqRqjQa+AHy86kIkSZI0FBtJkiSpWhGYAUyuuhBJkiQNxUaSJEmq1nzgK8DcqguRJEnSUGwkSZKk1mAjSZIkqeXZSJIkSdUbBXy+6iIkSZI0FBtJkiSpej3Ad6ouQpIkSUMprZEUQvhxCOHJEMK/a7ZNDiFcE0K4N/u+XLY9hBC+FUK4L4RwWwhh87LqkiRJLWg+cFzVRXQGM5gkSSpTmTOSzgF26bPtROD3Mcb1gN9nzwF2BdbLvo4EvldiXZIkqRUtIKUHNeoczGCSJKkkpTWSYozXk27mW2tv4Nzs8bnAPjXbfxKTvwGTQgirlFWbJElqQT3Ae6ouov2ZwSRJUpmavUbSSjHGx7LHjwMrZY+nANNq9nsk27aEEMKRIYSbQgg3lVemJEmqRAQuyb5UJDOYJEkqxKiqDhxjjCGEWMf7zgTOBKjn/ZIkqYUtAPYHAvCOimvpUGYwSZLUiGbPSHqid7p09v3JbPt0YLWa/aZm2yRJUjeKwK+rLqKjmMEkSVIhmt1IuhI4LHt8GHBFzfZDszuHbAPMqpl+LUmSutFuwHNVF9ExzGCSJKkQIcZyZiaHEC4AtgdWAJ4ATgYuBy4GVgceAvaPMc4IIQTgDNIdRl4CDo8xDnn9vdOqJUnqfDHGUHUN7cQMJkmSijBQBiutkdQMhhhJkjqfjaTWYwaTJKnzDZTBmn1pmyRJkiRJktqUjSRJkiRJkiTlYiNJkiRJkiRJudhIkiRJkiRJUi42kiRJkiRJkpSLjSRJkiRJkiTlYiNJkiRJkiRJudhIkiRJkiRJUi42kiRJkiRJkpSLjSRJkiRJkiTlYiNJkiRJkiRJudhIkiRJkiRJUi42kiRJkiRJkpSLjSRJkiRJkiTlYiNJkiRJkiRJudhIkiRJkiRJUi42kiRJkiRJkpSLjSRJkiRJkiTlYiNJkiRJkiRJudhIkiRJkiRJUi42kiRJkiRJkpSLjSRJkiRJkiTlYiNJkiRJkiRJudhIkiRJkiRJUi42kiRJkiRJkpSLjSRJkiRJkiTlYiNJkiRJkiRJudhIkiRJkiRJUi42kiRJkiRJkpSLjSRJkiRJkiTlYiNJkiRJkiRJudhIkiRJkiRJUi42kiRJkiRJkpSLjSRJkiRJkiTlYiNJkiRJkiRJudhIkiRJkiRJUi42kiRJkiRJkpSLjSRJkiRJkiTlYiNJkiRJkiRJudhIkiRJkiRJUi42kiRJkiRJkpSLjSRJkiRJkiTlMqrqAiSp1b0X2KvPtiuAsyqoRZIkqVuYwaTWZCNJkvpxNPCm7PGmwLp9Xt8QCMAPm1mUJElShzsaeNMHgR37ZLBLgNmw4XIQzjaDSVWykSRJ/dgK2G+Q19cFTiA1k85sSkWSJEmdbytgv81ZMohtAMyHdcfBCQHCj81gUlVCjLHqGuoWQmjf4iW1pKOBLYA3AGvn2P9e4GsYZKQyxRhD1TVocWYwSUU7+mjYYossg23LktPBa90N994IX/sLnGkIk0ozUAZzRpIkZY4GjgPWHMZ71gOOzx6bYyRJkobv6KPhuONgzTVzvmEDWG+DLIOdaQaTms1GkiQBHwY+AqxRx3t7m0kRr9eXJEkajg9/GD7yEVijjhBmBpOqMaLqAiSpFexJfU2kXusBuxRUiyRJUrfYc8/6mkgAbAnrHWoGk5rNRpKkrnckqREkSZKk5mk4g20M7FpMLZLys5Ekqau9B/gkw1sXaSAbAwcVMI4kSVKnew/wyXfDmus0PpYZTGou10iSJEmSJDXXu4CTgLUqrkPSsDkjSVLXOgSYDPQUNN584OWCxpIkSepUhwCT3ws9BTWRzGBSc9lIktS1VgCOAdYuaLxXAocXNJYkSVKnMoNJ7c1GkqSu9TAwu+oiJEmSuszDwOzLgOlVVyKpHjaSJHWtDYEJVRchSZLUZTYEJtwNvFDAYGsDOxYwjqTcbCRJ6lqnAPdXXYQkSVKXOQW4/xPABgUM9grSbdskNY2NJEmSJElSc/0OeLKAcW4Evl3AOJJys5Ekqav9AXiqoLGmAzcUNJYkSVIn+8PJ8NSdxYxlBpOay0aSpK72v8AdBY31EOnkmiRJkgZnBpPal40kSSrIOsAeVRchSZLULv4BPNPA+x8FbjODSc1mI0lS17sFmFHAOH8BTi5gHEmSpG5wy8dgxtnU30z6HfAlM5jUbDaSJHW944B/DucNU4DX13ytB6xAupetJEmScjkO+OcJwN/qePPjwD3F1iMpHxtJkgT8B5iZZ8dVgROBP9d8nQ7PvA/u/S6wfkkFSpIkdaD/ADPvAZ4d5huvAk5Jk5nuLbwqSYMJMcaqa6hbCKF9i5fUcr4JHApM6u/FVYGVgPcAH1z8pRnAmcAnAK4HjgWeB+4rp06p28QYQ9U1aHFmMElF+iZw6Ndh0vbAmsDkQXZ+EngE+BPMOA3OfCjLYJIKN1AGs5EkSTVOB3YAVqOmobQKcApw+JL7zwC+T7rzyGKuywaS1DAbSa3HDCapaKeTZbDzYNLBA+z0FPA14CswYz/4/hbwv0uEMElFsZEkScNwKnAkMGkl4CukqUo1ZpBOiF0KfKq/AWwkSYWxkdR6zGCSynLqqXDkXjUn9F4BLCBdw3YWzPj6EBlMUmFsJEnSMH1pBXjH6cBB6XnvLOsZwNmkSUoD+ivwdtJtaSU1xEZS6zGDSSrTl4B39D75OEx+AfhOzgwmqTA2kiRpuM4Cjlj09KPZ96/nff+fgTcWWpHUlWwktR4zmKRmGnYGayWrAiOBx0h3+X0GmFdpRVJuA2WwUc0uRJLaVVuGF0mSpDbXlhlsJVID6W+kxTffRFpY8zjgWmAZ0s1Z5lZVoFS/EVUXIEktaRlgfNVFFGd5UpaRJElSiVYAAnAzMB2WXy3LYH8ENgJ+A+wIXAXsAYypqE6pATaSJKmvicAPgf+pupD6LEfKL8uR1nWaDPwLeE3N894vp6VKkiQVYzkg3A3LrQeTRwySwX4No7YBfg7sDYyuqmKpPq6RJEl9XQgcUMA4Fa2R9FB22H+STooNZlfgd8D8souSGuAaSa3HDCZJS3oIeOMz8M/Jw8xg7yQ1lVw7SS1moAzmjCRJ6ms2bf0P+SbALQwdYAB+DbyVNAnLS98kSZLqt8lEuCXUkcEuAPbHqeJqGzaSJKmvdwEXFTDOSGBcAePkNB5YCphOmlqd11WktR53K6MoSZKkDjceWGo8TJ8Oyw0jhC2WwX4KHIRn9tQWbCRJUl+jKeb/HV8H/KGAcXKaDrxAWidckiRJzTF9LLzwJCzTaAg7Bzi4gIKkktlIkqS+zqGYNZIW0NTFh1YCXmzg/aNIi3RLkiQpv5VmwIsT63+/GUztxkaSJPV1EGl68YIGx7kBeEPj5eQ1k3Sdfb0uA96OYUaSJGk4Zs6HiQ3cgmCxDDYCg5hano0kSerP4cC5VRfRfBeT1hk/tOpCJEmS2sVyNDYtnJoM9mPgPY2XJJXJdeElqT/nAIdUXYQkSZJaWiCtmL1UwWNKLcwZSZLUn8OAHwP1TlNuYHqzJEmS2sTLFNtEAvgB8KGCx5QKZCNJkgbyHuD7db73D8AbiytFkiRJLShSzgnEM4DjShhXKoCNJEnSEs4B3lt1EZIkSa1uPGlWUkHOoSaDfR34eHFjS0WxkSRJg/kg8M1hvue3wA4l1DKECcALzT+sJElSVystg30UOLWMgaXG2EiSpKEcS1r08PM59r0K+AZwbZkFSZIkqeM5I0ktykaSJBXlMmBP4GoqmZEEsAwwq6CxfgD0kPpokiRJGlgpGcy7t6lFjaq6AElqGycDnwFOAFYB7gE2B47MXm+BO7UVWUKo+ZIkSVI/XgTGm8HUXUKMLfCbT51CCO1bvKT2Ffo8XlBVIf2bCSxb0FhHk24a4v/ZqkoxRrN0izGDSVJmJPA8ML6kDLaAdCIzzxILUsEGymBe2iZJwxVrvlqsiTQsfwLmANsDvwR2X3KX04DZpElYkiRJ6qOHdF1bwattL8xgI0iz4m0kqYV4aZskdZt/AK8CRpNOJ/yW9K/BW4H9gX2Ba4CfppNsI4EvAmOBL1RRryRJUiubX/yQi2WwkfAFp4CohdhIkqRucRvwDmAMqSvUa0z2fQRwISm1HEhqNJ2dXhoFfDr7/plm1CpJktROXgE8QZqdVCAzmFqRfU1J6jCrAs/298I44F/AJoO8eSwpqYwBvseihcRJfaXRBdUoSZLUUWbDqnGADNag0cDo44FTShhcqoONJEnqMC8xyOLYS5FmHOUxFjgd+FDjNUmSJHW6QTNYo74NfK6swaXhKa2RFEJYLYRwXQjhzhDCHSGEY7Ltk0MI14QQ7s2+L5dtDyGEb4UQ7gsh3BZC2Lys2iSp060DPF274WZg7ToGGg98BTgmPT0OeAb4ckPVSSqTGUySqrNEBivIcXPgmTlmMLWGMmckzQc+GmPcGNgG+FAIYWPgROD3Mcb1gN9nzwF2BdbLvo4kXVQhSarDTGrOiN1CWlw770ykvpYiTaU+AX4AbIA3DpFanBlMkqqwDsx8BcRnih/aDKZWUlojKcb4WIzxX9nj54G7gCnA3sC52W7nAvtkj/cGfhKTvwGTQgirlFWfJHWFW4BXUn8TqdcEYCl4D/BB4MUGh5NUHjOYJFXkGeBPwHLFD20GUytpyhpJIYQ1gdcANwIrxRgfy156HFgpezwFmFbztkeybX3HOjKEcFMI4abyKpak9rcZMHVFeKzRJlKNc4BvFjecpJKZwSSpyXaEzZ6FqcBjQ+6c3zmYwdQ6Sm8khRAmAj8Hjo0xPlf7WowxMsz1yGKMZ8YYt4wxbllgmZLUcR4Dpm8NPUWlmGPhueNhVkHDSSqXGUySKvAoPPYamP4k9BQ47HOYwdQ6Sm0khRBGkwLM+THGy7LNT/ROl86+P5ltnw6sVvP2qdk2SVK9LoLtVoD1gYcG2OVQ4Hrg/cBv+3n9S6T3rz8JvlHCVG1JxTODSVKFpgHbwHbrw/rrw0MDhLDcGezb8I1vlFSrVIcy79oWgB8Bd8UYa/+zvxI4LHt8GHBFzfZDszuHbAPMqpl+LUkarj8BW8N/R8O9wI7AxjVftwKHALsD3wfOB57vM8QXSP9QvIk0xjNHAic0pXpJdTKDSVIL+C/89164917YcS5svAdsfM8wM9ip8KaN4d7PwTMlLOAt1Sukmc0lDBzCG0i/xtwOLMg2f5J0jf7FwOqkE+T7xxhnZKHnDGAX4CXg8BjjoNfghxDKKV6SOsGTwCuW3HwJ8GlgHunyt5WAp0gBZg3S7ZpuIP1mOR0YQ7oF1FO1434DOLXE2qUaMcZQdQ3txAxWp6VJP7WXgddVXIukzrIBcBpcsj18evwwMtijMP8keOrsasqWBspgpTWSmqEjQ4wkFWUz4HfACotv3hh4kPTbYn/WIl2HP+iJryeB00lzrqWS2UhqPR2XwZYFriMtS94D3Nzn9QXAa5tdVIf5JbDyEPscC/yl/FKkSmwCG/8CHvwsvPRx4FXZ9g8DfwVOhbV27JPBTgO+Qro9glQBG0mS1I0eBcq4ifclwMeB/5YwttSHjaTW0zEZbAXSqf9RwFaD7LcAGEmaKrBtE+rqRA+R5sIN5i5gJnAS8AbSHLmvAr8otTKpeTYF7iedtTuPdBHy5sCGpNuy3dZn/4co9tZv0jDZSGoDywA/A/YAxrFowbX5pLVNAulk2fZVFNdBLgA+BnyedO3xfdWWI5WrrEbSdNJZsq+XMLbUh42k1tP2GWxVUiAYA2yTI4NFuO5PsP0bSZe/van5JbelK0mzvfYD/gUXrJ4zg90JLE+67ude0i/S3wMuLLtgqYk2IzWKpgLLkRqpT1dakbSEgTLYqGYXokVWIZ1kOYb07+Jo0kmua0gnvbbL9luQbQvZtmuy7Ts3s9gOcT6wF+mE2Kak5n/fhe1qfRl4NemE2HfxhJi00BRgvaqLkKT6rDIevrrdMDJYgO22yzLYG81gufwCeCupWfcLOH/FYWSwjWsy2Hrw3fXgF6uQPpzzyy5capJbs+8zqyxCqo8zkiq0DmkVzL8CO9Tx/l8DuxVaUef6GTCJNJtr/DDedyswmXRbm68x8O3TpZZ0AbAvMLak8X9Aul+tVDJnJLWets1g2Wq26ywFt7+pgQz26yyDGcQG9iwwqeAMdg9pKtN5xZXZEbYnLRDvuoWSCualbS1mCuky2J0aHOcKYJ9Gi+lg5wMTSDOKGvld+lbS4sTnAJc3WpTULLuQEvxyJY3/IOnubd8uaXwpYyOp9bRlBlsLOAum7FBgBruC1LBvv59GuS6A898GE0aXkMHuBr4I/KSxEjvGDqR7yF8AjCDdn7HWfcDxzS5KUqewkdRiNgD+U9BYlwAHYIbpz96kk1ZLFzTebcDngJ8XNJ5Uuj1I6Xv5ksZ/EJtJKp2NpNbTdhlsXdI16jsXnMH2hwMuMYMt5jzgANh7dIkZ7D+ka9/OLWjwdrUTqam2IvB30tpfr++zzwzg98AjwEeaWp2kDuAaSS1kKmltpKK8A5gHHEK6ll+L7E66NL8omwKfIq1XdWmB43aazUiz/Z1h3QKuAl4ucfw1SbfEtpEkqZVNBnYuIYNdagZbzNnAgcCokjPYhqQFlLqtkbRf9tVrXRbdbXCN9G2JDDaZ9MvCLFLDCdK95Y8ptVJJHc4ZSU22OvBNyrkcbRTQU8K47eq7wBEUG2J63UI6AXRJCWO3kx2Bg/rZPhXYhEV3vQGYBpzcjKK0pHeQ/kKsUNL45wMHlzS2hDOSWlFbZbB1gW/A6nuWkMHOBo6AUdEMBqSO2g/hu2ObkMHuIM3CuZ4087bTvZ102dpr0tOGMtiLwMWk1c5tKEkahDOSWsCapMUC96m2jK5xMOUEGEgnwXovQe/mZtJGwOGDvF772jOku/j2ehb4RBlFaUmXAF+hvEbS1sB7gR+WNL4kNWIlWHPPkjLY4cB78Nq2XlsDI5qUwTaBS56iO5pIAFuwsIkEDWawOfCJv1HujGVJHc1GUhOtCLytxPFPJ51UcGp1c7yadKvgbm4kDcfywPtqnj8PzCdNU1ebm0U63SlJrWY94NiSM9g34fSj4ZjY5RnsG6Q7eTbht4tXk2WwO0kzkjrdvjS0OvwSGWw0zH8lfOr5BuuS1LVsJHWQo0hr6HV1iFGpdgPeXPN88wbGWhp4N2mtzI8Dn25gLFXsbuA3VRchSX2sD3yedElQmR4ygwHwMM7MKsO+wCdhty0LzGBLw7s/DF9+CT7eA582hEkaphFVF9At1gA+UPIxTsLr82t9lvJn7G4L7FXyMVrFbqT/xo6v+dqhwTGXAb4AfLDBcZTDV0nXAjxZdSGS1CRrAvs3IYN9A05yjaQ0NT77IZjBCrIP8L+piVRKBlsAH5zR4ECSupKNpCaZCryr5GN8mS4/E9bH14E5JR9ja+CtJR+jFexKCjCvK3jcicCxBY+pAXwXeIF0uxtJ6nTrsfAmAKVnsE/Dl4MZDEhnh+Y2MYNtBexZ8oGq8DbS7dfmw649JWawQHmLWUnqaDaSpAZtTmq0dKpdSesYFR1gao0nXRKgko3GRpKk7jCC9P95zeAv4oucAsxr3uE23wp27bRpSfuSgterUxPpUz0lZrAJMP5jZjBJw2cjqQlWAw4scfyvk+5E4pmwJZ1B+VOrx5HO6nSiXUhrF5XZRAJYChfdbopvkO7e9lTVhUhSye4Gzm5CBvs6fO1/YYFrAy2hKRnsFph4dckHabalgTGwyy5p7aLXlRnCXoKlvgufGgN8uMTjSOo4LrbdBOuQFmEsywm4tuFAPkX6j/xoUrOiDDfQmXduewvwGeC1Fdehgn0VGAkcB/wBeAOwagPj3Ql0WoiX1N7WJi3st1ETMtgJEA1hi/s+8OHUnCg9g90Al3RaCPsJvCXAZz4Or92o5GP1kO66OgJYveRjSeoozkgq2VQ689LtdvIJ0vIwL5Y0/sbAdiWNXaXVgElNPN5Y4PAmHq+rfRn4JnAicD/wS1KQHI5/Az8ATgV+Umh1ktSYZYBXwtQ1Ss5gP8Azef05HvgeMMcMVq/VpsGkmU040NLAaTD2QDj8hCYcT1LHcEZSyTam+OuOI3BOn+ca3AmkSRiTSHcBXrrAsbcHDgGuL3DMVnAv6QqoDZp0vJHAVsDZTTpe1/ti9v3/gNmkuxutNoz3Pw6cC/y12LIkqWG3AF+AjY+Fj+xe7NCLZbC/mcEGdCzpt4wj4IRxZrBheTPcewo8tU2TMtgEGHkabHWOGUxSfjaSSrQq8OaSxn53SeN2st6G3kzgCNIJSw3sVTR2xdNwzCFNavlgk46nGl8mnTJennSp29osPr39H6RFwPpOr18t22YjSVILWnUFePNryhl7YQbbicXP7GlxR5F+0zgMPjIubZqJGWxIH4RXbdPEDDYHfnKJGUzS8NhIKtHKwLZFDBRItwElnfm69NIiBu1eHyFd0/kuYNnejW8m/SI9kH8BD/T/0n3ATUUV10K+A2xJ6iuUbS7wqyYcRwP4Eun/Z34LvAk4hjTj6CHSf/sTgfX7vOfXeOpSUsta+RbY9lRgh+LGjMDCCHYJcHBxY3es95OmIx0EjB8ggzXADNa4uXPhV4YwScNkI6lE04FrafDa7QAL9oGrLknX+V8RYf99+tnvykYO0n2OJQWZQ4FltyP9iz3YgobfA34D3EG6nfD6pPRyZ1pn+AdlFtsFliad1F2u4jq6Wu+lbv/Mvt8MXAccSVpQ+6QqipKk+hSSwUh3xL0K2DPCFVfC/r0v7D/we9THe0nNpMnATnDshJoM1uDQZrDGLb00nHMOLGcIkzQMNpJKtAZQ96X5I4GdgdHQcxnsS+oV7RuAK/rs+2tgHulU2W/xgv2cjn4DhFvh4K/BpKHuivGB7OtHpFuPvJN0NvLHwMOkX7Q70O2kNZiHs3TOcM0hhf2yFuJUHb5R8/jMyqqQpLo1lMF6LYCeq2HfXeDKCPvu03hdXav3esAfAfvD0RPTRNiDaeDGHv8lneDrNDcD28LtqzYpg82FF39X4kEkdaQQ2/iepSGEli/+LaTeTm6jSJeWjIOeq9KEgB5gl8Hes3N2kBHAGFJTSYN7HfAz4EtwxkfgoA3qCzIPA1+7Er59InBXkQW2jh+QJqWU5UlgpRLHl9T+Yoyh6hq0uLbIYG+B3w4rhC3SA1w3B3omwC6/IZ2ke0uBxXWzs0m3NQbO2BoOWqbODHY6fO04+HZxlbWOnwHvbFIGewZYocSDSGprA2UwZySVaDLp38lbgc3yvmkizP8d3EA6S5Ars1xT8/gNpAVzWz7eVWgr4ELSgsI/SGtBBuCVpAWmhzOz9zLg23uRmnjHA3cXW2rbWxZYi3QHHUizuTamMxc0kCS1jLoyGDB/PtzwN5jzhiyD9ZBO2Kk4hy96eNQPIBwAr1y2zgxWcGldZy7eNENSXWwklWhD0j+Kx5JmFg1qNMzbAm5aOl3iU3dmuZY0o2ku8HfSxf1a3AUsflcq4EPZ96vJ/7OfDjzY+2QP0uWIxwD3Nlhfi3mAtO7yysN94zKkSwAPJf1QX0VqKn0ZZh+QZqO/ah7c9M9BxpAkqQ7DymCkydw3zYcX/wA770nKU87wLt/7sgz2Krh6I9h5Ur63TZ8ODz5YXlmtou4MNojZZBlsLtz0O9IirJI0TDaSSnRD9rUFaQLGBPpfz3nuaLjrzTDltwXd5e2P2feJuPBMX68Cxg788j0sOhu2NumMZq1HSP+gQ7pd/WJnwnYlNZM+CNzfcKUt41TSWd2jhvOmZUh3aPkucBtpka/zFr386F/T2hXnPAW77wrMJ502liSpADcAN8yCLW7KMtgE2GgN0ro6myzaby5w1zyYcits+zyLziYVEsiUy/vSt3vOgOUOAiYNkcEehZ+cCt/u5OlIDwBPw6kr1JHBhvAoWQZ7GnZveCExSd3KNZKa6JXAT2uejwHWGQ1/fBO88xr4ObB9UQe7HXgt8HJRA3aAjYDfA6vk2/0rLHlp4WnAuUO98TLgbcMrrdWdCByYPV6Jfs6MLcfiq0G+gXQnvH7MIX0Mi2WXWSy6tc4C4N8NFCup47hGUutpuwy2Nfz0HNJl6F+HMRvCOnPhj3fBO5+Bn+9YYAZT/a4Cdh8kgz0OfJ50oqrTfQ74IJy4/BAZbBgWZrC5wPV42aakIQ2UwWwkVWg14Merw84PlTD45qR1adr6J1SwR4ApTTjONaRL3F4krcbdYU4GPlO7YRJpFtYpg79vNukjWABsMNBOPaQZTP8D/KeBIiV1FBtJraftMtirSOsj3gscB6tdBT9+FHb2F+nW8l3SEg0D+RZpBepu0We17SUyWE6zgUfmwIJpsMG6pNvBrT7EmyQJF9tuSdNGwc5rljT4v/DStqrsDNxJupPeoLfb6xD7M2QTCdJ6Uh8izW4a0GzgEOBS0hQ+SZKKcDuLXdI2bRMnY7SkD1ZdQIt5CphJfbe1q/Eg8KH74cTdSQuCPtbYeJKUu5EUQpgMEGOcUV45XWYNFq1nVNb4d+GsJEjTv2ybFmImKZAsA0xeGp5fPk2VHujOsbNJaxqMIzWRBr0T4QTgH6TLENcAypitJ0ltxgwmdamTSLO1D01PZ06CBydlGQx4nmFmsAeB9UusV1LXGDHYiyGE1UMIF4YQngJuBP4eQngy27ZmUyrsVKMo9hYM/bkDWLrkY7SLf5AuLFfDvkm6+donJsIjx8A5X4Qv1rz+MvB09ngO6W45r83eM2gTqdd4UqfqnqIqlqT2YwaTBKTr2dZKX9/8HKw1Cz5BWi7gHHJmsLnwlseRpMIM2kgCLgJ+AawcY1wvxrguaaniy0lXmqseo0j/r/7nqgvpIk+SzuioMGe+D1b7PPwdOIE0S/ox0v9pHJA9/j/gU6Q73EmShsUMJmlxpwGfgjMfg9We65PBXoaLnu0ng80D/gTsWFXRkjrRoItthxDujTGuN9zXmqXtFnrs9UrStfrNsCzwXJOO1eqatdh2r05dI2mp7PsHgK8Nvuto0tVqM+s5zlzS5/X0UDtK6nTduNi2GUzSoI4gNZZ6rz64EjgfOCM9XZjBbsUFwSTVrd7Ftv8ZQvgu6Y7n07JtqwGHATcXV14XGUm6sFnNN4t0OeHIJh1vNGnB8xeadLxmGE869TWadOH9ELYkzciuq582mrRo+Yr1vFmS2p4ZTNLAfgTMB76RPX8DadGkLDc1lMEkaQhDzUgaQ+p3782iuRzTST3vH8UY55Re4SDa8mzYpsD1pIsKy16/6DlgKukfFaWmzp2kGN4M84HfAHs26XjNcHL2PdQ8LstcYCyp8eqsOqmrdemMJDOYJEmq1EAZbNBGUqtr2xAzAtgM+FfJx5kIvFjyMdrJPcA6DL0yWFE69dK2k4HPNOE4c0ln1R4mXaIpqWt1YyOp1bVtBpMkSbnVdWlbCGEU6WzYPix+NuwK0tmweQXW2B0C8GrgLyUf5+WSx29H6wN3k2ZpjaN5DaVOMop0yVkzjCGtj+SloJK6kBlMkiS1qqEubbuAtE7uuaSliiH9Gn4YMDnGeEDZBQ6mLc+GbQacTfoJ3lbicZbCZtJg/gNsUPIxOnFG0peAE5t4vEj673hCE48pqeV044wkM5gkSapavYttbxFjXL/PtkeAv4UQ7imksm5zK3Aozbtrm/o3n3Tp1GjSLDG1ph5sIknqVmYwSZLUkoa6uGdGCOEdIYSF+4UQRoQQDgCeLbe0DhZJvyCXZX52DA3slaSFnG+tuhANqVl32ZOk1mIGkyRJLWmoRtKBwNuBJ0II94QQ7gUeB/bLXlM97gA2KWns3hkcOW7NLmAB6WfWg823VjQSL9GU1K3MYJIkqSUNemlbjPFB4ACAEMLy2bZnyi9LuSwgtQIX1GwbD7j8Zn5b1Dy+Dti+oHEjndmY6v1zNetywB7SotuS1GXMYJIkqVXlvm9VFl6WCSHsF0LYsMSalEckzdaIpHbgyOzLJlJjIsU0gX4F7Np4OS3nk8Bnm3SsTm3GSdIwmcEkSVIrGbSRFEK4vObx3sC1wJ7AlSGEd5VamfrnL9bleTPpb8QI4HL8WVdtHs5GktS1zGCSJKlVDXXXtjVqHn8c2CHG+N8QwgrA74FzyipM/VjAollIvXLPKdOw7AdcSHZRgZpuNukyTUnqXmYwSZLUkoZqQ9S2LEbFGP8LEGN8msVX5lEzRWwgNcOBpLWAfjzM911GOmcsSVL9zGCSJKklDTUjabMQwnOkX6fHhhBWiTE+FkIYgzflbp6+l/jYRGquI4C5wPtz7Hsh8M5yy6ncZ0k/jy+WNP5LpDsPSlJ3M4NJkqSWFGIc/kIwIYRJwEYxxr8WXtHw6mjvVWx673w1BfgHsCnwZM3rL7PoF+r2/pN2hm8CRw/y+nnAYXTPZ/VR4GsljGsjSVIfMcZm3Suy5ZnBJElSswyUwepqJLWKjgkxqwKPkJpImwAPAMtlrzl5vXUEBr/tfbfdZazIRtIMYF3gIWAS/ncvaTE2klpPx2QwSZI0oIEyWN0XSYUQbq+/HC3mUdKla1OAZ4DJpF+k/WW6tUQWfS79fXVbpP4GcGxBY0XgWRb9ty9JGpAZTJIkVWnQNZJCCPsN9BKwcvHldLH5NY97KqtCyi8C3wa+mz1/J3BuHeM8RWqiwuJ/DySpi5nBJElSqxpqse2LgPPpf67FuOLLkdRWamfO/ZS0MPzPhvH+x4C1svdJkmqZwSRJUksadI2kEMI/gcNijP/u57VpMcbVyixuKF6fL7WYkcBY4C3AL2q2/wv4IPC3Pvs/yqLZSJI0gG5cI8kMJkmSqlbXYtshhDcCD8UYH+7ntS1jjDcVV+LwGWKkFjUS2BX4Zfa8B7gJ2KnPfhF4sYl1SWpLXdpIMoNJkqRKedc2Sc01CphQ87wHeKGiWiS1tW5sJLU6M5gkSZ2v4bu2hRD2GOy5JC1mPjCr5ssmkiTVxQwmSZJaSe5GErBln+dbFVmIJEmS+mUGkyRJLWPIS9tCCBsCewNTSSuaTAeujDHeVX55g3NatSRJna9bL20zg0mSpCrVdWlbCOHjwIVAAG4E/p49viCEcGLRRUqSJMkMJkmSWtdQd227B9gkxjivz/YxwB0xxvVKrm9Qng2TJKnzdeOMJDOYJEmqWr2LbS8AVu1n+yrZa5IkSSqeGUySJLWkUUO8fizw+xDCvcC0bNvqwLrAUSXWJUmS1M2OxQwmSZJaUJ7FtkcAWwNTsk3TgX/EGHtKrm1ITquWJKnzdeOlbWAGkyRJ1Roogw3ZSGplhhhJkjpftzaSWpkZTJKkzlfvGkmSJEmSJEkSYCNJkiRJkiRJOdlIkiRJkiRJUi42kiRJkiRJkpSLjSRJkiRJkiTlYiNJkiRJkiRJudhIkiRJkiRJUi42kiRJkiRJkpSLjSRJkiRJkiTlYiNJkiRJkiRJudhIkiRJkiRJUi42kiRJkiRJkpSLjSRJkiRJkiTlYiNJkiRJkiRJudhIkiRJkiRJUi42kiRJkiRJkpSLjSRJkiRJkiTlYiNJkiRJkiRJudhIkiRJkiRJUi42kiRJkiRJkpSLjSRJkiRJkiTlYiNJkiRJkiRJudhIkiRJkiRJUi42kiRJkiRJkpSLjSRJkiRJkiTlYiNJkiRJkiRJudhIkiRJkiRJUi42kiRJkiRJkpRLaY2kEMK4EMLfQwi3hhDuCCF8Ntu+VgjhxhDCfSGEi0IIY7LtY7Pn92Wvr1lWbZIkSZ3KDCZJkspU5oykOcAOMcbNgFcDu4QQtgFOBU6LMa4LPAscke1/BPBstv20bD9JkiQNjxlMkiSVprRGUkxeyJ6Ozr4isANwabb9XGCf7PHe2XOy13cMIYSy6pMkSepEZjBJklSmUtdICiGMDCHcAjwJXAPcD8yMMc7PdnkEmJI9ngJMA8henwUs38+YR4YQbgoh3FRm7ZIkSe3KDCZJkspSaiMpxtgTY3w1MBXYGtiwgDHPjDFuGWPcstGxJEmSOpEZTJIklaUpd22LMc4ErgNeB0wKIYzKXpoKTM8eTwdWA8heXxZ4phn1SZIkdSIzmCRJKlqZd217RQhhUvZ4PLAzcBcpzLw92+0w4Irs8ZXZc7LXr40xxrLqkyRJ6kRmMEmSVKZQVk4IIWxKWrhxJKlhdXGM8XMhhLWBC4HJwM3AwTHGOSGEccB5wGuAGcCBMcYHhjiGIUeSpA4XY3Th52Ewg0mSpCIMlMFKayQ1gyFGkqTOZyOp9ZjBJEnqfANlsKaskSRJkiRJkqT2ZyNJkiRJkiRJudhIkiRJkiRJUi42kiRJkiRJkpSLjSRJkiRJkiTlYiNJkiRJkiRJudhIkio0BYg1X7OqLUeSJKkrmMEkqX42kqQm2wBYkH1N6/PaMsDzTa9IkiSp85nBJKkYNpKkJtiKRcHlLiDUfPU1IdtvdtOqkyRJ6kxLZLB1ICxIXwtfmJH2NYNJUj42kqQmCAzePOpvX0mSJDUmbAXhJgivhNAD3MOSwWwS0APhOTOYJOURYoxV11C3EEL7Fq+usD3wO1IoGW7Xdg4wruiCJKkNxRj93a7FmMHU6ranJoNtDvyDocPYi8BEM5gk9RoogzkjSSrZSOr7izYGeKngWiRJkrrFyG1hxDzg7+QLY0sBz5vBJGkoNpKkkuwEXNPA+wPpbJhBRpIkKb+ddoJr5gB/BEaRzurlEdK+ZjBJGpyNJKkkI0jZpREBGF1ALZIkSV1hVxhxHIzag/qC2DhglhlMkgZjI0kqwa7AlQWNNRJ4oaCxJEmSOtmuI+DKtwK/rHOA2cBK6aEZTJL6ZyNJKtiewM+BsQWNF4DxBY0lSZLUkd4Ge56fZbCR1B/ExgHT00MzmCT1r9ErbyTV2A84j2JDRw8wucDxJEmSOs1+V8B5I2H8RQ0OFIAJ6aEZTJL654wkqUCjSTf8KNII4P6Cx5QkSeoko/eEpX7EwiZQEcxgktQ/G0lSiwt4NkySJGlQY4CJxQ5pBpOk/tlIkgryduD7VRchSZLUZd7+dvi+IUySmsZGklSQ8cCkqouQJEnqMuPHw6RJVVchSd3DxbalBh0MfJ5CL8mXJEnSEA4GPn8ATPh61ZVIUndxRpLUoKWBNYFXlDR+D7BhSWNLkiS1q6WBNX8Nr/hqwQOPAe4wg0nSQJyRJNXp3cBHgeVKPk4E7i35GJIkSe1isQy2M3BMwQeYB7zdDCZJA7GRJNVpBWDjko8xH9ii5GNIkiS1k8Uy2CRgSsEHiDD/LjOYJA3ES9ukOhwBfLgJx4nAbU04jiRJUjtYLIO9HTip4APMAbY1g0nSYGwkSXVYBZhadRGSJEldZmEGOwA4lbRQZZEWAP8qeExJ6jA2kqQ6/AT4fsnHmAfsUPIxJEmS2snCDLYysHbBg78M7GwGk6ShuEaSVIc3AjuWMfBSwBVpVvWuwJ8hnRl7axkHkyRJai+lZTCABTDnLzUZTJLULxtJUh3WAtYrYqDlgW8BHwAuBEbC7J1gb+C63n0i8H/Z992LOKgkSVJ7WgtY7+3ASsB3gA8VNPBLMHvfPhlMktQvG0lSlV4EriQ1kXZNm3qAq2v3CdlrsbmlSZIkLeFyYJ+Ka7iFdEZvt4LGewE4AHqu6ZPBJEn9co0kqSorkC70P4iFTaQhXVxeOZIkSYwh5Y2LgQtIJ7Quqnl9r5rXQ9OrS+4DXgI2K2i8eaTZ35KkXJyRpJY3AjgHOLTiOmpdDmwMvLOeN68CfBWYwPDO6AXSbW4lSZLKMJYUut6Rno5YAOcsgEPfQWq2QMoj2evMY9GM6cNI06pLdjlZBvsNcDZweAODzSRdGje34bIkqavYSFJLGwmcBRwMzB9kvyNo7pVf/ybNqh52I2kK8AngNtIta+vxY+Dddb5XkiRpIKOBA9PDkcBZI+Dg/8ky2EH97P8/WQZ7D00LYgsz2N3AP2iskTQW2AY4uvG6JKmb2EhSSxkHnF7zfATwruzxYDmh9oTYUQzedCrKb4EtgP2H86blSEHs1/2//DI5sszfhnNASZKkoY0bB6d/c9HzYWWw1wJbZRnsaEqf4VNXBuvPeNJZwaNzZjBJEmAjSS1mAXAv8LVhvu/ImsfH0JxG0q3ADeQMMWsAB5Cuv59Ev1OZXgKOI004GtQrc5coSZI0uInAF2HBGLj33XVmsPemx8cA8z9C6Y2kYWWwHHJnMEkSYCNJLWQ88DnSGo+N+DKLGkmfAmY3OF4h5pDWFPjo4LucOdQ4gXS6z1NmkiSpCONh/IcLzGBfBD4On5rThAx2A2nR73o6SrOAzwMv58xgkqSFQozte0/xEEL7Fq/FTCD9W35cweNOIuWEsmwF7ARsCew32I5Tgc8y6NpGzwKT8xw04v0WJXWVGGNV94bSAMxgHWBZ4GMwYSn4/LEFZ7ClYdILTcpgm8F+/8uiBcDzmk7KZwwjg0lSlxkog/nrqCo3ETiZ4ptIACeSZjqV5R/Al8hxx9hVGbSJ9BL1r70tSZI0bEvDxE/Cyce2eQa7FfgCcGnON84APgN8Iz01g0nS8Hlpmyq3NHBCSWOfSLrM7UukoFCWW4FfAnv29+JU0i1x+5gFfCd7bIiRJEnNVloGOx5O/BLMn9OkDHYb7Pk54B5gXRZd6jYN+D1p1fAngV8AO8Ksz5rBJKkRNpJUqaWB95V8jJOAb1NuiLmJdOXa3cB6wN61L64BfHDx/WeSFrM8pcSaJEmSBlJqBjsZmAcnzYZvfw9eKjGELcxgt8N6t8PeGwEPZy8+CFwJPA08DvwEZj5oBpOkRrlGkio1lXSyqGwrkU5ENcOmwHt6n0whLbDdZ87406RFLeviGkmSuoxrJLUeM1j7mzoVpjUhhK20EjzZpBC2WAYbQEMZTJK6zEAZzBlJqszSwEFNOtYhwPcod1ZSr9uouanaBGDtEg5yJN5eRJIk1aUpGexsYC4c8nJFGUySVBrnNagyY4HNm3Ssr5F6LxOadLyFngPuKHjMAGxT8JiSJKlrjJ0Im+8/9H4NuQk4Br72fEUZTJJUGi9tU6WadWkbwGWkIPNMk44HwPbAdQWP6aVtkrqMl7a1HjNYm1sLpj7QhAz2PuBcuGxOBRlMktSwgTKYv46qa/yG5kyrXmgl4E0FjxnJf3tbSZKkKv0AWLbgDLYR8Krs8TrAFkUNLEnKyzWSVJmlgLc08XhvI9319eVmHGwl4P3AZ0oYu+yp6JIkqXNNgKV2alIG+zUwp+AM9jbgNcC5wK6kZtIZNa8/AdxYxIEkSQPx0jZVZiXgAuDNTTreFODRJh2L3YGrCh4zkgLZ7gWPK0ktzkvbWo8ZrI2tByvdAhfcCG8uO4QdAFwOU+YWkME2ANYC3gkcOsh+NwJHkdZokiQ1xEvb1HKeYPAcUKTrgTkljr8qsCPwyhKPAdhEkiRJjXkJnvg1HNqMu79+EK5fqsEMtgEpZH2OdEKtT3hcIoO9lnSbuB3xsjdJKomXtqkrvIPiFnhcA1i9z7Y9gI8BvweOXh7uXAW4nUXX8DcqkrphkiRJjZgOfJh0rVmZ/gLsBe94rsEMdgzwgfRwyAwG3AmwJfA7UhY7FLilkQIkSX15aZsqtSLwK9K/94V4LfB3UuOlxkrAkw0OvSawMumuI4cPst8f3gTvPwTuvh04vcGDQvqz3AC8oYCxJKkNeWlb6zGDtb8Vp8KvphWYwfqaBMxqMIOtDXwV1twvZwYjLVF5d+3GO0iXw91ebxGS1L0GymA2klS5VRjmdfO105T/mT3/Z/Z8DrA9MBe4GViQTkLtADzbYJ2nk06K5fGn3eGIq+DeBo8JwHxgdBEDSVJ7spHUesxgHWAqrDKtxPUjJ8EtsxrMYGcBRwwzg6W3LJ7B/gW8nrQw9x31FiNJ3WegDOalbWovm7H44omvJs1A2gJYAATS7B2AbYEX0yXyz95JasjUYSowGVhhGO9543NwzgPwgbXhReD+eg68gHT2rKeeN0uSJA1iHqmpskl5h9iRBppIq8HUyXVkMOAc0tVwCzPYeGA34GvAPqQTjv+ptzBJkjOSVLkVgD+T1lIcVCA1V4bhP8A2d8OsbYEZdRTHwpNhw/cW4LfpJNjbgMdJC0I+kPf9c0mNsrvqObgkdQ5nJLUeM1hnWGEV+POjOTLYcP0H/rM1bPM8zKp3jPPgrIPrzGCZATPYo6QZ7IVMHZekzuVd29SyniadPRrSuvnHvJ+UDbYGZh1MOp21LpXMwdscuJzUV/o1Sy4SOaAxwK3l1CRJkvT0fHjjvQzjLNcg7gMWwP33w71bwtaNNJFWBpZuvKQBM9iqwHWNjy9J3cpL29QSeoCHGaDJsgYwErgn/3ibAi/1PvlHzQurA9Pyj7MCheQYNgOuyB5fS5rq3QM8UsDYkiRJdXkKetaHhyfD6o3e3nZj4HbY9DXw0ssNjPMK4Iewwh4lZ7CRwGoMKxdKkhIvbVPLmAT8G5hSu3E06RKvHGqbMusD/WaY15IW5s657tB5wMH5dl3SdsClpEDUj0eArUjTrZcQSbfnnUe6Y4kkdTEvbWs9ZrDOMmkS/PvfMGXKkLsu7pGab2sD8wbJYHldBuzbYAYbwsIM9jjpri+SpH652LZa3kxgI9LEo5XJ/ifH3coey76vSY7+0I2kO3bknMI9kzSzaal8uy/uetKCjlcCyy/58lTSuuFT+3vvbNJZMkmSpJLNnAkbbQT33A0rB7IgNoTH4LHVgZgzgw23JhrIYEMYNINJkobkGklqKc+TljJ6ajw89RjwcJqc8xSLvtd+RdIMplUp5+ZmHwYuaGSAG0h3CRngliUjSMs3LbSA9Ad7upGDSpIkDc/zz8O6q8JTa8FTWciKT6XHEXgqpsdPkW1bFabEEjLYc8CcAjLYEJbIYJKk3JyRpJbzIrAiaa3px0jhZGXgSdJdRWp7MnX1W5YlpYdh3gGubn8HdgKuIa31tGy2vQdWeR7umASrLCCtSDmS9IeXJElqsheBFWfDmKnw2MPQszKsPAKefCrLYBuQFtVescRzXh8ClgH2LesAySrAHXhlmyTVwxlJallzSVeErUjq+bwC+G+ffVYgnSUbln8BryaFlBx/A2aTe5mmwY+5PPA6UsNoFnAXaWHKWcBDwFqks3CSJEkVmjsXll85y2AL4BXLw3+XJwWx7HL9ujJYHmcCe6SHhWSwwYygmBW9JanLuNi2WtN4am67VqJXA7cOvstY4AzgPSWV8ARpfcpm/HElqR252HbrMYOpVBcCBzQhgy2AtafBS2uWdABJanMDZTBnJKk7zSZ1bnLE4O9QXoCBdLbvPyWOL0mS1I5Kz2BPwn/WLPEAktShXCNJrWlsyeNvRbowPkcjaT5pnaaRJZXyJGmBcUmSJC1iBpOk1uSMJLWeiQx4l7PC3A5smm/X95Mu1y/LisCDJY4vSZLUVnqABWYwSWpVNpLUmuY34Ri3AJsVPGZg8b9VOf6GPUG6da4kSZKAg4CflX8YM5gk1cdGklrPC8CkJhxnQQlj7g9cSvqb9Vbgd9nj3i+Xi5UkSapepJwsKEldwDWS1HqWBmY24TibktZJKkrIvvYlTcnuVft4GJfUSZIkqSSPAqtVXYQktSdnJKm1LAtMJzWTynY7xV3adgTprNYFQ+z3KuCuRU+9d7IkSVLzxVVJzSRJ0rDZSFJrmUW6WP2Fko8TSU2dWxsc55hsrLPyHTICbJgexAdTz2yVBkuQJEnqKBeS1kkqyMIM1vs8wvTpsIohTJLqYiNJrecFYELJx9iQYi9ry+FmFr+q7eE1YLVpza1BkiSp5R0InF/ccEtksIdhNS9rk6S6uUaSNISzgHf398LHgC8P/f4/A2+seb7YettTgMdwWpIkSVIfA2awnPrNYPcB6zUwqCTJGUlSHgPebG2Iu7Bdx+IBZsD3Lw88M8yiJEmSOly9N7y9Dnjjn1l0M5TeL5tIktQwG0lqTS8DS1VdRPIe4Pt9N/4vuWYj5bIS8CSwHOkOb88VNK4kSVIb6zeD5XUdsF1xtUiSFrGRpNa1oKRx1wPuGd5blri7Wu9ZrUFcDeyUZ/BA+pvY+30i8Pzw6pMkSepE9dzh9uqrYaed6nyzJGlINpLUuuYAo4FxBYw1hXSL17VJ18YP01HAt3uffCb7GkKkzl5YAEbW80ZJkqQOcRgwGo76QU0Gy+NXEHeDBWWdkJQkudi2Wtx8ijmbNB9YE5hX39sXAMcB874EHzmB8hs940iXuC1T8nEkSZJa0YL0tWBBlsGAj+R5XyQtFSBJKo0zktT6eoAJdb53ZVJT5knqbiLVlvHxkfClIZpIl2eH3LORgwVgTCMDSJIkdYCjoec78HHgS0PsejkFZDBJ0pBsJKk9vERqJg2nobQS8ATpErmCzAdOBj43wOuXAgdmh2ywbyVJkqT5wHEw//QcGewKmPN2M5gklc1GktrHS9nXMqRpy8sAS7NoIaLlgLnZ4947oRXtUzDvVDglO/TSLJpmfQFwCIX2rSRJkjQPOHGQDHYBHLI0zOk9mydJKlWIsX1vZxBCaN/i1ZhlgVk1jyEtor0xqZk0q783FeCzwMdYbAHwMcB4UsZ5qchjzaGYhcYlqc3FGIe4T6aazQymSowDxmaPPwBj1oHxx8O8efBSoSFMkgQDZzAbSeock4EZJR9jKVIz6fiSjwNp1tUTwCpNOJYktTAbSa3HDKbKjSPd/OTFqguRpM41UAbzrm3qHGU3kSBNOXq5CceBNMVp0yYdS5IkqZ3MrroASeperpEktaoIPFV1EZIkSZIkLWIjSZIkSZIkSbnYSJKG6zTgG1UXIUmSJElS89lIkobr2exLkiRJkqQuYyNJakWzgS2qLkKSJEmSpMXZSJJa0RjgR1UXIUmSJEnS4mwkScN1JPC+ko8xDziq5GNIkiRJkjRMNpKk4foN8OsSx38ZeDNwU4nHkCRJkiSpDjaSpOF6OPsqw4vAW4G/ljS+JEmSJEkNGFV1AVLXegI4tM+2+cCfKqhFkiRJkqQcSp+RFEIYGUK4OYRwVfZ8rRDCjSGE+0IIF4UQxmTbx2bP78teX7Ps2qTKPAb8D3B1n69rqyxKktQpzF+SJKkszbi07RjgrprnpwKnxRjXBZ4Fjsi2HwE8m20/LdtP6jzTgXdh00iSVCbzlyRJKkWpjaQQwlRgd+Cs7HkAdgAuzXY5F9gne7x39pzs9R2z/aXWcxlwIIv+i83rYeC9pNlHkiSVwPwlSZLKVPaMpNOBjwELsufLAzNjjPOz548AU7LHU4BpANnrs7L9FxNCODKEcFMIwXtaqTp3AhcBNw/zfbMo945vkiSVkL/ADCZJkpLSGkkhhD2AJ2OM/yxy3BjjmTHGLWOMWxY5rlSX3wE/rboISZKSsvIXmMEkSVJS5l3bXg/sFULYDRgHLAN8E5gUQhiVnfWaSloxhuz7asAjIYRRwLLAMyXWJzXuDuDL2eODh9h3GnBKueVIkrqe+UuSJJWqtBlJMcZPxBinxhjXJK0mc22M8SDgOuDt2W6HAVdkj6/MnpO9fm2MMZZVn1SY3mbS+QO8/gBwLHAi6XI4SZJKYv6SJEllC83ICiGE7YHjY4x7hBDWBi4EJpNWmDk4xjgnhDAOOA94DTADODDG+MAQ4xp01Do2IS1tCrAhcDhwP/Bp4GdVFSVJ7S/G6OLPdSgrf2Vjm8EkSepwA2WwpjSSymKIUctal3Qe+EFcQ0mSGmQjqfWYwSRJ6nw2kiRJUluykdR6zGCSJHW+gTJYaWskSZIkSZIkqbPYSJIkSZIkSVIuNpIkSZIkSZKUi40kSZIkSZIk5WIjSZIkSZIkSbnYSJIkSZIkSVIuNpIkSZIkSZKUi40kSZIkSZIk5WIjSZIkSZIkSbnYSJIkSZIkSVIuNpIkSZIkSZKUi40kSZIkSZIk5WIjSZIkSZIkSbnYSJIkSZIkSVIuNpIkSZIkSZKUi40kSZIkSZIk5WIjSZIkSZIkSbnYSJIkSZIkSVIuNpIkSZIkSZKUi40kSZIkSZIk5WIjSZIkSZIkSbnYSJIkSZIkSVIuNpIkSZIkSZKUi40kSZIkSZIk5WIjSZIkSZIkSbnYSJIkSZIkSVIuNpIkSZIkSZKUi40kSZIkSZIk5WIjSZIkSZIkSbnYSJIkSZIkSVIuNpIkSZIkSZKUi40kSZIkSZIk5WIjSZIkSZIkSbnYSJIkSZIkSVIuNpIkSZIkSZKUi40kqU7HA18H3g/Mq/m6BNgN+L/qSpMkSepYZjBJqpaNJKkOnwBWAKYBZwCjar72A64A3gpcXVWBkiRJHWhhBvswnDEHRn0fRh0Eo+bAfnPgil/AW3eCq68GtgX+XG29ktSJQoyx6hrqFkJo3+LVtk4GTqp5PmqQfRcAfwR2KLUiSepsMcZQdQ1anBlMVTgZOOlE4DPASBg1CughBa7R2U496WvBX+CPO8IOo0jTlSRJwzZQBnNGkpTTV4EXgf9l8RlIgxkBbEdqJgFsDNxcVoGSJEkd6KtfhRdfhP99EUZ9DkaNzZpIACNZ1ETqfT4my2AR/pg1kcxgklScoX4Plrred4BDgTHZ13CNBF4PPA/8hzTLWpIkSYNbmMFGw5ilhvnm7WDk8/D6m+D5N5vBJKlIzkiShjAWmEh9TaReI7MxxgIvF1GUJElShxsLTPwGjHl/HW/OwtfIN8DE681gklQk10iSBvE94HBS+CjCfOA54L/AlgWNKUmdzjWSWo8ZTGX73vfg8P1h7FLAuAYHmw/zb4Tn3mAGk6ThcI0kaZh+CLyL4ppIkK4lnQwsW+CYkiRJneSHwLsmwtjJNN5EAhgFo14Lk28wg0lSEWwkSQNYmmKyS1/3Aa8rYVxJkqROsDQw7ijgggIHHQX3LWcGk6Qi2EiSmmwN4BdVFyFJktSqvgv8G9irwDH/A2tsbwaTpCLYSJKabDSwctVFSJIktarJwFRgQoFjrg2jzzeDSVIRbCRJkiRJ6mxjgClVFyFJncFGktSPbwG7lDj+6sA1JY4vSZLUjr71LdilxBBmBpOkxtlIkvoxlXLv6jEGWK/E8SVJktrR1KmwbIkhzAwmSY2zkSRVYBqwZ9VFSJIkdRkzmCQ1zkaSVIFXAJ+tughJkqQuYwaTpMbZSJL6cRJwbYnjjwM2L3F8SZKkdnTSSXBtiSHMDCZJjbORJPXjTuCpEsd/FDioxPElSZLa0Z37w1OblDe+GUySGmcjSerHScB2JY7/EvCXEseXJElqRyddCdvdXd74ZjBJapyNJKkfmwOrVF2EJElSl9n8LbDKOlVXIUkajI0kSZIkSa3hj8BDJYz7MHBsCeNKUheykSSV5VXAT4FPA2tnj38Kj34FProicDqwHPCdyiqUJElqLX8FphU85oPAB+DR38JHCx5akrrRqKoLkFrRGcAUYOt63rw5cFQ2wFuAJ4AtgT3hEeDoWXDlWsBWwCeB64uoWJIkqf01lMEGMgMe+T84GriyyHElqUuFGGPVNdQthNC+xavlXQgckHfn1wOHZI/XIjWQ+niYdBbs0t4NDwFr1l+fJHWLGGOougYtzgymMl24ExzwBeC1xYz38L/go1vUZDBJUi4DZTBnJElF2BB43+C7PE2fALMc8C3gOdJt4iRJkgS/Aw6jmEbSg/D0V2wiSVKRbCRJA/gpsC6wxVA7bgu8Y/BdHgZO67txGeDDwItAD3Dy8GuUJEnqNLkz2FD+Cw+fBKdd1HhNkqRFXGxbGsBVwF2D7fBG4BTSYtpvHXi3B7NdfjrQDhOA44FPDbtESZKkjjNkBsvpwenw6Z8NksEkSXVxRpJUr9eRFssewH2k4PIwcO5QYy0FfAKYB3y5kOokSZLa12XAxqSbmAzTfffBT/8MD6+bI4NJkobNGUnSIK4Ebh3oxb8C1wz83vuBzwJn5z3YeFIzSZIkqctd+Qu49eb63nt/hM/eD2efV2xNkqTEGUnSIC4BArAVsAM1J8W2Id2X9qUl33MvcAVwz3APNhv4dl1lSpIkdZRLgHA1bPUs7LADbJ5jZtK9wBX3wz1XkO5ycma5NUpSt7KRJA3hYuBJYMvajTsBn+9//9uBE+o50ALSdXCSJEni4ovhyYthy9MZ+hK3e+D238IJ/8YGkiSVzEaSlMOqwN3AcsDcbNtWffa5B7gOqHMWdlon6asYfiRJkjKrAnffAMuNg7lbAlv0yWB3wT3Xw3U3wc1nVVOjJHWbEGOsuoa6hRDat3i1pf2Bl18DnAT77rf4azcADeeX54BlGx1EkjpLjDFUXYMWZwZTs+0PvLwvsCfs+yrS2pJ/B66HG84pIINJkpYwUAazkSQN1+eAT5U0to0kSVqCjaTWYwZTpfYAJpFujytJKs1AGcxL2yRJkiS1j6uqLkCSutuIqguQJEmSJElSe7CRJA3XvaSVtSVJkiRJ6jI2kqTh+kv2JUmSJElSl7GRJA3X9sCbqy5CkiRJkqTms5EkDdc1wNVVFyFJkiRJUvPZSJKGa/fsq2jzgZtLGFeSJEmSpILYSJKG6/vZV5HmA38lXTYnSZIkSVKLspEktYKZwHZVFyFJkiRJ0uBsJEn1mAE8XXURkiRJkiQ1l40kqR7fBb5a0Fg9wMMFjSVJkiRJUolsJElV6gFuA7aouhBJkiRJkoZmI0mq0uPA5lUXIUmSJElSPjaSpHq9BLxQ53t7gGcKrEWSJEmSpCawkSTV6wzgpDrf+zCwPelubZIkSZIktQkbSVIj5gKzgJeH8Z4eYDJwAfDKMoqSJEmSJKkcIcZYdQ11CyG0b/HqLEcAZw2xTw+p4XQPLq4tScMQYwxV16DFmcEkSep8A2WwUc0uROpI80lrJo0CxvR5rQeYA9wKbNvkuiRJkiRJKpAzkqQivQM4j3TR6EhgHnADsEOVRUlSe3NGUusxg0mS1PkGymA2kqQy7AZ8CNi96kIkqf3ZSGo9ZjBJkjqfjSRJktSWbCS1HjOYJEmdb6AM5l3bJEmSJEmSlIuLbaulDXQK2tOgkiRJ5TGDSZIGYiNJLWssMLuf7RGn0kmSJJXFDCZJGoz/FqjlTCQFlf4CDKQzZJ4NkyRJKpYZTJKUR7vPSHoaeDH7rtaxAg18Ji8w8HRq1a2hz0Sl8XNpPX4mrWeNqgtQv14A7q66CC3BDNZ6/HelNfm5tB4/k9YzYAZr67u2AYQQbooxbll1HVrEz6T1+Jm0Jj+X1uNnIuXj35XW5OfSevxMWpOfS+vxM2kvXtomSZIkSZKkXGwkSZIkSZIkKZdOaCSdWXUBWoKfSevxM2lNfi6tx89Eyse/K63Jz6X1+Jm0Jj+X1uNn0kbafo0kSZIkSZIkNUcnzEiSJEmSJElSE7RtIymEsEsI4e4Qwn0hhBOrrqebhBB+HEJ4MoTw75ptk0MI14QQ7s2+L5dtDyGEb2Wf020hhM2rq7xzhRBWCyFcF0K4M4RwRwjhmGy7n0tFQgjjQgh/DyHcmn0mn822rxVCuDH72V8UQhiTbR+bPb8ve33NSv8AHSyEMDKEcHMI4arsuZ+JNAxmsGqYv1qTGaz1mMFalxmsc7RlIymEMBL4DrArsDHwzhDCxtVW1VXOAXbps+1E4PcxxvWA32fPIX1G62VfRwLfa1KN3WY+8NEY48bANsCHsr8Tfi7VmQPsEGPcDHg1sEsIYRvgVOC0GOO6wLPAEdn+RwDPZttPy/ZTOY4B7qp57mci5WQGq9Q5mL9akRms9ZjBWpcZrEO0ZSMJ2Bq4L8b4QIxxLnAhsHfFNXWNGOP1wIw+m/cGzs0enwvsU7P9JzH5GzAphLBKUwrtIjHGx2KM/8oeP0/6P+gp+LlUJvvZvpA9HZ19RWAH4NJse9/PpPezuhTYMYQQmlNt9wghTAV2B87Kngf8TKThMINVxPzVmsxgrccM1prMYJ2lXRtJU4BpNc8fybapOivFGB/LHj8OrJQ99rNqsmzq52uAG/FzqVQ2ffcW4EngGuB+YGaMcX62S+3PfeFnkr0+C1i+qQV3h9OBjwELsufL42ciDYf/frQW/51vIWaw1mEGa0mnYwbrGO3aSFILi+lWgN4OsAIhhInAz4FjY4zP1b7m59J8McaeGOOrgamks/gbVltRdwsh7AE8GWP8Z9W1SFLR/He+Wmaw1mIGay1msM7Tro2k6cBqNc+nZttUnSd6p+Vm35/MtvtZNUkIYTQpwJwfY7ws2+zn0gJijDOB64DXkaawj8peqv25L/xMsteXBZ5pbqUd7/XAXiGEB0mX4+wAfBM/E2k4/PejtfjvfAswg7UuM1jLMIN1mHZtJP0DWC9b5X0McCBwZcU1dbsrgcOyx4cBV9RsPzS7Q8U2wKyaab4qSHbN8I+Au2KM36h5yc+lIiGEV4QQJmWPxwM7k9ZNuA54e7Zb38+k97N6O3BtdgZTBYkxfiLGODXGuCbp341rY4wH4WciDYcZrLX473zFzGCtxwzWesxgnSe06+cRQtiNdJ3lSODHMcZTqq2oe4QQLgC2B1YAngBOBi4HLgZWBx4C9o8xzsj+cT2DdJeRl4DDY4w3VVB2RwshvAH4E3A7i647/iTpGn0/lwqEEDYlLRI4ktS0vzjG+LkQwtqkMzGTgZuBg2OMc0II44DzSGsrzAAOjDE+UE31nS+EsD1wfIxxDz8TaXjMYNUwf7UmM1jrMYO1NjNYZ2jbRpIkSZIkSZKaq10vbZMkSZIkSVKT2UiSJEmSJElSLjaSJEmSJEmSlIuNJEmSJEmSJOViI0mSJEmSJEm52EiS1FQhhF1CCHeHEO4LIZzYz+vbhRD+FUKYH0J4e5/XekIIt2RfV9ZsXyuEcGM25kUhhDHZ9rHZ8/uy19cs/Q8oSZLUgsxgkopiI0lS04QQRgLfAXYFNgbeGULYuM9uDwPvAn7WzxAvxxhfnX3tVbP9VOC0GOO6wLPAEdn2I4Bns+2nZftJkiR1FTOYpCLZSJLUTFsD98UYH4gxzgUuBPau3SHG+GCM8TZgQZ4BQwgB2AG4NNt0LrBP9njv7DnZ6ztm+0uSJHUTM5ikwthIktRMU4BpNc8fAaaEED4XQthrgPfUGhdCuCmE8LcQwj7ZtuWBmTHG+bVj9j1e9vqsbH9JkqRuYgaTVJhRVRcgSTHGT+fcdY0Y4/QQwtrAtSGE20nBRJIkScNkBpNUD2ckSWqm6cBqNc+nZttyiTFOz74/APwBeA3wDDAphNDbGK8dc+HxsteXzfaXJEnqJmYwSYWxkSSpmf4BrJfd4WMMcCBw5RDvASCEsFwIYWz2eAXg9cCdMcYIXAf03l3kMOCK7PGV2XOy16/N9pckSeomZjBJhQn+fZbUTCGE3YDTgZHAj2OMp4QQPgfcFGO8MoSwFfALYDlgNvB4jHGTEMK2wA9IC0COAE6PMf4oG3Nt0qKRk4GbgYNjjHNCCOOA80hnzWYAB2Zn0iRJkrqKGUxSUWwkSZIkSZIkKRcvbZMkSZIkSVIuNpIkSZIkSZKUi40kSZIkSZIk5WIjSZIkSZIkSbnYSJIkSZIkSVIuNpIkSZIkSZKUi40kSZIkSZIk5WIjSZIkSZIkSbn8P7qg7iRr4nMnAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 1440x4320 with 2 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "c='8'; fov = '001'\n",
    "ref = imread(f'tif/Cycle_0/Cycle_F{fov}.tif')\n",
    "mov = imread(f'tif/Cycle_{c}/Cycle_{c}_F{fov}.tif')\n",
    "\n",
    "ref_max = ref.max(0)   # max proj of z for each channel in reference\n",
    "ref_binary = ref_max[0] > threshold_otsu(ref_max[0])  # binarizing for channels - T/F\n",
    "\n",
    "mov_max = mov.max(0)  # max projection by each cycle in for loop (ref above)\n",
    "mov_binary = mov_max[0] > threshold_otsu(mov_max[0]) # binary of moved image\n",
    "\n",
    "out = imread(f'reg_bin_Cyc_{c}/Cycle_{c}_F{fov}_bin_reg.tif')\n",
    "\n",
    "fig, axs = plt.subplots(1,2,figsize = (20,60))\n",
    "axs = axs.ravel()\n",
    "\n",
    "im_reg = np.zeros((2048, 2048,3)) # empty color image\n",
    "im_reg[...,0] = ref_binary\n",
    "im_reg[...,1] = out.astype('bool')\n",
    "    \n",
    "im_orig = np.zeros((2048, 2048,3))\n",
    "im_orig[...,0] = ref_binary\n",
    "im_orig[...,1] = mov_binary\n",
    "    \n",
    "axs[0].imshow(im_orig[:500,:500])   #before reg\n",
    "    \n",
    "axs[1].imshow(im_reg[:500,:500])   #after reg\n",
    "\n",
    "axs[0].title.set_text('Before Registration')\n",
    "axs[1].title.set_text('After Registration')\n",
    "for ax in axs.flat:\n",
    "    ax.set(xlabel='0:1500', ylabel='0:1500')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "id": "cac9fd61",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAABJIAAAI4CAYAAAAxusOtAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/Il7ecAAAACXBIWXMAAAsTAAALEwEAmpwYAABZrUlEQVR4nO3dZ7gsVZm38fs5kcwhSVZUkjgCIiiYQEwgCogBDIjKiAkVRRRnDIgRXyOiIqKCyiiICowRBEyDoijBgAoiEiTDIZ+83g+r9j69c3V3dVeH+3ddfe3u6qrVz9kNnD9PrVoVKSUkSZIkSZKkmcyquwBJkiRJkiT1BxtJkiRJkiRJKsVGkiRJkiRJkkqxkSRJkiRJkqRSbCRJkiRJkiSpFBtJkiRJkiRJKsVGkqRREfH6iLglIu6LiPXqrqdVRf2PqOmznxIRf6vjsyVJUv+IiCdFxFVFbtm/7nqaUXfeiYj/ioiT6/p8adhFSqnuGiRVJCKuBTYElgNLgYuA16WUri9x7FzgHmDXlNLlnayz+LyfAbsCy4BFwC+AN6aUbur0ZxefvwfwjZTSZm2Ok4CtUkpXV1GXJEkaLEXm2QHYKKW0uGH7+cA5KaXPFK8rzxTFmA8ACbgbOB04KqW0vKrPmOHzjwG2TCm9vI0x9qCCzCapOs5IkgbP81JKawAbA7cAny153IbAKsCfm/3AyFr578nhRa1bAmsAH29hjI6JiDl11yBJkvpXRGwBPIXcyNl33NsPo4XcNcXnTJdZdijy1u7AgcCrq/jMKrSRISXVyH9ppQGVUloEnAlsN7ItIuZHxMcj4rriErYTI2LViNgaGJmevDAiLij2f2JE/C4i7i5+PrFhrJ9FxIci4v/IZ7oeERHbRsR5EXFnRPwtIl5cstaFwFnAjg3jTzlWRKwXEf8bEfcUdX0wIn7V8H6KiC2L58+JiL9ExL0RcWNEvD0iVgd+BGxSTCe/LyI2iYhjIuLMiPhGRNwDvDIiHh8Rv46IhRFxU0ScEBHzirF/UXzk5cUYB0bEHhFxQ0Mtjyp+Vwsj4s8RsW/De6dExOci4gdFfRdHxCPL/M4kSVJfeAXwG+AU4JCRjRHxD+ARwP8WGeLXxVujmaLY77kRcVmRIy6KiO0bxrg2It4ZEVcA9890AqyY6fR/jM1b042/U0RcWmSUb0fE6RHxweK98XnnnUXOurfIbU+PiL2A/wIOLP5Mlxf7TpYhXxURVxbHXxMRry32nS6zfaPh8/ctctbCYvxHjfs9vT0irigy7ekRscqM35ykKdlIkgZURKxGPuv0m4bNHwW2JgeILYFNgfemlP4OPLrYZ0FKac+IWBf4AXA8sB7wSeAHMXbtpIOBw4A1gduA84D/AR4CHAR8PiK2YwbFmAcAVxevV59hrM8B9wMbkUPZIePHbPBl4LUppTWB/wAuSCndD+wN/DultEbx+Hex/37kBtwC4DTyZYJvBdYHdgOeDrwBIKX01OKYHYoxTh/355oL/C9wbvHneBNwWkRs07DbQcD7gXWKP/+HZvp9SZKkvvEKcp44DXh2RGwIkFJ6JHAdxUzylNJuxf6jmSIiHgt8BXgtOYt9ETgnIuY3jP8SYB9yfls2XSERsS15dtRI3ppy/OKk2ffIDbB1gW8Cz59i3G2Aw4Fdirz1bODalNKPgQ8Dpxd/ph0aDmvMkP8CbgWeC6wFvAr4VETsNENmG/n8rYv6jgA2AH5IbtDNa9jtxcBewMOB7YFXTve7kjQ9G0nS4DkrIhaSr4N/JvD/IE8dJv+F/daU0p0ppXvJf7kfNMU4+wBXpZS+nlJallL6JvBX4HkN+5ySUvpzEVz2IoeGrxb7Xwp8B3jRNLUeHxF3A7eTGzVvKrY/d6qxImI28ALgfSmlB1JKfwFOneYzlgLbRcRaKaW7Ukp/mGZfgF+nlM5KKa1IKT2YUvp9Suk3RR3XkkPW7jOMMWJX8iV7H00pLUkpXQB8nxz6RnwvpfTb4nd4Gg1nCSVJUv+KiCeTL187I6X0e+AfwEubGOIw4IsppYtTSstTSqcCi8n5YsTxKaXrU0oPTjPOHyLifuBK4GfA50uMvyswpxh/aUrpu8Bvpxh/OTCfnLfmppSuTSn9Y4Y/22iGLMb/QUrpHyn7Ofkk3FNmGGPEgcAPUkrnpZSWkpdKWBV4YsM+x6eU/p1SupN8km/HkmNLmoSNJGnw7J9SWkBe7+hw4OcRsRH5DM1qwO+Lab8LgR8X2yezCfkMUaN/kWcxjWhcxPthwBNGxi7Gfxl51tBU3pxSWpt8ZmgdYGQRxenG2oAcbBo/e7rFxF8APAf4V0T8PCJ2m2bfCWNFxNYR8f2IuDny5W4fJje9ytgEuD6ltKJh2/jf4c0Nzx8gN54kSVL/OwQ4N6V0e/H6f5h+FvV4DwOOHJeHNifnixEz3lAF2ImcLw4EngCsXmL8TYAb09g7M036WcUlc0cAxwC3RsS3ImKTyfadaqyI2DsifhN5SYOF5OzWTN4azaxF7roe85bUMTaSpAFVnFn6Lvks0ZPJs34eBB6dUlpQPNYuFl+czL/JAaPRQ4EbGz+m4fn1wM8bxl5QTD9+fYla/wh8EPhcMXNqurFuI9/prfHOHZtPM/bvUkr7kS8tOws4Y5Laxxwy7vUXyDOxtkoprUW+1j9m+jMV/g1sHmMXkRz/O5QkSQMmIlYlX061e3Ey6mbypfI7RMQO0x896nrgQ+Py0GrFLPERpW7BXcz0OQP4NfDeEuPfBGxa5LIR0+Wt/0kpjczASsBxM9Q3ur24VO875JlEGxYnRH/Iyrw1059xTGYtat4c85bUMTaSpAEV2X7kmT5XFmdnvkS+5vwhxT6bRsSzpxjih8DWEfHSiJgTedHH7ciXZk3m+8X+B0fE3OKxS+NihzM4lXznuH2nGyvl29V+FzgmIlYrrvd/xRS/g3kR8bKIWLuY6nwPMDI76BZgvYhYe4a61iyOu6/4rPGNsVvIi2VO5mLyWa93FH+GPciXBn5rhs+UJEn9bX/yybztyJdR7Qg8CvglU+QWJmaKLwGvi4gnFLlu9YjYJyLWbKOujwKvKWarTzf+r4v6Dy9y4H7A4ycbMCK2iYg9i4bQIvKJy8a8tUVMf2e2eeRL424DlkXE3sCzGt6fKbOdAewTeYHvucCR5Ev0LprplyGpNTaSpMHzvxFxH7n58SHgkJTSyK1l30leYPE3xWVaPwW2mWyQlNId5LWKjgTuAN4BPLdhevb4/e8l/6V/EPnM0M3ks1HzJ9t/kuOXAJ8B3lNirMOBtYvtXycvsLh4iqEPBq4t/ryvI18iR0rpr8Vx1xTTuaeagv128noG95ID1+nj3j8GOLUYY8xd6oo/0/PIi0TeTl6T4BXFZ0uSpMF1CPDVlNJ1KaWbRx7ACcDLYvI7rB1DQ6ZIKV0CvKY45i5yhntlO0UVs8B/ARw13fhFhjkAOBRYCLycfKJvsrw1n9ygup2czR4CvKt479vFzzsiYtJ1Kovc92ZyQ+gucu46p+H9aTNbSulvRX2fLWp4HnkR8yUlfiWSWhBjL3uVpP4TEccBG6WUmll3QJIkSSVFxMXAiSmlr9Zdi6R6OSNJUt+JiG0jYvtiGvbjyWfLvld3XZIkSYMiInaPiI2KS9sOId8c5cd11yWpfpNNqZSkXrcmeYrzJuTr5j8BnF1rRZIkSYNlG/LlZqsD1wAvTCndVG9JknpBT13aFhF7kddImQ2cnFL6aM0lSZIkDTwzmCRJKqtnGkkRMRv4O/BM4Abgd8BLUkp/qbUwSZKkAWYGkyRJzeilNZIeD1ydUrqmWGH/W8B+NdckSZI06MxgkiSptF5aI2lT4PqG1zcATxi/U0QcBhxWvHxcF+qSJEk1SilF3TUMODOYJEmaYKoM1kuNpFJSSicBJwFERG9clydJkjTgzGCSJAl669K2G4HNG15vVmyTJElS55jBJElSab3USPodsFVEPDwi5gEHAefUXJMkSdKgM4NJkqTSeubStpTSsog4HPgJ+dazX0kp/bnmsiRJkgaaGUySJDUjUurfS9y9Pl+SpMHnYtu9xwwmSdLgmyqD9dKlbZIkSZIkSephNpIkSZIkSZJUio0kSZIkSZIklWIjSZIkSZIkSaXYSJIkSZIkSVIpNpIkSZIkSZJUio0kSZIkSZIklWIjSZIkSZIkSaXYSJIkSZIkSVIpNpIkSVJPmANcCfy57kIkSZKGSLMZbE4Ha5EkSSotgG2BBFxRbDuovnIkSZKGQgDbzoZ0aUMGmyaE2UiSJEm1mwdcXDwP4DHF81XrKUeSJGkozJsHF/8aCIjHNGSwaUKYjSRJklSr+cAvgR1rrkOSJGmYzAd+GbDjTs0d5xpJkiSpNqsCFwC71F2IJEnSEFkVuGA+7HJB88faSJIkSbVZAhxbdxGSJElDZsmqcOxPgCc2f6yNJEmSVJu5wNvrLkKSJGnIzJ0Nb9+9tWNtJEmSpNrMAZ5RdxGSJElDpp0MZiNJkiTVYjXgtLqLkCRJGjKrrQantRHCbCRJkqRaLAW+W3cRkiRJQ2bpUvhuGyHMRpIkSarFUuA7dRchSZI0ZJYuhe+0EcJsJEmSpNosAt5QdxGSJElDpp0MZiNJkiTVZhlwSt1FSJIkDZl2MpiNJEmSJEmSpGGyGHhHa4faSJIkSZIkSRomS4ETWjvURpIkSarVUuD9wHLgfTXXIkmSNCyWLoX3v7/5DGYjSZIk1WoZ8JG6i5AkSRoyy5bBRz4AfLi542wkSZKk2i0DPg6sUXchkiRJQ2TZcvj4Mc1lsEgpdaqejouI/i1ekiRNMJvcVBqxM3BJSlFTOZqCGUySpAEyG2Z/Apa9ZeWmnXeGSy6ZPIM5I0mSJPWMBHyh7iIkSZKGyXJIR5bPYDaSJElSz1gBvBE4ue5CJEmShsiKBG/8crkMNqfj1UiSJDUhAa8B5tVdiCRJ0rBYAek18Jp5MO/g6Xd1RpIkSepJh9RdgCRJ0jBJwCvgkG9Nv5uNJEmS1LMW1l2AJEnSsHkJLFw49dvetU2SJPW05F3beo4ZTJKkwTdVBnNGkiRJkiRJkkqxkSRJkiRJkqRSbCRJkiRJkiSpFBtJkiRJkiRJKsVGkiRJkiRJkkqxkSRJkiRJkqRSbCRJkiRJkiSpFBtJkiRJkiRJKsVGkiRJkiRJkkqxkSRJkiRJkqRSbCRJkiRJkiSpFBtJkiRJkiRJKsVGkiRJkiRJkkqxkSRJkiRJkqRSbCRJkiRJkiSpFBtJkiRJkiRJKsVGkiRJkiRJkkqxkSRJkiRJkqRSbCRJkiRJkiSpFBtJkiRJkiRJKsVGkiRJkiRJkkqxkSRJkiRJkqRSbCRJkiRJkiSpFBtJkiRJkiRJKsVGkiRJkiRJkkqxkSRJkiRJkqRSbCRJkiRJkiQJAlg2/S42kiRJUk+YCzxIDidLi8dOtVYkSZI04ObD3KXw4FKYtRSWLoGls2CnaUKYjSRJklSb1YGFwHzg/uLnImBO8YjaKpMkSRpcqwMLV4X598L9c2D+HFg0B+bMgTkBMU0Is5EkSZJqE8BawN3kGUlR/JQkSVLnxOqw1l1w99zmM5iNJEmSVIu1gFvJwWV+zbVIkiQNi7XWgltvg5jfWgazkSRJkrpqfeA6cgNp1ZprkSRJGhajGSxg1TZCmI0kSZLUVQFsSg4ykiRJ6o6qMpiNJEmS1HWzyJe2SZIkqXtmLYC1/tnmGJVUIkmSJEmSpN52N7Bze0PYSJIkSV11B/CououQJEkaMncAj0rFkzbYSJIkSV21DnB+3UVIkiQNmaoymI0kSZLUVXcBz6m7CEmSpCFTVQazkSRJkrpqAfCduouQJEkaMgsoMti9wI6tj2MjSZIkddXdwMF1FyFJkjRkRjPYCuAfrY9jI0mSJHXVcuBvdRchSZI0ZKrKYDaSJElS143cefapdRciSZI0REYz2APAk1sbY06F9UiSJJWyHPg9sEbdhUiSJA2R0Qy2Avg1sDvw8+bGcEaSJEmSJEnSsFkB/BJ4enOHOSNJkiTV5kHgGQ2vZwM/qakWSZKkYTGawRJwIcx+NvykZAizkSRJkmqzHDi/eD4bOBvYB/hBbRVJkiQNvjEZLMHZ58E++8APSoQwG0mSJKknrABOBl5bdyGSJElDZAVwcoLX/hDYD4jp94+UUhfK6oyI6N/iJUnSBAG8oPh5OrAecGdKM8QZdZsZTJKkwTIhg60Dd945eQZzRpIkSeoZCTizeD4XuKvGWiRJkobFhAw2TQhzRpIkSeppyRlJPccMJknS4Jsqg83qdiGSJEmSJEnqTzaSJEmSJEmSVIqNJEmSJEmSJJViI0mSJEmSJEml2EiSJEmSJElSKTaSJEmSJEmSVIqNJEmSJEmSJJViI0mSJEmSJEmldKyRFBFfiYhbI+JPDdvWjYjzIuKq4uc6xfaIiOMj4uqIuCIidupUXZIkSYPMDCZJkjqpkzOSTgH2GrftaOD8lNJWwPnFa4C9ga2Kx2HAFzpYlyRJ0iA7BTOYJEnqkI41klJKvwDuHLd5P+DU4vmpwP4N27+Wst8ACyJi407VJkmSNKjMYJIkqZO6vUbShimlm4rnNwMbFs83Ba5v2O+GYtsEEXFYRFwSEZd0rkxJkqSBYgaTJEmVmFPXB6eUUkSkFo47CTgJoJXjJUmShpkZTJIktaPbM5JuGZkuXfy8tdh+I7B5w36bFdskSZLUPjOYJEmqRLcbSecAhxTPDwHObtj+iuLOIbsCdzdMv5YkSVJ7zGCSJKkSkVJnZiZHxDeBPYD1gVuA9wFnAWcADwX+Bbw4pXRnRARwAvkOIw8Ar0opzXj9vdOqJUkafCmlqLuGfmIGkyRJVZgqg3WskdQNhhhJkgafjaTeYwaTJGnwTZXBun1pmyRJkiRJkvqUjSRJkiRJkiSVYiNJkiRJkiRJpdhIkiRJkiRJUik2kiRJkiRJklSKjSRJkiRJkiSVYiNJkiRJkiRJpdhIkiRJkiRJUik2kiRJkiRJklSKjSRJkiRJkiSVYiNJkiRJkiRJpdhIkiRJkiRJUik2kiRJkiRJklSKjSRJkiRJkiSVYiNJkiRJkiRJpdhIkiRJkiRJUik2kiRJkiRJklSKjSRJkiRJkiSVYiNJkiRJkiRJpdhIkiRJkiRJUik2kiRJkiRJklSKjSRJkiRJkiSVYiNJkiRJkiRJpdhIkiRJkiRJUik2kiRJkiRJklSKjSRJkiRJkiSVYiNJkiRJkiRJpdhIkiRJkiRJUik2kiRJkiRJklSKjSRJkiRJkiSVYiNJkiRJkiRJpdhIkiRJkiRJUik2kiRJkiRJklSKjSRJkiRJkiSVYiNJkiRJkiRJpdhIkiRJkiRJUik2kiRJkiRJklSKjSRJkiRJkiSVYiNJkiRJkiRJpdhIkiRJkiRJUik2kiRJkiRJklSKjSRJkiRJkiSVMqfuAiSp170G2HfctrOBk2uoRZIkaViYwaTeZCNJkibxZmD34vn2wJbj3t8WCOBL3SxKkiRpwJnBpN5nI0mSJrELcMA0728JHEUOMid1pSJJkqTBZwaTep+NJElq8GbgccATS+y7FfD24rlBRpIkqXVmMKl/2EiSpMKbgbcCWzRxjEFGkiSpPWYwqb/YSJIk4E3A24CHtXDsSJBJeL2+JElSM6bMYAcBT2l4/SPg+2N3MYNJ9YiUUt01tCwi+rd4ST3lXOCZbY7xXeAFFdQiaayUUtRdg8Yyg0mqyoQM9jJgZ2BP8mrbIy4Gfg38HDhr7BhmMKkzpspgzkiSNPQOI5/RkiRJUvdMyGAvBd4FPHqSnZ9QPGBCI0lSd82quwBJqtN/Av9Fc9fkT2U78kk0SZIkTW/SDLYvkzeRGu0GPHfsJjOY1F02kiRJkiRJ9XoRuSM0kycAe3e4FknTspEkaWgdDKwLLK9ovGXAgxWNJUmSNKgmZLAXAO8FHlNygMcCz1750gwmdZeNJElDa33gLcAjKhrvP4BXVTSWJEnSoJqQwV5FDlJl7Qa8cOVLM5jUXTaSJA2t64BFdRchSZI0ZMxgUn+zkSRpaG0LrF53EZIkSUOmkgy2DbB7+7VIap6NJElD60PAP+ouQpIkachMyGAXkKcpNeMpwKurqkhSM2wkSZIkSZLq80ng93UXIaksG0mShtrPgNsqGutG4KKKxpIkSRpkP6PNDHYtcEl+agaTustGkqSh9t/Anysa61/ATysaS5IkaZC1ncF+Cnw2PzWDSd01p+4CJGkgbASP3B6eC/wO4E5Gz5JJkiRpBpcBqwGPATaZYd/rGNOFeiQNGUxSx9lIkjT0LgO2B9Zt5eCNyXcNeRr833vhfSPbLwXeBtwD/KH9GiVJkgbNZTRksGOLjV8HXj7NQdeT11T6zMpN/0dDBpPUcTaSJA29twKPBp7Z7IEbAW8Gjp7kvccCFwIXA7u2U50kSdJgmjSD/R34VcPrLYDNyA2kfwFnMHpJm6R6uEaSJAF/BRY2e9DzGW0i3QFcNdk+a5ATkiRJkiaYkME+ADyl4fEx8om5DxavxzWRpsxgkjomUkp119CyiOjf4iX1nM8ArwAWlNl5A/JptHfl5ZBOyk8n9wfgce3XJw2rlFLUXYPGMoNJqlJTGazBjBlMUlumymBe2iZJhbcACdgT2Jxpwsz6wJHAO3OAOZF855EprQpsCVxdUaGSJEkDpHQGa1Aqg0nqCGckSdIkjgMOY4og8w648zi4FTgTeE+ZAf9IXk1SUtOckdR7zGCSOmXaDEZuIDWVwSS1zBlJktSEdwIrgBc1bBu5q9udwFeBD3W7KEmSpAFnBpN6nzOSJKmkI4ufn3gH+XRZM/5Knq99U6UlSUPBGUm9xwwmqZtGM1itVUjDZ6oMZiNJkpr1BuDDwNpNHvdX4FHVlyMNOhtJvccMJknS4Jsqg83qdiGS1Pc+D7yv7iKasx4wu+4iJEmShowZTIPINZIkqRuWAQu781HrFB+1ABg5hfAH4AXANeP2vacoTZIkSe0xg2lY2EiSpG64EtitOx91GfAU4PfA+g3bfzfJvnsDP8UgI0mS1K7LGJfB1gJikgx2H+y93Aym/uWlbZLUiqXAvcCiEvsuBx7obDmNHk0OMuvPsB/Aj4BnA2vgtGtJkqR2jGaw1cnh6i/kKUrjH0+FH802g6l/udi2JLXjUODkad5fDlwMPKnzpaxKnkZ9E/kEWLP2Bf630oqkarjYdu8xg0nSShMy2K3ABjMc9EzgImAR7LvCDKbe5GLbktQJy8mzkkYeS8Ztv4iuNJEAbgTuo7UmkiRJklozJoPNZ+UCSdM5D7if3FCS+oyNJElqxynk01Ajj+eTm0nnF6+f2r1SNiTnkVbNoVzukSRJ0kqjGWwueVpSmfUFRsyBOWEGU3/x0jZJGhD3A6u1OcaLge+RJ1T5H1j1Ci9t6z1mMElaaTSDLQTWbmGA58OLz4bvJTOYeouXtkmSZnQGeR3xV9RdiCRJ0rD4HpxxgBlM/cNGkiRJkiRJdQq8vk19w0aSJEmSJEl1+jawAjik7kKkmdlIkiRJkiRJUik2kiRJE5wCvKbuIiRJkobMKZjB1PtsJEnSgFgduK/uIiRJkoZMpRnsq9hJUs+zkSRJkiRJUrsWAHe2OYYLbqsP2EiSpAGyFnB3RWN9EVgOHFHReJIkSYNqLeDuu4B12h/LDKZeZyNJkgZIaufgs4BlKx+xDGYtgziigsIkSZIGWAJYD7ir/bHiRDOYetucuguQJNXsXOBp5FML404vvBk44RPAPOBj3S5MkiSpj6yoaJxZZjD1NmckSdKwm0U+rTDJ3wifAhbNgqM+DCwGPtDVyiRJkobSaAZzzST1IBtJkjTMzifPRprCbPKJsA/PhnfPw781JEmSumA0gwHvrrkWaTz/l0CShtXPgN0p9TfBHOC9wDFHkRONJEmSOm7OB+C9/wXH1F2I1MBGkiQNmE0ouc7jfPLprpLmAnPnAkcC/6+FwiRJkgbYJsBdmwO3Vzjo3Jy/5lY4pNQuG0mSNGAeoM27t81kHrkJJUmSpFEPAKkTQexo8tRwqUd0rJEUEZtHxIUR8ZeI+HNEvKXYvm5EnBcRVxU/1ym2R0QcHxFXR8QVEbFTp2qTpEH3SKo9GTbircAdwEdfUzz5bAc+RFJbzGCSVJ9HArdvBdxS4aCrwFvnFxmswmGlVnVyRtIy4MiU0nbArsAbI2I7cj/1/JTSVuRlXo8u9t8b2Kp4HAZ8oYO1SdJAW0hnZiV9EdgG+MAqwLrA6h34EEntMoNJUk0WAuluYEWFg34MvvihIoNVOKzUqo41klJKN6WU/lA8vxe4EtgU2A84tdjtVGD/4vl+wNdS9htgQURs3Kn6JEnN+0/gDcD9dRciaUpmMEnqATsA/65orAfgPx8wg6l3dGWNpIjYAngscDGwYUrppuKtm4ENi+ebAtc3HHZDsW38WIdFxCURcUnnKpak/rcDsBlw00w7NuEU4DONG14MHF/hB0iqlBlMkrpvB2Cz2+Cmx5HD2D/bGOwjwCcnyWBSjTreSIqINYDvAEeklO5pfC+llGjy6ouU0kkppZ1TSjtXWKYkDZybgBuB5VPtcADw2+bGvAe4u3HD6uRL3CT1HDOYJNVjNIPdXDx5CrAF8NcmBnlvccxxwL2TZDCpRh1tJEXEXHKAOS2l9N1i8y0j06WLn7cW228ENm84fLNimySpDU8Ftgb+Nf6NzwCPgVcAvwBeB/xkkuM/Uhy/NfDJzpUpqUJmMEmq32gGu5EcxJ7NylB1ybgM9ibg+8WB7waOh4/8C7a+2wym3tPJu7YF8GXgypRS4z/35wCHFM8PAc5u2P6K4s4huwJ3N0y/liS16J/AVcDTge0aHpe/Aw7+C+wDnAicBtw77tgPkv+i2L0Y447JPuC5wJ/xDm5SjzCDSVJvmJDBroPtroLLr4KDXwj7bAcnbgenbQf3nkruKG0HfA4+eHeJDCbVJPLM5g4MHPFk4JfAH1m5Zv1/ka/RPwN4KLkv++KU0p1F6DkB2At4AHhVSmnaa/AjojPFS9IA+zZ5tvRS4KZfwIZPgdvITaSHkW/XdBH5/yxvBOaRbwF120wD3w18Azi8M3VreKWUou4a+okZTJJ605gMRl6orpIMJnXIVBmsY42kbjDESFLztgOuJf/fIr8gX7ff4OHk6/BbOvP1I+A5bRQnTcJGUu8xg0lS88ZksEm0lcGkDpgqg83pdiGSpHr9ZYb327mxiCRJkiZnBtOg6Phd21TeWqxcX20V4OfF4/xiWwA/635ZA+eb5BVFTwG2rLcUafA8AVeDlNR3zGDdYQaTpMHgjKQabQz8P+AtwLeAucATgfOA2eRV/iEvbnAeOcQ8tXgO8MxuFjsgTgP2JS8OsT2wLRMXF270UWBH8qIRnwe+1+H6pK57K/Al4LEVjbcu8JiKxpKkDpk0g60B532vyGBLgH1gxWw474cQzzaDtcsMJkmDwzWSavRI8iqYvwb2bOF4lyIp73+ABcAewKpNHHc5+f+LzwY+ziS3T5cGwU/JtxOpyu3AV4B3VjimhpprJPWegclg68GeXy82zmFlh2g5cC75LN6zgR8DiXybS8xgzZg0g32BvLLwiHeSv5AGZjBJqp+LbfeYTcnTep/R5jhnA/u3W8wAOw1YnXw2a34b41xOXhjvFOCsdouSes0Tgc8CO1U45i/I96uVKmAjqfcMRAbbAPg6uVFURiLfSmn//NIMNr1JM9iJwEbA08jXE474FStXF/4I+SYQvwN+bgaTpDrZSOox2wB/rWisbwMHkvONxtqPnBHXrGi8K4Bjge9UNJ7UM35ETvtVuAJ4M3mBEakCNpJ6T99nsA2BLzM6w6i0BJwJvDi/NINNbUwG+yKwDrA3sMYMB15EXkjp38C7gAvzZjOYJHXfVBnMxbZrsBn5uvyqvAj4Bn6Zk9kHmFfheNsD7wFeWOGYg2gHcvZTH/kY8IeKxroVm0iSetJmwP/bmNzYaLaJBPlSt4YQYAab2pgMtg/5lzVTEwnyLNnNyTdv+DD5F/xsM1hZZjBJ3eCMpC57KPAZOjMVeg75kn5lnwcOpdpG0ojLyNnm2x0Yu588HXjZJNs3Ax4N/KRh2/XA+7pRlFr3dHJ4fy2wa4tjXAG8g7FfvtQmZyT1nr7NYJvC/scDB7Q52FeLn6/OP8xgY43JYJ8GXgOs1saAFwPHAD82g40Yk8F2Jd804wuw2WPg0W9t+Gv4Wrj+WDOYpNZ4aVsP2IK8WOALOjS+IWase6jukrbJXIZB5nDy0jpl3EG+GmDEXXjGrGftSb4v88vJ61SU8XvgJHLH8EcdqktDy0ZS7+nbDPYfTFjUuS2zgGQGG280g30CeD3N3elkKr8BPgD80AwGcPjT4LMHFi+2JP9D/lPyCaHGleBvgzu+W2SwW4H3msEklTdVBpvT7UKG2UPoXBMJ8gmftwArOvgZWmlH8uzrYQ4xzViPPNFlxL3AMvI0dfWYC4rHVeSpZZDv2dx4v+ufAzcBBwGXkL/IH3exRklqwkOAF2xGx/7v+dOYwSZ1CNU0kSDPuilCw44/NIPxaMYGK8i3IxxvA1jvtcWuC4HlcO/7zWCS2mMjaYAcDrwNQ4w65znkG62MaOcmX2uSrwj4KPmuv+9tYyx1yIWMLnLK5Yy9XO0P5GlmlwBXYhNJUu/bAHhpZ4Y2g3XJrsXjh3UX0n2jGWxPYG/YqZUQtgB4K6y5Brz6XvjosWYwSS1KKfXtg3yTjL54PAzSVyGlDj7+G9KsHviz9srjSEgPdPh3fjGkfXvgz9qNx3MgXVTx7+9eSJ+CdHsP/Pl8+PDRu4+684aPAclgj+3AryLy32dmsLGP0Qx2BCndW/Hv/Jek9OwhzGDPIqUPkdKF1fwe773RDObDh4+ZH1P9R8Q1krrkScCvOvwZXp8/0V3kky+d9HngjR3+jLrtTZ7+vFuHxr8DWL9DY0vqf8k1knpO32WwTcl/kY2/FKhdrpE0pdEMdjOwYcWDHwccPUQZ7Fmw23uAJ1c48L3AJ+COB2H9j1U4rqSBMlUG826lUpt2Iv8lP6g63USCvHzC2zo4viRpyG1K9U0klXM8uWnRAUORwTrRRIK8xsBRsOo8M5ik5tlI6oLNyevRdsonyHci8br8iU4AHuzwZ6wCrNHhz6jLXuTr5jvZRIJ8R2AXfJQkVW1z4KBN6MzaSJ+ATyQz2FRGM9iHybdx64ChyGDPo/om0oj7YbUPmsEkNc/FtrvgkeRFGCtXDHoUkE4otr0BOBETTeE95H/I30xuVnTCRQzmXUOeBRwDPKHmOiRJatUjgcMfTr6lWpU+C7y9yGAVDz0oxmSwr5K/gzWr/YyBz2BPA7bv4AetCrwG+FIHP0PSQLKR1GGbAc+rcsDDGp4fD4xcsbiUnGQ+U7xeBnyl+Dnk3kX+B/31wOqtDLA7sE3D60vId6wqbAc8FfhFqwX2qM3p/PpSjeYDrwK+2sXPlCQNrsozGOSsdRK5O6IZjWaw98Dqs8kLGq1V3fgDn8EOIv8BO2UOsK0ZTFLzbCR12HZUdN1xAK8EvpgzzCkNbyXIs5BGfK7hmK+Qm0xD7ihgNvkv5RdS8oTY7sDDyfeof0rD9u8C/1s8/wvscQ8cvC784qLKyu0JVwG3MbaH1kmzgV0wxEiSqlFZBhuRIJ0Cp7xuzCbNYDSD/Re8cDasuQGwP7BO+2PvARzM4DWSupbBlgKXmMEkNc9GUgdtAjytioECeBm5KVR4dZnjTiTPSPo6sKSKQvrbSJhcCBzKNCfE9iB/eW8Cdp3k/QOKB8CPyHcj2QL4b+DXVVTaGx5D/jV0w2Lga+QrMyVJatdoBnsI8MwKBkzAacCrS2YwjTGawd5ZZLBbybPsK2gmDaKuZbB5sHh3+No3zWCSmmMjqYM2Ap7Y7iCzgBeRm0HkHHNmM8efTJ5BYyNp1NvIv9ZXAmtPtsPbgX1KDrY3XE2+2o2PkFeWvJGBaCh9DtgZeEQXPmsJ8IMufI4kaTiMZrDtgPe1OdgK4ExIBzeZwTTBaAY7GtaeRT4R9zRg/dbGG81gA+ZzwM67wCM6HcLugSWvM4NJap53beugG4EL2hlgFvB8WPEtOIfcRDobeHGz4+wN7Eu+AFoAHEGeAXP3+Dd2AzZsbqxzgS9CvhTu2+Tb6LlC9dRWB55RPJ8P7Atr7gun7As8p76yJEmDo+0Mthw4G1acDeecDenAFjOYJjiCIoO9g/wL/TL5l3tr82ONZrBBtBuwZQfHXwT8OC/3cEoHP0bSYHJGUgc9jPITWyY1DzgzZ5nnk5tJz29lnNOKnxuTL8MSkNfJDODlFAsaPp48o2inNgfejbzo+eH0/WmyPwLXkxd9bMtarJyetynwQfKqjuvC4tNy2L8f8n2CRy4bTMBP2v1gSdIwajmDLQPOI19z/fwKMpgmNSaDHV1sPJdylyFeDfytQ4X1kD8eD9c/GjY/bOZ9W3InLD6kIYNJUhNsJHXQb4F30+L/C8+C5U+DC8khZgXw3AprU/YmVi5BteA4mm4iXQf8dbI3ngB8AXgtY+7w1m8+SV7ose0MsxV5PalGxeu7aZiItGrDfkuBvWjzlLIkaRiNZrC7yCd1dp5m58XAL4vnDwD75exlBuusMRkM4NJiw/bkta2m8k3gq9NksAFRWQabzCLgV+MymCQ1wUZSB61LvvXs5cAOTRy3DLhoFVj8Q3hWRypTo8PJueU/yIsbNrPu43eBz0715s7Ax4E926muz61NDoSt3HZkLvns5MiK9b+cZl9JkhqMZrDLYYdXkRedWR14XMNOi4DLyOv0FDNhlgEXkXtLZrDOG5PB3llksA+S75a7DWOXG7gGuAG4Nr+cNoMNiqvI12luWuGYD5IXRTqwwjElDR0bSR20LbkxcQT5rNZMijtwcj/V3GBE5b0R4C9w7o7wzAXljrmR0SwztbWARwN/brWy+l1DviJyo2YOWkD+F2BH8syscRaRfyWPYYar/2aT7+m7onguSVIJYzLYn8jrGG7N2MVgbiZP9zjFDFanNxY/R69se3ex4b3kmckjTgD+Jz8tlcEGwDUfh5sfBhsdXtGAD8KiH8KfX1Qig0nSNCKlVHcNLYuIvij+ccCJ5BNhj5rk/SXAleSTDRtAXmT7iVQ/A8M1kmZ0wrnwhCJBPoJ8RrPRDaz8FX6NkmfCfkU+s9bHPks+a1jK2sArgOOn3uUa8j/ip5DXgp/RCvIaVgC/L1uIpEGRUoq6a9BYA5vBVJsTWHmvksoy2AD47FFw+IvJ/5Bu3MZADwI/gmte0GQGkzTUpspgNpK66D+AbzS8ngc8Evg58BLgO8AekG+fcE8HCrCRNLMTyV/GWvAxJk5r/xRwarNjDkAj6WjgoOL5hkwzO2lt4BDyYuNTWAycTwuLoCbgCvIsJ0lDxUZS7xnYDKaeUFkGGwCjGewjsOHRTc4QH/Eg8GNYfECLGUzS0LKR1IM2B77CuCnUs8inz35b8Yf9DXgycHvF4w6i71Pt37AD0Ehq9D7gmKnefCHw7cnfWkQ+m7iC1pZMIpH//ZA0dGwk9Z6BzGBSL3srvO9dcEyzU+cehEU/hhsOaCODSRpaU2Uw/7esRtczSYBZi+qbSNeQb0lvE6mcfwP31V1EH1qDaU+TXUu+id2/2vmMLds5WJKkbNIMJvWyTwEfIC/AXSbTP1js+2O49oAKMpgkNSjdSIqIdSNi/OXK6nX/It/S/q66C+kjhwGnkVfc1AQLyU2hO4vX91LkmX2ZdMGCRcX+q5CnZ7d1F5xz2zlYkvqTGUwSwMLPwrVbw53vAq6Fe2+fpKf0QH6Pc2DR1rmJVEkGk6QG0zaSIuKhEfGtiLgNuBj4bUTcWmzboisVqnU3khcFuLvuQvrQ68grOd6ADaVxPgM8HHgX+ddzCvDhhvcfZGWoWUy+Y+ETimPaDjADdImgJE3HDCZpvNEMdjLc8HA45Qj48A3kQHYDPHgD3P6dvNPigyrOYJLUYKYZSacD3wM2SiltlVLakrxk81nAtzpc2/AJ8krGVbiZfJtbL9Fq3RvIiyicRD67ozFOIv96fgscBdxUPE4HDiye/xB4D7knV4nfVTWQJPU8M5ikSY1msNPgqM3hpuJx+uZw4Cs6lMEkqcG0i21HxFUppa2afa9b+n2hxwnWo5p1jG4DHkq+pkjV+BjwSvL9g1dr8tgBW2x7UvPJv58Tx26eS/6VLaziM1xsWxpaw7jYthlMUjsqzWCShlari23/PiI+HxFPiIhNiscTIuLzwKXVl6mW3dnw2AybSFV7B/AQ4OPk3/GD9ZbTcw4g/27G2ZmKTpsnXOdL0rAxg0lqWWUZTJImMdOMpHnAocB+wKbF5huBc4Avp5QWd7zCaQzc2bBWZiSNrH+0LvmenuqOd5PvnDGTZcD5wF6dLacnvJS8SHknrABmd2hsST1vSGckmcEkSVKtpspgc2Y4aAnwheKhTkvkhZ1XL7HPiAWdLEhTWkJef2ou+bKuqfyS4WgiddL4f+YlaQiYwSRJUq+a6a5tcyLitRHxo4i4onj8KCJeFxFzu1Xk0LiTfM7xASZeOpVYueDzmg0P1eNj5N//keTvZcm495cX22s9X9xly6n+zzvSRFqr4nElqceZwSRJUq+a6dK2b5LXaDuVfGNJyCvwHAKsm1I6sNMFTmegp1WvRm4sjVgOrEFuTqxaS0WazisZu9D0D8nrBg2bFwLfrmiskebpGhWNJ6lvDemlbWYwSZJUq6ky2EyNpL+nlLZu9r1uMcRIPeYA4JvkuY7TXjg7A5tIkhoMaSPJDCZJkmrV6l3b7oyIF0XE6H4RMSsiDsR7KEka77vkNaOe18YYNpEkCcxgkiSpR83USDqIfLHKLRHx94i4CriZPO/goE4XJ6mPtXIXwZE1kWwiSZIZTJIk9aRpL20bs2PEegAppTs6WlETnFYt9bCnA+cCUTymk4rHfcDaHa5LUt8ZxkvbGpnBJElSHVq9tK1xgDuAtSLigIjYtrLKJA2m84HZwBNZ2Sgab2T7ncW+NpEkaQIzmCRJ6iXTNpIi4qyG5/sBF5BXPzknIl7Z0cokDYbfkP9LsyMrG0cjj5uK99avqzhJ6k1mMEmS1KtmumvbpSmlxxbPLwJellL6Z0SsD5yfUtqhS3VOVZ/TqiVJGnDDeGmbGUySJNWt1UvbGkPCnJTSP4vBbqe1pXQlSZI0MzOYJEnqSXNmeH+HiLiHvFTu/IjYOKV0U0TMI69oIkmSpOqZwSRJUk+atpGUUpoqqKwGvLb6ciRJkmQGkyRJvWraNZJ6ndfnS5I0+IZxjaReZwaTJGnwtbpG0pQi4o+tlyNJkqRWmMEkSVKdpr20LSIOmOotYKPqy5EkSZIZTJIk9aqZFts+HTiNsXcOGbFK9eVIkiQJM5gkSepR066RFBG/Bw5JKf1pkveuTylt3sniZuL1+ZIkDb5hXCPJDCZJkurW6hpJRwD3TPHe89spSJIkSVM6AjOYJEnqQd61TZIk9bRhnJHU68xgkiQNvrbv2hYRz53utSRJkqpnBpMkSb2kdCMJ2Hnc612qLESSJEmTMoNJkqSeMeOlbRGxLbAfsBn5ziE3AueklK7sfHnTc1q1JEmDb1gvbTODSZKkOrV0aVtEvBP4FhDAxcBvi+ffjIijqy5SkiRJZjBJktS7pp2RFBF/Bx6dUlo6bvs84M8ppa06XN+0PBsmSdLgG8YZSWYwSZJUt1YX214BbDLJ9o2L9yRJklQ9M5gkSepJc2Z4/wjg/Ii4Cri+2PZQYEvg8A7WJUmSNMyOwAwmSZJ6UJnFtmcBjwc2LTbdCPwupbS8w7XNyGnVkiQNvmG8tA3MYJIkqV5TZbAZG0m9zBAjSdLgG9ZGUi8zg0mSNPhaXSNJkiRJkiRJAmwkSZIkSZIkqSQbSZIkSZIkSSrFRpIkSZIkSZJKsZEkSZIkSZKkUmwkSZIkSZIkqRQbSZIkSZIkSSrFRpIkSZIkSZJKsZEkSZIkSZKkUmwkSZIkSZIkqRQbSZIkSZIkSSrFRpIkSZIkSZJKsZEkSZIkSZKkUmwkSZIkSZIkqRQbSZIkSZIkSSrFRpIkSZIkSZJKsZEkSZIkSZKkUmwkSZIkSZIkqRQbSZIkSZIkSSrFRpIkSZIkSZJKsZEkSZIkSZKkUmwkSZIkSZIkqRQbSZIkSZIkSSrFRpIkSZIkSZJKsZEkSZIkSZKkUmwkSZIkSZIkqRQbSZIkSZIkSSrFRpIkSZIkSZJK6VgjKSJWiYjfRsTlEfHniHh/sf3hEXFxRFwdEadHxLxi+/zi9dXF+1t0qjZJkqRBZQaTJEmd1MkZSYuBPVNKOwA7AntFxK7AccCnUkpbAncBhxb7HwrcVWz/VLGfJEmSmmMGkyRJHdOxRlLK7itezi0eCdgTOLPYfiqwf/F8v+I1xftPj4joVH2SJEmDyAwmSZI6qaNrJEXE7Ii4DLgVOA/4B7AwpbSs2OUGYNPi+abA9QDF+3cD600y5mERcUlEXNLJ2iVJkvqVGUySJHVKRxtJKaXlKaUdgc2AxwPbVjDmSSmlnVNKO7c7liRJ0iAyg0mSpE7pyl3bUkoLgQuB3YAFETGneGsz4Mbi+Y3A5gDF+2sDd3SjPkmSpEFkBpMkSVXr5F3bNoiIBcXzVYFnAleSw8wLi90OAc4unp9TvKZ4/4KUUupUfZIkSYPIDCZJkjopOpUTImJ78sKNs8kNqzNSSsdGxCOAbwHrApcCL08pLY6IVYCvA48F7gQOSildM8NnGHIkSRpwKSUXfm6CGUySJFVhqgzWsUZSNxhiJEkafDaSeo8ZTJKkwTdVBuvKGkmSJEmSJEnqfzaSJEmSJEmSVIqNJEmSJEmSJJViI0mSJEmSJEml2EiSJEmSJElSKTaSJEmSJEmSVIqNJKlGmwKp4XF3veVIkiQNBTOYJLXORpLUZdsAK4rH9ePeWwu4t+sVSZIkDT4zmCRVw0aS1AW7sDK4XAlEw2O81Yv9FnWtOkmSpMFkBpOk6s2puwBpGEwVWKbaV5IkSe0zg0lS9ZyRJHXQHsAy4Nc11yFJkjRM9sAMJkmdYiNJ6rDZtPYv2jzggYprkSRJGhbTZrBNgKXAzRPfMoNJ0vRsJEkd8gzgvDaOD2AVDDKSJEnNmDKDbQ0sLh7Xkhf5eAhwy9jdzGCSND0bSVKHzKL9RcgCmFtBLZIkScNiQgbbgbyC9h/J043msTJgBbABcOvYMcxgkjQ1G0lSB+wNnFPRWLOB+yoaS5IkaZBNyGA7A78B5pMbSJMJYH3gtrGbzWCSNDkbSVLFngd8h5xXqhDAqhWNJUmSNKgmzWCzyNepzWSSwGUGk6TJ2UiSKnQA8C2qDR3LgXUqHE+SJGnQTJrBngCc3/qYZjBJmpyNJKlCc4HVKh5zFvCPiseUJEkaJBMy2JOBc4E1Wh/TDCZJk7ORJPW4ANatuwhJkqR+MgdYq8ljVgOuW/nSDCZJk7ORJFXkhcCJdRchSZI0ZCrLYCOLbkuSpmUjSarIqsCCuouQJEkaMhMy2FOBM2opRZKGwpy6C5D63cuBDwCr112IJEnSEJkyg60CbNDCgA8Aj2qzKEkaAjaSpDatCWzRwfGXY6aRJEkab9IMtgfw5RYHTIxZI8kMJkmTs5EktejVwJF0/rawCbiqw58hSZLULybNYHsDHyffpW2zFga9D9hl7CYzmCRNzjWSpBatD2wHbNzBz1gGPA7y4o+XFQ9JkqQhNmkGuwg4GXhoi4OuAP668uVoBpMkTeCMJKkFhwJv6vSHzIV0EVwx8noH8qmx3zHhjJkkSdIwmDSDPQd4P5Xe9STRkMEkSWPYSJJasDGtzZouZRXgAvIspJ3HvTey7SLgiZ0qQJIkqTdNmsHWY2Jmasa9wDPaOF6ShoyXtkkt+BpwYicGXg04D9gNlu4Ke061327AzztRgCRJUu+akMH2Af6rjQHvBp5NnvFdWMo0GUySZCNJasVTgKd3YuDZwJNhMTnT/Gq6fZ9KnqEkSZI0JCZksI2AbdsYcCnw65UvS2UwSRpyNpKkFjwc2KrqQVcHzoRFwL7AhVWPL0mS1OcqzWB3AQetfGkGk6RybCRJvWBN4AzgWbAcOLfmciRJkgbeYuD8lS/NYJJUjo0kqRfMI99xpFmn4+VtkiRpeF0AnNDCcXcAr6m4FkkaEjaS1PNmkRdW7CVnAd+suwiAF2EjSZIkdURfZLB/Ar9tYaBFwPcrKEiShtCcuguQpjMbOBl4ObBsmv0OBVJXKsr+BFwGvKSKwdYCjq9iIEmSpGr0VQa7CPgC8PouFiJJQ8xGknrKKsCnG17PAl5ZPH/VNMctZWWIOZzpA09VfgI8DnhxuwOtCrw0P30QeHMzx74OWNFuAZIkadiNyWBzYNYJRQZbAa96A3DiJAe9rkcy2D+An1OukXQr8B7g/rGbm85gkjTEIqVunkOoVkT0b/Ga1DzgTcDH2xhjFfLaid3wFsY2vpq2NvAJ4FB4AHgrcFIzx8+iu6cBJakGKSUv4u0xZrDBMw9403z4+P8jT0d6Q/HGCuBz5IA23meLn0fk/WrNYFsBzwb+A3jtNAdezYTbvrWUwSRpCEyVwZyRpJ6xKnAsOci046OsPBv2HvIl8D1nPeBo8t3aDs2bFtNkgHk7NpEkSVLbVgWOXQXmHcfEhtGsSbaNGNm+FHgHfHR5jRnsquLxiOLnIxk7Q+km8sm7hRMPbTqDSdKQc0aSesLqwAfIZ4OqtAC4u+IxG+0CPAPYGTigmQMXAO8C3rFy013Aus2M4fl5SUPCGUm9xww2OFYHPrAqvPUDwJFtDHQcsBw4BljaAxnsocDLGl7fRl70aRJNZzBJGhJTZTAbSardGsB7gaM6MPZHybOcHuzA2I0OZcpsMtF65EUE1ma0c/YAuc7jmvlQ/7dK0pCwkdR7zGCDYTSDrQvcUdGgHwKOhY8u6cEMNomWMpgkDQkbSepZGwP/7uD4HwQ+Qg4KnbIzOYg9r8zO2wB/zWfpPldseoCcu0r7MPDfzRwgSf3LRlLvMYMNho2Bf69Gvtz+PRUOvCZwXw9msEJbGUyShoiNJPWkNcmzqN/X4c/ZkHyTjk56HHAQef3G/abaaT3gKFj4zrygeNPBJZGv7+/E9C1J6lE2knqPGaz/jWawDag+JP03eeGhz8GGi3okgxUW0mIGk6QhZCNJPWkz4PoufE43Gkkjtgf+c6o3/wP4I9xOnkbdtERe9FKShoiNpN5jBut/oxmsE42kEevBhnf2SAYrtJzBJGkIedc29Zw1GbsGYicdDHyBzk6tHnEF8OZODJzwliKSJKltXclgXwUWD0gGkySN4dwG1WY+sFOXPuvjwGHkO5P0tdfVXYAkSep3oxlsVfI1YZ1wCbBsgDKYJGmUjSTV5nbau8tss54CrNLFz5vUQuB/WzguAadVW4okSRpOoxlsPrB7hz7kc4x2j3oig0mSKmMjSUPjx3RnWvW0bqC1hbITeW64JElSVRYCr+/8x/REBpMkVcZGkmqzGvCsLn7eC+iRadX3AT+tuwhJkjSsRjPYfGCvDn3ID4Gl+WnPZDBJUiVcbFu1WRN4eVWDNYagc4FnMKFN+uoL4fbFVX1gG24kr3V0dYl9VwA/KX5KkiRVYDSDrQYc0qEPeSlwb376avLldJKkwRAp9e/dW731bP8bvfVsq55e/PwJMLt4/hzgu4y5GP8XwAEHwh13ABcBD7bzoRNtAjwKuAX4U5kDNga+Rj4T+JRp9lsKzGu3Oknqb1Pdelb1MYP1v9EMtgFwawc+YAFwd5HBgDs68BHQQgaTJJU2VQazkaRatdVIegrwc6DE/15sSENGOgg4h5abSQ8DHjpu23OBdwDnk287+5eyg60PfIc8N/CJDduXA/8HLGNls0yShpSNpN5jBut/oxlsAfADxuaQdv0f8Gzg/nEZrE2VZjBJ0oxsJKknPYScXXZu9sBdyTOLSv6vxYQQcyBwHXAFpVd/3ALYiHwL21dNs9/PyFeu/a3csNnawI8aXi8GntbMAJI0uGwk9R4zWP8bk8HWpbopQxcDezKar6poJG1BkcE2h1dtOskOtwH/aDGDSZKmZCNJPWtj4N/NHrSC0k2ky8h55q7J3nwpcDalmkmfBt5S7iP5JXAocFXJ/SVJU7OR1HvMYINhNINV1Uj6A/AkYFF+eRnTZLAmfHpzeMtDyF2i/5xkhx8AbwL+aQaTpCpNlcG8a5v6z/bN7f50pgkw/wNsN/3xmxUfuX4Tn/kU4JTiuEc2cZwkSVLXLaP9BYauAHZjtIkEM2SwEkYz2JHAJUzeRALYB/g88DAzmCR1g40k1W4pTU5BvpzSs5H+Sl5uaFqPIK/SuNrkbx9TfOTLyn3kqCcWx51BnpK9SvFRkiRJvWA0g90DPLmNgf4CPA5YsnJTqQw2g2OAyzeCl61bYue9gC8Bm5vBJKnTbCSpdrcz/Y3LWvEP8pTmxwN3z7Tz6eQAtEvFRRR2As4CnkVeBmn8IpGSJEl1GJPBlpMDVFlXAQm4mjz9Z1ne3FQGK+O9wMEl930meTrSlsCmZjBJ6hQbSeoJy8lrX1dle2Br4N6yB9zEmKnYI9YH1qygnh3ISzFtDVxAvuvIZhWMK0mS1I7RDHYfOUBdWzymcy051PwT2IYxU4+azmBTWB9Ycz2aD2J7kjtZ3wYeBjtsbAaTpKq52LZ6xgLy5fmT3Yxj1Oi9aie6oeH51sCDzXz4/uRTVUvGbv468PJmxmnCDeRJUDd3aHxJGhQutt17zGCDZQHjMtg8xs5O2hS4seH15mOPbyuDTeHrwMtPJq+c3Y5LgecVz5cBt5jBJKmsqTLYnG4XIk1lIXmpor+Tb/E6QQC/nbj5puLnFrRxLf5ZwB7AzyfW9ABTLp/Uls3I60Z6VkySJNVpIeMy2BJWNotmkU/ibT7xuEoy2DQ1VZLBHsvKTtdfgUeZwSSpXV7app5yL/my9tuKB+TL728DUoLbNln53m3Fe5sCm1B9gIF8J9lvdmDcEbPId9yVJEmq05QZbAWkTcfmr65lsPuoZnrTiNnAOvmpGUySWueMJPWc+4GHkGdV30QOJxsBtwLbrDP2NrK3V/Wh95BvXdJlGwN/Ln5KkiTVadoMRocy2HSOKIp5fUXjbQX8HngsbHy3GUySWuWMJPWsJcB65ECzAtgg4J//HLvP+uQzYm3bC7ho8rcWMWHpJEmSpIE1IYOR19VuVFkGm8YiYMkiYHGFgz4cuLzC8SRpCNlIUt9ICRYs6P7nHgl8rYPjB51Zg0mSJKkKibwgd7cdCXztbcCJFQ66gjz1CjOYJLXKRpL6zwMVj7eIHCqm8DngPyv+yEYPIa/9KEmSpJVGM9hSqluC4FrgcfmpGUySWmMjSf1ndaqbS70EeDJw8dS7LKPFRSSDfF3/+MfcsbvdSl7cUpIkSSuNZrCjgE9XNGjDpW1mMElqjY0kDa+lwC7kRRen8TrgpGbGDXKzaC/yNf3jH39kzDL3DyGfHJMkSdJKoxlsFtX9X8s1wKPzUzOYJLXGRpL6U7v3mV0O7Ahc0X4pE+xPnun0wyne34Y8j3p2fnkL+da5kiRJmsQnyQsmVeERwN/yUzOYJLXGRpL601xav7xtBbAd8JfqyhmjzL9VjwT+QZ69JEmSpO74BzmHSZJaZiNJ/auZRtKKYv8V5Ivh/96RiuBA4MyS+z4MuJ7cTLKhJEmSNLlmc9JI5pvsUdU6m5I0xGwkqX/NplwgSMDGwL+BLYB/drasZqRNyX+GW+uuRJIkqQd9HnhzE/sn4FxyTpzssfXK3SRJrZkz8y5SDyvWGWIZK9uiwdh0sD5wJ7BZh2s5GPja1G+PlBQNr28ENodcYwLuAdbuVIGSJEkDLAE/AJ43cTNMkcEkSU1zRpIGwxxW3tFjEbBGw+s7a6yrwaXA9g2vr2OSALMWK6deL+pOXZIkSQPhe0xoIkHJDCZJKs0ZSRo8q1Y73MnAq9s4/lfAUxpez3iZv+slSZIkdT+DSZJKcUaSVEKp4PF14CVjN13I2ADTlHnA4lYPliRJ6n/xBuAzzR/XVgaTJE3LRpI0g/8ETiy7c5UrNwYwF1hS4ZiSJEl9YjSDHQF8fIadTwde2OGCJEmAjSSplNL9oYYQcy7wjHY/OHAetiRJGlppwpNpdiz2qSSDSZKmZCNJKuFw4LNld16Rf6SVT1u3lMrXfJIkSeoXoxnsaOBD5Y6pJINJkqZkI0kqYQXwVuCTZXY+C3h+RR88B7inorEkSZL6zGgGWwG8j7yG5Dzg7eTr3kZev7yuCiVp+HjXNqmk5cA7yetfv2ua/c5KcNDSis6ELQPWqWIgSZKk/jSawZbDu5YXGz9Nvvx/2cr9zgIOwtlIktRpzkiSmrCMfDLs2CnePxM4aF9YfGa+Kq1tc4BbqxhIkiSpf03IYMsZ00Q6k9xEWkxFGUySNCUbSVKTlpIv0V+zeLyt2P5N4GBg8WxglYo+bBmwaUVjSZIk9bEZM1hNdUnSsImUqrxfeXdFRP8Wr4Exj7we9lLgAcjrI323osGXFh8gSUMspeT9K3uMGUy9YEIGkyRVaqoM5hpJUpuWFA8A9gdObXKAW4FHVViQJEnSEBiTwSRJXWMjSarSfPJc67JuArYH7uxMOZIkSZIkVck1kqS63AA8Fri97kIkSZIkSSrHRpJUl+XALXUXIUmSJElSeTaSpCp9H3h93UVIkiRJktQZNpKkKt0PfAN43Qz7XQfs2flyJEmSJEmqko0kqWr3kdc/ms4S4Jou1CJJkiRJUoVsJEmd8EvgDXUXIUmSJElStSKlVHcNLYuI/i1eg28tYNvi+XbA0cAriteLgCvqKEqS+k9KKequQWOZwSRJGnxTZTAbSVI3rAE8Eri87kIkqf/YSOo9ZjBJkgafjSRJktSXbCT1HjOYJEmDb6oM5hpJkiRJkiRJKqXjjaSImB0Rl0bE94vXD4+IiyPi6og4PSLmFdvnF6+vLt7fotO1SZIkDSLzlyRJ6pRuzEh6C3Blw+vjgE+llLYE7gIOLbYfCtxVbP9UsZ8kSZKaZ/6SJEkd0dFGUkRsBuwDnFy8DmBP4Mxil1OB/Yvn+xWvKd5/erG/JEmSSjJ/SZKkTur0jKRPA+8AVhSv1wMWppSWFa9vADYtnm8KXA9QvH93sf8YEXFYRFwSEZd0sG5JkqR+9Wkqzl9gBpMkSVnHGkkR8Vzg1pTS76scN6V0Ukpp55TSzlWOK0mS1O86lb/ADCZJkrI5HRz7ScC+EfEcYBVgLeAzwIKImFOc9doMuLHY/0Zgc+CGiJgDrA3c0cH6JEmSBo35S5IkdVTHZiSllN6VUtospbQFcBBwQUrpZcCFwAuL3Q4Bzi6en1O8pnj/gpRS6lR9kiRJg8b8JUmSOq0bd20b753A2yLiavI1+F8utn8ZWK/Y/jbg6BpqkyRJGkTmL0mSVIno55NOEdG/xUuSpFJSSt5FrMeYwSRJGnxTZbA6ZiRJkiRJkiSpD9lIkiRJkiRJUik2kiRJkiRJklSKjSRJkiRJkiSVYiNJkiRJkiRJpdhIkiRJkiRJUik2kiRJkiRJklSKjSRJkiRJkiSVYiNJkiRJkiRJpdhIkiRJkiRJUik2kiRJkiRJklSKjSRJkiRJkiSVYiNJkiRJkiRJpdhIkiRJkiRJUik2kiRJkiRJklSKjSRJkiRJkiSVYiNJkiRJkiRJpdhIkiRJkiRJUik2kiRJkiRJklSKjSRJkiRJkiSVYiNJkiRJkiRJpdhIkiRJkiRJUik2kiRJkiRJklSKjSRJkiRJkiSVYiNJkiRJkiRJpdhIkiRJkiRJUik2kiRJkiRJklSKjSRJkiRJkiSVYiNJkiRJkiRJpdhIkiRJkiRJUik2kiRJkiRJklSKjSRJkiRJkiSVYiNJkiRJkiRJpdhIkiRJkiRJUik2kiRJkiRJklSKjSRJkiRJkiSVYiNJatHbgU8ArwOWNjy+DTwH+GF9pUmSJA0sM5gk1StSSnXX0LKI6N/i1dfeBawN3Ax8HJjd8N6K4jELOB94Vterk6TBklKKumvQWGYw1cUMJkndM1UGs5EkNel9wLsbXs+ZZt8VwM+BPTtakSQNNhtJvccMpjqYwSSpu6bKYF7aJpX0/4D7gf8mB5eRx3RmAU8lBxmA7YBLO1WgJEnSADKDSVJvcUaSNIPPAa8A5hWPViwHHgT+Sg41D1ZTmiQNBWck9R4zmLph2gz2BuC4SQ76MfAF4EhgHzOYJLVjqgw2UzNfGnrzgTXaHGN2McZ8DDCSJEllTJrBjiJPTZoHrDrJQfsBewO/zS/NYJJUPWckSdP4AvAqcviowjLgHuCfwM4VjSlJg84ZSb3HDKZOmzSDvQt4D5M3kMZbCtwL/AF4phlMklrhGklSk74EvJLqmkiQpwCuS77biCRJkiaaNIO9jzwTqUwTCWAuOXTtDpxvBpOkKnlpmzSFNYFVOjDu1cBuHRhXkiRpEEzIYMcCbwVWb2GwucCC/NQMJknVcEaS1GUPA75XdxGSJEn94MPAEbS/YCVmMEmqio0kqcvmAhvVXYQkSVI/WECeotSORwM/NoNJUlVsJEmSJEkaXPOBTeouQpIGh40kaRLHA3t1cPyHAud1cHxJkqR+NCaDfQA4sNrxzWCS1D4bSdIkNqOzd/WYB2zVwfElSZL60ZgMtjH5VmsVMoNJUvtsJEk1uB54Xt1FSJIkDRkzmCS1z0aSVIMNgPfXXYQkSdKQMYNJUvtsJEmTeDdwQQfHXwXYqYPjS5Ik9SMzmCT1vjl1FyD1or8At1U54DbACStf/ht4GcDtwBuAk4EXVPmBkiRJ/afyDDbOaAaTJLXMRpI0iXcDT61ioO2Bj5JXjXziys0PAP8HsIg8x3qHKj5MkiSpv1WWwRpdTT5xR0MGkyS1zEaSNImdyDcKacuO5FlIT5pmn1WApwOLgbPI6eal7X6wJElSf6okg413D/CrqgeVpOFlI0nqhJ2ATzF9E6nRfGA/YCmwHDi4Q3VJkiQNk6uAd9ZdhCQNFhfblqq2M/AJppyX/W/gyKmOnQs8vxNFSZIk9ZkvAT9t4/i/AW9aOca0GUySVJqNJGkSJwC/beXAxwPHAXtM/vYNwOHAOa2VJUmSNNDGZLCLyTOKWrUA2CU/NYNJUnVsJEmTuAD4ZysHPhzYc/K3rgPeCnxvpjHmAZ9t5cMlSZL6W8sZbLy/Ae8F/tREBpMklWIjSeqS24Ezy+w4Fzi0s7VIkiQNtH8DJwFnNZHBJEml2EiSpvAN4PfNHPA44OWTv3Udee1tSZIkTW9MBvsu8MsmB/gb8IX81AwmSdWzkSRN4fvAlc0csB3w3ImbryXPrP5GBTVJkiQNujEZ7KfAscAvSh78V+B9wLfNYJLUKXPqLkAaVFeTg8t1wKk11yJJktS3fgoE8N/A7lPs8yfy9Wt/h6tPN4NJUifZSJKmcQ7wGGCHFo79B/D+asuRJEkaChMy2HlAKp6Pbyb9EfgQcHp+aQaTpM6ykSRN49vkE2C7kG/GtlOJY64Czgb+3sG6JEmSBtmkGeynxZu/G7fz74HTzWCS1C02kqQZnAHcCuxccv8/Akd1rhxJkqShMGkG+ykrG0rjmMEkqTtsJEklbEK+Acg6wJJi2y7j9vk7cCFwabsfthT4cruDSJIk9b+uZjBJUimRUpp5rx4VEf1bvPrSi4EHi+fPH/feRQfDyV+r4EPuB9aoYBxJGhAppai7Bo1lBlO3TZvBgJO7W44kDYWpMpiNJKkqBwPtNpKWAl8HDm2/HEkaFDaSeo8ZTJKkwTdVBpvV7UKkgfUv4DdtHL8E+CY2kSRJkiRJPctGklSVXwCfbuP4JYzetlaSJEmSpF5kI0nqBUuAXwKvrrsQSZIkSZKm5l3bpF5wD/CcuouQJEmSJGl6zkiSqnQb8Ocmj1kC/LoDtUiSJEmSVDEbSVKVLgDeCPy15P5LgPOAfTtWkSRJkiRJlYmU+vfurd56Vj3rycCXga2neH8R8CfgVmCfbhUlSf1pqlvPqj5mMEmSBt9UGcxGktQpuwFfBx45bvsi4KfA87pekST1JRtJvccMJknS4LORJNXh8cBXx237J/DcGmqRpD5lI6n3mMEkSRp8NpIkSVJfspHUe8xgkiQNvqkymIttS5IkSZIkqRQbSZIkSZIkSSrFRpIkSZIkSZJKsZEkSZIkSZKkUmwkSZIkSZIkqRQbSZIkSZIkSSrFRpIkSZIkSZJKsZEkSZIkSZKkUmwkSZIkSZIkqRQbSZIkSZIkSSrFRpIkSZIkSZJKsZEkSZIkSZKkUubUXYA0nZhie+pqFZIkScPFDCZJmoqNJPWs+cCiSbYnnEonSZLUKWYwSdJ0/LtAPWcNclCZLMBAPkPm2TBJkqRqmcEkSWX0+4yk24H7i5/qHevTxndyH1NPp1bL2vpO1DF+L73H76T3PKzuAjSp+4C/1V2EJjCD9R7/XulNfi+9x++k90yZwSKl/j6vEBGXpJR2rrsOreR30nv8TnqT30vv8TuRyvHfld7k99J7/E56k99L7/E76S9e2iZJkiRJkqRSbCRJkiRJkiSplEFoJJ1UdwGawO+k9/id9Ca/l97jdyKV478rvcnvpff4nfQmv5fe43fSR/p+jSRJkiRJkiR1xyDMSJIkSZIkSVIX9G0jKSL2ioi/RcTVEXF03fUMk4j4SkTcGhF/ati2bkScFxFXFT/XKbZHRBxffE9XRMRO9VU+uCJi84i4MCL+EhF/joi3FNv9XmoSEatExG8j4vLiO3l/sf3hEXFx8bs/PSLmFdvnF6+vLt7fotY/wACLiNkRcWlEfL947XciNcEMVg/zV28yg/UeM1jvMoMNjr5sJEXEbOBzwN7AdsBLImK7eqsaKqcAe43bdjRwfkppK+D84jXk72ir4nEY8IUu1ThslgFHppS2A3YF3lj8O+H3Up/FwJ4ppR2AHYG9ImJX4DjgUymlLYG7gEOL/Q8F7iq2f6rYT53xFuDKhtd+J1JJZrBanYL5qxeZwXqPGax3mcEGRF82koDHA1enlK5JKS0BvgXsV3NNQyOl9AvgznGb9wNOLZ6fCuzfsP1rKfsNsCAiNu5KoUMkpXRTSukPxfN7yf+B3hS/l9oUv9v7ipdzi0cC9gTOLLaP/05GvqszgadHRHSn2uEREZsB+wAnF68DvxOpGWawmpi/epMZrPeYwXqTGWyw9GsjaVPg+obXNxTbVJ8NU0o3Fc9vBjYsnvtddVkx9fOxwMX4vdSqmL57GXArcB7wD2BhSmlZsUvj7330OynevxtYr6sFD4dPA+8AVhSv18PvRGqGf3/0Fv+e7yFmsN5hButJn8YMNjD6tZGkHpbyrQC9HWANImIN4DvAESmlexrf83vpvpTS8pTSjsBm5LP429Zb0XCLiOcCt6aUfl93LZJUNf+er5cZrLeYwXqLGWzw9Gsj6UZg84bXmxXbVJ9bRqblFj9vLbb7XXVJRMwlB5jTUkrfLTb7vfSAlNJC4EJgN/IU9jnFW42/99HvpHh/beCO7lY68J4E7BsR15Ivx9kT+Ax+J1Iz/Pujt/j3fA8wg/UuM1jPMIMNmH5tJP0O2KpY5X0ecBBwTs01DbtzgEOK54cAZzdsf0Vxh4pdgbsbpvmqIsU1w18GrkwpfbLhLb+XmkTEBhGxoHi+KvBM8roJFwIvLHYb/52MfFcvBC4ozmCqIimld6WUNkspbUH+e+OClNLL8DuRmmEG6y3+PV8zM1jvMYP1HjPY4Il+/T4i4jnk6yxnA19JKX2o3oqGR0R8E9gDWB+4BXgfcBZwBvBQ4F/Ai1NKdxZ/uZ5AvsvIA8CrUkqX1FD2QIuIJwO/BP7IyuuO/4t8jb7fSw0iYnvyIoGzyU37M1JKx0bEI8hnYtYFLgVenlJaHBGrAF8nr61wJ3BQSumaeqoffBGxB/D2lNJz/U6k5pjB6mH+6k1msN5jButtZrDB0LeNJEmSJEmSJHVXv17aJkmSJEmSpC6zkSRJkiRJkqRSbCRJkiRJkiSpFBtJkiRJkiRJKsVGkiRJkiRJkkqxkSSpqyJir4j4W0RcHRFHT/L+UyPiDxGxLCJeOO695RFxWfE4p2H7wyPi4mLM0yNiXrF9fvH66uL9LTr+B5QkSepBZjBJVbGRJKlrImI28Dlgb2A74CURsd243a4DXgn8zyRDPJhS2rF47Nuw/TjgUymlLYG7gEOL7YcCdxXbP1XsJ0mSNFTMYJKqZCNJUjc9Hrg6pXRNSmkJ8C1gv8YdUkrXppSuAFaUGTAiAtgTOLPYdCqwf/F8v+I1xftPL/aXJEkaJmYwSZWxkSSpmzYFrm94fQOwaUQcGxH7TnFMo1Ui4pKI+E1E7F9sWw9YmFJa1jjm+M8r3r+72F+SJGmYmMEkVWZO3QVIUkrpvSV3fVhK6caIeARwQUT8kRxMJEmS1CQzmKRWOCNJUjfdCGze8HqzYlspKaUbi5/XAD8DHgvcASyIiJHGeOOYo59XvL92sb8kSdIwMYNJqoyNJEnd9Dtgq+IOH/OAg4BzZjgGgIhYJyLmF8/XB54E/CWllIALgZG7ixwCnF08P6d4TfH+BcX+kiRJw8QMJqky4b/PkropIp4DfBqYDXwlpfShiDgWuCSldE5E7AJ8D1gHWATcnFJ6dEQ8EfgieQHIWcCnU0pfLsZ8BHnRyHWBS4GXp5QWR8QqwNfJZ83uBA4qzqRJkiQNFTOYpKrYSJIkSZIkSVIpXtomSZIkSZKkUmwkSZIkSZIkqRQbSZIkSZIkSSrFRpIkSZIkSZJKsZEkSZIkSZKkUmwkSZIkSZIkqRQbSZIkSZIkSSrFRpIkSZIkSZJK+f9tf2kbrN23jwAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 1440x4320 with 2 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "c='9'; fov = '001'\n",
    "ref = imread(f'tif/Cycle_0/Cycle_F{fov}.tif')\n",
    "mov = imread(f'tif/Cycle_{c}/Cycle_{c}_F{fov}.tif')\n",
    "\n",
    "ref_max = ref.max(0)   # max proj of z for each channel in reference\n",
    "ref_binary = ref_max[0] > threshold_otsu(ref_max[0])  # binarizing for channels - T/F\n",
    "\n",
    "mov_max = mov.max(0)  # max projection by each cycle in for loop (ref above)\n",
    "mov_binary = mov_max[0] > threshold_otsu(mov_max[0]) # binary of moved image\n",
    "\n",
    "out = imread(f'reg_bin_Cyc_{c}/Cycle_{c}_F{fov}_bin_reg.tif')\n",
    "\n",
    "fig, axs = plt.subplots(1,2,figsize = (20,60))\n",
    "axs = axs.ravel()\n",
    "\n",
    "im_reg = np.zeros((2048, 2048,3)) # empty color image\n",
    "im_reg[...,0] = ref_binary\n",
    "im_reg[...,1] = out.astype('bool')\n",
    "    \n",
    "im_orig = np.zeros((2048, 2048,3))\n",
    "im_orig[...,0] = ref_binary\n",
    "im_orig[...,1] = mov_binary\n",
    "    \n",
    "axs[0].imshow(im_orig[:500,:500])   #before reg\n",
    "    \n",
    "axs[1].imshow(im_reg[:500,:500])   #after reg\n",
    "\n",
    "axs[0].title.set_text('Before Registration')\n",
    "axs[1].title.set_text('After Registration')\n",
    "for ax in axs.flat:\n",
    "    ax.set(xlabel='0:1500', ylabel='0:1500')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "bbe6c7fc",
   "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"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}