{ "cells": [ { "cell_type": "markdown", "metadata": { "deletable": true, "editable": true }, "source": [ "# Point in Polygon & Intersect\n", "\n", "Finding out if a certain point is located inside or outside of an area,\n", "or finding out if a line intersects with another line or polygon are\n", "fundamental geospatial operations that are often used e.g. to select\n", "data based on location. Such spatial queries are one of the typical\n", "first steps of the workflow when doing spatial analysis. Performing a\n", "spatial join (will be introduced later) between two spatial datasets is\n", "one of the most typical applications where Point in Polygon (PIP) query\n", "is used.\n", "\n", "## How to check if point is inside a polygon?\n", "\n", "Computationally, detecting if a point is inside a polygon is most commonly done using a specific formula called [Ray Casting algorithm](https://en.wikipedia.org/wiki/Point_in_polygon#Ray_casting_algorithm).\n", "Luckily, we do not need to create such a function ourselves for\n", "conducting the Point in Polygon (PIP) query. Instead, we can take\n", "advantage of [Shapely's binary predicates](http://toblerity.org/shapely/manual.html#binary-predicates)\n", "that can evaluate the topolocical relationships between geographical\n", "objects, such as the PIP as we're interested here.\n", "\n", "There are basically two ways of conducting PIP in Shapely:\n", "\n", "1. using a function called\n", " [.within()](http://toblerity.org/shapely/manual.html#object.within)\n", " that checks if a point is within a polygon\n", "2. using a function called\n", " [.contains()](http://toblerity.org/shapely/manual.html#object.contains)\n", " that checks if a polygon contains a point\n", "\n", "Notice: even though we are talking here about **Point** in Polygon\n", "operation, it is also possible to check if a LineString or Polygon is\n", "inside another Polygon.\n", "\n", "- Let's first create a Polygon using a list of coordinate-tuples and a\n", " couple of Point objects" ] }, { "cell_type": "code", "execution_count": 48, "metadata": { "collapsed": false, "deletable": true, "editable": true }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "POINT (24.952242 60.1696017)\n", "POINT (24.976567 60.16125)\n", "POLYGON ((24.950899 60.169158, 24.953492 60.169158, 24.95351 60.170104, 24.950958 60.16999, 24.950899 60.169158))\n" ] } ], "source": [ "from shapely.geometry import Point, Polygon\n", "\n", "# Create Point objects\n", "p1 = Point(24.952242, 60.1696017)\n", "p2 = Point(24.976567, 60.1612500)\n", "\n", "# Create a Polygon\n", "coords = [(24.950899, 60.169158), (24.953492, 60.169158), (24.953510, 60.170104), (24.950958, 60.169990)]\n", "poly = Polygon(coords)\n", "\n", "# Let's check what we have\n", "print(p1)\n", "print(p2)\n", "print(poly)" ] }, { "cell_type": "markdown", "metadata": { "deletable": true, "editable": true }, "source": [ "Let's check if those points are ``within`` the polygon" ] }, { "cell_type": "code", "execution_count": 49, "metadata": { "collapsed": false, "deletable": true, "editable": true }, "outputs": [ { "data": { "text/plain": [ "False" ] }, "execution_count": 49, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# Check if p1 is within the polygon using the within function\n", "p1.within(poly)\n", "\n", "# Check if p2 is within the polygon\n", "p2.within(poly)" ] }, { "cell_type": "markdown", "metadata": { "deletable": true, "editable": true }, "source": [ "Okey, so we can see that the first point seems to be inside that polygon\n", "and the other one doesn't.\n", "\n", "- In fact, the first point is close to the center of the polygon as we\n", " can see if we compare the point location to the centroid of the polygon:" ] }, { "cell_type": "code", "execution_count": 50, "metadata": { "collapsed": false, "deletable": true, "editable": true }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "POINT (24.952242 60.1696017)\n", "POINT (24.95224242849236 60.16960179038188)\n" ] } ], "source": [ "# Our point\n", "print(p1)\n", "\n", "# The centroid\n", "print(poly.centroid)" ] }, { "cell_type": "markdown", "metadata": { "deletable": true, "editable": true }, "source": [ "- It is also possible to do PIP other way around, i.e. to check if\n", " polygon contains a point:" ] }, { "cell_type": "code", "execution_count": 51, "metadata": { "collapsed": false, "deletable": true, "editable": true }, "outputs": [ { "data": { "text/plain": [ "False" ] }, "execution_count": 51, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# Does polygon contain p1?\n", "poly.contains(p1)\n", "\n", "# Does polygon contain p2?\n", "poly.contains(p2)" ] }, { "cell_type": "markdown", "metadata": { "deletable": true, "editable": true }, "source": [ "Thus, both ways of checking the spatial relationship results in the same way.\n", "\n", "Which one should you use then? Well, it depends:\n", "\n", "- if you have many points and just one polygon and you try to find out\n", " which one of them is inside the polygon:\n", "\n", " - you need to iterate over the points and check one at a time if it\n", " is **within()** the polygon specified\n", "\n", "- if you have many polygons and just one point and you want to find out\n", " which polygon contains the point\n", "\n", " - you need to iterate over the polygons until you find a polygon that\n", " **contains()** the point specified (assuming there are no overlapping\n", " polygons)\n", " \n", "## Intersect\n", "\n", "\n", "Another typical geospatial operation is to see if a geometry\n", "[intersect](http://toblerity.org/shapely/manual.html#object.intersects)\n", "or [touches](http://toblerity.org/shapely/manual.html#object.touches)\n", "another one. The difference between these two is that:\n", "\n", "- if objects intersect, the boundary and interior of an object needs to\n", " intersect in any way with those of the other.\n", "\n", "- If an object touches the other one, it is only necessary to have (at\n", " least) a single point of their boundaries in common but their\n", " interiors shoud NOT intersect.\n", "\n", "Let's try these out.\n", "\n", "- Let's create two LineStrings" ] }, { "cell_type": "code", "execution_count": 52, "metadata": { "collapsed": true, "deletable": true, "editable": true }, "outputs": [], "source": [ "from shapely.geometry import LineString, MultiLineString\n", "\n", "# Create two lines\n", "line_a = LineString([(0, 0), (1, 1)])\n", "line_b = LineString([(1, 1), (0, 2)])" ] }, { "cell_type": "markdown", "metadata": { "deletable": true, "editable": true }, "source": [ "- Let's see if they intersect" ] }, { "cell_type": "code", "execution_count": 53, "metadata": { "collapsed": false, "deletable": true, "editable": true }, "outputs": [ { "data": { "text/plain": [ "True" ] }, "execution_count": 53, "metadata": {}, "output_type": "execute_result" } ], "source": [ "line_a.intersects(line_b)" ] }, { "cell_type": "markdown", "metadata": { "deletable": true, "editable": true }, "source": [ "- Do they also touch each other?" ] }, { "cell_type": "code", "execution_count": 54, "metadata": { "collapsed": false, "deletable": true, "editable": true }, "outputs": [ { "data": { "text/plain": [ "True" ] }, "execution_count": 54, "metadata": {}, "output_type": "execute_result" } ], "source": [ "line_a.touches(line_b)" ] }, { "cell_type": "markdown", "metadata": { "deletable": true, "editable": true }, "source": [ "Indeed, they do and we can see this by plotting the features together" ] }, { "cell_type": "code", "execution_count": 55, "metadata": { "collapsed": false, "deletable": true, "editable": true }, "outputs": [ { "data": { "image/svg+xml": [ "" ], "text/plain": [ "" ] }, "execution_count": 55, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# Create a MultiLineString from line_a and line_b\n", "multi_line = MultiLineString([line_a, line_b])\n", "multi_line" ] }, { "cell_type": "markdown", "metadata": { "deletable": true, "editable": true }, "source": [ "Thus, the ``line_b`` continues from the same node ( (1,1) ) where ``line_a`` ends.\n", "\n", "However, if the lines overlap fully, they don't touch due to the spatial relationship rule, as we can see:\n", "\n", "- Check if line_a touches itself" ] }, { "cell_type": "code", "execution_count": 56, "metadata": { "collapsed": false, "deletable": true, "editable": true }, "outputs": [ { "data": { "text/plain": [ "False" ] }, "execution_count": 56, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# Does the line touch with itself?\n", "line_a.touches(line_a)" ] }, { "cell_type": "markdown", "metadata": { "deletable": true, "editable": true }, "source": [ "- It does not. However, it does intersect" ] }, { "cell_type": "code", "execution_count": 57, "metadata": { "collapsed": false, "deletable": true, "editable": true }, "outputs": [ { "data": { "text/plain": [ "True" ] }, "execution_count": 57, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# Does the line intersect with itself?\n", "line_a.intersects(line_a)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Point in Polygon using Geopandas\n", "\n", "Next we will do a practical example where we check which of the addresses from [previous tutorial]() are located in Southern district of Helsinki. Let's start by reading a KML-file ``PKS_suuralue.kml`` that has the Polygons for districts of Helsinki Region (data openly available from [Helsinki Region Infoshare](http://www.hri.fi/fi/dataset/paakaupunkiseudun-aluejakokartat).\n", "\n", "- Let's start by reading the addresses from the Shapefile that we saved earlier." ] }, { "cell_type": "code", "execution_count": 58, "metadata": { "collapsed": true, "deletable": true, "editable": true }, "outputs": [ { "data": { "text/html": [ "
\n", "\n", "\n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", "
addressidaddrgeometry
0Ruoholahti, 14, Itämerenkatu, Ruoholahti, Läns...1000Itämerenkatu 14, 00101 Helsinki, FinlandPOINT (24.9155624 60.1632015)
1Kamppi, 1, Kampinkuja, Kamppi, Eteläinen suurp...1001Kampinkuja 1, 00100 Helsinki, FinlandPOINT (24.9316914 60.1690222)
2Bangkok9, 8, Kaivokatu, Keskusta, Kluuvi, Etel...1002Kaivokatu 8, 00101 Helsinki, FinlandPOINT (24.9416849 60.1699637)
31, Hermannin rantatie, Hermanninmäki, Hermanni...1003Hermannin rantatie 1, 00580 Helsinki, FinlandPOINT (24.9655355 60.2008878)
4Hesburger, 9, Tyynenmerenkatu, Jätkäsaari, Län...1005Tyynenmerenkatu 9, 00220 Helsinki, FinlandPOINT (24.9216003 60.1566475)
\n", "
" ], "text/plain": [ " address id \\\n", "0 Ruoholahti, 14, Itämerenkatu, Ruoholahti, Läns... 1000 \n", "1 Kamppi, 1, Kampinkuja, Kamppi, Eteläinen suurp... 1001 \n", "2 Bangkok9, 8, Kaivokatu, Keskusta, Kluuvi, Etel... 1002 \n", "3 1, Hermannin rantatie, Hermanninmäki, Hermanni... 1003 \n", "4 Hesburger, 9, Tyynenmerenkatu, Jätkäsaari, Län... 1005 \n", "\n", " addr \\\n", "0 Itämerenkatu 14, 00101 Helsinki, Finland \n", "1 Kampinkuja 1, 00100 Helsinki, Finland \n", "2 Kaivokatu 8, 00101 Helsinki, Finland \n", "3 Hermannin rantatie 1, 00580 Helsinki, Finland \n", "4 Tyynenmerenkatu 9, 00220 Helsinki, Finland \n", "\n", " geometry \n", "0 POINT (24.9155624 60.1632015) \n", "1 POINT (24.9316914 60.1690222) \n", "2 POINT (24.9416849 60.1699637) \n", "3 POINT (24.9655355 60.2008878) \n", "4 POINT (24.9216003 60.1566475) " ] }, "execution_count": 58, "metadata": {}, "output_type": "execute_result" } ], "source": [ "import geopandas as gpd\n", "\n", "fp = \"data/addresses.shp\"\n", "data = gpd.read_file(fp)\n", "\n", "data.head()" ] }, { "cell_type": "markdown", "metadata": { "deletable": true, "editable": true }, "source": [ "\n", "### Reading KML-files in Geopandas\n", "\n", "It is possible to read the data from KML-files with GeoPandas in a similar manner as Shapefiles. However, we need to first, enable the KML-driver which is not enabled by default (because KML-files can contain unsupported data structures, nested folders etc., hence be careful when reading KML-files). Supported drivers are managed with [`fiona.supported_drivers`](https://github.com/Toblerity/Fiona/blob/master/fiona/drvsupport.py), which is integrated in geopandas. Let's first check which formats are currently supported:" ] }, { "cell_type": "code", "execution_count": 59, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "{'AeronavFAA': 'r',\n", " 'ARCGEN': 'r',\n", " 'BNA': 'raw',\n", " 'DXF': 'raw',\n", " 'OpenFileGDB': 'r',\n", " 'ESRI Shapefile': 'raw',\n", " 'GeoJSON': 'rw',\n", " 'GPKG': 'rw',\n", " 'GPX': 'raw',\n", " 'GPSTrackMaker': 'raw',\n", " 'Idrisi': 'r',\n", " 'MapInfo File': 'raw',\n", " 'DGN': 'raw',\n", " 'PCIDSK': 'r',\n", " 'SEGY': 'r',\n", " 'SUA': 'r',\n", " 'KML': 'rw'}" ] }, "execution_count": 59, "metadata": {}, "output_type": "execute_result" } ], "source": [ "import geopandas as gpd\n", "gpd.io.file.fiona.drvsupport.supported_drivers" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "- Let's enable the read and write functionalities for KML-driver by passing ``'rw'`` to whitelist of fiona's supported drivers:" ] }, { "cell_type": "code", "execution_count": 60, "metadata": { "collapsed": true, "deletable": true, "editable": true }, "outputs": [], "source": [ "gpd.io.file.fiona.drvsupport.supported_drivers['KML'] = 'rw'" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Let's check again the supported drivers:" ] }, { "cell_type": "code", "execution_count": 61, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "{'AeronavFAA': 'r',\n", " 'ARCGEN': 'r',\n", " 'BNA': 'raw',\n", " 'DXF': 'raw',\n", " 'OpenFileGDB': 'r',\n", " 'ESRI Shapefile': 'raw',\n", " 'GeoJSON': 'rw',\n", " 'GPKG': 'rw',\n", " 'GPX': 'raw',\n", " 'GPSTrackMaker': 'raw',\n", " 'Idrisi': 'r',\n", " 'MapInfo File': 'raw',\n", " 'DGN': 'raw',\n", " 'PCIDSK': 'r',\n", " 'SEGY': 'r',\n", " 'SUA': 'r',\n", " 'KML': 'rw'}" ] }, "execution_count": 61, "metadata": {}, "output_type": "execute_result" } ], "source": [ "gpd.io.file.fiona.drvsupport.supported_drivers" ] }, { "cell_type": "markdown", "metadata": { "deletable": true, "editable": true }, "source": [ "Now we should be able to read a KML file with [geopandas.read_file](http://geopandas.org/reference/geopandas.read_file.html#geopandas.read_file).\n", "\n", "- Let's read the data from a following KML -file:" ] }, { "cell_type": "code", "execution_count": 62, "metadata": { "collapsed": true, "deletable": true, "editable": true }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Number of rows: 23\n" ] }, { "data": { "text/html": [ "
\n", "\n", "\n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", "
NameDescriptiongeometry
0Suur-EspoonlahtiPOLYGON Z ((24.775059677807 60.1090604462157 0...
1Suur-KauklahtiPOLYGON Z ((24.6157775254076 60.1725681273527 ...
2Vanha-EspooPOLYGON Z ((24.6757633262026 60.2120070032819 ...
3Pohjois-EspooPOLYGON Z ((24.767921197401 60.2691954732391 0...
4Suur-MatinkyläPOLYGON Z ((24.7536131356802 60.1663051341717 ...
5KauniainenPOLYGON Z ((24.6907528033566 60.2195779731868 ...
6Suur-LeppävaaraPOLYGON Z ((24.797472695835 60.2082651196077 0...
7Suur-TapiolaPOLYGON Z ((24.8443596422129 60.1659790707387 ...
8MyyrmäkiPOLYGON Z ((24.8245867448802 60.2902531157585 ...
9KivistöPOLYGON Z ((24.9430919106369 60.3384471629062 ...
10EteläinenPOLYGON Z ((24.7827651307035 60.09997268858 0,...
\n", "
" ], "text/plain": [ " Name Description \\\n", "0 Suur-Espoonlahti \n", "1 Suur-Kauklahti \n", "2 Vanha-Espoo \n", "3 Pohjois-Espoo \n", "4 Suur-Matinkylä \n", "5 Kauniainen \n", "6 Suur-Leppävaara \n", "7 Suur-Tapiola \n", "8 Myyrmäki \n", "9 Kivistö \n", "10 Eteläinen \n", "\n", " geometry \n", "0 POLYGON Z ((24.775059677807 60.1090604462157 0... \n", "1 POLYGON Z ((24.6157775254076 60.1725681273527 ... \n", "2 POLYGON Z ((24.6757633262026 60.2120070032819 ... \n", "3 POLYGON Z ((24.767921197401 60.2691954732391 0... \n", "4 POLYGON Z ((24.7536131356802 60.1663051341717 ... \n", "5 POLYGON Z ((24.6907528033566 60.2195779731868 ... \n", "6 POLYGON Z ((24.797472695835 60.2082651196077 0... \n", "7 POLYGON Z ((24.8443596422129 60.1659790707387 ... \n", "8 POLYGON Z ((24.8245867448802 60.2902531157585 ... \n", "9 POLYGON Z ((24.9430919106369 60.3384471629062 ... \n", "10 POLYGON Z ((24.7827651307035 60.09997268858 0,... " ] }, "execution_count": 62, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# Filepath to KML file\n", "fp = \"data/PKS_suuralue.kml\"\n", "\n", "polys = gpd.read_file(fp, driver='KML')\n", "\n", "#Check the data\n", "print(\"Number of rows:\",len(polys))\n", "polys.head(11)" ] }, { "cell_type": "markdown", "metadata": { "deletable": true, "editable": true }, "source": [ "Nice, now we can see that we have 23 districts in our area. We are interested in an area that is called ``Eteläinen`` (*'Southern'* in english).\n", "\n", "- Let's select the ``Eteläinen`` district and see where it is located on a map\n" ] }, { "cell_type": "code", "execution_count": 63, "metadata": {}, "outputs": [], "source": [ "# Select data and reset index\n", "southern = polys.loc[polys['Name']=='Eteläinen']\n", "southern.reset_index(drop=True, inplace=True)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Let's create a map which shows the location of the selected district, and let's also plot the address points on top of the map" ] }, { "cell_type": "code", "execution_count": 64, "metadata": { "collapsed": true, "deletable": true, "editable": true }, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAaIAAAEYCAYAAAAeWvJ8AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4wLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvqOYd8AAAIABJREFUeJzsvXlwnPWd5//69n2r1a1W65bv2xiw8RDbGIhjCMEJMDGEQAiQZBImm9TMb6t2JvXbrdqq3x41NTuTHaZmk6lcLAnMBMIMmIQkXDGJg218gG3A+MDGutWSWrda6vP7+0PqRrKultTdTx/fV5VKltz99Kdb3c/7+dxCSolCoVAoFFqh09oAhUKhUJQ2SogUCoVCoSlKiBQKhUKhKUqIFAqFQqEpSogUCoVCoSlKiBQKhUKhKUqIFAqFQqEpSogUCoVCoSlKiBQKhUKhKQatDZiJiooKuWzZMq3NUCgUCsUSOHnyZI+U0jff7fJSiJYtW8aJEye0NkOhUCgUS0AI0ZTO7VRoTqFQKBSaooRIoVAoFJqihEihUCgUmqKESKFQKBSaooRIoVAoFJqihEihUCgUmpKWEAkh3EKI54QQ54QQHwghPiGE8AghXhVCXJz4Xj7H/V1CiDYhxD9lznSFQqFQFAPpekSPA7+VUq4DtgAfAN8BXpdSrgZen/h5Nv4b8PulGKpQKBSK4mReIRJCuIDdwI8BpJQRKWU/cBfw5MTNngTunuX+WwE/8EomDFYoFApFcZHOZIUVQDfwhBBiC3AS+AvAL6XsAJBSdgghKq++oxBCB/w98BCwJ2NWKxR5TDwep6enh0AgQFdXF11dXfj9fvbsUR8BhWIm0hEiA3A98G0p5VtCiMeZOww3mW8Cv5ZStggh5ryhEOLrwNcBGhoa0jy8QqEdUkoGBgamCE4gECAYDJJIJFK38/v93HPPPRpaqlDkN+kIUSvQKqV8a+Ln5xgXooAQonrCG6oGuma47yeAm4QQ3wQcgEkIMSylnCZkUsofAD8A2LZtm1zEc1Eossbo6GhKaCZ/j0Qi8943EAjQ1NTEunXrcmCpQlF4zCtEUspOIUSLEGKtlPI84yG2sxNfDwN/M/H9wAz3fTD5byHEI8C2mURIochX2tvb+dWvfkVHR8eSjnPw4EFCoRBbtmxBr9dnyDqFojhId/r2t4GnhRAm4DLwKOOFDs8KIb4KNAP3AgghtgGPSSm/lgV7FYqcEI1GeeONNzhy5AhSLt1B7+rq4pe//CUnTpxg586dbNiwgfnC1QpFqSAy8SHLNNu2bZNqDYRCK5qamnjxxRfp7e3N2mPcc889XHPNNVk7vkKRDwghTkopt813u7zcR6RQaEE4HOa1117LyS6s/v7+rD+GQlEoKCFSKIAPP/yQX/7ylwwODmb9sUwmEytXrsz64ygUhYISIkVJMzo6yssvv8zp06dz9phCCK5cuUJVVZUqXFAoUEKkKGHOnj3Lr3/9a0ZGRnL6uMkQ4OnTp9m3b5/qm1OUPEqIFCXH8PAwv/71r/nggw80taO7u5snnniC6667jr1792K1WjW1R6HQCiVEipJBSsmZM2f47W9/y9jYmNbmpHjnnXc4f/48t99+O5s3b1Zl3YqSQwmRoiQYGBjgV7/6FR9++KHWpsxIKBTi+eef59SpU9x55514vV6tTcpLOjs7ee+992hoaGDVqlXodDPPbZZScvnyZd555x2EEGzdupXGxkYl8nmK6iNSFDVSSk6cOMFrr72W1jiefECv17Nz50527tyJyWTS2pysE41GOX/+PIFAgJ07d2I0GgkEAlRXV6eEQ0pJU1MTP/vZz1Jz/NxuN1u3buX666/HZrPR1NTE0NAQ4XCYd999l6amptRjeL1evvKVr2Cz2TR5jqVKun1ESogURUtvby8vvvjilBNSIeFwONizZw9btmwpuit5KSXNzc2cPn2as2fPEg6HgfHnHI/HGR0dZeXKlXzqU5+ira2N48ePEwgEZjyWXq+noqJi1v9PYrfb2b17NzfccEPRvZ75ihIiRcmSSCQ4evQoBw8eJBaLaW3OkqmqquL2229n2bJlWpuSMaSUHD9+nNdee41oNJqTx7RYLGzbtg2fz4fVakWv12MwGNDr9VO+Jv/OYDBgMKgMxmJRQqSYl3g8jpSyqD5o3d3dHDhwgLa2Nq1NyTgbNmxg3759RVVd19fXx4svvsiVK1e0NmVGHA4H3/rWtzCbzVqbUpCoET+KKUgpGRoaoqenJ/XV1NREd3c3FRUV7N+/n8rKabsNC4Z4PE4kEuHIkSNFKUIw3vfU3t7OvffeS01NjdbmZITy8nK+/OUvc/z4cV599dW882CHh4dpa2tjxYoVWptS1CiPqEBJJBLEYrFpyexYLEZvb+8UwUl+zRUCsVgs7N27lxUrVlBWVlZQMfTOzk4OHDhAdXU1er0+J7PitESv13P77bezbdu2gvo7zUcwGOSFF16gtbVVa1OmsHLlSsrKyojH48B4kUQikSAej7Np0yZcLhdSShKJBFLKaf8uLy8v2QkaKjRXxESjUf7t3/6N8+fPYzabcTqd2O12BgcH6e/vX/LagvLycqxWK+FwmJqaGuLxOMuWLcNut7Nq1aq8qeSKxWIcOnSIP/7xj6lKKr1enzphFDubN29m3759efP3yASJRILDhw9z8ODBKVtuCxGn08lnPvOZkl6IqISoSBkdHeXnP/85zc3Nmjy+w+Fg8+bN7N27V9Or8ba2Ng4cOEB3d7dmNuQDFRUV3Hffffh8Pq1NySiBQIDnn39+3kq4fEQIwfbt27n11ltLPrekhKgIGRoa4qmnnqKra6at7LmlpqaGG264Abvdjtfrpby8PCfClOmFdcWA0Wjks5/9LJs3b9balIwSj8f5/e9/zx//+MeC+VtXV1ezb9++osnhLRUlREVGMBjkZz/7GQMDA1qbMiNOp5Nly5bR2NhIY2MjXq8348L00Ucf8atf/SqrC+sKmW3btnH77bcXVRUkQGtrKy+88ALBYFBrU2bFZDLxyU9+khtuuGHWaQ+liBKiIiIWi/HDH/4wLzyhdLHb7SlRWrZsGT6fb9HCNDIywiuvvMKZM2cybGXxUV1dzb333kt5ebnWpmSUSCTCP//zP9PX16e1KdOoqanhC1/4Ai6XS2tT8g4lREXEq6++yuHDh7U2Y0lYLBYaGhpobGykoaEhVeE2F1JKTp06xauvvsro6GiOLC18zGYz+/btY9OmTVqbkhFCoRAnT57k6NGjhEIhrc2ZhsViSbU+mEwmbrnlFuLxOIlEgkQigcPhKOjWiKWghKhIaG5u5oknntDajIxjNBqpr69PiVNtbS1GozH1/93d3bz00ksFO54nH9i4cSN33HEHdrtda1MWRU9PD0ePHuX06dOp/qKysjLKyspSJdKxWIxIJMLY2FjeXqysX7+e++67T2szNEE1tBYBkUiEAwcOaG1GVohGo1y+fJnLly8D42XXNTU1NDY2otfrOXToUMGX72rN+++/z+XLl7n99tu55pprCqLnaHBwkMuXL3Pp0iWGhoaIx+N4vV6MRiN6vZ729vY586RJ76S1tTXj759169ah0+k4e/bsgu63du3ajNpRjCghymNef/31kknMx+NxWlpaaGlpweFwKBHKEKOjo7zwwgu899573Hnnnbjdbq1NmkI4HObKlSupi5Kenp4lHW9sbCwVFss0e/fuxePxcPz4cQ4dOsTQ0NC897FYLEWXr8sGSojylCtXrnDs2DGtzcg5er2e4eFhrc0oOgYGBvj3f/93ampquOmmmzQL18Xjcdrb27l06RKXL1+mra0t46KRjTzS1q1bcTqdwHjI8+TJk3MKUUNDA2vXruWaa67B4XBk3J5iQwlRHjIyMlK0Ibn5cLvdeV2mW4jU1NTQ3t4OQEtLC2+//TY33ngjO3bswGKxZPWxpZQEg8GUx/PRRx9ldS+UwWDA4XDgdDoZGxujr68v7enetbW11NfXT5vGbbFY2Lp1KzBeSv673/0Oq9XKqlWrUpO6J39fuXIla9asydpzLEZUsUKe0d3dzdNPP523/ULZpq6uLu9mjRU6Pp+PwcHB1M6fJBaLhR07drBmzZpU4j8ej0/7Sp4jkrPTZvq3TqdL3TaRSKRE4PLlywwODub2CV9FWVkZDocDvV6PEIKenh5MJhMulyv1miQLZRwOB9u3b1fbXDOEqporQFpbW3n66acZGxvT2hTNaGxsVJVyjAty8mSfPLknEgmEEOj1ekwmU0oI+vv7ZzzZ6/V6amtrczIOqqKiYsn5nVxiMpnm9Myqq6u58cYb2bhxY8kOLM0EqmquwLh06RLPPPNMzpaE5Sv5tgZAK4QQaXuGQgjsdjtGoxGj0YjBYECn0zEwMJCzmYSFtiNpvvBgR0cHzz//PK+99ho7d+5k8+bNas14FilKIYrFYly5coWqqqqCSBS+9957PP/886pSDPK2FyTX9Pb2YjAY0hJmKSUjIyM5sGo6ZWVlGI1GWlpaNHn8bDM0NMTZs2d59dVX2bx5M9u3b6e6ulprs4qOohQivV7Pyy+/TE9PD06nk6qqKqqqqqiurqa6ujov9u0kk7hnzpzh0KFDmtqSTyRzY+Xl5TgcDoQQjIyMYLVaEUIQjUbp7OzU2MrsMzIyQl1dHW1tbXk98NPlck0TIZ/Ph9lsZmxsjGg0mtrdk0gkMBgMOJ1ODAZDQYiXTqcjGAwSj8c5deoUp06dor6+nj/5kz9h3bp1KmyXIYpSiIQQbNq0iTfeeIOhoSGGhoa4ePFi6v+tVmtKlJJfuZgeHYvFaGpq4sKFC1y8eDEv52ZpiU6no6Kigt7eXvr6+qa8PkIILBYLo6Oj2Gw2KioqGB4eLuo+q0J4bjMN+DQajXOGFZPl+Q0NDZqtM0mXurq6aTYm+92cTifbtm1j69atBTu9Il8oSiGC8aVhb7zxxoz/Nzo6OqWrH8bnc13tOVVUVCx5km4oFOLChQtcuHCBS5cuZbV0tdBJJBKz7p+RUuLz+WhubiYUCqVODh6PB6fTSTAYLLr+o4qKirw/UV9diQekPf178n31ej0ulwubzYZer8+L5z1fBefQ0BAHDx7kD3/4AytXrmT58uWsWLFiSQN+S5Wirpr7wQ9+QEdHx6LvbzAY8Pv9U8SpsrJy3g9aMBjk/PnznD9/npaWlrwOrRQSTqdzziZCv9+P2Wyms7OzoAXfarWmRDffqa2tZWhoKDWxYXR0NO1lhWazGa/Xy/Dw8LSqP7PZjNlsToXRo9EoAwMDORt6uhRvzeFwpERp+fLllJWVZdi6wkGVbwNvvvkmr732WgYs+hidTkdlZeWUsF5lZSWBQCAlPoVUxlpoTG7OnA2dTpdKKHd0dBRMEUhVVRV6vb6gbHa73UQikZwIhM/nw2q10tbWltV18JkOGXq9XlasWMGKFStYtmxZ1puI8wklRIwnvv/hH/4hAxbNjRCC2tpawuFwya+uzjbV1dUL8nLNZjOVlZXEYjECgUBenuCdTidGo7EgckIzUVZWRjwez1lotLKyclG7uSorKzGbzUgpaW9vn/ZeSH6Os9lQLYSgpqYm5THV19cX3SLDySghmuCJJ57IaYjD6/Vis9lob2/P6lVbKeNyuRbVrW8ymfD7/cTjcQKBQF78fWw2GwaDQfPpA0vF6XQihMj68zAYDFRWVs7pFTudTmKxGFarFavVik6nQ6fTEQgEUs3ibrebWCyWEk+j0YjH45k1R5ktDAYD119/PbfddltRVuCphtYJNm7cmFMhCgaDBINBLBYLdXV1RZlE1xq3272oE14kEkmVDBsMhtT0gkAgoEkjbTIPUgzVkzabLSc9YLFYjPb2dmpra+nq6prWAD45rDaXPf39/VitViorKxkdHcVkMuVchGD8+Rw7doxAIMB9991Xsk2zaZWECSHcQojnhBDnhBAfCCE+IYTwCCFeFUJcnPg+bda5EKJRCHFSCHFKCPG+EOKxzD+FudmwYUOuHxIYH0ff1NTE8PAwtbW1VFVVaWJHMRIIBJYczojFYrS2ttLW1qZJg6LBYMDlchWFCAE59+ra2tpS6+gnn7wXEnodHR1lbGyMsbExzQftNjU18cMf/lATMcwH0q1Nfhz4rZRyHbAF+AD4DvC6lHI18PrEz1fTAeyQUl4L/AnwHSFEzdLNTh+Hw8GyZcty+ZDTaGtro7Ozk/Ly8tR038kkZ4fZbLbUcEbF7ITD4YyKR09PT04ncFRXV+N2u4smn+hwODRpMu7v76epqYlIJJLa9LvQOY02my1vxmr19/fz4x//mA8++EBrU3LOvJeVQggXsBt4BEBKGQEiQoi7gFsmbvYk8Abw15PvO3HbJGbSF76MsmnTJq5cuaLFQ08h2aRpMBgwm83Tph1PLjk2m83Y7XbMZnNqavBkhBBIKad8L4Ry30yRySqt0dFRPB4PRqMxqyelZEHLUloK8g2v10s4HNY03xaLxRb93jeZTBm2ZmlEo1GeffZZbrnlFnbv3l0y/UjpxDdWAN3AE0KILcBJ4C8Av5SyA0BK2SGEqJzpzkKIeuAlYBXwn6SUM2YZhRBfB74O43HeTLJ+/XpeeumlvOnnicVi8+YkwuHwjM2Cs1FTk1NHU3OCwWBGJz739vZSW1tLW1tbRo43mbq6OkZGRrJybC3xeDwMDQ0VbM+WECJvvdI33niDrq4u7rrrrrwTy2yQjodiAK4Hvi+lvA4YYeYw3IxIKVuklNcwLkQPCyH8s9zuB1LKbVLKbT6fL93Dp4XNZmPFihUZPWa+UYrhvEyfANva2mhsbMzoMevr62ltbS2aXFCS5OK5QhUhGJ9ckc9Dds+ePctPfvIT+vv7l3wsKSVvv/02v/nNb7h06dI0D1ZKycDAAFeuXOGdd97h/fffX/JjLoR0PKJWoFVK+dbEz88xLkQBIUT1hDdUDcxZ2C+lbBdCvA/cNHGMnLJx40YuXbqU64fNGUsdRVRoOByOrCTHm5qaFtyrNBmXy4XL5UIIkQqbejyeKT1Cer0+L0rHF4vZbEan08055aIQKITVFYFAgB/+8Ifcd999S7pIamlp4cyZMzQ1NXHs2DHMZjPLly9Hr9fT09NDT08P8XgcnU6Hz+ejurqa1atX58wbm1eIpJSdQogWIcRaKeV5YA9wduLrYeBvJr5P220thKgDglLK0Ymqup3AdzP5BNIlGZ4r5BPAXORjo2Y2cbvdWSuLTydPVF1dnZrM7HA4MJlMCCEIhULTGiJra2tpaGigu7ub0dFR4vE4brc71XszOjpKb29vwbw3vV7vvNMtCoF8CdXPRygU4qc//Smf+9zn2LJly4LuG4/HeeONNzh8+PCUc0Q4HObcuXNTbltXV8f+/fs1GUmUbg3st4GnhRAm4DLwKONhvWeFEF8FmoF7AYQQ24DHpJRfA9YDfy+EkIAA/k5K+W6Gn8OcDA0Nce7cOc6ePVswH/TFUGoL5bKZxO3p6Zl1zIvL5cJut0/xmOar1JopN9Tf3z8l5KLT6dIaX6QlyannxVIUo9UOp8WQSCTmvfCKx+MEg0ECgUDqq7W1Ne1KwtbWVr7//e9zyy23sH379pxGWdISIinlKWCm7tg9M9z2BPC1iX+/ClyzFAMXQygU4vTp03zwwQcFsfMkExRyrH4x5DosJISgvr6ejo6OrIQEE4lEXuf5ksK8GBFyu92Ew2Hcbjc6nQ4pJaFQCIfDQTgcntEbtFgsCy7FXiiZyL3kCqPRyPXXXz/ld6FQiFOnTqVEp7u7e8mRkXA4THNzMw0NDTktgCrKyQo2m41rrrmG6upqAoEAXV1ddHV10d3dvaBKtEIin5OumcZqtWb9JNLc3JzKFdXU1ExZPZENXC5XXntDzc3N1NfXp3Vh19DQgJQSnU5HIpFI3efq92jyb1hZOV5wOzQ0hMViYXBwECEEZWVlqUWJmSabx84GmzZtmpbTMplMDA8Pc+bMmSUfX6fTsXnzZnbt2kVFRcWSj7dQilKIAOx2O3a7fUozq5SSZ555hvPnz2tnWJbI1Xj8fMDj8eSkFLqzsxO73Z4TgVjs2KJcksyHzRUi8vl8BAKBBV3wTR5gGg6HSSQSjI6O4nK5CIVCWentcjqdBSVEly9f5oUXXmDv3r2MjY1x/Phxzpw5k7ELUKfTyZo1a/B6vRk53kIpWiGaCSEEK1euLDohMpvNRevpzUSuQlhSypzkESwWS0H0GIVCIcxm84z/l1zR0NrauujwUCIBw8NW7PYRhBivGMvWFtdCyhcnc4enT5/mwoUL6HS6jL8vBwYGOHDgAOvXr8/ocdOlZIQoHo/z5ptv8oc//EFrUzJOOBzO+0R3Jim2IbI+n68gcplmsxm//+M2wOREj8HBwSU3hiYS8OSTD9PSUk99fQsPP/wkOt24GNXV1dHW1paxKreysrK8n+mm0+moqKjAarVO8dyyFYI3mUx4vV7NJjmUhBC1t7fz4osv5v2bb7H4/f6iGhszFxaLpWD39sxGtquTDAYDHo8Hi8XCwMAAAwMDqXyDz+ejv78/rbBgOBymr6+PWCyW8RNiKGSnpaWeREJPS0s9oZAdh2OEcDhMa2srXq+XSCRCeXn5lK3HC+3JSq7dyLd2B7/fz8DAQGpVejAYXNTOpcUipaSsrIzjx4+zYcMG7HZ7zh4bilyIotEoBw8e5OjRowXTM7AY+vv7cbvdmEwmjEYjOp0OIQSJRIJEIkE0GiUSiRAOhwu+qMHj8RSV52ez2TIy+Vmn01FVVTVlKrmUkuHhYfr6+qad1JLvg+bm5rRXDzidTpxOZ1Zef7t9hPr6lpRHZLdPDT0lX6NkQUMyHD02NkZNTc20aexSytRnPhAIEI1G8Xg8jIyMaD5pW6fT4fV6sVqtxONxBgYGUhfJ2a4UnI1oNMq5c+c4d+4cv/vd7/jLv/zLWcOw2aCohai9vZ0jR45obUbWWehcOovFQmVlZUH2gxTbNkubzbaoeXlutxu73Z7Kl3V3dy9aIJJFAfOh0+mydhEgBDz88JOEQvZUjmg2kqsbksxnU3L3lFZbcO12O263G71ez+joKMFgUPMZd2azOeUZxuNxGhsb8fv9VFRUUFNTk1MRgiIXomL2gpaCVlddmSI5OqcY6Onpoa6ubs711DqdDo/Hg81mm+LlZKqEvbOzk9raWnp7e2f1mJcy9ihddDpwODJfHJLcPQXjUyGy7RF5vV7sdjuJRIL+/n6Gh4fzqnnWZDLxH/7Df8DpdGptSoqiFiLF7OTLDpaFkgwlWa1WzUMsmaK1tZWGhgZisRhGoxGDwZD6+4TDYbq6ujI2ZXw22tra8Pv9GI3GGfNFAwMDVFVVabJ3KJMEg8Gs9hAZDAZ6e3vz+r25Z8+evBIh0Gg/UK4olqvmbFDIXpHf78/rD/piaG5uToWNkhMM2tvbU5MIckEgEGBwcBCfz0dDQwOTp+B7vd6sidB42badXDxNs9mc9arL8vJpy6rzipUrV2ptwjSUR1SiFOK07uSYnUL15mbD4/GkGmebmppSv4/FYnR0dOR8CsDk/IXL5cLr9WZtseRsZdvZwu/3ZzU3qtfr834q+cDAgGaNq7OhhKhEKbTNj3a7HavVmjqJZKvRMZfU1tYSi8UIBAJzJtGTaxe0KDkeHBzE6XRmzSubrWw7W2R7Akk+9/StWrWKVatWUVdXp7Up0yi8y+IFkG+9AvlEIQmR3+8nkUhMyZMkBzPmEqvVSm1tLVVVVUvyKJM7X9ra2tLqbevq6tL05NHW1pa11zpZtq3TxWcs2840Foslq8eHcY9ycuNvvtDV1UVFRUVebnwtao9I5Yhmp1CEqKGhYUoD42QWMogzE0yecWcwGFJ9O319fQsKx9TW1i7Y5paWFqxWq2Z9YNkqllhI2fZSMRgMOcktRqNRAoEAer0ek8mEz+dL7aLSksHBQZ566iluueUW1q9fj9PpzJvlgCIfT9bbtm2TJ06cWPJxhoaGeP3113nvvfcKarZULmhoaCAejxOJRBgbG2NkZCSvPEi/3088Hp/3BGgwGLBarVmNywshUhMIZlu34Xa7cblcRCIRurq6pryWiQSpE21j4+JDio2NjVNySLlmKeHQya+BVtdAWoZzfT4fQoicTkuYD4fDwbe+9a2s9gwJIU5KKWdaITT1dsUsRElGRkY4efIkJ06cyPtEopYky6JhvNfAYDCkqraSc8UmM7l7/eqfpZQYjcYFf/A9Hg9ms3lBPSt+vz/j45vKy8txOp1Eo1G6u7txuVxpN0MaDAacTidms5n+/kF+9KMHuXy5isbGVh566IlFJ+N1Op2mU6Pn63eajVwXJMyExWIhkUhourcrF71YC2XHjh3s3bs3a8dPV4iKOjSXxG63s3v3bnbu3Mm5c+c4duxYwSe6s0EoFMp4Mre2tjatydJlZWU4nc5FnegyMaXZZrPh8XgA6O3tpa+vj76+vtT/O53OtIUoFoul7js8bOfyZT+JhI6mptolJeMTiURqXlyuqa+vX7TYZ7MgwWg0UlVVNe0iKTneKhaLEY1GKSsry1rlX7oEg0HNvdqrOXr0KDfddFNOcmdzURJClESv17Nx40Y2btxIZ2cnb731Fu+++64K22WRZKJ7NpFIrp9uaWlZ0gm2ra1tQWXOBoMBn8+H0WhkeHiY3t7eOUV4sZ70fDPUFkogEMhpXgxI/X0WS6Zfg8nU1NSkdWIfGhqioaGBrq4uzXro3G533m1STnqJWgtRSYTm5mJoaIjvfve7OXmsUqauri41fBLGK9AqKytpa2sjFotl5DEqKyvnjMFXVFRgt9sZGxtb1FrlxfbzZDo/YrFYiEQiOcvppROSm+85Xv3/S31NFjuA1Wg0UlNTQ09PT1bH7ng8nlR1mtlsZmxsTNPp/0KIVL9a8stms+FwOLjmmmuylidSobk0cTgcWptQErS2tuJ2uzGbzZjN5mnNm5mgq6trSujD5XLhdrtTpd/Jr8XidrsXJURXz1Bb6kk4OXE6F70qer1+3kkE6eSAJr8GS8kZVVZWYjKZ6OjoWJSXGo1GU+8Pq9VKWVkZZrOZRCJBKBSiv79/SRESn8/HyMhIXq0qWbZsGXfccUdqJXs+UvJCJITQrFmwlHAi5L81AAAgAElEQVQ6nbhcLtra2rIaCm1paaG+vj61YyeT67czkT/LVOI+V1PIk6tF5mKhOaCF3t5kMlFVVcXg4GBGqs6SFwJSjsxYUu12u7HZHIRCdiyWQQYHB9L62/t8vlTBz+DgYMaG0i4VIUTejx0qeSGC8as+JUSZx+PxYLVa0el0tLa25qRiMZFIZC1/0t3djcPhWNKsskwl7nOVaxgdHWV0dDRVRSilZGhoaMpJdqE5oPluX15ezsDAQEoEQ6FQxoqL0rkQ6O3t53//77um3MZuH/eeTCZTynvq6+tLVYmWl5fT29ubGo+UTw2tbW1tnD17li1btmhtyqwoIWJciIptfplW6PV6amtrGRwczKvwRKZYqhAtNnFvMBiorKxEp9MRDAZzPgX76ipCv9+PEILOzs4FN6XOd3ur1crY2BhCiIxfvKRzITDTbXS66d6T1WrF5XJhMpkYHByc4uknV5zH43H6+/s1aWb1eDxs376da6+9Nuf7hRaKEiJILRdTLJ6amhpgvGdLlcbPztUnYSlhZGT2E3hFRQUwXvqbT/PLkon3ZCHDQncJzXX7zs7OVITCZDJRWVnJ8PBwRkJd6VwIpHuxkPQWZyNZ4OF2u5FS5qxab/Xq1Wzfvp2VK1cWzAQVJUQoIVoKjY2NxGKxtHqFioFMzOlKnoTTCRP19PTQ0NBAX19fztoMFlJMkY3Ks8lh8kgkQmtrK0KIjExGSMd7y/TYof7+frxeL+Xl5cRiMQwGAyaTiXA4TG9vb0bCrGazmWuvvZYbbrgh7yZrp0NJC1EsFuP48eOYzWb8fj86nQ6dTpe6itDr9bS3t6uw3SzU1tbmVXNeLshkqCjdfFFzczM1NTV0dHRkfX7iQospnE7nlJBdtrh6isdSSMd7y/S22Llm3LndbhwOBzqdjrGxMXp7e+dtaTCZTFRUVGC1Wlm7di1btmzJy2Gm6VKSQiSl5N133+V3v/vdvOW4VquV6urqWQdvlhJCiFT57OjoaMl4QZPp6+vD6XTOK0jpeBULyRe1t7fnpJE1XXGsqKjAbDbnLAxrNpsXVHqfD7Pt0qW/v39a2LG8vByz2Zz6rPX19VFeXp4afNvd3U17ezsbN27khhtu0MjyzFFyQjQ2NsbTTz+d9iiZ0dFRmpubJ0o6bXkVp88lQgi8Xu+sTXkulwubzTblg5OpRtV8Y2hoCI/HM6UYw2QypUIss3kVV58cFxoCamlpyciImLlO0umKY39/f07/vgtZaJeJEnmthWwmL3PywsIkjY2NuTAn65SUEIVCIZ5++ulFiUnyqiXdqdDFhpQSq9U6bRVBsm/CZrNNuVr3eDzEYjGGhoaoqKhgaGiIyspKEokE0WgUg8FAT08P4XBYi6ezZBwOB729vVRXVxOPxzEajSkPcSavwmYbmfHkuNAQUFNT05KaWa8+SX/jG8/gdNoxm82p3qT/8l9+x8iIDbN5gIEBy7Qk+2KHny6WioqKBXleSy2Rz4chremSj2u/F0PJCFE4HObJJ59cckPc5GqhgYGBkprm3dLSgs1mY9myZalpBclBqVfHwHt7e9HpdOj1erq7uzGZTNNOJm63e4oQ2Wy2VD+G1rtb5iMUClFVVUVHR8e0JPpMXsXISPonx/muxru7uykvL19Ububqk3Rvr55IZO4dPclm5OQ09o8++mjBj5tLljrbLtdbYxeL2+1ODeotdEpCiJILoWZybRfLeMmqjoaGBgKBQMFe2adLspijr68v7SnGyQnIMHMDZn9/P/X19QghkFLS0tJCKBRKe2K3lkz2iK/OHc4Uckv35JjO1Xg0GiUWi1FbW8vY2BjhcBir1Trl/S2EwGKxoNPpMJvN2Gw29Ho9odDogk/SQ0NDqQuusrKydF+iJWMwGGhsbOTSpUsLut9Sq96yOaQ1kyxfvlxrEzJGUQ89TSQSdHR08Nxzz2V13IZOp8PhcGCxWDAajalycCklkUiE3t5ehBCUlZUxMjJCPB7H4XBgMpkYGhpaUoNktkh20sdiMXp7e3PSA+HxeHA6nVkfSJlJLBYLQoi0PLh08g7Dw3a++93/h0RCj04X5z/+x/8949V4VVXVlKbW8vJyhBBEIhHC4fCclZ5LyX+UlZWlhr9mcx1Fsil6aGgIu92e8/eD1jmidLDb7Xz5y1/O6xlyJT/0VErJ448/ntFZY7ORSCTSmms2+Yo16UGVlZVhMBjyJrFvNpsxGo3TOumzTfLEWmjTGCorK9POX6STD0rnatxqtU4TgYX8rZZSmpwUoGRoejEYDAbi8fiMVajl5eWYTKYpHrHX6825EGW6fDsbjIyM8OSTT/LQQw9RVVWltTlLomiFSAhBXV0dZ8+e1dqUORkYGNB0hfHVVFVV5bw3yOfz5XxkTb6STljJ7XZrvukzEolMqRRMF6/XSzQaJRwOp9ZnB4NB4vH4girj0qEQvJqlEgqFUmKUnG5SiBR1aE5KSU9PD+fPn+fChQu0trbmbS9QQ0MD8HH8PxqNajK9N7mzZ3h4eM4mvNkwGo14PJ4F7fvJdRXWXCRDbcm/wXzkekldvl20LOQCoqGhgdbW1hnfF3NFBXw+HwMDAwsSvVxXvmktemazmS996UvU1dXl/sHnoORDczDuFfl8Pnw+H7t27SIUCnHx4kUuXLjAhx9+mFfbEq8+uRiNRk2S9pN39thsNrxeL4lEgu7ubuLxOC6XC6vVCowXgSTzW8npFB0dHQQCASwWS1phq4qKiow9R6vVSnl5eSpHl5yQIYRIFU4k10cnT4bJPJ6UkpGRkSm5ML1ej8PhwGq1YjQaU8dJ9kklS9FzSTwez9kuovmYbz1EEpPJhM/nm/O9MFdoOll1WVdXl1rZkhx3lFzj0t3dPaVgKJeVb/lQ7h0Oh/nZz37Ggw8+mLqoLSTSEiIhhBv4EbAJkMBXgPPAM8Ay4Apwn5Sy76r7XQt8H3ABceB/SCmfyZDtC8Zms7Flyxa2bNlCPB6nqamJCxcuYLPZqK2tRUpJKBTirbfe0vyDHo1GU2u2tZrqkCzNnszVuSOXy4XT6WR4eHjKiUYIMSVnltwpA+Mn/3g8zuDg4JTqM51OR3l5OTabDZ1Ol8pH6PX61JfBYMBgMKDX61NCk0gkGBsbY3h4OKNl3/F4fNakvE6no6KiApvNlrHHS4ekaFdXV5NIJDTb+mk0GtPKESVDcUu92EjOnJuN8vLyKUKUy8q3fCn3jkQiPPXUU3zxi18suIq6tEJzQogngUNSyh8JIUyADfh/gV4p5d8IIb4DlEsp//qq+60BpJTyohCiBjgJrJdSzhlzyuWq8JmQUvK3f/u3mu22v5pczRnLJB6PB4fDwdjYGAaDgYGBgWkJZ4fDgdvtRggxYy9SIXD1hIVcYjabicVisw5DXWy4KN376XQ66uvrZ133PlcoLtN4vd7URVPyYiRX4TIp4f/+3489okceeVLTnJTBYOD+++/Pi2bXdENz8wqREMIFnAZWyEk3FkKcB26RUnYIIaqBN6SUa+c51mlgv5Ty4ly301qIAP7u7/4ur0qIc52LyATJstKkJ5P0YKLRKAMDU7deWq1WvF5v3uSK0kXr/NZsj7/YcNFi7jd5+y58vLpBi9fF6/VqckGjRY4oeaE300WAXq/nvvvuY82aNbkxZhYymSNaAXQDTwghtjDu1fwF4JdSdgBMiNGcxexCiO2ACZixO00I8XXg60BexDjzbWNrcgV2T09P3k8dSDJ5ikVylcFskyi0yLcslfLycs2Fc7YVJosNFy3mfsmG16qqKvR6PUNDQ5q9LsFgUJOCjlyXe0++MLXb7am+RCEE8Xic0dFRDh06pLkQpUs6QmQArge+LaV8SwjxOPCdhTzIhMf0M+BhKeWMZ3gp5Q+AH8C4R7SQ4y8GKSWxWIxIJDIlWR0Oh4lEInl5sm9paUGv19PQ0EBXV1fehA7TYXR0FJPJNKtnp0XZ+FKx2Ww57bWaibGxMYxG4zQRX2yOZCm5lbGxMUZHRxc0ZSQbnkRzczMVFRX09vbm3QVlJrh68O3IyMiM0Zs9e/bk0qwlkY4QtQKtUsq3Jn5+jnEhCgghqieF5mYc4jYR2nsJ+C9SyqOZMHqpxONxxsbGiMfjJBKJ1Ey0Y8eO8fbbb+f1lXk8Hqe5uRmTyURjYyOtra05W5i2WJxOJ8FgkEQiMWvYJN+fw0xkc7JAugQCAUwm0zQvYLFjbhZ7P5vNRjweX7AIZavaLLlQcGhoCKvVil6vJxAI5FWl7EJJ9kamc8Hm8/lYv359DqzKDPMKkZSyUwjRIoRYK6U8D+wBzk58PQz8zcT3A1ffd6Kw4Xngp1LKX2TU8iWg1+ux2+3Tfv/pT3+am2++mZMnT3Ls2LG8HmgaiURoamrC4XBQXl5OR0dHTqczJEu7o9FoqrTb4/FgMpno7Oykurqa0dFRDAYDNpstNXpmNpIeUyGdKPJlEVkkEqGrq4sqnw9XUxP1nZ142tpoWbeOo770jzPZO1lImMloNGK1WqddZMzn7WS72iwpzEmvNXnxNltxRT5jNBqpqKhIO0+8e/fuglkTDun3EX0beHpCWC4DjwI64FkhxFeBZuBeACHENuAxKeXXgPuA3YBXCPHIxLEekVKeytxTyCxWq5Vdu3Zx44038v7773PkyBHNSmTTYXh4mOHhYQwGA/X19YTD4SVPGL+a2traKbmIkZERYrHYlA+F1+ult7cXKSVms3nBnf/BYJDq6mrNJwYshOTaBK1oDIfxDw1hikTwX7nCij/+EdukAhD/e+9x9BvfIB03Y7HeiRCCioqKaX+3dI6X6+GiyYuhQhOh5KqVdD8bFRUVbNiwIctWZZa0PkkTwjFT5cO0IKSU8gTwtYl/PwU8tRQDtcJgMLBlyxauueYarly5wpEjR7h4cc5iP01JCoPf70/r9kII/H5/ygsZHR1N7QlyuVyplel9fX20tbVRV1fH2NjYrHuYJl8NL3YSeUdHx4K79bVEKyFaNzLC7n//d6rnmUrtDQTY0d3N4TTeE4v1Turq6ma8Sk/neEudkr0Y8sWLTRe32008Hl9QJeDu3bvR5esCpVko6skKmUAIwfLly1m+fDk9PT0cPXqU06dP5+1V1XzuuM/nS60MmO2En2xETa5nAFJVUH6/H71en7W+ptmqwPKRXH7Y/dEoq7u7MUUi7PjZz9CnmVNbffgwh++5Z97bLcY7mas6Ld3j5braLBgM0tjYSHd397Rm7XyjsrKSgYGBBV3Yeb1eNm7cmEWrskNRz5rLFp2dnbzyyis0NTXlXVWOXq+ntraWYDCI0+nEZDIhpSQcDjM4OJixSjubzUZFRQX9/f0ZnXCeT7PU5iMXtloSCW47fZpNv/41xgUW0fRWVvLTb3yDgRnEfab8TboVbEajEZ/PN+/0Ea3nr82FEILa2tpUfi3fqK2tpaOjY8Hnl7vvvpstW7ZkyaqFo2bNZZGqqiq+/OUvEw6H6e3tpampiY8++ojLly9r5iklPZ1IJEJnZyeRSCSrDbmhUIhIJEJZWRnDw8MZE+RgMEhdXR2BQCCvqxezjpTsbm9n6y9+gWuRw29dwSA2Kbm6tm+2/E063onD4cBsNqc1AiufVylIKVNevsfjwW6309nZmRfvucVe4Hg8HjZv3pwFi7KPEqIlYDabqa6uprq6mhtvvJFYLMYvf/lLzpw5kzMb0h0umkmSoY1s5HKSPREej4dQKFRQvVKZYv3ICDf9279Rffnyko5jiMf5+n//70SMRmJGIzGTiX/51re4NFa2qHxQWVkZkUikIEcxzUVvby+9vb2pMvhsL/2bi6V42TfddFPB5YaSKCHKIAaDgY8++ijrj+N0OvF4PIyOjtLd3Z1TEcpV42lvb2/eh+kyHdb2xWJ86uBB1rz5ZkaPa4pGMUWjEArxwP/5PwRWr+H1ijs5072GrbbT/OULP0CXiDNWVkbI52PI62WwvJyg00mH2czQxMnN5XIV3JiphRCJRFLvt+SA3lwVzuj1eqqqqhb9fne73QXrDYESooyzf/9+nn/++azsErLb7Xi9XgYHBzWbQtDZ2UllZSVjY2NZ336bb/m3q9GFw3jjcSJAVKcjLCVyEVekpkSCve+/z+Zf/hJzFvuoEghG+4ysOHaMl7gVgcQ/3IX4cO77jVksDFVUcOCRR0DjkvVc0dnZSWNjY04ey2w2U1ZWtqQJ5TfddFNBFfpcTWm8q3JIQ0MDf/qnf8pPfvKTRR+jtraWoaEhIpEIdrsds9lMPB6nq6srLwaxJpO7NTU1M07VzhT5OGZpMv6zZ3nkn/5pyu8SQtB03XW8/clP0m610qvTMVum3pJIcGNLC1sOHMCd4QneMXScYz0beB/dxM+7OcQxbsDBCCPY2cGbHOSTCOb27CxjY1haW7n23DnaNm3KqJ35TC6q6hwOBwaDYUkFE2VlZXlVoLAYlBBlmPb2dv7lX/5lQfcpKyvDarViMpkYHByccmWUzzmS9vZ2GhoaUiOTMonNZsNkMqHX6/N2/E98hitQnZQsf/ttlr/9NgC/+E//iU6LBXsshj0WwxqJYItEqPvgA5YdOYJlkT1XcxFDRwVBBiijjAG6qOBmfs9RPgEIBigDBIfZSTc+/DNP55rGssOHoQSEyGKxUFZWlvVGdo/Hw9jY2JKjJ4XuDYESoowipeT48eMLqvtP5kHyYW7ZYmhubk7NwOrv709tbF0qoVAoNb8sX4mm8eG/93/9rxxYMk4CQTc+uvClxGaAMg6zk+NsBwQgcTLIKDZ28CaVaYoQQEV7O3/5j/84j/80TtRm43tf/eqs3mA+4/f7CYVCmM3mRTdnz0dVVRXBYHDJVXoul4trr702Q1ZphxKiDCKE4K677uLmm2/m+PHjvP322ylPwWg04vF48Hq9XLlyBYfDgU6ny+tkfLokS2GFEBndm2SxWDJynGxh0Fgkk8JTSRcSwa38jsPsZAdv4mKAwQmP6Cb+wE7e5E12cgPHOMRNBCfut1CZKEs3hNjbS200SluBTTIApuRfrVYrTqcTs9mMTqdbUHtGcjoJfNwcHo/H0el0tLa2ZqTYZdeuXQXvDYESoqzgdrvZu3cvN998M4FAALfbjcPhSL0p29raeOKJJ/L6an8xSClpaWnB5/MRi8WWvCJheHgYp9OZt8Nn12rYdJ24Snh+zv0cZicxjBxmJ1dooA9vKkd0kE+mREtA2uG4pbC8p4e2mpqsP042GR0dzdtcpdPp5LrrrtPajIxQmEXnBUJy/47T6Zwyeqe2tpZ77723oKbjLoTu7m6sVuuSj9Pf3088HqexsRGj0ZgByzKHMx5n+bFjmj1+N76U8PyRnZxnLZ/gTQxE2cGb1NDJpgkRAtAh8S/CA1oKlUXg7eczu3bt0nzwbqZQI340ZHJ58vDwMB9++CEXLlzg0qVLeTvLLl30ej1+vz+tDvx0cLlcOJ3OJZW4ZpLbL17kxqef1uzxJXALB/kjOxFAHAMuBjjLemrozKngzEa/x0PrDTcQs1iIms2pr7DRSMRs5rLTSV8RhJW0wOFw8Bd/8Rd5L0RqxE8BMLkL2uVycf3113P99dcTjUa5cuUKFy5c4Ny5cxkrAMgl8Xic9vb2jK12GBwcZHBwEIvFkqo26s1wyfNCqJ2oitMKwXi47Q/cxK28AQgGKaMPL7Xkx/Ryd28v7pdfnvX/f/eNb3CoujqHFhUPO3fuzHsRWgjF80yKCKPRyOrVq1m9ejWf/vSneffdd3nzzTdnXcGQz3R3d+Pz+eju7s7I8cbGxmhvb6empkbT0m59Hswk0yHZzR8oYyBVqr2B97U2K23yMRpTCNjtdrZu3aq1GRlFCVGeo9frufbaa9myZQvnzp3jzTffzJvwVDrEYrGM5IuuxmAwaFrsocuTQhMd0IN3SvNqwVCkOdJss2PHjrzLmS4VJUQFghCC9evXs27dOq5cucJrr72WsfxLtslGUYaWYblbm5qoysFMwXQxkGBTAXlCSaQSogVjs9nYtm3elEvBoYSowEgu6nvggQf43ve+l/fLvWC8LyOT/UUwXtzhcrkQQkwRupl+vvr3VquVrq6uRU+DqDl7dgmWKxSL58Ybbyy4LbPpoISoQLHb7ezbt49nn31Wa1PSorW1NeNrwJcydDU58r+trW3BIb5yjQbOFh3KI1oQQoii6Ru6moIKKSumsn79eq655hqtzUgLKSXBYJCaPGlwTI78t1qt1NXVpX0/UzxOeR5u9FQUP42NjTgcDq3NyApKiAqcO+64A6fTqbUZaRGNRmlvb8fv9+Pz+bQ2BxgP8bW2tuL1elM7aOZi1dgYujxfT1EoVF26hF69lmmzfv16rU3IGio0V+BYLBbuuusunnrqKa1NSZvkVOPa2lqGh4fzYuBrcutoRUVFqiJpplxTw5EjuTeuSNn0yitUnz7NkS98gZMej9bm5D1KiBR5zcqVK/mrv/ormpubOXXqFOfOndPapLRIlqHX19cTDAbzovBivl6tnXlULVcMeAMBbv/e9zj5n/+zyhnNQUNDQ8FEPhaDCs0VCVarldWrV+P3+7U2ZcG0tLQQiURobGzM+4qg8itXtDah6DDGYlhVc+ucFLM3BEqIiorBwUHeffddrc1YFLFYjKamJhKJBDU1NdTX12M2m7U2awoikcDT2qq1GUVJucoVzYkSIkXB4Ha7+cpXvlLQM6hisRjt7e0pL6mhoUFrkwAoi8e599AhTHkw2qcYcUUiWpuQt9TW1lJWVqa1GVlFCVGRYbfb2bFjh9ZmZAQpJc3NzTQ2Nmpmg4jHue3iRb7x93/P+oMHNbOj2HFkaRNqMbBhwwatTcg6hXvprJiVXbt28fbbbxfk1O6ZaGpqSq1UzyVb+/r4xDPP4M1gE65iZmyhEJSXa21GXlLsYTlQHlFRYjQaufnmm7U2I6M0NzdTW1ubk8eqD4d5+IUX2Pf440qEcoR1ZERrE/KS6upqyktAoJVHVKRcd911HDlyRNPhoJmmo6ODyspKurI02cAej3P7iROsf+UVDHkyXbtUsBSJ955pSiEsB0qIiha9Xs8nP/lJnnvuOa1NyRiJRIL+/n7Ky8vp6+vL2HFFIsGtzc1c/9xz2NUJMePEdTrGrFYiNtv4l91OxOEgYrcTtdkI22y05ElRSr5RCmE5UEJU1GzYsEGT3Eo2iUQiGI1G7HY7IxkI52wZHGTHL35BZQYngys+JlhVxT994xuqWXUR+P1+vF6v1mbkBJUjKmKEEHzmM5/Jyj4gLRkZGcFisSxpOVh1JMJDv/oVd3/3u0qEssihBx5QIrRISiUsB8ojKnr8fj/bt2/nrbfe0tqUjBIMBqmurqazs3NBK6ct8TifPnWKDb/5DcZYLIsWKs7fdBOnXS6tzShYlBApiopbb72V999/v2jKuZN0dHSkH3pMJLi5rY2tv/gFziXsMVKkx4jDwa9vvVVrMwoWn89HRUWF1mbkDBWaKwHMZjO333671mZkhebm5nmnL2wcHuYbP/0pt/z4x0qEcsSxL36RQZ06vSyWUilSSJKWRySEcAM/AjYBEvgKcB54BlgGXAHuk1JOK2USQvwWuBH4o5RyX0asViyYdHbtFCrNzc3U1NTQ3t4+5feV0Sh7X3+dVUePamRZadKyYQN/yJMFiIVKKYXlIH2P6HHgt1LKdcAW4APgO8DrUsrVwOsTP8/E/wIeWqqhiqVRKKshFktXVxeeiZ02pnicz54+zVf/9m+VCOWYqNHIy3ffrQoUloDX66WyslJrM3LKvB6REMIF7AYeAZBSRoCIEOIu4JaJmz0JvAH89dX3l1K+LoS45erfK3JHJBLhzJkzWpuRVWKxGJFIhB19ffzJT3+KK4N9Ror0OX333bTl+SqPfGf9+vVFV+k6H+l4RCuAbuAJIcQ7QogfCSHsgF9K2QEw8X1JEi6E+LoQ4oQQ4kR3d/dSDqWYRDQa5V//9V8phdd0eHiYHpMJcx4s2CtFempq+G2J5TayQamF5SA9ITIA1wPfl1JeB4wwexhu0UgpfyCl3Cal3Obz+TJ9+JIkGo3y85//nCsltMztgt3Oy489Rlyv19qUkkICBx94gLgqUFgS5eXlRZ3PnY103jWtQKuUMtmI8hzjwhQQQlQDTHzPzgAwxaKIxWI8++yzXL58WWtTcs475eX8/qtf1dqMkuKDPXs463BobUbBs27dupILy0EaQiSl7ARahBBrJ361BzgLvAg8PPG7h4EDWbFQsSiampr48MMPtTZDMw7V1PDWF7+otRklwVBZGb/ZuVNrM4qCUitSSJKuH/1t4GkhxBngWuB/An8D7BVCXAT2TvyMEGKbEOJHyTsKIQ4BvwD2CCFahRDF2dCSZ6xcubJoe4fS5bdr13Lmjju0NqPoOfrAAwyrkFxGKPZNrLORVh+RlPIUsG2G/9ozw21PAF+b9PNNi7ZOsSRuvPFGQqEQhw4d0toUzXh++3asAwOsPnxYa1OKkitbtnC4RK/is4Hb7dbaBE1QlzFFzq233sq2bTNdQ5QIQvDzPXto3rRJa0uKjojRyCv79qmeoQziKtHZfEqIihwhBHfccQcbN27U2hTNSOj1/Os999C5bJnWphQVZ+66i44lTEBXTMXpdKIv0WpPJUQlgE6n45577mHlypVam6IZY3o9zz70EL0qjJQRgn4/L5dgv0s2KdWwHCghKhn0ej333XcftbW1WpuiGX16Pc997WsMl2j4I5O8+cUvElMFChmlVAsVQAlRSWEymXjggQdKarz81XSYTLz453/OmMWitSkFy6Xt23mnhK/es4USIkXJYLPZeOihh0r6TX/RauW3jz1GzKDWcS2UsNnMy7fdprUZRYkKzSlKCpfLxZe+9CVsNpvWpmjGabebg1/9KunvdlUAnLrnHrqVgGeFUr44VEJUolRUVPDggw9iKuFJyYerqzny4INam1Ew9FRX88qaNVqbUbQoj0hRktTU1E7TtlQAABO7SURBVHD//feXbMkowKurV3P6s5/V2oyC4NAXv0hCFShkDeURKUqW5cuX8/nPf74kBy0meeH66zl/kxoAMhcXd+zgjKo2zBpWq7WkoxNKiBSsX7+eO++8U2sztEMInrnlFq5s2aK1JXnJmMXCy3umTfNSZJBSDsuBEiLFBFu3bmVPCZ9spF7Pzz/3OTpWrNDalLzjnf37CZZw+DYXlHJYDpQQKSaxc+dObrzxRq3N0IywXs8zDz5IsAQXk81GV309rypxzjpKiBSKCYQQ3HbbbWwp4RDVgF7Pc1/9KkMlfmKA8a2rv//CF5CqQCHrqNCcQjEJIQSf+9znWFPCZbqdRiMvPPZYyU5fiBqNtK1Zw9EvfUltXc0Rpe4Rqc40xTR0Oh379+/nwIEDvP/++1qbowmXrVZ+/c1v8tl//EeMsZjW5mSVIZeL7tWr6V6zhqa6Oi5YrcSVF5RTSt0jUkKkmBGj0cj+/ftZt24dL730EmNjY1qblHPedblwfP3rfOr730cni2MGQ0IIemprCa5ZQ+fKlVz0+cZXOZRw+X4+oDwihWIONm3aRGNjIy+++CIffvih1ubknCOVldgeeohdP/2p1qYsijGzme6VK+lZs4a2xkbOlZUxorydvMJoNGK1WrU2Q1OUECnmxel08sADD3Dy5EleeeUVotGo1ibllNdXrMB2991c/8ILWpsyL/1eLz2rV9O1ejVXamr40GxWxQZ5jtvtLumGclBCpEgTIQTbtm1jxYoVvPDCC7S0tGhtUk755ZYt2Pr7WffGG1qbkiKu19NdX0/PmjV0rljBhYoKNZC0ACn1sBwoIVIsEI/HwyOPPMLhw4c5ePAgiURCa5NygxA8u3s3Dw0OsvzttzUxIWSz0b1qFT1r1tDS0MA5h4Ow8nYKHiVESogUi0Cn07Fr1y5Wr17N888/TyAQ0NqknCB1Ov71zjt5ZGiImosXs/54vZWV9KxZQ2DVKj6qruYjk0kVFRQhpV4xB0qIFEvA7/fzta99jTfeeIPDhw8ji6SybC6iej0/v/9+vvzjH1PR3p6x4444HASXLaNv2TI6Gxs57/HQp8bqlATKI1JCpFgiBoOBT33qU6xdu5bnn3+evr4+rU3KOkN6Pc8++ihf+t73cC3i+Q65XPQuW0bf8uUE6ur4qLycgF6vvJ0SRXlESogUGaK+vp7HHnuMV155hZMnT2ptTtbpnpi+sP/xx7GFQrPert/rpa+xkb7GRgK1tXxUXk638nQUk1AekRIiRQYxmUzs27ePtWvX8uKLLzI8PKy1SVnlI7OZ33zzm3zu8cfRJRL0VlfT19hIb0MDnVVVXHY6GVLFBIo50Ol0OJ1Orc3QHCVEioyzevVqvvnNb/LSSy8V/Yig9xwOQv/jf9AWChHW2hhFweFyuUq+hwjU0FNFlrBarezfv5/Pf/7zWIp8eOjlUIjK+nqtzVAUICo/NI4SIkVW2bRpE3/+53/OypUrtTYlq7S0tNDY2Ki1GYoCQ+WHxlFCpMg6LpeLBx98kDvvvBOj0ai1OVmjqamJhoYGrc1QFBAul0trE/ICJUSKnJAcEfTYY49RV1entTlZo7m5mdraWq3NUBQIyiMaRwmRIqd4PB4effRR9uzZg65IK8o6Ojrw+/1am6EoAJQQjVOcZwJFXpMcEfRnf/ZnVFZWam1OxkkkEvT19eHxeLQ2RZHnKCEaRwmRQjOqqqr4sz/7M3bs2KG1KRknEokwNjamekQUc6KEaBwlRApNMRgM7N27l0cffbToSllDoRA6na7oy9cVi8NqtWIymbQ2Iy9IS4iEEG4hxHNCiHNCiA+EEJ8QQniEEK8KIS5OfC+f5b4PT9zmohDi4cyarygWGhoaePTRR7Hb7VqbklEGBgZwOBwY1J4gxVUob+hj0vWIHgd+K6VcB2wBPgC+A7wupVwNvD7x8xSEEB7gvwJ/AmwH/utsgqVQuFwu7rvvvqIrYujp6cHn86kOesUUlBB9zLyfeCGEC9gN/BhAShmRUvYDdwFPTtzsSeDuGe5+O/CqlLJXStkHvAp8OhOGK4qThoYGPvOZz2htRsbp6Ogo6rJ1xcJRPUQfk86l5wqgG3hCCPGOEOJHQgg74JdSdgBMfJ+p/KkWmLxTunXidwrFrGzdupWtW7dqbUbGUdMXFJNRHtHHpCNEBuB64PtSyuuAEWYIw83CTLGIGbenCSG+LoQ4IYQ40d3dnebhFcXKHXfcUZRTCtT0BUUSJUQfk44QtQKtUsq3Jn5+jnFhCgghqgEmvnfNct/J0yDrgBnXWkopfyCl3Cal3Obz+dK1X1Gk6PV67r333qIsf25ublZhOoUSoknMK0RSyk6gRQixduJXe4CzwItAsgruYeDADHd/GbhNCFE+UaRw28TvFIp5cTgc3H///eiLcJFcW1sbVVVVWpuh0BAlRB+TbnnSt4GnhRBngGuB/wn8DbBXCHER2DvxM0KIbUKIHwFIKXuB/wYcn/j6/yZ+p1CkRU1NDZ/97Ge1NiPjSCkJBoN4vV6tTVFogE6nw+FwaG1G3pBWc4OU8hSwbYb/2jPDbU8AX5v080+AnyzWQIViy5YtdHR08NZbb81/4wIiGo0SCoVwuVwMDg5qbY4ih7hcrqJrU1gK6pVQFAS33XYby5cv19qMjDM6OgqAzWbT2BJFLlGl21NRQqQoCHQ6Hfv37y+6MUAAg4ODatxLiaHyQ1NRQqQoGGw2G1/4wheKcrleMBjE4/GocE2JoIRoKupdrygoqqqquOuuu7Q2Iyt0dnZSU1OjtRmKHKCEaCpKiBQFx8aNG9m1a5fWZmSF1tZW1fBaAighmooSIkVBcuutt7Jq1SqtzcgKzc3NahRQkaOEaCpKiBQFiU6n4/Of/3zRbkFtamqivr5+/hsqChIlRFNRQqQoWCwWC/fff3/RVpu1tLSonFERYjabMZvNWpuRVyghUhQ0Pp+PP/3TP9XajKwRCASorJxpsL2iUFHe0HSUECkKnrVr13LLLbdobUZWiMfjDAwMFGX/VKmihGg6SogURcHu3btZt26d1mZkhXA4TCwWK7o16qWKEqLpKCFSFAVCCO6++26KdYXI8PAwZrO5aPNhpYQSoukoIVIUDWazmfvvv///b+/eYqM4zzCO/18fwPbaJjZ2wcuxCWpFQCUBC1W0QkWIhEAKDcGYqupNquYmlaJKaVUprZSbXCSV2qi9aWMIqFWUXOSA0pzTKk1EFVpBAoRCGnLA1AYvDocSCFVweHsxY3ex1/iw3p2Z9fOTRp7ZnW/97OzsvJ5vxjNUVVVFHaUgzpw5Q0NDQ0neFmMyUSEaSoVISkpjYyN33nknZrluDpx8mUxG9zFKOBWioVSIpOQsWLCA1auH3KGkZHR3d+vqCwmmQjSUCpGUpBUrVrB48eKoYxSMrr6QTGZGXV1d1DFiR4VISpKZsWHDhpLuxurs7GT27NlRx5AxqKur0xXWc9ASkZJVWVlJe3t7Sd90rquri5aWlqhjyCipWy43FSIpaddddx1tbW0le/ICQG9vL01NTVHHkFFQIcpNhUhK3vz587n11lujjlEwfX19XLx4UbefTgB9RrlVRB1ApBiWL19OT08P+/fvjzpKQVy6dInKykqqq6u5dOlS1HEmrfLycmpraweGVCp11bS6UXNTIZJJwcxYv349vb29dHd3Rx2nIM6fP8/06dPp6+vj8uXLUccpGWY2bHEZXGimTp1a0t3AhaJCJJNGRUUFW7ZsoaOjgwsXLkQdpyBOnz5NS0sLmUyGK1euRB0n1lKp1DWLSv9QXV2t4lJgKkQyqdTX17NlyxZ27txZshvqkydPMnv2bLq6uqKOUnRVVVXX7BrrfyyVSuk06hhRIZJJZ86cOaxbt47nn38+6igF09XVxbx58+js7Iw6St6mTJkyqq6xVCpFRYU2aUmkT00mpWXLltHT08PevXujjlIwnZ2dzJ07l+PHj0cdZYiKiooRj7f0P6Yrjpc+FSKZtNauXcupU6diuaGeKMePH2fWrFlFOUGjrKxsVMdcUqmUDurLVVSIZNIqLy+nra2Njo4Ozp8/H3Wcgjl58iQzZswgk8mMq/3ggjJcsdFBfRkvFSKZ1Gpra2lvb2fHjh309fVFHacgrly5wtmzZ2lsbOTMmTMAVFdXj6prrKamRgf1peBUiGTSS6fT3H777ezatSvqKHkzM5qbm2lqarrmcRfdXE/iRIVIBFiyZAk9PT3s2bMn6ihj0tzcTDqdpqWlhXQ6zcyZM6msrIw6lsiYqBCJhNasWUMmk+Hjjz+OOkpO06dPJ51ODwwzZ87UGWVSElSIREJlZWVs3ryZjo4Ozp07F2mWxsbGq/Z0WlpamDp1aqSZRApFhUgkS01NDVu3bmX79u1Fu15bQ0PDQMHpLzpVVVVF+d0icaBCJDLIjBkz2LhxI0899dSEv/a0adOG7OmU8o37REZDhUgkh0WLFtHT08Pu3bvH/Rr19fVD9nRSqdQEphQpDaMqRGZ2DPgU+ALoc/dWM1sC/A6oBY4B33P3If8VaGb3Aj8EDOhw90cmJrpIYa1atYpMJsPRo0dHnLe2tvaqPZ10Ok1tbW0RUook31j2iFa5+ydZ09uA+9z9DTO7C/gJ8IvsBma2mKAILQc+B142sxfcfeRvtkjEysrK2LRpE9u2beP06dMDj9fU1DBr1qyrik5dXV2ESUWSLZ+uua8Cb4bjrwGvMKgQAQuBPe7+GYCZvQHcATycx+8VKZqqqira29s5ePDgQNGpr6/XpWxEJtBor93hwKtmts/M7g4fOwRsCMfbgDk52h0CVprZdDOrAdYNMx9mdreZ7TWzvb29vaN/ByIF1tzczOrVq1m4cCHTpk1TERKZYKMtRN9w96XAbcA9ZrYSuCsc3wfUEXS9XcXdjwAPEewxvQwcAHJe0MvdH3X3VndvbW5uHvs7ERGRRBpVIXL3E+HPU8CzwHJ3f8/db3H3ZcATwIfDtN3u7kvdfSVwBtDxIRERGTBiITKzlJnV9Y8DtwCHzOxL4WNlwM8JzqDL1b5/vrnAJoKiJSIiAoxuj2gGsNvMDgD/AF5w95eB75rZ+8B7wAlgB4CZpc3sxaz2T5vZYeBPwD3ufnZC34GIiCSauXvUGYZobW31Ur6Fs4jIZGBm+9y9daT5dMcrERGJlAqRiIhESoVIREQipUIkIiKRUiESEZFIxfKsOTPrBTrzfJkm4JMR54qXpGVOWl5Q5mJIWl5Q5kKZ5+4jXionloVoIpjZ3tGcNhgnScuctLygzMWQtLygzFFT15yIiERKhUhERCJVyoXo0agDjEPSMictLyhzMSQtLyhzpEr2GJGIiCRDKe8RiYhIAqgQiYhIpBJViMxsjpm9bmZHzOyfZnbvoOfvMzM3s6Zh2n9hZvvD4bmEZJ5rZq+G7Q+b2fw4ZzazVVnLeL+Z/dfMvhPnzOHzD4ftjpjZb6zA9wOfgLwPmdmhcGgvZNaRMpvZA2bWnfWZrxum/Voz+5eZfWBmP0tI5sfM7JSZHYp73pHWqVhz98QMQAuwNByvA94Hbgyn5wCvEPwjbNMw7S8kMPNfgTXheC1QE/fMWa/TSHBX3lhnBlYAfwPKw+Et4FsxzrseeA2oAFLAXqA+qmUMPADcN0LbcoK7OF8PTAEO9L/fuGYO26wElgKHCp11ApbxsOtU3IdE7RG5+0l3fzsc/xQ4AswKn/418FMgVmdf5JPZzG4EKtz9tbD9BXf/LM6ZB9kMvJSAzA5UEWwgpwKVQCbGeW8E3nD3Pne/SLBRX1vIvGHOa2UeyXLgA3f/yN0/B54ENhYm6f/lmRl3f5Pgj6miyCdvvu81SokqRNnCLqqbgb+b2Qag290PjNCsysz2mtmeYnQXDTaOzF8BzpnZM2b2jpn90szKixB1wDiXc7+tRHBr+LFmdve3gNeBk+HwirsfKUJUYFzL+ABwm5nVhF13qwj2ooomO3P40I/M7GDYldWQo8ks4N9Z010UeSM5jsyRyidvjrbxFvUu2XgGgi6qfcAmoIZgYU8LnzvG8N1c6fDn9eF8N8Q5M8EexX/CvBXA08AP4pw5q20L0AtUxn3dABYAL4Rtawm65lbGNW/43P3AfoIuuseBe6NYxuH0DIKutzLgQeCxHG3agG1Z098HfhvnzFlt51OkrrkJyntV2yQMidsjMrNKgg3y4+7+DHAD8GXggJkdA2YDb5vZzMFt3f1E+PMjgmMvN8c8cxfwjgfdGX3ALoL+6jhn7rcFeNbdLxcjL+SV+Q5gjwddnxeAl4Cvxzgv7v6gu9/k7msAA44WOu8wmXH3jLt/4e5XgA6CbrjBurh6r202cKLQeSGvzJHIJ2+utokQdSUc418JBvwBeOQa8xwj91+9DcDUcLyJ4ItbjIOl+WQuJ+iGaQ6ndwD3xDlz1vN7gFUJWTfagT8T7HVWAn8Bvh3jvOXA9HD8a8AhgmOJkSxjoCVr/MfAkznaVgAfERTa/pMVFsU5c9bz8yneyQr5LOMR16m4DpEHGOOH9E2CA7gHCbol9gPrBs0z8OUFWgm7AwjOjHo3/AK8S5G6uPLJHE6vCdu+C+wEpiQg83ygGyhLyLpRDvye4ODuYeBXMc9bFeY8TFDwb4pyGQN/DNfPg8Bz/RtNIA28mNV+HcGZXB8C9yck8xMExw0vE+zVFXS7kU/e0axTcR10iR8REYlU4o4RiYhIaVEhEhGRSKkQiYhIpFSIREQkUipEIiISKRUiERGJlAqRiIhE6n8CUE1rZKLM0QAAAABJRU5ErkJggg==\n", "text/plain": [ "
" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "import matplotlib.pyplot as plt\n", "\n", "# Create a figure with one subplot\n", "fig, ax = plt.subplots()\n", "\n", "# Plot polygons\n", "polys.plot(ax=ax, facecolor='gray')\n", "southern.plot(ax=ax, facecolor='red')\n", "\n", "# Plot points\n", "data.plot(ax=ax, color='blue', markersize=5)\n", "\n", "plt.tight_layout()" ] }, { "cell_type": "markdown", "metadata": { "deletable": true, "editable": true }, "source": [ "Okey, so we can see that, indeed, certain points are within the selected red Polygon.\n", "\n", "Let's find out which one of them are located within the Polygon. Hence, we are conducting a **Point in Polygon query**.\n", "\n", "- Let's first enable shapely.speedups which makes some of the spatial queries running faster." ] }, { "cell_type": "code", "execution_count": 65, "metadata": { "collapsed": true, "deletable": true, "editable": true }, "outputs": [], "source": [ "import shapely.speedups\n", "shapely.speedups.enable()" ] }, { "cell_type": "markdown", "metadata": { "deletable": true, "editable": true }, "source": [ "- Let's check which Points are within the ``southern`` Polygon. Notice, that here we check if the Points are ``within`` the **geometry**\n", " of the ``southern`` GeoDataFrame. \n", "- Hence, we use the ``.loc[0, 'geometry']`` to parse the actual Polygon geometry object from the GeoDataFrame." ] }, { "cell_type": "code", "execution_count": 66, "metadata": { "collapsed": true, "deletable": true, "editable": true }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "0 True\n", "1 True\n", "2 True\n", "3 False\n", "4 True\n", "5 False\n", "6 False\n", "7 False\n", "8 False\n", "9 False\n", "10 True\n", "11 False\n", "12 False\n", "13 False\n", "14 False\n", "15 False\n", "16 False\n", "17 False\n", "18 False\n", "19 False\n", "20 False\n", "21 False\n", "22 False\n", "23 False\n", "24 False\n", "25 False\n", "26 False\n", "27 False\n", "28 False\n", "29 False\n", "30 True\n", "31 True\n", "32 True\n", "33 True\n", "dtype: bool\n" ] } ], "source": [ "pip_mask = data.within(southern.loc[0, 'geometry'])\n", "print(pip_mask)" ] }, { "cell_type": "markdown", "metadata": { "deletable": true, "editable": true }, "source": [ "As we can see, we now have an array of boolean values for each row, where the result is ``True``\n", "if Point was inside the Polygon, and ``False`` if it was not.\n", "\n", "- We can now use this mask array to select the Points that are inside the Polygon. Selecting data with this kind of mask array (of boolean values) is\n", " easy by passing the array inside the ``loc`` indexing function of Pandas.\n" ] }, { "cell_type": "code", "execution_count": 67, "metadata": { "collapsed": true, "deletable": true, "editable": true }, "outputs": [ { "data": { "text/html": [ "
\n", "\n", "\n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", "
addressidaddrgeometry
0Ruoholahti, 14, Itämerenkatu, Ruoholahti, Läns...1000Itämerenkatu 14, 00101 Helsinki, FinlandPOINT (24.9155624 60.1632015)
1Kamppi, 1, Kampinkuja, Kamppi, Eteläinen suurp...1001Kampinkuja 1, 00100 Helsinki, FinlandPOINT (24.9316914 60.1690222)
2Bangkok9, 8, Kaivokatu, Keskusta, Kluuvi, Etel...1002Kaivokatu 8, 00101 Helsinki, FinlandPOINT (24.9416849 60.1699637)
4Hesburger, 9, Tyynenmerenkatu, Jätkäsaari, Län...1005Tyynenmerenkatu 9, 00220 Helsinki, FinlandPOINT (24.9216003 60.1566475)
10Rautatientori, Keskusta, Kluuvi, Eteläinen suu...1011Rautatientori 1, 00100 Helsinki, FinlandPOINT (24.9440942536239 60.17130125)
30Kampin keskus, 1, Urho Kekkosen katu, Kamppi, ...1031Urho Kekkosen katu 1, 00100 Helsinki, FinlandPOINT (24.9331155798105 60.1690911)
31Ruoholahdenkatu, Hietalahti, Kamppi, Eteläinen...1032Ruoholahdenkatu 17, 00101 Helsinki, FinlandPOINT (24.9252584 60.1648863)
323, Tyynenmerenkatu, Jätkäsaari, Länsisatama, E...1033Tyynenmerenkatu 3, 00220 Helsinki, FinlandPOINT (24.9212065 60.1587845)
334, Vilhonkatu, Keskusta, Kluuvi, Eteläinen suu...1034Vilhonkatu 4, 00101 Helsinki, FinlandPOINT (24.9473289 60.1718719)
\n", "
" ], "text/plain": [ " address id \\\n", "0 Ruoholahti, 14, Itämerenkatu, Ruoholahti, Läns... 1000 \n", "1 Kamppi, 1, Kampinkuja, Kamppi, Eteläinen suurp... 1001 \n", "2 Bangkok9, 8, Kaivokatu, Keskusta, Kluuvi, Etel... 1002 \n", "4 Hesburger, 9, Tyynenmerenkatu, Jätkäsaari, Län... 1005 \n", "10 Rautatientori, Keskusta, Kluuvi, Eteläinen suu... 1011 \n", "30 Kampin keskus, 1, Urho Kekkosen katu, Kamppi, ... 1031 \n", "31 Ruoholahdenkatu, Hietalahti, Kamppi, Eteläinen... 1032 \n", "32 3, Tyynenmerenkatu, Jätkäsaari, Länsisatama, E... 1033 \n", "33 4, Vilhonkatu, Keskusta, Kluuvi, Eteläinen suu... 1034 \n", "\n", " addr \\\n", "0 Itämerenkatu 14, 00101 Helsinki, Finland \n", "1 Kampinkuja 1, 00100 Helsinki, Finland \n", "2 Kaivokatu 8, 00101 Helsinki, Finland \n", "4 Tyynenmerenkatu 9, 00220 Helsinki, Finland \n", "10 Rautatientori 1, 00100 Helsinki, Finland \n", "30 Urho Kekkosen katu 1, 00100 Helsinki, Finland \n", "31 Ruoholahdenkatu 17, 00101 Helsinki, Finland \n", "32 Tyynenmerenkatu 3, 00220 Helsinki, Finland \n", "33 Vilhonkatu 4, 00101 Helsinki, Finland \n", "\n", " geometry \n", "0 POINT (24.9155624 60.1632015) \n", "1 POINT (24.9316914 60.1690222) \n", "2 POINT (24.9416849 60.1699637) \n", "4 POINT (24.9216003 60.1566475) \n", "10 POINT (24.9440942536239 60.17130125) \n", "30 POINT (24.9331155798105 60.1690911) \n", "31 POINT (24.9252584 60.1648863) \n", "32 POINT (24.9212065 60.1587845) \n", "33 POINT (24.9473289 60.1718719) " ] }, "execution_count": 67, "metadata": {}, "output_type": "execute_result" } ], "source": [ "pip_data = data.loc[pip_mask]\n", "pip_data" ] }, { "cell_type": "markdown", "metadata": { "deletable": true, "editable": true }, "source": [ "Let's finally confirm that our Point in Polygon query worked as it should by plotting the points that are within the southern district:" ] }, { "cell_type": "code", "execution_count": 68, "metadata": { "collapsed": true, "deletable": true, "editable": true }, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAaIAAAEYCAYAAAAeWvJ8AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4wLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvqOYd8AAAIABJREFUeJzs3WlwW9eVL/r/PhgOZmIkCIAEqXk2bYtWFEmeosiKbcWyY0me4tiOHcdxJ9X5ku7Uu111q96tvpXqfp2Ku7qdrjiO20nciYdEg4fYlh0pUaxZsiTbGi2JAAGSIAhOIDED+30gAZPiBJIHOAcH+1fFokgRwCIJYp2999prE0opGIZhGEYsnNgBMAzDMNWNJSKGYRhGVCwRMQzDMKJiiYhhGIYRFUtEDMMwjKhYImIYhmFExRIRwzAMIyqWiBiGYRhRsUTEMAzDiEopdgATsdvttKmpSewwGIZhmDk4ceJEN6XUMd3XSTIRNTU14fjx42KHwTAMw8wBIcRXzNexqTmGYRhGVCwRMQzDMKJiiYhhGIYRFUtEDMMwjKhYImIYhmFExRIRwzAMI6qiEhEhxEwIeYMQcp4Qco4Q8mVCiJUQspcQcmnkvWWK25sIIUFCyH8IFzrDMAwjB8WOiJ4D8C6ldCmAZgDnAPwYwIeU0kUAPhz5eDL/B8Bf5hIowzAMI0/TJiJCiAnALQBeBABKaYpS2gdgK4CXR77sZQD3TnL71QCcAN4XImCGYRhGXorprDAfQBjAS4SQZgAnAPw9ACeltAMAKKUdhJDaa29ICOEA/BuARwFsFCxqhpGwbDaL7u5uhEIhdHV1oaurC06nExs3sj8BhplIMYlICeBGAD+glB4hhDyHqafhRnsWwDuU0jZCyJRfSAh5GsDTAOD1eou8e4YRD6UU/f39YxJOKBRCJBJBLpcrfJ3T6cR9990nYqQMI23FJKIAgACl9MjIx29gOBGFCCGukdGQC0DXBLf9MoCbCSHPAjAAUBNCBiml4xIZpfQXAH4BAC0tLXQW3wvDlEw8Hi8kmtHvU6nUtLcNhULw+XxYunRpGSJlmMozbSKilHYSQtoIIUsopRcwPMV2duTtMQA/GXm/e4LbPpL/NyHkcQAtEyUhhpGq9vZ2vPXWW+jo6JjT/ezbtw+xWAzNzc1QKBQCRccw8lBs9+0fAHiFEKIGcAXAExgudHiNEPIkAD+A7QBACGkB8Ayl9KkSxMswZZFOp7F//34cOnQIlM59gN7V1YU333wTx48fx/r167F8+XJMN13NMNWCCPFHJrSWlhbKjoFgxOLz+bBnzx709PSU7DHuu+8+XHfddSW7f4aRAkLICUppy3RfJ8nziBhGDMlkEh988EFZzsLq6+sr+WMwTKVgiYhhAHz++ed48803MTAwUPLHUqvVWLBgQckfh2EqBUtETFWLx+N47733cPr06bI9JiEEra2tqKurY4ULDAOWiJgqdvbsWbzzzjsYGhoq6+PmpwBPnz6NLVu2sH1zTNVjiYipOoODg3jnnXdw7tw5UeMIh8N46aWXcMMNN2DTpk3QarWixsMwYmGJiKkalFKcOXMG7777LhKJhNjhFHz88ce4cOECNm/ejFWrVrGybqbqsETEVIX+/n689dZb+Pzzz8UOZUKxWAw7d+7EqVOncPfdd8Nms4kdkiR1dnbi008/hdfrxcKFC8FxE/dtppTiypUr+Pjjj0EIwerVq9HY2MiSvESxfUSMrFFKcfz4cXzwwQdFteORAoVCgfXr12P9+vVQq9Vih1Ny6XQaFy5cQCgUwvr166FSqRAKheByuQqJg1IKn8+H3/zmN4U+fmazGatXr8aNN94InU4Hn8+HaDSKZDKJTz75BD6fr/AYNpsN3/72t6HT6UT5HqtVsfuIWCJiZKunpwd79uwZ84JUSQwGAzZu3Ijm5mbZXclTSuH3+3H69GmcPXsWyWQSwPD3nM1mEY/HsWDBAnz1q19FMBjEsWPHEAqFJrwvhUIBu90+6f/n6fV63HLLLbjppptk9/OUKpaImKqVy+Vw+PBh7Nu3D5lMRuxw5qyurg6bN29GU1OT2KEIhlKKY8eO4YMPPkA6nS7LY2o0GrS0tMDhcECr1UKhUECpVEKhUIx5G/05pVIJpZKtYMwWS0TMtLLZLCilsvpDC4fD2L17N4LBoNihCG758uXYsmWLrKrrent7sWfPHrS2toodyoQMBgO+//3vg+d5sUOpSKzFDzMGpRTRaBTd3d2FN5/Ph3A4DLvdjm3btqG2dtzZhhUjm80ilUrh0KFDskxCwPC+p/b2dmzfvh1ut1vscARhsVjwrW99C8eOHcPevXslN4IdHBxEMBjE/PnzxQ5F1tiIqELlcjlkMplxi9mZTAY9PT1jEk7+baopEI1Gg02bNmH+/PmoqampqDn0zs5O7N69Gy6XCwqFoiy94sSkUCiwefNmtLS0VNTvaTqRSAS7du1CIBAQO5QxFixYgJqaGmSzWQDDRRK5XA7ZbBYrV66EyWQCpRS5XA6U0nH/tlgsVdtBg03NyVg6ncYf/vAHXLhwATzPw2g0Qq/XY2BgAH19fXM+tsBisUCr1SKZTMLtdiObzaKpqQl6vR4LFy6UTCVXJpPBgQMH8Le//a1QSaVQKAovGHK3atUqbNmyRTK/DyHkcjkcPHgQ+/btG3PKbSUyGo246667qvpARJaIZCoej+P3v/89/H6/KI9vMBiwatUqbNq0SdSr8WAwiN27dyMcDosWgxTY7Xbs2LEDDodD7FAEFQqFsHPnzmkr4aSIEII1a9bg9ttvr/q1JZaIZCgajeK3v/0turomOpW9vNxuN2666Sbo9XrYbDZYLJayJCahD6yTA5VKha9//etYtWqV2KEIKpvN4i9/+Qv+9re/Vczv2uVyYcuWLbJZw5srlohkJhKJ4De/+Q36+/vFDmVCRqMRTU1NaGxsRGNjI2w2m+CJ6erVq3jrrbdKemBdJWtpacHmzZtlVQUJAIFAALt27UIkEhE7lEmp1Wp85StfwU033TRpt4dqxBKRjGQyGbzwwguSGAkVS6/XF5JSU1MTHA7HrBPT0NAQ3n//fZw5c0bgKOXH5XJh+/btsFgsYociqFQqhf/6r/9Cb2+v2KGM43a78cADD8BkMokdiuSwRCQje/fuxcGDB8UOY040Gg28Xi8aGxvh9XoLFW5ToZTi1KlT2Lt3L+LxeJkirXw8z2PLli1YuXKl2KEIIhaL4cSJEzh8+DBisZjY4Yyj0WgKWx/UajVuu+02ZLNZ5HI55HI5GAyGit4aMRcsEcmE3+/HSy+9JHYYglOpVGhoaCgkJ4/HA5VKVfj/cDiMt99+u2Lb80jBihUrcOedd0Kv14sdyqx0d3fj8OHDOH36dGF/UU1NDWpqagol0plMBqlUColEQrIXK8uWLcOOHTvEDkMUbEOrDKRSKezevVvsMEoinU7jypUruHLlCoDhsmu3243GxkYoFAocOHCg4st3xfbZZ5/hypUr2Lx5M6677rqK2HM0MDCAK1eu4PLly4hGo8hms7DZbFCpVFAoFGhvb59ynTQ/OgkEAoI/f5YuXQqO43D27NkZ3W7JkiWCxiFHLBFJ2Icfflg1C/PZbBZtbW1oa2uDwWBgSUgg8Xgcu3btwqeffoq7774bZrNZ7JDGSCaTaG1tLVyUdHd3z+n+EolEYVpMaJs2bYLVasWxY8dw4MABRKPRaW+j0Whkt15XCiwRSVRrayuOHj0qdhhlp1AoMDg4KHYYstPf348//vGPcLvduPnmm0Wbrstms2hvb8fly5dx5coVBINBwZNGKdaRVq9eDaPRCGB4yvPEiRNTJiKv14slS5bguuuug8FgEDweuWGJSIKGhoZkOyU3HbPZLOky3UrkdrvR3t4OAGhra8PJkyexdu1arFu3DhqNpqSPTSlFJBIpjHiuXr1a0nOhlEolDAYDjEYjEokEent7i+7u7fF40NDQMK4bt0ajwerVqwEMl5L/+c9/hlarxcKFCwuduke/X7BgARYvXlyy71GOWLGCxITDYbzyyiuS3S9UavX19ZLrNVbpHA4HBgYGCmf+5Gk0Gqxbtw6LFy8uLPxns9lxb/nXiHzvtIn+zXFc4WtzuVwhCVy5cgUDAwPl/YavUVNTA4PBAIVCAUIIuru7oVarYTKZCj+TfKGMwWDAmjVr2GmuAmFVcxUoEAjglVdeQSKREDsU0TQ2NrJKOQwn5PyLff7FPZfLgRAChUIBtVpdSAR9fX0TvtgrFAp4PJ6ytIOy2+1zXt8pJ7VaPeXIzOVyYe3atVixYkXVNiwVAquaqzCXL1/Gq6++WrZDwqRKascAiIUQUvTIkBACvV4PlUoFlUoFpVIJjuPQ399ftp6ElXZG0nTTgx0dHdi5cyc++OADrF+/HqtWrWLHjJeQLBNRJpNBa2sr6urqKmKh8NNPP8XOnTtZpRgg2b0g5dbT0wOlUllUYqaUYmhoqAxRjVdTUwOVSoW2tjZRHr/UotEozp49i71792LVqlVYs2YNXC6X2GHJjiwTkUKhwHvvvYfu7m4YjUbU1dWhrq4OLpcLLpdLEuft5Bdxz5w5gwMHDogai5Tk18YsFgsMBgMIIRgaGoJWqwUhBOl0Gp2dnSJHWXpDQ0Oor69HMBiUdMNPk8k0Lgk5HA7wPI9EIoF0Ol04uyeXy0GpVMJoNEKpVFZE8uI4DpFIBNlsFqdOncKpU6fQ0NCAL33pS1i6dCmbthOILBMRIQQrV67E/v37EY1GEY1GcenSpcL/a7XaQlLKv5Wje3Qmk4HP58PFixdx6dIlSfbNEhPHcbDb7ejp6UFvb++Ynw8hBBqNBvF4HDqdDna7HYODg7LeZ1UJ39tEDT5VKtWU04r58nyv1yvacSbFqq+vHxdjfr+b0WhES0sLVq9eXbHdK6RClokIGD40bP/+/RP+XzweH7OrHxjuz3XtyMlut8+5k24sFsPFixdx8eJFXL58uaSlq5Uul8tNev4MpRQOhwN+vx+xWKzw4mC1WmE0GhGJRGS3/8hut0v+hfraSjwARXf/Hn1bhUIBk8kEnU4HhUIhie97ugrOaDSKffv24a9//SsWLFiAefPmYf78+XNq8FutZF0194tf/AIdHR2zvr1SqYTT6RyTnGpra6f9Q4tEIrhw4QIuXLiAtrY2SU+tVBKj0TjlJkKn0wme59HZ2VnRCV+r1RaSrtR5PB5Eo9FCx4Z4PF70YYU8z8Nms2FwcHBc1R/P8+B5vjCNnk6n0d/fX7amp3MZrRkMhkJSmjdvHmpqagSOrnKw8m0AH330ET744AMBIvoCx3Gora0dM61XW1uLUChUSD6VVMZaaUZvzpwMx3GFBeWOjo6KKQKpq6uDQqGoqJjNZjNSqVRZEoTD4YBWq0UwGCzpcfBCTxnabDbMnz8f8+fPR1NTU8k3EUsJS0QYXvj+2c9+JkBEUyOEwOPxIJlMVv3R1aXmcrlmNMrleR61tbXIZDIIhUKSfIE3Go1QqVQVsSY0kZqaGmSz2bJNjdbW1s7qbK7a2lrwPA9KKdrb28c9F/J/x6XcUE0IgdvtLoyYGhoaZHeQ4WgsEY146aWXyjrFYbPZoNPp0N7eXtKrtmpmMplmtVtfrVbD6XQim80iFApJ4vej0+mgVCpF7z4wV0ajEYSQkn8fSqUStbW1U46KjUYjMpkMtFottFotOI4Dx3EIhUKFzeJmsxmZTKaQPFUqFaxW66RrlKWiVCpx44034o477pBlBR7b0DpixYoVZU1EkUgEkUgEGo0G9fX1slxEF5vZbJ7VC14qlSqUDCuVykL3glAoJMpG2vw6iByqJ3U6XVn2gGUyGbS3t8Pj8aCrq2vcBvDR02pTxdPX1wetVova2lrE43Go1eqyJyFg+Ps5evQoQqEQduzYUbWbZosqCSOEmAkhbxBCzhNCzhFCvkwIsRJC9hJCLo28H9frnBDSSAg5QQg5RQj5jBDyjPDfwtSWL19e7ocEMNyO3ufzYXBwEB6PB3V1daLEIUehUGjO0xmZTAaBQADBYFCUDYpKpRImk0kWSQhA2Ud1wWCwcBz96BfvmUy9xuNxJBIJJBIJ0Rvt+nw+vPDCC6IkQykotjb5OQDvUkqXAmgGcA7AjwF8SCldBODDkY+v1QFgHaX0egBfAvBjQoh77mEXz2AwoKmpqZwPOU4wGERnZycsFkuhu+9o+d5hOp2u0JyRmVwymRQ0eXR3d5e1A4fL5YLZbJbNeqLBYBBlk3FfXx98Ph9SqVThpN+Z9mnU6XSSaavV19eHF198EefOnRM7lLKb9rKSEGICcAuAxwGAUpoCkCKEbAVw28iXvQxgP4B/HH3bka/N41F84hPUypUr0draKsZDj5HfpKlUKsHz/Lhux6NLjnmeh16vB8/zha7BoxFCQCkd874Syn2FImSVVjweh9VqhUqlKumLUr6gZS5bCqTGZrMhmUyKut6WyWRm/dxXq9UCRzM36XQar732Gm677TbccsstVbMfqZj5jfkAwgBeIoQ0AzgB4O8BOCmlHQBAKe0ghNROdGNCSAOAtwEsBPAjSumEq4yEkKcBPA0Mz/MKadmyZXj77bcls58nk8lMuyaRTCYn3Cw4Gbe7rANN0UUiEUE7Pvf09MDj8SAYDApyf6PV19djaGioJPctJqvVimg0WrF7tgghkh2V7t+/H11dXdi6davkkmUpFDNCUQK4EcDPKaU3ABjCxNNwE6KUtlFKr8NwInqMEOKc5Ot+QSltoZS2OByOYu++KDqdDvPnzxf0PqWmGqfzhH4BDAaDaGxsFPQ+GxoaEAgEZLMWlJc/eK5SkxAw3LlCyk12z549i1/96lfo6+ub831RSnHy5En86U9/wuXLl8eNYCml6O/vR2trKz7++GN89tlnc37MmShmRBQAEKCUHhn5+A0MJ6IQIcQ1MhpyAZiysJ9S2k4I+QzAzSP3UVYrVqzA5cuXy/2wZTPXVkSVxmAwlGRx3OfzzXiv0mgmkwkmkwmEkMK0qdVqHbNHSKFQSKJ0fLZ4ngfHcVN2uagElXB0RSgUwgsvvIAdO3bM6SKpra0NZ86cgc/nw9GjR8HzPObNmweFQoHu7m50d3cjm82C4zg4HA64XC4sWrSobKOxaRMRpbSTENJGCFlCKb0AYCOAsyNvjwH4ycj7cWdbE0LqAUQopfGRqrr1AH4q5DdQrPz0XCW/AExFihs1S8lsNpesLL6YdSKXy1XozGwwGKBWq0EIQSwWG7ch0uPxwOv1IhwOIx6PI5vNwmw2F/bexONx9PT0VMxz02azTdvdohJIZap+OrFYDL/+9a9xzz33oLm5eUa3zWaz2L9/Pw4ePDjmNSKZTOL8+fNjvra+vh7btm0TpSVRsTWwPwDwCiFEDeAKgCcwPK33GiHkSQB+ANsBgBDSAuAZSulTAJYB+DdCCAVAAPx/lNJPBP4ephSNRnH+/HmcPXu2Yv7QZ6PaDpQr5SJud3f3pG1eTCYT9Hr9mBHTdJVaE60N9fX1jZly4TiuqPZFYsp3PZdLUYxYZzjNRi6Xm/bCK5vNIhKJIBQKFd4CgUDRlYSBQAA///nPcdttt2HNmjVlnWUpKhFRSk8BmGh37MYJvvY4gKdG/r0XwHVzCXA2YrEYTp8+jXPnzlXEmSdCqOS5+tko97QQIQQNDQ3o6OgoyZRgLpeT9DpfPjHPJgmZzWYkk0mYzWZwHAdKKWKxGAwGA5LJ5ISjQY1GM+NS7JkSYu2lXFQqFW688cYxn4vFYjh16lQh6YTD4TnPjCSTSfj9fni93rIWQMmys4JOp8N1110Hl8uFUCiErq4udHV1IRwOz6gSrZJIedFVaFqttuQvIn6/v7BW5Ha7xxw9UQomk0nSoyG/34+GhoaiLuy8Xi8opeA4DrlcrnCba5+j+d9hbe1wwW00GoVGo8HAwAAIIaipqSkclCi0Ut53KaxcuXLcmpZarcbg4CDOnDkz5/vnOA6rVq3Chg0bYLfb53x/MyXLRAQAer0eer1+zGZWSileffVVXLhwQbzASqRc7fGlwGq1lqUUurOzE3q9viwJYrZti8opvx421RSRw+FAKBSa0QXf6AamyWQSuVwO8XgcJpMJsVisJHu7jEZjRSWiK1euYNeuXdi0aRMSiQSOHTuGM2fOCHYBajQasXjxYthsNkHub6Zkm4gmQgjBggULZJeIeJ6X7UhvIuWawqKUlmUdQaPRVMQeo1gsBp7nJ/y//BENgUBgTtNDo28bCoVKdoprJa0X59cOT58+jYsXL4LjOMGfl/39/di9ezeWLVsm6P0Wq2oSUTabxUcffYS//vWvYociuGQyKfmFbiHJrYmsw+GoiLVMnufhdH6xDTDf0WNgYKBkG0NDoRDq6+sRDAYFq3KrqamRfE83juNgt9uh1WrHjNxKNQWvVqths9lE6+RQFYmovb0de/bskfyTb7acTqes2sZMRaPRVOy5PZMpdXWSUqmE1WqFRqNBf38/+vv7C+sNDocDfX19RU0LJpNJ9Pb2IpPJlG1NMplMIhAIwGazIZVKwWKxjDn1eKZ7svLHbkhtu4PT6UR/f3/hqPRIJDKrM5dmi1KKmpoaHDt2DMuXL4dery/bYwMyP48onU5j3759OHz4cMXsGZgNnueh0+mgVquhUqnAcRwIIcjlcsjlckin00ilUkgmkxVf1CC3kZ9OpwPHcXMe5XEch7q6ujFdySmlGBwcnLarg06nK2qN0Wg0wmg0iv7z12g0henoRCIBt9s9rhs7pbTwNx8KhZBOp2G1WjE0NCT6NDbHcbDZbNBqtchms+jv75fUKF+j0eCHP/zhpNOwM8HOI8LwSOjQoUNih1FyM+1Lp9FoUFtbW5H7QeR2mqVOp5tVvzyz2Qy9Xl9YLwuHw7NOEPmigOlwHCd6EgJQOLohb7qY8mdPiXUKrl6vh9lshkKhQDweRyQSEb3HHc/zhZFhNptFY2MjnE4n7HY73G63IEloJuT1V30NOY+C5qLU+zNKLd86Rw66u7tRX18/5fHUHMfBarVCp9ONGeUIVcLe2dkJj8eDnp6eSUfMc2l7JLb82VPAcFeIUp89ZLPZoNfrkcvl0NfXh8HBQUltnlWr1fi7v/s7GI1GsUMpkHUiYiYnlTNYZsrv90On00Gr1Yp+mJlQAoEAvF4vMpkMVCoVlEpl4feTTCbR1dUlWJfxyQSDQTidTqhUqgnXi/r7+1FXVyfKuUNCikQiJd1DpFQq0dPTI+nn5saNGyWVhACZJyK5XDWXQiWPipxOJ65evSp2GILy+/3QaDRwOp3w+/1Ip9NQKpVwOBxlex7ni3nypdjxeLwwhWSz2Sqism86PM+XfD3GYrFIuqBmwYIFYocwjqwTETO5SuzWnW+zU6mjuclYrdbCxlmfz1f4fCaTQUdHR9m7AIxevzCZTLDZbJI4WFII+URfKgqFQvJdyfv7+0XbuDoZloiqVKWd/KjX66HVagsvIqXa6FhOHo8HmUwGoVBoyivo/LELYpQcDwwMwGg0ymZ2odQdSKS8p2/hwoVYuHAh6uvrxQ5lnMq7LJ4Bqe0VkJJKSkROpxO5XG7MOkm+MWM5abVaeDwe1NXVzWlEmT/zJRgMFrW3raurS9QXj2AwWPafdaloNJqSP0Y4HB6z8Vcqurq6YLfbJXniq6wTkVyu4kqhUhKR1+tFV1fXhNVc+Uac5ZLvcdfZ2Vk4tsHr9c544dfj8cy4fLetrU3Ug9xKXSxRDkqlsixFBOl0GqFQCAqFAlqtFl6vVxKH8A0MDOC3v/0t/vKXv0z6NyUWWW9ojUaj+PDDD/Hpp59WVG+pcvB6vchms0ilUkgkEhgaGpLUCNLpdCKbzU77AqhUKqHVaks6L08IKXQgmOy4DbPZDJPJhFQqha6urkl/lnOZUmxsbByzhlRulT4dKmb8DocDhJCydkuYjsFgwPe///2S7hkqdkOrrBNR3tDQEE6cOIHjx49LfiFRTPmyaGB4r4FSqSycH5PvKzba6N3r135MKYVKpZrxH77VagXP8zPas+J0OgVv32SxWGA0GpFOpxEOh2EymYquhFIqlTAajeB5HgMDAzAYDOB5Hl1dXXPa1c9xnKhdo6fb7yRlGo0GuVxO1HO7pLgXa926ddi0aVPJ7p8loglks1mcP38eR48eregru0ri8XiK6ixdU1MDo9E46xe6uV7t6nQ6WK1WAEBPT8+4RW2xRyN5pUi6xWhoaEAoFJLcAYwqlQp1dXXjLpLy7a0ymQzS6TRqampEr/xTq9VwuVySeB7lcRyHH/3oRyVbO2MtfiagUCiwYsUKrFixAp2dnThy5Ag++eQTNm1XQvmF7smSRP746ba2tjld6QeDwRmVOef36KhUKgwODk6YfEaTykg6FAoVfUCdUPK/Hylyu91FvbBHo9HCeqNYe+jMZrPkEnl+lFiOIo6pVNWIaCLRaBQ//elPy/JY1ay+vr7QfBIYrkCrra1FMBhEJpMR5DFqa2unnIO32+3Q6/VIJBKzOlZZKqd6ajQapFKpsq3pSXFKbrYNWFUqFdxuN7q7u0vadsdqtRaq03ieRyKRELX7PyGksF8t/6bT6WAwGHDdddeVbJ2IjYiKZDAYxA6hKgQCAZjNZvA8D57nx23eFEJXV9eYKTSTyQSz2Vwo/c6/zZbZbJZEIsp3nC7HXhWFQiGpztC1tbVQq9Xo6OiY1Sg1nU4Xnh9arRY1NTXgeR65XA6xWAx9fX1zmiFxOBwYGhqSVGeFpqYm3HnnnYUj2aWo6hMRIUS0zYLVxGg0wmQyIRgMlnQqtK2tDQ0NDYUzdoQ8fltKx7GXqwt5/mgRManVatTV1WFgYEDQqrN4PD5hCbPZbIbBYIBCoUAqlUJ/f39Rv3uHw1Eo+BkYGBCsKe1cEUJgsVjEDmNKVZ+IgOGrPpaIhGe1WqHVasFxHAKBQFnWWXK5XMnWM8LhMAwGgyRGCOVaa8i/WOerCCmliEajJX2RtVgs6O/vLyTBWCxW1uKivr6+cd9ffvSkVqsLo6fe3t5ClWi+v1x+f5iUNrQGg0GcPXsWzc3NYocyKZaIMJyI5Na/TCwKhQIejwcDAwOSmp4QiliJSKlUora2FhzHIRKJlL0Ldm9v75ifFEijAAAgAElEQVQD9pxOJwghJYlDq9UikUiAECKZIpGJRk9arRYmkwlqtRoDAwNjRvr5I86z2Sz6+vpE2TxqtVqxZs0aXH/99WU/X2imWCICCoeLMbPndrsBDO/ZYqXxwrHb7QCGjy+QUv+y/MJ7KQoZOjs7CzMUarUatbW1GBwclMxUV95kU3t5+Z+L2WwGpbRs1XqLFi3CmjVrsGDBgorpoMISEVgimovGxkZkMpmi9grJQbn7dHV3d8Pr9aK3t1eS2wxKUXk2epo8lUohEAiAEFKxnR36+vpgs9lgsViQyWSgVCqhVquRTCbR09MjyDQrz/O4/vrrcdNNN0mus3YxqjoRZTIZHDt2DDzPw+l0guM4cBxXuIpQKBRob29n03aT8Hg8ktqcVw5iTBX5/X643W50dHRIrn+i0WgcM2VXKtd28ag0U/W4yxdHcByHRCKBnp6eabc0qNVq2O12aLVaLFmyBM3NzZJsZlqsqkxElFJ88skn+POf/zxtOa5Wq4XL5UJbW1tF/yEIgRBSKJ+Nx+NVMwoarbe3F0ajsewJqb29vewbWadit9vB83zZRig8z8ui8epEJiqOsFgs4Hm+8LfW29sLi8Uy5sDC9vZ2rFixAjfddJNIkQun6hJRIpHAK6+8UvS8djweh9/vh9lshk6nk9Q8fTkRQmCz2SbdlGcymaDT6cb84Qi1UVVqotEorFbrmGIMtVpd8kq2trY2ybQa6uvrK+vvt9QH2knNRKPMiTq2NzY2liOckquqRBSLxfDKK6/MKpnkr1qK7QotN5RSaLXawhVZXn7fhE6nG3O1brVakclkEI1GYbfbEY1GUVtbi1wuVzgGu7u7e05NQMVkMBjQ09MDl8uFbDYLlUpVlhGiz+cTdDOrRqOBXq8Hz/OFvUmUUnAch2Qyib6+vnGL7OXutGC326sqCc2EFI/9no2qSUTJZBIvv/zynDfEja4W6u/vl0x5aTm0tbVBp9Ohqamp0K0gFoshFouNmwPv6ekBx3FQKBQIh8NQq9XjXkzMZvOYRKTT6Qr7MaR0VspEYrEY6urq0NHRUfZF9HA4DIvFMuu1mZqaGtTU1KC9vR2JRGLaaq78ZuR8N/arV6/O6nEZYZnN5kKj3kpXFYkofyDUTA8jm0ogEADHcfB6vQiFQhV7ZV+sfDFHb29v0V2M8x2QgYk3YPb19aGhoQGEEFBK0dbWhlgsVnTHbjGNHhGXe+0wnU4jk8nA4/EgkUggmUxCq9WOeX4TQqDRaMBxHHieh06ng0KhQDweR1dX14xaFUWj0cIFV01NjeDfz2SUSiUaGxtx+fLlsj1mJZk3b57YIQhG1okol8uho6MDb7zxRkn2IORyOfj9fnAcB5PJBI1GA5VKVSgHp5QilUqhp6cHhBDU1NRgaGgI2WwWBoMBarUa0WhUEjv1r5XfSZ/JZNDT01Oyho3XTucZjcaKmvbUaDSixKvX68cka5VKBavVilQqhWQyiXQ6XRhVCt37zOv1or+/v6R99/Kboi9fvgy9Xl/SBqWV6uLFi+jq6pJ0D7liybb7NqUUP/vZzwTtNVYK+eQklYV9nucLRyOUU11dXdm7BQhBjL0t+cMLxZy+nMs6kVKpRDabnXAkabFYoFarx1z4VOr+oXLQ6XR49NFHUVdXJ3YoEyq2+zZXjmDEQAhBfX292GFMq7+/v9CVQArq6urKnoQcDkdFJiGxmM1m0dfQUqnUrPat2Gy2QnVlfX09GhoaoNPpwPN8YeOumMclVJpYLIaXX3654qt5ZTsiAoZHRd3d3bhw4QIuXryIQCAg2b1AXq8XwBfz/+l0WpSWJvkzewYHB6fchDeZ/BTRTM77kdJ5NxqNBoSQwu9gOuXe2yOl0cFMR7FerxeBQGDC54VSqZx0VsDhcKC/v19yh8pJCc/z+OY3vym5i292VPgEYrEYLl26hIsXL+Lzzz+X9BNbpVIVDo4Ti06ng81mQy6XQzgcRjabhclkKkwNDQwMFEZP+e4UHR0dSKfT0Gg0qK2tnfZF0263IxKJCHKBoNVqYbFYCmt0+Q4ZhJBC4UT++Oj8i2F+HY9SOm4dQqFQwGAwQKvVQqVSFe4nv08ql8uVfUrR4/GAUiqJK+Bik7BarYbD4ZjTcznfcy5/ZEu+3VH+GJdwOCz7gqHpqNVqPPLII4WLWikQNBERQswAfglgJQAK4NsALgB4FUATgFYAOyilvdfc7noAPwdgApAF8M+U0lene7xynNCazWbh8/lw8eJF6HS6wh94LBbDkSNHJPGHDgxfRUq5q4PJZILRaMTg4OCYxWutVgue5wujuvyZMsDwi382mx2TyACA4zhYLBbodDpwHFdYEFcoFIU3pVIJpVIJhUJRSDS5XA6JRKKsU4ocx8FqtcJkMuHKlStle9w8l8uFXC4n2jSWSqUqnLszFZvNhnQ6XfK12rmUs8uJSqXCQw89JJmKOqET0csADlBKf0kIUQPQAfh/APRQSn9CCPkxAAul9B+vud1iAJRSeokQ4gZwAsAySumUc07lPCp8IpRS/Mu//ItoZ9tfS6p9xqZitVphMBiQSCSgVCrR398/bsRhMBhgNptBCJlwL1IluLbDQjnxPI9MJiNaM1SO49DQ0DDpce9TTcUJzWazFQ6vE3v9TGxKpRIPPvigJDa7ClasQAgxAbgFwIsAQClNjSSSrQBeHvmylwHce+1tKaUXKaWXRv7dDqALgKPYb0IshBBJdeRub2+X3NzvdHp6egqJnFIKq9UKr9cLr9cLl8sFnU6HwcFBBAKBwv6hSvsegeHpS7Ekk0m4XC7RHj+Xy8Hn80Gr1cLj8RQ+ny9E8Pv9ZTtwMhKJIB6Pi/r7KDeDwTDhSb2ZTAa/+93vcPHiRRGimp1i9hHNBxAG8BIhpBnDo5q/B+CklHYAAKW0gxAyZTE7IWQNADWACXenEUKeBvA0AEnMcUrtxNb8Edjd3d0Vc8U3uotFviJqsk4U+dY/lcRisYheZCGFC6b8hte6ujooFApEo1HRfi6RSERSBR2lMnp9Tq/XF/YlEkKQzWYRj8dx4MABLF68WORIi1NMIlICuBHADyilRwghzwH48UwehBDiAvAbAI9RSid8haeU/gLAL4DhqbmZ3P9sUEqRyWSQSqXGLFYnk0mkUilJvti3tbVBoVDA6/Wiq6tLMlOHxYjH41Cr1ZMucNfV1UmimedM6HQ60dclEokEVCqVJJJ4IpFAPB4XvWjA7/fDbrejp6dHcheUQri28e3Q0NCEG343btxYzrDmpJhEFAAQoJQeGfn4DQwnohAhxDUyGnJheNptnJGpvbcB/BOl9LAQQc9VNptFIpFANptFLpcr9EQ7evQoTp48KYk/6slks1n4/X6o1Wo0NjYiEAhI8sC00YxGIyKRCHK53KTrQFL/HiZSys4CxQqFQlCr1aKPAnQ6HbLZrOhJKC9/oGA0GoVWq4VCoUAoFJJ0pex08nsji7lgczgcWLZsWRmiEsa0iYhS2kkIaSOELKGUXgCwEcDZkbfHAPxk5P3ua287UtiwE8CvKaWvCxr5HCgUCuj1+nGf/9rXvoZbb70VJ06cwNGjRyXd0DSVSsHn88FgMMBisaCjo6Os3Rnypd3pdLpQ2m21WqFWq9HZ2QmXy4V4PA6lUgmdTldoPTOZ/Iipkl4opHIQWSqVQldXF+ocDph8PjR0dsIaDKJt6VIcdpR+STZfQSe1YpN8Ys6PWvMXb5MVV0iZSqWC3W4ves/aLbfcUjHHhAPF95r7AYBXRhLLFQBPYLjQ4TVCyJMA/AC2AwAhpAXAM5TSpwDswHChg40Q8vjIfT1OKT0l3LcgLK1Wiw0bNmDt2rX47LPPcOjQIUnv9B4cHMTg4CCUSiUaGhqQTCbn3GH8Wh6PZ8xaRL4l0eg/CpvNhp6eHlBKwfM8Ojo6ZvQYkUgELpdrxrcT00QLxeXUmEzCGY1CnUrB2dqK+X/7G3QjlWMA4Pz0Uxz+7ncBrnQNVAghsNvtFfF7y18MVVoSyh+1UuzP2G63Y/ny5SWOSlhVtaF1NiilaG1txaFDh3Dp0iWxw5mW0+ksKnESQuB0OgujkHg8XjgnyGQyFY5M7+3txcDAAOrr65FIJEre4LOSes4JeS7QTCwdGsItf/wjXEV0pd77ve/hoNNZslikdGpsMcSewpwps9mMbDY7o9mZb3zjG1i1alUJoypeseXbsu6+LQRCCObNm4d58+ahu7sbhw8fxunTpyV7VTXdcNzhcBSODJjsBT+/+TB/PAOAQhWU0+mEQqEo2b4mKVSBFYsr4UjjWs50GovCYahTKaz7zW+gKHJNbdHBgzh4330lianSXtSB4ZF3Y2MjwuFwYd+RVNXW1qK/v39G6242mw0rVqwoYVSlwRLRDNjtdmzZsgUtLS14//334fP5JFeVEw6H4fV6EYlEYDQaoVarQSlFMpnEwMDAjM5kmijR5EdbOp0OdrsdfX19gu6aZ4loLE0uhztOn8bKd96BaoZFND21tdh1zz2Cx6RSqeBwOCouCQFfVJjlF/7z62tS4/F40NHRMePXl5tvvrmsF0hCYYloFurq6vCtb30LyWQSPT098Pl8uHr1Kq5cuSLaSCk/0kmlUujs7EQqlSrpGS6xWAypVAo1NTUYHBwULCFHIhHU19cjFApJunqx5CjFLe3tWP366zDNsvmtKRKBjlIIWdtnMBjA87xkWmDNFqW0MMq3Wq3Q6/Xo7OyUxHNutiNNq9UqmSm5mWJrRALKZDJ48803cebMmbI9ZrHNRYVUjqkNq9WKWCwm6b1SpZqaWjY0hJv/8Ae4BOphl1KpkFGpkFGr8T/f/z5CsyyyqKmpkeweOyHk+yGW+tC/qczlObV161Zcf/31Akc0N2yNSARKpRJXr14t+eMYjUZYrVbE43GEw+GyJqFybTzt6emR/BqE0BdxjkwGX923D4s/+kjQ+1Wn01Cn00Ashof/8z/RtXgxuhcuRFahgOvcOWj6+kCyWSRqahBzOBC12TBgsSBiNKKD5xEdmeoxmUwVVZgwU6lUqvB8yzfoLVfhjEKhQF1d3ayf72azuWJHQwBLRILbtm0bdu7cWZKzhPR6PWw2GwYGBkTrQtDZ2Yna2lokEomSd1SW2vrbtbhkErZsFikAaY5DklLQWczPq3M5bPrsM6x6803wJd5HZerthenIESw8cmT6Lx6R0GgQtdux+/HHAZFL1suls7MTjY2NZXksnudRU1Mzp2Mybr755opaX71WdTyrysjr9eIb3/gGfvWrX836PjweD6LRKFKpFPR6PXieRzabRVdXV0nXfYqVX9x1u90TdtUWitSngJxnz+Lx//iPMZ/LEQLfDTfg5Fe+gnatFj0cB0xSyajJ5bC2rQ3Nu3fDXI4O3tcDeBbAXgCbADwPoIgdfZpEAppAANefP4/gypUlDVFKylFVl29cOpeCiZqaGjQ3NwsYVfmxRCSw9vZ2/M///M+MblNTUwOtVgu1Wo2BgYExV0ZSXiNpb2+H1+sttEwSUv44aYVCIdn2P9kJrkA5SjHv5EnMO3kSAPD6j36ETo0G+kwG+kwG2lQKulQK9efOoenQIWhK1RJnO4AfAvgZgNcxnIT+E4AVwIqR98BIm+HiNB08CFRBItJoNKipqSn5Rnar1YpEIjHn2ZNKHw0BLBEJilKKY8eOzajuP78OIoW+ZbPh9/sLpbB9fX2CHU4Xi8UK/cukKl3EH//2f/3XMkQygR9iONn8EMOJ6NmRj3swnJzyI6IZsLe344f//u8oZmUsrdPh+SefnHQ0KGVOpxOxWAw8z5esd15dXR0ikcicq/RMJpPkChRmgyUiARFCsHXrVtx66604duwYTp48WRgpqFQqWK1W2Gw2tLa2wmAwgOM4SS/GFytfCksIEXSnvUajEeR+SkUptSSZn3p7HsPJJj8iAr5IOvnpuFl2fqwpdgqxpweedBpBifTjm4nR669arRZGoxE8z4PjuBltz8h3JwG+2ByezWbBcRwCgYAgxS4bNmyo+NEQwBJRSZjNZmzatAm33norQqEQzGYzDAZD4UkZDAbx0ksvSfpqfzYopWhra4PD4UAmk5nzEQmDg4MwGo2SbT67ROwtBqMTz6mRf9888n9PY2yyOYUZTcMJYV53N4Jud3kfVGDxeFyya5VGoxE33HCD2GEIovK24FaQ/Pk7RqNxTOsdj8eD7du3V1R33JkIh8PQarVzvp++vj5ks1k0NjZCpVIJEJlwjNks5h09Km4Q+cTz7MjHzwM4gBlPuZVKrQxG+1K2YcMG0RvvCkUe30UFWrJkCf7pn/6p8PHg4CA+//xzXLx4EZcvX5ZsL7tihUIhQZqCxmIx+Hw+mEwm1NbWzqnEVUjrrlyZccsdweUTznkAH2F4Gq7Mo56pNBw9ivspRUajQZrnC29JlQopnscVoxG9MphWEoPBYMCNN94odhiCYZ0VJCidTqO1tRUXL17E+fPnBSsAEIPQRztoNJpCtVFPOUqeJ/HtV19Fw7lzoj3+GB/hi0KE9SLHMgN//u53ccDlEjuMirR582asXbtW7DCmVWxnBTY1J0EqlQqLFi3C3XffjR/+8IfYunUr7Ha72GHNSjgchkPAw9kSiQTa29uh0WhEXaRViD0aGu1n+KIaroJI8SK4Euj1eqxevVrsMATFpuYkTqFQ4Prrr0dzczPOnz+Pjz76SDLTU8XIZDKCrBddS6lUilrswUmp0OR1zLoKTlQyXSMttXXr1kluzXSuWCKqEIQQLFu2DEuXLkVrays++OCDiumAXIqiDDGn5W73+VBXhp6CckdZIpoxnU6HlpZpZ7oqDktEFSZ/UN/DDz+M559/XvKHewHD+zKEPslzcHAQJpMJhJAxiW6ij6/9vFarRVdX16y7QbjPnp1D5Awze2vXroW6AvdmTYclogql1+uxZcsWvPbaa2KHUpRAICD4MeBzabqqVqvh9XoRDAZnPMVnEanhrOywEdGMEEJks2/oWqxYoYItW7YM1113ndhhFIVSikgkArdENjjmW/5rtVrU19cXfTt1NguLBE/0ZOSvsbERBoNB7DBKgiWiCnfnnXfCaDSKHUZR0uk02tvb4XQ6Ba2km4vBwUEEAgHYbLbCGTRTWZhIgJP48RSVou7yZSjYz7Joy5YtEzuEkmFTcxVOo9Fg69at+O1vfyt2KEXLdzX2eDwYHByURMPXSCQCALDb7YWKpInWmryHDpU/OJla+f77cJ0+jUMPPIATVuv0N6hyLBExkrZgwQL8wz/8A/x+P06dOoXz58+LHVJR8mXoDQ0NiEQikii86O7unvL/17NqOUHZQiFsfv55nPhf/4utGU3B6/VWzMzHbLCpOZnQarVYtGgRnE6n2KHMWFtbG1KpFBobGyVfEWRpbRU7BNlRZTLQss2tU5LzaAhgiUhWBgYG8Mknn4gdxqxkMhn4fD7kcjm43W40NDSA53mxwxqD5HKwBgJihyFLFrZWNCWWiJiKYTab8e1vf7uiO/JmMhm0t7cXRkler1fskAAANdksth84ALWUWvvIiCmVEjsEyfJ4PKipqRE7jJJiiUhm9Ho91q1bJ3YYgqCUwu/3o7GxUbQYSDaLOy5dwnf/7d+wbN8+0eKQO0OpjkyXgeXLl4sdQslV7qUzM6kNGzbg5MmTFd21ezSfz1c4Ur2cVvf24suvvgqbgJtwmYnpYjHAYhE7DEmS+7QcwEZEsqRSqXDrrbeKHYag/H4/PB5PWR6rIZnEY7t2Yctzz7EkVCbaoSGxQ5Akl8sFSxUkaDYikqkbbrgBhw4dErU5qNA6OjpQW1uLrhJ1NtBns9h8/DiWvf8+lFLqrl0FNDIZvQutGqblAJaIZEuhUOArX/kK3njjDbFDEUwul0NfXx8sFgt6e3sFu1+Sy+F2vx83vvEG9OwFUXBZjkNCq0VKpxt+0+uRMhiQ0uuR1umQ1OnQJpGiFKmphmk5gCUiWVu+fLkoayullEqloFKpoNfrMSTAdE7zwADWvf46agXsDM58IVJXh//47nfZZtVZcDqdsNlsYodRFmyNSMYIIbjrrrtKch6QmIaGhqDRaOZ0OJgrlcKjb72Fe3/6U5aESujAww+zJDRL1TItB7ARkew5nU6sWbMGR44cETsUQUUiEbhcLnR2ds7oyGlNNouvnTqF5X/6E1SZTAkjZC7cfDNOm0xih1GxWCJiZOX222/HZ599Jpty7ryOjo7ipx5zOdwaDGL166/DOIdzjJjiDBkMeOf228UOo2I5HA7Y7XaxwygbNjVXBXiex+bNm8UOoyT8fv+03RdWDA7iu7/+NW578UWWhMrk6EMPYYBjLy+zVS1FCnlFjYgIIWYAvwSwEgAF8G0AFwC8CqAJQCuAHZTScaVMhJB3AawF8DdK6RZBomZmrJizdiqV3++H2+1Ge3v7mM/XptPY9OGHWHj4sEiRVae25cvxV4kcgFipqmlaDih+RPQcgHcppUsBNAM4B+DHAD6klC4C8OHIxxP5VwCPzjVQZm4q5WiI2erq6oJ15EwbdTaLr58+jSf/5V9YEiqztEqF9+69lxUozIHNZkNtba3YYZTVtCMiQogJwC0AHgcASmkKQIoQshXAbSNf9jKA/QD+8drbU0o/JITcdu3nmfJJpVI4c+aM2GGUVCaTQSqVwrreXnzp17+GScB9RkzxTt97L4ISP8pD6pYtWya7StfpFDMimg8gDOAlQsjHhJBfEkL0AJyU0g4AGHk/pxROCHmaEHKcEHI8HA7P5a6YUdLpNH73u9+hGn6mg4OD6FarwUvggL1q1O12490qW9sohWqblgOKS0RKADcC+Dml9AYAQ5h8Gm7WKKW/oJS2UEpbHA6H0HdfldLpNH7/+9+jtYoOc7uo1+O9Z55BVqEQO5SqQgHse/hhZFmBwpxYLBZZr+dOpphnTQBAgFKa34jyBoYTU4gQ4gKAkfelaQDGzEomk8Frr72GK1euiB1K2X1sseAvTz4pdhhV5dzGjThrMIgdRsVbunRp1U3LAUUkIkppJ4A2QsiSkU9tBHAWwB4Aj4187jEAu0sSITMrPp8Pn3/+udhhiOaA240jDz0kdhhVIVpTgz+tXy92GLJQbUUKecWOo38A4BVCyBkA1wP4vwB+AmATIeQSgE0jH4MQ0kII+WX+hoSQAwBeB7CREBIghMhzQ4vELFiwQLZ7h4r17pIlOHPnnWKHIXuHH34Yg2xKThByP4l1MkXtI6KUngLQMsF/bZzga48DeGrUxzfPOjpmTtauXYtYLIYDBw6IHYpodq5ZA21/PxYdPCh2KLLU2tyMg1V6FV8KZrNZ7BBEwS5jZO72229HS8tE1xBVghD8fuNG+FeuFDsS2UmpVHh/yxa2Z0hApirtzccSkcwRQnDnnXdixYoVYocimpxCgd/ddx86m5rEDkVWzmzdio45dEBnxjIajVBUabUnS0RVgOM43HfffViwYIHYoYgmoVDgtUcfRQ+bRhJExOnEe1W436WUqnVaDmCJqGooFArs2LEDHo9H7FBE06tQ4I2nnsJglU5/COmjhx5ChhUoCKpaCxUAloiqilqtxsMPP1xV7eWv1aFWY8/3voeERiN2KBXr8po1+LiKr95LhSUipmrodDo8+uijVf2kv6TV4t1nnkFGyY7jmqkkz+O9O+4QOwxZYlNzTFUxmUz45je/CZ1OJ3YoojltNmPfk0+i+LNdGQA4dd99CLMEXhLVfHHIElGVstvteOSRR6Cu4k7JB10uHHrkEbHDqBjdLhfeX7xY7DBki42ImKrkdrvx4IMPVm3JKADsXbQIp7/+dbHDqAgHHnoIOVagUDJsRMRUrXnz5uH++++vykaLebtuvBEXbmYNQKZyad06nGHVhiWj1WqrenaCJSIGy5Ytw9133y12GOIhBK/edhtam5vFjkSSEhoN3ts4rpsXI6BqnpYDWCJiRqxevRobq/jFhioU+P0996Bj/nyxQ5Gcj7dtQ6SKp2/LoZqn5QCWiJhR1q9fj7Vr14odhmiSCgVefeQRRKrwYLLJdDU0YC9LziXHEhHDjCCE4I477kBzFU9R9SsUeOPJJxGt8hcGYPjU1b888AAoK1AoOTY1xzCjEEJwzz33YHEVl+l2qlTY9cwzVdt9Ia1SIbh4MQ5/85vs1NUyqfYREduZxozDcRy2bduG3bt347PPPhM7HFFc0WrxzrPP4uv//u9QZTJih1NSUZMJ4UWLEF68GL76elzUapFlo6CyqvYREUtEzIRUKhW2bduGpUuX4u2330YikRA7pLL7xGSC4emn8dWf/xwclUcPhhwh6PZ4EFm8GJ0LFuCSwzF8lEMVl+9LARsRMcwUVq5cicbGRuzZsweff/652OGU3aHaWugefRQbfv1rsUOZlQTPI7xgAboXL0awsRHna2owxEY7kqJSqaDVasUOQ1QsETHTMhqNePjhh3HixAm8//77SKfTYodUVh/Onw/dvffixl27xA5lWn02G7oXLULXokVodbvxOc+zYgOJM5vNVb2hHGCJiCkSIQQtLS2YP38+du3ahba2NrFDKqs3m5uh6+vD0v37xQ6lIKtQINzQgO7Fi9E5fz4u2u2sIWkFqvZpOYAlImaGrFYrHn/8cRw8eBD79u1DLpcTO6TyIASv3XILHh0YwLyTJ0UJIabTIbxwIboXL0ab14vzBgOSbLRT8VgiYomImQWO47BhwwYsWrQIO3fuRCgUEjuksqAch9/dfTcej0bhvnSp5I/XU1uL7sWLEVq4EFddLlxVq1lRgQxVe8UcwBIRMwdOpxNPPfUU9u/fj4MHD4LKpLJsKmmFAr9/8EF868UXYW9vF+x+hwwGRJqa0NvUhM7GRlywWtHL2upUBTYiYomImSOlUomvfvWrWLJkCXbu3Ine3l6xQyq5qEKB1554At98/nmYZvH9Rk0m9DQ1oXfePITq63HVYkFIoWCjnSrFRkQsETECaWhowCWljDcAABApSURBVDPPPIP3338fJ06cEDuckguPdF/Y9txz0MVik35dn82G3sZG9DY2IuTx4KrFgjAb6TCjsBERS0SMgNRqNbZs2YIlS5Zgz549GBwcFDukkrrK8/jTs8/inueeA5fLocflQm9jI3q8XnTW1eGK0YgoKyZgpsBxHIxGo9hhiI4lIkZwixYtwrPPPou3335b9i2CPjUYEPvnf0YwFkNS7GCYimMymap+DxHAmp4yJaLVarFt2zbcf//90Mi8eeiVWAy1DQ1ih8FUILY+NIwlIqakVq5cie9973tYsGCB2KGUVFtbGxobG8UOg6kwbH1oGEtETMmZTCY88sgjuPvuu6FSqcQOp2R8Ph+8Xq/YYTAVxGQyiR2CJLBExJRFvkXQM888g/r6erHDKRm/3w+PxyN2GEyFYCOiYSwRMWVltVrxxBNPYOPGjeBkWlHW0dEBp9MpdhhMBWCJaJg8XwkYScu3CPrOd76D2tpascMRXC6XQ29vL6xWq9ihMBLHEtEwlogY0dTV1eE73/kO1q1bJ3YogkulUkgkEmyPCDMlloiGsUTEiEqpVGLTpk144oknZFfKGovFwHGc7MvXmdnRarVQq9VihyEJRSUiQoiZEPIGIeQ8IeQcIeTLhBArIWQvIeTSyHvLJLd9bORrLhFCHhM2fEYuvF4vnnjiCej1erFDEVR/fz8MBgOU7Jwg5hpsNPSFYkdEzwF4l1K6FEAzgHMAfgzgQ0rpIgAfjnw8BiHECuB/A/gSgDUA/vdkCYthTCYTduzYIbsihu7ubjgcDraDnhmDJaIvTPsXTwgxAbgFwIsAQClNUUr7AGwF8PLIl70M4N4Jbr4ZwF5KaQ+ltBfAXgBfEyJwRp68Xi/uuususcMQXEdHh6zL1pmZY3uIvlDMped8AGEALxFCPiaE/JIQogfgpJR2AMDI+4nKnzwARp8pHRj5HMNMavXq1Vi9erXYYQiOdV9gRmMjoi8Uk4iUAG4E8HNK6Q0AhjDBNNwkJpqLmPD0NELI04SQ44SQ4+FwuMi7Z+TqzjvvlGWXAtZ9gcljiegLxSSiAIAApfTIyMdvYDgxhQghLgAYed81yW1Hd4OsBzDhsZaU0l9QSlsopS0Oh6PY+BmZUigU2L59uyzLn/1+P5umY1giGmXaREQp7QTQRghZMvKpjQDOAtgDIF8F9xiA3RPc/D0AdxBCLCNFCneMfI5hpmUwGPDggw9CIcOD5ILBIOrq6sQOgxERS0RfKLY86QcAXiGEnAFwPYD/C+AnADYRQi4B2DTyMQghLYSQXwIApbQHwP8BcGzk7f8d+RzDFMXtduPrX/+62GEIjlKKSCQCm80mdiiMCDiOg8FgEDsMyShqcwOl9BSAlgn+a+MEX3scwFOjPv4VgF/NNkCGaW5uRkdHB44cOTL9F1eQdDqNWCwGk8mEgYEBscNhyshkMslum8JcsJ8EUxHuuOMOzJs3T+wwBBePxwEAOp1O5EiYcmKl22OxRMRUBI7jsG3bNtm1AQKAgYEB1u6lyrD1obFYImIqhk6nwwMPPCDLw/UikQisViubrqkSLBGNxZ71TEWpq6vD1q1bxQ6jJDo7O+F2u8UOgykDlojGYomIqTgrVqzAhg0bxA6jJAKBANvwWgVYIhqLJSKmIt1+++1YuHCh2GGUhN/vZ62AZI4lorFYImIqEsdxuP/++2V7CqrP50NDQ8P0X8hUJJaIxmKJiKlYGo0GDz74oGyrzdra2tiakQzxPA+e58UOQ1JYImIqmsPhwDe+8Q2xwyiZUCiE2tqJGtszlYqNhsZjiYipeEuWLMFtt90mdhglkc1m0d/fL8v9U9WKJaLxWCJiZOGWW27B0qVLxQ6jJJLJJDKZjOyOUa9WLBGNxxIRIwuEENx7772Q6xEig4OD4Hletuth1YQlovFYImJkg+d5PPjgg9BoNGKHUhI9PT2wWCyyPBajmrBENB5LRIysWK1W3H///SBkosOBK18oFGLnGFU4lojGY4mIkZ2FCxdi48ZxJ5TIRjAYZN0XKhhLROOxRMTI0rp167By5UqxwygZ1n2hMhFCYDQaxQ5DclgiYmSJEIJ77rlH1tNYPp8P9fX1YofBzIDRaGQd1ifAfiKMbKlUKjzwwAOyPnQuEAjA5XKJHQZTJDYtNzGWiBhZM5vN2L59u2yLFwAgHA7DbreLHQZTBJaIJsYSESN7TU1N2Lx5s9hhlEwmk8HQ0BA7froCsN/RxJRiB8Aw5bBmzRp0dnbi1KlTYodSEvF4HCqVClqtFvF4XOxwqpZCoYDBYCi86fX6MR+zadSJsUTEVAVCCO6++26Ew2EEg0GxwymJgYEB2Gw2ZDIZpNNpscORDULIpMnl2kTD87ysp4FLhSUipmoolUrs2LEDL7zwAgYHB8UOpyQikQhcLhdCoRByuZzY4UiaXq+fMqnk37RaLUsuJcYSEVNVTCYTduzYgf/+7/+W7Qt1R0cH6uvrEQgExA6l7DQazZRTY/nP6fV6VkYtISwRMVWnoaEBd911F9566y2xQymZQCCAxsZG+Hw+sUOZM7VaXdTUmF6vh1LJXtIqEfutMVVp9erV6OzsxPHjx8UOpWR8Ph+8Xi/8fr/YoYyjVCqnXW/Jf451HJc/loiYqvW1r30NXV1dknyhForf74fH4ylLgQbHcUWtuej1eraoz4zBEhFTtRQKBbZv344XXngBAwMDYodTMh0dHXA6nQiFQrO6/bUJZbJkwxb1mdliiYipagaDAQ888ABeeuklZDIZscMpiVwuh97eXlitVvT09AAAtFptUVNjOp2OLeozJccSEVP13G43tmzZgl27dokdypwRQuBwOGC326dcd2GH6zFSwhIRwwBobm5GZ2cnDh8+LHYoM+JwOOB2u+FyueB2u1FXVweVSiV2WAwzIywRMcyITZs2IRQK4erVq2KHMiGbzQa32114q6urYxVljCywRMQwIziOw7Zt2/DCCy+gr69P1FisVuuYkY7L5QLP86LGxDClwhIRw4yi0+nw4IMP4sUXXyxbvzaLxVJIOPmk8/+3d7+xddV1HMffn3aDpuufIKt1HcMJRrPNKMxmMWoWm2UEZmTSiNMYn2DkyUwWk2lM0IQnPAATJfpEoYLREHgwwGAGG9PgyAjDdLB2dUUGS43b2CjzT5jEwObXB+fXetu13K63955z6ueVnPSce87v9nNPz73f+/ud03tbWloa8rvNisCFyGya7u5utm7dyq5duxb8vjs7Oy/q6SzmL+4zmwsXIrMZrFu3jtOnT3PgwIF530dHR8dFPZ1ly5YtYEqzxWFOhUjSGPAWcAE4HxG9kj4B/AxoA8aAr0XERf8VKGkH8E1AwP0Rce/CRDerr76+Ps6cOcOxY8eqbtvW1jalp9PT00NbW1sDUpqV36X0iPoi4s2K5QFgZ0Tsl3Qb8B3gB5UNJH2MrAhtAN4B9kjaHRHVn9lmOWtqaqK/v5+BgQHOnj07eXtraysrV66cUnTa29tzTGpWbrUMzX0UeDbN7wP2Mq0QAWuAgxHxNoCk/cAtwD01/F6zhmlpaWHbtm0MDw9PFp2Ojg5/lI3ZAprrZ3cE8LSkQ5JuT7eNADen+VuBVTO0GwE2SrpSUiuwZZbtkHS7pEFJg+Pj43N/BGZ11tXVxaZNm1izZg2dnZ0uQmYLbK6F6DMRsR64CdguaSNwW5o/BLSTDb1NERGjwN1kPaY9wBAw4wd6RcR9EdEbEb1dXV2X/kjMzKyU5lSIIuJU+vkG8DiwISJejogbIuKTwMPAa7O0/UVErI+IjcDfAJ8fMjOzSVULkaRlkton5oEbgBFJ70+3NQHfJ7uCbqb2E9tdDfSTFS0zMzNgbj2ibuCApCHgj8DuiNgDfFXSK8DLwCngQQBJPZKerGj/qKSjwG+B7RHx9wV9BGZmVmqKiLwzXKS3tzcW81c4m5n9P5B0KCJ6q23nb7wyM7NcuRCZmVmuXIjMzCxXLkRmZpYrFyIzM8tVIa+akzQO/KXGu1kOvFl1q2IpW+ay5QVnboSy5QVnrpcPRkTVj8opZCFaCJIG53LZYJGULXPZ8oIzN0LZ8oIz581Dc2ZmlisXIjMzy9ViLkT35R1gHsqWuWx5wZkboWx5wZlztWjPEZmZWTks5h6RmZmVgAuRmZnlqlSFSNIqSc9IGpX0J0k7pq3fKSkkLZ+l/QVJh9P0REkyXy3p6dT+qKTVRc4sqa9iHx+W9G9JXyxy5rT+ntRuVNJPVOfvA1+AvHdLGknTtnpmrZZZ0p2STlb8zbfM0v5GSX+W9Kqk75Uk8wOS3pA0UvS81Y6pQouI0kzACmB9mm8HXgHWpuVVwF6yf4RdPkv7cyXM/Adgc5pvA1qLnrnift5H9q28hc4MfBp4DmhO0/PA5wqc9/PAPmAJsAwYBDry2sfAncDOKm2byb7F+RrgMmBo4vEWNXNqsxFYD4zUO+sC7ONZj6miT6XqEUXE6xHxYpp/CxgFVqbVPwa+CxTq6otaMktaCyyJiH2p/bmIeLvImaf5EvBUCTIH0EL2Ank5sBQ4U+C8a4H9EXE+Iv5F9qJ+Yz3zppzvlbmaDcCrEXE8It4BHgG21ifp/9SYmYh4luzNVEPUkrfWx5qnUhWiSmmI6nrgBUk3AycjYqhKsxZJg5IONmK4aLp5ZP4I8A9Jj0l6SdIPJTU3IOqkee7nCV8hh6+Gv9TMEfE88Azwepr2RsRoA6IC89rHQ8BNklrT0F0fWS+qYSozp5u+JWk4DWVdMUOTlcBfK5ZP0OAXyXlkzlUteWdoW2x5d8nmM5ENUR0C+oFWsp3dmdaNMfswV0/6eU3a7toiZybrUfwz5V0CPAp8o8iZK9quAMaBpUU/NoAPA7tT2zayobmNRc2b1t0BHCYbonsI2JHHPk7L3WRDb03AXcADM7S5FRioWP468NMiZ65ou5oGDc0tUN4pbcswla5HJGkp2QvyQxHxGHAt8CFgSNIYcBXwoqQPTG8bEafSz+Nk516uL3jmE8BLkQ1nnAd+QzZeXeTME74MPB4R7zYiL9SU+RbgYGRDn+eAp4BPFTgvEXFXRFwXEZsBAcfqnXeWzETEmYi4EBH/Ae4nG4ab7gRTe21XAafqnRdqypyLWvLO1LYU8q6El/guQcCvgHvfY5sxZn7XewVweZpfTvbEbcTJ0loyN5MNw3Sl5QeB7UXOXLH+INBXkmNjG/A7sl7nUuD3wBcKnLcZuDLNfxwYITuXmMs+BlZUzH8beGSGtkuA42SFduJihXVFzlyxfjWNu1ihln1c9Zgq6pR7gEv8I32W7ATuMNmwxGFgy7RtJp+8QC9pOIDsyqgj6QlwhAYNcdWSOS1vTm2PAL8ELitB5tXASaCpJMdGM/BzspO7R4EfFTxvS8p5lKzgX5fnPgZ+nY7PYeCJiRdNoAd4sqL9FrIruV4D7ihJ5ofJzhu+S9arq+vrRi1553JMFXXyR/yYmVmuSneOyMzMFhcXIjMzy5ULkZmZ5cqFyMzMcuVCZGZmuXIhMjOzXLkQmZlZrv4LgT1XH4z7SE4AAAAASUVORK5CYII=\n", "text/plain": [ "
" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "# Create a figure with one subplot\n", "fig, ax = plt.subplots()\n", "\n", "# Plot polygons\n", "polys.plot(ax=ax, facecolor='gray')\n", "southern.plot(ax=ax, facecolor='red')\n", "\n", "# Plot points\n", "pip_data.plot(ax=ax, color='gold', markersize=2)\n", "\n", "plt.tight_layout()" ] }, { "cell_type": "markdown", "metadata": { "deletable": true, "editable": true }, "source": [ "Perfect! Now we only have the (golden) points that, indeed, are inside the red Polygon which is exactly what we wanted!" ] } ], "metadata": { "anaconda-cloud": {}, "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.5" } }, "nbformat": 4, "nbformat_minor": 2 }