Wednesday 11 September 2019

openstreetmap - Using PostGis BBox and ST_intersects Queries?



CartoDb Recommend here that the following query is faster than an St_intersects Query


SELECT * FROM {table_name} WHERE the_geom && ST_SetSRID(ST_MakeBox2D(ST_Point(-73.9980, 40.726), ST_Point(-73.995, 40.723)), 4326)

I have the following Query that works just fine and returns GeoJson as needed,


WITH c AS (SELECT -73.976440429688, 40.743095232182, -73.970947265625, 40.738933241136, 4326 As srid)  ,
bbox AS (SELECT ST_Transform ( ST_MakeEnvelope( -73.976440429688, 40.743095232182 ,-73.970947265625 , 40.738933241136 , 4326) , 900913 ) As geom FROM c )
SELECT row_to_json(fc)
FROM ( SELECT 'FeatureCollection' As type, array_to_json(array_agg(f)) As features
FROM (SELECT 'Feature' As type, osm_id As id,

( SELECT row_to_json(t)
FROM (SELECT osm_id, "name" ) t )As properties,
ST_AsGeoJSON(ST_Transform (lg.way ,4326 ) )::json As geometry
FROM planet_osm_polygon As lg, bbox
WHERE ( building is NOT NULL) AND ST_Intersects(lg.way, bbox.geom)
) As f ) fc;

I have been trying to create the same JSon result using ST_MakeBox2D in a similar fasion to the CartoDb example and failing miserably.


My objective is to have the fastest possible query for any geometry that any part of itself WITHIN the bounding box.


My Questions are:




  1. Which method will return the fastest possible result, assuming that the Columns in the WHERE clause are indexed properly?

  2. How do I change my query to use St_makeBox2D with the geom && syntax as per the CartoDb Example?



Answer



The use of && will return quicker than ST_Intersects, as it is only checking bounding boxes, which is the first stage of any spatial intersection query, via the bounding boxes stored as part of R-Tree indexing (technically a GIST index in Postgis, but the bounding box is still used in index construction). As, the second stage, the actual intersection is being skipped, you will save time, though this will depend on the complexity of the actual polygons.


You can simply change this line,


WHERE building is NOT NULL AND ST_Intersects(lg.way, bbox.geom)

to



WHERE building is NOT NULL AND lg.way && bbox.geom;

in order to just do a bounding box check instead of a full intersects. You can check the relative speed by using the \timing switch in psql. Note, however, that because an ST_Intersects query checks the bounding box as part of the calculation, this part of the query will get cached, so you need to clear the query cache to get reliable testing results of ST_Intersects vs &&. This post has information on clearing the cache.


You can simplify the initial creation of you bbox to one line, as


WITH bbox AS (SELECT ST_Transform(ST_MakeEnvelope(-73.976440429688, 40.743095232182,
-73.970947265625, 40.738933241136, 4326), 900913) As geom)

which will be exactly the same as using ST_MakeBox2D(ST_MakePoint(x1, y1), ST_MakePoint(x2, y2).


Putting this all together yields,


WITH bbox AS 

(SELECT ST_Transform(ST_MakeEnvelope(-73.976440429688, 40.743095232182,
-73.970947265625, 40.738933241136, 4326), 900913) As geom)
SELECT row_to_json(fc)
FROM ( SELECT 'FeatureCollection' As type, array_to_json(array_agg(f)) As features
FROM (SELECT 'Feature' As type, osm_id As id,
(SELECT row_to_json(t)
FROM (SELECT osm_id, "name" ) t )As properties,
ST_AsGeoJSON(ST_Transform (lg.way ,4326 ) )::json As geometry
FROM planet_osm_polygon As lg, bbox
WHERE ( building is NOT NULL) AND lg.way && bbox.geom

) As f
) fc;

No comments:

Post a Comment

arcpy - Changing output name when exporting data driven pages to JPG?

Is there a way to save the output JPG, changing the output file name to the page name, instead of page number? I mean changing the script fo...