Creating a GEOJSON Linestring from a polyline
In a previous post, we had set out to take a pair of latitudes and longitudes and eventually create a map of directions between the two points. We had discussed the process of geocoding the locations and using APIs to produce a polyline which represented the directions between the two points. In this point we will proceed to translate the polyline and produce GEOJSON code representing the directions.
Recalling from the previous post, we were looking to map the walking directions from the Bucktown library to the Map Room in Chicago. The polyline produced in the previous post had the value of “{~x~FpadvOTa@??Ky@??w`@X”.
The Python code below will decode the polyline into coordinates and produce a GEOJSON code snippet representing the directions. Normally I tend to deploy the resulting GEOJSON to the web, ofttimes using D3.js and Leaflet within a webpage. For here, I will merely show the GEOJSON as it is shown at geojsonlint.com
Here is the python code which creates the GEOJSON
import json
from shapely.geometry import Point,mapping
#### Define code function ... copied from http://seewah.blogspot.com/2009/11/gpolyline-decoding-in-python.html
def decode_line(encoded):
"""Decodes a polyline that was encoded using the Google Maps method.
See http://code.google.com/apis/maps/documentation/polylinealgorithm.html
This is a straightforward Python port of Mark McClure's JavaScript polyline decoder
(http://facstaff.unca.edu/mcmcclur/GoogleMaps/EncodePolyline/decode.js)
and Peter Chng's PHP polyline decode
(http://unitstep.net/blog/2008/08/02/decoding-google-maps-encoded-polylines-using-php/)
"""
encoded_len = len(encoded)
index = 0
array = []
lat = 0
lng = 0
while index < encoded_len:
b = 0
shift = 0
result = 0
while True:
b = ord(encoded[index]) - 63
index = index + 1
result |= (b & 0x1f) << shift
shift += 5
if b < 0x20: break dlat = ~(result >> 1) if result & 1 else result >> 1
lat += dlat
shift = 0
result = 0
while True:
b = ord(encoded[index]) - 63
index = index + 1
result |= (b & 0x1f) << shift
shift += 5
if b < 0x20: break dlng = ~(result >> 1) if result & 1 else result >> 1
lng += dlng
array.append((lat * 1e-5, lng * 1e-5))
return array
# This is our polyline between the two locations
test_polyline = '{~x~FpadvOTa@??Ky@??w`@X'
gjson_dict={}
gjson_dict["type"]= "FeatureCollection"
feat_list = []
type_dict = {}
ls_dict = {}
prop_dict = {}
coord_list = []
type_dict["type"]= "Feature"
ls_dict["type"]="LineString"
coords = decode_line(test_polyline) # Translate the polyline
for c in coords:
coord_list.append([c[1],c[0]]) # GEOJSON looks for long,lat so reverse order
ls_dict["coordinates"]=coord_list
type_dict["geometry"]=ls_dict
prop_dict["Property1"]= 'Map Room Directions'
prop_dict["Property2"]= 'Walking from Library'
type_dict["properties"]=prop_dict
feat_list.append(type_dict)
gjson_dict["features"] = feat_list
print gjson_dict
with open('poly_for_blog_geojson.json', 'w') as outfile:
json.dump(gjson_dict, outfile, sort_keys = True, indent = 4,
ensure_ascii=False)