Author: Regina Obe 2014/03/14 in tips ( newbie, geometry, ST_Intersects, ST_Intersection, ST_CoveredBy )

Doing an `ST_Intersection`

is much slower than relation checks such as
`ST_Intersects`

, `ST_CoveredBy`

, and , `ST_Within`

.
In many situations you know the intersection of 2 geometries
without actually computing an intersection. In these cases, you can skip the costly `ST_Intersection`

call. Cases like this:

- geometry a is covered by geometry b -> intersection is geometry a
- geometry b is covered by geometry a -> intersection is geometry b
- geometry a does not intersect geometry b -> intersection is empty geometry

This kind of question comes up a lot: As discussed in stackexchange Acquiring ArcGIS speed in PostGIS

**Examples**

For this exercise, we’ll grab the portion of each parcel that falls in a neighborhood.
Unlike stack exchange one, we’ll use `ST_CoveredBy`

instead of `ST_Within`

check. Both constructs are similar except that `ST_CoveredBy`

allows
for geometries to be wholly within the boundary of another and tends to be a bit faster to compute than `ST_Within`

.
This subtlety becomes more important if you are comparing linestrings such as what portion of a road falls in a county
as detailed in Subtleties OGC Covers Spatial.

SELECT p.parcel_id, n.nei_name , CASE WHEN ST_CoveredBy(p.geom, n.geom) THEN p.geom ELSE ST_Multi( ST_Intersection(p.geom,n.geom) ) END AS geom FROM parcels AS p INNER JOIN neighborhoods AS n ON ST_Intersects(p.geom, n.geom);

In this particular case, we’d probably want to exclude the case where a parcel just borders a neighborhood (in spatial speak `touches`

).
We wouldn’t really consider a bordering parcel as having any part in the neighborhood, though it intersects the neighborhood.
So a slightly more complicated but more accurate statement would be to exclude from consideration the case where a parcel borders a neighborhood using
`ST_Touches`

.

SELECT p.parcel_id, n.nei_name , CASE WHEN ST_CoveredBy(p.geom, n.geom) THEN p.geom ELSE ST_Multi( ST_Intersection(p.geom,n.geom) ) END AS geom FROM parcels AS p INNER JOIN neighborhoods AS n ON (ST_Intersects(p.geom, n.geom) AND NOT ST_Touches(p.geom, n.geom) );

This raster question comes up quite a bit on PostGIS mailing lists and stack overflow and the best answer often involves
the often forgotten `ST_Reclass`

function that has existed since PostGIS 2.0.

People often resort to the much slower though more flexible `ST_MapAlgebra`

or dumping out
their rasters as Pixel valued polygons they then filter
with WHERE val > 90,
where `ST_Reclass`

does the same thing but orders of magnitude faster.

PostGIS 2.2.0 came out this month, and the SFCGAL extension that offers advanced 3D and volumetric support, in addition to some extended 2D functions like `ST_ApproximateMedialAxis`

became a standard PostgreSQL extension
and seems to be a fairly popular extension.

I’ve seen several reports on GIS Stack Exchange of people trying to install PostGIS SFCGAL and getting things such as

ERROR: could not open extension control file /Applications/Postgres.app/Contents/Versions/9.5/share/postgresql/extension/postgis_sfcgal.control

as in this question: SFCGAL in PostGIS problem.

There are two main causes for this:

I’ve been asked this question in some shape or form at least 3 times, mostly from people puzzled why they get this error. The last iteration went something like this:

*I can’t use ST_AsPNG when doing something like*

SELECT ST_AsPNG(rast) FROM sometable;

Gives error: Warning: pg_query(): Query failed: ERROR: rt_raster_to_gdal: Could not load the output GDAL driver.