Monday, 25 February 2019

Better/Faster ST_PointOnSurface for PostGIS

And now for the final chapter in the saga of improving InteriorPoint / PointOnSurface.  For those who missed the first two episodes, the series began with a new approach for the venerable JTS Geometry.interiorPoint() for polygons algorithm.  Episode 2 travelled deep into the wilds of C++ with a port to GEOS.  The series finale shows how this results in greatly improved performance of PostGIS ST_PointOnSurface.

The BC Voting Area dataset is a convenient test case, since it has lots of large polygons (shown here with interior points computed).
The query is about as simple as it gets:

   select ST_PointOnSurface(geom) from ebc.voting_area;

Here's the query timings comparison, using the improved GEOS code and the previous implementation:

Data size Time Time
OLD		 Improvement Time
ST_Centroid
5,658 polygons
(2,171,676 vertices)		 341 ms 4,613 ms x 13 369 ms

As expected, there is a dramatic improvement in performance.  The improved ST_PointOnSurface runs 13 times faster than the old code.  And it's now as fast as ST_Centroid.  It's also more robust and tolerant of invalid input (although this test doesn't show it).

This should show up in PostGIS in the fall release (PostGIS 3 / GEOS 3.8).

On to the next improvement... (and also gotta update the docs and the tutorial!)

at

2 comments:

  1. Bruce Rindahl26 February 2019 at 14:04

    Just curious, does ST_PointOnSurface(geom) return a point on the interior of the surface? In your documentation link, example 1 returns the point (obviously), and example 3 returns the center of a square. However, examples 2 and 4 returns an edge vertex of the surface.
    The tutorial link specifies "ST_PointOnSurface(geometry) returns a point that is guaranteed to be inside the input argument."

    Thanks!

    ReplyDelete
  2. Dr JTS26 February 2019 at 14:26

    Yes, ST_PointOnSurface returns a point in the interior of the surface - when it is given a polygonal input (e.g. examples 1 and 3 in the doc).

    JTS/GEOS and hence PostGIS extend this to handle other input geometry types too, such as lines and points. In those cases, an interior vertex/point is returned. This may or may not be useful, but it's nice to have spatial functions handle all geometry types where possible.

    The docs could be clearer on this aspect - hope to improve them sometime soon.

    ReplyDelete