Elijah Robison | GIS Blog

A scrapbook of GIS tricks, with emphasis on FOSS4G.

PostGREsql/PostGIS Implementation of Google’s Encoded Polyline Algorithm

with 8 comments

[Edit 30 Jan, 2014]

I added an additional PostGREsql method to perform Polygon encoding by concatenating polygon geometries (delimiter: †) and their inner rings (delimiter: ‡) together into one massive encoded block of ring features. I also provided an example JavaScript method demonstrating how to bring the amalgamated polygon feature encodings into your Google Map.

By some uncanny twist of the fates, I’ve elected to use, had to use, and/or been asked to develop applications that use Google Maps ASCII Encoded Polyline expressions. In previous encounters, I’ve used a PHP class to handle the encoding task, and most recently I wrote a Python method to decode these expressions so that I could return a 3rd-party’s encoded geometries to WKT and import them into a spatially aware database.

So far so good.

However one thing has always bugged me about using the PHP solution–I don’t like using a piece of middleware to handle what I consider to be a responsibility of the data layer. Mark McClure’s page, which is basically the seminal authority on this topic, provides external links to implementations for Perl, Ruby, PHP (note: I prefer the PHP class linked, above), Java, and Mathematica. Also, by searching Stack Overflow, you can find implementations of the algorithm in both C# and VB.Net. But for all my efforts searching, I could never dredge up an implementation for either MySQL or PostGREsql/PostGIS. Bummer.

Looking up, it seems version 2.2 of PostGIS might include some built-in Google encoding conversion methods. While this is cool enough for a hat tip, unfortunately, it’s too inconvenient to wait that long, and even then, there’s no guarantee the implementation will work the way I expect with complex Polygon geometries; for instance, maybe it will encode only the exterior ring of Polygons, ignoring MultiPolygons completely, etc. For that matter, it’s equally possible there could be some bugs. So with this said, and even though the previously-mentioned PHP implementation does the job, my boss was cool-enough to let me take a crack at implementing the algorithm as a PostGREsql/PostGIS function, and then share the results with the world. Since some initial testing confirms my PostGIS implementation works, I’ll just post the code parts and hope others find it useful.

For what it’s worth, if anyone finds a bug or has recommendations for improvements, please don’t hesitate to drop me a line.

 

Sample query calling the first encoding function on the EXTERIOR RING of Polygon geometries:
(Also works on single-part LINESTRING features.)

/************************************************************************
 * Note that the encoding method can accept a LINESTRING only, which
 * is the geometry type used to represent the ring parts of a Polygon.
 * To help understand this, and why, please see the trailing discussion
 * section, which elaborates on this situation further.
 ************************************************************************/
SELECT
  GoogleEncodeLine(ST_ExteriorRing(wkb_geometry)) as Google
  FROM polygons_wgs84
  WHERE ST_GeometryType(wkb_geometry) = 'ST_Polygon'
  LIMIT 10 ;

 

[Added 30 Jan, 2014]

Sample query calling the second encoding function on Polygon and MultiPolygon geometries:
(Preserves multi-part polygons and their inner-ring parts, a.k.a. “holes”.)

/************************************************************************
 * This encoding method will accept Polygon and MultiPolygon geom types.
 * The output returned is an amalgamation of Polyline encodings, where
 * individual geometries and their interior rings are concatenated
 * together using string delimiters, †, and ‡, respectively.
 ************************************************************************/
SELECT
  GoogleEncodePolygon(wkb_geometry) as GooglePolygon
  FROM polygons_wgs84
  LIMIT 10 ;

 

Implementation functions to execute/save in your PostGREsql instance:

[Added 30 Jan, 2014]

/*************************************************************
 * Pass in either a Polygon or MultiPolygon geometry. Returns
 * an array of ASCII-encoded Polygon feature parts, including
 * multi-part geometries and their interior rings.
 ************************************************************/
