PostGis Encode Polyline

Functions

-- Based on the word done by elrobis
-- http://cartometric.com/blog/2014/01/27/postgresql-postgis-implementation-of-google-encoded-polyline-algorithm

CREATE OR REPLACE FUNCTION pim.encode_polyline(
    in_wkt character varying
)
RETURNS character varying[] AS $$
DECLARE
  g GEOMETRY;
BEGIN
  g = ST_GeomFromText(in_wkt, 4326);
  RETURN ARRAY[pim.google_encode_polygon(g)];


END
$$ LANGUAGE plpgsql;

ALTER FUNCTION encode_polyline(character varying)
  OWNER TO postgres;


CREATE OR REPLACE FUNCTION pim.google_encode_line
(
  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 = COALESCE(ST_PointN(g, p), ST_GeometryN(g, p), g);
    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 = pim.google_encode_signed_integer(deltaX);
    enY = pim.google_encode_signed_integer(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 pim.google_encode_signed_integer(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;



CREATE OR REPLACE FUNCTION pim.google_encode_polygon
(
  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 = COALESCE(ST_NumGeometries(g1), 1);
 g = 1;

 IF (GeometryType(g1) = 'POLYGON' OR GeometryType(g1) = 'MULTIPOLYGON'  OR GeometryType(g1) = 'MULTILINESTRING') THEN
     FOR g IN 1..ng BY 1 LOOP

         g2 = COALESCE(ST_GeometryN(g1, g), g1);

         IF (GeometryType(g1) = 'POLYGON' OR GeometryType(g1) = 'MULTIPOLYGON') THEN
            -- Get ExteriorRing now; if there are any holes, get them later in the loop..
            r1 = ST_ExteriorRing(g2);
            gEncoded = gEncoded || pim.google_encode_line(r1);
        ELSE
            gEncoded = gEncoded || pim.google_encode_line(g2);
        END IF;
     END LOOP;
  ELSE
    gEncoded = gEncoded || pim.google_encode_line(g1);

  END IF;

 RETURN gEncoded;
End
$$ LANGUAGE plpgsql;

END
$dontcollide$;

Tests

import unittest

from django.db import connection

class TestCityStats(unittest.TestCase):

    def _run_query(self, query):
        cursor = connection.cursor()
        cursor.execute(query)
        row = cursor.fetchone()
        return row[0]

    def test_point(self):
        result = self._run_query("select pim.encode_polyline(st_astext('POINT(0 1)'::geometry))");
        self.assertEqual(result, ['_ibE?'])

    def test_point_one(self):
        result = self._run_query("select pim.encode_polyline(st_astext('MULTIPOINT((0 0),(1 2))'::geometry))");
        self.assertEqual(result, ['??_seK_ibE'])

    def test_point_two(self):
        result = self._run_query("select pim.encode_polyline(st_astext('MULTIPOINT(0 0,1 2)'::geometry))");
        self.assertEqual(result, ['??_seK_ibE'])

    def test_linestring(self):
        result = self._run_query("select pim.encode_polyline(st_astext('LINESTRING(2 3,3 2,5 4)'::geometry))");
        self.assertEqual(result, ['_}hQ_seK~hbE_ibE_seK_seK'])

    def test_multilinestring(self):
        result = self._run_query("select pim.encode_polyline(st_astext('MULTILINESTRING((0 0,1 1,1 2),(2 3,3 2,5 4))'::geometry))");
        self.assertEqual(result, ['??_ibE_ibE_ibE?_}hQ_seK~hbE_ibE_seK_seK'])

    def test_polygon(self):
        result = self._run_query("select pim.encode_polyline(st_astext('POLYGON((0 0,4 0,4 4,0 4,0 0),(1 1, 2 1, 2 2, 1 2,1 1))'::geometry))");
        self.assertEqual(result, ['???_glW_glW??~flW~flW?'])

    def test_multipolygon_one(self):
        result = self._run_query("select pim.encode_polyline(st_astext(" \
            "'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))");
        self.assertEqual(result, ['???_glW_glW??~flW~flW?~hbE~hbE~hbE??~hbE_ibE??_ibE'])

    def test_multipolygon_two(self):
        result = self._run_query("select pim.encode_polyline(st_astext('MULTIPOLYGON(((-1 -1,-1 -2,-2 -2,-2 -1,-1 -1)))'::geometry))");
        self.assertEqual(result, ['~hbE~hbE~hbE??~hbE_ibE??_ibE'])