1

I have three spatial datasets - mitch_ibra has about 2700 polygons, hev has ~1,200,000 polygons and prop_p1 with about 180,000. For each polygon in mitch_ibra I want to know how many polygons of prop_p1 are inside, with a minimum size of 5ha or overlaps the mitch_ibra polygon by 5ha. Finally I want the number of these polygons intersect polygons from the hev table.

WITH region AS
  (SELECT geom
   FROM mitch_ibra
   WHERE sub_name_7 = 'Tingha Plateau'
     AND lscape_nam = 'Dumaresq Channels')

SELECT COUNT(props.id)
FROM
  (SELECT DISTINCT initprop.id
   FROM
     (SELECT id,geom
      FROM prop_p1
      WHERE st_area(ST_transform(geom,3308))>50000) AS initprop
   JOIN region ON st_intersects(initprop.geom,region.geom),
     (SELECT hev.geom FROM hev JOIN region ON ST_Intersects(hev.geom,region.geom)) AS hevsel
   WHERE CASE
             WHEN ST_Within(initprop.geom,region.geom) THEN ST_Area(ST_Transform(initprop.geom,3308))
             ELSE ST_Area(ST_Transform(ST_Intersection(initprop.geom,region.geom),3308))
         END > 50000
         AND ST_Intersects(initprop.geom,hevsel.geom)) AS props

Above is my query so far, which gives me a processing time of about 658s for a region with 1300 polygons selected in the prop_p1 subquery. All have tables have GIST indexes. How can I optimize my indexes/query to reduce processing time. Is this as good as I'm going to get? (I have to run this over 300 different selections in the WITH statement)

  • Indexes cannot be applied to a CTE. You can use a materialized view instead, or a temp table, or replace references to region by the complete query (i.e. remove the CTE) – JGH Sep 21 '17 at 11:58
  • Thanks, that actually got me ~ 3.5% decrease in processing time – Liam G Sep 21 '17 at 22:52
1

Indexes cannot be applied to a CTE. You can use a materialized view instead, or a temp table (both with an index), or replace references to region by the complete query (i.e. remove the CTE)

3

You could check if the polygons are at least close to eachother, before checking if they intersect:

...JOIN region ON initprop.geom&&region.geom and st_intersects(initprop.geom,region.geom)...

http://postgis.net/docs/manual-2.3/geometry_overlaps.html

And can you tune your PostGres configuration? Adding work_mem (at least for this session) can do miracles.

  • +1 for Postgres config. Note that ST_Intersects includes a bounding box comparison that makes use of the index. – JGH Sep 21 '17 at 11:50
  • Adding && bounding box intersection actually decreased performance unfortunately, but i'll have a play around with the tuning - thanks. – Liam G Sep 21 '17 at 22:58
  • Actually by replacing the st_intersect in the hev, region join with the && bounding box join, I saved 28% processing time – Liam G Sep 21 '17 at 23:15
1

I'm not sure I got it right, but I'd try this query:

SELECT
  props.mi_id,
  COUNT(1)
FROM (
  SELECT
    a.id,
    a.geom
  FROM (
    SELECT
      mi.id,
      mi.geom,
      ST_Area(ST_Transform(ST_CollectionExtract(ST_Intersection(mi.geom, pp1.geom)), 3308)) area
    FROM mitch_ibra mi
    JOIN prop_p1 pp1 ON (ST_Intersects(mi.geom, pp1.geom))
  ) a
  GROUP BY
    a.id,
    a.geom
  HAVING SUM(area) > 50000
) props
JOIN hev ON (ST_Intersects(props.geom, hev.geom));

This gives you count of hev ids that intersect with mitch_ibra ids that intersect with prop_p1 and that intersection is > 50000.

If that's not what you're after, be a little more specific please. Sorry for potential typos, can't test the query without the tables.

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

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