CREATE OR REPLACE FUNCTION GoogleEncodePolygon
(
  g1 GEOMETRY
)
RETURNS TEXT AS $$
DECLARE
 ng INT;        -- Store number of Geometries in the Polygon.
 g INT;         -- Counter for the current geometry number during outer loop.
 g2 GEOMETRY;   -- Current geometry feature isolated by the outer loop.
 nr INT;        -- Store number of internal ring parts in the Polygon.
 r INT;         -- Counter for the current inner-ring part.
 r1 GEOMETRY;   -- Exterior ring part isolated BEFORE the inner loop.
 r2 GEOMETRY;   -- Inner-ring part isolated within the inner loop.
 gEncoded TEXT; -- Completed Google Encoding.
BEGIN
 gEncoded = '';
 ng = ST_NumGeometries(g1);
 g = 1;
 FOR g IN 1..ng BY 1 LOOP
     g2 = ST_GeometryN(g1, g);
     if g > 1 then gEncoded = gEncoded || chr(8224); END IF;
     -- Get ExteriorRing now; if there are any holes, get them later in the loop..
     r1 = ST_ExteriorRing(g2);
     gEncoded = gEncoded || GoogleEncodeLine(r1);
     nr = ST_NRings(g2);
     if nr > 1 then
       -- One (1) is because interior rings is one-based.
       -- And nr-1 is because ring count includes the boundary.
       FOR r IN 1..(nr-1) BY 1 LOOP
         r2 = ST_InteriorRingN(g2, r);
         gEncoded = gEncoded || chr(8225) || GoogleEncodeLine(r2);
       END LOOP;
     END IF;
 END LOOP;
 RETURN gEncoded;
End
$$ LANGUAGE plpgsql;

 

/*************************************************************
 * First of two methods. Pass in a geometry (LINESTRING only).
 * Returns ASCII-encoded point array for use in Google Maps.
 ************************************************************/
CREATE OR REPLACE FUNCTION GoogleEncodeLine
(
  g GEOMETRY
)
RETURNS TEXT AS $$
DECLARE
  pt1 GEOMETRY;
  pt2 GEOMETRY;
  p INT; np INT;
  deltaX INT;
  deltaY INT;
  enX VARCHAR(255);
  enY VARCHAR(255);
  gEncoded TEXT;
BEGIN
  gEncoded = '';
  np = ST_NPoints(g);

  IF np > 3 THEN
    g = ST_SimplifyPreserveTopology(g, 0.00001);
    np = ST_NPoints(g);
  END IF;

  pt1 = ST_SetSRID(ST_MakePoint(0, 0),4326);

  FOR p IN 1..np BY 1 LOOP
    pt2 = ST_PointN(g, p);
    deltaX = (floor(ST_X(pt2)*1e5)-floor(ST_X(pt1)*1e5))::INT;
    deltaY = (floor(ST_Y(pt2)*1e5)-floor(ST_Y(pt1)*1e5))::INT;
    enX = GoogleEncodeSignedInteger(deltaX);
    enY = GoogleEncodeSignedInteger(deltaY);
    gEncoded = gEncoded || enY || enX;

    pt1 = ST_SetSRID(ST_MakePoint(ST_X(pt2), ST_Y(pt2)),4326);
  END LOOP;
RETURN gEncoded;
End
$$ LANGUAGE plpgsql;

 

/**************************************************************
 * Second of two methods. Accepts a signed integer (LON or LAT
 * by 1e5) and returns an ASCII-encoded coordinate expression.
 *************************************************************/
CREATE OR REPLACE FUNCTION GoogleEncodeSignedInteger(c INT)
RETURNS VARCHAR(255) AS $$
DECLARE
  e VARCHAR(255);
  s BIT(32);
  b BIT(6);
  n INT;
