Exploring Geopandas Basics in Python.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Exploring Geopandas Basics in Python\n",
    "--------------------\n",
    "http://geopandas.org/\n",
    "GeoPandas is an open source project to make working with geospatial data in python easier. GeoPandas extends the datatypes used by pandas to allow spatial operations on geometric types. Geometric operations are performed by shapely. Geopandas further depends on fiona for file access and descartes and matplotlib for plotting.\n",
    "\n",
    "The goal of GeoPandas is to make working with geospatial data in python easier. It combines the capabilities of pandas and shapely, providing geospatial operations in pandas and a high-level interface to multiple geometries to shapely. GeoPandas enables you to easily do operations in python that would otherwise require a spatial database such as PostGIS.\n",
    "\n",
    "How to install geopandas: http://geopandas.org/install.html#installing-geopandas\n",
    "\n",
    "## Dependencies\n",
    "Installation via conda should also install all dependencies, but a complete list is as follows:\n",
    "* numpy - Fundamental package for scientific computing with Python http://www.numpy.org/\n",
    "* pandas (version 0.15.2 or later) - High-performance, easy-to-use data structures and data analysis tools http://pandas.pydata.org/\n",
    "* shapely - Python package for manipulation and analysis of planar geometric objects (based on widely deployed GEOS) https://pypi.python.org/pypi/Shapely\n",
    "* fiona - Reading and writing spatial data (alternative for geopandas) https://pypi.python.org/pypi/Fiona/1.7.9.post1\n",
    "* six - provides utility functions for compatiblity on both Python versions https://pypi.python.org/pypi/six\n",
    "* pyproj - Performs cartographic transformations and geodetic computations https://pypi.python.org/pypi/pyproj\n",
    "\n",
    "Further, optional dependencies are:\n",
    "* geopy 0.99 (optional; for geocoding) -  Geocoding library: coordinates to address <-> address to coordinates. https://pypi.python.org/pypi/geopy/0.99\n",
    "* psycopg2 (optional; for PostGIS connection) https://pypi.python.org/pypi/psycopg2\n",
    "* rtree (optional; spatial index to improve performance) https://pypi.python.org/pypi/Rtree/\n",
    "\n",
    "## Plotting Packages\n",
    "* matplotlib - 2D plotting library which produces publication quality figures in a variety of hardcopy formats and interactive environments https://matplotlib.org/\n",
    "* descartes - Use geometric objects as matplotlib paths and patches https://pypi.python.org/pypi/descartes\n",
    "* pysal - library of spatial analysis functions to support the development of high level applications http://pysal.readthedocs.io/en/latest/\n",
    "* shapely - manipulation of geometrict shapes https://pypi.python.org/pypi/Shapely\n",
    "* bokeh - interactive mapping https://bokeh.pydata.org/en/latest/\n",
    "\n",
    "There was one absolutely amazing site with lessons from the University of Helsinki on GeoPandas, much of which was adopted and modified for this notebook.  Take a look at https://automating-gis-processes.github.io/2016/index.html for a more detailed walk through."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "%matplotlib inline\n",
    "import geopandas as gpd\n",
    "import pandas as pd\n",
    "\n",
    "# Visualisation\n",
    "import matplotlib \n",
    "import matplotlib.pyplot as plt\n",
    "import matplotlib.pylab as pylab"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Importing Shapefiles\n",
    "---------------\n",
    "Shapefiles are a common Geographic Information System (GIS) file type used throughout the industry. Shapefiles are technically made of several files.\n",
    "\n",
    "Mandatory files \n",
    "* .shp — shape format; the feature geometry itself\n",
    "* .shx — shape index format; a positional index of the feature geometry to allow seeking forwards and backwards quickly\n",
    "* .dbf — attribute format; columnar attributes for each shape, in dBase IV format\n",
    "\n",
    "To import a shapefile you can use the following:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "pntHistMark = gpd.read_file('VT_Roadside_Historic_Markers.shp')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Because shapefiles can be large, downloads are usually zipped up. In order to conserve space and not have all those subfiles showing up, you can actually import it directly from a zip. VFS stands for Virtual File System which is an abstraction layer on top of a more concrete file system. The purpose of a VFS is to allow client applications to access different types of concrete file systems in a uniform way."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "vctCntyBound = gpd.read_file('/VT_Boundaries__RPC_polygons.shp', vfs='zip://VT_Boundaries__RPC_polygons.zip')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Saving a File\n",
    "-------\n",
    "Like in Pandas, it is not hard to save to a shapefile.  Adjust the name and location to best fit your case and see what happens."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# Select first 50 rows\n",
    "selection = pntHistMark[0:50]\n",
    "\n",
    "# Write those rows into a new Shapefile\n",
    "#Change the filename to the correct path, unlik the read_file my Windows computer required the complete path\n",
    "selection.to_file(driver='ESRI Shapefile',filename=r'C:/Users/kendall.fortney/Desktop/PythonGIS/VTHistoric_Markers_Redone.shp')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Pandas\n",
    "\n",
    "Like Pandas, Geopandas contains all the functionality that can be applied to non-GIS dataframes. Most of the functions used in Numpy and Pandas can be applied here to non-geometry columns"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style>\n",
       "    .dataframe thead tr:only-child th {\n",
       "        text-align: right;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: left;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>OBJECTID</th>\n",
       "      <th>RPC</th>\n",
       "      <th>INITIALS</th>\n",
       "      <th>SHORTNAME</th>\n",
       "      <th>LONGNAME</th>\n",
       "      <th>ShapeSTAre</th>\n",
       "      <th>ShapeSTLen</th>\n",
       "      <th>geometry</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>1</td>\n",
       "      <td>CC</td>\n",
       "      <td>CCRPC</td>\n",
       "      <td>Chittenden County RPC</td>\n",
       "      <td>Chittenden County Regional Planning Commission</td>\n",
       "      <td>1.609237e+09</td>\n",
       "      <td>201814.709240</td>\n",
       "      <td>POLYGON ((452203.9672731007 241948.2255608736,...</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>2</td>\n",
       "      <td>SW</td>\n",
       "      <td>SWCRPC</td>\n",
       "      <td>Southern Windsor County RPC</td>\n",
       "      <td>Southern Windsor County Regional Planning Comm...</td>\n",
       "      <td>8.908414e+08</td>\n",
       "      <td>128848.691548</td>\n",
       "      <td>POLYGON ((497499.9973051278 114802.4954725726,...</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>3</td>\n",
       "      <td>BC</td>\n",
       "      <td>BCRC</td>\n",
       "      <td>Bennington County RC</td>\n",
       "      <td>Bennington County Regional Commission</td>\n",
       "      <td>1.488553e+09</td>\n",
       "      <td>203533.897896</td>\n",
       "      <td>POLYGON ((451638.0272736274 89835.76545525473,...</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "   OBJECTID RPC INITIALS                    SHORTNAME  \\\n",
       "0         1  CC    CCRPC        Chittenden County RPC   \n",
       "1         2  SW   SWCRPC  Southern Windsor County RPC   \n",
       "2         3  BC     BCRC         Bennington County RC   \n",
       "\n",
       "                                            LONGNAME    ShapeSTAre  \\\n",
       "0     Chittenden County Regional Planning Commission  1.609237e+09   \n",
       "1  Southern Windsor County Regional Planning Comm...  8.908414e+08   \n",
       "2              Bennington County Regional Commission  1.488553e+09   \n",
       "\n",
       "      ShapeSTLen                                           geometry  \n",
       "0  201814.709240  POLYGON ((452203.9672731007 241948.2255608736,...  \n",
       "1  128848.691548  POLYGON ((497499.9973051278 114802.4954725726,...  \n",
       "2  203533.897896  POLYGON ((451638.0272736274 89835.76545525473,...  "
      ]
     },
     "execution_count": 5,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "vctCntyBound.head(3)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style>\n",
       "    .dataframe thead tr:only-child th {\n",
       "        text-align: right;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: left;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>OBJECTID</th>\n",
       "      <th>name</th>\n",
       "      <th>descriptio</th>\n",
       "      <th>website</th>\n",
       "      <th>town</th>\n",
       "      <th>county</th>\n",
       "      <th>address</th>\n",
       "      <th>video_url</th>\n",
       "      <th>coordinate</th>\n",
       "      <th>coordina_1</th>\n",
       "      <th>Organizati</th>\n",
       "      <th>YearCast</th>\n",
       "      <th>FiscalYear</th>\n",
       "      <th>geometry</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>1</td>\n",
       "      <td>\"Snowflake\" Bentley -- Jericho's world famous ...</td>\n",
       "      <td>For fifty years Wilson A. Bentley, a farmer an...</td>\n",
       "      <td>None</td>\n",
       "      <td>Jericho</td>\n",
       "      <td>CHITTENDEN</td>\n",
       "      <td>Jericho Center Cir., in village</td>\n",
       "      <td>None</td>\n",
       "      <td>-72.973052</td>\n",
       "      <td>44.469608</td>\n",
       "      <td>Vermont Division for Historic Preservation</td>\n",
       "      <td>2003</td>\n",
       "      <td>2016.0</td>\n",
       "      <td>POINT (-72.97305209009137 44.46960766018835)</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>2</td>\n",
       "      <td>\"Vermont's First College\"</td>\n",
       "      <td>On this site Castleton State College, Vermonts...</td>\n",
       "      <td>None</td>\n",
       "      <td>Castleton</td>\n",
       "      <td>RUTLAND</td>\n",
       "      <td>VT Route 4A in front of the town library</td>\n",
       "      <td>None</td>\n",
       "      <td>-73.178299</td>\n",
       "      <td>43.611718</td>\n",
       "      <td>Vermont Division for Historic Preservation</td>\n",
       "      <td>1986</td>\n",
       "      <td>NaN</td>\n",
       "      <td>POINT (-73.17829915991881 43.61171830002038)</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>3</td>\n",
       "      <td>1st Normal School -- Pioneer in Teacher Training</td>\n",
       "      <td>The first recognized school in America for the...</td>\n",
       "      <td>None</td>\n",
       "      <td>Concord</td>\n",
       "      <td>ESSEX</td>\n",
       "      <td>U.S. Route 2, at Shadow Lake Road</td>\n",
       "      <td>None</td>\n",
       "      <td>-71.889163</td>\n",
       "      <td>44.428617</td>\n",
       "      <td>Vermont Historic Sites Commission</td>\n",
       "      <td>1954</td>\n",
       "      <td>NaN</td>\n",
       "      <td>POINT (-71.88916283008687 44.42861686006097)</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "   OBJECTID                                               name  \\\n",
       "0         1  \"Snowflake\" Bentley -- Jericho's world famous ...   \n",
       "1         2                          \"Vermont's First College\"   \n",
       "2         3   1st Normal School -- Pioneer in Teacher Training   \n",
       "\n",
       "                                          descriptio website       town  \\\n",
       "0  For fifty years Wilson A. Bentley, a farmer an...    None    Jericho   \n",
       "1  On this site Castleton State College, Vermonts...    None  Castleton   \n",
       "2  The first recognized school in America for the...    None    Concord   \n",
       "\n",
       "       county                                   address video_url  coordinate  \\\n",
       "0  CHITTENDEN           Jericho Center Cir., in village      None  -72.973052   \n",
       "1     RUTLAND  VT Route 4A in front of the town library      None  -73.178299   \n",
       "2       ESSEX         U.S. Route 2, at Shadow Lake Road      None  -71.889163   \n",
       "\n",
       "   coordina_1                                  Organizati  YearCast  \\\n",
       "0   44.469608  Vermont Division for Historic Preservation      2003   \n",
       "1   43.611718  Vermont Division for Historic Preservation      1986   \n",
       "2   44.428617           Vermont Historic Sites Commission      1954   \n",
       "\n",
       "   FiscalYear                                      geometry  \n",
       "0      2016.0  POINT (-72.97305209009137 44.46960766018835)  \n",
       "1         NaN  POINT (-73.17829915991881 43.61171830002038)  \n",
       "2         NaN  POINT (-71.88916283008687 44.42861686006097)  "
      ]
     },
     "execution_count": 6,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "pntHistMark.head(3)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "OBJECTID        0\n",
       "name            0\n",
       "descriptio      0\n",
       "website       236\n",
       "town            0\n",
       "county          0\n",
       "address         2\n",
       "video_url     189\n",
       "coordinate      2\n",
       "coordina_1      2\n",
       "Organizati      0\n",
       "YearCast        0\n",
       "FiscalYear    197\n",
       "geometry        0\n",
       "dtype: int64"
      ]
     },
     "execution_count": 7,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# A quick way to get a count of Nulls in the data\n",
    "pntHistMark.isnull().sum()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style>\n",
       "    .dataframe thead tr:only-child th {\n",
       "        text-align: right;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: left;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>YearCast</th>\n",
       "      <th>name</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>45</th>\n",
       "      <td>2017</td>\n",
       "      <td>2</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>44</th>\n",
       "      <td>2016</td>\n",
       "      <td>14</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>43</th>\n",
       "      <td>2015</td>\n",
       "      <td>12</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>42</th>\n",
       "      <td>2014</td>\n",
       "      <td>3</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>41</th>\n",
       "      <td>2013</td>\n",
       "      <td>10</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>40</th>\n",
       "      <td>2012</td>\n",
       "      <td>12</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>39</th>\n",
       "      <td>2011</td>\n",
       "      <td>10</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>38</th>\n",
       "      <td>2010</td>\n",
       "      <td>4</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>37</th>\n",
       "      <td>2009</td>\n",
       "      <td>1</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>36</th>\n",
       "      <td>2008</td>\n",
       "      <td>7</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>35</th>\n",
       "      <td>2007</td>\n",
       "      <td>10</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>34</th>\n",
       "      <td>2006</td>\n",
       "      <td>5</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>33</th>\n",
       "      <td>2005</td>\n",
       "      <td>5</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>32</th>\n",
       "      <td>2004</td>\n",
       "      <td>3</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>31</th>\n",
       "      <td>2003</td>\n",
       "      <td>9</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "    YearCast  name\n",
       "45      2017     2\n",
       "44      2016    14\n",
       "43      2015    12\n",
       "42      2014     3\n",
       "41      2013    10\n",
       "40      2012    12\n",
       "39      2011    10\n",
       "38      2010     4\n",
       "37      2009     1\n",
       "36      2008     7\n",
       "35      2007    10\n",
       "34      2006     5\n",
       "33      2005     5\n",
       "32      2004     3\n",
       "31      2003     9"
      ]
     },
     "execution_count": 8,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Or create counts or other mathmatical functions like in this case a count of new Roadside Marks by year\n",
    "yeargroup = pntHistMark[['name', 'YearCast']].groupby(['YearCast'], as_index=False).count().sort_values(by='YearCast',ascending=0)\n",
    "\n",
    "# View of abbreviated list\n",
    "yeargroup.head(15)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "<class 'geopandas.geodataframe.GeoDataFrame'>\n",
      "RangeIndex: 236 entries, 0 to 235\n",
      "Data columns (total 14 columns):\n",
      "OBJECTID      236 non-null int64\n",
      "name          236 non-null object\n",
      "descriptio    236 non-null object\n",
      "website       0 non-null object\n",
      "town          236 non-null object\n",
      "county        236 non-null object\n",
      "address       234 non-null object\n",
      "video_url     47 non-null object\n",
      "coordinate    234 non-null float64\n",
      "coordina_1    234 non-null float64\n",
      "Organizati    236 non-null object\n",
      "YearCast      236 non-null int64\n",
      "FiscalYear    39 non-null float64\n",
      "geometry      236 non-null object\n",
      "dtypes: float64(3), int64(2), object(9)\n",
      "memory usage: 25.9+ KB\n"
     ]
    }
   ],
   "source": [
    "pntHistMark.info()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "1991.6101694915253\n",
      "2002.0\n"
     ]
    }
   ],
   "source": [
    "# The mathmatical function in Pandas still apply here\n",
    "print(pntHistMark.YearCast.mean())\n",
    "print(pntHistMark.YearCast.median())"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Plotting\n",
    "-----------------\n",
    "First a quick plot of what is in both of these shape files, one is a series of lines, and the other is a set of polygon vectors. Shapely holds the power to create a variety of points, lines and polygons"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAMUAAAD8CAYAAADHTWCVAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJztnXl4Y2d97z+vdlmWZcvbeMb22DOeNZPZMjMJScgGTQI0\nhbQkpBu5dAmFFNrS3hKgfcK9wL20lEJpKSWUlBRuSyg0TR5Cmoa07M1syWRmMplkFnsWj2dsj3db\n1nL03j90PCNrZEmWziJL7+d59Fg+OjrnlXS+591+7+8rpJQoFIrLOOwugEJRbihRKBQZKFEoFBko\nUSgUGShRKBQZKFEoFBkoUSgUGShRKBQZKFEoFBm47C6A0TQ1Ncmuri67i6EoM/bv3z8spWwuZN+K\nE0VXVxf79u2zuxiKMkMIcarQfVXzSaHIQIlCochAiUKhyECJQqHIQIlCochAiUKhyECJQqHIQIlC\nocig4ibviuFHrw8xMB4hkZQkNEkiKUkm9b9SEteSJDRJPJn6qyUvb0voz7WkxONy4HE68LpTf/0e\nJ26ng7iWZDauEUskiSclmn6suCZxCvC4HDgdDlwOgdMhcDkESQkOAU6HwO104HKm/nqcDoI+F+GA\nh7aQn+7mALXe8vsZ59b+CyGQUhJNJJmJaczEEkRiGtMxjZlogkhcYzaeJKZpJLTU952UoCUlsUSS\nSFwjEtOIJjTimkRKiQTiWpJoPEk0kaQ56OXs6AzbVzZw/xu6CJT4fZTft2kxkZjG+76xn+mYZndR\niub61Y2cG4sggWg8SUxLktCSOB0Cr8tJjddJ0Ocm4HHq2xzUeFwEvC7qa9xE40kAklISiWnEk0li\nidQjrqXEK5E4hEAIQTIpL50jrkmiCY1oIkkkpjEbTz2PJlLHdDsFWjJ1oZvBw3dt5L6dnZwbj1Dn\nc5csCFCi4IevDy5pQQCcujhN/9jsAq/GLS3LFWfXzMsWs6o5wK9dtxK308Hq5lrDjlv1fYqVjQGE\nsLsUpTEwPovTscQ/RBGcujjDCycvGn7cqhfFhrY6HnjjKruLURJJCcvqfHYXw3K0pOThJ19hcHKh\nWrI4ql4UAI21HruLUDLhgNvuItjCyeFpRqZjhh5TiQKYiCTsLkLJ+D3V2T30OFOjdkaiRAG4nEu/\nPV6t6U9jWpKTQ9OGHlOJAkiYOEJiFY6lPlpQAn/65GHGI8aNslW9KJ56+Rxf/tEJu4uhKIELE1H+\n8+gFw45X9aL4/HOvmzqWbhVL/xMUz7t2dHD3tnbDjlfVohiZjnFy2Nj2qF0kq7RPUedz8ds3dRt6\nzKoWxSvnxu0ugmFoZsVRlDFBn4u/uGcLPS1BQ49btaLQkpJ/+Gmf3cUwjJgea1RNfOadm7n9qmWG\nH7dqRTE8FeU/jw7aXQzDiCaWdvzWYgl4nLxxTUFpnBZN1YqizuemJei1uxiGEYlVV03hEIKpqDmT\nrlUrCr/HyefftZXGwNIP8QBoCVbG5ygEp0Pwu7f10GpSvFfVigLg+p4m/uE9O/G5l/bX4Hc7ODMa\nsbsYlvGGVY289+bVph1/aV8NBrC5vZ4/un2d3cW4gjq/iw1tQXpaAnQ0+GkJeqnzu3BnifO5ur2e\n4Sljg+LKmf2nRomYuAamOqPIMljbauyQXqn43Q6aAl5eHZjM+no44MHtEFyYjNLVWMO+vhGLS2gv\nq1sC+D1O046ft6YQQnQIIf5LCHFECPGKEOL39O0fF0L0CyEO6I+3pr3nI0KI40KI14QQd6Rtv0YI\ncUh/7QtCpAJ2hBBeIcTj+vbdQoiutPfcL4Q4pj/uN/LDQyqQ7pnDA0YftiSW1/tzTiqOTMcYn42z\nsa2O+hqPaUs9y5Ff2t7Ot977BlPPUUhNkQD+UEr5ohAiCOwXQjynv/Y5KeVfpO8shNgI3AdcBSwH\nvi+EWCul1IAvAb8N7Aa+B9wJPAP8JjAqpewRQtwH/BnwLiFEGHgY2EEqkmG/EOIpKeVoaR/7Mi+e\nHuOf95wx6nCGUMgqutl4kiMDE3hdgl1dYY4NTjI6Y+/SU7NZ21rLw7+wkRqTw+Tz1hRSygEp5Yv6\n80ngVWBFjre8HfimlDIqpewFjgO7hBBtQJ2U8gWZinP+R+Adae95TH/+beBNei1yB/CclHJEF8Jz\npIRkGE8fLK9aYmdXA8cGpwreP5qQ7OkbIRLTuLY7TMhfmYuNNrbV8Zl3bqHOZ/7nW1RHW2/WbCN1\npwf4gBDioBDiUSFEg75tBZB+6z2rb1uhP8/cPu89UsoEMA405jhWZrkeEELsE0LsGxoaKvjznBia\n4hu7C7YtMB2ngHNjsxQTxjSbSLK7d4R4QmNXd5g6X2V1F29e18yWjnpLzlWwKIQQtcB3gN+XUk6Q\nagqtArYCA8BnTSlhAUgpH5FS7pBS7mhuLnyW88zITFmFR2ztqKd/rLSh1Zl4kj29I2hJybXdYYIV\nIo6+4WnL4rsKEoUQwk1KEP9PSvmvAFLKC1JKTUqZBL4C7NJ37wc60t7erm/r159nbp/3HiGECwgB\nF3McyxAKiSztaTEudUo+nA7jRsinYxq7e0eQujiWOs8duWDZKFsho08C+CrwqpTyL9O2t6Xtdjdw\nWH/+FHCfPqLUDawB9kgpB4AJIcR1+jHfDTyZ9p65kaV3Av+p9zueBW4XQjTozbPb9W2GUEgtMTAW\nYdOKOrqbaow67SWaa+eHmUzMGt9RnoppFRFBm0hKPvydg6Z8R5kUUrfeAPw6cEgIcUDf9lHgl4UQ\nW0mNCvUB7wWQUr4ihPgWcITUyNWD+sgTwPuBrwF+UqNOz+jbvwp8XQhxHBghNXqFlHJECPEJYK++\n3/+WUhp2uyhkFMPndnK4fwKA7Z31nB+f5dx46SlVmmo9TMzG2bS8jmgiybHBKc6aMCvtdghOjcwY\nflw76Ls4w9MHB/jlXZ2mnifvVSGl/AmQbYzwezne8yngU1m27wM2Zdk+C9yzwLEeBR7NV85iuDCR\n/+LW0ppYL54ewyFSI0S9w9MlzSJ3NQbYd2qUw+dSglu/LMjR89kn60phW2cDeypocs/ItdgLUdVh\nHj89Ppzz9WV1PsYyxv6TEvb2jTIRiRfdkW1v8PPi6flTLWYIosbt4PVB449rJ08dOMcZk2u+qhVF\nMinzLkWt8y98wcc0ye7e1B14V1cY9yLS5DTVWjMLfXV7/RWiXuocGZhgxuTcv1Urir19Ixw8m3s5\naiHzBZOzCfb0jdAY8LC9sx6/O3dMzpb2EAfOmL8MNuhzcbi/cpbbpvPJp4/wo9cLn49aLJUxiF0E\nRwYmcr6+rjXITLzwRSznJ6Kcn4jiENDTHCAc8BBPSvpHIwxORtm0vA4JHF/EbHUpbGirY09v5fQl\n0vnxsWEmInFuWmvOyruqFcXPTuTOVl3rdfLahcW3x5MSjg9NQ1rWuuZaL26XA6cQlqX9T1bAMGwu\nDvaP8+LpUbZ3NuTfeZFUrSgGJ6M5XxcGZtwbmooyNJU639aOeg6cGTPs2Aux79QoK8M11AfcvGxB\nc81qpImZ1qu2T/EbN3TlfN2sNEonBidptmht+KmRGV4+M053U4DN7SFLzmkV9TVu0ybyqlYUneGF\nZ6g3t4d4+axh0enzmIxqNFmc+r93eJqDZ8fpaQ6waUWdpec2i/FInMf3nmF4KneNXwxVK4r6misv\nzJ6WWjavCHHw7Dhmxgm+OjDJrm7j28L5OD40zWy8fAIgS0FK2Nc3yuCE8aKo2j5FV2MNN69t5uj5\nCVY2BhiajFo2MgTw4qlRelpqLT0ngNdVGfdBIeB9t6xm43Lja76qFYUQgs/es5m3f/GntgxdJpIw\nGYlT53MxMWuNacxVy+t45Vzuoehy5/aNrfz1r2yjd3iadSatra+M20aRNAV9PPHgDYRtyv10YTJK\ne4PflAjcbJiZAcMqdveO8LMTF1m/rM7QEcJ0qloUAC1BH29Y3Wjb+Y8MTNISNN/EcXtnfUVkWB+P\nxHnPP+xlYNy8PFdVLwqAh+5cT1ejNXfrTAIeZ97Z9VJxO0TJK/rKBb/byXtvXkVbyG/aOZQogI5w\nDdd221NbdIRrTF8EtK2zgQsmjNLYQSSusc3ktdpKFDofun3toiJdjeLo+Ul8bqdpqTvrfM6Kykge\n8rtNF7gShU5rnY9NK+yZ9V3dHDBt/mB9W4iXz46nUnA2B0w5h1V4XQ5u7Gni/uu7TD1P1Q7JZmNX\nV5iXTpsfl5SOxymKsrx1OwVSypyTjI0BDy/rcVZzKTi3doQYnoqZsvTVLJpqPfzKtSt55/Z2OsLm\n9SXmUDVFGp02dLa3djZwcXrxy1q3dtSzsjHAsrqF46i6mwJEM1Rz4Mw458Yi7FjZcEXihHKlJejj\nd2/tobOxxrRh2HSUKNKYiVrf9haw6FioudrlxNA001GNLVmC/TrCVy55nSMpU1G0Y5EYu5ZAbqi/\n+7Vr8Fg4E69Ekcb1PdaPQO3uHSGW0NjYVseurjCtBUTQbu24XLtMRhO8fHacnV0NpF83Ib8775LX\nuCZTs/kSdnWH8brKz6B+echneQ1e3rcIi0naFCu3flnoUsYNp0Owo6uBk0PTjGRpVvlc2ZMR7O0b\nZVVzgEhMI+C5nJanECajCfb0ppbUbm4PsK9vtCx8ubd21PPlX7/G8vOqmiKNx/67z/JztoV885o5\nWlKyry9lSpLt7r25Y+FkBCeHppmIxEkUOe9xcTrG3r5RVrfUEq6xN1HzutYgf3XfVtMsvHKhagqd\n44NTfHv/2fw7GsyyOh8DWZKrReIae3pHaA566Wjw8+LpMQIeJ0fzzH5PxzSmL5aWAub44BQb24KM\nWJgJpKuxhnNjs/zxnevwuZ3cvW0FAa89l6cShY5dKfldToHTIRac1R6ajDI0GWVtay1Bn5v9p8xZ\n/JROOOAuapi4GH5pezs3rW3irs3L6bs4zapm63L3LoQSBalF/vkSo5nF3r5R2hv8jE7HciY1eP2C\ndesuVjfXsrfPfPE11Xr5ozvWXopjKgdBgOpTAPBXzx+zNbVkU8BrWZaPfLSFfLxoQW0E8Cdv22Bq\nYF+xKFGQP32mGVy9oo41rbV4XYLzBeS0tYq2kA/NoqGnQmzM7ECJAkxZ0piPsUicYxemSErKRhTd\nTQFeNCHMJXPizSHgXTs6WF0mzaVMVJ8C+MXt7fzjf1tn85W+Njtu1W25AAIG2fCG/G48LgfRuMYv\nXdPOn75tIy/0XuSL/3Wc1qCPe3Z02LqwKx95RSGE6CBl2thKyoviESnlX+nOpY8DXaT8Ke6dcy0V\nQnyElOOpBnxQSvmsvv0aLvtTfA/4PSmlFEJ49XNcQ8rB6F1Syj79PfcDf6IX55NSyjnDSMOwKg/T\nHA02zwFkY1d3uOS16i6H4P23rOZDt68DUgMYDr2JdP3qJq5f3VRyOa2gFMvg/wE8L6X8tBDiIeAh\n4MNLzTIYQLPwbu12Co4uYBpvF1ctr2NvCYKor3Hzt7+6nVVNtSwLXZ5sc5RpnyEfpVgGp9v8PsZ8\n+98lYxkM0BAw7s7tcznIFbu2aXmIyag12TsKoS3k4/TITElhHW6ng+tXN80TxFKmFMvgVt3HDuA8\nqeYVLCHL4DmCPjfdTcYswFnR4MftdHL1ilDWFCwxrbySkbXUeZksIcWOx+ngf+rNpUqhFMvgS+h3\nftt6jMVaBqfzjq1XaC0v4YCHTcvr2NZZz9Ur6tjZ1UBDwEMkrnGof5z6jL5Dnd/FayY4FpVCKdnJ\nhYDH33sd9+7syL/zEqKg0adslsHABSFEm5RyQG8aDerbS7EMPpvFMviWjPf8oKBPtkja6heu+gMe\nJ8vr/YRq3DgQTEUT9I9FGJmOZY1knSOeUStMRBJcvSLEoTIyUyk2eBBg84pQ2cxCG0kho09ZLYO5\nbPP7af1vuv3vPwkh/pJUR3vOMlgTQkwIIa4j1fx6N/DXGcf6b9Isg4UQzwL/R7cLhpRl8EeK/rQ5\n+OFrQ/Q0B6jzu3E4BLFEksnZBMNTUSZnExwrIr3ldMaipfoat+npbBZDV2PNpWWqxfDy2XGODkxw\n7aryHV4thlIsgz8NfEsI8ZvAKeBeWFqWwXMc7h/nh68PMWVwB3h0Zn4tsrY1WFbuQkZEof7w9aGK\nE4WQZhkx2MSOHTvkvn37Ct7/O/vP8uHvHCypGbEQTocgmZSXOlt2JFReiM3tobyef4XQEvSy52Nv\nNqBE5iKE2C+l3FHIvlUd5nFxKsrDT71iiiAgtWBobmKwvd5fNoJwCBgtIllCNnZ0WW8pYDZVLYqv\nv3DK8CZTJnPJm1c0lE806PbOBs4YlOKmI4f5zVKlakXx1Z/08vnvHzP9PAFvKp7otMmG6IvByCwx\nNyyR0I3FUJUBgf/12iCffPqIJedyOhysba21dJFQPuYWNrWFfBw4MzYvKDHodRIOeKj1ufG5nbgc\ngqSUHB+cYjRjeWrQ5+LGHiWKiuDHrw+bZvSYiZaUhPzlFwB4djTC2dEIzbUeGmu9TMzGuTgVYzKq\nMRmNAPObV06HYHN7CEFqtE6TMDmb4MmX+7l7W3vWcyxVqlIU+09ZNyw6HU1wrozT4Ad9bo4WMMuu\nJeWl0aqQ38261iDDU1H+Zd+ZihNF1fUppqMJSyfQXr8wSTSRJOg1Zq2C0TQW4dQ6Homzp2+Ek8PT\nnBub5ZlD9iR9MIuqE8Wzr5y3dGFPUkI0kaTFhvxFhVDqaHTfxZmyivo1gqoTxf995qgt5w36yq9f\nAakUOqViVyYUs6g6UeQylTcTKxMEF4rP5eDsaGlDxY0BD79z82qDSlQelN8vZTK5olrNpBzDabqa\nAiU3n25c08SGNusTP5hJVYnicP84vTY5hJZTggJIzTHUGdCke+n0GB/61gEmZ61LsWk2VSWKL/3g\nhG3nnogYc9Fsbg+xflnppuotQS8vnSl9qfvcTH2tTXlfzaByPkkepJT8+Njil6oaxcDELGuaA2iS\nov2s3Q7B2dEZRqbjLA/56AjXcGJoiuGpxTUJg14nJ4enDZnAbKr18JZNbWhJicsGI00zqJqa4vzE\nLBMlrEUulUhMI1zr5eTwNCsba9jV1UDIv7h70paOekamUzXOufFZdveOMDIdY/OKEJsXYWLZ01Jr\n2Iz+hrY63ryhBZezci6lqqkprEgYnI8ZPV/sqYsznLo4g9sp2N5Zz0xMK2hWOVsmwaSEg/ry1s5w\nDc1BDwdOj+VMfel0GHMBr24O8MVf3W6JD52VVI6883Bdd5ieFnvXE2emQYprkhdPj3H0/CSd4Rp2\ndYWpXWDme0t7KK+j6emRGfafGqOx1suuroYFfcGNitj9y3u3GtJZLzeqRhQtdT6eeP/1WU0TrSJX\ncrDTIzPs6RshrkmuWdnA2tb5Ao7EC89KPjgZZU/fKEGfm11dYXxpcyQdDX4GDZiw87ocl8LiK42q\nEQWkZpVvXd9i2/kdBTQzookk+0+N8vqFKTrDNezsauDqFXVFhZ6PTMfY0zeC1+3k2u4wtR6nIanv\nV9T7+df3X09PS+mjYOVI1fQp5jifxUrLKhbb8j49MmNIU2c8Emd37wi1Xhd9F0ufp7lry3LDkseV\nI1UnCqtsq7JRSE1hJlPRRMnLb5tqPTx462pqPJV76VRV8wmgv4zXNiwFhqdiC/rzVQpVJ4prV4Ut\nPV+Nx8mOrgY2Vkh8kNMhONw/wUQFhXVkUnWiaAlau65h04oQ+/pGOTIwYauvnlFoSclHnziEp4Im\n6zKp3E+2ADMx62a1PS4Hx8soYYERrKj386c/vxGfuzKHY6EKRdFvUL6jQtjSHmJkxp5QdbNYFvLx\ncxtb8++4hKk6Udy4xpqULC6H4NTF8sn1ZBT+Cq4h5qg6Udxm0eTdlo56Q2aOy4lVzQH+19uvsrsY\nplN1ougM1/CL2xdv0LJYBsvEBtgobuhp5In331C2Nr9GklcUQohHhRCDQojDads+LoToF0Ic0B9v\nTXvtI0KI40KI14QQd6Rtv0YIcUh/7Qu67wVCCK8Q4nF9+27dQmzuPfcLIY7pj/uN+MBCCD7ylg14\nTVwz3RH2V1SGi5DfzUffuqEsk7qZQSFXxtfIbr74OSnlVv3xPYAMZ9Q7gb8VQsw1QuecUdfoj7lj\nXnJGBT5HyhmVNGfUa4FdwMNp5i0l0Rz08sY1xdmAFcKZkQjRRHl525XCw3dt5Krl9gVSWk0h7qg/\nImWkUghLwhk1riVNzxKY0JJ0NS79jNwrG2sqOqQjG6W0IT4ghDioN6/m7uCWO6MWw0xUM31hTFyT\nl9Lwp9PR4Lc1fH2x3Lujgzs3LbO7GJZS7C3gS8AnSDmifgL4LPAbRhVqsQghHgAeAOjs7My7v8eV\nygT+wsn8tUXI78bnduB1OfG4HDTqF7qUEk2mZniTSUlM0/A4nbx6fuJS5o70RUGdYT9NtV5eOjO2\npPymD54ds7sIllOUKKSUF+aeCyG+AnxX/9cWZ1Qp5SPAI5Cy98pX/o/926GCBAGpWJ8LE5eHVo/n\n2X9Vc+BSJO7gZJR1rbXUet28eGaU0yNLLxix0oP/slFU80nvI8xxNzA3MvUUcJ8+otTNZWfUAWBC\nCHGd3l94N/PdVOdGli45owLPArcLIRr05tnt+raScS6i6dSYpQmUi4aa1P5NtR52djVwbHCK/adH\n5yUKKMO8aAtyfQWasuSjkCHZfyZl5btOCHFWd0P9c3149SBwK/AHkHJGBeacUf+dK51R/57UzfYE\n851RG3Vn1A8BD+nHGiHVNNurPwxzRv3k3ZsIeAqbmV1sPqOklOzsamA8Emdv32jWDHzaElLF/tP2\nJ3ywmqp1R/27H57g0wUkW97VHTbc5jfkdxNLaETi5T9s+/QHb6yI4VjljloAv3Pzau7asjzvfmbc\nNMYjca5urzf8uEbz3ptWVYQgFkvVigLgtvX5J/BmF5FFYzG8OjBBbYFNOLvYuLwyFkYtluqalcng\nuy8P4HM7aAx4qfO7qHG7cDoFUkqiiSQTkTi9Jq3pnpxNcG13mN0GN82M4kM/t5a7NuevSSuRqhXF\nvr4R9vReZDaepH8sQr8Nw/FHz0/iczmYLbOQkFXNAd53y+qceaoqmaptPv3k+DCTUXOaRoWS6luU\nT5vd7RT8/pvX8OSDN+Cu4OWm+ajamsIu85ZMMr2p7aKp1st//MFNWUNTqo2qvB08+8p5vvHCKbuL\ngRBQ73ezqQw6tF/8lW1KEDpVV1PEtSQff+qVkm2tjEBK2HcqNTnWVOthVXMtA2MRzli4jhzgt27s\n5tpVjZaes5ypOlF86ulXGbAxdeZCDE/FGJ5KjUT1tNTSUJMyfZ+0wFPj3p0d+XeqIqpOFMcG8/tA\n2M3xwVRaHLdTsLa1tqjkyoshWWFRDaVSVX2Kidk4e3uXTixPnc/NiUHz80ZNRCpn6awRVJUonj18\nnphWXnMCuehpqc3pSGQE160Ks7PLkFW+FUPViGJ0OsZn/+N1u4tRMC6H4JgFtcTGtlDF2XOVStWI\nIppIMhbJPjdRjhO3V68ImT6Xsrk9xPtuWW3qOZYiVdPRXhby8cY1zTx3JLVo8MN3rmd7Zz1dTQEC\nXhdPHujnY08cznMU65iyIOftHVctoznoNf08S42qEQXAS2kLZu7etmLeWuldXWG8LkdZpKZZvyxY\nkFtqqTTXKkFko2qaT/1jEcYjqZCKHSsbaKydP3u7pjXIlo7yWOPgsqg9d2rEPlencqZqaopXz03o\n9lqSO65adkXAm5SSUwb4wZVCe4Of5fV+w1f6ZVLrdfE7N6/iwVt7TD3PUqVqaoo3b2zlK+/eQTjg\nYX3bla6eUsLn7t3KjpX2DU+2BL2mC+INqxr56Ydv43dvW6NGnRagamoKgJvWNrP7o2/KGhbtcAiu\n72liV3eYD37zJb536Lzl5XNZEK5997YVhGqqIydssVRNTTFHvnUCLqeDP7x9nUWlmY+WZ2JxdXOA\nTcvr2NHVQEORF7bbpWqHfFRVTVEoy0N+lod8nLM4cHAiT/BfyO/mxdOXlwiuX1ZLnd/D6ZGZgv3B\no0sgg4jdKFFkwe9x0t5QY7kozo/nDhnPrOWOnr88493dFKC51gMIogmN8Uic/rHIpRSekApPv3md\nednWKwUliixIKTlwxtpF2yG/+9KQ8ULkimbtHZ6md/jK0bPWoJfGWi81Hif37GinLeQvuayVjhJF\nFoQQlodTL6vz5RVFMU2fC5NRLug2Y9f3VF8KzGKouo52IWhJScLipXkxTaOnOZBzn+kSQz++8Pwx\nfvj6UEnHqAaUKLJgVgK0XPQOz+Bx5U6OZkSA4JmRynNsNRrVfMpCxAZRABw9P8H6ZUHq/G76R2fo\nH7vc0fe6HCVn/qjxOLluVbjUYlY8ShRZiMTsEUVScikQcHVzgM6wn7aQn6HJKDNxreBh14XY1R2m\np+XK2XzFfJQoslAOkbJ+t5MTQ9OXjF7qDZiFnq4gx1YzKdYyOCyEeE638n0u3bW03C2DC8GOPkUm\nmeaLYwYkTdvVrZpOhVCsZfBDwPNSyjXA8/r/S8YyOB/lIIoLJpjThwNq/UQhFGsZnG7z+xjz7X/L\n3jI4H7M2h0K0BL2cMmGUaMMy1Z8ohGKHZFt1HzuA80Cr/nxJWAbnY8aCpaC56G7KPV9RLK1LyJXV\nTkqep9Dv/LZm0xJCPCCE2CeE2Dc0VPrklF1DsnOYsdipMeChu9EcsVUaxYriwpxDqv53UN9eimUw\nWSyDsx3rCqSUj0gpd0gpdzQ3lx7wNmPBkKzH5aC1zsuqpgAh/+WRpY1tdZxPsyg2ik++Y1PV+k0s\nlmJFkW7zez/z7X/L3jI4H1aIYm1rLRcmopwcniaZTLKzqwEhoNZnzii5z13eVmLlRN5fQLcMvgVo\nEkKcJTUi9GngW7p98CngXkhZBgsh5iyDE1xpGfw1wE/KLjjdMvjrumXwCKnRK6SUI0KIOctgMNAy\nOB9WjD550sLAJ6Mae/tGWdUUwG/SxfuJp48wPBXlnh0qmXI+8opCSvnLC7z0pgX2/xTwqSzb9wGb\nsmyfBe5Z4FiPAo/mK6PRmC0Kv8fJubEr106cHJ6mKWiOR8TJoWk6wjWmHLvSUAGBWTB7RrvO51qw\n3yAwp93GlJk9AAAJmklEQVR/67pmrlMeFAWhRJEFs2OfcjWRzFrH8dr5ScbLxEqs3FGxT1lIJJM0\nBjzUeJx4XQ68bgcupwO3w4HTkVqAlEhKovHkpYs4qdsMz8Q0IjGNqRxxRrk6vWYt4zg3Pst3D53j\nV69dac4JKggliixMRzUuTscoZbqgOeglmtDwu5z43A48Licel8DtdOB2OlKjTaQENmfxBZAwySpA\nCLhrS3X6Yi8WJYoslNrRFqT6DSeGokyQe3Z8Xev80Iu4SaJY1xqkzqfyPRWC6lNkodSOtgRODBVW\nzciMYACHSVn73rHNkgiZikCJIgtWRcl2NdZcER+Tqy9SLOuXBXn3G1RfolBU8ykL+ZKSGUGN20Hf\nxRnCgflNmnwZDBdLR9jP0x98I04V4lEwqqbIwkSeVDNG4NcXEY1Mzz9XXEsSNjDXa9DrVoJYJKqm\nyIIVoeMLzVX0XZxBCNjQFqTW6+LE4BQjM3ECHifNQS91PjdetwOBIKppTEQSDE9GiWlJPE4HbpfA\n5XDgdTnwuFJDyNG4hlfFPhWMEkUWrAgd97gWrqSlhFcHUgkMhEgtOhqcjDJ9MffCo2giCRkT5W/e\n0KoEsUiUKDKQUpa88m57Zz2jM/HU/IQzNfHnEFzyg5BS4hCCk1nSXF5ZHhicLD6UfGTa+DD0SkeJ\nIgMjwsYnZxNZ87rawf3Xd9ldhCWH6mhnYMSQaDmkyJlDuRUtHiWKDCZnSx95snuN9xw1Hic3rVFJ\nlReLEkUGkwbMUVixcq8QklJSX2PO+oxKRokiAyMu6HIRxablIbuLsCRRosig1LUUNe7y+UpfuzDJ\n9w4N5N9RMY/y+QXLhFI9IHye8hnQm5xN8KUfnCBWRh3/pYASRQaljj5NzsbZ1llvUGlK51D/OKeV\nJ8WiUKLIoNTM3HFN4naUz9d6384Oelpq7S7GkqJ8fr0yYcqA0SfNYr+8XOzuHTFtNV+lokSRgRFh\n4+XkAxHXkqat+65UlCgyMGKeQiuTqzDoc/HV+3fmDD5UXIn6tjKYipY+o23WktLFIAR87T27WKfS\n7y8aJYoMjIh9Kod+dqpbUx411lKjDH6+8sKI5pP99USKs6NXpuZU5EeJIgMjaorhqdL9ro1AiaI4\nlCgymDUgbmncgjXehfCZZ1/jyQNZLT0UOVCiyGDGgKWo5bKe4oNvWqOSKhdBSaIQQvTpNsAHhBD7\n9G2W2AmbxUy0NFH4yygg8NZ1zbTWKZ+7xWLEL3irlHKrlHKH/r/pdsJmEUskiZU4++t32x8Q2FDj\n5h1bl7O2VQ3HFoMZv+DbSTkfQcoK+AfAh0mzEwZ6deeiXUKIPnQ7YQAhxJyd8DP6ez6uH+vbwN8I\nIYRu/2U4RqTg93scYFP83bbOeu7b2cGt61poUTVE0ZQqCgl8XwihAV+WUj5CbjvhF9LeO2cBHKdA\nO2EhxJyd8HCJ5c7KTLz0kSePy550Ml2NNXzzgevw2nT+SqJUUdwopewXQrQAzwkhjqa/KKWUQgjT\nZ5CEEA8ADwB0dnYWfRwjVsx5XfbMUvRdnOHIuQm2dTbk31mRk5L6FFLKfv3vIPAEsAtr7IQzy2GI\nZXCpnWywN3tGqfmqFCmKFoUQIiCECM49J2Xpexhr7IRNwYgsHHY2X8plfmSpU0rzqRV4Qr8zuoB/\nklL+uxBiLybbCZuFEeky7Upm7HM7uOOq1vw7KvJStCiklCeBLVm2X8QCO2EzMKL5YVeC7xt7mlTi\nM4Mon5mmMiCaWJp9iru3reDLv74j/46KglCiSMOI8Aw77tWbVoSUB4WBKFGkYYQoymh5tqJIlCjS\niBrQ0TbLHD4X58ZUiLiRKFGk8fyrg/l3ykM8af1cweN7z1h+zkpGiUJHS0pe6L1iXnDRTBswAbjo\nc8YSljm6VgNKFDpxLUnAgJSXZpnD50LKVCZAhTEoUegc6h83ZCmqHaIA+Ox/vGZIlK9CieISawxK\nLRnX7Bl+euHkCA98fZ9KpmwAShQ69TUe2kKlr0Gw86L88bFhfnJ8SPUvSkSJIo0NbXUlH8NOUfjc\nDs6MRHjbF37M8JRyRS0WJYo0NrSVvnxz1oBQkaLPHU/y8FOvcGJomsd+1mdbOZY69i8oLiN+701r\n+a0bVxFNJJmNa8wmNKLxuef637hGNJEkGteYjSeJJlJ/Z+MaMS3JTEy7/P745efRtH0jcY1oIvWa\nWXN9zx25wAduW6PyyBaBEkUaHpcDj8s640QpZZqALgsxEtOI6KKKxJLM6PMQkbjGTCxt30v/p0Qb\niSUuiW4qmuDfXurn3p0d+QuimIcShY0IIfC5nfjcal11OaHqVoUiAyUKhSIDJQqFIgMlCoUiAyUK\nhSIDJQqFIgMlCoUiAyUKhSIDJQqFIgNhYhZKWxBCDJHKTGg2TZiU/dwAVNmuZKWUsqBEwxUnCqsQ\nQuxLM6opK1TZSkM1nxSKDJQoFIoMlCiK5xG7C5ADVbYSUH0KhSIDVVMoFBlUpSiEEE4hxEtCiO/q\n/39cCNGv+4EfEEK8NW1fw7y/hRD36/7ix4QQ95MFO73Jc5VvgXKVzfdmKFLKqnsAHwL+Cfiu/v/H\ngT/Kst9G4GXAC3QDJwCn/toe4DpS2fefAd6ib38/8Hf68/uAx/XnYeCk/rdBf96Q5Zx9QFPGtj8H\nHtKfPwT8mdXlW6BcZfO9GfmouppCCNEOvA34+wJ2v+T9LaXsBea8v9vQvb9l6peb8/6ee89j+vNv\nA2/S74Z3AM9JKUeklKPAc8CdBRY7/ZiPZZyrHMqXrbzlWK6CqDpRAJ8H/hjITND0ASHEQSHEo2nN\nk0s+3jpzHt8rKND7G5jz/l7oWJnMeZPv162QIbc3uVXly1YuKJ/vzTCqShRCiJ8HBqWU+zNe+hKw\nCtgKDACftbpsadwopdwKvAV4UAhxU/qL+h3WjiHDbOUqp+/NMKpKFMANwC8IIfqAbwK3CSG+IaW8\nIKXUpJRJ4Cuk/MDBWO/vhY41D2mfN3nO8mUrVzl9b4ZiV2fX7gdwC5c72m1p2/+AVHsY4CrmdxhP\nsnCH8a369geZ32H8lrzcYewl1Vls0J+HM8oUAIJpz39Gqv38GeZ3tP/cyvLlKFdZfG+GXxt2X5xl\nIoqvA4eAg6QM7dN/7I+RGj15DX2kRN++Azisv/Y3XJ4I9QH/QqpzuQdYlfae39C3Hwfek6VMq/SL\n6WXgFeBj+vZG4HngGPD99IvCivLlKFdZfG9GP9SMtkKRQbX1KRSKvChRKBQZKFEoFBkoUSgUGShR\nKBQZKFEoFBkoUSgUGShRKBQZ/H/bfJDz13RVYAAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x2362719bfd0>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "vctCntyBound.plot()\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAARgAAAD8CAYAAACsLLusAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAGQlJREFUeJzt3X+Q3PVdx/Hn+y4L2QulF8q1wh3Xq9oGBwO59ixxYlVQ\nmtpQeoVCQNLR6YzR0TpQajApWOi0NdGzJc7o6CCjw0yQBjG9tsFOik2jTjSMd94laZpEBgmhG2yu\nTa+15CiXu7d/7O6xt7ff/X73x2f3Nn09ZjLcfve7u2+O7IvP9/Pra+6OiEgIbc0uQETOXwoYEQlG\nASMiwShgRCQYBYyIBKOAEZFgFDAiEowCRkSCUcCISDBLGvlhl156qff19TXyI0UkgNHR0e+4e1fc\neQ0NmL6+PkZGRhr5kSISgJm9kOQ8XSKJSDAKGBEJRgEjIsEoYEQkGAWMiASTOGDMrN3Mxsxsd+7x\ng2aWMbPx3J/3hitTRFpRJcPUdwFHgYsLjj3k7n9W35JEFpfhsQxDe45zanKKyzvTbFq7gsH+7tjn\nml3bYpAoYMysB1gHfAa4J2hFIovI8FiGLbsOMzU9A0Bmcoq7d47zyS8fYd3Vl/GPo5l5z23ZdRig\nIV/yUrU18vOTSHqJtB24F5gtOv77ZnbIzP7WzJbXtzSR5hvac3zuC1zoe2en2XHg5ILnpqZnGNpz\nvGm1NfLzk4htwZjZjcBpdx81s18ueOqvgE8BnvvnZ4EPl3j9RmAjQG9vbx1KFmmM4bEMmcmpil+X\nmZziLZufmrtkAYJcxpyKqC3qeDMkuURaA9yU68RdClxsZjvcfUP+BDP7G2B3qRe7+8PAwwADAwO6\nhYG0hOGxDJv+4WDVr3eyQfPRneMsaTOmZ7N/9et5GXN5Z7pkAF7ema7pfesp9hLJ3be4e4+79wG3\nA3vdfYOZXVZw2geAbwSqUaThtuw6NBcKtXBY8D71uIwZHsvw8o/OLTieTrXPtZoWg1oWO/6pma0i\n+zs8Afx2XSoSabLhsQxT08XdjfVVy2VMcedu3vKOFA+876pFM8IFFQaMu+8D9uV+/lCAekSaLknr\nwsj+n7VatVzGRHU8d1ywZF64FI8wfXTnOHfvHKe7gWGjmbwiRcq1LpZ3pDixbR0PrV9V9fsb2S/8\nmm17GR7LVPz6qI7nwrpLhVA+EPP9QNV8dqUUMCJFyrUu3LOtg1r+71/LF314LINFPFdYd9wlWKOG\nsxUwIkU2rV1BOtVe8rnJqem5UFjekarofUsFQ6Vf9KE9x0temhnM69xNcgnWiOFsBYxIkcH+brbe\nvJLuiC9pPhQeeN9VpNqj2hMLRfXZVPJFjzrXmT/sXS4k8xoxnK2AESlhsL+77HBvZnKKwf5uhj54\nDd2daQzo7kyz5qcuod2yodNmyb5gr08nbwlFnVschsUhWRyDjRrObuievCKtpNyli/FaX0xUf8ya\nbXsTzQR++dVzifp1hscyvPzqwrkv8Fqn8XVXdvH1YxMLhqabNWStgBGJUO7SxYG7d44ztOd45Jc1\n6aXP9IzzsSeys4bLfemH9hxneiZ6cDwzOcWOAyfnPS6eNZwPmXx4hg4ZBYxIhKip+IXKTf1P8vq8\nGXe27DrMyAtneOrQS3zv7DQAnekUD96UnTxXTadsYSdyM1Zem3vjlgcNDAy4blsirSJqxmwpHak2\nli+7cN4lCLDg9dVM0Eu1GUO3XsPQnuNVLb40osOuM51i2YVLKr50MrNRdx+IO0+dvCIR8h2l6VT8\n1+Ts9CyZyal5ixxHXjgz19Ga7wS+c3Vv7OhOselZn7sUq/S1kA2XqNbP5NT0vLrrPQFPASNSxmB/\nN5csu7Di1znwWK4/ZP/m63l+2zr2b76egTdfwtIEgVXsVG7UqjCwOtOp2GHy/GhR0iHpek/AUx+M\nSIxqJ6Q52U7VqPVBxcpdPuUDonDUangsw4NfOsLk1HRkDVPTM9y9c5x0qo32NmMmwQrxek7AUwtG\nJEYtE9IK+z2iFinCa5dPpVokqTZbMGdleCzDpicPlg2XQlPTs4nCBeo7AU8BIxKjlglp+fkyEN0y\nMJi7fFp2wfyLis50iqFbr1nQ8frJLx8pO2RdSX2F6j0BTwEjEmOwvztygWGc/GXS/cOHy17+5C+f\nClsk6VT73BB1sfwwdq0c5nVCb715ZV2HrdUHI5LAQ+tXcffO8apeWzwBrlC+xRC1gfeDXzqy4Atf\nz1Ge7s40+zdfX7f3K6YWjEiMfOuiUDrVTmcFa4ii9Cxfyt07xyPnt0xOTS8IlHqN8jRiPZICRiRG\nVOsiSQdr3KXVs6dfTvT5haqZbNeZTrFhdW/Qy6FSdIkkEqOaLzRkw+XO1b08/syLzNQwY764c7jd\nLPH7ndi2rurPrQcFjEiMcvNT0qn2kksBive9jeqDSaJ42DhpuOS3jWgmBYxIjHJf5603r4zdBuHT\ngysB5loyZpBe0sbZBHcuKNVP0p1wEeUd114Re05oWuwoEqNv81ORz1V6CZJkAWW+3REVWKXeo82y\nQeiebbncce0Vc8EWQtLFjmrBiARUvNHTyz86F7s624ENq3sjA6J4b5eOC9p5+dXX3vOCJcbAmy+p\n279DLTSKJFJkeCzDmm17ecvmp1izbS8dEYsTC2fpRr3Pll2H561WTjq1//FnXiz7/GB/N/s3X8+d\nq3vnhQtklwXcs3O8IbcliZM4YMys3czGzGx30fGPmZmb2aX1L0+ksfJrfApD4ZVzpftK8rva3fk3\n/1Hy+XJrj+Ik6cgdHstEdh7PUr/5MrWopAVzF3C08ICZXQG8G6i+i1xkESm1xidujeD+585w//Dh\nBcdrWZUcNwI0PJbhnifKzyxuxG1J4iTqgzGzHmAd8BngnoKnHgLuBb5Y/9JEGq/aNT47DpxcsNl2\nJVtmFosbAfr4rkOxwdeI25LESdqC2U42SObaimb2fiDj7gdDFCbSaop3hrvuyq6qdqAr18Gb7x+K\nG+Juo7ZV4PUSGzBmdiNw2t1HC451AB8HPpHg9RvNbMTMRiYmJmoqViS0eqwvguxSgq8fmyh7b6JS\nlrQZjx04WfK+1YWdxuWkU218bv2qhtyWJE7sPBgz2wp8CDgHLAUuBr4CvAs4mzutBzgFvNPd/zfq\nvTQPRha74bFM1aumixnwfME8mcIh6ySzz4pnEKfaIG5unhk8vzX88oC6bfrt7lvcvcfd+4Dbgb3u\nfou7v9Hd+3LHvwW8vVy4iLSCwf5uNqzurXr/l0LFfSD5oeXnE07OKw6hBBN/ufPa3oTVNYbmwYgU\n+fTgSh5av6qm90i1L9zmMi/k/JSQs3erUVHAuPs+d7+xxPE+d/9O/coSaa7B/u6a+mOWXbAksg8k\n1PyU4vtTLwZqwYhEqGUx8vfLzNgNMT+lUTezr5QCRiTCZA373pa791E95qdsX7+q4ZtHVUOLHUUi\nRE2Uazdj1r3sSNDU9Cz3Dx8u2Seyae0KNj15sKa7Amz9p2/yzH03VP36RlHAiETYtHbFgm0RUm3G\nBUvaFiwwLGXHgZPsOHCSjlQbF6bamTw7PTfLd9kFSxIvfCzl2//3atWvbSQFjEiEwm0RMpNTGNn7\nRE8nCJdCZ6dn52be5mf5VrsIstWoD0akjMH+bjatXUGqzRJNjktianpmUWxn2QgKGJEYD37pCNMJ\nb7ua1Ix7VeuUWo0CRiRGLX0l5by99/VVt2QW45yXUhQwIk2y/7kzrP7J5RW/brHOeSlFASMSY3lH\nfVZYl7L/uTO86XUXxJ6Xb+ks5jkvpShgRGI88L6raG8L1yl7w1U/UTJkCj8x32dT6i4Di5kCRiTG\nYH83n731mmDv//fPnOQHr8wfti51s7ep6ZlFsc9uJRQwIgnkt3EIYdZZMC8masxqMeyzWwkFjEhC\nnx5cGSxkkloM++xWQgEjUoFPD65k+/pVTZnDkmqL3mNmsVLAiFRosL+bW97RHTuHpd7dwhctjd5j\nZrFSwIhUaHgswz+OZmJvjlbvu77Xsn1Es2ixo0gChRt2t5kluvNid2eaU9+fIsGpibRa/wuoBSMS\nq/ge00nCBeC6K7vKhkslywRaafZuIQWMSIxq7zH99WMTkWuGujvTzFbQtGml2buFFDAiMaqde3Jq\ncopNa1csGHHKt0aSXvIs70i1ZLiAAkYkVrV9H21mfHTnOBcuaWN5R2rB/rmlwqeUH75yLuitTkJS\nwIjEyG84VamZ3L69k1PTvDI9y0PrV3HdlV187ImD9G1+io89cZC3974+duuF6VlvuSUCeQoYkRiD\n/d1ctLS2Adep6Rk+vusQOw6cnOsknnFn/3NnuO7KLravX0W5DGu1JQJ5iQPGzNrNbMzMducef8rM\nDpnZuJl91cwuD1emSHPVYw7K2Yh7vz72zEkG+7vLzptpxSFqqKwFcxdwtODxkLtf7e6rgN3AJ+pa\nmcgikvQLvuanLql4tzn37FB4uUGlVhyihoQBY2Y9wDrgkfwxd/9BwSnLqP/ERZFFo+8N8aGRTrXx\n2G/9fOTIUTmbnjwY+Vybcd6PIm0H7gXmtfHM7DNm9iJwJ2rByHns3//nTOw5W2++GsiGwdabVy64\n82JHmbs9lrsJ269f29wV3LWI7bkysxuB0+4+ama/XPicu98H3GdmW4CPAA+UeP1GYCNAb2/r/qLk\nx1fc5UuqDYZuXTWvlTHY312y1XHPE+MU3qCgzaDcDQs2rO4teXfIVpGkBbMGuMnMTgCfB643sx1F\n5zwG3FLqxe7+sLsPuPtAV1dXTcWKNEPcEHFxuEQZ7O/mc7fNv6d0XOuklcMFErRg3H0LsAUg14L5\nA3ffYGZvdfdnc6e9HzgWrEqRJoobIh7aczxxH0lxy2bNtr2R53amw2023ii1DO5vM7MVZPtlXgB+\npz4liSwunR0pvldmmLraOSr3D2cXUEZ58KarqnrfxaSigHH3fcC+3M8lL4lEzifDYxm+HzMHptQQ\nduH2Dvkb3he2XO4fPsyOAycj37OVR44KaT8YkTKG9hyn9PS4rDZbOEclv71DfgV2/ob38FpoPP7M\ni2U/t853qm0aLRUQKSPu8udzty3s4C21vUPxLUfi9pRplVvDxlHAiJRRbgZvd2e65GVMVCgVHi+3\n2VSrbi5VigJGpIyoldSp9ugd/qNCqfD4HddeUfKcdKqtZTeXKkUBI1LGYH83Q7deM2/IeHlHiqEP\nXhMZAuU2mcrL32Mp35JpN2PD6l6OfurXzptwATCv147ECQwMDPjIyEjDPk+kWeJGkVqdmY26+0Dc\neRpFEgkgaqnAjxtdIolIMAoYEQlGASMiwShgRCQYBYyIBKOAEZFgFDAiEowCRkSCUcCISDAKGBEJ\nRgEjIsEoYEQkGAWMiASjgBGRYBQwIhKMAkZEglHAiEgwiQPGzNrNbMzMduceD5nZMTM7ZGZfMLPO\ncGWKSCuqpAVzF3C04PHTwM+6+9XAf5O7f7WISF6igDGzHmAd8Ej+mLt/1d3P5R4eAHrqX56ItLKk\nLZjtwL0QeRfNDwNfKfWEmW00sxEzG5mYmKiiRBFpVbEBY2Y3AqfdfTTi+fuAc8BjpZ5394fdfcDd\nB7q6umoqVkRaS5LblqwBbjKz9wJLgYvNbIe7bzCz3wRuBH7FG3mDJRFpCbEtGHff4u497t4H3A7s\nzYXLe8heNt3k7mcD1ykiLaiWeTB/AbwOeNrMxs3sr+tUk4icJyq6s6O77wP25X7+6QD1iMh5RDN5\nRSQYBYyIBKOAEZFgFDAiEowCRkSCUcCISDAKGBEJRgEjIsEoYEQkGAWMiASjgBGRYBQwIhKMAkZE\nglHAiEgwChgRCUYBIyLBKGBEJBgFjIgEo4ARkWAUMCISjAJGRIJRwIhIMAoYEQlGASMiwSQOGDNr\nN7MxM9ude3yrmR0xs1kzGwhXooi0qkpaMHcBRwsefwO4GfjXulYkIueNRAFjZj3AOuCR/DF3P+ru\nx0MVJiKtL2kLZjtwLzBb6QeY2UYzGzGzkYmJiUpfLiItLDZgzOxG4LS7j1bzAe7+sLsPuPtAV1dX\nNW8hIi0qSQtmDXCTmZ0APg9cb2Y7glYlIueF2IBx9y3u3uPufcDtwF533xC8MhFpeVXPgzGzD5jZ\nt4CfB54ysz31K0tEzgdLKjnZ3fcB+3I/fwH4Qv1LEpHzhWbyikgwChgRCUYBIyLBKGBEJBgFjIgE\no4ARkWAUMCISjAJGRIJRwIhIMAoYEQlGASMiwShgRCQYBYyIBKOAEZFgFDAiEowCRkSCUcCISDAK\nGBEJRgEjIsEoYEQkGAWMiASjgBGRYBQwIhKMAkZEgkkcMGbWbmZjZrY79/gSM3vazJ7N/XN5uDJF\npBVV0oK5Czha8Hgz8DV3fyvwtdxjEZE5iQLGzHqAdcAjBYffDzya+/lRYLC+pYlIq0vagtkO3AvM\nFhx7k7u/lPv5f4E3lXqhmW00sxEzG5mYmKi+UhFpObEBY2Y3AqfdfTTqHHd3wCOee9jdB9x9oKur\nq/pKRaTlLElwzhrgJjN7L7AUuNjMdgDfNrPL3P0lM7sMOB2yUBFpPbEtGHff4u497t4H3A7sdfcN\nwJeA38id9hvAF4NVKSItqZZ5MNuAG8zsWeBXc49FROYkuUSa4+77gH25n78L/Er9SxKR84Vm8opI\nMAoYEQlGASMiwShgRCQYBYyIBKOAEZFgFDAiEowCRkSCUcCISDAKGBEJRgEjIsEoYEQkGAWMiASj\ngBGRYBQwIhKMAkZEglHAiEgwChgRCUYBIyLBKGBEJBgFjIgEo4ARkWAUMCISTOx9kcxsKfCvwIW5\n85909wfM7Brgr4GLgBPAne7+g1oLGh7LMLTnOKcmp7i8M82mtSsY7O+u9W1FpAmStGB+BFzv7tcA\nq4D3mNlq4BFgs7uvBL4AbKq1mOGxDFt2HSYzOYUDmckptuw6zPBYpta3FpEmSHJvanf3H+YepnJ/\nHHgb2ZYNwNPALbUWM7TnOFPTM/OOTU3PMLTneK1vLSJNkKgPxszazWwcOA087e7PAEeA9+dOuRW4\notZiTk1OVXRcRBa3RAHj7jPuvgroAd5pZj8LfBj4XTMbBV4HvFrqtWa20cxGzGxkYmKi7Odc3pmu\n6LiILG4VjSK5+yTwdeA97n7M3d/t7u8AHgeei3jNw+4+4O4DXV1dZd9/09oVpFPt846lU+1sWrui\nkjJFZJGIDRgz6zKzztzPaeAG4JiZvTF3rA24n+yIUk0G+7vZevNKujvTGNDdmWbrzSs1iiTSomKH\nqYHLgEfNrJ1sID3h7rvN7C4z+73cObuAv6tHQYP93QoUkfNEbMC4+yGgv8TxPwf+PERRInJ+0Exe\nEQlGASMiwShgRCQYBYyIBKOAEZFgzN0b92FmE8ALdX7bS4Hv1Pk9VUNr1tDsz/9xquHN7l5+5iwN\nDpgQzGzE3QdUg2po9uerhoV0iSQiwShgRCSY8yFgHm52AaiGvGbX0OzPB9UwT8v3wYjI4nU+tGBE\nZJFqyYAxs51mNp77cyK32x5mljKzR83ssJkdNbMtTajhzoLj42Y2a2arGllD7rmrzew/zOxI7vex\ntJE1mFmfmU0VPFfzdh6V1lDwfK+Z/dDM/qDRNZjZOwuOHzSzDzShhhvMbDT392DUzK4PVUOxJNs1\nLDruvj7/s5l9Fvh+7uGtwIXuvtLMOoBvmtnj7n6iUTW4+2PAY7njK4Fhdx8v+SaBajCzJcAO4EPu\nftDM3gBMN7KGnOdyOyEGFVMDwOeArzSphm8AA+5+zswuAw6a2Zfd/VwDa/gO8D53P5XbjXIP0JA9\nUVoyYPLMzIDbgHwiO7As9wVLk93Gs+ZbqVRYQ6E7gM+H/PyIGt4NHHL3gwDu/t0m1NBwpWows0Hg\neeDlZtTg7mcLnl5K9u9oo2sYK3j6CJA2swvd/Ueha2nJS6QC7wK+7e7P5h4/SfYv0kvASeDP3P1M\ng2sotJ7sdqKhFdfwNsDNbI+Z/ZeZ3duEGgDekmuu/4uZvavRNZjZRcAfAp9swGeXrCFXx7VmdgQ4\nDPxOiNZLXA0FbgH+qxHhAou4BWNm/wz8RImn7nP3L+Z+voP5X+B3AjPA5cBy4N/M7J/d/X8aWEP+\ntdcCZ939G9V8do01LAF+Afg54CzwNTMbdfevNbCGl4Bed/+umb0DGDazq6q9OV+VNTwIPOTuP8z+\nT7021f59yN2F4yoz+xmyu0N+xd1faWQNuddeBfwJ2RZuY7h7S/4h+yX6NtBTcOwvyfY75B//LXBb\nI2soeO4h4ONN+j3cDjxa8PiPgE3N+D0UnLOPbF9EI38P/0b2rqMngEngDPCRJv8e9jb695A73gP8\nN7Am1GeX+tPKl0i/Chxz928VHDtJ7rrTzJYBq4FjDa4hvxH6bTSg/yWihj3ASjPryPVH/RLwzUbW\nYNnN4ttzP/8k8FagqpZktTW4+7vcvc/d+4DtwB+7+180sgYze0vuvwFm9mbgSrKB18gaOoGnyN6J\ndX/Az16glQPmdhY2A/8SuCh3vfufwN95dk/hRtYA8IvAi17lpVmtNbj798iOnPwnME72mvupRtZA\n9ndwKDdU+iTZvoeQ/WFR/y0aqVQNv0B25Gic7C2Wf9fdQ650LlXDR4CfBj5RMIz9xoA1zNFMXhEJ\nppVbMCKyyClgRCQYBYyIBKOAEZFgFDAiEowCRkSCUcCISDAKGBEJ5v8BXYTA0ZJFnuIAAAAASUVO\nRK5CYII=\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x23627676828>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "pntHistMark.plot()\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Interestingly this quickly showed an outlier in the historical marker data. In that dataset there is both coordinates and the geometry so it may be useful to use one to correct the other. Most shapefiles will only have the geometry.\n",
    "\n",
    "Coordinate Reference System (crs)\n",
    "-----------------\n",
    "Each geospatial dataset should have a coordinate system which allow you to know what datum was used. Datums are a frame of reference that relate to mathmatical models of the Earth. Knowing if two datasets are in using the same datum can be the difference between everything lining up or everything completetly off.\n",
    "\n",
    "For a fill list that can be use to decode the CRS, take a look here:http://spatialreference.org/ref/epsg/"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "{'init': 'epsg:4326'}"
      ]
     },
     "execution_count": 13,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# current coordinate reference system\n",
    "pntHistMark.crs"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Coordinate Transformation\n",
    "\n",
    "Sometime a CRS needs to be assigned or be changed.  Here is a simple way to do just that that, the transformation of the points happen automatically."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "{'init': 'epsg:4327', 'no_defs': True}"
      ]
     },
     "execution_count": 14,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "#create a new dataset to modify\n",
    "newdata = pntHistMark\n",
    "\n",
    "from fiona.crs import from_epsg\n",
    "\n",
    "# Set the GeoDataFrame's coordinate system to WGS84\n",
    "newdata.crs = from_epsg(4327)\n",
    "\n",
    "# verify crs definition change\n",
    "newdata.crs"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Table joins\n",
    "---------------------\n",
    "Table joins are something that is done in data science in general and the same methods used in Pandas can be applied here. In this case the two datasets don't have a common index so in this example I am going to merge based on the OBJECTID but ideally this would be a more relevant field."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style>\n",
       "    .dataframe thead tr:only-child th {\n",
       "        text-align: right;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: left;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>OBJECTID</th>\n",
       "      <th>name</th>\n",
       "      <th>descriptio</th>\n",
       "      <th>website</th>\n",
       "      <th>town</th>\n",
       "      <th>county</th>\n",
       "      <th>address</th>\n",
       "      <th>video_url</th>\n",
       "      <th>coordinate</th>\n",
       "      <th>coordina_1</th>\n",
       "      <th>...</th>\n",
       "      <th>YearCast</th>\n",
       "      <th>FiscalYear</th>\n",
       "      <th>geometry_x</th>\n",
       "      <th>RPC</th>\n",
       "      <th>INITIALS</th>\n",
       "      <th>SHORTNAME</th>\n",
       "      <th>LONGNAME</th>\n",
       "      <th>ShapeSTAre</th>\n",
       "      <th>ShapeSTLen</th>\n",
       "      <th>geometry_y</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>1</td>\n",
       "      <td>\"Snowflake\" Bentley -- Jericho's world famous ...</td>\n",
       "      <td>For fifty years Wilson A. Bentley, a farmer an...</td>\n",
       "      <td>None</td>\n",
       "      <td>Jericho</td>\n",
       "      <td>CHITTENDEN</td>\n",
       "      <td>Jericho Center Cir., in village</td>\n",
       "      <td>None</td>\n",
       "      <td>-72.973052</td>\n",
       "      <td>44.469608</td>\n",
       "      <td>...</td>\n",
       "      <td>2003</td>\n",
       "      <td>2016.0</td>\n",
       "      <td>POINT (-72.97305209009137 44.46960766018835)</td>\n",
       "      <td>CC</td>\n",
       "      <td>CCRPC</td>\n",
       "      <td>Chittenden County RPC</td>\n",
       "      <td>Chittenden County Regional Planning Commission</td>\n",
       "      <td>1.609237e+09</td>\n",
       "      <td>201814.709240</td>\n",
       "      <td>POLYGON ((452203.9672731007 241948.2255608736,...</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>2</td>\n",
       "      <td>\"Vermont's First College\"</td>\n",
       "      <td>On this site Castleton State College, Vermonts...</td>\n",
       "      <td>None</td>\n",
       "      <td>Castleton</td>\n",
       "      <td>RUTLAND</td>\n",
       "      <td>VT Route 4A in front of the town library</td>\n",
       "      <td>None</td>\n",
       "      <td>-73.178299</td>\n",
       "      <td>43.611718</td>\n",
       "      <td>...</td>\n",
       "      <td>1986</td>\n",
       "      <td>NaN</td>\n",
       "      <td>POINT (-73.17829915991881 43.61171830002038)</td>\n",
       "      <td>SW</td>\n",
       "      <td>SWCRPC</td>\n",
       "      <td>Southern Windsor County RPC</td>\n",
       "      <td>Southern Windsor County Regional Planning Comm...</td>\n",
       "      <td>8.908414e+08</td>\n",
       "      <td>128848.691548</td>\n",
       "      <td>POLYGON ((497499.9973051278 114802.4954725726,...</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>3</td>\n",
       "      <td>1st Normal School -- Pioneer in Teacher Training</td>\n",
       "      <td>The first recognized school in America for the...</td>\n",
       "      <td>None</td>\n",
       "      <td>Concord</td>\n",
       "      <td>ESSEX</td>\n",
       "      <td>U.S. Route 2, at Shadow Lake Road</td>\n",
       "      <td>None</td>\n",
       "      <td>-71.889163</td>\n",
       "      <td>44.428617</td>\n",
       "      <td>...</td>\n",
       "      <td>1954</td>\n",
       "      <td>NaN</td>\n",
       "      <td>POINT (-71.88916283008687 44.42861686006097)</td>\n",
       "      <td>BC</td>\n",
       "      <td>BCRC</td>\n",
       "      <td>Bennington County RC</td>\n",
       "      <td>Bennington County Regional Commission</td>\n",
       "      <td>1.488553e+09</td>\n",
       "      <td>203533.897896</td>\n",
       "      <td>POLYGON ((451638.0272736274 89835.76545525473,...</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>4</td>\n",
       "      <td>2nd Vermont Volunteer Infantry Regiment</td>\n",
       "      <td>On June 24, 1861, the first Vermont soldiers w...</td>\n",
       "      <td>None</td>\n",
       "      <td>Burlington</td>\n",
       "      <td>CHITTENDEN</td>\n",
       "      <td>North Avenue in Battery Park</td>\n",
       "      <td>None</td>\n",
       "      <td>-73.220763</td>\n",
       "      <td>44.481325</td>\n",
       "      <td>...</td>\n",
       "      <td>2011</td>\n",
       "      <td>NaN</td>\n",
       "      <td>POINT (-73.22076252968733 44.48132500982513)</td>\n",
       "      <td>NW</td>\n",
       "      <td>NWRPC</td>\n",
       "      <td>Northwest RPC</td>\n",
       "      <td>Northwest Regional Planning Commission</td>\n",
       "      <td>2.285635e+09</td>\n",
       "      <td>242793.916924</td>\n",
       "      <td>POLYGON ((460860.0944790148 279521.6205874127,...</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <td>5</td>\n",
       "      <td>A.W.Gray &amp; Sons</td>\n",
       "      <td>Near this site Albert W. Gray manufactured his...</td>\n",
       "      <td>None</td>\n",
       "      <td>Middletown Springs</td>\n",
       "      <td>RUTLAND</td>\n",
       "      <td>Route 133A and Montvert Rd.</td>\n",
       "      <td>None</td>\n",
       "      <td>-73.118337</td>\n",
       "      <td>43.482350</td>\n",
       "      <td>...</td>\n",
       "      <td>1998</td>\n",
       "      <td>NaN</td>\n",
       "      <td>POINT (-73.11833662009371 43.482349910201)</td>\n",
       "      <td>LC</td>\n",
       "      <td>LCPC</td>\n",
       "      <td>Lamoille County PC</td>\n",
       "      <td>Lamoille County Planning Commission</td>\n",
       "      <td>1.207382e+09</td>\n",
       "      <td>166207.247550</td>\n",
       "      <td>POLYGON ((491679.967300891 254475.715569835, 4...</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "<p>5 rows × 21 columns</p>\n",
       "</div>"
      ],
      "text/plain": [
       "   OBJECTID                                               name  \\\n",
       "0         1  \"Snowflake\" Bentley -- Jericho's world famous ...   \n",
       "1         2                          \"Vermont's First College\"   \n",
       "2         3   1st Normal School -- Pioneer in Teacher Training   \n",
       "3         4            2nd Vermont Volunteer Infantry Regiment   \n",
       "4         5                                    A.W.Gray & Sons   \n",
       "\n",
       "                                          descriptio website  \\\n",
       "0  For fifty years Wilson A. Bentley, a farmer an...    None   \n",
       "1  On this site Castleton State College, Vermonts...    None   \n",
       "2  The first recognized school in America for the...    None   \n",
       "3  On June 24, 1861, the first Vermont soldiers w...    None   \n",
       "4  Near this site Albert W. Gray manufactured his...    None   \n",
       "\n",
       "                 town      county                                   address  \\\n",
       "0             Jericho  CHITTENDEN           Jericho Center Cir., in village   \n",
       "1           Castleton     RUTLAND  VT Route 4A in front of the town library   \n",
       "2             Concord       ESSEX         U.S. Route 2, at Shadow Lake Road   \n",
       "3          Burlington  CHITTENDEN              North Avenue in Battery Park   \n",
       "4  Middletown Springs     RUTLAND               Route 133A and Montvert Rd.   \n",
       "\n",
       "  video_url  coordinate  coordina_1  \\\n",
       "0      None  -72.973052   44.469608   \n",
       "1      None  -73.178299   43.611718   \n",
       "2      None  -71.889163   44.428617   \n",
       "3      None  -73.220763   44.481325   \n",
       "4      None  -73.118337   43.482350   \n",
       "\n",
       "                         ...                         YearCast  FiscalYear  \\\n",
       "0                        ...                             2003      2016.0   \n",
       "1                        ...                             1986         NaN   \n",
       "2                        ...                             1954         NaN   \n",
       "3                        ...                             2011         NaN   \n",
       "4                        ...                             1998         NaN   \n",
       "\n",
       "                                     geometry_x RPC INITIALS  \\\n",
       "0  POINT (-72.97305209009137 44.46960766018835)  CC    CCRPC   \n",
       "1  POINT (-73.17829915991881 43.61171830002038)  SW   SWCRPC   \n",
       "2  POINT (-71.88916283008687 44.42861686006097)  BC     BCRC   \n",
       "3  POINT (-73.22076252968733 44.48132500982513)  NW    NWRPC   \n",
       "4    POINT (-73.11833662009371 43.482349910201)  LC     LCPC   \n",
       "\n",
       "                     SHORTNAME  \\\n",
       "0        Chittenden County RPC   \n",
       "1  Southern Windsor County RPC   \n",
       "2         Bennington County RC   \n",
       "3                Northwest RPC   \n",
       "4           Lamoille County PC   \n",
       "\n",
       "                                            LONGNAME    ShapeSTAre  \\\n",
       "0     Chittenden County Regional Planning Commission  1.609237e+09   \n",
       "1  Southern Windsor County Regional Planning Comm...  8.908414e+08   \n",
       "2              Bennington County Regional Commission  1.488553e+09   \n",
       "3             Northwest Regional Planning Commission  2.285635e+09   \n",
       "4                Lamoille County Planning Commission  1.207382e+09   \n",
       "\n",
       "      ShapeSTLen                                         geometry_y  \n",
       "0  201814.709240  POLYGON ((452203.9672731007 241948.2255608736,...  \n",
       "1  128848.691548  POLYGON ((497499.9973051278 114802.4954725726,...  \n",
       "2  203533.897896  POLYGON ((451638.0272736274 89835.76545525473,...  \n",
       "3  242793.916924  POLYGON ((460860.0944790148 279521.6205874127,...  \n",
       "4  166207.247550  POLYGON ((491679.967300891 254475.715569835, 4...  \n",
       "\n",
       "[5 rows x 21 columns]"
      ]
     },
     "execution_count": 15,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "join = pntHistMark.merge(vctCntyBound, on='OBJECTID')\n",
    "join.head(5)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Point in Polygon & Intersect\n",
    "\n",
    "Often GIS will have a mx of points, lines and shapes. Determining if a point, i.e. a road marker, is within a shape, i.e. town boundries, can be quite useful and there are some built in functions in Shapely to do just that. Performing a spatial join (this will be a seperate notebook) between two spatial datasets is one of the most typical applications where Point in Polygon (PIP) query is used\n",
    "\n",
    "There are basically two ways of conducting PIP in Shapely:\n",
    "* using a function called .within() that checks if a point is within a polygon\n",
    "* using a function called .contains() that checks if a polygon contains a point\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "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)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "True\n",
      "False\n"
     ]
    }
   ],
   "source": [
    "# Is p1 in the polygon\n",
    "print(p1.within(poly))\n",
    "# Is p1 in the polygon\n",
    "print(p2.within(poly))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "metadata": {
    "scrolled": true
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "True\n",
      "False\n"
     ]
    }
   ],
   "source": [
    "# Does polygon contain p1?\n",
    "print(poly.contains(p1))\n",
    "# Does polygon contain p2?\n",
    "print(poly.contains(p2))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Intersect\n",
    "\n",
    "Another typical geospatial operation is to see if a geometry intersect or touches another one. The difference between these two is that:\n",
    "* if objects intersect, the boundary and interior of an object needs to intersect in any way with those of the other.\n",
    "* if an object touches the other one, it is only necessary to have (at least) a single point of their boundaries in common but their interiors shoud NOT intersect."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "metadata": {
    "collapsed": 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": "code",
   "execution_count": 20,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "True"
      ]
     },
     "execution_count": 20,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "line_a.intersects(line_b)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "True"
      ]
     },
     "execution_count": 21,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "line_a.touches(line_b)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 22,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "True"
      ]
     },
     "execution_count": 22,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Does the line intersect with itself?\n",
    "line_a.intersects(line_a)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Summary\n",
    "\n",
    "This is just an introduction to Geopandas, there are a lot of applications and libraries that can be explored. I again will strongly recommend this link for a more comprehensive walkthrough: https://automating-gis-processes.github.io/2016/Lesson1-Intro-Python-GIS.html. The power in using Python and R is the scalability with large datasets and by their nature, most of GIS is made of big data. Explore and share!"
   ]
  }
 ],
 "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.1"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}