[postgis-users] Spatial join optimisation

[Bill Binko]

On Fri, 13 May 2005, Rhys Ickeringill wrote:

> Hi List,
> 
> I have two large tables I wish to spatially join - one which contains
> polygons, the other points of latitude and longitude. The join criteria is
> on polygons containing points.
> 
> At present I have the query:
> 
> SELECT pntid, polyid
> INTO PolygonsContainingPoints
> FROM Polygons INNER JOIN Points
>     ON contains( Polygons.wkb_geometry, SetSRID( MakePoint(
> Points.longitude, Points.latitude ), 32767 ) )
> 

As a non-PostGIS expert, but one who's had some success in this recently, 
here's what I'd do:

A Quick and dirty way to make it use that index is this:

SELECT pntid, polyid INTO PolygonsContainingPoints FROM
  Polygons JOIN Points ON contains(wkb_geometry, 
  SetSRID(MakePoint(longitude,latitude), 32767)) AND wkb_geometry ~ 
  SetSRID(MakePoint(longitude, latitude), 32767));

The extra (superfluous) criteria should cause it to use the index.


However, I would really suggest something a bit more useful.  Create a 
column for the Points...

SELECT AddGeometryColumn('','Points','the_point_geom','32767','POINT',2);

update Points set the_point_geom =  SetSRID(MakePoint(longitude,latitude),32767);

create index point_index on Points USING GIST (the_point_geom);

vacuum analyze Points;

select update_geometry_stats(); -- depending on your Postgresql version


Then you could write this:

select pntid, polyid into PolygonsContainingPoints FROM 
  Polygons JOIN Points on contains(wkb_geometry, the_point_geom) AND
  wkb_geometry ~ the_point_geom;

This would use both index and result in a query measured in seconds, not 
hours.

Bill