Finding the Optimal Location for a New Park

This example notebook will show how to find the next potential location for a new park in San Fransisco. To accomplish this, three factors will be taken into consideration when deciding on a possible spot: existing parks, schools, and Bay Area Regional Transit (BART) stops. By calculating Euclidean Distance for these three factors and then weighing them together, we will be able to produce a visual representation of where is and is not a good location for a new park.

Importing the Libraries


In [ ]:
import geopyspark as gps
import fiona

from pyspark import SparkContext, StorageLevel
from shapely.geometry import MultiPoint, MultiPolygon, shape

import folium

Setup the SparkContext


In [ ]:
conf = gps.geopyspark_conf(appName="park-siting", master="local[*]")
sc = SparkContext.getOrCreate(conf=conf)

Set map display parameters


In [ ]:
center = [37.8, -122.2]
zoom_start = 9.5

Download the Geometries as GeoJsons


In [ ]:
!curl -o /tmp/bart.geojson https://s3.amazonaws.com/geopyspark-demo/bayarea/bart.geojson
!curl -o /tmp/school.geojson https://s3.amazonaws.com/geopyspark-demo/bayarea/school.geojson
!curl -o /tmp/parks.geojson https://s3.amazonaws.com/geopyspark-demo/bayarea/parks.geojson

Read in the GeoJsons as Shapely Geometries


In [ ]:
with fiona.open("/tmp/bart.geojson") as source:
    bart_crs = source.crs['init']
    bart = MultiPoint([shape(f['geometry']) for f in source])

with fiona.open("/tmp/school.geojson") as source:
    schools_crs = source.crs['init']
    schools = MultiPoint([shape(f['geometry']) for f in source])

with fiona.open("/tmp/parks.geojson") as source:
    parks_crs = source.crs['init']
    parks = MultiPolygon([shape(f['geometry']) for f in source])

Calculate Euclidean Distance for Each Geometry

Three new TiledRasterLayers will be produced from the Euclidean Distance calculations for each geometry. All resulting layers will have a zoom_level of 12.


In [ ]:
bart_layer = gps.euclidean_distance(geometry=bart,
                                    source_crs=bart_crs,
                                    zoom=12)

schools_layer = gps.euclidean_distance(geometry=schools,
                                       source_crs=schools_crs,
                                       zoom=12)

parks_layer = gps.euclidean_distance(geometry=parks,
                                     source_crs=parks_crs,
                                     zoom=12)

# Persists each layer to memory and disk
bart_layer.persist(StorageLevel.MEMORY_AND_DISK)
schools_layer.persist(StorageLevel.MEMORY_AND_DISK)
parks_layer.persist(StorageLevel.MEMORY_AND_DISK)

Weighing the Layers Together


In [ ]:
weighted_layer = -1 * bart_layer - schools_layer + 3 * parks_layer

# Persists the weighted layer to memory and disk
weighted_layer.persist(StorageLevel.MEMORY_AND_DISK)

Reprojecting, Pyramiding, and Calculating the Histogram


In [ ]:
# The following code may take awhile to complete
reprojected = weighted_layer.tile_to_layout(layout=gps.GlobalLayout(),
                                            target_crs="EPSG:3857")
pyramid = reprojected.pyramid(resample_method=gps.ResampleMethod.AVERAGE)
histogram = pyramid.get_histogram()

Creating the ColorMap

The below code creates a ColorMap instance using the Histogram from pyramid for its breaks. For the color, the matplotlib color palette, viridus will be used.


In [ ]:
color_map = gps.ColorMap.build(breaks=histogram,
                               colors='viridis')

Running the Server


In [ ]:
tms = gps.TMS.build(source=pyramid,
                    display=color_map)

tms.bind('0.0.0.0')

In [ ]:
m = folium.Map(tiles='OpenStreetMap', location=center, zoom_start=zoom_start)
folium.TileLayer(tiles=tms.url_pattern, overlay=True, attr='GeoPySpark tiles').add_to(m)
folium.GeoJson(data='/tmp/bart.geojson', name='BART stops').add_to(m)
folium.GeoJson(data='/tmp/parks.geojson', name='Parks').add_to(m)
folium.LayerControl().add_to(m)
m

Cleaning up


In [ ]:
tms.unbind()