BEGIN
 e = '';
 s = (c::BIT(32))<<1;

 IF s::INT < 0 THEN
   s = ~s;
   END IF;

 WHILE s::INT >= B'100000'::INT LOOP
   b = B'100000' | (('0'||substring(s, 28, 5))::BIT(6));
   n = b::INT + 63;
   e = e || chr(n);
   s = s >> 5;
 END LOOP;
 e = e || chr(s::INT+63);

RETURN e;
End
$$ LANGUAGE plpgsql;

 

[Added 30 Jan, 2014]

JavaScript method demonstrating how to add Polygon encodings to a Google Map object:
(This client implementation works for either the single or the multi-part polygons.)

/*************************************************************
 * JavaScript! Pass-in an encoded text block created by either
 * PostGIS method, GoogleEncodePolygon() or GoogleEncodeLine(),
 * and render it in your Google Map object. If you don't want
 * the map to zoom to each rendering, just remove the "bounds"
 * variable and any references to it.
 ************************************************************/
function renderEncoded(encoded_path)
{
   var bounds = new google.maps.LatLngBounds();
   var $encodedGeoms = encoded_path.split("†");
   for (var i=0; i<$encodedGeoms.length; i++)
   {
       var encodedGeom = $encodedGeoms[i];
       var $encodedRings = encodedGeom.split("‡");
       var polyPaths = [];
       for (var j=0; j<$encodedRings.length; j++)
       {
           var ptarray = google.maps.geometry.encoding.decodePath($encodedRings[j]);
           polyPaths.push(ptarray);
       }
       var polygonObject = new google.maps.Polygon(
       {
         paths: polyPaths,
         strokeColor: '#890000',
         strokeOpacity: 1.0,
         strokeWeight: 2
       });
       polygonObject.setMap(map);
       polygonObject.getPath().forEach(function(e)
       {
           bounds.extend(e);
       });
   }
   map.fitBounds(bounds);
}

 

And some additional discussion..

There are two “gotchas” when it comes to implementing the encoding algorithm with respect to Polygons:

1) Polygons, as geometries, can be composed of many rings. The outer ring is considered to be the boundary, and various inner rings are often called “holes”. So this is a specified, understood, and accepted built-in many-to-one relationship between polygons and their internal ring geometries.

And 2) It’s not rare to find polygon tables containing both Polygon and MultiPolygon data types. I think this happens because ESRI allows it, and so in an effort to play well with others, other GIS systems have accommodated it. At least, I know this is true for MySQL and PostGIS.

Here’s why this makes trouble–Google’s encoding algorithm is only intended to represent individual point arrays as a singular geometry. Basically, as long as your first point equals your last point, it’s considered to be a closed geometry, and you can add it and render it in a Google Map as a polygon. The algorithm itself isn’t designed to represent nested arrays, which would be necessary to render either a Polygon with “holes” or a MultiPolygon, which could potentially define many polygons with holes of their own! As such, I suspect there could be considerable disagreement as to how a Polygon-to-Google-Encoded method should actually handle Polygons..

The only solutions I can imagine for this problem would require “faking” a one-to-many relationship by perhaps delimiting together several encodings to account for MultiPolygons and/or single feature Polygons with interior rings. But this starts to get weird. So to keep things somewhat simple for the sake of the post, I chose to stay true to the algorithm’s intent and return a single encoded geometry expression from my method. And the sample query demonstrates this by calling the method against the outermost ring (i.e. the boundary) of a Polygon geometry type, which PostGREsql regards as a LineString, anyway.

[Added 30 Jan, 2014]

