1

I want to select some points based on a intersection with a polygon coming from an other database as a WKT geometry.

For the moment, everything is working fine:

SELECT p.name, p.geom, p.class
FROM point_table p
WHERE ST_Intersects(
        ST_Transform(p.geom, 2056), 
        ST_Transform(ST_GeomFromText('POLYGON(((6.12 46.22,...,6.12 46.22)))', 4326), 2056)
      )
AND (...other conditions...);

But let's say the point have a 'class' attribute depicting the object they represent. When class='city' I not only want to select the point which are inside the polygon, but also those which are in a 1000m buffer and for class='Town' I'd like to apply a 500m buffer.

I can do this in two different ways; either buffering the polygon itself, or the points, where the class attribute of the points matches city or town using the appropriate distance for each.

As the polygon is quite large, I prefer to buffer the points themselves.

I was thinking of something like:

SELECT p.name, p.geom, p.class
FROM point_table p
WHERE 
  CASE
    WHEN p.class = 'City'
      THEN ST_Intersects(
             ST_Buffer(ST_Transform(p.geom, 2056), 1000),
             ST_Transform(ST_GeomFromText('POLYGON(((6.12 46.22,...,6.12 46.22)))', 4326), 2056)
           )
    WHEN p.class = 'Town'
      THEN ST_Intersects(
             ST_Buffer(ST_Transform(p.geom, 2056), 500),
             ST_Transform(ST_GeomFromText('POLYGON(((6.12 46.22,...,6.12 46.22)))', 4326), 2056) 
           )
    ELSE ST_Intersects(
           ST_Transform(p.geom, 2056),
           ST_Transform(ST_GeomFromText('POLYGON(((6.12 46.22,...,6.12 46.22)))', 4326), 2056)
         )
  END
AND (...other conditions...);

I think this actually works, but I'm wondering if there is a better approach?

For example, one which would not repeat the massive polygon 3x.

I always want to "pre-compute" the geometry in my SELECT statement, using an alias, and then using this, but it's not interpreted in the rest of the query as if it doesn't exist at that moment, but only at the very end of the query:

SELECT p.name, p.geom, p.class,
CASE 
  WHEN p.class = 'City'
    THEN ST_Buffer(ST_Transform(p.geom, 2056), 1000)
  WHEN p.class = 'Town'
    THEN ST_Buffer(ST_Transform(p.geom, 2056), 500)
  ELSE p.geom
END new_geom --<---- That would be a nice solution, but I cannot use 'new_geom' hereafter.
FROM point_table p
WHERE ST_Intersects(
        new_geom,
        ST_Transform(ST_GeomFromText('POLYGON(((6.12 46.22,...,6.12 46.22)))', 4326), 2056)
      )
AND (...other conditions...);
1

1 Answer 1

3

Some notes about your queries and ideas:

  • in 9 out of 10 cases if you could use a buffer, don't!
  • computational complexity increases with M * n computation steps (with M & n being the amount of vertices included per geometry); with M = 100000 (input Polygon), the difference of n = 1 (Point) and n = 32 (default buffer Polygon) per geometry comparison is massive
  • whatever expression you use, in any direct or indirect filter conditions (WHERE, JOIN, ORDER BY, ...), needs to be covered explicitly by an index to gain index driven performance; this includes the ST_Transform result
  • depending on volatility, and in non-fenced queries, PostgreSQL' s query planner is perfectly able to cache simple expressions and function results; this should be the case with your repeated polygon creation - in either way, this is the least of your concerns, performance-wise

Some considerations:

  • you absolutely want to use ST_DWithin; geometric proximity searches can usually truthy out early, while intersections and the likes always need a full scan over the Cartesian product of vertices
  • PostgreSQL cannot resolve dynamic index conditions; with a CASE statement as condition, PostgreSQL will fall back to sequentially scanning for matches

With all this in mind, you generally want to run sth. like this:

SELECT *
FROM   <points> AS pt
WHERE  ST_Transform(pt.geom, <SRID>) && ST_Expand(ST_Transform(<POLYGON>, <SRID>), 1000)
 AND   ST_DWithin(
         ST_Transform(pt.geom, <SRID>),
         ST_Transform(<POLYGON>, <SRID>),
         CASE
           WHEN p.class = 'City' THEN 1000
           WHEN p.class = 'Town' THEN 500
           ELSE 0
         END
       )
;

or, with a join on a CTE

WITH
  poly AS MATERIALIZED (
    SELECT ST_Transform(<PLYGON>, <SRID>) AS geom
  )
SELECT *
FROM   <points> AS pt
JOIN   poly
  ON   ST_Transform(pt.geom, <SRID>) && ST_Expand(poly.geom, 1000)
WHERE  ST_DWithin(
         ST_Transform(pt.geom, <SRID>),
         poly.geom,
         CASE
           WHEN p.class = 'City' THEN 1000
           WHEN p.class = 'Town' THEN 500
           ELSE 0
         END
       )
;

where the index gets utilized for all matches in the largest proximity via an explicit && bbox filter, so that calculations are only done on that results set.

Note that

  • a functional index on the transformed result is mandatory for index driven performance, i.e.

    CREATE INDEX ON <points> USING GIST ( ST_Transform(<points>.geom, <SRID>) );
    
  • the same could be achieved using the GEOGRAPHY type, for when a single projection is out-scoped

  • this would ultimately be improved by ST_SubDivide'ing your Polygon into a(n indexed) table


Depending on a multitude of factors (with some of the more obvious being table size, category size and approx. result set size) it may be more performant to run different queries, e.g. a UNION ALL of two separate ST_DWithin calls on category filters - but that is some later stage optimization and probably not necessary with the general optimization in the above query.

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.