1

I need to do a point (summit) in polygon (contour) check on a very complex (contour) polygon. I only need the highest point within the polygon (ORDER BY elevation DESC LIMIT 1) . ST_Within takes forever (2-3 minutes) with the SORT clause on the result, but is very quick if I forget the order and just apply a limit of 1.

I thought maybe I could sort the point input data to ST_Within (very fast to do) and then just use LIMIT 1 on the result of ST_Within, and get only the highest summit. But that does not reliably return the highest point, so I guess the WHERE ST_Within does not check the points in the order they are supplied But does that reliably return the highest peak?

How can I do this efficiently & reliably? Can I force the query to test the input data in the order supplied? Or is there another trick I'm missing?

Original query is as follows:

SELECT * FROM peaks p 
  WHERE 
    p.elevation>#{ele} 
  AND 
    ST_Within(p.wkb_geometry, 
     (SELECT wkb_geometry FROM contour WHERE fid=#{contour.fid})) 
  ORDER BY p.elevation DESC LIMIT 1

The following appears to reliably return the highest peak, and do so much faster. But it relies on ST_Within in the outer query checking the input points in the order that they are presented by the sub-query. Is this approach valid, will this consistently return the same result as the above slower query?

Get all peaks within bounding box of contour and sort descending (~1 sec), then test that list for ST_Within the contour, returning only the first found (typically 3-10 secs)

SELECT * FROM (
  SELECT * FROM peaks p 
  WHERE 
    p.elevation>#{ele} 
  AND 
    ST_Within(
      wkb_geometry, 
      (
        SELECT ST_Envelope(wkb_geometry) 
        FROM contour 
        WHERE fid=#{contour.fid}
      )
    )     
  ORDER BY elevation DESC
) AS p2 
WHERE 
  ST_Within(
    p2.wkb_geometry, 
    (
       SELECT wkb_geometry 
       FROM contour 
       WHERE fid=#{contour.fid}
    )
  )  
LIMIT 1;

===

Further information as requested:

Table sizes:

 Schema |        Name         |   Type   |  Owner   |    Size    | Description 
--------+---------------------+----------+----------+------------+-------------
 public | contour             | table    | mbriggs  | 68 GB      | 
 public | contour_fid_seq     | sequence | mbriggs  | 8192 bytes | 
 public | peaks               | table    | mbriggs  | 7720 kB    | 
 public | peaks_ogc_fid_seq   | sequence | mbriggs  | 8192 bytes | 

Individual contour size (vertices) - max for example elevations:

# select ele, max(ST_NPoints(wkb_geometry)) as vert from contour  where ele in (1,100,1000,2000) group by ele;
 ele  |  vert   
------+---------
    1 | 1803355
  100 | 2097301
 1000 | 1019721
 2000 |   16959

Original query:

                                       QUERY PLAN                                       
----------------------------------------------------------------------------------------
 Limit  (cost=15637.61..15637.62 rows=1 width=204)
   InitPlan 1 (returns $0)
     ->  Index Scan using contour_pkey on contour  (cost=0.43..8.45 rows=1 width=44763)
           Index Cond: (fid = 13873)
   ->  Sort  (cost=15629.17..15641.21 rows=4816 width=204)
         Sort Key: p.elevation DESC
         ->  Seq Scan on peaks p  (cost=0.00..15605.09 rows=4816 width=204)
               Filter: ((elevation > '1000'::numeric) AND st_within(wkb_geometry, $0))

Faster query (but will this reliably return max altitude peak?):

                                                        QUERY PLAN                                                        
--------------------------------------------------------------------------------------------------------------------------
 Limit  (cost=230.97..230.98 rows=1 width=204)
   InitPlan 1 (returns $0)
     ->  Index Scan using contour_pkey on contour  (cost=0.43..8.45 rows=1 width=44763)
           Index Cond: (fid = 13873)
   ->  Sort  (cost=222.52..222.53 rows=2 width=204)
         Sort Key: p.elevation DESC
         InitPlan 2 (returns $1)
           ->  Index Scan using contour_pkey on contour contour_1  (cost=0.43..8.45 rows=1 width=32)
                 Index Cond: (fid = 13873)
         ->  Bitmap Heap Scan on peaks p  (cost=4.70..214.06 rows=2 width=204)
               Recheck Cond: ($0 ~ wkb_geometry)
               Filter: ((elevation > '1000'::numeric) AND st_within(wkb_geometry, $1) AND _st_contains($0, wkb_geometry))
               ->  Bitmap Index Scan on peaks_wkb_geometry_geom_idx  (cost=0.00..4.70 rows=56 width=0)
                     Index Cond: ($0 ~ wkb_geometry)

Indexing & fields used above:

dem1=# \d contour
                                       Table "public.contour"
    Column    |          Type          | Collation | Nullable |               Default                
--------------+------------------------+-----------+----------+--------------------------------------
 fid          | integer                |           | not null | nextval('contour_fid_seq'::regclass)
 ele          | double precision       |           |          | 
 wkb_geometry | geometry(Polygon,4326) |           |          | 
Indexes:
    "contour_pkey" PRIMARY KEY, btree (fid)
    "contour_wkb_geometry_geom_idx" gist (wkb_geometry)

Peaks:

dem1=# \d peaks
                                         Table "public.peaks"
    Column     |         Type         | Collation | Nullable |                Default                 
---------------+----------------------+-----------+----------+----------------------------------------
 ogc_fid       | integer              |           | not null | nextval('peaks_ogc_fid_seq'::regclass)
 elevation     | numeric(24,15)       |           |          | 
 wkb_geometry  | geometry(Point,4326) |           |          | 
Indexes:
    "peaks_pkey" PRIMARY KEY, btree (ogc_fid)
    "peaks_wkb_geometry_geom_idx" gist (wkb_geometry)
8
  • What about index on peaks.elevation like CREATE INDEX ... ON peaks USING btree (elevation)? Feb 29 at 20:20
  • I can try that. But won't that just make the final sort more efficient (which is already sub-second so little to gain). The issue as I see it is that instead of doing the sort first, then checking (starting at the highest) only until it finds the first peak within the contour, it is checking all 3000+ peaks AND THEN doing the sort.
    – madpom
    Feb 29 at 21:33
  • 1
    Please post the output of EXPLAIN ANALYZE <your_query>; and add info about table sizes, contour complexity (vertex counts) and current indexation.
    – geozelot
    Mar 1 at 10:21
  • 1
    Table sizes and EXPLAIN added. EXPLAIN ANALYZE would be pretty meaningless at present due to other activity on the database. Will add when current data processing ends.
    – madpom
    Mar 2 at 19:58
  • The primary question remains however: will my 'faster query' above always return the highest peak? i.e.: Can I rely on the outer query to test the points returned by the inner query in the order they are presented by the inner query?
    – madpom
    Mar 2 at 20:04

1 Answer 1

0

You should try using ST_DWithin, possibly as an AND filter on existing. See What is difference between ST_DWithin and ST_Distance for proximity search in PostGIS?

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.