Since I wanted to handle the more complex geometries, I wrote the wrapper method GoogleEncodePolygon(), to first iterate over ST_NumGeometries() and gain access to any multi-part features, then second, iterate over ST_NRings() using ST_InteriorRingN()–you could also do this using ST_DumpRings()–and gain access to any interior rings of the Polygon geometries, themselves. Then, for each ring part, I call GoogleEncodeLine(), and concatenate together all those expressions into one massive block of “compound” expressions. I chose to delimit each geometry encoding using an extra-special character that would never be used by Google’s algorithm; for example chr(8224), which corresponds to “†”. I then further delimit the internal ring parts using another special character, chr(8225), which corresponds to “‡”, and return all these concatenated together as a compound encoding expression. Then, on the client-side (a JavaScript example is provided above), I merely split the compound expression against my delimiters, loop over the various expressions, and add them to the map individually. Note if you are attaching attributes to your features, you’ll need to remember to include them explicitly to each unique Polygon added to your map.

Written by elrobis

January 27th, 2014 at 12:20 pm

8 Responses to 'PostGREsql/PostGIS Implementation of Google’s Encoded Polyline Algorithm'

Subscribe to comments with RSS

  1. Very useful article, I tried you PlPgsql scripts but my setup (Pg.9 + PostGIS 1.5) give on the “GoogleEncodeSignedInteger”, row:10:
    s = (c::BIT(32))<= B'100000'::INT LOOP
    I see you forgot a semicolon after the int cast. I fixed it but still is giving me errors, when run in a query. Something like this:

    "ERROR: operator does not exist: bit <= integer
    LINE 1: SELECT (c::BIT(32))<= B'100000'::INT
    ^
    HINT: No operator matches the given name and argument type(s). You might need to add explicit type casts."

    Can you help me?
    Thank you!

    Charlieboy

    10 Feb 14 at 10:27 am

  2. Hey there, Charlieboy, I appreciate your feedback. Looking at the code in the blog post vs. a regular ASCII text file I kept, and alongside the code that was ultimately stored in my database instance, something’s definitely amiss. Fortunately, the code in my ASCII file and my database are agree with each other. If you’re willing, I’d like to email you the code blocks, perhaps as the text file itself, and you can tell me if those work for you. If they do, I’ll fix-up the blog. ….I’m not yet sure how to explain this. Either the blog “ate” some of the code when I saved it, or I got in a hurry and mistakenly chopped some lines out of it. For example, it’s totally missing the “WHILE” condition of that loop. My vote is human error. Please hit me back with a comment reply if you don’t receive an email from me shortly. Best, elrobis

    elrobis

    10 Feb 14 at 12:41 pm

  3. Hi elrobis the code you mail me is soooo good!
    Now everything works as expected.

    Thank you very much.
    Ciao, grazie!

    Charlieboy

    10 Feb 14 at 1:26 pm

  4. Thanks for the great function. It works for me with basic polylines but when I try and encode a function such as
    ‘MULTIPOLYGON(((0 0,4 0,4 4,0 4,0 0),(1 1,2 1,2 2,1 2,1 1)), ((-1 -1,-1 -2,-2 -2,-2 -1,-1 -1)))’::geometry

    I get an invalid geometry
    ???_glW_glW??~flW~flW?‡_ibE_ibE?_ibE_ibE??~hbE~hbE?†~hbE~hbE~hbE??~hbE_ibE??_ibE

    What I expect:
    ???_glW_glW??~flW~flW?,~hbE~hbE~hbE??~hbE_ibE??_ibE

    As tested against https://developers.google.com/maps/documentation/utilities/polylineutility

    Keith

    1 Aug 14 at 1:38 pm

  5. Keith, hi–thanks for the kudos. Based on the geometry encodings you provided, I think my implementation is working as it should. I’ll try to explain. (And perhaps I’ll find out later that I am wrong!?)

    First, Google’s algorithm doesn’t actually support MultiPolygon geometries. If I go to the URL you provided (https://developers.google.com/maps/documentation/utilities/polylineutility) and plug in your expected geometry: ???_glW_glW??~flW~flW?,~hbE~hbE~hbE??~hbE_ibE??_ibE

    It gives me this ..which doesn’t look right.

    On the other hand, I wanted to accommodate MultiPolygon geometries. So I wrote these routines to concatenate together multiple outer and inner ring parts using the symbols ‡, and †, which identify inner ring parts (‡) and subsequent/new geometries (†), respectively. In order to actually use this contrived way of handling polygons with holes and MultiPolygons, I also had to write a JavaScript function to iterate over the encodings and add them to the map. If you double-check the post, the JavaScript method is the last code block I provided.

    When I glance at the expression my routines created for your MultiPolygon, I can see it representing both an initial polygon with an inner ring, but also a second polygon without any inner rings, which is pretty much what I expected. So I decided to check it out. Fortunately I still had a test page floating around, so I uploaded it to a public location and gave the expression my routines provided you a whirl..

    http://www.elrobis.com/site/apps/polyline/google_maps_encoded_line_test.php

    And when I paste the encoding you provided and click “render”, it gives me this ..which looks about like what I’d expect.

    In summary, Google didn’t provide a way to encode polygons as much as they provided a way to encode individual polylines and/or rings. (Which are basically the same things, only a ring has a final point that equals its initial point.) But I, personally, needed to encode features in a real world polygon dataset, which included inner rings/holes, as well as MultiPolygons, which could include any variety of crazy polygons.

    It’s possible, maybe likely, that alot of people would protest my solutions to the multigeometry/rings problems, particularly because the encodings created by my routines are subject to client-side dependency upon my JavaScript routine. And that would be a fair complaint. However, Google (that I am aware of) doesn’t handle the more complex shapes, either in their routine or in their client. So my solution was more or less a best effort. Anyway, I’m pretty sure it’s working as advertised.

    On another note…… where did you get your expected geometry encoding from, and how did you test that it was a legit encoding, rendering in the map as you expected? Because the expected encoding you provided certainly didn’t look the way I expected it would when I checked it in the test utility you provided.

    Hopefully this provides some helpful insight.

    Best, Elijah

    elrobis

    1 Aug 14 at 2:58 pm

  6. Elrobis,

    Thank you very much. That was my fault for not looking closely at your javascript function. Re-reading it, it makes sense. I apologize for the snap judgement.

    I’ve been tasked with migrating our database away from its plpython dependency. The current tests and functions are based on the old cgpolyencode v0.1, so they themselves might be out of date. However, I’m trying to make the update from the encode exactly match the original python version. I think I have a working version based off your code here: https://gist.github.com/keithhackbarth/44d44be1507de484dbb0

    Anyways, thanks again for making this public. You saved me hours of time. I know this function is slated for a future PostGis release. Keeping my fingers crossed!

    Keith

    1 Aug 14 at 5:44 pm

  7. Thank you for providing such a useful solutions. I appreciate very much the assumption that this conversion belongs to the data layer. Anyway I couldn’t use them on my pg 8.4/PostGis 1.5 database. This is the pl/pgsql error I got.
    Sorry for this stupid question, I’m sure a better understanding of pl/pgsql would let me fix it, but …

    CREATE FUNCTION
    psql:cartometric_gencodepoly.sql:126: ERROR: “$1” is declared CONSTANT
    CONTEXT: compilation of PL/pgSQL function “googleencodeline” near line 15
    psql:originale.sql:159: ERROR: syntax error at or near “‘100000′”
    LINE 1: SELECT $1 ::INT >= $2 ‘100000’::INT
    ^
    QUERY: SELECT $1 ::INT >= $2 ‘100000’::INT
    CONTEXT: SQL statement in PL/PgSQL function “googleencodesignedinteger” near line 14

    carlo

    4 Jul 15 at 12:10 pm

  8. Hey Elrobis,

    Was wondering if you’d be kind enough to shoot me an email. Hoping to pick your brain and so bootstrap some thinking I’m doing on some Google route data / PostGIS / searching goodness I’ve got in the works.

    Cheers. 🙂

    L

    Luke

    6 Oct 15 at 2:19 pm

Leave a Reply