01_geohash_seattle.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2019-10-24T13:58:07.035103Z",
     "start_time": "2019-10-24T13:58:05.761858Z"
    }
   },
   "outputs": [],
   "source": [
    "import pandas as pd\n",
    "from math import sin, cos, sqrt, atan2, radians\n",
    "import numpy as np\n",
    "import os\n",
    "import time\n",
    "import gc\n",
    "import geopandas as gpd\n",
    "from geolib import geohash\n",
    "import matplotlib.pyplot as plt\n",
    "# geohash.encode(latitude, longitude, precision)\n",
    "# >>> geohash.encode('70.2995', '-27.9993', 7)\n",
    "# >>> gkkpfve\n",
    "# >>> geohash.decode('gkkpfve')\n",
    "# >>> (70.2995, -27.9993)\n",
    "# >>> neighbours = geohash.neighbours('gcpuyph')\n",
    "# >>> neighbours\n",
    "%matplotlib inline"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2019-10-24T13:58:07.046257Z",
     "start_time": "2019-10-24T13:58:07.039079Z"
    }
   },
   "outputs": [],
   "source": [
    "from shapely.geometry import Point, Polygon"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2019-10-24T13:58:07.789960Z",
     "start_time": "2019-10-24T13:58:07.050744Z"
    }
   },
   "outputs": [],
   "source": [
    "city_df = gpd.read_file(\"king_county_admin_mapfiles/city.shp\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2019-10-24T13:58:08.510468Z",
     "start_time": "2019-10-24T13:58:07.793318Z"
    }
   },
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAWMAAAHfCAYAAACFyrnvAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8QZhcZAAAgAElEQVR4nOzdeViU1QIG8PcIwyaYqARugIGamhvgvmsu1yXLJTXNXa+WlnUzs0zLsjTbNCvLMDVTc08rNUtTQQXB3BXFVEpRSSzBYVjP/YOBkFhmhhnOLO/veeZp+Ga+b96pfD2cbxNSShARkVoVVAcgIiKWMRGRVWAZExFZAZYxEZEVYBkTEVkBZ9UBzK1atWoyMDBQdQwisjOxsbF/Sil9LLV9uyvjwMBAxMTEqI5BRHZGCHHFktvnNAURkRVgGRMRWQGWMRGRFWAZExFZAZYxEZEVYBkTEVkBljERkRVgGRMRWQGWMRGRFWAZExFZAZYxEZEVYBkTEVkBljERkRVgGRMRWQGWMRGRFWAZExFZAZYxEZEVYBkTEVkBu7vtki07ePAgFixYAGdnZ7i4uNzz0Gg0Jf5s6LKS3uPs7AwhhOp/DUQOiWVsRV566SW0bNkSrVu3RkZGBjIzM5GRkZH/KPizVqst8j2lrVfSsuzsbLOUuiHvKeq1wstL+wyNRoMKFfjLHdkHlrGVOHToEBISErBnzx44O6v5z5KTk3NPSZtS6sW9JzU1tcj3ZWZmIj09vcj1C7+n4D/zXiv4W4RGoyn2Udrrzs7O+dvKe174Lw43N7d/PVxdXYtc7ubmBnd3d7i7u0Oj0fA3DioVy9hKvP3225gxY4ayIgaAChUqwNXVFa6ursoyGENKiaysrCJLPisrK7+0C5Z34WV5j6ysrPxt5T3PWyc9PR0pKSnQ6XTQ6XRIT09HWlpa/vOCP+c98n5OS0tDTk7OPQVduKyLenh4eBT7c2l/CeQ9OO1kW1jGVuDUqVOIiYnB+vXrVUexKUKI/JFtxYoVVccpVlZWVrFFXdxDq9UiLS0Nf//9N65fv460tDTcvXsX6enpRf6lUNQjJyfnX6P/kn47KO7h4+ODjz76iFNCFsYytgLLli3DhAkT4ObmpjoKWYCzszM8PT3h6elZrp+bnZ39r9F/cb8ZlPR47rnnMG3aNNStW7dc8zsalrFiUkps3LgRe/fuVR2F7IyTkxOcnJzK/Jf8N998g9jYWJaxhfH3DsWuXbuGjIwM/o9OVis0NBRHjx5VHcPusYwVO3r0KEJDQ7mjhaxWSEgIYmNjVceweyxjxWJjYxESEqI6BlGxQkJCcPToUUgpVUexayxjxWJjYxEaGqo6BlGxfH19UbFiRVy6dEl1FLvGMlbs6NGjHBmT1eO8seWxjBW6fv06dDodAgMDVUchKhHnjS2PZaxQ3nwxd96RtQsNDWUZWxjLWCHuvCNbwZ14lscyVujkyZNo2rSp6hhEpapRowY0Gg0SEhJUR7FbLGOF4uPjebIH2QzuxLMslrEiUkpcvHgRwcHBqqMQGYQ78SyLZaxIUlISNBoNvL29VUchMghHxpbFMlYkPj4eQUFBqmMQGSxvZMydeJbBMlaEUxRka2rVqgUpJa5du6Y6il1iGSvCkTHZGiEE540tiGWsSHx8PEfGZHM4b2w5LGNFOE1BtogjY8thGSvCaQqyRTwt2nJYxgr89ddf0Ol08PX1VR2FyCgBAQFIT09HYmKi6ih2h2WswMWLFxEUFMQLBJHNyduJx3lj82MZK8Cdd2TLWMaWwTJWgDvvyJZx3tgyWMYKcOcd2TKOjC2DZawAR8Zky4KCgnDnzh0kJSWpjmJXWMbl7Pbt2zh58iQaNGigOgqRSYQQaN68OUfHZsYyLmfz5s3D4MGDUb16ddVRiEwWGhqKmJgY1THsCsu4HP3222/48ssv8frrr6uOQlQmXbp0wY4dO1THsCss43L08ssvY9q0afDz81MdhahMunfvjjNnzuDq1auqo9gNlnE5WbFiBaKiovD888+rjkJUZi4uLnjkkUewadMm1VHsBsu4HGzatAkzZ87Ezp07UbFiRdVxiMxi8ODBWL9+veoYdoNlbGG7du3C5MmT8cMPP6B+/fqq4xCZDacqzItlbEEREREYMWIEtmzZgubNm6uOQ2RWnKowL5axhRw9ehQDBgzA119/jXbt2qmOQ2QRnKowH5axBZw7dw59+vTBp59+ih49eqiOQ2QxnKowH5axmV2+fBk9evTA22+/jYEDB6qOQ2RRnKowH5axGd28eRPdu3fHCy+8gNGjR6uOQ1QuOFVhHixjM3rllVfQp08fPPPMM6qjEJUbTlWYB8vYTC5cuICtW7dizpw5qqMQlStOVZgHy9hMPvzwQ0yaNAne3t6qoxCVO05VlB3L2Azu3LmDtWvXYtKkSaqjECnBqYqyYxmbwapVq/Dwww+jZs2aqqMQKcGpirJjGZeRlBJLlizBlClTVEchUopTFWVTahkLIWoLIfYKIc4KIU4LIZ4t8NpUIUScfvk7BZbPFELE61/rWWB5qBDipP61xUJ/r3ohhKsQ4hv98ighRGCBdUYJIS7oH6PM9cXN5eeff4aLiws6dOigOgqRUpyqKBtDRsZZAP4npWwAoDWAp4UQDYUQXQD0B9BEStkIwLsAIIRoCGAogEYAegH4RAjhpN/WpwAmAqirf/TSLx8H4LaUMhjABwAW6LdVBcAcAK0AtAQwRwhhVXvI8kbF+r9XiBwWpyrKptQyllImSimP6p+nADgLoCaAyQDmSynT9a/d1K/SH8A6KWW6lPISgHgALYUQ1QFUklIeklJKAKsAPFpgnZX65xsBdNOPmnsC2C2lTJZS3gawG/8UuHJXrlxBREQEhg8frjoKkVXgVIXpjJoz1k8fNAcQBaAegA76aYV9QogW+rfVBPB7gdX+0C+rqX9eePk960gpswD8DaBqCdsqnGuiECJGCBFTnnesXbp0KZ588kleo5hIj1MVpjO4jIUQngA2AZgmpbwDwBmAN3KnLqYDWK8fzRb1+7osYTlMXOefBVJ+LqUMk1KG+fj4lPpdzEGn0yE8PBxPPfVUuXwekS3gVIXpDCpjIYQGuUX8tZRys37xHwA2y1zRAHIAVNMvr11g9VoArumX1ypiOQquI4RwBnAfgOQStqXchg0bEBISgrp166qOQmRVOFVhGkOOphAAwgGclVK+X+ClrQC66t9TD4ALgD8BbAMwVH+ERB3k7qiLllImAkgRQrTWb3MkgG/129oGIO9IiUEA9ujnlXcB6CGE8NbvuOuhX6bcF198gYkTJ6qOQWR1unfvjsTERKxbt051FJvibMB72gF4EsBJIcQx/bKXASwHsFwIcQpABoBR+gI9LYRYD+AMco/EeFpKma1fbzKAFQDcAezQP4Dcsv9KCBGP3BHxUACQUiYLId4AcET/vrlSymRTv6y5XLhwAefOnUPfvn1VRyGyOi4uLti4cSN69OiBJk2aoGHDhqoj2QSR25/2IywsTMbExFj0M15++WWkp6fjvffes+jnENmyL7/8Eu+88w6io6Ph5eWlOk6ZCSFipZRhlto+z8AzUlZWFlauXIlx48apjkJk1caMGYP27dtj/PjxsLdBnyWwjI20a9cu+Pv781cvIgN89NFHiI+Px6JFi1RHsXqGzBlTAeHh4Rg7dqzqGEQ2wc3NDRs3bkSrVq0QFhaG9u3bq45ktTgyNsKNGzewd+9eDBkyRHUUIptRp04dfPnllxg6dCiuX7+uOo7VYhkbYfXq1Xj00UdRqVIl1VGIbEqfPn0wZswYDBs2DFlZWarjWCWWsYGklJyiICqD1157DRqNBrNmzVIdxSqxjA10+PBhZGdnc86LyEROTk5Ys2YN1q5di61bt6qOY3VYxgZavnw5xo4dy0tlEpVBtWrVsH79ekycOBHR0dGq41gVlrEBUlNTsXHjRowcOVJ1FCKb16pVKyxevBiDBg1CmzZtsGLFCmi1WtWxlGMZG2DDhg3o0KEDqlevrjoKkV0YOnQofvvtN7z88svYsGED/P39MW3aNJw9e1Z1NGVYxgbgjjsi83N2dka/fv3w/fffIyYmBp6enujatSs6d+6MtWvXIj09XXXEcsUyLkVcXBzi4+PRp08f1VGI7FZgYCDefPNNJCQkYMqUKVi+fDn8/f0xY8YMXLx4UXW8csEyLsXy5csxcuRIaDQa1VGI7J5Go8GgQYOwe/duREREICcnB23atEHPnj2xefNmuz5GmVdtK0FmZib8/f2xd+9ePPjgg2bZJhEZR6fTYfPmzVi0aBECAgLwzTffKDmqiVdtU2jHjh144IEHWMRECrm5ueGJJ57Avn37cO7cOYSHh6uOZBEs4xKEh4fzUplEVsLNzQ3r1q3DzJkz7fKoC5ZxMRITE7F//34MHjxYdRQi0mvYsCHmzZuHYcOGQafTqY5jVizjYnz11VcYMGCAXdyhgMieTJgwAUFBQXjppZdURzErlnER8i4KxCkKIusjhMCyZcuwefNmfP/996rjmA3LuAiRkZEQQqBNmzaqoxBREapUqYLVq1dj/PjxSExMVB3HLFjGRVi+fDnGjRvHiwIRWbGOHTti4sSJGDlyJHJyclTHKTOWcSGpqanYsmULLwpEZANeffVVpKWl4d1331UdpcxYxoWcOXMGderUga+vr+ooRFQKZ2dnfP3113j33Xdx5MgR1XHKhGVcyOXLl1GnTh3VMYjIQAEBAViyZAnGjBmjOkqZsIwLuXLlCgICAlTHICIjDB48GNevX8fVq1dVRzEZy7iQy5cvIzAwUHUMIjJC3tFPhw4dUh3FZCzjQjgyJrJNbdu2xcGDB1XHMBnLuJArV65wZExkg1jGdkRKicuXL3NkTGSDWrRogZMnT9rsNStYxgXcvn0bTk5OqFy5suooRGQkDw8PNGjQALGxsaqjmIRlXABHxUS2zZanKljGBXDnHZFtYxnbiYSEBPj7+6uOQUQmyitjW7ydHMu4AJ1Oh4oVK6qOQUQmql27NjQaDX777TfVUYzGMi4gOzsbFSrwXwmRrbLlkz/YPAXk5OTAyclJdQwiKoNmzZrh1KlTqmMYjWVcQE5ODkfGRDZOo9HY5PWN2TwFsIyJSBU2TwEsYyJShc1TAHfgEZEqbJ4CuAOPiFRhGRfAaQoiUoXNUwDLmIhUYfMU4ObmBq1WqzoGETkglnEBfn5+uHHjhuoYROSAWMYF+Pr6soyJbFxWVpbqCCZhGRfg5+eH69evq45BRCaSUmLjxo3o0KGD6ihGYxkX4Ofnh8TERNUxiMhEP/30E9LT09GnTx/VUYzGMi6gdu3aSExMREZGhuooRGSCd955By+++KJNHhVle4ktyNXVFQEBAbhw4YLqKERkpNjYWJw7dw7Dhg1THcUkLONCGjZsiDNnzqiOQURGyM7OxuzZs/Hcc8/BxcVFdRyTsIwLYRkT2ZacnBxMmDABWq0WkydPVh3HZCzjQho0aMAyJrIROTk5mDRpEuLj47F9+3a4u7urjmQylnEhHBkT2QYpJaZMmYLTp0/j+++/h6enp+pIZcIyLqR+/fqIj4+32QPHiRyBlBLTpk1DbGwsduzYAS8vL9WRyoxlXIiHhwceeOABnDhxQnUUIiqClBLTp09HZGQkdu3ahUqVKqmOZBYs4yK0a9cOkZGRqmMQUSGZmZn43//+h59//hk//vgjKleurDqS2bCMi8AyJrI+Z8+eRZs2bXDmzBns3r0bVapUUR3JrFjGRWAZE1mPnJwcfPDBB+jQoQPGjx+PHTt2oFq1aqpjmZ2z6gDWKCgoCBkZGUhISIC/v7/qOEQO6/Llyxg9ejSysrIQFRWFoKAg1ZEshiPjIgghODomUkhKifDwcLRo0QK9e/fGvn377LqIAY6Mi9WuXTtERETY7HnuRLYqMTEREyZMwNWrV7Fnzx40btxYdaRywZFxMTgyJip/GzZsQLNmzdC8eXNERUU5TBEDHBkXKyQkBPHx8bhz547dHMdIZK2Sk5MxZcoUxMbGYtu2bWjVqpXqSOWOI+NiuLi4ICQkBIcPH1Ydhciu7dy5E02aNIGPjw9+/fVXhyxigGVcIk5VEFnWBx98gEmTJmHVqlVYtGgRPDw8VEdShmVcApYxkeUcPXoUb7/9Nvbv34+uXbuqjqMcy7gEbdu2RXR0NC8aRGRmWq0WTzzxBBYtWsRj+fVYxiWoUqUKfH19eRsmIjObPn06QkNDeehoATyaohT16tXDhQsX0KBBA9VRiOzCDz/8gO+++w7Hjx9XHcWqsIxLERwczJExkZncvHkT48ePx9q1a+3qimvmwGmKUtStW5dlTGQGUkqMHz8eI0eORKdOnVTHsTos41LUrVsX8fHxqmMQ2bxly5bhjz/+wNy5c1VHsUqcpigFpymIyu7SpUt45ZVXsH//fri4uKiOY5U4Mi5FQEAAbty4AZ1OpzoKkc36+OOPMXbsWO4ILwFHxqVwdnaGj48Pbt68yeMhiUyQkZGBr776CgcOHFAdxaqVOjIWQtQWQuwVQpwVQpwWQjxb6PUXhBBSCFGtwLKZQoh4IUScEKJngeWhQoiT+tcWCyGEfrmrEOIb/fIoIURggXVGCSEu6B+jzPGljVWlShUkJyer+Ggim/fdd9+hfv36qFevnuooVs2QaYosAP+TUjYA0BrA00KIhkBuUQPoDiAh783614YCaASgF4BPhBBO+pc/BTARQF39o5d++TgAt6WUwQA+ALBAv60qAOYAaAWgJYA5Qghvk7+tiVjGRKYLDw/HuHHjVMeweqWWsZQyUUp5VP88BcBZADX1L38A4EUAssAq/QGsk1KmSykvAYgH0FIIUR1AJSnlISmlBLAKwKMF1lmpf74RQDf9qLkngN1SymQp5W0Au/FPgZcbljGRaa5evYpDhw5h0KBBqqNYPaN24OmnD5oDiBJCPALgqpSy8Gk0NQH8XuDnP/TLauqfF15+zzpSyiwAfwOoWsK2CueaKISIEULEJCUlGfOVDFK1alXcunXL7NslsncrVqzA4MGDUbFiRdVRrJ7BZSyE8ASwCcA05E5dvAJgdlFvLWKZLGG5qev8s0DKz6WUYVLKMB8fnyJWKRuOjImMl5OTg+XLl3OKwkAGlbEQQoPcIv5aSrkZQBCAOgCOCyEuA6gF4KgQwg+5o9faBVavBeCafnmtIpaj4DpCCGcA9wFILmFb5YplTGS8ffv2wcPDAy1atFAdxSYYcjSFABAO4KyU8n0AkFKelFLeL6UMlFIGIrc0Q6SU1wFsAzBUf4REHeTuqIuWUiYCSBFCtNZvcySAb/Ufsw1A3pESgwDs0c8r7wLQQwjhrd9x10O/rFyxjImMl7fjTn/QFJXCkJFxOwBPAugqhDimf/Qu7s1SytMA1gM4A2AngKellNn6lycD+AK5O/UuAtihXx4OoKoQIh7A8wBe0m8rGcAbAI7oH3P1y8oVy5jIeFFRUejVq9z3t9usUk/6kFJGoOi524LvCSz08zwA84p4XwyAh4pYrgMwuJhtLwewvLSclsQdeETG02q18PT0VB3DZvB0aANwZExkvLS0NLi7u6uOYTNYxgZwcXFBRkaG6hhENoVlbByWsQH4PxWRcaSUSE9Ph5ubm+ooNoNlbACWMZFxdDodXFxcUKECK8ZQ/DdlAJYxkXH4Z8Z4LGMD8H8sIuOkpaXBw8NDdQybwjI2AMuYyDj8M2M8lrEBtFot/8ciMgLL2HgsYwPwfywi43AAYzyWsQFYxkTG4Z8Z47GMDZCQkIAaNWqojkFkM1jGxmMZG+DYsWNo1qyZ6hhENoNlbDyWcSmklDh+/DjLmMgILGPjsYxLcfXqVTg7O8PPz091FCKbweOMjccyLgWnKIiMl5qayjI2Esu4FMePH0fTpk1VxyCyKb/99hsCAwNVx7ApLONScGRMZLzz58+jfv36qmPYFJZxKbjzjsh4cXFxLGMjsYxLkJqaiqtXr6JevXqqoxDZjLS0NCQmJnKawkgs4xKcPHkSDRs2hLNzqbcKJCK9+Ph41KlTh39ujMQyLgF33hEZj1MUpmEZlyAuLg4NGjRQHYPIpnDnnWlYxiW4efMmfH19VccgsilxcXHcz2IClnEJkpKSUK1aNdUxiGwKpylMwzIuQVJSEnx8fFTHILIZUkqWsYlYxiW4desWKleurDoGkc1ISkqCEIK/UZqAZVyC5s2bIzIyUnUMIptx/vx51KtXD0II1VFsDsu4BEOGDME333yjOgaRzeAUhelYxiXo168fIiIikJycrDoKkU04cuQILx9gIpZxCby8vNC9e3ds2bJFdRQimxAREYEOHTqojmGTWMal4FQFkWGSk5ORkJDAkbGJWMal6NOnD6Kjo5GUlKQ6CpFVi4yMRKtWrXhNChOxjEvh4eGB//znP9i0aZPqKA5h3rx5aNSoEZo0aYJmzZohKirK6G1s3boVZ86cyf95xYoVuHbtWv7PnTt3RkxMjFny0j8iIiLQvn171TFsFsvYAEOHDsXq1atVx7B7hw4dwnfffYejR4/ixIkT+Omnn1C7dm2jt1NaGZNlcL64bFjGBujduzcuX76M48ePq45i1xITE1GtWjW4uroCAKpVq4YaNWogNjYWnTp1QmhoKHr27InExEQAwLJly9CiRQs0bdoUAwcOhFarxcGDB7Ft2zZMnz4dzZo1w4IFCxATE4Phw4ejWbNmSEtLu+czf/zxR7Rp0wYhISEYPHgwUlNTy/1724O0tDQcP34crVq1Uh3Fdkkp7eoRGhoqLWHevHly3LhxFtk25UpJSZFNmzaVdevWlZMnT5a//PKLzMjIkG3atJE3b96UUkq5bt06OWbMGCmllH/++Wf+uq+88opcvHixlFLKUaNGyQ0bNuS/1qlTJ3nkyJF//ZyUlCQ7dOggU1NTpZRSzp8/X77++usW/572aN++fbJly5aqY1gUgBhpwe7iTLuBJkyYgHr16mHBggWoWrWq6jh2ydPTE7GxsThw4AD27t2LIUOGYNasWTh16hS6d+8OAMjOzkb16tUBAKdOncKsWbPw119/ITU1FT179jTq8w4fPowzZ86gXbt2AICMjAy0adPGvF/KQXC+uOxYxgby8fFB//798cUXX2DGjBmq49gtJycndO7cGZ07d0bjxo3x8ccfo1GjRjh06NC/3jt69Ghs3boVTZs2xYoVK/DLL78Y9VlSSnTv3h1r1641U3rHFRERgQkTJqiOYdM4Z2yEqVOn4pNPPkFWVpbqKHYpLi4OFy5cyP/52LFjaNCgAZKSkvLLODMzE6dPnwYApKSkoHr16sjMzMTXX3+dv56XlxdSUlKK/TlP69atERkZifj4eACAVqvF+fPnLfLd7Fl2djYOHjzIkXEZsYyNEBoaipo1a2L79u2qo9il1NRUjBo1Cg0bNkSTJk1w5swZzJ07Fxs3bsSMGTPQtGlTNGvWDAcPHgQAvPHGG2jVqhW6d++OBx98MH87Q4cOxcKFC9G8eXNcvHgRo0ePxqRJk/61A8/HxwcrVqzAsGHD0KRJE7Ru3Rrnzp0r9+9t606dOgU/Pz9ebraMRO68tP0ICwuTljyGdO3atfj888+xd+9ei30GkS35+OOP8euvv+KLL75QHcWihBCxUsowS22fI2MjDRw4EHFxcTh58qTqKERW4cCBAzy+2AxYxkZycXHBpEmTsGTJEtVRiJSTUuLAgQOcLzYDlrEJJk6ciPXr1+P27duqoxApdeXKFUgp8cADD6iOYvNYxibw8/NDt27deGlNcnh5o2Le2aPsWMYmGjhwIDZv3qw6BpFSPNnDfFjGJurTpw/279+PO3fuqI5CpExsbCxatmypOoZdYBmbqFKlSujQoQN27typOgqREllZWTh79iwaN26sOopdYBmXQYsWLXiIGzms8+fPo2bNmqhYsaLqKHaBZVwGgYGBuHz5suoYREocP34cTZs2VR3DbrCMy6BOnTq4dOmS6hhEShw/fhxNmjRRHcNusIzLwM3NDRkZGapjECnBkbF5sYzLIDk5GVWqVFEdg0gJlrF5sYzL4Pbt2/D29lYdg6jcJSUlIS0tDf7+/qqj2A2WcRmwjMlRnThxAk2aNOGZd2bEMi4DljE5Kk5RmB/LuAySk5NZxuSQWMbmxzIug9u3b3MHHjkkHtZmfizjMuA0BTmijIwMnD9/Hg899JDqKHaFZVwGnKYgR3Tu3DkEBATA3d1ddRS7wjIuA46MyRFxvtgyWMZlwDljckQsY8tgGZcBR8bkiE6cOMEytgCWsYnS09ORkZHByweSw+HI2DJYxibKGxXzDCRyJNevX0dWVhZq1KihOordYRmbiPPF5IjyRsUchJgfy9hEnC8mR8QpCsthGZuIxxiTI2IZWw7L2EQcGZMjYhlbDsvYRH/++SfnjMmh6HQ6XLx4EQ0bNlQdxS6xjE0UFRWF0NBQ1TGIys3Zs2cRHBwMV1dX1VHsEsvYBDk5OdizZw+6dOmiOgpRuUlISEBgYKDqGHaLZWyC06dPo1KlSggICFAdhajcJCcno2rVqqpj2C2WsQn27NmDrl27qo5BVK5u3brFMrYglrEJWMbkiHg3dMtiGRspKysL+/fv53wxORyWsWWxjI109OhR1KpVC76+vqqjEJUrlrFllVrGQojaQoi9QoizQojTQohn9csXCiHOCSFOCCG2CCEqF1hnphAiXggRJ4ToWWB5qBDipP61xUJ/grsQwlUI8Y1+eZQQIrDAOqOEEBf0j1Hm/PKm4BQFOSrOGVuWISPjLAD/k1I2ANAawNNCiIYAdgN4SErZBMB5ADMBQP/aUACNAPQC8IkQwkm/rU8BTARQV//opV8+DsBtKWUwgA8ALNBvqwqAOQBaAWgJYI4QQulpbyxjclQcGVtWqWUspUyUUh7VP08BcBZATSnlj1LKLP3bDgOopX/eH8A6KWW6lPISgHgALYUQ1QFUklIeklJKAKsAPFpgnZX65xsBdNOPmnsC2C2lTJZS3kbuXwB5BV7u0tPTcejQIXTq1ElVBCJlWMaWZdScsX76oDmAqEIvjQWwQ/+8JoDfC7z2h35ZTf3zwsvvWUdf8H8DqFrCtgrnmiiEiBFCxCQlJRnzlX8B5DkAACAASURBVIwSFRWFBg0aoHLlyqW/mcjO3Lp1i2VsQQaXsRDCE8AmANOklHcKLH8FuVMZX+ctKmJ1WcJyU9f5Z4GUn0spw6SUYT4+PsV/iTL6+eefOUVBDik9PR3p6enw8vJSHcVuGVTGQggNcov4aynl5gLLRwHoC2C4fuoByB291i6wei0A1/TLaxWx/J51hBDOAO4DkFzCtpTgfDE5qrybKfCi8pZjyNEUAkA4gLNSyvcLLO8FYAaAR6SU2gKrbAMwVH+ERB3k7qiLllImAkgRQrTWb3MkgG8LrJN3pMQgAHv05b4LQA8hhLd+x10P/bJyd/fuXfz6669o166dio8nUorzxZbnbMB72gF4EsBJIcQx/bKXASwG4Apgt/5vy8NSyklSytNCiPUAziB3+uJpKWW2fr3JAFYAcEfuHHPePHM4gK+EEPHIHREPBQApZbIQ4g0AR/TvmyulTDb1y5ZFREQEmjdvzhuQkkPifLHllVrGUsoIFD13+0MJ68wDMK+I5TEAHipiuQ7A4GK2tRzA8tJyWlpERAQ6duyoOgaRErxIkOXxDDwDHTx4kFMU5LA4TWF5LGMDZGVlITo6Gq1bt1YdhUgJlrHlsYwNcPz4cfj7+/N/RnJYnDO2PJaxASIjIzlFQQ6Nc8aWxzI2AOeLydH99ddfqFSpkuoYdo1lbACOjMnRBQYGIj4+XnUMu8YyLkVCQgLS09MRFBSkOgqRMi1atEBMTIzqGHaNZVyKvFExTwMlR9aiRQscOXIE/1z1gMyNZVwKzhcTAbVr10ZOTg7++OOP0t9MJmEZl4LzxUSAEIJTFRbGMi5BSkoK4uLiEBISojoKkXJ5UxVkGSzjEkRFRaF58+ZwdXVVHYVIOZaxZbGMS8ApCqJ/hIWFISYmhjvxLIRlXALuvCP6h6+vL7y8vHDx4kXVUewSy7gY2dnZOHz4MNq2bas6CpHV4FSF5bCMi3Hq1Cn4+fmhWrVqqqMQWQ2WseWwjIvB+WKif2vZsiV2796NjIwM1VHsDsu4GJwvJvq3Tp06oU6dOpg6dSp35JkZy7gYHBkT/ZuTkxO+/vprHDx4EB999JHqOHbFkBuSOpxr164hJSUF9evXVx2FyOp4eXlh+/btaNOmDerVq4devXqpjmQXODIuQmRkJNq0acOLAxEVIzAwEBs2bMDIkSNx9uxZ1XHsAsu4CJwvJipd+/bt8c4776Bfv364deuW6jg2j2VchISEBDzwwAOqYxBZvdGjR2PAgAEYNGgQMjMzVcexaSzjIvj6+uL69euqYxDZhLfffhuenp6YMmUKj7AoA5ZxEZo2bYro6GjVMYhsgpOTE9asWYODBw9i0aJFquPYLJZxEf7zn/9g586dyM7OVh2FyCbkHWHx0UcfYcqUKUhPT1cdyeawjIvg7++PGjVqcHRMZITAwEDExsbi+vXraNu2LS8oZCSWcTF69+6N77//XnUMIptSuXJlbNiwAWPGjEGbNm2wceNG1ZFsBsu4GEOGDMHKlSuh0+lURyGyKUIITJkyBT/88ANmzJjBaQsDsYyL0bx5c4SEhODzzz9XHYXIJoWFhXHawggs4xK8/vrrmD9/PrRareooRDaJ0xaGYxmXoFmzZmjbti0++eQT1VGIbBanLQzDMi7Fa6+9hoULFyI1NVV1FCKbVnDaomPHjkhKSlIdyaqwjEvx0EMPoXPnzli2bJnqKEQ2L2/a4uGHH0b79u1x5coV1ZGsBsvYAM888ww+/fRT5OTkqI5CZPOEEJg3bx6eeuoptG/fHqdPn1YdySqwjA3Qtm1buLu7Y8+ePaqjENmNZ599FvPnz0fXrl1x8OBB1XGUYxkbQAiByZMnc0cekZkNHz4cK1euRP/+/R3+JCuWsYGGDx+OX375BVevXlUdhciu9OrVC9u3b8e4ceOwatUq1XGUYRkbyMvLC8OGDeOOPCILaN26Nfbu3YtZs2bh/fffVx1HCZaxESZPnoxly5bxItpEFtCgQQNERERg2bJlmDNnjuo45Y5lbISHHnoIQUFB2LZtm+ooRHbJ398f+/fvx0cffYSEhATVccoVy9hIkydPxqeffqo6BpHd8vHxwciRIx1uhznL2EgDBgzA8ePHcfnyZdVRiOzWlClTEB4ejrS0NNVRyg3L2Eiurq4YOnQoVq5cqToKkd0KDg5Gq1atsGbNGtVRyg3L2ASjR4/GypUreUYekQVNnToVS5YsUR2j3LCMTRASEgJPT0/s379fdRQiu9W+fXvExcWpjlFuWMYmEEJg9OjRWLFiheooRHZLCKE6QrliGZvoiSeewNatW3kHaSIyC5axifz8/FCjRg2cOHFCdRQisgMs4zLo0KEDDhw4oDoGEdkBlnEZ1KxZE3/++afqGERkB1jGZSClVB2BiOwEy7iMHG2PLxFZBsu4DDgyJiJzYRmXgZSSI2MiCzlx4gT8/f1Vxyg3LGMiskorVqzAyJEjVccoN86qA9iynJwcODvzXyGRuaWlpWH9+vU4fvy46ijlhiPjMvjrr7/g7e2tOgaR3fn2228RFhaG2rVrq45SbljGZZCUlAQfHx/VMYjszooVKzB69GjVMcoVy7gMbt68yTImMrOEhARER0fj0UcfVR2lXLGMy4AjYyLze/fddzFu3Dh4eHiojlKuuPepDJKSknD//ferjkFkN27cuIHVq1fjzJkzqqOUO46MTZSTk4Nbt26hWrVqqqMQ2Y33338fw4cPh5+fn+oo5Y4jYxPdvn0bnp6e0Gg0qqMQ2YVbt25h2bJlOHbsmOooSnBkbCJOURCZ1+LFizFgwACHOuuuII6MTcSdd0Tmc+fOHXz88cc4fPiw6ijKcGRsIh7WRmQ+n3zyCXr27Ing4GDVUZThyNhEv/32m0OdHURkKVqtFh9++CF++ukn1VGU4sjYRHv27EHnzp1VxyCyecuWLUPbtm3x0EMPqY6iFEfGJsjIyEBkZCS+/vpr1VGIbNrdu3excOFCfPvtt6qjKMeRsQmio6NRt25dVKlSRXUUIpuVk5ODkSNHonv37ggNDVUdRzmOjE3w/fffo1u3bqpjENm0V199FTdu3MCaNWtUR7EKLGMjHThwAOHh4Th06JDqKEQ266uvvsLatWsRFRUFV1dX1XGsAsvYCJcuXcKQIUOwatUqBAUFqY5DZJMiIyPxv//9D3v37uXhoQVwzthA586dQ8eOHfHqq6+iV69equMQ2aTLly9j0KBBWLVqFRo1aqQ6jlXhyNgAJ06cQK9evfD2229j1KhRquMQ2aQ7d+6gb9++mDlzJgc0RSh1ZCyEqC2E2CuEOCuEOC2EeFa/vIoQYrcQ4oL+n94F1pkphIgXQsQJIXoWWB4qhDipf22x0N9aWQjhKoT4Rr88SggRWGCdUfrPuCCEKPcmjI6ORvfu3fHhhx+yiIlMlJ2djWHDhqFDhw6YOnWq6jhWyZBpiiwA/5NSNgDQGsDTQoiGAF4C8LOUsi6An/U/Q//aUACNAPQC8IkQwkm/rU8BTARQV//I++txHIDbUspgAB8AWKDfVhUAcwC0AtASwJyCpW9p+/fvR9++fREeHo7HH3+8vD6WyO688MILSE9Px+LFi6Efg1EhpZaxlDJRSnlU/zwFwFkANQH0B7BS/7aVAPLukdIfwDopZbqU8hKAeAAthRDVAVSSUh6SUkoAqwqtk7etjQC66UfNPQHsllImSylvA9iNfwrconbt2oWBAwdi7dq16Nu3b3l8JJFd2rBhA3744Qds2LCBl5wtgVE78PTTB80BRAHwlVImArmFDSDvepI1AfxeYLU/9Mtq6p8XXn7POlLKLAB/A6hawrYK55oohIgRQsQkJSUZ85WKtGXLFjz55JPYunUrjycmKqMjR45g3LhxvJN6KQwuYyGEJ4BNAKZJKe+U9NYilskSlpu6zj8LpPxcShkmpQwr66Eyq1evxuTJk7Fz5060a9euTNsiIiAlJQVeXl6qY1g9g8pYCKFBbhF/LaXcrF98Qz/1AP0/b+qX/wGg4OXMagG4pl9eq4jl96wjhHAGcB+A5BK2ZXY5OTmYM2cOZs2ahZ9//hkhISGW+Bgih5OamgpPT0/VMayeIUdTCADhAM5KKd8v8NI2AHmHF4wC8G2B5UP1R0jUQe6Oumj9VEaKEKK1fpsjC62Tt61BAPbo55V3AeghhPDW77jroV9mVlqtFkOGDMFPP/2EqKgoHv9IZEYsY8MYMjJuB+BJAF2FEMf0j94A5gPoLoS4AKC7/mdIKU8DWA/gDICdAJ6WUmbrtzUZwBfI3al3EcAO/fJwAFWFEPEAnof+yAwpZTKANwAc0T/m6peZza1bt1C5cmU4OTnh559/hq+vrzk3T+TwWMaGKfWkDyllBIqeuwWAIvduSSnnAZhXxPIYAP+6aKmUUgdgcDHbWg5geWk5TfX333/D2dkZP/30E6ZPn46JEyeicePGlvo4IoeTkpLCMjaAw58O/cADD0Cr1SImJgZVqlRBr1690KVLF2zfvh05OTmq4xHZPI6MDePwZZwnMDAQr7/+Oi5fvoyJEydizpw5aNq0KTZu3MhSJioDjUYDnU6nOobVYxkXotFoMGzYMMTGxmL+/PmYP38+WrRogb1796qORmSTHnjgAfz222+qY1g9lnExhBDo06cPoqOj8eKLL2Ls2LEYOHAgfv/999JXJqJ89evXx/nz51XHsHq8alspKlSogCFDhuCRRx7BggUL0Lx5c7z22mt46qmnUKEC/y6jf2RnZ2P69OnIzMyEm5sb3N3d4ebmds/zopa5ubmhcuXKCAgIUP0VLKJevXr4+eefVceweiL3cF77ERYWJmNiYiy2/bi4OIwdOxYuLi5YuXIl/P39LfZZZFuSk5NRtWpVLF68GFqtFunp6UhLS0NaWlr+c51Ol79Mp9PlLz958iQSEhJQu3bt0j/Ixhw6dAjPPvssoqOjVUcpEyFErJQyzGLbZxkbLzs7GwsXLsT777+Pzz77DI899phFP49sw7Vr1xAaGorExESj133++ecRHx+Pb775Bu7u7hZIp86tW7cQFBSE27dv2/QV2yxdxvw92wROTk546aWXsH37djz33HOYMWMGsrOzS1+R7JpOpzO5SN9++21UrlwZTZo0wZo1a2BPg6SqVavCyckJN2/eLP3NDoxlXAatWrVCTEwMYmJi0K9fP9y5U9L1k8jeCSFw6dIlDB5c5PlLJXJ1dcWqVauwdOlSLFy4EGPHjrWrQm7RogV27NhR+hsdGMu4jKpVq4adO3fC398fHTt2xLVrFrmOEdmAOnXqYMeOHThx4oTJ2+jWrRsiIyMRGxuLrVu3mjGdWi+99BLmzZuHrKws1VGsFsvYDDQaDT799FM8/vjjaN++PY+pdGDBwcG4fPkylixZAq1Wa9I2PDw88Morr2Dp0qVmTqdO586dUbNmTaxZs0Z1FKvFMjYTIQRefvllTJ8+HZ07d8aFCxdURyIFgoODsX//fvz000+oU6cO3nzzTZN26PXq1QsRERF2dfbna6+9hjfeeIOj42KwjM1s8uTJmD17Nrp164bLly+rjkMKtGrVClu3bsWePXtw+fJlNGzYEP379zfqLM777rsP9913HxISEiyYtHxxdFwylrEFjB8/HtOnT0f37t25B9mBNWrUCF988QV+//13PPLII3j88cdx7Ngxg9fv2LEjfvzxRwsmLH8cHRePZWwhU6dOxdChQ9GvXz+T5w7JPnh6emLcuHGYPXs2Zs+ebfB6vXv3xp49eyyYrPxxdFw8lrEFzZ07F/Xq1cPo0aPt6jAlMs348eNx/vx5fPXVV6W+V6fTYfv27XZ5ijRHx0XjGXgWptPp0KlTJzz22GN46aWXVMchxX799Vd06dIFf/31V7HvuXnzZv4dZ3r27AkfHx94enrC09MTFStW/NfzTp06oUqVKuX1Fcyic+fOGDt2LEaOHKk6isEsfQYeLxRkYW5ubti0aRPCwsLQpk0bdOrUSXUkUsjX1xdOTk4lvqdatWrYvXs3srOzkZqairt37yI1NTX/kZiYmL/88OHDmDJlCqZNm1ZO38A85s6diyFDhiA4OBht27ZVHccqsIzLQa1atfDll19ixIgROH78uM2NYsh8qlevDq1WW+Lt6ytUqICHH37YoO1NnTq11HK3Rh07dkR4eDgeffRRvPPOOxg9erTqSMqxjMvJf/7zHwwcOBBTpkzhzgs7dv78eXz88cf5d7eoVq0aXFxc8h8ZGRnQ6XRmuw1RZmYmnJ1t849x7969sW/fPjzyyCM4ceIE3nnnHZv9LubAHXjl6K233kJMTIxdneZK99q/fz+ioqLg5+cHDw8PSCmRkpKCxMRExMfHIy4uDrNmzTLb1cuysrKg0WjMsi0VGjRogOjoaJw8eRJ9+/YtcS7d3jnuX0MKeHh4YNmyZRgxYgS6du2KSpUqqY5ExVi/fj3efPNNODs7Q6PRwNnZGc7OzujUqRPmzp1b7HqZmZlo1qwZXnjhhXLJacsj4zze3t7YsWMHXnjhBbRq1Qrbtm1D/fr1Vccqd7b9X9EGderUCd26dcMbb7yBhQsXqo5Dxci7xsTOnTvh5OSErKwsHD16FOvWrSuxjLOyssq1HDMzM216ZJzH2dkZH374IRo3bowOHTpg1apV6NWrl+pY5YrTFAosWLAAX375JeLj41VHoWKEhISgffv2iIuLQ6tWrdCuXTs0bty41KIt73Is7/K3tHHjxmHz5s0YO3Ys3n//fYc6Pt9+/ivaEF9fX0ybNg2zZs3CunXrVMehYvTo0QO//vorxowZAyC3aO/cuYNffvklf4ecEAKenp75PycmJnJkXEbt27fH4cOH83fsLV26FG5ubqpjWRzLWJHnnnsOwcHBOHHiBJo0aaI6TrnIzMxETEwMXFxc4O7uDg8PD7i7u8Pd3R2urq7QaDRWdZPXoKAgfP/99/k/16lTB/fffz/mzJmDjIwM3Lp1Czk5OahQoQIyMjKQkZGBxMREfPTRR+WW0R7LGAD8/f0RGRmJ0aNHo2XLlnj33XfRo0cP1bEsimWsSMWKFfHCCy9g3rx5+Oabb1THKRf79u1D9+7dERISAq1Wi7S0NGi1Wmi1WmRkZCAzMxNOTk7QaDTQaDTIyMjA33//jTFjxuDw4cNwcXGBRqOBi4sLXF1d8x9ubm73/Fx4WVHPXVxc4Obmhk6dOhV7mFnnzp3x2GOPIS0tDe7u7ggODsauXbvK+d9ayey1jIHcPyPr16/Hli1bMGXKFAQEBOCdd95B8+bNVUezCJaxQv/9738xf/58xMfHIzg4WHUci3N2dkbHjh2xb9++Il+XUiIrKwuZmZnIzMyEt7c3MjMzcezYMXz88ceoU6dO/gg0PT09/5F3l+XilqWmpt7zc97zXbt24YMPPsB///vfIvN4eXkhOzsbmzZtwogRIyz5r8Zktn5oW2mEEBgwYAD69euHL774Ar1798bDDz+MN9980+6u28EyVsjT0xMTJkzA4sWLsXjxYtVxLM7NzQ06na7Y14UQ+aNiAHBxcYFOp0NWVhYCAwPx4IMPmjXPa6+9hqtXr5b4nvDwcKxcudJqy9geDm0zhEajweTJkzFixAi89957CAkJwZgxY/Dyyy/bzRmt1jNB56AmT56M1atX4+7du6qjWJyrq2uJZVxYYGAgAgICkJiYiPvuu8/seWrUqFHqPQufeOIJHDlyBDdu3DD755tDRkYGXFxcVMcoN15eXnjttddw6tQppKamon79+li4cKFR/19ZK5axYrVr10b79u2xfv161VEszs3NDenp6Qa//9y5c7h79y5SUlJQo0YNs+epXr36PbdEunr1KpYvX44xY8agQYMG+OWXX+Dm5obevXtj8+bNZv98c0hPT4erq6vqGOWuevXqWLp0KQ4cOICDBw+ifv36WLVqFbKzs1VHMxnL2AqMGTMGK1asUB3D4kqbpihvlSpVwo0bN3Dy5EksWrQItWrVwosvvoiwsDCMHz8ekydPhk6nw2OPPWa1O1kdtYzzPPjgg9iyZQvWrFmDpUuXIiQkBNZ0CV2jSCnt6hEaGiptjU6nk97e3vL3339XHcWirl27Jn19fVXHyHfp0iXZrFkz2ahRI+nn5yefeeYZmZOTI6WUMicnRzZt2lRWqFBBApAA8l+zJkFBQfLChQuqY1iFnJwcuWbNGunj4yNXrlxp9u0DiJEW7C77n/m3Aa6urnjkkUewZcsWTJ06VXUci3F1dbWqufHAwED8+uuvRb4mhDDqfnWqpKenO9SccUmEEBg2bBiaNm2KRx99FDExMXjvvfds5mgTTlNYif79+2Pbtm2qY1iUm5sbUlNTVcewK44+TVGUhg0bIjo6GhcvXsTDDz9sMzcFZhlbiYcffhiHDx+267Jyc3ODEMKhrjdgaSzjolWuXBnbt29Hx44dERYWhqNHj6qOVCqWsZXw8vJCaGgoDhw4oDqKxVSoUAHOzs7IyMhQHcVuZGRksIyLUaFCBbzxxht477330KdPH6u/MBfL2Ip07drV7m7NXpirq6tRh7dR8aSUnDM2wODBgzFnzhz06dMHycnJquMUi2VsRTp06ICIiAjVMSwq78I6VHZ5t7o/duwYLly4gOvXr0Or1XIaqAiTJk1Cv379MHDgQKv9zUzY23+4sLAwaavHGd69exf3338/kpOT7fZXTycnJ+h0OpvZw23NsrOz0b9/f1y7dg0pKSm4c+cOUlJSkJGRAS8vL3h5eaFSpUr3/DPv0bhxY0yYMAHp6enYunUr3Nzc/vVePz8/s90eyhpkZ2dj0KBBuO+++/Dll18a/d2EELFSyjALxWMZW5smTZpg+fLlCAuz2H9zZXQ6Hdzd3ZGTk2NXf8itTVZWFlJSUu4p6Dt37uQ/P3XqFJYtW4aUlBQcPHgQjzzyCNq0aXPPOvHx8QgPD8fYsWMtkvHcuXN49tln4ePjAw8Pj2IfFStWLPF1d3f3Uq/NIaXMv8DU7du30bZtW/Tp0wdLly416v9DS5cxjzO2Ms2bN8exY8fssoxdXV0hhMDVq1dRq1Yt1XHslrOzM7y9veHt7V3k65cuXco/vVun06Fx48bYvn37Pe+ZOHEiMjMzLZYxKSkJP/74I1atWgWtVou7d+9Cp9NBq9Xi1q1b+cvyLrWa97zwpVe1Wi2cnZ3zy9nJySm/ePMeeZcZzbt0anJyMj7//HPMnj0bNWvWtNh3NBbL2Mo0btwYp06dUh3DIoQQaNSoEW7cuMEyVkij0eDKlSt49tln0a1btyLvoqHT6Sx6d428S8YOHz68TPsQ8ka9ecWclZWVX7p5D41GYxO/iXFPipV58MEHce7cOdUxLKZx48Y4c+aM6hgOrVatWjh27BjOnDmDxx57rMgytPTxy87OzrjvvvuQkJBQpu0IIeDq6gpvb2/UrFkTAQEB8PPzQ5UqVfJvh2ULRQywjK1O3bp1rf54yLJwc3NDWlqa6hgOr2nTpti+fTtmzJiBxx9//F+v63Q6i5axj48Pnn/+ebz++usW+wxbw2kKKxMQEIDff//dbg8Bu3jxotVeqN3RuLm54a233irytbS0NHh4eJj9M7VaLT788EO0bdsWLVu2RGRkpNk/w1bZ3592G5d3iFFSUpLqKBZx9+7dYu85R9Yj775/5hYXF4dXXnkFixYtQlpaGk9YKYBlbIV8fX2t9s4SZSGlRHJyMsvYBlhqB17e2ZceHh44ceIEGjZsaPbPsFUsYyvk6elplxcMWrp0KapUqWL2e9mR+VmqjPNuLuDm5pZ/dTXKxTK2Qn/99RcqV66sOoZZ/frrr5g9ezbWrFljl3Ph9ibvBB1zyxsZ513D+8SJE7hy5YrZP8cWcQeeFUpNTYWXl5fqGGZz4cIF9OvXD59++inq1aunOg4ZwBwj45ycnPyTN1JTU+Hr65s/Mt64cSMCAgKQnJyMCxcuoEqVKgadTWfPHPebWzF7ukbtjRs30LNnT8yZMweDBg1SHYcMZMyhbeHh4Vi+fDlSU1ORmpqaX75arRbu7u6oWLEisrKy8OCDD+Kzzz4DkHsG3ksvvQQAGDBgAKSU+WfTVa5cGSdPnsT9999vse9njVjGVkhKaTMHqpdm8eLF6NWrFyZMmKA6Chno0qVL+PPPPw0eGe/duxfdunXDgAEDULFiRXh6esLT0zP/9GQAyMzMROXKlREYGFjsVeXyzqYLCgpyyGPROXlnhZycnGz6luMFXbx4EW3atFEdg4ywbNkyADD4OOO861s0a9YMdevWRfXq1eHl5ZVfxEDuKdiNGjXC8ePHi91O3tl0mZmZdvOboTFYxlbIw8PDLkYGUkpERkaiZcuWqqOQERo2bIiuXbsaPH9r6Pxys2bNSizjPI569xKWsRWqVKkS/v77b9UxyiwuLg5CCO60szE9evQw6o4zhpZx06ZNDSpjR717CeeMrVDVqlXx559/qo5RZrt370b37t3tZv7bUbi4uKBSpUoGv9+YkfG0adPg4uJS4nWK7WkHtjFYxlaoevXquHbtmuoYZbZu3Tq88MILqmOQkYwtQ0NPnW7ZsiU++eST/OsWa7Va3Lx5855rE2u1Wjz++OMOeYib431jGxAQEGDzB8Lv378fiYmJ6Nevn+ooZCRjy9jQkbFGo+FRNSXgnLEVCg4OtunLaJ44cQJDhgzBe++955AjHFtn7OUzLX0hekfBMrZCDRo0wOnTp1XHMEl8fDx69OiBRYsW4bHHHlMdh0xgqZExlYxlbIUeeughnD17Nv9W7LYiKysLTz75JGbOnFnkBcvJNuSVsaE3K7bUdSwcDcvYCnl5eSEgIAAnT55UHcUo8+fPh6enJ6ZOnao6CpWBh4cHTp48iQoVKkAIASEEKlSoACcnJzg7O0Oj0cDFxQWurq75d26pWLGi6tg2jxN6Vqpt27aIiIhA8+bNVUcxSExMDD766CMcPXqUV2WzcQ0b3EMsOgAADFFJREFUNsy/M7SU8l+PwsudnJyg0WhURrYL/FNjpTp37oy9e/eqjmEQrVaLESNGYPHixVZ163Mqu6JGxYVHxixi82AZW6mHH34Ye/fuzR+hWLMZM2YgNDQUQ4YMUR2FyGaxjK2Un58fgoODceDAAdVRSrRr1y58++23WLJkieooRDaNZWzFHnvsMWzatEl1jGLdunUL48aNw5dffglvb2/VcYhsGsvYig0ZMgQbNmywyqkKKSUmT56MwYMHo1u3bqrjENk8lrEVCwoKQr169fDDDz+ojvIva9aswenTp/HWW2+pjkJkF1jGVm7ChAlYunSp6hj3SEhIwHPPPYfVq1fzYH8iM2EZW7nHH38csbGxiIuLUx0FQO5NJkePHo1p06bZzDHQRLaAZWzl3N3dMWnSJLz//vuqowAAPvzwQ2RkZGDGjBmqoxDZFWHo+ee2IiwsTMbExKiOYVZ//vkn6tevj+PHj6NWrVrKcpw6dQpdunRBVFQUHnjgAWU5iFQQQsRKKcMstX2OjG1AtWrVMH78eLzxxhvKMqSnp2PEiBFYsGABi5jIAljGNmLGjBnYvHkzzp49q+TzZ8+ejTp16mDMmDFKPp/I3rGMbUSVKlUwc+ZMPPfccwZf2tBc9u/fj6+++gqff/4572dHZCGllrEQYrkQ4qYQ4lSBZc2EEIeFEMeEEDFCiJYFXpsphIgXQsQJIXoWWB4qhDipf22x0P+pFkK4CiG+0S+PEkIEFlhnlBDigv4xylxf2lZNnToVv//+e7melXfnzh2MGjUKn3/+OXx8fMrtc4kcjSEj4xUAehVa9g6A16WUzQDM1v8MIURDAEMBNNKv84kQwkm/zqcAJgKoq3/kbXMcgNtSymAAHwBYoN9WFQBzALQC0BLAHCGEQ59zq9Fo8Nlnn+GZZ57BrVu3yuUzn332WfTo0QN9+/Ytl88jclSllrGUcj+A5MKLAeTdy/s+AHm3Mu4PYJ2UMl1KeQlAPICWQojqACpJKQ/J3N+xVwF4tMA6K/XPNwLoph819wSwW0qZLKW8DWA3/v2XgsNp3749nnjiCYwbN87i0xWbN29GREQE3nvvPYt+DhGZPmc8DcBCIcTvAN4FMFO/vCaA3wu87w/9spr654WX37OOlDILwN8AqpawrX8RQkzUT5fEJCUlmfiVbMdbb72FGzduYN68eRb7jOvXr+Opp57CV199BU9PT4t9DhHlMrWMJwN4TkpZG8BzAML1y4vauyNLWG7qOvculPJzKWWYlDLMEeY1XVxcsHnzZixbtgwrV64sfQUjSSkxduxYTJw4Ea1btzb79ono30wt41EANuufb0DunC6QO3qtXeB9tZA7hfGH/nnh5fesI4RwRu60R3IJ2yIA1atXx86dOzFz5kysXr3arNv+7LPPcPPmTbz66qtm3S4RFc/UMr4GoJP+eVcAF/TPtwEYqj9Cog5yd9RFSykTAaQIIVrr54NHAvi2wDp5R0oMArBHP6+8C0APIYS3fsddD/0y0mvQoAF2796NmTNn4r333jPLHPL58+cxa9YsrF69mrfTISpHpd6QVAixFkBnANWEEH8g9wiHCQAW6UeyOuQeJQEp5WkhxHoAZwBkAXhaSpmt39Rk5B6Z4Q5gh/4B5E5xfCWEiEfuiHioflvJQog3ABzRv2+ulLLwjkSH16hRI0RGRqJfv36IiorCmDFj0KVLF7i5uRm9raysLDz55JN47bXX8OCDD1ogLREVh9emsBNarRYffPABdu3ahVOnTmH48OGYPn06/P39Dd7G3LlzcfDgQezYsYMndxAVwmtTkEE8PDzwyiuvYP/+/Thx4gQqVqyI5s2bY+7cuQbdKSQ6Ohoff/wxli9fziImUoBlbIdq1aqF+fPn49ixYzh06BC6dOmCkg75u3v3LkaMGIElS5agRo0a5ZiUiPKwjO1Y7dq18f3336NTp07o0KEDrl+/XuT7XnzxRbRq1QqDBw8u54RElKfUHXhk2ypUqIB58+bBxcUFffv2xYEDB+65VdKOHTvw3Xff4fjx4wpTEhFHxg5i9uzZCAoKuucOHbdu3cL48eOxYsUKVK5cWWE6ImIZOwghBJYuXYoNGzYgNjYWUkr897//xbBhw9ClSxfV8YgcHqcpHIi3tzdeffVVzJkzB48//jji4uLMfvYeEZmGxxk7GJ1Olz9nfOzYMTRt2lRxIiLbYOnjjDkydjBubm54+umnUa1aNRYxkRVhGTugJUuWqI5ARIVwBx4RkRVgGRMRWQGWMRGRFWAZExFZAZYxEZEVYBkTEVkBljERkRVgGRMRWQGWMRGRFWAZExFZAZYxEZEVYBkT/b+9cwuxqgoD8PfbpKJRjmaQmpeBhJJM0jShwpduVlj0kKkg2VNqUdCFaaCUMNCCIIq0ByXIyhsFQTWhDxbd1dTG1JzMcmrARqkwI7T+HtY/zjqbubhzztnrTP8Hm1ln3fY366z9733W3jPHcRLAg7HjOE4CeDB2HMdJAA/GjuM4CeDB2HEcJwE8GDuO4ySAB2PHcZwE8GDsOI6TAH3u26FF5Bfgh7Ps5kKgrRd0KkG1uFaLJ7hruah21zGqOrxcO+xzwbg3EJFt5fxK7t6kWlyrxRPctVy4a/f4MoXjOE4CeDB2HMdJAA/GnfNK0QI5qBbXavEEdy0X7toNvmbsOI6TAH5l7DiOkwAejB3HcVJAVat6A1YDR4CmKO9ZYB+wG3gLGJJpMxo4DjwS5U0GvgaagRfoWMIZAKyz/M+BsVGb+cAB2+ZH+eOs7gFr2z+vKzAW+BPYadvKSrnmHVNgIvApsMe8BiY6pnOj8dwJ/ANMStT1XOBVc9oL1Cc8V/sDa8xpFzCj4Ln6tDnuBD4ARkRl9bbP/cBNlRzTHmNZ0cH0bDfgeuCqzJtxI1Bj6eXA8kybTcAGSoPxF8B0QID3gFssfyEWCIHZwDpLDwUO2s9aS9da2XpgtqVXAvfndSUE46Yufueyuub0rLGJf6W9Hgack+KYZsbwCuBgwu//HOBNSw8CDmGBIEHXRcAaS18EbAf6FThXz4/SD0b7uZxwshhACJjfUcG52mMsKyKA9vZG94HrTmBt9PoOwhl+CRaMgYuBfVGde4BVlm4EpkeBp83esNN1rGyV5YnVaZ+004HGvK5d1auUaw7PmcBrRXnmff+j/GeAZam6Wl/v2D6HAd8SDvwUXV8C5kVlW4CplXLtwbMeeDlKx58wGq2fio1pd9v/Yc14AeFMh4gMBh4HlmbqjARaotctltdedhhAVU8BvxEOjtP5mTbDgF+tbravM3Y1xonIVyKyVUSuS8g19hwPqIg0isgOEXksIc+sa8zdwBsJu24E/gBagR+B51T1WKKuu4BZIlIjIuMIH/kvKdJVRJaJyGHC0tST2X1m+khiTGvOpFK1IiINwClgrWUtBZ5X1eMiUlK1k+baQ1ne/LyurcBoVT0qIpOBt0VkQtGunXjWANcCVwMngC0ish34vUjPLlzb86cBJ1S1qQefIl2nAn8DIwgfgz8Skc2Juq4GLgO2Ef4vzCdWXpirqjYADSJSDywGnvqPfVdkTKEPP00hIvOB24C5ap8XgGnAChE5BDwEPCEiiwlnr1FR81HAz5ZuIZzlEZEa4ALgWJyfadMGDLG62b7O2FVV/1LVo5beTljfGl+kaxdj2gJsVdU2VT0BvEtYw0tuTCNm03FV3O6Tmusc4H1VPamqR4CPgSkpuqrqKVV9WFUnqeosYAjh5lWhrsbrwF3ZfWb6SMGzb64ZAzcD3wDDu2mzhNIbeF8C19CxgD/T8hdRuoC/3tJDge8JVy21lh5qZRsoXcBfmNcVGE7HzYU64Keo/7K75vCsBXYQbjLVAJuBW1McUyvrRziQ6jL5SbkSltPWmM9gqzMxUddBwGBL3wB8WMlx7cTz0ij9ALDR0hMovYF3kI5jrCJj2m0cK1eArNRGuMJpBU4SDrL7CI+hHKaTx8KidksoDcZTgCbCFeiLdDzaMtAGt5lwx7UuarPA8puBe6P8OqvbbG0H5HUlnM332OTZAdxeKde8YwrMM9cmYEWqY2r1ZwCfdTIfknIFzrO2ewgB8NGEXccSHhXbSzgZjyl4rm6yfe4m3AQdGfXRYC77sScmKjWmPW3+59CO4zgJ0GfXjB3HcaoJD8aO4zgJ4MHYcRwnATwYO47jJIAHY8dxnATwYOw4jpMAHowdx3ES4F9icUZGTVZJPwAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 360x864 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "fig, ax = plt.subplots(figsize=(5,12))\n",
    "city_of_interest = 'Seattle'\n",
    "single_city = city_df[city_df['CITYNAME']==city_of_interest].copy()\n",
    "\n",
    "single_city.plot(color='none', \n",
    "                   edgecolor='black', \n",
    "                   legend=True,\n",
    "                  ax=ax)\n",
    "\n",
    "_ = single_city.apply(lambda x: ax.annotate(s=x['CITYNAME'], \n",
    "                                           xy=(x.geometry.centroid.x, \n",
    "                                               x.geometry.centroid.y), ha='center'),axis=1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2019-10-24T13:58:08.759881Z",
     "start_time": "2019-10-24T13:58:08.514775Z"
    }
   },
   "outputs": [],
   "source": [
    "single_city.to_crs({'init': 'epsg:4326'}, inplace=True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2019-10-24T13:58:08.878453Z",
     "start_time": "2019-10-24T13:58:08.763991Z"
    }
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "[-122.45968008362588,\n",
       " 47.47334752062687,\n",
       " -122.21723116846032,\n",
       " 47.73424598428268]"
      ]
     },
     "execution_count": 6,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "seattle_bounds = single_city.bounds.values[0].tolist()\n",
    "seattle_bounds"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2019-10-24T13:58:09.087542Z",
     "start_time": "2019-10-24T13:58:08.882317Z"
    }
   },
   "outputs": [],
   "source": [
    "upperleft_lon = single_city.bounds.minx.values[0]\n",
    "upperleft_lat = single_city.bounds.maxy.values[0]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2019-10-24T13:58:09.361287Z",
     "start_time": "2019-10-24T13:58:09.091886Z"
    }
   },
   "outputs": [],
   "source": [
    "upperleft_geohash = geohash.encode(upperleft_lat, upperleft_lon, precision=9)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## test the function to draw the bounding box"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2019-10-24T13:58:09.515042Z",
     "start_time": "2019-10-24T13:58:09.365071Z"
    }
   },
   "outputs": [
    {
     "data": {
      "image/svg+xml": [
       "<svg xmlns=\"http://www.w3.org/2000/svg\" xmlns:xlink=\"http://www.w3.org/1999/xlink\" width=\"100.0\" height=\"100.0\" viewBox=\"-122.4597089767456 47.73422069549561 4.634857177165941e-05 4.634857177165941e-05\" preserveAspectRatio=\"xMinYMin meet\"><g transform=\"matrix(1,0,0,-1,0,95.46848773956299)\"><path fill-rule=\"evenodd\" fill=\"#66cc99\" stroke=\"#555555\" stroke-width=\"9.269714354331881e-07\" opacity=\"0.6\" d=\"M -122.45970726013184,47.73426532745361 L -122.45970726013184,47.734222412109375 L -122.4596643447876,47.734222412109375 L -122.4596643447876,47.73426532745361 L -122.45970726013184,47.73426532745361 z\" /></g></svg>"
      ],
      "text/plain": [
       "<shapely.geometry.polygon.Polygon at 0x7faf5df32160>"
      ]
     },
     "execution_count": 9,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "hash_bounds = geohash.bounds(upperleft_geohash)\n",
    "n_y = hash_bounds.ne.lat\n",
    "w_x = hash_bounds.sw.lon\n",
    "s_y = hash_bounds.sw.lat\n",
    "e_x = hash_bounds.ne.lon\n",
    "nw_pt = [w_x, n_y]\n",
    "sw_pt = [w_x, s_y]\n",
    "ne_pt = [e_x, n_y]\n",
    "se_pt = [e_x, s_y]\n",
    "Polygon([nw_pt, sw_pt, se_pt, ne_pt, nw_pt])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2019-10-24T13:58:09.789901Z",
     "start_time": "2019-10-24T13:58:09.518969Z"
    }
   },
   "outputs": [],
   "source": [
    "#modified from source: https://blog.tafkas.net/2018/09/28/creating-a-grid-based-on-geohashes/\n",
    "def is_geohash_in_bounding_box(current_geohash, bbox_coordinates):\n",
    "    \"\"\"Checks if the box of a geohash is inside the bounding box\n",
    "\n",
    "    :param current_geohash: a geohash\n",
    "    :param bbox_coordinates: bounding box coordinates\n",
    "    :return: true if the the center of the geohash is in the bounding box\n",
    "    \"\"\"\n",
    "\n",
    "    coordinates = geohash.decode(current_geohash)\n",
    "    geohash_in_bounding_box = (bbox_coordinates[1] < coordinates[0] < bbox_coordinates[3]) and (\n",
    "            bbox_coordinates[0] < coordinates[1] < bbox_coordinates[2])\n",
    "    return geohash_in_bounding_box\n",
    "\n",
    "def build_geohash_box(current_geohash):\n",
    "    \"\"\"Returns a shapely Polygon for a given geohash\n",
    "\n",
    "    :param current_geohash: a geohash\n",
    "    :return: a list representation of th polygon\n",
    "    \"\"\"\n",
    "    hash_bounds = geohash.bounds(current_geohash)\n",
    "    n_y = hash_bounds.ne.lat\n",
    "    w_x = hash_bounds.sw.lon\n",
    "    s_y = hash_bounds.sw.lat\n",
    "    e_x = hash_bounds.ne.lon\n",
    "    nw_pt = [w_x, n_y]\n",
    "    sw_pt = [w_x, s_y]\n",
    "    ne_pt = [e_x, n_y]\n",
    "    se_pt = [e_x, s_y]\n",
    "    return Polygon([nw_pt, sw_pt, se_pt, ne_pt, nw_pt])\n",
    "\n",
    "def compute_geohash_tiles(bbox_coordinates, geohash_precision):\n",
    "    \"\"\"Computes all geohash tile in the given bounding box\n",
    "\n",
    "    :param bbox_coordinates: the bounding box coordinates of the geohashes\n",
    "    :return: a list of geohashes\n",
    "    Geohash length\tCell width\tCell height\n",
    "    1\t≤ 5,000km\t×\t5,000km\n",
    "    2\t≤ 1,250km\t×\t625km\n",
    "    3\t≤ 156km\t×\t156km\n",
    "    4\t≤ 39.1km\t×\t19.5km\n",
    "    5\t≤ 4.89km\t×\t4.89km\n",
    "    6\t≤ 1.22km\t×\t0.61km\n",
    "    7\t≤ 153m\t×\t153m\n",
    "    8\t≤ 38.2m\t×\t19.1m\n",
    "    9\t≤ 4.77m\t×\t4.77m\n",
    "    10\t≤ 1.19m\t×\t0.596m\n",
    "    11\t≤ 149mm\t×\t149mm\n",
    "    12\t≤ 37.2mm\t×\t18.6mm\n",
    "    \"\"\"\n",
    "\n",
    "    checked_geohashes = set()\n",
    "    geohash_stack = set()\n",
    "    geohashes = []\n",
    "    #geohash_precision = 7\n",
    "    # get center of bounding box, assuming the earth is flat ;)\n",
    "    center_latitude = (bbox_coordinates[1] + bbox_coordinates[3]) / 2\n",
    "    center_longitude = (bbox_coordinates[0] + bbox_coordinates[2]) / 2\n",
    "\n",
    "    center_geohash = geohash.encode(center_latitude, center_longitude, precision=geohash_precision)\n",
    "    geohashes.append(center_geohash)\n",
    "    geohash_stack.add(center_geohash)\n",
    "    checked_geohashes.add(center_geohash)\n",
    "    while len(geohash_stack) > 0:\n",
    "        current_geohash = geohash_stack.pop()\n",
    "        neighbors = geohash.neighbours(current_geohash)\n",
    "        for neighbor in neighbors:\n",
    "            if neighbor not in checked_geohashes and is_geohash_in_bounding_box(neighbor, bbox_coordinates):\n",
    "                geohashes.append(neighbor)\n",
    "                geohash_stack.add(neighbor)\n",
    "                checked_geohashes.add(neighbor)\n",
    "    return geohashes"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2019-10-24T13:58:15.264093Z",
     "start_time": "2019-10-24T13:58:10.014225Z"
    }
   },
   "outputs": [],
   "source": [
    "geohash_precision = 6\n",
    "geohash_list = compute_geohash_tiles(seattle_bounds, geohash_precision)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2019-10-24T13:58:15.274237Z",
     "start_time": "2019-10-24T13:58:15.268238Z"
    }
   },
   "outputs": [],
   "source": [
    "seattle_geohash_df = pd.DataFrame(geohash_list, columns=['geohash'])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2019-10-24T13:58:18.709814Z",
     "start_time": "2019-10-24T13:58:15.294117Z"
    }
   },
   "outputs": [],
   "source": [
    "seattle_geohash_df['geometry'] = seattle_geohash_df.apply(lambda x: build_geohash_box(x['geohash']), axis=1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2019-10-24T13:58:18.725271Z",
     "start_time": "2019-10-24T13:58:18.714124Z"
    }
   },
   "outputs": [],
   "source": [
    "seattle_geohash_gpd = gpd.GeoDataFrame(seattle_geohash_df[['geohash']], crs={'init':'4326'}, geometry=seattle_geohash_df['geometry'].tolist())"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2019-10-24T13:58:20.995531Z",
     "start_time": "2019-10-24T13:58:18.742219Z"
    }
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<matplotlib.axes._subplots.AxesSubplot at 0x7faf5bccea58>"
      ]
     },
     "execution_count": 15,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAc0AAAG2CAYAAAAdoN0fAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8QZhcZAAAgAElEQVR4nO3dd5xU1fnH8c+hSREFBBUERVHQ2CirqERiQWzYFRu2mKDGFiyo2HuLsf4STYgNxa7BhoqFYncRO9gQFTWAIkqRtnt+fzwz7rLOstPP3Lnf9+s1rztlzzzPXoZ95t57ivPeIyIiIg1rFDoBERGRqFDRFBERSZOKpoiISJpUNEVERNKkoikiIpKmJqETyET79u19165dQ6chIiIlZPLkyd977zsUI1akimbXrl2prKwMnYaIiJQQ59yXxYql07MiIiJpUtEUERFJk4qmiIhImlQ0RURE0qSiKSIikiYVTRERkTSpaIqIiKRJRVNERCRNKpoiIiJpUtEUERFJk4qmiIhImlQ0RURE0qSiKSIikiYVTRERkTSpaIqIiKRJRVNERCRNKpoiIiJpahI6gRCcc6FTEBGJvbZt2zJ37tzQaWQklkUTwHufVTvnXNZtc2n/6qvQr58Dso2dS9tc2yt2tGKvvK1z0KgRNG78223jxvDDD46OHf0KzzVpsuLP1L01amS3SZMcO+7of33cuHFNvPputb8D33+/45BDGv69U31vvu8+x6GHrti2USMYMgR22y2NvZbD34ZQbUshdtS4XH7hYquoqPCVlZU5v0/oD0mm7ZcuhYoKeP99R+fOnj/8wZ6v/TYN3X/gAcfgwT7jdsn7jz3m2Hdfn3a72s89+aRjzz19gz9X33PjxjkGDPht7HTaTpjg6N//t/u7ujq93+O11xx9+9Yfu74cvIcpUxy9emX/WUm2r/t3ZWV/Z7y33+2ddxybbWb7vKqKX7fV1TU/k3xc9zZ3rmPVVf0Kz9Vu27Dy+qLRsSPMnGkFdKWtVTSL3r7W+0z23lfk/EbpxFLRLF7bbNuffz5cdhmAY+FCT8uWxYmbr/aKHa3YK2ubLLh1i2lVVc1tzTUd33zjf328fHnqn6t9Sxb1HXZwPP+8X+H9kzFrF/7at9oOO8wxevTKf+/6dsvhhzvuvXfFF88+G77+2s70bLtt9vutIXH8nOWjfa33UdFMJY5F8+23Yeutk380opO3Ykc3dlTzLkTs006D66+HM86Aa68tbuxitI1y7DrvU7Siqd6zJcx7GDbMvl2fckrobETiZ//9bfvoo+mempZyp6JZwl5/HSZOhLZt4eKLQ2cjEj/bbgtrrQXTp8O774bORkpBbHvP5tJrK9ceX5m2//FHaNMm99jFzluxoxs7qnkXMnavXuFiF7JtyNhR7D0b26IZhXP4/fpZB4THH4e99orutQfFjlbsqOZdqNjjxsHAgbDppvDBB8WNXei2pRA7anR6tkTNnw9vvmnj1JJDTESk+HbYwc70fPghfPxx6GwkNBXNEjVxonXV32orWG210NmIxFfTprD33nb/0UfD5iLhqWiWqBdesO3OO4fNQ0RW7EUr8aaiWaJUNEVKx8CB0KoVVFbCV1+FzkZCUtEsQbNnw3vvQfPmDc9CIiKF16IF7LGH3dfRZrypaJagiRNt26+fFU4RCU+naAViPOQkCmOiXnjht5NzRyFvxY5+7KjmXYzYkybVP2m+/r2K1zaU2BbNUh4Ttc8+NjbzvvvgkEPyEzvqY7kUOxptyz32oEHw1FNw661w3HHFjV2ItqUQO2p0erYEvf++bTffPGweIrKiAw6wrU7RxpeKZolZsAC++MLGhnXvHjobEaltr71swpEXX7TpLSV+VDRLzIcf2naTTaxwikjpaN/eZuhavhyefDJ0NhKCimaJSZ6a3WyzsHmISGrJXrSPPBI2DwlDRbPE6HqmSGnbbz/bPvusXU6ReFHRLDHJVRRUNEVKU6dONunI4sXwzDOhs5Fii+2Qk1IfEzVoUP5jR3Usl2JHq22cYh90ULjY+WobMnYUh5zEtmiW4pioWbNg7bVtVZN581JPbFCKeSt2ecWOat7FjD19OnTrBquuCnPm2Mxd+vfKrn3U6PRsCandCSiCnyWR2NhgA+jZ065pJhdXkHhQ0Swhup4pEh2aizaeVDRLiHrOikRHsmg+/LAmOogTFc0SoqIpEh2bbmrr3f78M1x7behspFhUNEvEY4/BW2/ZLEBbbBE6GxFJx+WX2/bGG8PmIcWjolkCPv4YjjrK7l95JbRpEzYfEUlP3762KtGiRaEzkWJxuXQXLraKigpfWVmZ8/tEsZuziEi5cc5RXV2dj/eZ7L2vyENKDdI4zQzlc0zU7Nlw6KG2YsImm8Cbb9q4r0LEjvpYLsWORtu4xj7iCLjnHsdRR3nuvLO4saO6z5Lto0anZwOZMMHGeb34InToYNc0V1YwRaR0XXSRbUeNgo8+CpqKFJiKZgCXXQY77QTffQfbbw9TpkCPHqGzEpFsdetm2+pquOCCsLlIYaloFtGsWbY9/3zwHs49144011knbF4ikh/Nm9uSYXnoeiElSkWzSF56yU7Hgi1kO3asHXE2ie1VZZHyc/LJtj333LB5SOGoaBZYVRVccgkMGAD/+5899847sOuuYfMSkfw76yxbcOG552D8+NDZSCHE9jgn5HI2nTtHcxkexY5P7KjmHTp2+/Y17XfcsXixo7rPoth7NrZFs9BdrGfMgO22s84+HTrAPffAwIEaRqDYpR87qnmXSuz5861j0Jw58MQT9a+Nm6/YUd9nUaPTswWwbJmNv/zuOyuc77xjBVNEyl/r1jBihN0/91zrUSvlQ0WzAK6+Gl5/Hbp0sW+anTqFzkhEiun446FzZ3jvPXjwwdDZSD6paObZl1/CFVfY/bvugnbtwuYjIsXXvDlceKHdP/98O/sk5UFFM8/OPBN++QUOPjjzjgAiUj6OPho22gg++4ysptaT0qSimUcvvQQPPQQtW2p9PZG4a9IELr3U7l98MSxeHDYfyQ8VzTxZvhxOOcXujxhh1zNFJN4OOsgmNfnmG/jnP0NnI/mgpcFERCSIsl4azDnXGKgEvvHeD3LOPQAkpxlvA8zz3ves06YH8ECtpzYALvDe3+Cca5d4rSswAxjsvf8x218kU/kclzRnDnTvDvPmwX//a4vSZtI+l9jFaKvY8Yod1bxLNbb3tjDDK6/A5ptbj/r11stf7Kjvs6jJ5PTsqcDU5APv/cHe+56JQvkI8GjdBt77j2v9TB9gEfBY4uWzgRe89xsBLyQeR9K551rB3HVX2Hvv0NmISClxzpYM23hjeP992HprePXV0FlJttIqms65zsCewMgUrzlgMHBfA2+zM/C59/7LxON9gLsS9+8C9k0nl1IzeTKMHGkX/W+4wf6DiIjUtv768NprsMsutvj8jjvCvfeGzkqyke6R5g3AcCDVyeftgVne+08beI9DWLGwruW9/w4gsV0zVSPn3FDnXKVzrnLOnDlpplsc3sPpp9v21FPtm6SISCpt2sDTT8Nf/gJLl8KQIXDeeZoxKGoaLJrOuUHAbO/95Hp+5FAaOMp0zjUD9gYeyjRB7/2/vPcV3vuKDh06ZNq8oJ55BiZMgLZt7cMvIrIyTZrA//0f3HwzNGoEl19uY7olOtI50uwH7O2cmwHcD+zknLsHwDnXBNifFTv7pLI78Lb3flat52Y55zom3qcjMDvD3IOqroZzzrH7I0bYt0gRkXScdJIdda62Gjz8sD33zTdhc5L0NFg0vffneO87e++7YqdYX/TeD0m8PACY5r2f2cDbpDoafRw4KnH/KGBM2lmXgPvug3fftfGYJ50UOhsRiZpdd7XrnBtsYI+33tr6SEhpy3VpsLrXKXHOdQJGeu/3SDxuCewCHFen7VXAg865Y4GvgINyzCUj+VoD7uuvoUWLMLGL2Vax4xU7qnlHOfa33zoqshxpGNV9FsUhJ7Gd3CCXcUk33eQ55RT43e9sFYPGjYsXO65juRQ7Gm0VO/vYRx/tf52j9rLL7LJPOjUldN6hYtd5n6JNbqBp9LKQnE/yiisyK5giIvW5/Xa45horlOedZ8NSrrvO1uNVD9vSkevp2ViaM8cWl9ZEBiKSL87ZKkndu8Phh1vP/AkT7LU11rAiutNOsPPOtnpKBM9slgUVzQzMrtW/96qr9KEVkfzbZx+YPh2eew5eeMFuX39tvWyTPW3XWaemgEpx6ZpmBk49FW66ybHnnp4nnyxu7JBtFTtesaOad7nG9h4+/xxefNEK6Isvwvffr9CajTbyvxbRHXaATIa0R3Wf1Xmfol3TVNFM0xdfQI8esGyZ4733PJtvXrzYodsqdrxiRzXvuMSuroYPP6wpoE884YAV2265pR2J7rKLDW1ptJLeK1HdZ3XeR0UzFS0NJiJSPlw5Lw1WbjL5svDuu9CrFzRtCkuXRuPbaD7bKna8Ykc177jGrtt28WKbNOGFF+CWW+Cnn6B/fxgzJvXMZaH3WdRoyEkazjnHriuccELoTEREVq55c+tpe9llMGkSdOoEEyfamp4zG5q7TRqkotmACRNg7Fho3drWzRQRiYrNN7ejzk02gQ8+gG23hY8+Cp1VtKloroT3cNZZdv+MMzLrkSYiUgrWXRdeftnGls+cCb//PbzySuisoktFcyX++1944w1Yc0047bTQ2YiIZKddO3j+eZuQ5ccfYcAA+/smmVPRrMfy5Tb3I8D558Oqq4bNR0QkFy1awCOPwNCh1lnogAPgtttCZxU9Kpr1uOsumDbNlu0ZOjR0NiIiuWvSBG69FS6+2MZ7Hn+8PR+hkYfBxXbISbpdnadPh1VWya5trrFLqa1ixyt2VPOOa+xc827UKJp5hxDbormysUXXXgvDh0PPnrYobO3ZNMplLJdiK3Yh2ip2tP69Hn8c9tnHZhS64AI7Ai1W7GT7qNHp2TrmzYMrr7T7V1218umnRESirPZKTSNHQlVVuFyiQiWhjquvtt5lO+4IAweGzkZEpPA22AC+/RbGjw+dSelT0azF+5reZFdeqaW/RCQehgyx7ahRYfOIAhXNWmbPtqPMNm1g661DZyMiUhzJovnII7BoUdhcSp2KZi2ffGLb7t11lCki8bHRRrDNNrBggSY9aIiKZi0ff2zb7t3D5iEiUmxHHGFbnaJdOa2nKSIiQWg9zQhJ9WVhn31s3NIDD8DgwanbxXEsl2LHK3ZU845r7Hzmve++tu7m3/8Ow4YVJ3bU6PRsLbWvaYqIxI1O0TZMRTNh+XL4/HO7v9FGYXMREQlhzz1t9MCUKfDhh6GzKU0qmgkzZsCyZdC5M7RqFTobEZHia97cCidozc36qGgmJE/N9ugRNg8RkZDatrXt4sVh8yhVKpoJs2bZtlOnsHmIiISUXNVpyZKweZQqFc2EJol+xJqwWETiTEVz5TROU0REgtA4zQip+2XhkUfgwANhv/3g0Ufrb6exXIpd7rGjmndcY+c776uugnPOgTPPhGuuKXzsqNHp2YTk6dlly8LmISISkk7PrpyKZkLTprZdvjxsHiIiISWL5tKlYfMoVSqaCcmiqSNNEYkzHWmunIpmgk7PiojAzz/bVn8LU1PRTNDpWRGJu9mz4fLL7f7AgWFzKVUaciIiIkFoyEmE1P2yMHcurLEGtGwJ8+dDo3qOwdUtXbHLPXZU845r7HzlPWaMLQ3WqhV88AF07Vqc2FGj07MJ7drB2mvDokXw5ZehsxERKZ4vvoATTrD7V16ZXsGMKxXNWjbd1LYffBA2DxGRYtpyS/juO/j97+HEE0NnU9pUNGtJFs333w+bh4hIof30Exx+uN2fPx8OOADGjKn/0pQY7Z5attnGtuPHB01DRKSgXnnFji5Hj7bHI0fCQw/ZZSpZORXNWgYMsO3EifDLL2FzERHJt+XL4cILoX9/67tRkehveuyxEME+OUGoaNbSoQP07m0zYUyaFDobEZH8mT7diuUll4D3cPbZdsQpmdE4TRERCULjNCOkvi8L48fDjjvCZpul7hAU97Fcil3+saOad1xj19f2iy/g6qvhjjtqJl8/+GD4xz9WvHYZep9FjU7P1rHddjWDe7/5JnQ2IiKZmTYNjjoKNtoIbrvN5pA98EB4+224/3519smVimYdzZrVdAh6+OGwuYiIpGvKFDjoIPjd7+Duu+25I4+EDz+0nrG9eoXNr1yoaKYwZIht77orbB4iIunYc0/rxPjww7b4xPHHw6ef2t+wTTYJnV15UdFMYa+9oG1b++b27ruhsxERWZH38Pzz1v8C4Omnbd7s006za5n//Cesv37YHMuVimYKq6wChx1m93W0KSKlwnt44gmbiGWXXWomYjn3XBt3ed110KlT0BTLnopmPY46yrb33KPFWEUkrKoqeOAB6NkT9t4b3nwT2revWfvyssvssRRebIecpNvVec4c6xyUTdtcY5dSW8WOV+yo5h2n2N9/b0eYucbNtX3I/R1CbItmOmOLrr0Whg+H/faDRx+158ptLJdiK3Y+2yp2fmP//DPsthu89hp07GhT4B19tF1CykfcXNtrnKasYMgQm/H/iSfsiFNEpFh+/hl23dUK5rrr2pR3xx23YsGU4lPRXImOHe1b3vLlcN99obMRkbj46ScrmK+/DuutZx1+1Bu2NKhoNuDoo22rXrQiUgwqmKVNRbMBe+0FbdrYFFTvvRc6GxEpZ8mC+cYbNQWza9fQWUltKpoNaN4cDj3U7utoU0QKpXbB7NpVBbNUqWimITnRwdixYfMQkfI1cKAKZhTEdshJNl2dp07Nvm2usUO3Vex4xY5q3lGO/eab1nbGjMwLZlT3WRSHnKRdNJ1zjYFK4Bvv/SDn3ANAj8TLbYB53vueKdq1AUYCmwEe+KP3/jXn3EXAn4HkYI4R3vuns/5NMpTp2KIBA+CFFwCiOSYq6mO5FDsabRU78/ZffAEbbOBwzjNlCmy5ZXHi5qN9HMdpZnKkeSowFVgNwHt/cPIF59x1wE/1tLsReMZ7f6BzrhnQstZr13vv/5ZZymH0758smiIi+XPeebYdMiTzginFl9Y1TedcZ2BP7Iix7msOGAz8ZiSjc241oD/wHwDv/VLv/bxcEg7lD38InYGIlJspU2D0aLt/6aVhc5H0pNsR6AZgOFCd4rXtgVne+09TvLYBdvr1DufcFOfcSOdcq1qvn+Sce885d7tzrm2qwM65oc65Sudc5ZyA0/JsvXXNHLQ//hgsDREpI2edVXN/vfXC5SHpa7BoOucGAbO995Pr+ZFDSXGUmdAE6A3803vfC1gInJ147Z9AN6An8B1wXao38N7/y3tf4b2v6NChQ0PpFkyLFrbCAMAHHwRLQ0TKxLhxdlt99dCZSCbSOdLsB+ztnJsB3A/s5Jy7B8A51wTYH3ignrYzgZne+zcSjx/Giije+1ne+yrvfTXwb2DrrH+LIunWzbaffx42DxGJturqmqPMs89e+c9KifHep30DdgCerPV4N2BCA20mAT0S9y8Crk3c71jrZ4YB9zcUv0+fPj4fsF68uummm266Bbw55/L1N73SZ1DLcrnlOk7zEOqcmnXOdQJGeu/3SDx1MnBvoufsdOCYxPPXOOd6JnbeDOC4HHPJiM+im/Sdd8IxxzgOPdT/evE+UxpGoNilHjuqeUcp9qabwkcfwciRcOyx0cm7ELGjJqOi6b0fD4yv9fjoFD/zLbBHrcfvABUpfu6ITGKXAp2eFZF8mDHDtoMHB01DsqBp9DKgoiki+bBkiW2bNw+bh2RORTMDHTva9ocfbHJlEZFMVVXZzTloEtuJTKNLRTMDtU+/62hTRLKxbJltmzVb8W+KRIOKZpZUNEUkG8lTs8nJUiRaVDSzpKIpItlYutS2q6wSNg/Jjsulu3CxVVRU+MrKypzfJ4rdnEVEyo1zjurqVLOzZvw+k733vxmlUQixvQydy7gk8IwYAZdfnl17jb1T7FKOHdW8oxJ7+nTrib/++nY/19hR32dRo9OzWdL1CBHJRvL0rP6GRJOKZpaaNg2dgYhEUbIjkK5pRpOKZpb0LVFEsqEjzWhT0cySjjRFJBsachJtKppZ0gdeRLKhISfRpqKZJR1pikg2dKQZbRqnKSIiQWicZoRk82Vh2TJo1szGaX73Hay9duZxNfZOsUs9dlTzjkrsyy6D88+H006D667LPXbU91nU6PRsBqZMsW337tkVTBGRt9+2be/eYfOQ7KhoZmDiRNv27x82DxGJrsmTbauiGU0qmhlQ0RSRXHz/PXz1FbRsaWesJHpUNNNUXQ2TJtl9FU0RyUbyEk/PntC4cdhcJDsqmmn64AOYN8/ur7de2FxEJJqS1zP79Ambh2RPQ05ERCQIDTmJkEy/LAweDA89BBDN7t1R75au2NFoq9grb7/hhraA/bvvwhZb5Cd21PdZ1Oj0bBq8hwkTQmchIlE2b54VzFVWgU02CZ2NZEtFMw3ffguzZ0O7dqEzEZGoeucd226xhabhjDIVzTTMmmXbLl3C5iEi0aVJDcqDimYakkVzrbXC5iEi0fXmm7bt1StsHpIbFc00qGiKSC6WLoWxY+3+zjuHzUVyo6KZhsWLbTt7dtg8RCSaXnwRfv7ZrmduuGHobCQXsR1ykk1X52efzb5trrFDt1XseMWOat6lHvu996C+HynlvEuxbSixLZqZji266CK4+GIAR3W1r/eD3xCNvVPsUo8d1bxLNXZVFXTsCHPm/HZ8Zj5iR32fRY1Oz6bp9NOhQwe7P2ZM2FxEJDpeecUKZrdusPnmobORXKlopql1a7jgArt/zjmwfHnYfEQkGh591Lb771//qVmJDhXNDAwdattp0+DOO4OmIiIR4D089pjd33//sLlIfqhoZqBZs5r7F14IixaFy0VESt/bb9v6mZ06wdZbh85G8kFFMwu9etnUesnTtSIiqSRPze63HzTSX9uyoH/GLNx6q/0H+Pvf4dVXQ2cjIqWqdtGU8hDbISe5dHXu27embb9+xY0dx7Fcih2ttor9WwMGFDZ2VPdZFIecxLZo5jouackSm3j5o4/gtNPguusya59L7GK3Vex4xY5q3qUW+/LL4bzz4Oij4Y47Chc76vssanR6NkurrGI9aBs3huuv12laEVlR7aEmUj5UNHOw1VYwfLh1Kz/mGPjll9AZiUgpmDHDes62agW77BI6G8knFc0cXXghbLopfPKJnYoREUmOzdxzT2jePGwukl8qmjmqe5r2lVdCZyQioSWLpnrNlh8VzTyoqICzzqo5TatJD0Tia+FCeO01G5a2++6hs5F8U9HMkwsugM02g08/hSuuCJ2NiITy+us2N3WvXrD66qGzkXxzuXQXLraKigpfWVmZ8/tEsZuziEi5cc5RXV2dj/eZ7L2vyENKDdI4zQw1NC5pzz3h6aetV+3VV2fePpfYhWqr2PGKHdW8SyX2jjvC+PF2XXPffQsfO+r7LGp0ejbPbKFquOUWW0NPROJjyRI7PQvw+9+HzUUKQ0UzzyoqbMqsRYtg3LjQ2YhIMb31FixebMPQ2rcPnY0UgopmAey4o22T3zhFJB4mTrRt//5h85DCUdEsgL59bauiKRIvKprlT0WzANZe27Y//RQ2DxEpruTkJttvHzYPKZzY9p4txnI2n3wCqX40qsvwKHZ8Ykc179CxFyyw9p07Fzd2VPdZFHvPxrZoFrKL9RNPwN57wx57wFNPZd4+l9iFaKvY8Yod1bxLITZ4jjkGbr+9eLGjv8+iRadnC2DuXNu2axc2DxEpPl3PLG8qmgXw44+2bds2bB4iUhy1J7VR0SxvKpoFoCNNkXj56CPbdu4M668fNhcpLBXNAtCRpki81B5qEsHLdJIBFc0C0JGmSLxofGZ8qGgWgI40ReLDexXNOIntkJNijC3aa69wsfPdVrHjFTuqeYeO/bvf6d+rWG1DiW3RLOS4pB49bGKDqVNh440zb59L7EK0Vex4xY5q3qFi3347HHssgP69smkfNWmfnnXONXbOTXHOPZl4/IBz7p3EbYZz7p162rVxzj3snJvmnJvqnNs28Xw759w459yniW3ZnMxMXtPU6VmR8pc8NSvxkMk1zVOBqckH3vuDvfc9vfc9gUeAR+tpdyPwjPd+Y2DLWu9xNvCC934j4IXE48jzXtc0ReJERTNe0iqazrnOwJ7AyBSvOWAwcF+K11YD+gP/AfDeL/Xez0u8vA9wV+L+XUCaa5yXtvnzoaoKWrWCZs1CZyMihfT11/DFF7DaaqEzkWJJ90jzBmA4UJ3ite2BWd77T1O8tgEwB7gjcWp3pHOuVeK1tbz33wEktmumCuycG+qcq3TOVc6ZMyfNdMNJHmVquIlI+Usu/7fddmHzkOJpsGg65wYBs733k+v5kUNJcZSZ0AToDfzTe98LWEiGp2G99//y3ld47ys6dOiQSdMgdD1TJD7+9z/bdusWNg8pnnSONPsBezvnZgD3Azs55+4BcM41AfYHHqin7Uxgpvf+jcTjh7EiCjDLOdcx8T4dgdlZ/QYl5vvvbauiKVL+kv/f11gjbB5SPA0OOfHenwOcA+Cc2wE4w3s/JPHyAGCa935mPW3/55z72jnXw3v/MbAzkJilkceBo4CrEtsxufwimSr02KIJE+qfTiuqY6IUOz6xo5p3qNiXXBIudq5tQ8aO4pCTXMdpHkKdU7POuU7ASO/9HomnTgbudc41A6YDxySevwp40Dl3LPAVcFCOuWSkEOOSqqpg3XXh22/h1Vdh220za59L7EK2Vex4xY5q3iFiH3wwPPggjB4Nhx2mf69s2kdNRkXTez8eGF/r8dEpfuZbYI9aj98BKlL83A/YkWfZmDTJCmbXrrDNNqGzEZFCS56ebd8+bB5SPJp7No/uSxxzH3qoVjoQiQMVzfhR0cyTpUvhoYfs/qGHhs1FRIpDRTN+VDTz5LnnbIzmZpvB5puHzkZECs179Z6NIxXNPBk92rY6yhSJhwUL7AxTy5Z2k3hwufR8KraKigpfWVmZ8/tEsceWiEi5cc5RXZ1qormM32ey9/43HU4LQUuDZShVF+v77oPDDrMes6+9lnn7XGIXo61ixyt2VPMuduy33nu7fYQAACAASURBVIKtt4bevWHyZP17Zds+anR6Ng+SvWYPOyxsHiJSPOoEFE8qmjmaOxeeeQYaNYLBg0NnIyLFoqIZTyqaOXrkEVi2DHbeGdZaK3Q2IlIsKprxpKKZI/WaFYknFc14UtHMwaxZNjF7s2aw336hsxGRYlLRjCcVzRw89ZQNcN55Z2jTJnQ2IlJMKprxFNshJ/lczmbs2Mzmmo3qMjyKHZ/YUc07ROzaHQD171W8tqHEtmjmOi5p0SL7hvnLLzBzJqyzTmbtc4ld7LaKHa/YUc272LH/9jc480xbHuz++/XvlW37qNHp2Sw9/7wVzIqK9AumiJSPAw+07ZNPwqJFYXOR4lHRzNKYMbbdZ5+weYhIGF272oxACxfaJRqJBxXNLFRVwRNP2P299w6bi4iEc9BBtk0uCyjlT0UzC2+8AXPm2DdNLQMmEl+1T9FKPKhoZuHxx227zz6Z9ZoVkfJS+xStxIOKZhaS1zN1alZEkqdoJR5iO+Qkl67O06ZZ2513Ln7sOI7lUuxotVXsaLUNGTuKQ05iWzRzGZcEnsMOg3vvza69xnIpdinHjmreIWP37Qtvvul4+GHPAQcUN3ZU91myfdTo9GyWNNRERJKSp2hHjoTly8PmIoWlopmB5FyTTZvCbruFzUVESsfBB9v2mWdg4ECYPTtsPlI4KpoZeOop2+64I6y2WthcRKR0dOli27XWgpdegj59bGialB8VzQyo16yIrMzbb8N229l81NtvD7feaishSflQ0UzTL7/As8/afRVNEUmlUyc70jzlFFi2DE44AY45xv5+SHlQ0UzTiy/WTMqcPBUjIlJXs2Zw443Wu75lS7jrLjv6nD49dGaSD7EdcqIxUYqt2IVpq9i/9c470K1b4WJHdZ9FcchJbItmJmOLqqtt+a///Q8gmmOioj6WS7Gj0Vaxf9t+3jw48khb5ME5uPBCOPdcaNKk4ba5xC1Ge43TlJTeessK5rrrhs5ERKKmTRv473/h8svt8UUXwSabwB132HVPiRYVzTQ8/7xt99wzbB4iEk2NGsGIEdaZsFs3+Owz+OMfoXt3uO220NlJJlQ007D66rZdsiRsHiISbbvsAtOmwahRsPHGMGMGHH+8vXbzzeplGwUqmmno08e2yckNRESy1aQJDBkCH3wADzxQsybvKafA+uvDdddpqbFSpqKZhm22gS23hFmzQmciIuWicWMYPNh61gL07m1/Y844w9bpvPJK+PnnoClKCi6Xnk/FVlFR4SsrK3N+nyj22BIRKTfOOaqrq/PxPpO99xV5SKlBGnKSpuXLYcMN4csvHaNGeYYMyS6uhhEodqnHjmre5RLbe3jhBbj0Upg40V5v3dqOQM88E1q0KM28s20fNTo9m6YmTWxsFdjUWFOnhs1HRMqTczBgAEyYYLcBA2D+fBvfuemm8Pjjms82JBXNDPzpT7ZdsAD220/XG0SksPr3h3HjYPx46zD0xRe2lu+gQTZsRYpPRTMDyTMJm28OH39ss3zk4XS8iMhK/eEPtoLKjTfasoRPP21HneedFzqz+FHRzMJjj9ksH2PGwBVXhM5GROKgSRMblvLJJ3D00bB0ac0sQ488olO2xaKimYVu3WD0aDvyvOACePnl0BmJSFystZZNwffqq9Crlz134IEwcKD6WhSDimaWdt8dzjnHvt396U+weHHojEQkTrbd1ubFBmjb1qb73GILGD7c+l1IYWicpoiIBKFxmhGSr3FJr7wC229vs3tMnmzf9DJpn0vsYrVV7HjFjmrecY1du21lJey1l63KdMMNcOqpxYudbfuo0enZHPXrB3/5i01+cOyxthURCaGiombVlHPOgU8/DZtPOVLRzIMrr4QuXexb3o03hs5GROJs771tQvhffrHlxzQsLr9UNPOgdWu49Va7f/758PnnYfMRkXi78UbrZfvyy1qvM99UNPNkjz3gsMPs291xx2nMlIiE065dzRjOsWPD5lJuVDTz6IYbYI01bLLlBx8MnY2IxFnnzrbVcLj8UtHMow4d7PomwOmna6yUiISzyiq2XbIkbB7lJrZDTnLp6pxO22++sWudIWIXoq1ixyt2VPOOa+yVtZ04sWbe7GLHLmTbUGJbNAs5LumNN2CbbaBpU/jgA+jePbP2ucQuRFvFjlfsqOYd19j1tZ082Yag9OwJU6YUN3Ym7aNGp2cLoG9f6+q9bBmcdlrobEQkjnR6tjBUNAvkyiuhZUt46ikNMBaR4lPRLAwVzQJZc00YPNju33572FxEJH5UNAtDRbOAjjzStk8/HTYPEYmfzz6zbZPY9lwpDBXNAlpnHdsuXBg2DxGJF+9trV+AY44Jm0u50dJgIiIShJYGi5BidLGeNQvWXtsmPZg9O/P2ucTOZ1vFjlfsqOYd19h121ZXQ58+8M47cN11DffgD73PokanZwtIF+JFpNgeftgK5jrrwAknhM6m/KhoFpCKpogU088/w/Dhdv+CC6BFi7D5lKO0T8865xoDlcA33vtBzrkHgB6Jl9sA87z3PVO0mwHMB6qA5cnzzs65i4A/A3MSPzrCe19W/UybNrXtsmVh8xCReDj9dPjySzs9qw5AhZHJNc1TganAagDe+4OTLzjnrgN+WknbHb3336d4/nrv/d8yyCFS5s61bdu2YfMQkfL31FMwcqSd4brrrpov7ZJfaZ2edc51BvYERqZ4zQGDgfvym1r0/e9/tl177bB5iEh5++EH+NOf7P5ll8Gmm4bNp5yle03zBmA4kKpv8PbALO99fZPFeeA559xk59zQOq+d5Jx7zzl3u3Mu5fGYc26oc67SOVc5Z86cVD9SspJFc621wuYhIuXt5JPt702/fjBsWOhsypz3fqU3YBDwj8T9HYAn67z+T+D0lbTvlNiuCbwL9E88XgtojBXuy4HbG8qlT58+Ph+wQq6bbrrpplvAm3MuX3/TK30D9SNft3SuafYD9nbO7QE0B1Zzzt3jvR/inGsC7A/0qa+x9/7bxHa2c+4xYGtgovd+VvJnnHP/Bp5MI5e88UUYl3TttdaTbdgw+PvfM2+fS+x8tlXseMWOat5xjG1jwR3g+b//g7/8pXixc22bbB81DZ6e9d6f473v7L3vChwCvOi9H5J4eQAwzXs/M1Vb51wr51zr5H1gIPBB4nHHWj+6X/L5cjIr8bVA1zRFJN+8h6GJC14DBsDxx4fNJy5ynRHoEOp0AHLOdQJGeu/3wE7BPpb4NtEEGO29fybxo9c453pih+kzgONyzKXk6JqmiBTKqFHw+ON2/z//gUYadV8UGRVN7/14YHytx0en+JlvgT0S96cDW9bzXkdkEjuKdKQpIoXw3HNw0kk1j9ddN1wucaPvJgWkI00RybdbboE99oD58+Hww0NnEz8qmgXy/fcwbRo0bgxdu4bORkSibvlyOPFEG15SVQUjRsDdd4fOKn5iu8pJoT38sH3Id9sN2rQJnY2IRNmPP8LgwfD889CsmV3DHDKk4XaSf1pPU0REgtB6mhFSyHFJn30GG21kKwzMmgWtW2fWPpfYhWir2PGKHdW8yzH2Sy/BAQfYkebmm8MTT8B66+Uvbq7tNU5Tcvbhh7DTTnZ///1XLJgiIun6979h4EArmHvtBa+88tuCKcWnoplHjz1mcz9+/bVtb7wxdEYiEjVVVTaL2NCh1i/izDPtb4u+gJcGFc08mDvXun7vvz/89JOdThk3DtZYI3RmIhI1F1wAN9xgS3vdfjtcc431wpfSoKKZoyeftGV4Ro+Gli3h5pvhwQe1YrqIZGfSJNvecYcWki5Fse0IlKt58+Cvf7XFXsFOx955J2y4YdC0RCTiFi60bffuYfOQ1HSkmYVnn4XNNrOC2bw5XHcdTJiggikiuVu0yLatWoXNQ1LTOE0REQlC4zQjJNMvCy+8AH/8I3z1laNZM88ll8AZZ2R+gV5j7xS71GNHNe9yib3GGta5cM4caN++sHFzbR/HcZqxLZrp+uYbOOssuPfemufefts6/4iI5FvymqZOz5YmXdOsx5IlcNVV0KOHFcxVVoHLLrPXVDBFpBC8r1nmS+tjlib9s6QwcaJ18DnnHPvWt99+MHUqnHtu6MxEpJw5ZxOyA7z1VthcJDUVzRQOOWTFx6usYpMVfPFFmHxEJD522MG248eHzELqo6KZwj/+YRfjV1/dHt9/Pxx3HGywgT0+7TT7QC9fHixFESlTKpqlTUNOREQkCA05iZBsviwsXw5NmzqGD/f897/wySc1r3XqZIvC/vGP1nmoPhpGoNilHrtu29VWg/nzM3oHoP7YrVrZBCH9+jUcO1Pl8u+12Wa2YtKkSfD73xcubq7t4zjkJLZHmrl+SLyHDz6weWZHj4bp02t+5g9/gL/8xToQNW2a/9jFbqvY8Ypdt22LFrB4MXz8sV3fr662W1VVzdb7mse9ezveeMP/+rj27aST4KOP7H132w0OPtgWOlhttdzzzvfvHTL2SSfB//2fXQq67rrCxc21fcjYdd6naEeaiQIQjVufPn18Ptivnb+21dXev/yy98ce632rVt7bnxDv11nH+6uu8v7HHwsXuxhtFTteseu2bdzYPs+XX+79kiW5xV682PsRI7xv0aLm/0nr1t7ffLP9P4rqPst37EmTvHfObuPGFS5uru1Dxq7zPpW+SHVIHYHywDk71TRyJHz7LdxyC2y8sU2McPbZtnDsiBHw/fehMxXJXHIdx3PPtZmxcrHKKnD55fDVV9bhbvvt7dTvySc3fEQVJ7//PZx/vn2tOPRQmDkzdEaSpNOzBWpbXW3Xba65pqYX3KqrwoIFjnnz/K89cwsRO99tFTteseu2ff112HZbu7/VVnDrrdCrl31ZzEfse+6BI46w8YlLl0ZznxUidlUV7LEHPPccbLONLQqRHMOZr7i5to/j6VkdaRZIo0aw++7w0kvw6qt2/WbBAnutWzdbd3PZsrA5iqRjm21sog+wAfd9+kCXLnD88fDIIzZHai6GDLFrm0uX5p5rOWnc2GYj69LFvriccUbojARUNIti221h7Fh45RV7/MMPcMopsOWW8OKLYXMTScell8JTT9l45U6d7NLDbbfBgQfCmmtab89hw+C777J7/wMOyG++5aJ9e3joIetQePPNNmZcwort6VkREQkriuM0Y3ukmW3PqVza1m6/eLHnsss8zZt7wNO+vefBBwsXO195K3b5x8607eLFngkTPJttZp9lgIULM497zDHWdvjw6O2zQseurvYcfrjt34039vz8czTyTvc9oiS2R5rZ/t75vvD9+ecwdGjNadohQ2x8VnLcWr5ix7WzgGIXr+0XX0BFBcyd6+jb1/Pcc6k/x3U984xd/09EZ4MNPH/8IzRpYqclmzatuV/3ua5doXfv3HPPtW0xYi9cCH372qQHO+1kvY979Cj9vAvZvtb7FO1IU0WziG3ra+89/POfdqH/l19sjtuHHlrxj0GuseP6H0uxi9v2o49g001tRqB77oHDD2+4zRZbwPvv/xqdlc0mlMp778Hmmydal/m/18cfW8esefOso9Cxx8K//lX6eReqfa330TR6ceKczSC08862wso778B221khPeaY0NmJpO93v6u5n+wt3pBRo+CSS6B7d1vD9oILrGf5smU2deXy5TWPaz83frz13P3uu5qiWe569LCZyC6+GP7zH/jXv+z5ESNg+HBo0yZsfnGgI80itk2n/eLFcOqpNf8Z/vpX+Nvf7FtlKeet2OUTO1XbefPg6adrhoXMn29jBlu2XHG8pnMwZIgdLV54IVx0UeHyHjjQlux79lm7n2n7XGLnu302badNswknHn3U9nfbtlY8TzrJ1gMuZOx8tM1H+1rvoyPNuGre3Lryb7WVHX3ecINdKxo9OnRmEmfnnmvX0DKRXEpPCmPjjW2crHM23/WECXDmmXDjjXbkfuSR9mVb8ktFs0T96U+w0UY26fuYMTXfpEXS9dNP9qXr++9X/OOZ/GLfrJlNUNCtW8PvlZwCsn9/m9lq6VJYe+2a90q+r/dw3312tLPnnvn7XVJJHuHmYcRC5L30knWoOuccePddW23pb3+DK66Avfeuf/YmyVxsT8+KiEhYGqcZIaHGJWXT/ssvPRtuaG179/bMmxeNvBU7bOyff06Om7Tb9dfX3JJj/g4/PL24Bx9sPz96dOnss113tZzGji2Pf698tl2yxHPTTZ4OHWr+/Xff3fPww55ffimNvGu/R5To9GwErLuunX7p0gXefhsGDbJJnFu0CJ2ZlLLWra0ndp8+dtr0z3+2BaDBJkm/9167JYc2OWe9UgFuuskeJ0/KPPCAbUvxb1wp5hRas2a2cszRR9vqMdddZ1N5jh1rY2cPOMCGA+2wQ+hMo0dFMyI6d7btOuvAyy/bJAgPPWQTw4vUZ8stbTjGO+/Aa6/BgAH2fLJ4Apx++m/bnXpq6vfLpFdmoekqS8Nat7YezCeeaEN7Ro+GyZPhjjvs1rGj/dzkyfblSfu0YSqaEfPss7Z256OPWo/GK68MnVE0eW89DD/80DrJ1L01alTTeebEE1O/x5572uo1f/sbfPll6p855RR7L+dW3KZ6ru4W7N+39vO173fqBIMHN/yHbq+9rGiecor9vs5Z3ldeCbNn1+yPqir47DM7Gjn55JoONskjufbtrZ1ET4cOcNppdps2zYrn6NE2IxnYTE49esBhh9ltww3D5lvScj0fXcxbnz59fD4Q0ZXKk23HjfO+cWPvwfuHHy583Fzbl2Lsr77yvqa/58pu1PvaWmt5/+ab2bXNNXby9sILDf/etX/XceOy32fpKtZnZbfd7Hd66qnix853+xBtq6u9f/11a7/mmit+rrbe2vsbbvD+u+8KEztf7Wu9T6UvUh3SkWYEDRhgRzfDhtmMQb16aUxcplq2rLl/7712lJW8VVfX3D/pJLjllhXb/vyzDSJfuBCWLLHnunb97WnOk0+2IR/e23vW3a7sOe/h8svh7LN/+3Pew/XXW4y33rJ5SFemSxdYYw1bku6vf7UhCeUwfs8njoB1iSI7ztlctmBLvb34ov1fePRRePNNu512ms1UdvjhsP/+dro39opVnfNx05FmTdvqau/339/+hG67rffLlhUubq7tSzF2VZX3ztn+W9m+S9V+wYLfHvFtvXX6sdO1svb//rfFPeyw3762dKm1/cc/vL/5Zu9/+MH7X37xft11rc2oUdnHzTXvfLYdMMB+n2efLX7sfLcvpbwXLfL+gQe832cf75s2rfmMt2vn/TXXeL9wYeFi5/A+RTvS1DhNEREJQuM0IyTbbxm5tM21faq2Y8fapa/mzT3Tp0cn71KI3b277buPPore7/3NN5Z748aeDTZIXP5c4QaDBnlatqx5rmvXmvv/+U/0/r3q3vr3t9/lpZdK/9+rVNum27662vPMM57evet+nmD58uxjJ+NHSWyLZrnYbTc49FCb6H3YsNDZRMsaa9j2hx/C5pGNNdeE9daz667Tp6/42scf2/aJJ2DSJFt7EmDGjJqfmT+/KGkWVFWVbcvh+mypcw523RUqK+Hhh2GTTWo+T5tvDk8+GTS9olLRLAPXXWfj7saMsUmbJT3JZZTmzAmbRzaaNLG1Kz/+GD791ArnV19ZZ6Hu3Wt+rndv+PZb+7lPPrEhJTNn1j8OM0qSZ/VUNIvHOZsY4f334a677LmpU21Y06mn1qyCU85UNMtAx462lh7Y2M0InvEIYpVVbPvuu2HzyFbLllYgN9wQ1l/fesmmulzfvr393EYb2eTs66xT/FwLIXmkqd6zxde4sa2iAtaTv2lTm0Wqf3/78lbO9HErE8OGQbt28MorOtpM18KFtq09O45Eh07PlobTT7fLAOuuC2+8YUPgxo4NnVXhaJxmmWjd2mZ8uegiuPZazSmZjj/8wRYxnjs3dCaSjWTRfPRRO1vQrFnN42bN7OinWbMVb/U9lzzrINnp29fmxT7ySFusfI89bCzzxRfbpYRyoiEnIiISRBSHnJTZd4D0ZftlwTmXUzfpXNqn0/aYY+DOO+0a59VX5yduru1LNfa//w1Dh8Kxx8LIkcWNXej2ccj7jTdswevFi2HZMruNGuXYbz/P0qX8elu2jBUe133txx9/jY73nuuvt5lwGvLjjzWdyTLNva5M2u67r3X6q9Wa5DCj7OTSvuG2U6fCxhvX0zqCBzCxLZrl6s9/tqI5apSt2q7rPfVLdljQ6dlo6tu3Zhq4pFGj7PRsJsaNg4ED7f5jj1khBTjhBDvFuHSpTZeYLLJbbWWvh/p7nxwq1bgxvPqq7YMJE1acCnL58prtsmV2P7lNvpa8P2wYXHVVzePabZP3wU5hN29u2+Tt5JPh7rt/+1rz5nZatksXWGutMPupUFQ0y8y221oPyc8/h4kTYccdQ2dU+h57LHQGEtKaa9of+OXLbX7VpNVXr1mSr7YWLeCXX2rGvxbb738Pt99uxWzVVe25/v2zf79hw+Css7Jre/LJcMQR2ceOIvWeLTPOwUEH2f1Mv3HHTfIPZLl9E5bMbLGFTQQBNUUI6l87NLlQd7HP4ixbZmdHBg+2hegh9VqoUlgqmmVon31s+9RTGrO5MhtvbKvYz5oFX38dOhsJJbm+KNiEEcOG2ZqSqY6gvK85fZvsrVsM3ttQjvXWs3Uvzz/fnp85s3g5iNHp2TK01VbQti188YXNFNOtW+iMSlOLFjYd2Btv2L7q0iV0RhJaly7w97/X/3pyxpumTYt7TXPpUltAHGwZr/fes/u1l7iT4tCRZhlq3NjGIIINOpb6ffaZbbt2DZqGRETyKLPY1zOTcZO+/NK2Wt+y+DROU0REgtA4zQgp9zFsL70EO+1k3dFffz06eRcz9pQp0KePnWZbvDj10UM5/t6FbFvusefMsd62a6wB33+f39iLF3sWLrTpHRcs4Nf7bdrYfMEdOtT8/Prr2yWFRGvGjPE0aWK9gBs3tm23bql7/6aKHfLfK2piWzTL3aab2nbqVHUGSuXVV22qL+/hqKPCDR+QaKl9TTMbI0bYULBkQUwWSKi/ty7A6NErPq4pmCbZ+a+2Zs3s+mf79tnlKqmpaJap9u3tCOrnn2sGJ4sVyTvvhJNOgkWLbJmj224LnZVERXINyWzmU/3hB7jyyvpfb9LEhry0alVzmz4dfvrJPrP/+Id1Ukpeh69t0KAVJyV4/XUbSzprlopmvqX9T++cawxUAt947wc55x4AeiRebgPM8973TNFuBjAfqAKWJ887O+faAQ8AXYEZwGDv/Y9120t25s+3AtG6dflNmJytxYvhxBNtYDjY9Hm33qr9I+nbd1/bJo8OM7F4sW3bt4dnn12xOLZr99vOPmCngNdayy63PPSQzVJUl3M140yTevSw9VM1I1j+ZdJ79lRgavKB9/5g733PRKF8BFjZUPodEz9b+0Lt2cAL3vuNgBcSjyVPliyxrVZvqHHKKVYwW7Swb+4jR6pgSmaSc6imKnANSZ7abdXKFgfv0cOuObZtW3+b9u2hXz+L9+yz6cdKnl3S5zv/0tqlzrnOwJ7A5cBpdV5zwGBgpwxj7wPskLh/FzAeyHIyJ5GGTZxo2zFjYJddwuYi0dSnD7z8cnZjn7MdrrLXXjZ07Iwz7MxR06YrdvgB69SWfNykiU2jCTrSLIR0v4fcAAwHUo0K2h6Y5b3/tJ62HnjOOeeB27z3/0o8v5b3/jsA7/13zrk1UzV2zg0FhgKsu+66aaYryUHPixaFzaNUfPWVna5aZRWbu1MkG0ceCTfemF0/gVyK5vDh9hnefffUP9O7d+rn1cEt/xosms65QcBs7/1k59wOKX7kUOC+lbxFP+/9t4miOM45N817PzHdBBNF9l9g4zTTbdeQXLo659pNupixFy2qmbkkSnkXKvaSJZnPolIOv3cx28Yh9ocfpp4RKJ32U6dm37Z+qdumO8tVVD9nIaRzTbMfsHeiQ8/9wE7OuXsAnHNNgP2xDj0pee+/TWxnA48BWydemuWc65h4n47A7Cx/h6x477O65dI21/aZtl1vPQ94PvkkWnnnO3Z1tWfjjW1fPPlkfH7vuOVdjNhvvWWfo969M2//5pvWtk+feO2zdN4jShosmt77c7z3nb33XYFDgBe990MSLw8ApnnvU04b7Jxr5ZxrnbwPDAQ+SLz8OHBU4v5RwJjfvoPkInnd5ZNPwuYR2j33wLRpsPbaNesmimQjucJJNh1sQk3BJ/mV69yzh1Dn1KxzrpNz7unEw7WAl51z7wJvAk95759JvHYVsItz7lNgl8RjyaPNNrNtcnLnuDrpJNteeaX+YEluktcyM+lg4721Sw450Wcw2jL6vuS9H4/1ck0+PjrFz3wL7JG4Px3Ysp73+gHYOZP4kpmKxACf118Pm0cIv/xiSzyBTfCw3342849ILpJHmq+9ZlPbVVVBdXVNMV1tNbtf+/m6U6uqaEabRvGUse23t+3EtLtdlYfFi20Q+nPP2eNrr4XTTivuUk5Snjbc0DqSLVoE337729fnz6+/bePGNrXd3nsXLj8pPBXNMta1q/0nTzXtVrlasgT2398K5pprwuzZNr5NJB/WWQf+9z/48Udo1MgKYePGdr9DB5g3b8Xnkvf1ha18xLZolnuX+Hy3jWLs2bPDxc5H25Cxo5p36Nht2ujfq1htQ4lt0cy2q3OpL11UV2UlbLUVgGPJEk+zZlmFDr58UEPtly2DwYPhv/+1eTxffBG23DJ6/16lEDuqecc1dlTzTraPmlx7z0qJ69OnphftoyubHTjCli2DQw+1gtmmDTz/vBVMEZF8U9Esc87Zyh5gywpFcCzxSi1fDkccAY88AquvDuPGQa9eobMSkXKlohkDRx5p27fequlRWg6qqmwYyQMP2ETWzz5bM8xGRKQQVDRjoPZcqyNG/HbcWBRVVcExx9iK9quuCs88A337hs5KRMqdA+9SmwAACL1JREFUimaMdOoEb78Nd90VOpPcVFfDn/8Mo0bZ2oRjx8J224XOSkTiQEUzRq6+2rbDh9uK8FFUXQ3HHw933GGLST/1lJb6EpHicVGaZb6iosJXVlbm/D5R7OYsIlJunHNU5+F6kXNusve+KD0aNE4zQ1EdE5Vs+/nnNhxj4UK4+27reVqs2Nmy/1ieU0+Fm2+2haQffzy9FUui+u8VMnZU845r7KjmnWwfNTo9GzPdutnK8wAnnGCL6UbB6adbwWzWDB57TEt8iUgYKpox9Mc/wuGH29HmvvvC3LmhM6pf8kvs9dfb6hCPPAK77x42JxGJLxXNGHIObrsNeva0ydz33bdmrb9S4j2cd57db9IEHnoIBg0Km5OIxJuKZky1agVPPGGrNkyaZNPQJVeWLxUXXwxXXGH3778f9tknbD4iIiqaMda5s41xbNPG5m0dMqR0Cuell1rRbJT4hB5wQNh8RERARTP2Nt/cpp9r3RoefNCO5hYsCJvTVVfBBRdYwRw1KmwuIiK1xXbIidauS23sWCugIWLXVV1tHZZCxM5X+6jGjmrecY0d17xDiG3RjNuYqHTafvop7LUXfPyxzed6001w9NHWcagYed9wAwwbZvf/8x/r5ZtJ+1xiF6J9VGNHNe+4xo5q3sn2UaPTs/KrjTaCN9+0xZwXLLCitdtuVkwL7ZZbagrmbbfVFEwRkVKioikrWG0166l6993Qtq0tJbbppvbarFmFiXnrrXDyyXb/H/+AoUMLE0dEJFcqmvIbztn0etOm2RHf8uX2/Nprw/rrw377wbXXwpQpuS9qPXKkzUwEdjo4eV9EpBTFdsL2uF0/yKXtu+9Cz56Oxo09VVUrvrbeelZghw6FLl0yi33nnVaUvYe//73m9Gw+c4/qv1fI2FHNO66xo5p3PtrXep+iTdiuI01p0JZb2vaXX2yu2rvvhmOPhY4d4csv4bLLYIMNrADOnJnee95zT03BvPrq+gumiEgpUdGUtDVtCr/7nR1ZjhxpBfKll+Dgg6343XEH9OhhvWBXttrPfffBUUdZm8svt/U9RUSiILanZ0VEJCytpxkhUTyHH4VrD2PG2GnXuXOhTx8YNw7atbO2jz5qw1mqqmwi9ksvLXzuUf33Chk7qnnHNXZU8062jxqdnpW82mcfG+vZrRtMnlwzyfrjj9tp3KoqOPtsuOSSsHmKiGRDRVPyrls3u9bZqZOtoAJw4IE2dOX0023lkgh+wRQRUdGUwujSxSaAT1q2DE491cZ3qmCKSFSpaErB9OsHxx1n9088Ea6/XgVTRKJNRVMK6qabbHvzzSqYIhJ9sew927Zt26x7bTnncl4KJ0Ts0Hk3ahS9fRbX2FHNO66xo5o32N/iqIll0Zw7d27oFEREJIJ0elZERCRNKpoiIiJpUtEUERFJk4qmiIhImlQ0RURE0qSiKSIikiYVTRERkTSpaIqIiKRJRVNERCRNKpoiIiJpUtEUERFJk4qmiIhImlQ0RURE0qSiKSIikiYVTRERkTSpaIqIiKRJRVNERCRNznsfOoe0OefmAF+GziOF9sD3oZMoQdov9dO+SU37JTXtl/q1B1p57zsUI1ikimapcs5Veu8rQudRarRf6qd9k5r2S2raL/Ur9r7R6VkREZE0qWiKiIikSUUzP/4VOoESpf1SP+2b1LRfUtN+qV9R942uaYqIiKRJR5oiIiJpUtEUERFJk4pmHc65g5xzHzrnqp1zFbWe38U5N9k5935iu1Pi+ZbOuaecc9MS7a5q4P3Xdc4tcM6dUeu58c65j51z7yRuaxbuN8xOoP3SJ/G+nznnbnLOucL9htkp1H5xzm1d6/PwrnNuv1qvxfbz0sB+KfnPCxR036Rsn3gtzp+Zle2XzD8z3nvdat2ATYAewHigotbzvYBOifubAd8k7rcEdkzcbwZMAnZfyfs/AjwEnFHruRVileIt0H55E9gWcMDYlbUvt/2S+Lkmifsdgdm1Hsf289LAfin5z0uB903K9vrMrHS/ZPyZaYKswHs/FaDuFw7v/ZRaDz8EmjvnVvHeLwJeSvzMUufc20DnVO/tnNsXmA4sLEDqBVXs/eKc6wis5r1/LfH4bmBf7INdMgq1XxI/l9QciFSPvWLvl6h8XqCg+6a+9kvy/CsURLH3C9COLD4zOj2bnQOAKXU/jM65NsBewAt1GzjnWgFnARfX8553JE6bnF+qp5XSkM/9sg4ws9bjmYnnoijj/ZJ4va9z7kPgfeB47/3yWi/H8vOSeD3VfimnzwtkuW8aaB/bz0w97bP6zMTySNM59zywdoqXzvXej2mg7abA1cDAOs83Ae4DbvLeT0/R9GLgeu/9ghSf18O9998451pjpymPAO5O65fJoxLbL6n+Uwc52gq0X/DevwFs6pzbBLjLOTfWe7+YeH9eUu4XSujzAuH2zUrax/ozU0/7rD4zsSya3vsB2bRzznUGHgOO9N5/XuflfwGfeu9vqKd5X+BA59w1QBug2jm32Ht/i/f+m0Re851zo4GtCfCBLqX9gv3Hrn2qpTPwbTb55SrQfqkdf6pzbiF2PaYy5p+X2vFr75eZlMjnBcLtm/rax/0zU0/7rD4zOj2bpsTh/1PAOd77V+q8dhmwOvDX+tp777f33nf13ncFbgCu8N7f4pxr4pxrn3ifpsAg4IMC/Rp5V6j94r3/DpjvnNsmcSrpSGCl30RLSa77xTm3fuIbNM659bAOEjPi/nmpb79E/fMCedk3KdvrM5O6fdafmYZ6CsXtBuyHfQNZAswCnk08fx7WUeWdWrc1sW8nHpha6/k/JdrsDVySIsZFJHqJAq2AycB72EXqG4HGofdD6P2SeFyB/ef+HLiFxAxWpXQr1H7BTp99mHj9bWBffV7q3y9R+bwUeN/U1z7un5mU7bP9zGgaPRERkTTp9KyIiEiaVDRFRETSpKIpIiKSJhVNERGRNKloioiIpElFU0REJE0qmiIiImn6f9qZ9K9Mf/QPAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 504x1440 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "fig, ax = plt.subplots(figsize=(7,20))\n",
    "single_city.plot(ax=ax, color='None',edgecolor='blue',linewidth=2 )\n",
    "seattle_geohash_gpd.plot(ax=ax, color='None',edgecolor='black')\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# save to json if you want to view on a web map or QGIS etc"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2019-10-24T13:58:21.398513Z",
     "start_time": "2019-10-24T13:58:20.999800Z"
    }
   },
   "outputs": [],
   "source": [
    "seattle_geohash_json = seattle_geohash_gpd.to_json()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2019-10-24T13:58:21.534973Z",
     "start_time": "2019-10-24T13:58:21.402676Z"
    }
   },
   "outputs": [],
   "source": [
    "with open(f\"seattle_geohash_boxes_level{geohash_precision}.geojson\", 'w') as f:\n",
    "    f.write(seattle_geohash_json)"
   ]
  },
  {
   "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.6.7"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}