From some millions of points I am trying to select those points that are within a polygon and save them into a new table. Both datasets with spatial index and both on projected coordinates.
The query I am trying with is:
CREATE TABLE result AS SELECT t1.point_id, t1.geom FROM points t1 INNER JOIN boundaries t2 ON ST_WITHIN(t1.geom, t2.geom) WHERE t2.country LIKE "ENGLAND"
The query has been running for hours.
On the same machine I have done exactly the same thing using ArcMap in 25 minutes. I really don't understand this. I was under the impression that doing geoprocessing tasks directly on a RDBS (PostgreSQL or SQL Server) was the fastest option (provided the query is the right one).
Is there anything wrong with the query? Is there any way to optimise it? Am I missing something here?
PS: I tried with ST_Intersect instead of ST_WITHIN with the same result.
As suggested in 'Comments', I have tried with FROM being the smaller table (table with 3 polygons).
CREATE TABLE result AS SELECT t1.point_id, t1.geom FROM boundaries t2 INNER JOIN points t1 ON ST_WITHIN(t1.geom, t2.geom) WHERE t2.country LIKE "ENGLAND"
CREATE TABLE result AS SELECT t1.point_id, t1.geom FROM boundaries t2 INNER JOIN points t1 ON ST_Contains(t2.geom, t1.geom) WHERE t2.country LIKE "ENGLAND"
Unfortunately, these two options seem to be pretty slow too. 2 hours running and still processing. I must admit that the polygon geometry is pretty complex (England boundary). But still, as indicated in the descriptions, other GIS tools take 25 minutes in doing the same job.
Explain analysis: (from last query written above, the one with "ST_Contains")
"Nested Loop (cost=111.22..26768.81 rows=3007 width=44)" " -> Seq Scan on gb_bdry t2 (cost=0.00..17.88 rows=3 width=32)" " Filter: ((country)::text ~~ 'ENGLAND'::text)" " -> Bitmap Heap Scan on final_def_sgm_centroids t1 (cost=111.22..8906.96 rows=1002 width=44)" " Recheck Cond: (t2.geom ~ geom)" " Filter: _st_contains(t2.geom, geom)" " -> Bitmap Index Scan on final_def_sgm_centroids_geom_idx (cost=0.00..110.96 rows=3007 width=0)" " Index Cond: (t2.geom ~ geom)"
Different approach following suggestion from losbaltica:
First we limit the number of points with &&:
CREATE TABLE temp1 AS SELECT t1.point_id, t1.geom FROM points t1 WHERE t1.geom && (SELECT geom FROM boundaries WHERE country = 'ENGLAND')
Then we do the actual intersection:
CREATE TABLE temp2 AS SELECT t1.point_id, t1.geom FROM points t1, (SELECT geom FROM boundaries WHERE country = 'ENGLAND') AS t2 WHERE ST_Intersects(t1.geom, t2.geom)