DEMO 1.3 - Update a stratified design - option 2.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<div class=\"alert alert-block alert-info\">\n",
    "<b>Note:</b> To run this notebook you will first need to generate a stratified design, using the <b>demo-stratified-design</b> notebook.\n",
    "</div>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Update Stratified Design Option 2: Exclude radius around specific sites\n",
    "**In this notebook we show how a stratified design can be adapted, based on new information gained in the field. Sites which have already been sampled are left in place, while unsampled locations will be adjusted based on the new constraints imposed by the user.**\n",
    "\n",
    "In this option we mask a radius around a specific site (or sites) which cannot be accessed. This information can be coded directly into the sample design csv, to avoid updating the mask manually (as in option 1)."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Import required libraries"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "import os\n",
    "from osgeo import gdal, osr\n",
    "import numpy as np\n",
    "from random import randint\n",
    "from scipy import ndimage\n",
    "from copy import copy\n",
    "import time\n",
    "import pandas as pd\n",
    "from matplotlib import pyplot as plt\n",
    "%matplotlib inline\n",
    "\n",
    "def get_file_info(file_path):\n",
    "    '''\n",
    "    Function which extracts a geotiff file as numpy array, and saves geographic projection information\n",
    "    input: \n",
    "        file_path: (str) Path to the file\n",
    "    output: \n",
    "        file_map: (np.array) The extracted map\n",
    "        nbins: (int) The number of categories in the map\n",
    "        res: (float) Resolution of the map in meters\n",
    "        GeoT: The geographic transform used to project the map\n",
    "        auth_code: The authority code matching the projection\n",
    "    '''\n",
    "    file_raw = gdal.Open(file_path)\n",
    "    prj_info = file_raw.GetProjection()\n",
    "    GeoT = file_raw.GetGeoTransform(); res = GeoT[1]\n",
    "    file_map = file_raw.ReadAsArray()\n",
    "    nbins = len(np.unique(file_map))\n",
    "    return file_map, nbins, res, GeoT, prj_info"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# User Inputs"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<div class=\"alert alert-block alert-info\">\n",
    "<b>Note:</b> Sections where users can enter their own values are highlighted with #'s. \n",
    "    \n",
    "<b>Example files currently entered link to datasets in the raw data folder, which can be used for testing. </b>\n",
    "</div>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Specify a unique name for the project as ```save_folder```. A folder with your chosen name will be created in the ```results``` directory, and all outputs for the design will be saved there. E.g:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "############################################\n",
    "save_folder = 'Adapted_Stratified_Demo_Opt2'\n",
    "############################################"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Results will be saved in the folder results/Adapted_Stratified_Demo_Opt2\n",
      "Files will start with the time stamp 2019_10_30_140158\n"
     ]
    }
   ],
   "source": [
    "directory = 'results/{}'.format(save_folder)\n",
    "if not os.path.exists(directory):\n",
    "    os.makedirs(directory)\n",
    "    \n",
    "# Create unique time stamp to add to file names (to avoid overwriting)\n",
    "ts = time.gmtime()\n",
    "ts = time.strftime(\"%Y_%m_%d_%H%M%S\", ts)\n",
    "\n",
    "print('Results will be saved in the folder {}'.format(directory))\n",
    "print('Files will start with the time stamp {}'.format(ts))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Enter path to **original** invalid areas mask (`mask_path`) used to generate the stratified design . This should be in GeoTiff format and saved in the 'raw' data folder. E.g:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [],
   "source": [
    "######################################\n",
    "mask_path = 'raw/InvalidAreasMask.tif'\n",
    "######################################"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Load the .csv file output by the original stratified design. In the \"sampled\" column, **unsampled** sites should be labelled with a **zero**. Sites which have been **visited** and are known to be accessible should be tagged with a **one**. The site which requires masking shoud be labelled with a **two**. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "###############################################################################################\n",
    "csv_path = 'tagged_stratified_output2.csv'\n",
    "###############################################################################################"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Enter the radius (in meters) which you would like to exclude around your selected sites (tagged with 2's). E.g:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "#############\n",
    "radius = 3000\n",
    "#############"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Load data"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "To project the final sites into longitude and latitude, we extract geographic information from the invalid areas mask."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Loading raw/InvalidAreasMask.tif ...\n",
      "\n",
      "Number of categories in mask (should only be two): 2\n",
      "Pixel resolution (m): 30.0\n",
      "GeoT info: (604343.5084861958, 30.0, 0.0, 5302852.190465175, 0.0, -30.0)\n"
     ]
    }
   ],
   "source": [
    "print('Loading {} ...'.format(mask_path))\n",
    "\n",
    "mask, nbins, res, GeoT, prj_info = get_file_info(mask_path)\n",
    "\n",
    "print('\\nNumber of categories in mask (should only be two): {}'.format(nbins))\n",
    "print('Pixel resolution (m): {}'.format(res))\n",
    "print('GeoT info: {}'.format(GeoT))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Load the tagged csv file"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {},
   "outputs": [],
   "source": [
    "sampled_csv = pd.read_csv(csv_path)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Update the mask based on information in the csv**"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {},
   "outputs": [],
   "source": [
    "def update_mask(site_df, mask, radius, res):\n",
    "    \n",
    "    # Create array same dimensions as input mask\n",
    "    imheight, imwidth = mask.shape\n",
    "    mask_update = np.ones((imheight, imwidth))\n",
    "    \n",
    "    # Set the point to mask as a zero\n",
    "    center_pixel = site_df.loc[site_df['sampled']==2]\n",
    "    x = center_pixel['row'].values; y = center_pixel['col'].values\n",
    "    mask_update[x, y] = 0\n",
    "    \n",
    "    # Threshold distance transform and add to original mask\n",
    "    dist_im = ndimage.distance_transform_edt(mask_update)*res\n",
    "    dist_im[dist_im<radius] = 0; dist_im[dist_im>=radius] = 1  \n",
    "    \n",
    "    new_mask = dist_im * mask\n",
    "    \n",
    "    return new_mask"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "metadata": {},
   "outputs": [],
   "source": [
    "updated_mask = update_mask(sampled_csv, mask, radius, res)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Plot updated mask**"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAuYAAAGUCAYAAAB5m+cEAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvnQurowAAIABJREFUeJzsvXvYLkdZ5ns/kg8QlCRkBVxJENBECDoXiNmA4iASUSJgmIgIoglO3NmOqKjoiM4oxNOg25HDwEaXBAnIOUbJZiPCcJCtAwwJIqfAEBHIIguSBTkAEVjBmj+6O6u++qq669hV3X3/ruu7vvftQ9XT1VV3VT/1VL+ilAIhhBBCCCGkLl9T2wBCCCGEEEIIB+aEEEIIIYQ0AQfmhBBCCCGENAAH5oQQQgghhDQAB+aEEEIIIYQ0AAfmhBBCCCGENAAH5oQQoiEiLxKRa0XkAx7H3lFEPisixxrb/0pEHjdy3kNF5HX95x8Skac5jvtCqP2EELIl1qbZHJgTQshuXgzgET4HKqW+COCNAB4zbOsF/7sBvM4zjcuUUs8MN5MQQghWptkcmBNCiIZS6u0APhdwyisAPF77/u8AvEEpdbOIPEBE/oeI/EP//17mySLyJBF5Xv/5niLyDhF5t4j8dtKFEELIBlibZh+TIxFCCCnJqd95qrr5hpuzpHXow4c+COBL2qYDSqkDU+eJyE8DgFLqj41dbwDwQhE5QSn1WXSC/9/6fR8G8BCl1C0i8n0Afg/AD49k8xwAL1BKvUREnux3RYQQ0hbU7Hg4MCeENM/NN9yMCy6+IEtaFz7wwi8ppc4IPc8i7sP2r4jIZQAeKyJ/AeB+6KZKAeBYABeLyGkAFICdiWwejKOdwEsB/H6onYQQUhtqdjwMZSGEkHSGqdHHAnitUupIv/23AbxVKfVtAB4N4PYeaakyJhJCCOlpVrM5MCeEkHTeCuA0AE9GJ/gDxwL4VP/5SR7p/D2Oxj4+MZdxhBBCdtGsZnNg3gAi8usi8sLcx3qkpUTk1Bxp5bBBRP5YRH7D59itISIvFpHfqW3HFhCRVwB4B4B7ichBETm/3/7TQ8yiiVLqXwH8BYATALxd2/UHAP6LiPw9gNt4ZP8UAE8WkXej6yDIAugXg/3dTHndo9fCaqGopg0i8tcicp7PsVtDRN4mIj9V2441szbN3mRDKYmIPAnAUwF8M4CbAPwlgF9TSt3gOkcp9Xu+6Yccm4KIvA3AnyulsjwE+KCUsjagEETkxQB+HMA3KqWuSTbKP9+3AfgeAPdTSv2jtv2vAJwN4HuVUm+byx4Sj1LqCY7t1nhFbf9T0Im0vu0dAL5F2/Qb/fa3AXhb//nF6F73BaXUPwP4Tu14vkaxACKiAJymlLpK2/YMAKcqpX68cN5F8xGRjwP4KaXUfy+Rvg2l1FmpafQael8A36CU+nKyUf75fhzASQBOUkod1ra/t7fnnkqpj89lDwlnbZpNj3lGROSp6AL/fwXdk9ODANwdwJtE5LaOc/hwlAkRuSO6RRg3YmJKqVC5/y8A52p5nICuDlxXIC9CCFkFInIPAP8WXazuD00c6+PFDOWfAdw6uBORfwPgawvkQ8gkHJhnQkTuBOBCAD+nlHqDUupI/5T9OHSD8x/vj3uGiFwiIn8uIjcBeFK/7c+1tM4VkU/0v071GyLy8f61PdCP1aYIzxORT4rIYRH5T1o6D+jfr3mDiBwSkee5HhAmru2h/fTQU6X7da1DIvKT/b4HicindbEUkX8nIu8LtcEM1xCRX+nPuUZE/r2HqT8M4AYAvwVg17Sqo9y/RkSeJiL/1Jf1q0Xkzto5r+mv7UYRebuIfOtE/i8D8KNaWTwB3YzJV7Q0neUhHc/qy/hGEXmfiHybpZy+XkTeKiLPFRHxKBdCSACa5v16r6sfF5EnavtPEJHLROQmEfmf6GZI9fOfIyJX9/uvEJF/229/BIBfR6cTXxCRf+y3HysiF/Wa8CkR+Z1BR0TkNiLyh70dHwPwyIDreJKI/F1//vUi8s8icla/7/Eicrlx/C9K97YKiMgjpXuX8039tTxjJJ9bwzUi7T0XwDvReSFN7X6xiLxARF4vIl8E8L0icrs+j0+KyGekC4P82v7440XkdSJyXX/NrxORUybyfyk0p0pvw0sMO5zlISK37/uWz/ba/m4RuaulnPb3uv7LHmVCNgoH5vn4LnSrdy/VNyqlvgDgrwE8XNt8NoBLAByHbjB3KyJyHwD/DzqP7350nveTJ/L+bgD3AnAmgN8UkdP77V8F8IsA9qGbajkTwM8EXtfAN2i2nA/g+SJyvFLqnQC+COBh2rE/BuDlKTb0Hdgvoyu30wB8n4eN56FbxPFKAPcWkfsb+81y/3l0v/71PeimMq8H8Hzt+L/u874LgPfAuFcWrgHwIQDf338/F4a4Y7w8vh/AQ9BNox0H4EcBfFY/WTov/JsB/L1S6ueVUnyDByFl+AZ07fRkdNpyQI7+2Mjz0b1XeT+Af9//6bwb3SvY7oxOC18jIrdXSr0B3XuRX6WU+jql1H374y8GcAuAUwF8OzotGOKS/08Aj+q3n4HuLRIhPBDAR/pr+QMAF/UP9Jehi8k9TTtW1+4votOw49ANrv+DiDwG08TYey46fX0ZgB+wDGp/DMDvAvh6AH+Hbmb6W9CV8ano7tFv9sd+DYA/Q+cQ+0YA/wLgeRP5vxPAnUTk9P6B6EcB/LlxzFh5nIeuf7wbupjln+7zvRXpZgX+FsDzlFJ/OGEP2TAcmOdjH4DDSqlbLPsO9fsH3qGU+iul1L8qpf7FOPaxAP5fpdTfKaW+gk5spgZfFyql/qWPbf5HdHFxUEpdoZR6p1Lqlt57/yfoBqExHAHwW/1MwOsBfAHdwwDQDYafAHTeXAA/2G9LseFxAP5MKfWB/id0nzF2sIh8I4DvBfBypdRn0A1ezcVIZrn/XwD+k1LqYB/T+Ax07zU9prf9RUqpz2v77ivdT/eO8RIA5/Yd+HF9vNqtTJTHEXQdz70BiFLqSqXUIe30k9AJ+2uUUv95wg5CSDq/oZT6slLqbwH8fwAe1w/cfhjAbyqlvqiU+gC6gfWtKKX+XCn12b6d/1cAt8NRvdxFPwg9C8Av9OldC+BZOPqmh8cBeLZS6mql1OcA/JfAa/iEUupPlVJf7e3cD+CuSqmbAbwWR7X7NHTac1l/DW9TSr2/18v3odN0X+32tldEvhvdIPrVSqkrAPwTuoG4zmuVUn/fL9j7MrrB/y8qpT6nlPo8uoedx/d2f1Yp9RdKqZv7fb/raffgNX84uh+Z+ZS+c6I8jqAbkJ+qlPpqr/M3aaffB1188tOVxw/jkG3DgXk+DgPYJ/bY5f39/oGrR9I5Sd/fi+dn3YcDAD6tfb4ZwNcBgIh8Sz+N92npwjd+D7sfEEL4rPHQcWs+6Dws54jI7QCcA+A9SqlPJNqwqxwAfGLi+J8AcKVS6r3995cB+DER0X8cwCz3uwP4y37q8QYAV6LzaN+1n459pnRhLjcB+Hh/zpTtl6KbPfg5dEK/i7HyUEq9BZ1n5/kAPiMiB6QLkRp4JLq4x9EFLYSQSb6KvT8csoNugDVwfe8UGPgEOl06Ed2LE5z6JF3Y35XShaTdgM6b6tKOu/d5H9K06E/QzdQB4Vpocmv/0PcnwG7tHmKrfwzAXw3HiMgDpQuZu05EbkTnBS6h3ecBeKO28PLl2OtU0dM7EcAdAFyhldcb+u0QkTuIyJ9IFw56E7o3bhwn07HpL0VXBk/C3pnOqfJ4KYC/AfBK6UIv/8Doe56IbqB/yYQNhHBgnpF3oHuSP0ffKN2CxLPQeXAHxjzghwDcGg/Xx82dEGnTC9A9+Z+mlLoTutjG7DHJSqkPoRPfs7B7KjTFhkPopgUHvnHi+HMBfFM/4P00gD9CJ5r62wLMcr8awFlKqeO0v9srpT7VX8fZ6EJojgVwj/6cUdv7Tu2vAfwHWAbmmCgPpdRzlVLfAeBb0U3V/op27p+i64Be39crQkgcn8TRNj1wT+weRB5vtLNvRBeudh26sBOrPkkXT/6r6DzHxyuljkO3IH1o5zYd+jKAfZoO3UkpNaxpCdXCEN6IzqF0P3QDdF27X47Oe343pdSx6BwCWbW7798eB+B7NO3+RXSzk/fVDtXL7DC6MJFv1crrWKXU8LDxVHSzEw/sNfYhQ3ZjRvfOpH9GN+N7qeUQZ3n0M8kXKqXugy6s9VHYHbP+jN7ul3s8IJCNw4F5JpRSN6Jb/PnfROQRIrLTx5S9BsBB2AdpNi4B8GgR+S7pFgVeiPjB9Neje2XjF0Tk3ugGi6V4ObqY7Yegu+ZUG16NboHmfUTkDgCe7jpQRL4T3eKrB6CLObwfgG+D3fOi88cAfldE7t6nc6KInK3Z/WV0sxV3QOfZ9uXXAXyPsr9iy1keIvJ/9F6ZHXTxjF9C59nT+Vl08aKvGxY7EUKCeRWA/ywip0i3CPz70P3Kn+nRvFBEbtsPth+FLozsq+gGbs/ovbP3wW6d+Xp0A/frABwjIr8JQJ/5+gyAe4jI1wBAH672RgD/VUTu1NvzzSIyhEm8GsDP97YeD+BpuQqhnwW9BMD/jS4e/k3GdXxOKfUlEXkA9oaXuAix9zHoNO4+OKrdpwP4/7F7YKvb/K/onBTPEpG7AICInCwiP6DZ/S8AbpBuMb+z77BwPoCHGTMlA87yEJHvFZF/0w+6b0I386Jr9xEAPwLgjgBeOtx7QmywcmREKfUH6AZlf4iucb4LnTfkTOX5Xlal1AfRhUG8Ep3n4fMArkU3SAzll9GJx+fRCdmrItLw5RUAHgrgLdqUZLQNSqm/BvBsAG8BcFX/38V56GIQ36+U+vTwB+A5AB4l2ptWDJ6DzgPyRhH5PLoFQA/s970EnffsU+gWdL7Tx+7e9muUUq4fGxkrjzv1267v8/4surqkp60AXICuXr1WRHx+LpgQspvfAvA/0C0kvB7dosgn9vHiA5/u912DLjTup5VSH+73/Sy6cJBPo3uTyJ9p5/0Nulmz/4WuHX8Ju0MxBsfFZ0XkPf3ncwHcFp3WXI9usLy/3/enfZr/iG4Rus2bm8LL0c0MvsYIV/wZAL/Va+Nvohtw+xBi73no1hJ90tDu5wF4oiM0FOhmJK4C8M4+XOW/42gM/7PRhfwdRqfbb/C0G0qpf1JKXe7YPVYe34Dunt2ELiTyb2EsHlXdmrFz0IUovYiDc+JCFF/q0DQi8nXoXgF4mupeZE/I5jjp9JPUBRdfkCWtCx944RVKqTOyJEZWiYg8FN0PrE29Zo8QYoGaHQ+f2BpERB7dT5HeEZ3H9P04uviQEEIIIYSskNkH5n389UdE5CoRyRYrtzLORjd9eg2692g/XnFqgxBSAWo2IYTMx6wD835hxPPRvSnjPgCe0C+cIRpKqZ/SVpqfqZT6SG2bCCHbY6ua3b+zmmEshJDZmdtj/gAAVymlPtYvhHglOu8wIYSQ9qBmE0LIjLhWPJfiZOxenX4QR9+CAQAQkQvQvXUCO1+78x377r4PRw7rv/mwm5195m9E2BlLowQ7+3aS8xzS8LnGIa8c+Zo2LInS9zm0PHzv31QaS8J1vYc+fOiwUurEmc0haUxqNkDdtqVB3faHul0f2/XecOgG3HzDzdl/+4SMM/fA3HaDd8VO9z9XewDoVvU++iGPBgDsP38/Dl10aO/Z/b4pXOfmZLDDzMu1PSRNX3JeZ2jetSh9b/W6F3s/UsrSdn1j7aEFbNd74QMvDP3FQlKfSc0GqNuuNH2hbueHuh2Oeb0HzjtQyZJtM3coy0Hs/kWwU9AtcBzF1Tj07YcuOuSs8HM1hCkbQht5ToGNSatlARmYw8bh3tXo8JYo7mRVRGk2QN3OAXU7LQ/qNlkicw/M3w3gNBG5Z/+rlo9H9wMvTnSvxVQj82mAJRvpkLYtj7kapev6WhaEFNti7mfIObk8JjHX6BJ3239CChGs2QB1O9YGE+p23DnUbbJkZh2Y978q9rPofhXsSgCv7n/pchLXtFRoRZ9j+mz4r3dOofnGNtwYD8/WKHUvXFPhIWmMpWWm0bIHplW7SBgpmg1Qt32hbk9D3S5Pq3ZtjdnfY66Uer1S6luUUt+slPrdkHNTGswcuIQ8ZjottoHYOpiWySEEta9Tv+8+0/I+1zzcP9d9TK0fpal9T0g+UjQboG77QN2eH+p2vXzIOIv95U/fCj7nE6AZO6l/1r8vpfIvxU5gd4eWs3Obqj/mfZ1zWj7XNGsJ6HkhNqjb5VmKnQB1u+Q5MVC322BxA/OQ6cUcK6tjGYspC7E/taG07oHJbZs+DZ0jr6lzasUM1rqnvvm2XOfI/FC3w+1ouQ1Rt+OgbhMfFjcwD10AMmc8lytO0LSBlb8srji/nPWg5v2s6dWgR4XEQN0mU1C315k3CWdxA3NfWquIMZ6X2ON9zq85LWdSWiBzeLD0dKbiEuegtfptkqvMybZorc5Qt91Qt+NsaRnqdhvM/QNDUYyt6ncxNjVWi1iRTxVAm9chJtZzSMcs21beKGDzuJToHEMIKZepKfxcdXnOttFK3SDzQ92mbvtA3Q6zibq9fpofmA8VJOZpd25x963INTqdmBXirk7Bdl7thuyyqTY+5TI2veoS/Tmn+lOpXTfI/FC380DdrgN1u437sFWaDmU5cvgIgHFhaYnW7dNZkq0pLGFAqC/0ck0lxsRG2jqFFj2SZF1Qt8uxJFtToG678yHrp3mP+RIoKSJLEKjajIlVC7GUvuKc+5VdrXqjCGkB6nZdqNvT+Y5tI+uFA/NMTMWZ6ceUzqcGNe2xTSPqU4Z6jKXvPWilfH09JT5xsFz0RMhuqNvU7RKsRbeH2S8yL4scmPs0vjkrs6/ouuLQptLOLTam8On2rWEAleMdtbnKPVeZul7nZmLz8uS+p7Z0OdVKpqBup0Hdnoa6PW0Ldbt9FjMwb+VJOBehi6FMoUjxwIS8JWENhHpdcuddA/Mel7LD9PhQ6InO2vSFuj0f1G3q9lZpevEnsHvhQ2vp6WmlxpCVOCclj5DryVGeJa8v1r4lipStfue8jrHBwdYGDsQNdTv+nJQ8qNvUbVf6ts9T+0gdmh+Yx+KqYK7V0yXyynX8GLmuxfWKrbkoLQix17NEkdepab8++KHgEx+o22npULfTzmuFVnR7Z99ONTu2TNOhLKmVwmdhRWq6oQzTRqnTVK5ri7UtZSqr1YVNA7kWb6VcX6tlE0qJ2FmyLqjbfjZQt8ehbueDur0sVusxB/YKX+2KaQpNqhjr31OvbfBI5bKpJXKUzRrLJYaaXkayDajb/lC33ay1XGKgbi+Lpj3mqZSoXKlP34PIl5iqymFbytRhi4259pTmlhYuEZID6nbYudTtdeRf+5pJO6zaY26SWvFzPMGn2lFKRHNNG7ZC7pjU0tdWKuY0d31p6R6TbUDddkPdnk6vJNRtUoJNDcxjiFm9b6PVhjEmhEud/ipV1i2+gQAYtyvW5lbuJSExULf9aaWtU7f99o3Ryr0kaXBgPkKL4l4qLds15licNDdzdaQl3hKRknduW0p0GoTMAXWbuj2WD3WbtM4mBuY1G6MNm1iWfp+uy46xVf16B9dS+bmYw0Z9sZXPPYvpGEPEVf9eMv51uN4cb6YgxAfqttsO6nZ4HtTtZdQHsvKBeYqwp7yGyrTBlq65LSSfubwgS/K8zFF+usD5EhOf6vKClfa2uOxIHZAQEgJ1Ow3q9t7zqNtwbiPtscqBeQ5Pi/6EnZLGgCkMsU+xqU/yJb1QNZ/GfcslpLz1e1bztVs1p14JmQvqttsO6jZ1u/V8ST5WOTBvjVxTSDlERv8capOP56G2yOf0COTsDFNfiTYnuWJ0CVky1O15oG7ngbq9Hlb1HvPWnhRdwqpTw2YfkR+bEm6tnE1struux6fzLbXyfo7zl5In2S6t1Tfqdh2o22m0fn+JP4vymK9txXHM1Nwceft0TEvAx4Pgcw9S6lZq+dUs/8HztMS2RdqBuj1P3tTt8GOmbKh1fgrU7eWzqIE54PYELE2UfOMga1xLSU/E3NS2NWYRkXl8iQ4+FAo9SYG6XR7qdj6o26QmqwhlMSv/EipizVcX5RDwXG8/mItcU5Nj08lj5RYq0EP9qFlPbKQurCNkgLodBnU7/nzqNnV7SSxyYN5ixQ+lpnfDFI/YsrSd02rjz9UhxSy8iqXVOu4SeVsnSMgAdTsN6nY81G3q9pJYXCjLQIsV35eWphy31BiXXGeWwpbqEwlnyW2Qul2HJdeZpbCl+rQEFjEwrxknNWe+c8akmfnE5NtCDJ0vIdOXoZQqh5anH9lZkimo22XyGvvum0arumJC3c4LdXsZND0wP3L4SG0TZqOFhhxqw5IbeU7bp1bBxw5QWqgTviy5LpC8ULfbtmHJbZW6nZcl14U103yMec6KExuvlvsJ2Od9tLUIuc4lCVDMvY+JFSxxP5e0YKtlbxGZD+r2vFC3d59D3Q6Dut0WzQ/McxLbQEpXWDaIcqR6PELF1fWmidR73Lq4651hyelnsj2o29uDuj0P1O022czAXK+AoU/UJZ8ml175zfe21qaUELa60h6oZ5tvnq2WG2kf6nYZqNv1WYJubyksrSWajjHPSay4m+cTN2spI9d1TC0Wsu3LtcAod+fZQmdMyBTU7fKspYyo22QtrNpjbnpbbNtJPloQjdS4Pp9rqHWdrnrLAQtZE9TteaFul4W6TUJZtcdcr8Q5KjQbxV5aXhAVQ4lXvE2t/vdNA9jtyeF7k8kaoW6Xh7rtnyZ1m8xN8wPz1F/gGv77vP+VlTiOlsqt5ruTbeSK5bSdO+d1tlSmpH2o2+3TUrlRt8vQUpkSf5ofmKc+rbpWW8fmxYpupyWRb4kSXhLWQdI61O1lQN22Q90mNWl+YB6La/FGbCxbjWmopVG7bLYgfFu4RrJdqNvzU7tstqBpW7hGko9FDcxDBCTVw2KLCastYGScFr0brnfE1oJ1mMwNdZuMQd2ehnV4WyxqYB5CbEU2PTZsEMsi12uugDK/XhiaR44FSCH5EVIT6vY2oW6PQ93eFosYmKe+SsmWnv6nbyPLp0URG3vHru0Ym6jnuC7WcTIX1G0SAnXbDev4tljEe8z1xUAlGu9WKr250rzkL+PVJGTh2BS5yqjme3Zd6eZ8LzTfMU1MqNt5oG6HQ932z4O63R6LGJgP2CoQK9Y4Ztlsbbo35QcratWr0vmWuO9DR8i2SEyo2+FQt6nbJtTt7bCoUBZz+tL1nez9YQOWTzhzxD3O4UmMySPEU8S6RWxQt8OhbqdD3fY7hnWrXRYxMNeZqkhLrWy5f2BhiWVQC9uAAUgXX59YQ9t9Sq3DeufuymMM3+PpZSG+ULf9WGIZ1IK6PW3TWD6kXZoOZdnZt7Nqocr162K5WFvsYqrHIVdZxKaT0tH4invKNa6prpB8ULfnhbpN3Y7Jg7RL0wPzNdOimLZmTwq+orimawbcnULodY6V39rKjBBfqNtloW5Tt8kCQ1nWAhtJWVqa1pvrdVmm1yh2SjTFDtcULyFrgHW5LNRt6jahx5yskFZi6HIIW+ov0OX0uKTkuf/8/cCB5KQJISuFuh1vQynd3tm3k5wuCYcec7IaYqcBS3YIKb8AF3JOzs6EEELmgrqdBnV7fXBgTlaF+StsPuR4TdVU2jU8IKEs9c0YhJBlQ92Oh7q9Pjgw1+CT526WUh4pryyziXsuL8bwF5Nei0Lbok2ELEWn5mIp5UHdnocWbSLjcGCusf/8/dnfS0vKkvte5XzV1thrs1LtLvkDJBRysiSo28uDuk3dJm44MDfgtNBuWu/sUu+V3qHPea1Tg4ma5d76PSfEhLq9m9bbMHU7P63fc+IP38pCnJgitOaOr6SoucotZaq01PuUW/vxFEJIGNTtPFC3SS3oMSej5H6nak4hzTV9XUrIYj06ruNL/LqdKx+KOyHLhbodD3Wb1CZ6YC4idxORt4rIlSLyQRF5Sr/9ziLyJhH5aP//+H67iMhzReQqEXmfiNw/10WQsgwemBRKxIDm+Fni0qKZSyznElwK+3qhZm8L6nY81G1SkxSP+S0AnqqUOh3AgwA8WUTuA+BpAN6slDoNwJv77wBwFoDT+r8LALwgIe/VkSKAY+flENUc4j5QYsFLTJrmlOLcwsafTiYVoGZnhrodB3WbEDfRA3Ol1CGl1Hv6z58HcCWAkwGcDeDi/rCLATym/3w2gJeojncCOE5EWJN7XELlK65TC1JSRHoJ3oNYkS/JmE1j+8wOkYJPckDNzg91e/60qdtk7WRZ/Cki9wDw7QDeBeCuSqlDQNcRiMhd+sNOBnC1dtrBftuuViYiF6DzzuDYbzg2h3mLZ/AS6It5xsS8JUEwbS+VRwwtlBNX0pMa5NTsPj3qtgF1ezqPGFooJ+o2KUnywFxEvg7AXwD4BaXUTSLiPNSyTe3ZoNQBAAcA4KTTT9qzP4TWxM4HlyCG/DKabX/ucghZ8V96JXpMB9JKvQixnSvvSQ5yazZA3aZux9lA3SZkL0kDcxHZQSfwL1NKXdpv/oyI7O89L/sBXNtvPwjgbtrppwC4JjRPszGspaKHCHhqPjneIWt+dk3nznF/YsqrpVeJmSKvf9fLsAVbybKpodkAdTtHPtRt6jbZBtEDc+ncLBcBuFIp9UfarssAnAfgmf3/12rbf1ZEXgnggQBuHKZPffFpyCGNfY5G7hK5WlNhOWMOTaGfW4ByvXVgSKsmtfMfWKK3kvhRQ7MB6nYOqNu7oW7vhrq9LlLeyvJgAD8B4GEi8t7+7wfRifvDReSjAB7efweA1wP4GICrAPwpgJ/xyWRYAOOzQMZ2jG+nUOK1UC01lhLXN5VfqRjFOdKcu7xManlbWqmvpAizaDZA3c4FdTssTeo2WQPRHnOl1N/BHoMIAGdajlcAnhySx5HDRyIs24truikXupC3uihkbQ13KR1HCmu7Z6Quc2g2QN3Oydo0gLpNyDRZ3sqyJHI14rFFPnPkH0qNuMESebYmwkulJa8gIVNQt8vnUzJP6nYeqNvbICWUpSlyNPyQCh9ybKkp15rETj8vjTVeEyGtQN2eF+o2Ie2zmoF5Dkp4T0oLxJZ2Pl2/AAAgAElEQVSenmt6rWoLfe78t1RvyLqhbrcNdTsfW6o3W4YD857QX+qq3eAHQjubXHb7lBVFJB8sS0L2Qt0Og7o9LyxLEgMH5j25X+lVAjbyerDsCWkP6jYZg2VPlsjmFn+a+Dbckq+RMnG9KaB2B1OTJb8/mBCSF+r2MqBuExLOpj3mIa/KmkvczTznPC8XodPLhBDiC3W7DNRtQtpg0wPzlhYDjeU1CGZrwtmSLYSQbUDdTqMlWwghe9l8KAvQllC5bJnqYHTP0NSxw/5cP/OcM71UGwgh24C6HQ91m5B2adpjvrNvZ/Y8TZFs6ZVLNgE1t4WK7NzvES6RR+37Qwg5CnUbu2yhbtupfX8IaZWmPea5fto5hlpeBPOXvaY8Kr7eFldeS8K8J3PYb5ueTkmrJS8fISWgblO3dajbhITRtMd8TuYUu9BGPhWjqO8LvY4lC47L9pLXlPJrgEsua0JahLq9PKjbhIyzqIF56UU0tunQkvnE7ncdH1o+OctybvEaE9mYTi7G/qV5rgipAXV7/Hjq9tF9IVC3yVpZxMDcbIC5xGTstVs1G2/oopgWVv0v5ZVkNnw6yBbKmJAlQd2ePr62plC3CWmPpmPMgbTpPt90a68OT8m/Rvyejx2lyH19Ywuz6FkhJA7q9vS5OtTtMKjbZM0swmMO5F2F7rNKPiX9lPOnbAyJWzTTNMkpYHOJYc53A/t0jqW8frlp4S0UhJhQt6nbAHXbBXWb2GjaY26+divVCzO2ulqfFgv9RbmYJ3WbEIcQIvJj6DYP1xQjYnMKX6n37/qm15rI134fsQ/sfLYDddv//Ni0qNtHoW6Xo+YblrZM0x7zI4eP7BHelA7epwHEro4PFcaUVeKhTOUx5fHxSb/GwCs1T59yqQU9KWSpULfzQN2OO5+6TZZO0wNzYG9MIbDXWxCajkmKxyTUjlzn5cqjhalgX3KLXsueipgp39Y6hDkHMaQtqNvpULftULfLQt2uz6JCWYD601a6h8UW61Z7QVIMsdOgc1Hy9WetXHcuW1q5pqW1AZIP6vY8tNLWXVC3508nhx2kPk0PzH1JiWHcSkV0lVEuQZhKIyWfmFjQkHTXRAvXtJU2RdKgbk9D3XanuyZauKattKkl0HwoCxA+ZRlSyXM2CNvCoilqTBvpnqO5BCE1nxLlswYhMq+h5hQkpz+JDnU7L9TtcmnODXWbjLGIgXlJEcrdeaR4CWqIfI58S05ZbkE0Yq/PdZ5eB1tZqEa2B3U7P9TtdqBuk1I0PTDXV/eXIGb6dOycpYpRC9NoNmwLyHKk5ZPmEu+jyRyetak632rdIuWgbs9Dq22Lup0GdZs0PTDf2bdTdDFQKHNU1jUIS+4V+DkWZrVerjGvbAvZn/v6fYSd4r5NqNvLhLodDnWblKDpgbmNkIrqW8l9KqLP+2JzVeileHBc5Zu7YdcoizlfJWb7bDsutF6MHR9bx8bOG0Sdwk5MqNvtQN1Oh7pNStL8W1lslSlGSMYa0lRFZ4V1M5TdoYvyv+4pVdhNj02ptwTkYo7pS59tuVniq+hIGtTttqFu52PNuo1Li2dDLDTtMU/9OVizwfh4T2yYT5tzCUSrQmSjNXHXac0b0MJ9zeVtDEmPbAPqdv327Qt1258W7it1exs07zF3EVqhUt/FmiOtHHm3gu5paUGwBgZvy5g41fLgxg4SStjqm+ZUWa19IRbJC3W7LtTtcKjbZG6a9piPMaf3IySvXKLcsrgDRxt+qfsQc/0tewpCyykmNjEHPgOamHJuvT6TeaBu14W6HcbWddv2K76kPJvxmOvELp5oXXRrYivTHK/KWvLTe23vSihjrzmb+j5g83ot/T6SfFC324K6vRfq9tH/Ry5NC0sjcSx2YG5bHOJLrHDUFPmUNxOUxFaWOYXZJ63aZaDTglCH4ivevmnZ2sn+8/cDB6KTJSuBut1RW7Oo27uhbtt1mx7zOiwylCV1UUhMI5x7IUro65jmZqzDyWWPz7W1JO5z4LrenLG4qWztnhA/qNu7j6NubwfqNglhUQNzXWRTRKR0BczVGdQSbx9c9yH3qvwpWiqfOWyxebn0z3rdm4opjK2jLZU5aR/qdjtQt/dC3SatsdhQllj4yqyyzPnUrXcyNZ/2a95rn7xTy2fs/NplT7YBdbss1O328qZub5fFDMxzxVPl8NyQ6WnaGkI/N7Xr0NR1D/bZjvO5R+b5sa8NI9uFut0W1O36dYi6TaZYVChLLUIFpIVXJZXA51VQua8993ReLlp43VjqlHnoubXLnJAQqNsd1O2jULfJEljEwLx2xWr5KbN22bgoEbNoxuLVomR9GPOW+JxnfvZNx+wszM9Le4UYqc+a22kqtcvGBXU7LW3qNslB86EsuRvz2ipfrCCk4rOKP+fUaCvTnnNP9bru79j0ZO43Xwx2uO63fp+HYxjDuG2o2+NQt8tC3aZuL5lFeMxLETr9lCPdUpW+lY5r7Y16rmlEXdhLvRprSrhNW8bS0b8TUhLqdn7W3m6p2/Z09O+kHZr3mOfG9Bj4eBBiKq35JOrbmFIo6SUw0x57Cl86JeItl5DmWF6Af/2ao66TbUHdzpM2ddsf6jZw5DB/+bMGzQ/Mc06tmGnN0UDMp9I5BbDWdOmSyHU/XPc2dZFP6XsXUidZj4gv1O14qNvTULep22um+VCWnOKu/98Sqde8xTIbwzZNGVtGsdPnc71FwQdXfCPZLtTtdKjbeaFu74a63S7ND8xzkho3lsuGlMqfMj07B0tp2FPT02MxeD5xhCGr6mP26286KF2nhrLSy2xswMQwFpIT6nZ5qNt7j5+yI2Y/dZv4sJmB+VCJWxD50KmuqWncqWsr8ZQ+V14l8LFxKjbTFDwfQuteyyvk9anUVm0ky4e6nQ/qNnWbur0Mmo4x39m3UyTd2JjBXA3ON/+pVdM2e/QFH2MinDoFZrPNXDDVIqYojXk2fLGlleP6YzqEmPP0/GLtbvV+k/mhblO3c0PdHs+Pur0umh6Yl1gRXCKmrJQNU0I81gGY3pqpzsIHlwfIJuytPo2HTFXayjR0wU1ofUv1DM5xD2La0P7z9wMHChhDmoO6Td3ODXU7nZg2VOohm4zTfChL7qnI2oTaUOMJ3sXUNGAL5RtCSJzgWBkOZWI7Rp+unhLe2p0ivS4kF9Rt6nYpqNu7oW6vj6Y95rmZ+2k4B7UbvYtBsFq1z4XvVLQtNtHnWodzh+PNdEIEfwzfzieEklPZ7ARILNTtfFC3x8+lbu+G7zGvQ/Me89zUFvfY+LMc6J6TlKlhU6RS0itFik02UfYRv7E4SNfnlOnp2AELIUuDuk3d9jlXh7pNlsqmPOZbJ/UJP1eauZmK0zO9IKEehuF8X5GfOjaX12qqo4jx7JTyvjDGnJA4qNvU7SlK6TZjzOvQtMd8Z9/OKqcj9afmrZFLPGzpjKVtCmqsByyHKOcqg9x1qEUPGlke1O31Qd2mbpP5aHpgnpuxJ2GfbbkJWWmu/18yuQTSjAk0/4ZjzXP1/6Wmp/X0pxYMuToq37zMPEI7PptNerqtLngi24G6XR/qNnblTd0mJdnUwNx3EUjOypvT0zB3o3JNL8akkQuzDFy2ueINS9ll4hOPOCXGc5V1TLzpGgYbZBlQt8PzBKjbMVC3SQtsamA+YHtKr/EkWVssffMKbfi+i2RSrsfW4cV4OUoyJpwur5B5DT42jx0TMmBxea1ytA12CCQV6nZYXtTtOKjbR+FbWeqwucWfY5W1hMhP5WfzEAD5YuPmJmbBSsg55vG28/Wy841frDlwnPIYhaZhbp9Kw9cj6Zun8zwu/iSRULfLQt0OZwu6vfN2Lv6swSY95rVxxR6agtSSuJfyUIWmFyJGIaId2jGVwNfLEhJzGWPrlLcmZPo0ZrqVkBahbu9NN+fx1O2jtKLb9JjXYXMe81q4psDMbS2J+kBpQYtZzKN3OGanOPZ9Ki3drjGb56BmXXDVxZjrvzUdeszJwqBu74W6Pc6adJse8zoke8xF5DYi8g8i8rr++z1F5F0i8lEReZWI3Lbffrv++1X9/nuk5l2bEo19eFJtSehzezz19GI9A2NTgHq6PmUZOg1Z896YndDYvUmxc2yamSybLWs2QN3OkR51OwzqNgkhh8f8KQCuBHCn/vvvA3iWUuqVIvLHAM4H8IL+//VKqVNF5PH9cT+aIf9qzB3bOBelvAtTcZmubeY+12dXzN8g8iFCONWJTHl0cuAS1TnyNkmplwxjaY7NajZA3U5Nl7o9zpp0m6EsdUjymIvIKQAeCeCF/XcB8DAAl/SHXAzgMf3ns/vv6Pef2R9PeloQ95IMMXB6LJxratiGjwfIVYZjHpgYz8WcAuuKba3BUI5T8YwDvseReaBm52ft9Zu6HQd1m8SSGsrybAD/EcC/9t9PAHCDUuqW/vtBACf3n08GcDUA9Ptv7I/fhYhcICKXi8jlN99wc6J57bCUSj4mIqnX4LsQMHZRjB5z6DpuKffBxGcKdA5cXqqY80gVsms2QN2uDXW7TZau2zv7GGNeg+iBuYg8CsC1Sqkr9M2WQ5XHvqMblDqglDpDKXXGHY67Q6x5JALf1eWxmIt1Ys93Ya58H4u3C41fXCI+HWrqdfp0qqQNSmk2QN2uCXV777FLhrpNUjzmDwbwQyLycQCvRDcd+mwAx4nIELt+CoBr+s8HAdwNAPr9xwL4XGimtZ88U4j1KMyBb5mWttOMz9Pz87XRdwrRTNusW0utZyGUvJ81YirJKFU0G6Bul4K6Td0ukfYWyrBlohd/KqV+DcCvAYCIPBTALyulnigirwHwWHTCfx6A1/anXNZ/f0e//y1KKav3ZeDI4SOTq7gHagukL0OljxGuEtTIe2rFvR4PB+wWCp/V+uY5ocfVFqXQmMpWWZq9a2cOzQao23NA3aZul2Jp9q6REj8w9KsAfklErkIXj3hRv/0iACf0238JwNNyZuobB9cCttXac7OUsvLBdi1T19bytU8tauL0I8lMFc0GqNuhLKWsfKBuE2Iny8BcKfU2pdSj+s8fU0o9QCl1qlLqR5RSX+63f6n/fmq//2NT6e7s29mzInwNFbu2pyUm/1zlrq/q94mhc63+z1WGLQu9jdjBQY3rHHuDA6lLKc0GqNul8qZu781nKVC3SQiL/OXP0IqjH9dCB1HLhtgpQnNbqv0hMWy2adGBGOFYQixiyAKn1GvI3TbMtkmRJwPU7Tio29TtsbRK6vaRS/ke8xqUCGWZDd9XAK3Ja1Mal6djYK4ydK3MXzJmPXT9z+FZ0svPrPvm/czRNkIWbJFtQ93OD3W7HFvWbb4usQ6LHpgP2FaCxzCXiNQSK98OUd9Xq2NcUyyljTGRzyW4rjRK3M8W4m/JsqBu+0HdbgfqNpmDxQ/Mc4n7nLQ0Jeqz2Kak0PpM/80l9KXvS67rCIlPXFK7INuBuu0PdTvNllSo22RuFj8wJ9OkiGQpgdCn/1oRoSV5eqbKrJUyJYTEQd32g7pN1sYiF3+WQl+ospYGEipaodOmsbSwmCdkMVOOvEqk6bJ/6Dhb6kAJKQF1e1u6PSfUbVIDeswNci7kqM3S7dfJLVI1BuVz5jenqK+pnpFlQt0muaFuk1pwYG5BbyBLrbw133k7RerrspbGnN68Got5hocc/Y+QuaFul2WpZRoLdZvUggPzEZY4lRTbwGqKe623B6yBVl4nZ3tNGsWe1KCF9hAKdXtbULfJGIwx9yD0RxVqkCMmsTRjcXUmc4aa6Hm67MmRbklaEXkKOmkF6nYeQnR7TVC3SS3oMZ8g5Ml2KU+bLYjBFHOWYekp8NbrQ074XlzSAtRtkkrr9SEn1O224MDck9Cpp7nF3se2EtNnaxCv4RpKroRfQzn5Yr6jekvXTtqCum2HbdKPLZUTdbsdGMoSSMy0T4kBX4gNJZ9+zfKw5VXz1U++90u3kVN76dhiFwmpBXV7b9ot6zapg3m/jxw+UsmSbUOPeQShHoyaA7254uSGP93bNNc7b8eucc3lvhRYFqQFqNt786ip2y1CrToKy6Ie9Jgn4Dvlk0vcWvG2uLAJuu51KSXypucnNp9cXvOtCZqPZ23/+fuBAzMZRMgI1O3d1NLt1qBu72Vn385M1hAdDswr4jOgDBWLqeNDpyd9j5+rk0u1w+f8VIF2pbHGqeFc9YOQpUDdXifU7b0wlKUOHJhnICSO2XX81CupvDySHuh5u86x2ZYjXGTuqeFU73loXj7b1sYWrpGsE+p2m7o9J9Rt0iKrizGvGROoh0GMNfhBMF2vKHI9tU/l64NP+ZRa8DT3vYnJb+weTJ2zFabeXLHWTpyUg7o9zpZ0m5SBur0cFu0xb3GKKTR+UT8+JL4v9LrNdFortyXTYj0shctjqF//VsqCxNFie6Fub5MW62IJqNvLYtEDc8A+zVi7sYU8kfpMhYaK81jnYq7AjymnGuUbOp265unXWrA8SS6o2/756+ktTbdJfajby2PxA/OBliqfTwxgaHyjLz4en9RBa+03B0wRa59tirqlelWLqTJgZ09iaal9Ube3yVofWKjby6X5GHM9LkqvaCkCVbszmIppjCHkmmKu3xWbZnZaLTX2HHHivmnUrlOlCBmEMHaRDFC3/aBut0HtupUb6vayWZzH3LeyuJ6CbZVxLlEKWRwUk/YcDWms3FpryHOW7Ro7ttBBw9jDTWt1g8wLddudNnWb5CSnbuPSXFaREBY3MAfyTCe2MJAai0d04fuqrBgbQtMqKei5OqyQujJWJ8bSaaEujZESl+qDWXecb7bgDwxtGuq2O60UG0LT4kDcTmshLS3oNn9gqA5Nh7K4Xm6/NmEZGkULC59aKdtYO8xV5jmvx3WPWikzF75vnIilhheTtAt1e14b9P9kPVC3t0vTA/NStFgJXbGAtuNCyP2u3FYZewexvi/HvXe9x7hl5rLRtx4TEkqL7Yy6TUpC3d4mqx+YtzyIyultcb2P1Kex1S6TOfPPGRta21sWSujr2vRrHHu4MeveksqEtAl1u33dXiMtDk6p29tjkTHmPjjjXDORMx5Nn7KKiRk00xk+2xqsnrarU4jJO5ac4SatxQgugdDXsOkiz5hWkhvq9jJ0m9RlLt12haWRsqzWY14yfiqH+E29ZaCEwLYYU5Y7Bnxq4djoCvQNogu2TyevT3lOeVl8jyNkgLrtztPMi2yXuXSbiz/rsNqBuUmJhYCxmE+7tvgun4Yztc+cwvL1rizZ6+JzX5a4gLMGYzH7NoZ6XPptAmQ7ULe3odstsfRypG4vn1UNzKcW+ZWIC4zFJsT6tqkpV9dijTGPsE8nMpZvy3GfJJ2pBUC6gE91Xr7HEULdpm6TeKjb62PxA/PSMYlzEhIvGHKNQ0Pz9SYvzSOT6363en058BFk27axgUSpN1KQ9UPdnmbtuk2moW5vk1Us/pxD0EsuLNSnRn3zCD0nxHZT5G2eoZxlkXu6OjbPpQ4MXIxNTdrqjWsxm56WyVxlduTwEXYUK4O6vWzdJmVYi25Ts+NZtMc8tnL5xvrVwNdLMjb160rX9nnqWPMNAGas2loW9q3hGlyYMbG271OCnuLhXHPZknCo29RtMg11e7ssemAey9giCN84vRL5h6A3Ut+8zMZtbjcbu88iyZB4NNuTfmhHVQIK0HR9Gttee1BEtgF1u45uk3ahbq+TxQ7MSwym5hyg+XgwcnUEU6uxzbx8YyRDG3aJKVUOqvNg69RdXjV9ypyeNxICdduPlnSbtAt1e50sMsbcFJklV7AcgunjMTLjDWOmuMx8cop9iWufQl8Is+Q6lIOQjt4nTrbGugHSNtTt3axBt0ldqNvrZJEec5fHYC34eGB8F3no+3wa3FRMY45G29L08po9Bz7XFVKPbGVlzoLYwpMo9ASgbtu2LUm3yTxQt8kiPeZrxzbdZD7J6sf4xiv6xB7azjOPsz1x+z5pu2yIeVJfY+eeg5DyjFkf4OOh8z2WkLWwZt0m5WlSt9/ulQ3JzCI95ltlLHZsiqnFHuZT83BOjI0hhHiFUvLZErZYQtv9DWHKI2fmTwjpWKNuk/xQt8nAogfmsYO6VvFtGHpD1WOkx0R8OMfc75rm0ss25Ol8sMVXTGJFh6vK3UyJcY7wHZY9iYW63aZuk7pQt8nAogfmOmupcGONzyXOIfFmLqbCD1x22TrZXB2vrwit5d7nwLcsYuuIeZ9Z9iSFtdSfNeg2qQd1m+gsdmBuW6AwVdlcFdLcNhxn/oWkmUrI1KfufdG32Y4bi/F2nTu2GMT0ttgWkZh5TRGy+MR1TWugpICmeF/GPHU50ifrhbq914bWdJukQd0mqTS/+NM1zTYmuFMDQv2YkAY0dqweG5aDkAGnPh2ql5d5nSkemrFjY+LKx6ZPQ0VtzYNzPTxIx+d+2AYuY/XB1x5XPq6BAdke1O1plqbbxB/qNkmh+YF5qSfPXIPCEsTY4BJ5fb+evs2rPdbozTR9BKKFslwiPkI5do9dg43cU+dj6VLctw11248WdZvEQd0muVhsKEsJcoj71BSqD7HTSVPTxGPeJ3Na0yTGi22zrQRrEhPXfXAJq62j9imPlDo2/Pl4IlsYMJF1Q93enY4O2988ULdJTpr3mOfAp5LFTpOG5lt6EDkm8mOx3/qxoWLhysN3mji13Nf4pJ96PSXLYyhvn4GAa0qXkCmo2/PrNkkri7XpNqkDPeYGpStjqlcmFL0hhsYbjg2ubd9tnhxXGrb9saypU3EJZ82wId3bMny3eRhdXqM13R/SJtTt3dtidZvEQd0mOVntwNwcHPqeMwc+9uQerI5NgfkKypTHxnW8a/DNp/JxfD2GIcebxNZFU/QJyQF1u75ukzSo2ySVVQ7MzUrfmvhM2WPbn7PxuuLIzOOn4uOGY2wxkD6emqm8fFizwMxxbb7l3lobIuuDut3Rkm5vjRzXT90mqaw+xrzFaTrfmOuQc0LzH4R4LDbRVnY+3huf40ldQtuDWVds97jFtkaWSYt1aUu6TdqEur0NVukx12nR81ILM77M5qGy7dOnv2zeGrNhj3VGplcmtuOq5d3Ry6mFehU782CK8hRjA4EpDx4hobTSvlqgBd3eAnNeP3WbjLF6j3mLmI1lrjx99o15vk3xt4WnhCzqTOl4a3nlWxevEPt8jh3zpuj7KPJk7WxJt7n+Z16o20Rn9R7zkqRU4lqiNzWlOTRWW0yiyVTDnmPRUUsdZY10SmLrsG0zBhRzsiSo2/G6vVZCPNCtQ91ePhyYJ5DSSH0aRc7p3JDQD70x+xybi9C0liwssbaX6BhyDFSWfC/ItqBuh3ldyVGo22QOODBPRI/jG77bjklJOxdjoq17W2KmbG0xjmOkxJrXXqySKwQn5typ2ETf9M04VR1XfXbdIx+RX4KniWwH6vbR8/X/W8T3YSuWJev2kcNHvOwieRGlVG0bnJx44onqnHPOqW3GHmIXXtQmZNAc4nWZimObGrD5DrR90wzFx0aXPSHXV4rQePuQTsaWnkvQze1mPgcOHLhCKXWGd+YaObUgxQ4yDXU7L7V1e+2DdpfGUbeBSy+9FNddd514Z6xBzY6HHvMIQldQt0Jr9ujeGp9FKPpxOTuLktOTpco8dTp+zC7dm+jqvG3iPraoiJDaULeXmX8pxnSQuk1qkjQwF5HjROQSEfmwiFwpIt8pIncWkTeJyEf7/8f3x4qIPFdErhKR94nI/fNcQh2WWoHNhmkLKUm9tpBpNBc+NuS4Bz5iNCXSsdOSsYw90OQqk4HQ+Na1duJrYcuaDVC3ffOw5bc2qNvU7VZJ9Zg/B8AblFL3BnBfAFcCeBqANyulTgPw5v47AJwF4LT+7wIAL0jMuzprEXnXNh9MMdDTGYuLNPNOyd9cbR6bTgq66OboJH3J5e03Y259H4z0cp/zukk0m9ZsgLoNxOn2GqFuU7dbJHpgLiJ3AvAQABcBgFLqK0qpGwCcDeDi/rCLATym/3w2gJeojncCOE5EFv+YtvQKbTbS0Kdn22ITc4HKmJjkflLPMUA3PRuhC3N0kS9ZP3RvR2hZmvdtrKP2ScMUenpg2oOafRTqdpxurxHqNnW7NVI85t8E4DoAfyYi/yAiLxSROwK4q1LqEAD0/+/SH38ygKu18w/223YhIheIyOUicvmXvvSlBPPmYU0VWW+YU9eVO7wjNfbO9n8qvzEbQ+JRQxba5MC1eCe0DGPKyzZD4ZpCXVPbWAlFNBugbtdkTt1eS7lRt6nbLZMyMD8GwP0BvEAp9e0AvoijU6A2bCt797wSRil1QCl1hlLqjNvf/vYJ5s3D0jwvY8Jj8xCnTJ/qAjJ2Tq4wEl/bcsVi+nhXcnhfbOfbyjbV42HzmI3NeLjSML0wpBmKaDZA3S5NS7q91DZN3aZuL4WUgflBAAeVUu/qv1+CTvQ/M0x39v+v1Y6/m3b+KQCuScg/CXMaKVfDWAOhXhXz+JB4N9f5MWIVm9/Usa4YTNs0Yu7pwFivSighdvt2WhT65li0ZgPU7THm1O0lhj1Qt8N1e2ffTpJ9JI7ogblS6tMArhaRe/WbzgTwIQCXATiv33YegNf2ny8DcG6/0v9BAG4cpk/nxiVI+udQ4V+SyE8NlvT9ZuyiLt5jDV0XjJiwkFiPha/w2PJ2HRsiZsNf6tSgjwcsB1PhO7H5cmq0PZas2QB1u1XdbhmfhxXqtv086nY9Ut/K8nMAXiYi7wNwPwC/B+CZAB4uIh8F8PD+OwC8HsDHAFwF4E8B/IxvJnNVkFhxDz02B2NTVmOYT9iuzs0VizaI11h8npmfa5/ruBwdQmgaPp7dqYGBzZtnnjN3ZzcVM6qj3/fUznkNnfpKmUWzAeq2jTXqdusDOOp2fL6kDseknKyUei8A268xnWk5VgF4ckw+5hNojtivUI/KWL5zCpPLhtCnZttTvW+5jInGmLfAln6uUIeQgRKhOnkAACAASURBVP2YeMU8QJhhLiGD+9TjfLF1zMO2HO3Ldp6ZF6nPXJoNULdt9vhs35JuzwV1286Ubh85fCTCapJK07/8ubNvx1rRcw3ifI/1mQ5rkamnbPPYHJ2nKy/bYLiUuPt6cVLydeVRo6Py8abYyt/0EIV0VK68bJ/JtqBup7FV3Z4D6rY9L9tnUo8kj/ncpE7TLD1/X0Kf7s2pz1xib3sS16fcTI+Pj5drzPaQadXQaWRXZ+Tj3fLxHJtCHMKUMPt41nRbx67X975Ytx8YPZWslNq6WTt/X7ag23NC3d57nI8tOlz8WYemPeY2cjT82uLhE27hooQnyJZeiTIai4EzRb/W1KEtzZCp4rGp0rEytXWyU/dgbGraduyUuNs+2zppelVIKNRt6nZpqNtum6nby0K6MMI2Oen0k9QFF1+wa1tO4YmpqGPegTGPQaxohNg0la9vIy/J2NRciE226wr1Cpjl0UL5uHB1wq4OM3R60+f6UuvrhQ+88AqllC2+eZITTzxRnXPOOUn5Dxw4cCDaDjINdTvMpql8W9ClUrodmv9wXmvl42Lpun3gvAO45sprbL9nMAk1O57FecxbfurzFXfXtlhGwwcMbE/Tc4pXyDShT1rmNG4IPjakeMlyYrPBFRsYMm2rH2P+mdtbKAeyTKjbfum7tq9Vt800bYzpGXWbur1GFhVj3grmtF3rFX9KBOdsvK6YvDEviO2zeU0pnYQ5zWfrhFrAxzPkO3Xrul7XVLV5rIvWyoyQAep2PNTteKjbJJTFecxzEiMKLoFssWLbrq92ZzTW2bj26cJj8yq4RN/0IvjYpAt+bFmVrAuuWMaQOjhcm5mGqwMJ8T7p6daua2SdULfnh7qdBnWbhLA4j3nLlaamyI8JkrkvJuwjN2NiojM23WsTuakyKOV5yjHFG5KPef0+ZelD6Dmm99H83qoXi8wLddudN3Wbuq3Tkm6TOixqYM7KkpfaQm/zmNiOMY9PSdu23XZMq+hl4BL1sWnnXNjCAfR8TY9Yy2VKysJ7nxfq9t48bN9bYqm6zR8YqsOiBua+TMUQjj2Fb4lWyiHE82E7LlTQcgv4WH3L4dWZwtVBjU0/x+DjVbEdc2vefI85GYG67Ucr5UDdTmMJur3zdr7HvAaLiTHPMWVlxr+tjbF4Nde1jz3Fz8HgTU3Nv5b9U3WqVGcytS03uudbrzN6nXMJPz3m24W6PQ11e36o2366TY95HRbhMQ9tJPq0jHn+GoXdZOoaW/C2mKQI4ZyDvjFvlc8Ub2reel4hpMYg+tYZDsDJAHU7DOp2OajbefMhZVnEwDwUsyKWFrMWBdOGj/jM4T1YKmPlN8eimRr3wbdel4yLJNuAum2Hup0Gddv/ONaZNlhMKItJrakhG0sQ94ExWynu0yytjGxT3lP1NcXD4zpnSW2ElIO6HQd1O42llVErur2zjzHmNVicx1x/+o2ZflqKl6QkPlOmZC9jU+0tMyXyY9dhtrcxxtJZSlmRMlC306Fux0Hdjq83jDGvQ9Me8yOHj1jjDG3C7rsYZuviPgYX6e3GtfBML6fW65PP/TQXA9n2j6XlU29aLyeSD+r2vFC3d0PdPrp/LC2fekOPeR0W5TE3Vxbr/0k8FPVpXCLXctnpNsfGHOrb9Ws1B16ucmD7JNTtMrSsPa1A3U7U7bd7ZU8y07THnJQnVqDW3rGa4rhkr1TsAh+boNtevWXzQo15cgghaVC37VC3qdtrYHED86nV1cSPVMFae3nr5bPEa7V5isxryn3/XZ4ZQqjbeaBuj0PdHoe6vQwWMTDXn+Km4qmIm6FR52h8ay7vJVybz4IeHzH3qQ++izr56i2iQ93OA3XbjyVcG3Wb+ND0wHxn307Q64KInalGXGr6qqX7lWtg0Mo15RTPXGXjWthHtgV1Ow/Uber2GNTt9bKoxZ/EH99V3SXT12Paaj+F+06lT73KrfZ16MS8di42H8C/rCjshMRB3d5ri8926rY9H4C6vUQ4MF8ZvsJuemNKNcqQleW1ySmYc4vvQGieIe+7NYk9jxCyG+p2PNRt6vbaaDqUhYQx1bjHpj5jhCGXXa1g8yzEitcc12zmEXsPfeqFy8PiWz5LqQOEzA11Ow3qtjuPVN3mDwzVgR7zhRPaqPWn5RgR0r02qWmEHA/M31EsoWOKtdFHmF3iHtPp7T9/P3Ag+DRCVgl1uxzU7Xy6zR8YqgM95hHkWiGfw46Sx7vS8HnbQmk71k704NfjmJA3A6TCe01agbpN3S4NdZvkgB7zCGrHZqU0GlucYmz+Ng+K7iUx9/t4XFzvcQ09Z+mUupap2NHci4JqtxVCBmrXRer29DlLh7pNckCP+cJoScRstuidR2gsXahQ6+LRUrnUImXBlqv8zDJmORMSTkvthrrdFtRtYkKPeWFsK6B9hc6V1tg5UzGEORpoK8JqdiRzrahvGd/rNzsDfYrbR+wJWTPU7XJQt/dC3SY69JgXQH9C1YXX98k15Fg9D/38kpjTnDkbvk8nFnMuOUqKh2YqXd4DslSo23nSNqFu54G6vR04MC9IrkU7toZT2/sxNvUZg8+iJIpHOrHvvPXdliMWlpCaULf9oW7PA3V7WzCUJQFfL8oc1GpQucR9aiqX5MFHsG0/XjF2f2wLvjh9SlqFuk3dXhrU7W1Bj3kktYVozBMTQ+kG6eNRKV2mWxad0Fel2Tx7PjGMej612wghJrXrJHU7rw1rh7q9TTgwH6HkQpxc5BCtLQlfjWutVb6D2Jrekyl8vCemt8XMp1Q8JCFTULfXB3Wbur0lODA3MJ8yp+IEWyDkidp1buz5ZJqa09X6nyu+1Xbfx2Jkhz/bVKgrPUJKQt0muaFuk1psfmDu2/jGBL8Uvk/Juk2uc3KvjLdNg42t9PdJfyqNFPsGG1rymtXAVk9CysRXyCnypCTU7aPpxNpG3V4O1O1tsaqBeWjjtYn2HAtZfCu/PtWUOpU1Nr0bEyvoWs2d6qUqMZ3G6bmj2MphqmymOnDdE8NOlIRC3aZuu9KkbndQt7fFqt7KMlSyEIGzHaOLvE3wp6aWbHnp6cQKqXldMQ0ppAPzWQySszG7pthS0xqYo/OeC72uu+qg7XpDO7vUY9dS3qQc1O1pqNvr0JGl6faRw0e80yH5WJXHfCBUwGznu6bUbOI+tc3XprHGFPpEmzq16CPuZl62Y11COzdLEPaQcnENFsa8YHr9nZriH7unUza5vhMyBnWbuu1jR2tQt0lumvaYhzytmU+f+pOn2XB8Ku6Yl8RVkUM8BzYBDW0QU52Ir+COpelDqN1zN3x9ui6k43KllXJ+aF4ur5+rvpuYHhizvk3VoVJeGrJeqNvjULf9oG5Tt7dK8x7zoWFO/emeBtsUkc07MjU9NHbu2PG+1xV6zhSpnpaYOLYStpRArxu5xd3nWn3Lw6xreh00PYC6d8S3jpqY7SgmDds59LpsG+q2P9RtN9RtO9TtddP8wNyHsek4Hw/LlMiXqqypjSgHqVOnLlps4CmdlW2KfMD3WkPzd4n2mCdR/z9Vt31tDO0wCPGBuh0PddvvPOq2Hy3e9y2z6IG5GU81Vrl0T0eMhyT3sS3gKgezPMfErXX06dA58horH9NDaNtv++5zzlS+Ptg6rrGp1tT8yDahbqdB3c6fF3V7ee1gzSxmYO6anvSZotQrqctDYyOmosYIYCnR9JnuNTGn4Wz75yY1Xs4cCIQyVcd8Bg22QYhP+fp6D23nhnbMPu1grENZWudPykPdDoe6fXQbdZu6vUWaH5jHTj+ZmJU0x1Sg6/wWKvqYMPiKke2clGuLPTfkfqV0QFMeDlOcY+//lKfQ7DTMzsl2b6c6oamOxdd+l3eG3haiQ92Og7odDnV7ej91e1k0PzAH0oRhbF9MLJdpk+vp1ifdUlONqZ1gjum1UueGllOIsPnYp4uzj9fPB13QzfNyDHCmvGm+aeifY+8H2Q7U7TCo2/Y8qdvU7a3R9MB8Z99O9ie7XJ6cmLx896XgK1JjT/+up/rajdkUllr2mAMDs46mdEauuh7qPdHTNP9PDXxC0redS7YNdTsc6vY8dlC33eeSdmh6YH7k8JHsjXhqGsq2fSqdGLExG1quqd7Yc4HxKd7Wpr1q2GOKur592Jcq8mb6uuib3pmp64/pEFzY8mytTpA2oG77pRObBnU7DOo2dXtpND0wB8o8zdkqu6tBuqZNfRvMlPclVQhydIA2L0xK463toUnB1ekN98qnM9fF0Ly/tvN9jonFtzOYgmJOQqBux6fvC3X7KNRtO9TtZdL8wDyGXE+ZA7aG7VvZp44b83aEpJmjAdvyDWnYpWIvUxnzkrliBPVzzTg/W0dobrOdo3uwbPGJPtOWIffZnEatSUv1gbQJdTsc6jZ1uyQ7+3Zqm7BJVjkwB6afXn2fnqe2pzQgXRzMdHymu2xP9rHo5ZWaViuiMobLI2Hzfrjuk08n7PLc6elOpTEnpjcwV0e9hDpB6kPdDoO6vXs7dTuvbh85fCRLOiSM1Q7Mbdieum2fx7b5nDcwVwPNPV0Vk5atc8jtfUkRHFsHOqRnE26fMg31vrnE3icd04tSum5xCpS0AnXbP71QqNvTx1G3ydwcU9uA0oxNC/k2LBdjT9W5SYk5mzo3xzXYpgFzkmN6z8drUlLUYgYDY4OM2Low5e0xPU6p3HoNB5KTIhuBuk3dHqBu19PtnbczlKUGm/CYu57YbXFi5uexNOdgzFuUi1xP9KWe1kvF3OW2d8ozpOcX4m1x2elTF2Km8HNNsxOSAnV7HOp2HqjbpDUWOzB3CV9IpXTFZMWkYW6zPeWHNhjfeMkpphqs2cmlNOxcNuuUEJwUD9ZYmnNMTY+JvovYOpMDTq2SAeq2P9Rtd5ohULfjoG7XQ5RStW1wctLpJ6lHP+TR1n16hY2Zmsw1nTk2zWYKU2pFD03H5/jQaTofL4h5b9bA3Nfkune+daCUtyqFAwcOXKGUOiPm3BNPPFGdc8451e0g01C3d0Pdrgd1O41LL70U1113ncScS82Op3mPuWsqs0S6Ifg8+Q4enVxP5KHpuDwWPteeEt9pO3ZqSprsptT0ci14n7cFdTssT/P4wQ7TLt9zfbdPHUvdDmNtus3XJdah6cWf5qt6Sky3hRDT4AYPT+3GOtjgmp4tbZ+5OGWOPFPQy6s1sQ2JNxwwr2fuOtlS+ZGyULfzQd0Og7qdF74usQ5JHnMR+UUR+aCIfEBEXiEitxeRe4rIu0TkoyLyKhG5bX/s7frvV/X77zGV/vC01soTeg0bcnmIck0B60x5xWydiEtUYmwreT/GOsMS+ZrlljuPmuJO2qG0ZgPU7Rx5UrfjoG6TNRA9MBeRkwH8PIAzlFLfBuA2AB4P4PcBPEspdRqA6wGc359yPoDrlVKnAnhWf1yz1G4AJUTEJbBzLhjxuS7fa5/DW2Qrs9z52gQ3NL7U7BRsx8zh8dLtMD+TuqxdswHqdkoeY1C397IF3WYoSx1SY8yPAfC1InIMgDsAOATgYQAu6fdfDOAx/eez++/o958pIlGLCmoSuuhmbPsYU3F+KeRs4GNCMpa3zSvTkodtTjtCxH1qejZXHG9sOvr91T1+tQdM5FY2p9kAdduEup3OFnSboSx1iB6YK6U+BeAPAXwSnbjfCOAKADcopW7pDzsI4OT+88kAru7PvaU//gQzXRG5QEQuF5HLb77h5ljzkjCfUPX/JaYWfWwplW5o+vq1xwj52Lmh07dzCXELnc4ULi/RUEZj12AekxqbaZ7PQXkblNJsgLrtsqVUutTtaajbcbbo30kdUkJZjkfnUbkngJMA3BHAWZZDh/cx2jwte97VqJQ6oJQ6Qyl1xh2Ou0OseV7YRM71FJyjwvo0Nts5pUW+FCHpxyxqCu0MQtEFr6TIzxFDOmz3ibtcQodGwiml2QB123YOddudPnXbPz3q9vZICWX5PgD/rJS6Til1BMClAL4LwHH9NCkAnALgmv7zQQB3A4B+/7EAPpeQfxIur0qJPFwNy7chlWpwMZ6L0DLyPT71SR/IV04l6kJpXB2k7Tr0gYZ5rQw7WTWL1myAuh1qwwB1u02o28RGysD8kwAeJCJ36OMOzwTwIQBvBfDY/pjzALy2/3xZ/x39/reoSr9uVKoC+zzVTp1jklu8YmxYCj5TtAO+nZvufSlFjjjNsXri8h7avlPcV81iNRugbofasBSo29RtspuUGPN3oVsQ9B4A7+/TOgDgVwH8kohchS4e8aL+lIsAnNBv/yUATwvJL2fFKyFoY9NQPue6pqZsT8e5iI1XtKWxNKY6zjmEvURerjTMe72WTp34M7dmA9Rt6nZeqNvU7S2Q9FYWpdTTlVL3Vkp9m1LqJ5RSX1ZKfUwp9QCl1KlKqR9RSn25P/ZL/fdT+/0f881nKk4wJK7N93hfci1kMRveXHFyW0KfCtTLtPbUNJCnTurX5ePRY73aHnNpNkDdZvvKA3V7d36sV+sn9XWJs2A+mdoEf46YuDnT0wW/ZGOs4TlpSVimvFq6rTmmLksRYpOtjqVcU4vlQepD3aZul4K6Td1eM4sYmAP2RRIh4lRyIU6pdFqNHUu9Zt8OK/dCK1dZjsX3ubbnrE96x5EarzhVX0zvU+p1UNzJGNTtdqBuU7f19Ei7LGZgnoOQyug7tTQnrdiS2umETNvltMk2SBj+j3WwphiW6HhzT4kO32352Oy3lYFPJ9DqIISsB+p2Hqjb1G0zPdImx0wf0g7D9KAPQ4X1Pd48tkYslzn1q7OF2DLb/fIVbh+PQ2i6pXDZG3KPQ7xFLkwvkn7+WF0kJATqNnXb9zzbMaHploK6TeZiVR5zMwYrpnLmnoaLYSzd2Dxb6Rym7Aidrovx4qR0ljmmEcfqpc+Ua+57acYs5hL3VuocaRvqdlyac0Ldpm6T+ViMx9y3UYZWSnMKybZQRP9cIlZtyi7bvtin7NqYi6PMfbZpOf28UFzlGOLhyXXPQ64lJK4yB7Y6Pnasj5dGv5+t1D8yL9Tt3fuo235Qt/0ordu4NJOhJIhFeMxTPCm+aQ/4egZc58/FmAcgh3egBC7x9ilPXy+MnmapMoi953rsYwl7cqVrlqO5z5a/rUx04W+xPpKyULen7fDdVxPq9nZ1e2ffThbbSBiLGJiXJkWkawtpyNRZbtELmboE7OVsE+QpwXB9n8prDF+vXki6tvKJ8ZpN2aT/D70nvunbvo/dDw7KSWmo2/F5U7fH06RuU7drsphQFrIbnyk82/YWGpvvdG9IKISvB20sHR/xjRXoUlO7OkN5jS0OyoXpmTE7v1bqGiEtQd3eDXWbuk32soiBecmKuhTMxpuSTm5K35eQDivVGzNV11KuNdc91DHT8/FYmceF2DSW1pCOzdNEsd8e1G3qtms7dXsZuk3qsIhQljkqyBIqYas26tNwY1Nyvt6Dse9Tdvimm5JWCwyejjGPk3mMK+Z3Kj50qiympj9bL0tSBup2R6s2Urfnh7pNfGjaY37k8JFbP8c+KS6VJVxfiCfUfCrPNW03NTXsss82deibbiwxQhdq25SXasozM3VP6f0mU1C324a6HQZ1m8xN8x5zW4UsJX62J9mWyb1gJCTfAd/yipm2i7XJlm8ropSjftkW68TmbYs1NAdTtv+2dMy0ltSWSF6o226o29PnU7fH86Zur5vmB+a50StiK43eJLRh1GpIsUJfovx9Og1fe3OXZ06xM9NKSdfmiQnxhFHAyVxQt/NB3faDuk1qsbmBObB3WqdVoXcxdFK5n259y2HsyXwqPf14m7ikeJNCRd5FrvoQE3c5xwAkNp7Q1kG64iVreATJuqFuu9P1gbodZw91m8xN0zHmuRl70mwl/nHKhpKNxleAcuVjE4zU2MWpPHVCYvTmprQttk7EzHOqTErdK0J0qNtpeefKn7o9DXWb5GBTA/OpuLWaDbymsM+Jyxuhe5LGntrHyilUoFqLY/T1fswlqD7iPtBaWZL1QN2uD3XbDXWb5GZToSz6VGJLhIr7WhZomDF3qdOVYwvOfNIusfjF57qmpg9Nm0pPN5pT7mNTu1MdMiGpULfbgrrdQd0mpdjUwFx/sl8qLQl7iThJ33vjirfTBckUHpcQlS7TqfTHOhW9vpqCm9pxTBHa4S65XZF2oW7nhbrtB3Wb1GJToSwt4ut1aUnYB3TbQqcjTcZEa8yj4tpmeipaLD9fXNc6Vb45BzNDWr7luPSBFCFjULePHq+nadsXso26Td0mHJhXxbWwxUf05xQs22IbX4EJTX9sf0x6rU4fm2UYY2PsdfmKr63TdXXALQ9ECMkJdXtv+mP7qdu7oW6TKTY3MG/hadCM7zIbkk1QB2qJu2mzaUesTa7rcU15puZXA59rLFUvxwYRKXna4hdt3jc9P0JioW6HQd1Oh7q9f9ev+JL52FSMOVBPGMZizKZi6HIIaQ5sXoxSC0hs167bUCrfUtS2d8jfNYiwfa4V20mICXU7Hup2PLXtpW5vk80NzGsR07jNBlZDIFydUs5Ox3VdY+n6TnXaOtO50L1Wesc0bNP3Ddg60By4yso1iNDtColPDNlOSOtQt91Qt6nbpAyrG5i3+DTuEkmf81oKB3A9uYdeV2gaodc9Jkxzl6Hp8bCJ+ZTo567TpoCP2e06x2WTr1eREJ0W6wh1e286IWlQt6nbJI5VDcxtwmGjRiMPxTYVWIsxMfD1gKQQ05mEeg1yk1PUclyDjz2++YQONmoPTEjbULfLQN0Oh7odni7Jz6oG5iFPqy1XOpdtrT3BmsJrehnGzhmbYjVxxUe2/FQf0kFPTTunXqMrLrEE+rTvHPmR5UPdnhfqthvq9jz5kXFWNTAH1iPyrTE1BeYzzTZ2nM+9ck1xtirygL9nquQUqG6HT9qxxwwdLtsVCYW6XQbqdhzUbVKT1Q3MB0IaVksVMsQDMSdTcW3D/5QGHhPLGXtuS9jiGAdS66evhyvEW1JzepysG+p2Xqjb5aBuk1KsdmDeKmttAIOwh8awmef4TKeaU4U1hSXH4qfBIzIm9LHkTMt34dDUdkKWBnX7KNTto8dRt0kJVjkwD63Uc1dEfYGL7XNJSk235UrHJy3du1NbRHJeeyq2xVY+njBXGZr3wxUzOuaJqX1/yHKgbruhbueFuk3dbplVDsxbQxdxc5v5uVV8PCI+x7qImb5e2xRc6rWYHhzft124yBHru6b7Q7YFdXsa6jZ1m+SHA3OUEQpzum/pFT0kfq30tS69LG3k8lCYnjzdQ5XTC7KWek2WC3V7Gup2WajbpAQcmPfEVv4lVvTUa43Nb63TYznEM9f0oa2sbTGh5vHmvjliHAlJhbo9DXXbDnXbzyYyP8fUNmDp6LFaQyNtWez1aTP9qTzVZtsq8hKLYqbyLpWXfi0+i2Zi7Mh9D2qlZdanltsD2SbUbdyahi2f4XNJqNt780+Fur0OVu8xT3mVk8+xJZ9YS2EuwPFduBPSEdjKvVTZlPZ+2Rbn+Jwzp8dhrIOdwlZuOTobelxILNTtvVC3w6Bux9lD3a7P6j3moYK9NXxF2+wUQvOwpZWDXNOILpuWIFKu8s2xKCkHrXsjSXtQt8ehblO3XeTUbVyaJSkSyOoH5r74PqmusRPwFfmSXqwYXG8VCM13bfe09kA4V9wlIVNQt6nba4G6TQZWH8rii2sKb22N30WJzs3n9U+xnhzXeUv2oOSitrfRnELeShsi80Pdpm6vBeo2GeDAXGNsKmlNYm+bQislCjaRL7noZ0g/h8gsoVNwdZot2d6SLWR9ULf9oG63A3WbjMFQFtgFoNbq5MGWOVbGx6ZdOj4wV56pdrY+tWe+oaFm/jr6NnpeSCmo22FQt9uAuk2moMcc9cRcR189P1djjcknV/nMXc4tC3Us5kDALNOcHpiUDnaNZU/qQ932h7rdDtRtMgUH5g3QemyXaV+uNH3TmqNMSolQibLTmfIY5uyQbVPprU2/EjIX1O3pY6nbdqjbZIxVDczXUtH0p+lS1xQbm5jrlU4h5CoHmxiVFmCXDaF5hRyfa5pUt9WcftW36ffHjENtbbBC2oO67Q91++g22+dSULfJnKxqYN56bFkoeqOxfY5FF+kQ70eJp2zf9HIKxZjIl8KWhymgPmn4dkilVu2bAj4VU7um9kjKQN32g7pN3Y6Bur08Fj0wtzWMJT3pjTVs23W00lhKrMhvjZydWc7rmyuu1jYNbeZt87aY35fUHsk8ULfrQN0OTysX1G0SwqIG5lNPqkuLmxqr/HNdh08DtDV2Wzo+adlEo4YItFhPpjwpoeWUWq6u+z41+KCoEx3q9rw2DFC354G6TXKzqIG5jVKLJuagRZGxMTxlu57EY4TH916VLKOxQULNtzy0hK0sxgZY5sKi1q6HtAF1uzzU7florU5Qt5cN32O+QGqIj62Rx9gSI55bFolanY2eZ4hYb/leETJG7XYMULfngrpNUliEx9wVk5hrFXNufBeGxDbclOvNuQDJd7osZL/rnLlFrpU6VdOTqHfiS/NqkvpQt/emHwt1249W6hR1m6TQvMfc1dDMJ/+a01gm5lNrabtC8ygdV5b7LQv6vW5FeEvhez/miA0cS7flRW6kPtTtaajb64G6TXLS9MB8Z9+OdXsrQj6FXuFLCr3NA+WawhyjpU4SqPN6rCVh3ve57p0rH7MD5v3aJtRtP6jb24S6TaZoemDuwhSilkQJGBfKEpU+Nb8SAhErPhQFP2z1f44OurW2RpYDdTtvftTt5UHdJj4sIsbcxlgsYA2RMO2Zy+tiy89mj895vueEEBrjFvvGANJGmTGmkYxB3d6b/5g9PudRt5dNC2VG3W6LxQ7Mp+Ko5hZ5nwVNJeL3SnlyantAaudfg7XE/y3RZjIP1G3q9tqgbpPcLHZgDrTxpJmD2KfVkg0pZ9n6dBh6h7WW+5qDpZUFxZ1MsbQ6mYbDVAAADCVJREFU7YK6Td12sbSyoG63xWIH5uarrcyK1VrDmPLItNYwcthkS8OV5uAt28IKfhut1VcfXNPp5rUs8dpIGajbZaFuz0tr9dUH6nb7LHLx50BLC4laE6USQumzICh1Sri1cpyLpXqcpu7XEq+JlIW67Ya6vSyo26QEkx5zEXmRiFwrIh/Qtt1ZRN4kIh/t/x/fbxcRea6IXCUi7xOR+2vnnNcf/1EROa/M5WybEgs4Qjwwsce16HnaKiH3Ys43WJAwqNvLgbpNUqFurwufUJYXA3iEse1pAN6slDoNwJv77wBwFoDT+r8LALwA6DoEAE8H8EAADwDw9KFTiKW1J7oW7BmbIs6Rtu3z2DYXtikz/Y8sH59FdaQoLwZ1e5IW7KFuk1agbrfB5MBcKfV2AJ8zNp8N4OL+88UAHqNtf4nqeCeA40RkP4AfAPAmpdTnlFLXA3gT9nYak1AA2sHsSHLcF3pg6qML81RHrn/WX5fGNlof6jaxQd1eJ9TtdRG7+POuSqlDAND/v0u//WQAV2vHHey3ubbvQUQuEJHLReTym2+4ede+nGJC0tFFPocw8762gc8r7Wz324yPpdg3B3WbULdXCnV7PeRe/CmWbWpk+96NSh0AcAAATjr9JOsxrUFvAUnFRwjnXGhkm9J0LUzT38xAFgl1m5AIqNukBLEe88/0U53o/1/bbz8I4G7acacAuGZk+ypYcsVesu1rwfcezPErhGMxqS5xt/0nTULd1lhyXV2y7WuBuk1KETswvwzAsEL/PACv1baf26/yfxCAG/sp078B8P0icny/eOj7+22jHDl8hBWmILXePUtP1W5ylUdKOq532wLjU5u8l4uCur0CqNttQN0mpZgMZRGRVwB4KIB9InIQ3Sr9ZwJ4tYicD+CTAH6kP/z1AH4QwFUAbgbwkwCglPqciPw2gHf3x/2WUspcmOSkZZFfcgWvabtrCm2rP1RRe0rRNa1p3oup/aQNqNvjLLneUrfbgbpNSjA5MFdKPcGx60zLsQrAkx3pvAjAi4Ksm5najYx08JVNcaTWXVeHO8VYu+HCvzpQt8ncULfjaFm3cWmSaSSS2FCWVWB7d6zvSnW+IioNH+FgZ9sWrvsx9TYA3keSE+p2Pajby4O6vTw2PTDXSRVrVuIwQn6lbCtl2+KAYcqmrdwb0ibU7Xmhbu+Fuk1yI90sZpuIyOcBfKS2HQD2AThMGwC0YUcLNgBt2NGCDYCfHXdXSp0Yk/iJJ56ozjnnnJhT93DgwIErlFJnZEmM7KER3V5Su9iCDUAbdrRgA9CGHdTsRsn9HvPcfKSFmyEil9e2owUbWrGjBRtasaMFG1qygzRBdd1upT62YEcLNrRiRws2tGJHCzYQOwxlIYQQQgghpAE4MCeEEEIIIaQBWh+YH6htQE8LdrRgA9CGHS3YALRhRws2AO3YQerTQl1owQagDTtasAFow44WbADasKMFG4iFphd/EkIIwIVEhBCyJKjZ8bTuMSeEEEIIIWQTcGBOCCGEEEJIAzQ7MBeRR4jIR0TkKhF5WuG8XiQi14rIB7RtdxaRN4nIR/v/x/fbRUSe29v1PhG5fyYb7iYibxWRK0XkgyLylLntEJHbi8j/FJF/7G24sN9+TxF5V2/Dq0Tktv322/Xfr+r33yPVBsOe24jIP4jI62rYISIfF5H3i8h7ReTyftus9aJP+zgRuUREPtzXj++cuV7cqy+D4e8mEfmFGmVB2oWaPb9m9+k2o9u1NbtPu7pu19bsPl3q9kJpcmAuIrcB8HwAZwG4D4AniMh9Cmb5YgCPMLY9DcCblVKnAXhz/x29Taf1fxcAeEEmG24B8FSl1OkAHgTgyf01z2nHlwE8TCl1XwD3A/AIEXkQgN8H8KzehusBnN8ffz6A65VSpwJ4Vn9cTp4C4Ertew07vlcpdT8tvm3uegEAzwHwBqXUvQHcF12ZzGaHUuojfRncD8B3ALgZwF/OaQNpG2p2Nc0G2tLtFjQbqK/bVTUboG4vmSYH5gAeAOAqpdTHlFJfAfBKAGeXykwp9XYAnzM2nw3g4v7zxQAeo21/iep4J4DjRCT5N3mVUoeUUu/pP38eXUM+eU47+rS+0H/d6f8UgIcBuMRhw2DbJQDOFBFJsWFARE4B8EgAL+y/Sw07LMxaL0TkTgAeAuAiAFBKfUUpdcPcdmicCeCflFKfqGgDaQ9qdgXN7vNuQrcb1mwzv6L3pEHNBqjbi6LVgfnJAK7Wvh/st83JXZVSh4BOgAHcZS7b+mm9bwfwrrnt6Kci3wvgWgBvAvBPAG5QSt1iyedWG/r9NwI4IdWGnmcD+I8A/rX/fkIFOxSAN4rIFSJyQb9t7nrxTQCuA/Bn/RTxC0XkjhXsGHg8gFf0n6u1EdIcLdzzTWp2n38Lut2CZgP1dbs1zQao24ui1YG57cm5lfc6FrVNRL4OwF8A+AWl1E1z26GU+mo/9XUKOi/Y6SP5FLFBRB4F4Fql1BX65rntAPBgpdT90U3xPVlEHjJybCkbjgFwfwAvUEp9O4Av4ujU45x2oI8P/SEAr5k6tJQNpFlavuer1mygvm43pNlAfd1uRrMB6vYSaXVgfhDA3bTvpwC4ZmYbPjNM4/T/ry1tm4jsoBP4lymlLq1lBwD0U29vQxc7eZyIHGPJ51Yb+v3HYu/0cgwPBvBDIvJxdFPiD0PnjZnVDqXUNf3/a9HF5j0A89+PgwAOKqXe1X+/BJ3o16gXZwF4j1LqM/33KnWTNEkL93zTmg1U1e0mNBtoQrdb0myAur04Wh2YvxvAadKt6L4tummYy2a24TIA5/WfzwPwWm37uf0K5gcBuHGYFkqhj6+7CMCVSqk/qmGHiJwoIsf1n78WwPehi5t8K4DHOmwYbHssgLcolf6LVUqpX1NKnaKUuge6e/8WpdQT57RDRO4oIl8/fAbw/QA+gJnrhVLq0wCuFpF79ZvOBPChue3oeQKOTocOec1tA2kTanY9O6rrdguaDbSh241pNkDdXhzHTB8yP0qpW0TkZwH8DYDbAHiRUuqDpfITkVcAeCiAfSJyEMDTATwTwKtF5HwAnwTwI/3hrwfwgwCuQrfK+SczmfFgAD8B4P19rCAA/PrMduwHcLF0b1j4GgCvVkq9TkQ+BOCVIvI7AP4B/aKW/v9LReQqdN6Ox2ewYYxfndGOuwL4y67vxTEAXq6UeoOIvBvz1gvgf7d3xy5yVVEcx7+/WChpRCIEuzQRC0GLhYCYIIQlMQELFSxETRcEFfwHhKRKLSJksLAJRLCQYBSDhawigst2rl0IaURYsBIjqMfivuyGZTe7yc7O3J18PzAwPIb7bjNnztx37znwLnBpSHiuD2Pvm+Q8kuwH5oGzd1ye9HdEnTJmTy1mQ99xe5IxG/qJ21OP2WDc3qsyhj+pkrSrJtneOclJWrmzh4BPqurC3cZLcgaYq6p3xjJBVg8TfllVTyeZA96sqvc2+NyN4d4r47q3JO2UMfv+Y3aXK+aSNA1Zq8c9T9tz+XOSK1W1PK05VdUisDit+0tSr2YxZve6x1ySpmFH9biTfJrWPe/HJNeTvDpc/yzJqXWfeyXJoSTfJ1kaXs9tMOYLWeukeCDJtbQybBfZuJKCJD0oZi5mu2IuqXsrKyvfjEajx8c03CMZWnUPRlU1Gt5vVMv3CECS88BiVW11qPEJ4HngKdqBqs9pPxavAV8N+06PA2/TgvR8Vd1Kcph2SGvTR7a0vdQ/VNX5JKdpHfokqSvG7FX3HLNNzCV1r6rWt1/fLZvW8q2qD7Y5xhdV9R+wnOTgcO1r4MMkD9NayS9U1V9JHgU+SvIs8C/w5BZjHwNeHuZzNckf25yTJE2MMXvVPcdsE3NJWjOOWr5/3/E+AMPqynfACdoqzO3yZe8DvwPP0LYW3trG+J7Yl6Rm5mK2e8wlac1u1uO+TCtBdpRWVhBaY5XfhtWaN2hVBe5mAXgdIMmLwGNjmpsk7UUzF7NNzCVpUFX/ALfrcf9Kqwn9C7T9ikle2sHw12iPNb8dDikBfAy8leQn2iPRP7cY4xxwLMkSrXnKzR3MR5L2tFmM2dYxlyRJkjrgirkkSZLUARNzSZIkqQMm5pIkSVIHTMwlSZKkDpiYS5IkSR0wMZckSZI6YGIuSZIkdeB/PuKRRKYqHqoAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 864x504 with 4 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "plt.figure(figsize=(12, 7))\n",
    "\n",
    "plt.subplot(121)\n",
    "plt.imshow(mask, cmap=plt.cm.get_cmap('Accent_r', 2))\n",
    "plt.title('Original Invalid Area Mask')\n",
    "cbar = plt.colorbar()\n",
    "cbar.set_ticks([0, 1])\n",
    "cbar.set_ticklabels(['0: Invalid', '1: Valid'])\n",
    "\n",
    "plt.subplot(122)\n",
    "plt.imshow(updated_mask, cmap=plt.cm.get_cmap('Accent_r', 2))\n",
    "plt.title('Updated Invalid Area Mask')\n",
    "cbar = plt.colorbar()\n",
    "cbar.set_ticks([0, 1])\n",
    "cbar.set_ticklabels(['0: Invalid', '1: Valid'])\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Generate updated stratified design"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now we have our inputs, we are ready to generate the design. \n",
    "\n",
    "This function keeps the sites which have already been sampled fixed, and then adjusts the location of all unsampled sites based on the new mask constraints."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "metadata": {},
   "outputs": [],
   "source": [
    "def update_stratified_design(mask, sampled_csv):\n",
    "    \n",
    "    nsp = len(sampled_csv) \n",
    "    sampled_df = sampled_csv.loc[sampled_csv['sampled']==1]\n",
    "    nsampled = len(sampled_df)\n",
    "    print('{} sites already sampled and will not be moved'.format(nsampled))\n",
    "    print('{} sites to be adjusted based on mask update'.format(nsp-nsampled))\n",
    "    \n",
    "    imheight, imwidth = mask.shape\n",
    "    sampled_x = sampled_df['row'].values\n",
    "    sampled_y = sampled_df['col'].values\n",
    "\n",
    "    sites = np.ones((imheight, imwidth))\n",
    "    sites[sampled_x, sampled_y] = 0\n",
    "    mask_aux = copy(mask)\n",
    "    x_vals = [sampled_x]\n",
    "    y_vals = [sampled_y]\n",
    "\n",
    "    for i in range(nsp - nsampled):\n",
    "        \n",
    "        print('Plotting site {} of {}'.format(i+1, nsp-nsampled))\n",
    "\n",
    "        # Generate EDT image of all sampled/selected sites\n",
    "        dist_im = ndimage.distance_transform_edt(sites) * mask\n",
    "        \n",
    "        # Extract coordinates with the maximum distance value\n",
    "        dist_mx = list(zip(*np.where(dist_im == dist_im.max())))\n",
    "        \n",
    "        # Choose one maximum coord pair at random\n",
    "        idx = randint(0, len(dist_mx) - 1)\n",
    "        x, y = dist_mx[idx]\n",
    "\n",
    "        # Make all elements of EDT map in invalid region 0\n",
    "        dist_im = dist_im * mask_aux\n",
    "        \n",
    "        # Extract coords of pixels with maximum distance value\n",
    "        dist_mx = list(zip(*np.where(dist_im == dist_im.max())))\n",
    "\n",
    "        # Choose one max coord pair at random\n",
    "        idx = randint(0, len(dist_mx) - 1)\n",
    "        x, y = dist_mx[idx]\n",
    "\n",
    "        # Save x and y coords\n",
    "        x_vals = np.append(x_vals, x)\n",
    "        y_vals = np.append(y_vals, y)\n",
    "\n",
    "        # Code chosen site to be zero (a feature pixel) in site array\n",
    "        sites[x, y] = 0\n",
    "        \n",
    "    print('Adapted stratified design complete!')\n",
    "\n",
    "    return x_vals, y_vals\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 22,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "9 sites already sampled and will not be moved\n",
      "21 sites to be adjusted based on mask update\n",
      "Plotting site 1 of 21\n",
      "Plotting site 2 of 21\n",
      "Plotting site 3 of 21\n",
      "Plotting site 4 of 21\n",
      "Plotting site 5 of 21\n",
      "Plotting site 6 of 21\n",
      "Plotting site 7 of 21\n",
      "Plotting site 8 of 21\n",
      "Plotting site 9 of 21\n",
      "Plotting site 10 of 21\n",
      "Plotting site 11 of 21\n",
      "Plotting site 12 of 21\n",
      "Plotting site 13 of 21\n",
      "Plotting site 14 of 21\n",
      "Plotting site 15 of 21\n",
      "Plotting site 16 of 21\n",
      "Plotting site 17 of 21\n",
      "Plotting site 18 of 21\n",
      "Plotting site 19 of 21\n",
      "Plotting site 20 of 21\n",
      "Plotting site 21 of 21\n",
      "Adapted stratified design complete!\n"
     ]
    }
   ],
   "source": [
    "x_adpt, y_adpt = update_stratified_design(updated_mask, sampled_csv)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Results"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Plot Design**\n",
    "\n",
    "We can now check how our design looks over the landscape."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 28,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYAAAAIYCAYAAABkGjsMAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvnQurowAAIABJREFUeJzsvX+YZsdV3/ktu7vHM2PJYDmhW5ZbTqSE2O4ILWB1d4Jib4ZsDwlKggKamJAwlkarySZPlF5nA5uwohEsIdn1NmZJdgb9mgRju0WigIcN0yQTMGLd0yaxJ0oT2KwMuEXUbWLZII9mNP021P5x721VV1fdW1W3ft17z+d55pl+3/e+VXXrxzmnzjl1X8Y5B0EQBDE8Xpe6AQRBEEQaSAEQBEEMFFIABEEQA4UUAEEQxEAhBUAQBDFQSAEQBEEMlEEoAMbYScbYL6VuBwAwxpYYYx9O3Q4RxtjPMsa+U3j9A4yxLzDGthlj04yxK4yx1zuU+3bGGGeMjfltsVHdv8IYe2+EerzcY9nHf9hXuxzqP8MY+19qPs9u3saEMfYLjLFTnsvct+5S0GkFUA7KlxhjhyLW+ZuMsW8MWP7fY4z9RikQfosxtiJ81noSqhYy5/ybOOf/tPz8bQA+AOCdnPNJzvkm5/yNnPPfa1Ovpi2/yRi7xhj7MmPsdxhjn2SMnWaMtZ6XnPN3cc5/wUMzo1D28a8nrP805/z7AYAx9l7G2G+1KY8x9vOMsf/KGHuZMfYfGGN/Qfr82xljn2OMvcIY+ynG2Jvb1NdFpHWXxEjtrAJgjL0dwN0AOIA/n7Qxniitgb8K4Bs5528E8PUALlp834elfSuAlzjnv+2hLBPu4ZzfUNb7QwC+C8ATkeomwvEwgCnO+Y0A/nsAH2aMTQEAY+xdAM6imOtfBeAqgH+SqqGDhnPeyX8AHgHw/wD4PwD8jPTZTQA+DuBlAJ8C8P0Afkn4/EMAXig///cA7hY+WwLwzwGsAPgygE8D+Jrysx8H8PsArgG4AuDvlu/PAfgkgN8B8B8AvFco7w8B+ERZ1r8G8KMAPqy5px8F8MOaz/5XAL8H4NWy7h8t3+cA/gaA/w/Ab9TdH4DjAHYAjMoy/kP5/i8AOAXgG8t7+/3y83MA3l7WMVZe+yYUAnoLwH8B8AMAXl9+9noA/zuALwD49bJde99V3NNvolB24nt3lfXPlK8PlWVuAvg8gDMADpefvQXAz5T9/kUAzwJ4nVw2gMMA/imALwH4VQB/F8BvSe34OwCeA/C75di/QdPm2nts6J/by7nwu+X3V4RyOYDbhfl7vhy/Xy7L+CXp2tPlmH8JwD8GwBRtfUM5nm8pX38PgF0AN5avfwDlfCvH+gcAHJXmwBUAN6NYF08D+Gco5vKvAPh6w7V6F4p5e1f5+gcBfET4/DYU8/IGxXcZgGUAv13223PC3PhzAD5T9tMLAJaE77297Kf3l599qeyzd5dl/A7KNVRefxKFPPk/y3p+DcAx4fNfAHBKeH1/OZe+BGAVwK2ae38DgA8DeKms85cBfJW07t5R9s/vlf39O23mvpUcDSGcY/wD8DyA/wHA16EQaF8lfPaxcrIeBTCDYiGKC+g7UCyyMRTujm2UC76c6CMA3wpgHIVg+A0A4yqhBeCt5eD+WRQ7qj9Tvv4D5edrKJTUIQB/CsXi0SmA7ygH839CYf2/Xvp83yQUhMG/BvBmYXI03d+HdeUCeC/2C8e3Y7+A+ykU1ttRAH8QhYJ9qPzsNIqF87ayPT8PSwVQvr8J4K+Xf/8wCmX+ZgA3oBCM/6D87B+Ui2K8/Hc3SkGI/Qrgh1AI3q8EcAsKASArgE+hEHRvRrGwT2vaXHuPDf3zUQB/v5wnbwDwDdI4VgrgY+W/IwDeiUKAyQrgZwB8BYBpAP8VwHFNe38RwF8q//45AJ8F8E3CZ99S/n0OwA+o5oAwb15FMc9fX/b9pYY1+jPldziAC3hNOf80gO+Srr0C4OsUZSygMGK+AoUyeAeKnUXVzj9e9ucdKITkX5Tm7Zmyr/+7si0/VY7LW1EolfeU159EoRwXUcylEygUwZsVa+QvopA/70Cxxr4HwCc1ffAQijl7pOy3r8NrClgs86Q4xm3mvpUcbSOEU/0D8A0ohHRl2fwagMXy79eXn/0x4foflDtXKu9LeM3KXxIndjm5tvCaFf2b2K8AvgvAj0vlrQL4ThSLcxfAUeGzj0CjAMrP/wqAfwPgFRSK5LuFz/YmjPAeB/CnG/pLvj8nBYBiu34dpaIpP38fgJ8v//63EAQnikXnogAuoRCUrOyH24TP5vHaTudRFMLk9rqyUVjqC8Jnp3BQAXyH8PofATijabP2Hg36558B+DEAtyjK5Sh2CNX8/WrhM9UOQFQeT4vzRCr3+wH8SNm+bRSumR/Cwd3BOTQrgH8jvH4ngGsGa3UcwDehXJ/lexchKVgURtp7Fd//0wD+M4pddq2Fi0JgLkvz9q3C5y8BOCG8/hcA/nb590kAL0IQoiiU919VrJGfBfCAcN3rULixblW06X4U3oE7FJ+JZZ6Uxth57tv862oM4DsB/Bzn/Avl64+U7wHAH0Ax2V8Qrv+c+GXG2AcYY7/KGPtdxtjvoNi2v0W4ZO+7nPPfB/BbKKxDFbcC+LYyiPk7ZXnfAGCq/M6XOOev6Noiwzn/Cc75N6KweE4DeJQxtlD3Hey/V5P7c+VWFAt6S7jXsygsKqC4X22/W/BWFDuhP4DCcvr3Qn0XyvcB4H9DYYn9HGPs1xlj360pT27XC4prtoW/rwJ4o2FZ4j029c/fRbGwP1VmKd2vKF81f9u09xMoBPrXAviPKHaL70EhUJ8X1pAJcp1vaIo7cc5HnPOfBbDAGKtidVcA3ChdeiOK3bH8/X+LwjX6jwF8njH2Y4yxGwGAMTYrBJt/F8V6kef554W/rylei/32X3gpXUs+B/W6vxXAh4Qx/iKKcX2r4tofR2EQfowx9iJj7B8xxsYV18n4mvu1dE4BMMYOA7gPwHvKNMVtFNu2r2GMfQ2K7fAuii16xbTw/btRWO33AfhKzvlXoNjqMeH6twnXvw6F2+DF8i1xggDF4vxxzvlXCP+Ocs5/CMXO4SsZY0dVbamjXDg/idLnqakb8vsG96crw4QXUFi4bxHu9UbO+bvKz7eg6XdTGGPvRrGQfgmFn/wagHcJ9b2JFwFycM6/zDn/AOf8DwO4B8D/yBg7pih2C8UYVrxNcY0pdfdY2z+c823O+YOc85tRuAb+CWPsdqn8av76au8nAXw1gG8B8AnO+X8q2/znUCgHFW3miI4xFL5+oIgffE31QZn+egiFpX+wMZz/COf86wC8C8AfReEiBQrD7+MA3sY5fxMKlwhTlWHIWxlj4ven8dq6F3kBhVtPXPOHOeefVLR9xDn/Ps75OwH8CQDfDOCvKcqU+9zX3K+lcwoAhf/t91BsQe8s/70DRRDkr/EiXfEZAEuMsSOMsXfitd0BUPjSdlEstDHG2CM4aI18HWPs3tK6+dsoFvWl8rPPAxDztT8M4B7G2AJj7PWMsTeUaXS3cM4/B+DfAfg+xtgEY+wbUAyWkjIV7M8xxm5gjL2OMfZNKCb9uqZuFU3393kAb3dJteScb6HwI3+QMXZj2cbbGGPvKS95GsDfYozdwhj7SgDGVklZ3jej8H1/mHP+H8vd12MAlhljf7C87q3Vjogx9s2MsdvLRfsyinmhSld9GsD/zBj7SsbYWwH8Tdt7l8pS3mNT/zDGvo0xVgn2L6FY9Pvaq5i/fwxqgWEE5/wqCh/638BrAv+TKBSQTgF8HsBNjLE3udTJGPtjjLFvYowdZoyNM8a+A0X8q6rvJ1CsmbtL4+hRAM9wzg/sABhj7y4t/XEULpEqWAoUc/2LnPNXGWN3Afh2l/YK/EEUYzvOGPs2FHLlXymuO4NiPr2rbOObyusPwBj7bxljf5wV52heRuHeU83RzwO4hTE2Aex5HnzM/Vq6qAC+E8BTvMhP367+odgm/pVSaP9NFFu7bRS+zaeE76+i8OH9ZxRbvFdxcIv90yiCQF9Ckap2L+d8VH72DwB8T7kt+zuc8xcA/AUAfw+F0H0BhYVS9e23A5hFsU38XhR+YB0vl+Vsooju/yMUwdAqP/hDAL6VFWcffkRTRtP9/WT5/0uMsU/XtEXHXwMwAeA/oeiff47C3QUUE3YVRSbUp1EIsibOM8a+XLbx76MImL9f+Py7UGx1LzHGXkYRH/nq8rM/Ur6+giLY/k+4Ovf/URRuvN8or//nKJS6C033WNc/7wawzhi7gsJyfZhz/huKOv4mCrfdNgoXwkdbtBcoBO84Cp929foGFEHgA3DOf62s89fLea5zf+pgKGIGv41iTTyMwvf+6bL8X0HhrvmJ8pobUCR0qLgRRZ9/CcV8fglFZgzK7zxazp9HUCjnNqyjmFNfQJF1962c85fkizjn/xLAP0Th1nkZwAaKOIeKSRRz4GUUyQWfQGE0yvxbFDujbcZY5ZbzMfdrqTImiBLG2BKKwMp3pG4LEQbG2F8H8Jc55+9pvDgDGGP/EMAk5/w7Gy8mnGCMnUQRkP2G1G2JSRd3AARhBWNsijH2J0uXzFejSI39l6nbpaN0odzBCu4C8AAybi/RXaI/o4UgEjCBIhvnD6FwrX0MeZ88vQGFC+ZmFC6SD6JwSxKEV8gFRBAEMVDIBUQQBDFQSAEQBEEMFFIABEEQA4UUAEEQxEAhBUAQBDFQSAEQBEEMFFIABEEQA4UUAEEQxEAhBUAQBDFQSAEQBEEMFFIABEEQA4UUAEEQxEAhBUAQBDFQSAEQBEEMFFIABEEQA4UUAEEQxEAhBUAQBDFQSAEQBEEMFFIABEEQA4UUAEEQxEAhBUAQBDFQSAEQBEEMlN4qAMbYGxljTPMZY4y9MXabiIHD2BuhmZNgjGEIc3LgfZCbXOqlAig78RMAluXOLl8vA/gEKQEiGsKcPCAAhTnZawE48D7IUS71UgEAeAXAswAehtDZQic/XH7+SrIWEkNj35zcE4DDmpND74P85BLnvJf/ADAAPwyAl/8feJ26jfRvYP8AxoEf5gAv/z/4OnUbqQ+C/stNLrGyUb1E0qwVHwKwyPt84xlTbm9fUfV/OV5HOedX4rcsEjVzEkOZkwPvg5zkUq8VALDX2b9fvX7wwQe1MSiRqQemsPXEFqYemAKAfX9Xr6vr5PdiI7YxdB11qPpLbNNoNML58+cxOTmJ+fn5fePAOcfa2hq2t7dxzz33YHx8XFm/qt9TIc8JYzjH0tyjey+XLj0CGMzJmPPLtX/lNmrLMeyDVGsqNJxzPPbYY3uvH7n0yAG5tHTXUvOkaElfYwAA9mnaPZ57+TljI6NOsE49MLXv85QTVaw/lGA0uT+xblH4V++PjY1hcnISGxsbWFtb2xuHSvhvbGxgcnISY2NjjW1I3eeAY19zjuPLq/veOr68WjhANKS4z7Z1iuvjABZ9kIOi900130VWl1eTbH5ev7S0FL3SGEjbrA89cumRuZcuvYTP/PxnsLOzg3d/77tx5TN6T8MNX3vDvr/F16prb/jaG2rLi0WoNrgsxKpfKmHCGMMtt9yCnZ0dbGxsYGdnB7fccsue8J+ZmTmwM6jQ3deVz1ypHZusKAXf3Mo6Lp2YxeNPPIA3fPlVzK2s4w1ffhXPz922ZwVvPbGFK5+5snffoefW1ANT++pwFbzVOtCOiUUfVNTdu9zu3BGNnZmZGZz++Gm8+uVXsb6yjle//Cpum7ttb/6/963v/b7Q7VGbWh1HFv4AFhljD9/7+L3AKWBjYwNHl4/ijhvv0LqDZBdG04JI6f7R1S23ObXFDBRKYH5+HkAxDhsbGwBQK/x7gST4LiwuAIwV/wOYW1kHADx14x1G7qAQ1M0lW5TuMcM+qN43aW+XkIX//Pw8tp/cxkJ5/+vl/S8sLkRbB311AR0FcDekwEolfGZmZrB5eRO7u7sHvli7ddWQg2BVUbWrrbvE90ITlUBFG+Gfa/+LTFwbYfry5j7BB2BPAF46MYupi5/FIcWcjIHYh23Gu85gMumD6cubmLg2UrZLJEbcyze7u7vY3t4+YOwwxrCwuIDZE7PYvLyJkXD/oeltEFjONln61NLejXLO8cLZFzB9eloZtDIOZEEd9IodrJODr/Lnqvdt62iLWL9oCVW03QGksgZtAsETV3ewc3hcbd1yjolrI+wcmXAeK5e5p0tycEEXKxMx7QNdeaoyVQkHuTIajTA2NrZvnlf3wznH6NoIE+X9xwgC99IFBAB1qYSMMaXwr/BhAcVCXri+LSNfgrVaoKptsKgMXJRAaleAqfCsBJsSxuo/N2yHKy673qbvqD437QMTZSK+3wXhD+BAdpsIY2xP+Meiry4gI6oJVbl9UgsSV+TJ73Mx+Cxr8v5JPPfycwcCvpVbTs4OIvJDdCc27ZTlrC2bOsQ1qVubXRH6TaS8j97uAEyR8/i7rgRyVmSjayNsXt7E7InZfQF4MSawvb2N3d3dWkspF0ItXJVbpi/CzgTT+dsly7+OlOt18AqgDblMPpUSSx2bUDFxZAInz5zE+OFxbD+5ve+zSgnYCv+UiydGn8rBWd/1uZ4fCX39kKg7cBqa3gaBZcQgsIhNwFf1vdQWWlOqqs8AWYiJ6ePAUWpiLNrQuw2bNtgmRYQmtWHTFl0f0UngTMlF+KvaIiML/y4FuIm8Hn0B5NMOwg+DdwG5ps5VwbAuLAjxHk0f6RBD2PdFoYScA6HHQveMK11b2pbvk77MH18Gmgu0AxCwcf9UwdauTEKbTCebFLy2bbJ53wc5PEPIhq4r4hzdhjkT+94GtwNQaVtbQZ6D66cNLoorpCXn2y2ly+zq2njFbG9TsNkmM8e2TtM4lvgd332T0xyJ6VkYlAJQTTjXAc9hogDu2/JciGEhqu43Z/dd6vGJVb/oSm1qS+w2pSTmvByMC0gnBJquMS0rFbYujZzanpo+9UVOysxlNy1+N2bGk+gW7Zp70AeD2gFUtPVx52AliO3IafH7IIe+TUXqHansCvHxfCCXR0bYXO9676bBb9ey24xlrDU9mB2Ab1IL3VyUUJdJPYYyKc+SyAkCbU+U+7ambQ00Fbodfw5nLFIxyB2AD3ISvqYWg89DYVV5gP+JHjogl+vCtG1X2wB3rv0g0tRG3dw33SWEmmN1mYIma3DriS3griBN2wftABzIQfiLwtc2Q8PFuqubtD77Q7TKTBZ2mzp8kmpO5DAXm/DpnpLHv+lUsrgTsUmFbkuKFGcXBrsDaEo96wIx2y0LZt122qeA7mN8IwRdUAJDRbUuVNmIqRisAgDCH25KPbg62rYr5cEhH3X7HvdUiiqn3PWKnNpSEXNsTHYocqwFSNdvg1IAKSxmX+RgCZv6Lun0Zzfo+w7LNBjsk6716WBiAKEClSHKVpGDAHRpg+/0RB/0RTn7yLQJnXcfy+duSuh15HqvcnyO0kA7gu4gi88B9DkhxOwEn8gTN5Rg8dkPXaZPgew+0bU+JAXgiCzkZaGa40RQpaP5oOmEpU98KdecrFIXRIvRB9Sv6UnRd4OKAXSRUJNCpwy64B6pa2MOmRUhcX18yZChPtJDOwAHcvNr2qJ7LlLOh4l8L+KuCoWQ/e2jT3Lq15SnkbsCKQBPhJhsOQc+gXjB71APu8tJWDUhH2oKWU+XaWp/28dbmJTZpT4kBZARvv26FV2akLa0UZR9tOiIekKcI+ny+iIF0BJfEyqF8Pf1wK+Qz1OpaDry31diWP26OsV6TdsQ0yUotlH+O3TdKuQsuC5ACsADXRnsCl2qZm4WscmCErNXuh6bkclpXuXUryb9EqPvdIZbl5QAZQE54vt5HiFOEKraVfeQKtf7SC0curDQbMjhfnw908mlHN33cugXFaEeWhgDUgAW6LJnQpTrs8y6yVgJ/tjCX7XIYx+j9/kAO5+kSmXNpR9sz9Pkqhi6ALmAMiFUxk+dxR8q7mBC3XNacnzgW2xs29/mfnM+AyK79bo+rrlBOwADYlgYIaxfk/Jycvv4OkktK7auWog22SZt7rELDzBru1Ml1JACyIQUwr/r2Ci4OrdKFwSgDp2S67ryq+h6+3Nn0C4gE4si1gRMcWLR9eRvaGFp0i7bNoQcx1hzRHbbNaXG9uFQXM6ZXX3YkQxaAdQhWo4xJmDsSZ7zxI3dF7kKGB2xng+VkroHF+Y2Xl1WBKQAUJ8uGWtwu5D10bXn/lSundDtTi2QfPWXqUvN9ACgz3EU25ersO2iIqAYQEklLHI+GOUT21TDrgl/QP3zfF1boE34vB/d3K+zxnW0nS/VWFUKvEvj1qWYEimAGro06WzIUfCL9cmBWx+WX9eEiA0h7q0p46jpbIkP+jpmORmZ5AIS6ONka0PKyen7jEKfx7bP99ZVQiQyhGCwOwBaNGpymJQh29BXqzImMeZIyDGKMQdyWEcmDE4B5Lz4u+Q7TIGP5xXlPP6A3RzQJS+EFp4xsDkEZ0vspI6c1/RgXECqCH3OAxOSusdD5I5rVk/ugr+i7RiEvk+bDKAUZ1tyRuyTXO5ncDsAkZDPPnEtO9YuoCsWsS9yWXA+aWvtt/1+7LnTdfed2P5cdvuDVgC+6eJDzHKYhCnp4v23nWdyppVr2aq+C2Xhdlnwi8hKQHdNLAbjAopF2+yVvkz0HFH1bVf629T1ouPQaARwvvd63/zkHBNXd9o2MYu+zFGhi+Nm0kcxD5SRAgiAyVa1zpedw0LqI6qDYTkKDBEfwuDQaIQPnD+P+9bWAM73l8c5ji+v4uTpc5i4umPcJ6prqu+GcNXkPk4m5OL2ESEF4AEXQRJCy09c3dln5e3Dk5XXdXw9crpLXB8bw/OTkzi2sbGnBADsCf+5lXVs3jmNncPje9/xcZI3BE3l5j6mubWPYgCOqKwc+YFVMQd74uoOTp4+h807p3FhcQFg7LUPy4U+fXkT586cxM6RiWjtyo3cFmATXhQWY3h6fh4AcGxjA0ffdRQX+NSe8L90YvbgnEF+Vndu7ekDtANwRFyMIYSKbZk7h8exeec05lbWcXx51cjKGyJdESLy+LeeY6USuDgzg7mVdSzNPVor/NvSNmZRh+hqql6nJoc2uDBIBeDD95vdaUjGcGFxAZdOzO5TAk1W3hDpwmIN8jRaYSdQ0aU5oeuHHHZ1ObTBhUG6gHwMVl0anU8BY1VeqQQAYG5lHXMr6wBAwl9B13PKneC8iAEIHF9e9T43Qp2CHWL8JjSD3AG0pW7raXrAyvRpitYLSFACFST81eSuBLy2rRT+xzY2cHFmBkuXHjmwWySGBykAB2xyek1QCXxnN1Xp9hGhBa4n18dKtDEgDiAJ/6fn57H15LbSZeiTUK62ECmmXXALhmDQCqDNoLedNLILycsElHz+ZOWFIdbOoW5O2NR/aHcXt29v7wn/vd0gY3jqxjtw6cQspi9vYuLaqG2TndtIpGGQMYC2+BDWde4jJzQBXzEmAJA7yBehLUbZMGgjTK+Pj+OD99yD62Nj+8Z+64ktTJ26GRf4FCaujbymB/tUkqGfrDlkRTVIBVBNzlADX1duKMExcW2E6cubBwO+ghKorDzXhZ7bkwyHgo95en28Jv2XMe9nQ2zbLM+tuoyfoT3IMCSMD8Qt8NBDDylv1MVSyVUATlzdKfL8VRY+596svJBH2knJ7Ed32DDXRy27WuvynIr9gy05KpOzZ88G36oPcgdQ0ST8XZRDSgFWK9w9WnlNzzByeSyGSflDJNZ5k1CuFdNyY1r1Te5X38pVVoo5KZvBBoFNBsJloIYsvGwC4/JJUfl0J6EnRPZZCNoaAqHwmmFlWVduD4QbpALQpVrqjpjbDlhOA5wbNkI/J0upb1QGUE5KN5d2AM3ndFzbmtM9AgOMAeQ2ADmSm5VCvIbLD7bUxRFCWrtiXeI1ovKxKcsXJvdcdy+ubbN1M1EMwCMk0OoJnWpH+MHGh6yy8GP52MXDkrqDjipCPTzOtnz5O/I99CVxZDAKYOjUPbtI9X5udKWdOdEU7AzZlz5dJG2UgrzbcWlXDHfPTe+7CV/4yBfAFBl8rHjzKOf8ilNDaiAFMBB0W/GuCNSutDM0Q42L+LC6fe1yTdphU/5oNMJj3/YYJicnMT8/v08JlMJ/GcDdjLH3+FYCgwwCEyRQu4hrSnJscvh9DHl+q5I7QmFSvnjN2NgYJicnsbGxgedefg6T908C2Cf8HwbwLIBXfLeVdgAEkTmuwdpUSj50YLkNbdvmM4VUvHYexe80rJePbFlYXAB+bE/4fwjAIg+QsTOYLKClTy0N40aJXiLGQMRMmia6vtPzkXETux2udXPOsbq8uqcESoIJf4BcQASRLXVCJ+dYgK+fg3Q9iNlF4Q8AjLHC8t9PMOEPkAuIILJE95sTvn7NLscMIJGcFFysx1RUOwCJZcYY7QD6RE6Tu+v0rS99PqVWlQUT8rk/vix/IE/XVVOb2rRZdP/MnpjFI5ceAQr3z8MolECQQ2GkABJQ+W/7JrxS0Zd+jPEEzFiWf6ynecqPE0k5F1zrloX/wuJClQq6iMBKgBRAIrqUg0/kSdP8EQVS7Mcr6/Bl+Oge1dDFNTW6NsLm5U1Z+KN0+1RK4G4AR33XTTEAohf04flFtlasyXUpLOKmB6m1bXfdoyZS0aYNE0cmcPLMSYwfHj9wEphzzhljiwh0Eph2AESnkd0OfXEHhSCXvrHZubQpxzeqh+r52nVMHJlQPgYCKJRACOEP0A6A6AGiRRhKKIQsP/STL2M9R0l2ObnUZ3rQKtapXl39faG/OwDG3lijUjFxdSdyg4iQVJaYb2Eq/35BCFRWpI+6Ugir3NNLu1x/CPq5A2DsjQA+AeBZMLYIMYeWcxxfXsX05U2cO3PS+49hE+nIxcXhio8smqYngMai7nHPPmIEhB/6qQCKhyY9iyJ9CmBsEevfuyf851bWcenEbPED6sQ+fATVTJ66GCJ451t4pPyNBNe0Rt2Pv+SCrf8/t/b3jX66gPanTz0MYBmc4+5Tz+wJ/wuLC0CYsxWdxvcPdTQ9ziC3h31V5GCFtm1DF/3jJPDj0tcdQKEEivQpAHh4ae5RACDh74h77nnVAAAgAElEQVT4EDLgYFC0SdDXuTdySeWTSeGOCHUKOEdyULJDp587gIrXdgJ72Ah/3QQdWrqheL/yvZv2he7ZNkAcYeXjV6VC4ivQ3NXDUDYMbf2FpN+Pg97/gwoA1DsAW19vX3+eMNWiCtmPunsyHedYYzwky1/G5zP2+8TSXUv0o/DO7Bf+H1q69MjDVQAYKHYCW09uH/ha3aLP1VVRUeeiqftOn3E9XZuLJW3rhupiEFXVxr7Py1zopwtIEv4AFsEYLiwu4NKJWcytrOP48iqg2f34ziSJNZnFhWQi+Ie8yHK7d9/Bd1/fJ/pNPxVA8dCku1EJ/9LPtfXkNp668Q5cnJnB1MXP4tDurrYAnd/bxv+tujbUglQ9+CuUUPFNjPbYuHxiBX9jKWFdPV0wArqwg+ky/XQBcX4FjL0HwCt7wr+a6Izh6fl5HNrdxfVxP+cAKsEh/h+TugBr7tACjzdudUkNuY1Dbu3pK/1UAEChBHQwZiz8bX3HMSz8HMrxRQ4B9RRpqame0pkqC6sr5KgMQ9JfBdAzXHcWuQl8Gyau7hSntVVpu5xj4trI+FEeOheIyXt9oUs7ACIOfY0B9AoXK7kL/t06Jq7u4OTpc+pgfflIj5OnzwV5qF8frX8Tcm1XTIamCEkBdIRcBbrPRy+I97dzeBybd04fzNgSnue0eec0Pc/JEznOLZkutLFrkAsoICq3Td8sjGCLskzbBbDv7Ib4ML8Qj/To2/j0CRob/5ACCIiPHO5cJ32UbCdJCVSKgJ7nFIZc5xoRDlIAEbH54ezcsXkYXCtKJVAJf8DueU45keu4kuAfLqQAIiAKSZPFJv/odQpMD5GpXFxe21z6/EWOL6+S+4cgPEBB4ICIz5PJ2Z1jikk8o+7ha9b3L/2Az9KlR4we5eFSv49f4yKIrkE7gIC4HnLKQQDp2h5NiUnCv7L4VYHhtjuBZPdIEIkhBRCJ0ELFp+vF1xNEXVxZ1XW3vu8mTF/ePBjwFZTA9OVNq8Ngcpuq1wQxVEgBZEQbAR561yC6sEzcJW0F686RCZw7c1J9ErhUArbCX2xbqkcxAHns8AgCGFgMIMVWvytP5BRRPblUFF51TzZVCTnXfv7cR19S/mYDAIAxJ+GfA1U8IsZvDpg+BTXXg4ZEWAalAGRcJ7ztIxna1pcK1X02Ca0u328fsR0DUgbDYlAKwNep3LaKoyuLq00721q24rP5QxDa+raNeYSizT2SMug/g1IAKmIGZ+XFpPuhmT7Q9jHPXQ/O5tb+3NpD5MHgFIAcwAyZn2+aW+5D+KsCtKnIoQ1N9E3h1mEzHrn8FjIRh0FmAckT3Gf6ZNMPboT+wRjf92KrIG0yhVKSi5ALmZHUlPRgU2+sH8nJZVyGwuB2ACK+F15qQedj8YgWYJ8tx9RjpcJn/9mW1aQs2rZt4uqO/uQ255i4uuPl/ilmYcegdgAxrdIQwlD3eGlb67+uHNkK092HrXIQy5X/VjGkw1o+kwNc+kqn8H31e/XjPpt3Th88uV2e+J6+vFmc+2iR2kuC357BKABdSmOMeiraPEHT5+LUbbVDCloTpaJD7Le2wWXXNoTC1/24lBfrQJz44z6A8PgO6XEf9OM+r7H1xBZwV/h6BqMA6ixbXxaY6SL2VU8bn2kOFrbNjkX1d9t6c1ECgF/r1SVmE5REP+5DNDOoGIBqkfnMda5zZ7gK6zq3VU6ZP6GIvUuLRcgMtCxdIaUSqJ7mujT3qFb4+0hoyLIPGkgRvxiMAlAtNNESbBv8rHvkQ9MOo6k+sV119+FCSmGYgyCOyeAPVgk7gQqV5e/LsOlSP6c6QT8YBQCo/eg633rd5NMJfZ1fvWlA254R6LIgTdX2FPXqlPdgdgGKH/e5+9QzVr/r0EdSjtWgFIAtTdZ9jPqbUFmUNu3NQVDU9XPq9vmuP9Y9ZmcUaH7c59jGhvWP++gIlcUUihx2g4wPRPsufWop6Y26DLRLtkzbgHYuiyb1oaA+WOa5jKXux3207ztQ16/Z9EOJ6Rw4e/Zs8Kg47QAyxkVYtLUqUlskFW1OrYao3wc5ZF6lYOLaSPvjPk/deAcuzszs/biPCybCPwdru2qHCbHmB+0AImA78Vzz5H3Rd+EUc3eRUujkNI4TV3eUP+6z9cQWwDlu/fa3OB8Ca+pjVSZdir6xnQu0A/BICgvAtU4xDa6ujFD3lIu11AdCH65rc1I75hjvHJnQunemTt0cTPjL16TYEZjWlSLeOJiDYLFOPfqmKRupi/eUIyF3BaHGSXVCunrdBVK7KkP3Uxcy+AazAwD8ZmDoDpXZ1N10relx/lD0aSegu5ccTwXbIt9Xbn5mE0wteXEcbbPdmuazj3MaJmWkzi4UGcwOQDVpdCcOXQamTa5/27oJM0RLOeb23yd93PWZ7L5sDS7Tel0+q/D1CBnben0yGAVQYRMwCllP3fds2xBaKKROyfRFbDdJjMVc96iQ3NClKXd1bslGpS+3FJ0EjoTsamnzEC2fWzqXbWhI/zWQv3Cpo+sKTDVPXU+xp0CezznPJdO+E8fAZr3WPTImRb8MbgcgIroDclg0ObRBpJrcubXLFlEw5ix8VHSxzSKu8Ym233ElVF2ma4hcQInIIRiYQxtE+iD8RVIcJmtTp4+daSpC+udTII+luDZUWViua5nSQBPis/PbuE76JnhzhPo4DLlb+KFo+xyiVHNx0DEAmZwmYg4pmHWT0vSwWi7oskhkqy5E1k4XLXgbujIHXBDHT47B+CB1n9EOQMKnG8Z1cHMWAnXBvNxcWDJ1u7LQWVQ25Np/MqmFV2hU49CVsTGFdgACvlLqKmuo65PF1arLVTCY3I/vMbPpi67Pl1g0PQKjbT/mckgrBoPeAciD7EtwpUqdDFVvrgLdFyHvzyYQ3CWjIeWcUOXf1z0OQ7xedxZB9b22bUx9kt+EQSuAPpKLsO6SMMupnV3qt9TYnC2QFUJIKz+XNWjCYBWAzr/nawGmmARdmni5EUrw9uERIH2YVyb3YDsHxJ2Halch15njGFMMQCLHQTIlp7bn1JamxS9mNPmmC26APuAr08o27tX1HRspgIxouy3tg6XmC9PURLG/Qz9Oo4vEnFM+5n7sg3ddZzC/CPbQQw9xIO7gxk4DzUEB5LB42vZDaN9wDn1kQqr55OPxF7qMPpOzLarvm2QINrVbl3Sia9PSXUv0i2BdpEuHYnwHw7p070BcYZxzemFuh7lM+srmwW0m16s+lw881n3Pte9S9vlgg8Ch8DGYrn5FF8sp1OTLydqt6xfdobAc2h2D2IfiTJCfs+MDG2VR505K3Te+oR1AQFwtvhDCJ6X1mZNl2UW61Hc573JM8eGCksuzPVUca8wHFQOINTF9BKJCkIsgCT0OPu6z60KsIkUQN3Y2VagDXbYpvE0uIts6zp49GzwGMCgXUGi3hHjARA4o5SJ8cyDkOFA/v0aovtDN5xRptKlPcbfNIkstGwbnAgqd4aFyd5BQigP1cxxi9XMOLhLfCRK5MagdQBcszpABSPHUoi83SZtyfO4EbNuR6nlNofHtuw5VdlW+6SE9uS1ykDhV4L6N+0l+RlEKBrcD8Ino8hH/5UyMwy66slWLmALT/nBRgr6sbNMxtynfxkDIbd2J86tprpELKAKhMmtC1BVSOLmUrVJwOotM9V2x7qbrXWmqP8aJ3z4ing9IYeR0cax0sUDx81yMj8EoAKKgbmLalCGXp/rM5pqKEC418X9bBdYVYggUVyvWtW9dv+f7UJvPswi5CP6KQcUAYpDbAKtwOSxmuw03zXtuEswuNJ3orK7pU7ZWU9/6Sml0wbZun3XFjg1UfS23I1cDY1DnAKq/c8+1z7l9vrMicsiyyHVxtsFG6KpccyYWvc6lp/Ld26RT+kCV7NCmjhTnKehZQIRXchL+IcpzrberVr8Om/txEZA2Al7ecanwbQSEOhgWi5jzkVxAnnF1I8SaoLm2L9ThsDrXSIj6QtG0W3LJ2GlzXkUMdOq+K7rdVCnIvvo+hsCM7R7cemILuCt8PYPbAcQSZDaMRiPoXHGcc+xc3WnbLADtt8CxAo0u9ZhYfWLZXRH8gPqQYdv8c5WP2rZPVArEZ2qnKXVlyQqqb7u9ttAOIDGj0Qjnz5/HbVdvw8LiAhh7ze3HOcfq8io2L2/i5JmTmDgykbClems6VODWFNOAm2vmj3iPoYKKch0yba3/JgtWFswuqcIm7Qnh6zdtR5eEfywDZTA7gJSHjuoYGxvDbcduw/rKOlaXV/d2ApXwX19Zx/Sd0xg/PO6tTtd87q74Vn3vVlIL/7rvyW3ymbJoS12/hxK+MWILOeXt+4Z2AIFpWtiMMSwsLgAA1lfWAQALiwt7wn/2xOyBnUEqZOsqR+HfFZqsU5u+1Z11EMs0FWC24ype7yuF1GVuNcV6XMuvM3r6oBRIAQRGtyjlCSgqgUoR5CT8gTyDpzHz+H37rU2Coi6Koa1yDq3c6+aRTpmlRDevQs23mPc8GBdQSkwGVNwJVLQR/j63rTkswj6jc2OImTa63PsQrs2m8uqC603fzyGRoE1/HRqNsPX4i+ryOMeh0ci57IqYOwtSAAEQF6ap26Ty+YuIMQGXNvgk1xiKiOvCaaMs2yzWJgtSd7hK3Dno2h5SiNS5Vmz70jRryIZQh8pufd9N+MD587hvbQ0o1+XeGHKO+9bW8N2fvIBb33eTlzpjQC4gz5ikpMnXiAHfyu1TvQbcdwKuC0GVOqdrex9wvac2i1TOoxff81V3V3zVsitMxEeiQhvE+ncOj+P5yUkc29jA0eWjuLC4ADCGqfsncXx5FXMbG7h0Yhaf+8gXMHXqZuv2pFhbtANIjEr4V+6g2ROzB7KDfNCUEy37X5smZk4CJvYiaqtkY+C6e4vZRpfsJ9U8tnVfWcEYnn38Xlw6MYu5lXUcX14FOC+E/8o6Lp2YxYXFBUydutlpVyneF+0ABsLo2giblzcPBHzFmMDm5U2Mro2czwHYugh8ZGDEIKXiUe2O5M9NFKfJd+VsHlPLXpUiatJnuSh0E2EePdefscLyBzC3so65cpdeCX8w1roNMdcTKYDETByZwMkzJzF+ePyAm6dSAqNrI7z00ZcA9NMF44MU/WLi7mtCdq817c5U9di6GZquD51N1ZbUu66tJ7dxYXFhT/gDwIXFBWw9ub2328q5/0RIAWRAnWXPGMPEkYlk2/HQh6DqmLi6g53D44Aq/lFmXFwfH99rn0ioA1tNn9f59Nse2hIDvwBw6/tuwucef1HfP7u7ynGzEVCphZlpimzMNlY+f5Hjy6t7O4MuQQqAOIBqsaUQ/idPn8PmndN7W+s9OMfdp57B7dvb+OA99+wpgZDoApW2WSxNilTn+pFfHxqNcPL0Ofzq+E14en7+QP/ct7aG27e38RM/+SBUT5KqE+ypBGtdO1QkUU4Kn3/1GnhtJ9CGvXuih8ERscnFxbRzeBybd07vC7YB2FuAxzY28PzkJK6PFTZM6HaLbhrTw1h1n9lY4Cquj43hV8dvwrGNDbz/5ef29c99a2t7/fO5j3yhtuxcxrsJW9dYEBTCv4oJVIHhu089A3DudEYj1LmOOkgBDADTCZWVMJAWlpxxcXFmZp/lG1ow2Bx28vG5eI1SKDCGp+fncXFmBnMr63tKoBL+cv/o6IqvusIlTlJ9ry0T10aYvry5T/gD2DdX3zF6Cbd++1ucyreZX74gFxCxRwo/fy2ajAtT4ZYrTSdldZ8fcHmUSuDou44W/YN+9I+ONq4pH3N758gEzp05qY5LlXN14toIO1JMTzy8p2t7qnVHO4CBkJVgt0FQAhUphVtoq1MUEnX57nso+ufZx+890D+qeEJdgLWP+Mix3zkyoZ97jB0Q/hVNwf/Y+f8VpABKurYVdqEp+JilICjdPiLiUfxY2ByKU2Uk1QV3TdB+pwyIi+yLmWC//9z2TAgRh1TrbzAuINcDMFkKxRbUWR/Z3Wsp3OYEn3bl4wb0O4GQ92FzEKuuLbbWtzLjRfD5H8hIOaHPSMlZ4DcdjMu57W1IdV+DUQCupAjMxCbH+9p6/EVlQPPp+XkAaFQC3tvjYYG2OYlbJ/wvzszgWSEjBcBevOSpG+/oVCxAd25B/FuXktu2niEyKAWgO9TUV6tCh2ry57Ag9uW57+7i9u3tgwFNQQncvr2NQ7u7Uc4ByL75JmSr1fdjGOT+mRL656kb7wBOANOXN3HrmT+Nz5WnyPtCnTuubXrt0BiMAqgLwvR9e9lFro+PF4e8xsaUGRfPPn4vfvrsCweEf+iTyy6WZ53V6mrN6vpn64mtPSWQo/BvulfbMXNds7kFwOU5EqtdFAQWMAnwqd6ri+DnqFRymfQVur67Pq5+DMTUA1PYenJbKfwbs2YS4cNqlRH750AfMnZA+Oc27m1pO75N6bjV/yHmkSrom8IrwXw+Zjhnlj61tHejTRrWRzA4B5dKl/CxmFXZNyFwDf7mpJBiYXLmQeUua6JNXzYFmkNi4n6tXi/dtRQ8eDPIHYDL8eyQdRB6t5zNd3Pq81R53V1CNW6h3LE+TmaHakebud+WwcQAbAidRhi6ji6i8o0D7pkz1L/dxkRg2ygJH5lDsdhrIz0Mjug7TVZy04JtyrMPQdOpTlVbuiB4QhBrHEzoils25qEwUgBEMnISDkQ3sTlIRxyEXECRkTNVqveGhu8gXsw+HKo1b0PM8TBx7aRIsewCtANIgI9nw3SZPt8vCZf8qNZb6nmXun4VpAASUk1M1a6gr7S1/FMLWBNLk8ijH1TZNbnMnxz6ByAFkBWpJ2dIbBRcH/qhD/fQhlzuX5Vzn5rQh8xsGORBsByxEY5d82OGPLFZV6fvPtLdRw7uhdzo0vxMgcnjMGIcBKMgcCbY5rt3IYjsKhR93E8s4d/0GUHkDLmAMsJVaOWwlZSJLfxD3n9ufdsFqM+6Ae0AMqNNdlBOp2Ft3SJt22r6fVvXEAkyIgSqR16kWK+0A8iYthMipfCKKfxtsW1bri623CHlaUbKtHAKAncE3wenfOPSvpSC1eRgEAmw9pDydGPriS2cPXuWngZKFLgupFwzVFILBtniyik1r2vQLqm7UAygQ9g+KjeG8A/h6kmV5loXQyGloIYEfxjoYXBEa0Sh5SLA6qxhW0vZ5fn+uUDC3x3qu4K2CR2hoB1AB/Fh2ds+S72q1xaXjJs2WTo2uwyf1xH029o6chX+AAWBe4GcSmbyDH3bH1yxxUc+f9PPCbrW63LqmoTaa+S2Q+syOncnBYEzJFchYPKTiKrfHNVd3yao1+a7NjuAmJlNuY47kSdyQkHOOwBSABJ1Pu+cMdl+6zJfVHnIus/qynZ9Zk9TW4m8yH0tpERnNOgUQep5TgpAgTxYXZnwLo8qbrJSTKxym59FlP+p6vCZitmVses61M95+/p1DFYB7FzdgS7+wTnHaDRSCqIcBs0E8bcGgHbtbiuQU2UL0e4iHDk9diQ1XT47MjgFsPXEFnau7uCxb3sMz5x65oAS4JxjbW0N58+fx2g00pbRVEcMVItQ5YrR/TBGTieEQ1j9qnuU6xmy4GpLV4VeKmRjrC4GF4vBKQAA+MJHvoDJyUlsbGxgdXkVnHNsPbGFFx9/EWtra9jY2MDk5CTGxtRZsrkIDV9pmbncD9BeqNSd6KUdAeEbl3Tquu/F3k0MKg1U7NjK0t/Y2MDMzAzm5+cPvGZMn4Vl+/yYXIVM28nm6obxUb6ta4ksVv/kOq9j4DqfTGNmlAbqGbHjGWOYn5/HzMwMNjY28NhjjxkLf6A5yh/TzeJKFwSir6ysLtxrF+liv6Zss6k8oEdBeEa1taqUgIiJ8BfLFP9XfQYM10qiTB4iV2ImNlTk+GypwSgAFZUbSGRtbU2bHaSiay6fJny12/cEV6WNEnnQtfHwlR3Xtg0VKftvsApAjgE8+OCDe+4gWyXQJ2JNxhx/+pEYFr6NNJt4mHzYMhWDCQI/9NBDezeqCgAzxrTv61AFFrtq+bcRrCYB8TZB2NTWGmFGV+e+T8Rn+7SVDUt3LVEQOAS7u7vY3t4+IOTFwPD29jZ2d3e1ZeSeTmlLiLMCchlNzx/Sfd9nm5rqJAgZ14OMXdjlDnIHAACj0QhjY2NKC59zjt3dXYyPjyvL0j29T/eZD0TLwld5wP4nXuquacK2XVV9bS0k2hHkxZAUatN6lHcCLn1DO4CAjI+Pa907jDGt8AfiC54QwU8TS9rHs/V12RYk/IkuY7N2fD6i3DeDVQC+CW39+HR9qE4k6h4K17bOOlfZkCzGIUBK2QzVAxFTrQVSABpkv7WJRg81iCHKNZl81TU2T/tUoXoej/iZbYBYbJNtTIHoJxNXdwCdO5vz4nMBl4OEJt9RCXeb78SGfhLSgr4KEN2htUowm0xMGyGuKtf2KaBNbrG+ZGcRzUxc3cHJ0+eweec0LiwuAKJrl3McX17F9OVNnDtzEjtHJqzK1s1R1fuqOWc7x33H+pqgHYCC1C6KlFtpV0uk7iFsTfh8qJ3OVUbuif6yc3gcm3dOY25lHceXV1/bCZTCf25lHZt3TuNzH/lC7fysMyqaUpHr5pet+zbmXKUdgARZigVylpApsfqvycKXs5tkdxYphDDEtmABAIwVlj+AuZV1AMCFxYU94X/pxGzx+ZPbB9oqo1ICpu5f3+deYkAKIAFJFkkkXO/NxGWj23abpOP1tb+JEkkJVIqgEv5bkvA3xeSRDTpFUjfnRAXRRlm0hVxACfCVfhmifpMDWb7bULewdAGypgWjE/ypgm1EBAQlUFHFBEIcyqrbDYhxrqb5m3I+kgIQGJKVqBKOIdJM665pstxNMM0CIsE/AEqfv8i+mIADLid75Qw7GwMqtgwiFxAObsH67KKRsd0N+K5XtRNoW6+uHJNdB9FRhIBv5fapXgM4sDMwpS6Lx4ehU5VLLqCEyIGfoQh/W0Iqi9ALgMa0xyiEf+UOunRidi87aOr+SStLXCXwfbhvxV2BbAiFOPVfB+0AiH00TTxdFo3ute3W19TXb3O8Xrb8XQ6fEfkycW2E6cub+4Q/AIAxPHXjHcAJYPryJr549gVcr3nEi4hO2JvMR5eUz1Du2CZoB0DsQxaWOt+579PQJgrDx6KgnUD/2DkygXNnTh48BAbsKYFzZ07izaen933UNBfanIdx+V6KbDVSAEQjdScd2z5yweRgTl2d8vddDqGRUug+O0cmDgp/lHORMXzuoy/tE7Ki68XH+MsGTN0Zg5wY7OOg65BzdIeG62MUdC4X0zJsgu+2OdQ+Mo4IM3JbM7bzuc18UCkA07knX0ePg86AIQoH29Q1+XuAmUXUto26cuvaroox5CawCD+kSP3V5fy3iWWFhBRAiUp4VZhOpDqfuen1uSkcU0Gqel91b3XfcwmeyW1tUjiqz3Prc8L9CZk5j2WOSoAUgEDIyWNyclVMAas7QRjTsrGdtKpsBpMtd52iCSEEZB9wakusL4TsR99zPpRAbtqB2iRVhGZwaaA6n3HbQx0614IPF0OTtZ2r8KrrK5s0TmK4iGvINC1T5fP3dWjLFt3hQ9c4m28GpwB8BHh81NUF4da2r+QAl4nFVbcwyF/ff2ytY928Mpk3OoXiOseqtsiCXj7sJX4mvp9ifpMLKAKq3YHv8kO5hdpMyDrBrntteqCGFMGwqJvb8u6gjQBvM69UaaZN9YnfTcHgdgAhToDmZM03LRQbfB68MnVj1aWSyp+r0nVzGouh4VMp+xSOMY0FeX6aXCta/rF3AYPbAaQUEKl99invXRUc1l1nGh9oE0zLKRBH1OPT9Rry+6bWf05poYNTAD5oM8FS5CbL9Ye83hbZkq+7TvU902wkFaqsKyItujG1dam0bYPt9fKORc7qs6kv5jwkBWDB0HLHfd+jq08/tOuKLP88sZl/vs/TtI19iUoh5/lICsCSPiiB1K4gk1RalUVumgZYV59YriruQMrAjVCpzqnKdFE+XZQHgwsChyBEYDkkpovVd5qcSX1NqXF1sQSb3GrT8xxEHFzdPj7xmfYslmkzZ2NDOwAPpB7EmIRclHI2hLx4TE/42mZgVN9LHaTvKrn2l0u7fO4C5RhTjtDTQAdITgvWVKjrrHud0JYPoTXV37VdXE7kMp98Z3bZzE2bHaVpm2I8DZQUwADJZcECbrsnl5z/2KmEQyHnueSjbabnappcoy5pn/Q4aGIfvremXWLqgSnc+r6bAM73uYf2+oRzHBqNtN+33YbnJNiIZlLMaVVCgfi+7vOcIAXQIXxNchdrOPSjJpp8r188s4mTp8/hvrU1QNi1bj2xBXCO+9bW8IHz52uVgG2+da6LNhdSHmhUWfsx26OqT7dG6rLMUhtjpAA6RFNeccjsnBjl1y2G62Nj2LxzGsc2NnDf2hqm7p/cJ/yPbWzg+clJXB/TJ7aFcBEQaah873Vp2TGEa5uT5zlAMYAOYXO6Vby+C4LPaIFwjve//BzmVtZx6cQsLiwu4PjyKuZW1nFxZgZPz88rfxe2wjVlNNfFm5Lc5pCPRyuYjrPtOqy+I2e5qcoSv3f27NngMQA6B9ABbCaxKisht8XqDGN46sY78MrMKzi2so65lXUAMBL+Iib52iHPO/SZVP1Wlylmgu3BL9PyVeuwzk0UG3IBZY5LemKvBRdjhbAXsBH+JnRhx5SK3E5LN6X72gb+m+5NPGfSlCGkixGIn6U+I0AKIHOqyfHFM5v7gp/74BwTV3citiohpc9fRA4Mu2Cy8Al3P3vIRz3okhRsdiMxfPku6cuhIQXQAQ6NRvjA+fO4+9Qz2Hr8xf0fco7jy6s4efpcp5WAzvqSUz3f//JzOLaxgYszM3jowQdxcWZmLzDsqgRUrp8cHk3QBWRLOFYfxRaistVevddE6mylJkgB1JDLQF0fG8Pzk5MHBV0p/OdW1rF55zR2Do+nbagH5AWytyITMD0AACAASURBVNjun9wLAIs+/6fn562UgOkD5nwEFfuKyuKWH60hE7rfVOnErunOLtR9T9VXstI8NBrp5y5jDIy90alhDZACqCGbrZok6I4vr+4T/lVGjE8/eE5sPbGFL559AVMXP7sn/KdO3Vx8KPTN7dvbOLS7qyxDFgpaRVOTKUQcJNWjNOTgqirxoU3srC7V2mZ3qJtn4ucTV3fwgfPn1QYMYwzAMoBPhFAClAbaJYSc94q+Cn/V4j00GuH62Nhrwl+8jnMc2t3F9fGDuyBVMM5G0Ns856WP6FIYddemoG5MXXYCPh7zUNVdtyvaemILU/dPKtOZz/7Yj70OhfB/GMCHACzCs8CmNNAuUVq7ogLok/BvWmjXx8e1PuetJ7aMhL8NQ3f3mKYwioRIA63LmdcJatWuwGdbRFRzsskltu8zxop1DOBYmdpcZroFFf4A7QCCIk5CL9Ziz3cANpa5KJCarKymekw/09HHncCh0Qhvfuht++aVyW6rwpfgTRFg1rXDVAm6uBGnHpgCOMfdp57Zt74RUPgDpACC4d1tIAj/apson4rtuhKwVQAqYe3qztGVESM9MDeqrLPnJycPnrEoM7GmLn4WH7znHrz59DQA/2cncg7CBx1rznH2scfEd14XSvgXhRPZoJ3cCuFfnYq9dGIWcyvre4HhrmMqfOueAVP3/T4J6hBMPTBVm3V239oa5lbW9567ZJMS6aNtMYh9OEtMc5bPuABYLgPBQaAYQADqJqrLJD60u4vbt7cPPPJg6tTNuMCL8qYvb2Li2gg7RybcGt0CX9aaLrtClWlS/W2aflh3/F6uV/bhmtxHjod8bNnzn5+6GU+jOG1duSOeffzePffEpROzeHZxAVOSXPLpVpPL9dWvdbtGU4PDd8bT1hP7H2pY7eiX5h79EIoYAMAYxQDaoHMByUIFSPvYZR1VBsye8BfL5jyJ8Hdxt8i+1FCugqZ26YJ5urJNhVuXFQAg+a8ff7E25hQiyBqaOl++qTuxjZxQzj3FQw7BGJZmvy94FtCgXEB1Akm2KFNQd0rw+vi4OiAHAIwlsfybkC10Gz+xfNhId+jI5WSlboxV7WvaQXRNANowdepmPPv4vfvee+rGO5LHmkzHr+n7TXOn+lz8Z1OPqjwVt377WzB9eXNP+G89uV3tCjiARRTC/24AR50qrmEwO4ClTy1xoF7rp7L868gx57rCpr9crH7TLCDfNMUOXHPCu0SVlVLlp1dcnJkplIKkBFyzp3yvFZv8fdd21NVjul7l+TRxdac4yS/0695PQhYxgKPg/IpxIw0ZdAygC35b3/7G1NgsNBP/fQh0O4E6yzOkKzEGB+5PcdL8+PIqjq2s4+jy0QNZZ7rzGRUmfeeiDOTc/zbft7mmbb3yd2p38IWV7l34AwPfAYReqD6smxA7gBCZGyZ9mMLf7wt5rqjmjm1KairqgqBVQFKZYiwohQMJCZ7WUpMikdHNZdU4qSxwm92caexAdU1TXSpi/Cj84BSATMhFGloBVHXoJpvOugmZuucSGPZRbgzqhH5TX9eVKe8gfLZPrEO8RkV1TdM5gPvW1nD79jY+eM89tYfBfGHS5qZ53maH4dJGH2uMFIBHdAoAyHcH4CP1zGfGjQ2+dy45KQDA7rCYixK3Ge+2fm6VD1vOOtuHwUngUJjeq+l3XAW1bm36XGMxFMCgsoB0hBKMviw63WubNqQOGNchZvzE9vfbYJt9Ypr5lDLDSJelJWed7YOxJMJfRjdfbGiTRSZ+VzeuPtoYEtoBIH83kAqf1mEI2riuXMqLgexKsLE0ZQ743RV1mO4KTCx73fs59KsLdTtbH8ozZr/o2ko7gB4Q0r3UNMlzFf5dpem+TD5XWY8VdeMlWpo2Vqvcpj4If0B9LqSiTZxJZ7G77BRs603BoNNA+4CrXzokbQRM7FiFC03+dBHxOtGClxWBTV81ZbjI9fsWLro5FwtVn9q4Ol12p23us2mNquZELAbtAoo1eXMXaKEwDQSbuDW6bK3K2KQsNil1ldBo+r7P2FSqcXHZAblkabkgB4VNjBrVmJ09eza4C4h2ABa4WFNDFf4+yD2A7SpQ64LGKrdQnbtD937IYHqTYk+trHUxlBi7S3l3Uv1fp4jEHUDsuT7YGIDrJDX1BYb0GXaFttvxpiyL1FSuB59tk90Bqs9SIccfch0T8X/V56F87qI7SjYOXGM2oaEdgCXi4Kryf7vgw45Jk5XsulD6TBeCtCmztXLLFKubw7p4RQ47JWCgCsC142Xhr4IEmJ6mLBfX1NHU+F7QbQOOPn38FT6Dom2I5cc3bYNrADeXeTzIIHCIBUK0p2kR1X2u2okRbqh2tiKphX/qsVUpgDYyQXc/MYLAg4sBuEweOQebCIOJu8gkm8J0jGks9ejGIqXw9eG7DxG/8CH8bc6E+GSQLiAZVapWaitjqMTcQtMYq6nWQh/XQ5Mbt8md56NPTAL8W09sAXc5FW/FoBRA0yEdFbaamILA4emLMMqVIZ3B0H1uIuhtExhs+pF2AJ5p4/qR/zb9DuGXPgmiHKjzZfe1r0UL32ZN111fpyxyTVyoGIwCUGE6qERacl5AfcJWMHaNJoOurTxQ5f/bKNYU/T64ILANJHjS42MM+irQfCH7vfs4721cOiaYBMfFOIoYwFbVk2qODnoHAJDPPld8CqE+CjTfdLGPTN0rsvD3rQya2ia3M6e+HrwCqIMUQ3xyWhx9p+t93Sau1/Rd1/z+HB/hUcegFAAJ9DwJfQI490BcG5pSFnWuia6eunYh1HN/XMrPLbV2MAqAhH9+uCyc0IumS4KP5nRYTAO/NuQ2tygITETH5ESnTT62SX2213bh1HdTP7pa/7bk3k9tsJmrLqSeZ4PZARBpcbGUmp5JE5K+JAe4BiBtlWbKnZNunrR1G6r+1l1n6lbLDdoBEMGJ6SPNcZH5xNZiFPswtbUZEp/35nIA1KbMnBjk00CJeNg+mC2UNelaZpu2dCGe0IU2muDz8RW+nkfV9iFx9JOQROdos/B8H0SShZuLsGtztF/37KmcBG5ObWlDbgcGfZwjiAEpgA6Qyg9uiw/h7xP5aL6P/qMspG5g04++n0DrOs9SrHOKAXSALvhuc7PATMsWg3iqzCPbMuvGyiSjxKStvsjJ4g1dXh0mY9L0DJ+uxgpIARDOVAItZ+Evtq/ulGbddXUC3eZ9GZXy0V1jU64pKR+30TTebdsWahflMk9D7xzaQApgoJgIxdjEcK+4YLOAXc4cyGmaumfI9IkY99Xm4W66suqyqmLt3HxCCqDH6KxzWfCIlrxpINbWlRHa4qtrQ8wF6EPJ5CwwukTsA29dPDtCaaA9QSVAQi4AF4spFZUlHdKdYnq/OfULsR8XI8V2l2FzuHHpriX6UXgiH0x9/rlaQb7TTFVl53rvOZOiz1KNU27zhBRAh5H9xSmtS9m6Tt2GUNa+vHBdTuaqxixnQrczxVwxdbk1XdcmEcBnEoUr5ALqEOQ2SENdMDZ1oFZXf+p2dRETAa9y4di4hmxiPeQCIvbo2mLucl64Tb2pxqXpgFvX5ktsql2NqfVffUf8v66Pm7LrfD66og10EjgDmvyCXVrMsoWU8xmBLtPmERWhyTXQrVOUtieBTV1iqnWtc02mmuO0A8gInT+wywLQ92nTrvjNh0xugr9CbldTZpjJie4612DTPM1hLlMMIAIqK0NOCRMJma7ogq11aeuSMO0Hn5Zlbn1sQk5WftfRzS0Zl/5uKstU6FMMoOPURfhNfI8hFrurxeGS6WLS/kOjEaAxQrYefxFfPLN5oOw2hLS6Qu5QVMK/bV2prU9fuD6eQRXM9TG/mj5PnfkjQgrAAB8Dpto2diWfXvZ9uhyxV33n0GiED5w/j7tPPXNQCXCO+9bW8IHz5/eUgK9FE+Isg08hEouutLMJX7501dw2KU+1LkxcuTn0P7mADPEdqDWdLDr3UV05srUYynUiY7Lb2XdNKeSPbWzg4swMnn38Xmw9uX3g/afn5wGm3g23ceXknN5pQhfdWL6pxkm3Q6qoc7nK14rENtLE+ugHYTLEJmugTvA2TayqnrqsAdNFb5OV03SNaaCsrnyxrKfn5wEAxzY2cHT5KJ668Y4owr+JLghUm7nYd0zWWJPwt4kzievE1revgtJAM8XEklYNfBsXgk27bJSLaqKaLhITTLe4+95n7DUlsLKOOawDgFb4+xB4bfywvncGruWJ91AZCjm0q019gJmRY6r4bQwEcT2p0jZVckD+v4sMRgE0WcG6QVRp+rrrxfrECWJq8ava7BsfFospRnWUSuDYxsbeW3WWf1va9GuOOw3fbbKZt3VltN1tNl1jsxuVrXzxuyojqM4F2ycGFwSum0xy9oqJoGxyl5hamimsidCT27js0ucvct/a2oHAcJ2f1tfOoEuEyjaqsDGUmspyoS7rSdz51O2AXHd6dd/zkRDis7w2DE4B2NBkfTQFluSydO+3ySDxsciS+i6lgO9DDz6IizMzOLaxgfvW1jB1/2Rt5pHPxdMFv7+Mr7RT0wwYl7hT0zW6OuuMIt17ppa9Dap4nPy5S5k5zDdSAAqaJpKID+szJ5dEVDTZPk/Pz+8pgePLq9pzAjJds+DbkEqAiELZZ3+r7qcuztXUPvnvEK67Lgd/K0gBKLCdeKrAkak//4tnNgHO1Z9zXhyUaqBrgq9q76HdXdy+vX0w4MsYnn38Xlw6MYvpy5uYuHawD8TFp7MGu9YvMjphJrs9QsYRdG2oe08s01ZIut6LSnnoYgQ2sTzbenXvqdqRA4NVAG0HwceWsDoIpfJ3iwehKiVgEuTqEtfHx/HBe+45EPCdemAKYAwXFhdw7sxJ7ByZ2Pc9k/71aSGHcrE0YboLtWmDaX2m5bi4LE36QeXXFw2zuh1DU5t15btiE4zOjcEcBHvooYc4EPfwTOME0x140ryvsmCqeip8xBBCxCFchEtdOU2Wp880Rl9l+Zx7LmXp5k9T+dV3fCkQk/pM63G9J9syqjmgmws2fWM6n2IcBBvUDsCnUGiqBzCYlJK/u9oJ1Al/07p9X9s1cgwO+/Ibx8D37kOF6MqybZNLXWI5qrqb2mET0M4lyNvEYBSATWDXZz3i+0okJXD2sccaT8G28cv6xqf1X11rMk422/3ckBVBW/dQbEziBa7lmljhbcrWuZFsBLbrbqXuvVQMRgHEEP51ZdduoYXTsBWy8Fe5fUIIurZ9Y+OXbSrHNb2ua7jcq43lKtZjUq6clqzyn4sC1HVnoxO+8t+6/rFx/5i0Q1e+SR226a65MBgFkANaJWB4ECpEe1wsILkMkRButjqlqrOwct8FVPjKevGBauxE37dJm0T3SvV98X+TNtQlV7jumMTyTXFdDzkLfJnBPAoiJk2Wyb7PFU/EPL68imMrxTNxXB+JYCKI66zIWAE9U+rao7I+XSzqVAvXR92m966zpHWJBaKBoBtXnXtVfF8MoOp2FLkr7bo2dknoi9AOwAIfbpd92QD3Tx4I+G49uY0Liwu4dGJ2X2BY/m4oZAvG9X5DttXXIrQNAoYidZBVFGx1isTEHSP+L+8G5PJV7TCxoF36S5X2mbPCidU22gFY4NtKmbg2wjtGLykPQl1YXAAAvOPyJg7t7uL6+Pi+dgDt/Y4hrpEFgE9kQRLLDZIrNvNRJ/xEC71NO1R/113flFZZp2Rc29h2l6tSaqHmytYTW8BdQYreB+0AHKizqmwsyp0jEzh35qTazSMchBKFf1P9QPoMEd9+UNX9+nSZ+CgvFbqgu86F4yq8fGQp2SoK03piW/Jd2UWYQAqgBbKwN5kM8jU7RyYwdepm9XcZKz73PMlymrSuQbmQKYeh+8d3BpfcL3XKzVVBp1SQovvGdLegKsMnIQydFJALqCVNLgMbl4JppkUM14+OEMFi0z4Kte2OpRBDug18uQV9E3q+tnXZ2YyJfC8hYzexoB2ABSZuH5vPQ2zBQ04eVdkxJmtIwRXLiquzYEPgq666OW869q6uH59zvk5B2rSp6xa/zOCeBVSHyhpok6LWNOna+FVNJrIvi9NmgbmW32TlyXX4DNbW9WWoYHaugiSXtoVMIqgLPtuWF3IXsHTXEj0LyBdNFktd8LaNleu6a2j6ji5LQmelt82gqP539cHWtadOGariLHLKYV09bclBGMbCtF9zwKaNqhiCj3Htw9wYjAKwxcciCL2QTPKlfWZbNP1ti60/NUbgMoabJlfBkVt2iy/jKVWw35YY8QWZwSsA05TOWNRZ2a7B5BBCrU2ZdVaYyU7MdwZGUzmqOdIla1mFbt7HVk668TYZ37YxgtwUcYr2DCYLqEnY5IIvK0CMXeQ20St0bqCUgsm0v3y6ElKhi3N05Z5M21m39sV5l/t6CcFggsBLn1rau1HXoG5sfOU3+wyWpkx9Cxm0a5Pe2jbYHjI91KTeilBJBCZtULmfTMbKZlfcdH+5CX/6QZhAdEH4A2EPO8VuSxN1bW17H22D4FUZutdtd5cx00Pr6m6jBF2oE/5ty61zJVbI/Z6L+ysmg3EB5UpdqqiPctsgL9CcrCMbbHd8PgLSXegrE+EX8j5cM8pSump8p1arlEDM+xqkCygnxMnsGqwK7TNPkZ0QYjHYuBZs6+5KoBFIM54xODQa4frYmPonVDnHF8++gOvj485zy3e/6cqr3o9xDoB2AAlRWdZyUApI75sMtUtpqrPt+Qu5z5rKbFOnzr9sQg4COfUca8uh0QgfOH8ez09OFo9VFwXq/ZO4+9QzuH17Gz/xkw9ix6LckLvgxl1XhKeB0g4gICqhEGL72rQD8GGdxhYQIeqzsfxSuxiAOMFXsZ5OKwHph5X2nrDLOd7/8nOYW1nfe3/q1M0Hvm4ahwjhWtURIwhMO4CIxPKlm05m2wWfQgkAcayuXJDHrq69uv4xGSfVbieVwvNSr/C72sc2NgAUv6Z339oa5mSloGmHSMhdb05zkBRAIFIOss2ishGyKvdUDFJYpjks0jaBZts5kDILyduuR1IClSKQhb/O1ecalO4yg0wDrYi1xQ5ZR90pZls/dF1ZIjEWRYj0vL4vZhnT8fedhqlqR9Pc8pYKKyiBCpXlX9VXlwbbhq6kmg9OAYiD3pVBEqk7LdsWkz4ZmhDtAi5zIeZOrkm4e21DGQsQef/Lz+39rnYsurJOBucCEre7dZrf1iqKNeChFVeKjB8dPnKtU+aMx6TOwraNI8RSDt7nmBAIvnRiFhcWF3B8eRVzK+vACeCpG+/QxgD6Pj90DG4HAMR5WFpIYlptqgBjDsrBBFHR+2pzV+4daD5d3fSdWPfqZVerEP7V72pfOjGLuZX1YmdgsRPo0lx3ZXA7gApdnrjqutwI1Sadtag7p9AFqnvyGUdQ9dNoNMLY2BiYZGFOPTAFzjleOPsCxsfHvbRB1yb5vap+l3uvygy1exLL9SFoD+3u4vbtbVycmcGziwvYenJ777MLiwsAgNsvfhaHdndxXRqHprrb9kEO5zx0DHIHEIIuCcU6TO+jL/dri8paHY1GOH/+PNbW1iCfq3nx8RexuryK8+fPYzQaGZfvE5NzICohHGpnHCKp4Pr4OD54zz17Ad/qfqYemNrbCfzETz54QPhX9evutc4VbKq46spPDSmABmKkhg1hq5mKGItubGwMk5OT2NjY2KcEOOdYW1vD+so6JicnMTbmf8PdNHea3DmyYIoxD0MJw+vj4/sOee2rgzHsHJkwapvJZ31Zr4M/CWx6KtTm+qZy6lxNTeXHfhxDHTHcA7mhmwuT909idXkV6yvrmD0xiztuvANra2vY2NjAzMwM5ufnD7iHTMs2uV5EFbdpKtM0L74NIc4m+DjlbkNMdw49DjoT5KwhH7hY/bGtjqZTxL7qSH0i1eReVNeI720/uY2FxQXMnpjF+so6HnvsMWxsbGD2xKyR8Hdpk0kMqw2+51us9oXaoYcU/ql2FIPcAbgGyFwtU511ZZNzn/OWM1drPRR182Dy/kk8Ovfo3utHLj2yJ/xddnc6a16Fbl7bWt519dsi7z587O5inlVx6Q/dPTbtxGQFQzuAwNhOFNeJ5ZpdFMPy8IFuN9OX2IZ4D3Upi5xzPHPqmX3vrS6v7sUEXA9s2aZJttkZ+J5fcqaPj/KbUlt9zrk6BeyyQzYR/k1l+2SwaaCx6Vr6pAvyAq/bweg+y1HB6dJgxXZXwl/0+T/38nNYX1kHACwsLtS6gULODdP4kkjbna6qz3yNbd1OSUxfldvhC5MdiDxH5PMVKmUV+/wFMNAdQCpB7GMy5q5EmrJNVBkgsdLkVDsq2++r4JxjdXn1QMD3jhvv2IsJVDsBm7p9xZ5sBYsPY8XUUg6BKr3VddxNs6d039X1fUgFZQPtACKSu/Bugy/lFiN7w3c5o2sjbF7exOyJWSwsLmC7PIRUKQGcADYvb2J0bYSXPvqSUZmy8G/j904173zv7GyVp0oAp8wwq+qW25TSOzDIIHBFjtkmTeTgShLbkNqCsaFtm+v6fefqDsYPj4MxdqD8yfsntcK/qX2m7gYXBRAjkyWk28e2Llsl6jJfbIP0dQkiMX4ScpAuoKZtWeh6u45syXQB1VjbuAeaPp84MrHn45evZYxhwuAQkm2dofEZxG9bTtN80yVayLsQ1XW271dl26LKiFK1Pea4D2YH8NBDD/G6Do8pzFIvbFu6JOjraGtFmlh3bWmyDl3qTLkD8FFfqJ2Fyc7JNJFB/p6OpvLE92gH4Jmm9DHiILZ9k3Pqp2kOd93rqpxQVluT8Ff5kGVMArAh5rythW5CyLnU1CZVRpFrSrd8fS5rZFAKAGhvNcVsQ9fwnV3jAx+uh7rXqUiRMthEzLbIQtRn3bJbsM7Sl9M9bV2JsTLgdAxOAZjS1mIxnRA5YxOM1H0vxPUuVIu1KbWybSAxNKp5ZeqaMP3cB6JQ812fLpXYx9jI8S1R0boEhXMxGHSQAvCMmJUhLtbUgsMWG8tXFK6xcvlN0blq6twnplacqmzfrqC6+pves/m+D1QC38WFqPrbpv62qNpgsh5MAs6m9caCzgF4pEkIdiV10iS1EEgzYV1S8nz1d53fv62lqKvPNXAdI2AdihTzSly7IdxLurGsjdHc1braRkgBeMRUcHaBHIW/K75cA4Cbq8HWRaP7noqmwKtO6OQ4R1NZ/Db1m/ZdnWA3rScGg0oDjVlflwRkRRvB1AeBImK707CNl7gIEVX+uEmgUtW2UGmUruWnFv5t2mD6PZs5svUEPQ20c8jZA12ib8I/Jiq3QRM+lJNrUDLkeNkqTtN+k9udep3ZxDlU7sEmYt0bKQDPDF0Y5kIqReWycG2EgkqIpPL3tz1XY6KMdBk/FS6K16dLUEedtZ+TgUgxAI90Ufi7Wv6536uLr75Nbr3Jd3XBWhNBKl4nW/5yimvXMBmruvRMl/pywVf6qiu0AyCsqcuGyYW6jAvbFEvTulIJo5QHiWJRpyRT37vr/A91TsIG2gF4JkdhKNOUraSyMJvKS40qJVP+W3fPrnnrtqmCbQOwsc5ahEK1iwHymD+AX2u8bsesmoepdgKkAAZOnbDX5bnL30+9gOsWjnh/bQ4o1WGbh99G+OXQ3ybYCLmYgq/OrSaeA7Apx+Vz3Xf29VuEcwDkAhowsrCvPZSi+DsXYaSz6m0Fs0+alI186KgLO0dbcr43laGjamvdLtnnvdnuJn1BCsAzOQhEFbZ56iK5HybSWdSpsy9sslx8xyVS0LSzkeeM6VmIEPff5N6MvWMMUY8JpAAMsZmEuQhGkbaLSLZQcrpHldCvE0a+BUpsAZ2rVS0jKjdxzuhSO1XKMGTMo64NqutUWUiqctq2JebaopPADbimP3ZhgZqSk7CvwyQFs20gtq6upjLbBj5VPutcMHEf6r5n6qrzeb9t4zDid1X3YPqequzqGjoJnAGu1mNOi9OVrmWc1LVXFwy2IZTP1/TaVFZiakLvAtqMaRvPQA5GImUBGSBbj7lZXynoYh/o0hAr2uzuTL4rCxvTcrrQz7IlbDo/xD5RKThR6frqB7E8H/M4VnwgBLQDsKBNILVrmAbouorO59z0HV/16j7rcp+K+BZkocpr4/7xcU3q8SYF4ICp0PC1dZ24ugPoYjWcF59HpmvuoQpd6quNElDde5uFLJeXWijEpikjp6vo7ienQ3CkACLQRlhOXN3BydPncHx59aAS4BzHl1dx8vS5JEogV2yCsarXrviIEXRVsQJuOxhVbEYWkCGCv6ZWujymJum6Te213XmGhBSABSnyyncOj2PzzmnMrazvVwKl8J9bWcfmndPYOTzeWFYuky40povb1Q0g4yOLpAvWr48TsqrrZUEbco7alG26y3cpWybVoTkKAhuiChbZDLjzwDKGC4sLAIC5lXUAwIXFhT3hf+nEbPE5c8sYyymtMEYbfKaDhkgtTT0GdeSaiaNCFfDXBZzryqhDl76qcueZ7ApSGAGkAAxJujAlJVApAhvh32Sp5Cx42mCSeeOjDtvrTXLJh4bvbBpffeqya3AR6CmUAB0Ei0yrAeYcS3OP7r1cuvRIrfAfqkDx5aogodx9dBY5YB4r0gX8bT0AdWdUVJ/FOAhGO4CuUPr8RY4vrx7YAXRZYPlypfg43DOUeEkOhPT9yzsBU2Xg2xKvu7eU84yCwF1ACPheOjGLhx58EJdOzO4FhrcefzF1C3uHGJDrQoDWF3X3GipQGSoN1OTQnYuiD522GnO+kQsoIk4DKwn/yuLfevxF3Le2hmMbG3vvbz257b/RkVFZZa6umLY+2FSWma2bwvcJ2areNo85cK3fVzkqXJI2dIFe2/J0ddSVQc8CypwYmnri2gjTlzf3C/8ntgDG8PT8PC7OzGDq4mfxxbMvBG9LDFwyKGQqK9XWustB+Iv1m1iaIYR/U52hCCn8mz5LReoUYNoBONJGD77U+gAADitJREFUONlaVxNXd7BzeBxbT24ftBo4x6HdXVwfbz4H0CXaCmOXfm4itFKoCzrKitF13tXV64uUylO+F9e2xDiTINajqot2ABnjGiR0+c7OkYm9QO+B7zPWO+GfG6F9vjJNVmGMeZd6B+QLlzGLJfzr6o4FZQEloM5HmZMbIhdixQAqUpyTULW1i4IoJ+vflVSHslJAO4DI6LaoJOzribkg5eyfUH5ak6yaunkxFCGVgi4p3TaQAkhAk9Afcg566qCY2A7x/+rvVG3rQ8phSEzSV32W2RdIAUSmyb0wZMEvIypK137JtT9NhYvuPELfTym7pPD67o/QWUlVHa7nEXxACiASuVsTqkMzMQ9C1e2G2pQXynXT5rspvx8an354mzrFeIlOmKYI4pveR6pxJQUQCXki5Ga9VYEvVT54rLbmkItuQtvUQl/k2D8p5rWNsg/pSpOVkGudMceVFADRaH3EmpA+TlaK/4cidytcRUzrN0W9Kqtf9TqUgjI5uGfSF7EVKCmABrq42G1pmnSmT01MjS64ntNuy3dbTJReDuOTah21scSb0MVnXA7txT5rUkEngQmrSSc/XTHHYKQq1TbUwnI5kRuKpidaxvKLu/S3rZEhC0zT74d6dIauXbo66059V+/TSeAE6CyVHKyoUNRtjeuESug8ebk+F1Ic4KojZHtUY6P7PCShlIrubxsL2wemKaeiMmjjFgrJoE8C21q+faSpD5ombuh+sdlpuIyn6XfqrNqcdkHiuLiOba607eMQaaKypW/iiqtL+459Cpl2AAb0OQ7gK/AaGtN22lqDNvevOpdgYoGmmDu5Wpw6urjGRMWv+7/pe6khBWBA3SLv2qSV8dH+XPqgzjXgS9H14dBeLuPVByphrsp+Mk1LjXneRmbQLiAdLlZhalJZFbncv4hN0M2GNn2cm9DNrT0Vch/HcjXaoIuL2brfdIojJoNUAKa+UZd0rlSIk0nO1PFdDxDuUFPIPO0cvp+r4M2VnNeciEnAN0cGpQDapDuGFKo+MfFBtkFlsagsItf+akpd7MIYiKTc3neRLo1tRV2b63ajOcyJwSiAus5WbTlV1+c+OXVb51TbaNO6m3KqQxN6V6OqL4fF3wdMzjrExHW9pZoPFATGwUMbprnFuaMKWIaYaE2HdKq/bf2dqu93Kc7R18SBULhkU4nrNQfFaiL8c5IvpABKbCZO6kmmwyQHufo79eEtm21wzDMHqrxsG7rmA+4SNim+qcbAVLjnMkdIAQioovmqa3LR3k253qpTiTGp83+aIiuuGLTZBXY1GNg3clmjuUMKQEAOXOYeB5CFjW4LHMMqssmDbtMOcfdCArb/9HWMc5EjpABK5Gwf8b2cSWEh64gllENnOlW43ksX5k1O5OQTHxqDyQJqQhT6Xcr/NyG1YvCNiavOR/nV37Y+3a7Pl5zoY1/mtB5pByBBroVhoxp7k4C1agdJEE2ZbyoFF1Pp0Q5AQVee90KCxh82Kak+riOGh8nciC1vSAEIpBL2fXAz6RAD0KGsY9n9YuoacmmLPFYk8IeDi8uxac6nnj/kAhJImTvcV1K41ExSMX0FeMWYkXwoiTDDVRHHJFR9ugfLxYIUQEkXF23X2htjcpsqHN99J+8KUlt2fSf13G87vrm4mQetAESh39UFm3ohdIGYiqer8yglXZjDprKibg7IRmYO9z2YH4Vf+tSS8kZzGIQ2kMCJS4i4AtGtddiUyVOnAMRrmp50Sz8K3wLG2BsZY8oO5Jxj5+oOgO4uWLI2ib7Q9pHhsVGdFfLRphRKsJcKgDH2RgCfALAsKwHOOVaXV3Hu9Lk9JUAQRHzaPm8pNSoLvsn6V31+aDQCdJ4YxhgKeRaEXioAAK8AeBbAwxCUwIuPv4jV5VWsr6xj+s5pjB8ez2Yy2dLVdned1NZnn+hKX8qCXY73NK1Fld+/+v/QaIQPnD+P+9bWDiqBQm4tA/hEKCXQy3MAnHPOGFssXz4MFMJ/bW0NGxsbmD0xi4XFBWg8RIRHdq7uYPzwuLKvOecYXRth4sjEgc/oRG0aJq7uYOfwOKBaG5xj4toIO4rxciXEGRifjwgxfdZWnYVf157rY2N4fnISxzY2AABPz88XfV8og2UU8utDKIxa7/Q6CMxe06APV+/NzMxgfn6+s8K/S0Jx5+oOzp0+h+k7pw8o3MoVt3l5EyfPnFQqgZypEy7igu/SeE1c3cHJ0+eweec0Liwu7FcCnOP48iqmL2/i3JmTXpVAzgqgqR6T54Y1todz3Le2hmMbG7g4M4On5+f3XqMQ/osIJKh7uQOoEHYCewrARfjLj4hO6X7pklAZPzyO6Tunsb6yDgB7SqAS/usr65g9MYvxw+OJW+qXrqYW7xwex+ad05grx2tPCZTCf25lHZdOzBY7hIzJ1T2qVQSMFZY/gGMbG3u7AQQW/gDw+qWlpVBlJ0fYAcxV7+3s7OCWW26xUgJXPnNF+bcLUw9M4YavvQFXPnMFUw9MNZZXXVN974avvaFV/TFhjOG2udvw6pdfxfrKOl798qu4be62fcK/q644k3lw5TNXko3X1hNb9vUzhufnbsMbvvwq5lbW8YYvv4rn527bJ/wP7Aw8UK2HriG3WXcf1RhUY6KEMfzKLbfgnk9/Wnz3T4QU/kCPXUCS++dDDz744MNVDCC2G6jJEnT9Cbmu7AZEi7+iy8IfcLMyuzBWAPZZ/BWhhH9FrlZ7HV69AYIbSCD4DqCXWUCy8AewyBjD/Pw8ZmZmsLGxgbW1tdDKFUB4l1EXFg5jDAuLC/ve67Lw7z2MFcJeIKTwBzqkHEMgxQAeevBBXJyZAcosxpALpZc7AOEcwLMAFrlwk4JyuBvAezjn3dt7dgxVMB6lYuZ9nIBdp2a8olhNQ0JhrIJzrn3fM70MAnPOrzDG3gPgFVnACIHhoyT8w6PajQmvwRgjJZATDeMFxkgJ+OUoCmN0v5AvlECVyn53eZ13edVLBQAUSqDmM44AnUnsRyX8VWc0SAlkgt4a3TdepAQ8wvkVlMbqgT59re+PIpCx2ksXEJEH5IrrGNJ47RNI0niFEkhEXEgBEEEplcABV1z5GQO54vKiHC+lhV+OFwn//kAKgCAIYqD0Mg2UIAiCaIYUAEEQxEAhBUAQBDFQSAEQBEEMFFIABEEQA4UUAEEQxEAhBUAQBDFQSAEQBEEMFFIABEEQA4UUAEEQxEAhBUAQBDFQSAEQBEEMFFIABEEQA4UUAEEQxEAhBUAQBDFQSAEQBEEMFFIABEEQA4UUAEEQxEAhBUAQBDFQSAEQBEEMFFIABEEQA4UUAEEQxEAhBUAQBDFQSAEQBEEMlLHUDSDceNvb3sZfffXV1M0YPONvGU/dBALA1q9trXLOj6duR9cgBdBRXn31Vdx7772pmzF4ph6YSt0EAsD3zX7fW1K3oYuQC4ggCGKgkAIgCIIYKKQACIIgBgopAIIgiIFCCoAgCGKgkAIgCIIYKKQACIIgBgopAIIgiIHCOOep20A4wBjbAEBHgQmi4A2c85nUjegadBK4u7zKOf/61I0giBxgjP271G3oIuQCIgiCGCikAAiCIAYKKYDu8mOpG0AQGUHrwQEKAhMEQQwU2gEQBEEMFFIAgWGMHWeM/b+MsecZY99tcP1JxtiPem7D28u0UTDGvp4x9iOa636TMUbPVSeiwRh7kjH229X8bLj2KGPsJcbYm6T3f4oxdl/N997LGPuZ8u8/r1uHjLErtu3vOqQAAsIYez2AfwzgmwC8E8D7GGPvTNkmzvm/45z/rZRtIAiBcwCMfsmLc/4KgJ8D8Ber90pl8A0AfsawjI9zzn/Ivpn9hBRAWO4C8Dzn/Nc55zsAPgbgL5h+mTF2jjH2I4yxTzLGfp0x9q3l+yuMsT8rXfeXSkv/WcbYp8t/f0JRpmgN3cQY+znG2GcYY2cBsLY3TBA2cM5/EcAXLb7yUQB/WXj9LQAucM6vMsbuKtfKZ8r/v1r+srjDZoz9IcbYGmPslxlj39/qRjoKKYCwvBXAC8Lr3yrfA2PsUcbYnzcoYwqFhfPNACrL5WMATpTlTAA4BuBfAfhtAH+Gc/615edKV4/A9wL4Jc75fwPg4wCmDdpDEMFhjJ1mjJ1WfHQBwNcxxm4qX/9lFEoBAH4NwJ8q5/MjAH6woZoPAfi/OOfvBrDtodmdg04Ch0VlUXMA4Jw/YljGT3HOfx/Af2KMfVX53s8C+BHG2CEU2+df5JxfK7fDP8oYuxPA7wH4ow1l/ykA95bt+b8ZY18ybBNBBIVzfkbz/g5j7OMAvpUx9i8A3InCLQQAbwLwTxljfwTFOhtvqOZPAvhL5d8/DuAftm54xyAFEJbfAvA24fUtAF60LOO68DcDAM75q4yxXwCwgMLSryygRQCfB/A1KHZ3Js8Kojxgomt8FMD3oFgPP805H5Xvfz+An+ecfwtj7O0AfsGgrEHPf3IBheWXAfyR0tc4gWK7+nFPZX8MwPsB3A1gtXzvTQC2yh3DXwXw+oYyfhHAXwH+//btGKWhIIrC8H+WIK5B3IWl4BosFCwjuAYrwQVkC0LADYhgFRQEu2AvugQ7YSzmhVQhCSaF3P+rhlcMt3i8c2cuD5KcAHtbqk3apSfgABixaH6gv/9fw/p8jX2mLOYJp9sq7j8xAHaotfYDXNI/0O/ApLU2g41mAMs80K9wHocBM8AYOEvyQr/++V6xxzVwlOQNOAY+/lCPtLEkd8AzcJjkM8nF8HzZDIChwbkH9ulNzNwtcJNkyurmB+AKGCV5pYdHOf4JLElFeQKQpKIMAEkqygCQpKIMAEkqygCQpKIMAEkqygCQpKIMAEkq6hdDhr3Xq8m4kgAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 504x648 with 2 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "plt.figure(figsize=(7, 9))\n",
    "plt.imshow(updated_mask, cmap=plt.cm.get_cmap('Accent_r', 2))\n",
    "plt.title('{} design with {} sample sites'.format('Adapted Stratified Design', len(x_adpt)))\n",
    "plt.axis('off')\n",
    "cbar = plt.colorbar(fraction=0.02, orientation='horizontal', pad=0.01)\n",
    "cbar.set_ticks([0, 1])\n",
    "cbar.set_ticklabels(['0: Invalid', '1: Valid'])\n",
    "nsampled = sum(sampled_csv.sampled==1)\n",
    "nsp = len(sampled_csv)\n",
    "plt.scatter(y_adpt[:nsampled], x_adpt[:nsampled], c='black', marker='x', linewidth=1.5, s=70)\n",
    "plt.scatter(y_adpt[nsampled:], x_adpt[nsampled:], c='red', marker='x', linewidth=1.5, s=70)\n",
    "plt.savefig('{}/{}_{}Site_Adapted_Stratified_Plot.tif'.format(directory, ts, nsp))\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Black crosses represent sites which have already been visited, while red show the sites which have been moved based on the new mask constraints."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Save the design"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 29,
   "metadata": {},
   "outputs": [],
   "source": [
    "def LongLatConvert(x, y, projection_info, GeoT):\n",
    "    \n",
    "    # Project row and column (will be for bottom left corner)\n",
    "    x_proj = x * GeoT[1] + GeoT[0] \n",
    "    y_proj = y * GeoT[5] + GeoT[3]\n",
    "    \n",
    "    # Shift from bottom left to centre of the pixel\n",
    "    x_proj += GeoT[1] / 2.0\n",
    "    y_proj += GeoT[5] / 2.0 \n",
    "    \n",
    "    xy_proj = np.stack((x_proj, y_proj), axis=-1)\n",
    "    \n",
    "    # Make spatial coordinate system\n",
    "    srs = osr.SpatialReference()\n",
    "    if srs.ImportFromWkt(projection_info) != 0:\n",
    "        print(\"Error: cannot import projection '%s'\" % projection_info)\n",
    "        sys.exit(1)\n",
    "        \n",
    "    srsLatLong = srs.CloneGeogCS()\n",
    "    ct = osr.CoordinateTransformation(srs, srsLatLong)\n",
    "    long, lat, height = list(zip(*ct.TransformPoints(xy_proj)))\n",
    "    return long, lat"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 30,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Convert from row/col to projected\n",
    "long, lat = LongLatConvert(x_adpt, y_adpt, prj_info, GeoT)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now save the results to a csv file..."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 31,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Design results saved as a csv in results/Adapted_Stratified_Demo_Opt2 directory \n",
      "File name: 2019_10_30_140158_30Site_Adapted_Stratified.csv\n"
     ]
    }
   ],
   "source": [
    "# Record the sites which have already been sampled\n",
    "nsp = len(sampled_csv)\n",
    "nsampled = sum(sampled_csv.sampled==1)\n",
    "\n",
    "# Reformat and save to csv\n",
    "result = pd.DataFrame()\n",
    "result['longitude'] = long; result['latitude'] = lat\n",
    "result['row'] = x_adpt; result['col'] = y_adpt\n",
    "result['sampled'] = [1]*nsampled + [0]*(nsp-nsampled)\n",
    "result.index += 1\n",
    "csv_filename = '{}_{}Site_Adapted_Stratified.csv'.format(ts, nsp)\n",
    "result.to_csv('{}/{}'.format(directory, csv_filename), index_label='site')\n",
    "\n",
    "print('Design results saved as a csv in {} directory \\nFile name: {}'.format(directory, csv_filename))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Example output:"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "You should now be able to navigate to the results folder (specified above) and open the csv using your preferred software. We can load and print the csv here as a quick check. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 32,
   "metadata": {},
   "outputs": [],
   "source": [
    "view_csv = pd.read_csv('{}/{}'.format(directory, csv_filename))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 33,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>site</th>\n",
       "      <th>longitude</th>\n",
       "      <th>latitude</th>\n",
       "      <th>row</th>\n",
       "      <th>col</th>\n",
       "      <th>sampled</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>1</td>\n",
       "      <td>172.575026</td>\n",
       "      <td>-42.517812</td>\n",
       "      <td>834</td>\n",
       "      <td>377</td>\n",
       "      <td>1</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>2</td>\n",
       "      <td>172.268461</td>\n",
       "      <td>-42.419761</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>1</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>3</td>\n",
       "      <td>172.336065</td>\n",
       "      <td>-42.615136</td>\n",
       "      <td>174</td>\n",
       "      <td>726</td>\n",
       "      <td>1</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>4</td>\n",
       "      <td>172.447568</td>\n",
       "      <td>-42.422772</td>\n",
       "      <td>491</td>\n",
       "      <td>19</td>\n",
       "      <td>1</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <td>5</td>\n",
       "      <td>172.638446</td>\n",
       "      <td>-42.415069</td>\n",
       "      <td>1015</td>\n",
       "      <td>0</td>\n",
       "      <td>1</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>5</th>\n",
       "      <td>6</td>\n",
       "      <td>172.273458</td>\n",
       "      <td>-42.486438</td>\n",
       "      <td>10</td>\n",
       "      <td>247</td>\n",
       "      <td>1</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>6</th>\n",
       "      <td>7</td>\n",
       "      <td>172.643392</td>\n",
       "      <td>-42.478758</td>\n",
       "      <td>1024</td>\n",
       "      <td>236</td>\n",
       "      <td>1</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>7</th>\n",
       "      <td>8</td>\n",
       "      <td>172.272357</td>\n",
       "      <td>-42.612079</td>\n",
       "      <td>0</td>\n",
       "      <td>712</td>\n",
       "      <td>1</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>8</th>\n",
       "      <td>9</td>\n",
       "      <td>172.407246</td>\n",
       "      <td>-42.445969</td>\n",
       "      <td>379</td>\n",
       "      <td>103</td>\n",
       "      <td>1</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>9</th>\n",
       "      <td>10</td>\n",
       "      <td>172.494765</td>\n",
       "      <td>-42.613175</td>\n",
       "      <td>608</td>\n",
       "      <td>726</td>\n",
       "      <td>0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>10</th>\n",
       "      <td>11</td>\n",
       "      <td>172.645018</td>\n",
       "      <td>-42.610305</td>\n",
       "      <td>1019</td>\n",
       "      <td>723</td>\n",
       "      <td>0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>11</th>\n",
       "      <td>12</td>\n",
       "      <td>172.418827</td>\n",
       "      <td>-42.538764</td>\n",
       "      <td>405</td>\n",
       "      <td>447</td>\n",
       "      <td>0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>12</th>\n",
       "      <td>13</td>\n",
       "      <td>172.545802</td>\n",
       "      <td>-42.443373</td>\n",
       "      <td>759</td>\n",
       "      <td>100</td>\n",
       "      <td>0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>13</th>\n",
       "      <td>14</td>\n",
       "      <td>172.318540</td>\n",
       "      <td>-42.547257</td>\n",
       "      <td>130</td>\n",
       "      <td>474</td>\n",
       "      <td>0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>14</th>\n",
       "      <td>15</td>\n",
       "      <td>172.485336</td>\n",
       "      <td>-42.493344</td>\n",
       "      <td>590</td>\n",
       "      <td>282</td>\n",
       "      <td>0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>15</th>\n",
       "      <td>16</td>\n",
       "      <td>172.568583</td>\n",
       "      <td>-42.581389</td>\n",
       "      <td>812</td>\n",
       "      <td>612</td>\n",
       "      <td>0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>16</th>\n",
       "      <td>17</td>\n",
       "      <td>172.645247</td>\n",
       "      <td>-42.549515</td>\n",
       "      <td>1024</td>\n",
       "      <td>498</td>\n",
       "      <td>0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>17</th>\n",
       "      <td>18</td>\n",
       "      <td>172.354996</td>\n",
       "      <td>-42.492798</td>\n",
       "      <td>233</td>\n",
       "      <td>274</td>\n",
       "      <td>0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>18</th>\n",
       "      <td>19</td>\n",
       "      <td>172.415081</td>\n",
       "      <td>-42.599328</td>\n",
       "      <td>391</td>\n",
       "      <td>671</td>\n",
       "      <td>0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>19</th>\n",
       "      <td>20</td>\n",
       "      <td>172.497704</td>\n",
       "      <td>-42.552619</td>\n",
       "      <td>620</td>\n",
       "      <td>502</td>\n",
       "      <td>0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>20</th>\n",
       "      <td>21</td>\n",
       "      <td>172.342830</td>\n",
       "      <td>-42.418914</td>\n",
       "      <td>204</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>21</th>\n",
       "      <td>22</td>\n",
       "      <td>172.417770</td>\n",
       "      <td>-42.492038</td>\n",
       "      <td>405</td>\n",
       "      <td>274</td>\n",
       "      <td>0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>22</th>\n",
       "      <td>23</td>\n",
       "      <td>172.372333</td>\n",
       "      <td>-42.567698</td>\n",
       "      <td>276</td>\n",
       "      <td>552</td>\n",
       "      <td>0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>23</th>\n",
       "      <td>24</td>\n",
       "      <td>172.318503</td>\n",
       "      <td>-42.458641</td>\n",
       "      <td>135</td>\n",
       "      <td>146</td>\n",
       "      <td>0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>24</th>\n",
       "      <td>25</td>\n",
       "      <td>172.602838</td>\n",
       "      <td>-42.448533</td>\n",
       "      <td>915</td>\n",
       "      <td>122</td>\n",
       "      <td>0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>25</th>\n",
       "      <td>26</td>\n",
       "      <td>172.489811</td>\n",
       "      <td>-42.451139</td>\n",
       "      <td>605</td>\n",
       "      <td>126</td>\n",
       "      <td>0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>26</th>\n",
       "      <td>27</td>\n",
       "      <td>172.271506</td>\n",
       "      <td>-42.570212</td>\n",
       "      <td>0</td>\n",
       "      <td>557</td>\n",
       "      <td>0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>27</th>\n",
       "      <td>28</td>\n",
       "      <td>172.270635</td>\n",
       "      <td>-42.527265</td>\n",
       "      <td>0</td>\n",
       "      <td>398</td>\n",
       "      <td>0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>28</th>\n",
       "      <td>29</td>\n",
       "      <td>172.538805</td>\n",
       "      <td>-42.485343</td>\n",
       "      <td>737</td>\n",
       "      <td>255</td>\n",
       "      <td>0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>29</th>\n",
       "      <td>30</td>\n",
       "      <td>172.462579</td>\n",
       "      <td>-42.581981</td>\n",
       "      <td>522</td>\n",
       "      <td>609</td>\n",
       "      <td>0</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "    site   longitude   latitude   row  col  sampled\n",
       "0      1  172.575026 -42.517812   834  377        1\n",
       "1      2  172.268461 -42.419761     0    0        1\n",
       "2      3  172.336065 -42.615136   174  726        1\n",
       "3      4  172.447568 -42.422772   491   19        1\n",
       "4      5  172.638446 -42.415069  1015    0        1\n",
       "5      6  172.273458 -42.486438    10  247        1\n",
       "6      7  172.643392 -42.478758  1024  236        1\n",
       "7      8  172.272357 -42.612079     0  712        1\n",
       "8      9  172.407246 -42.445969   379  103        1\n",
       "9     10  172.494765 -42.613175   608  726        0\n",
       "10    11  172.645018 -42.610305  1019  723        0\n",
       "11    12  172.418827 -42.538764   405  447        0\n",
       "12    13  172.545802 -42.443373   759  100        0\n",
       "13    14  172.318540 -42.547257   130  474        0\n",
       "14    15  172.485336 -42.493344   590  282        0\n",
       "15    16  172.568583 -42.581389   812  612        0\n",
       "16    17  172.645247 -42.549515  1024  498        0\n",
       "17    18  172.354996 -42.492798   233  274        0\n",
       "18    19  172.415081 -42.599328   391  671        0\n",
       "19    20  172.497704 -42.552619   620  502        0\n",
       "20    21  172.342830 -42.418914   204    0        0\n",
       "21    22  172.417770 -42.492038   405  274        0\n",
       "22    23  172.372333 -42.567698   276  552        0\n",
       "23    24  172.318503 -42.458641   135  146        0\n",
       "24    25  172.602838 -42.448533   915  122        0\n",
       "25    26  172.489811 -42.451139   605  126        0\n",
       "26    27  172.271506 -42.570212     0  557        0\n",
       "27    28  172.270635 -42.527265     0  398        0\n",
       "28    29  172.538805 -42.485343   737  255        0\n",
       "29    30  172.462579 -42.581981   522  609        0"
      ]
     },
     "execution_count": 33,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "view_csv"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Columns consist of:\n",
    "* A **site** ID\n",
    "* The coordinates in **longitude** and **latitude**\n",
    "* The coordinates in **row** and **column** (corresponding to pixels in image)\n",
    "* A **sampled** column: There should be ones in the sites which have already been sampled, and zeros in all the sites which have been adjusted based on the new constraints. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.7.3"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 1
}