ST_difference resulting in multiple overlapping polygons and no difference

Question

I have a dataset with two layers of polygons:

Layer A: Areas of interest from the corine land cover inventory dataset

Layer B: Road network resulting from buffered linestrings that were aggregated by st_union (in combination with st_dump to expand the multipolygon)

Below are layer A and B (excerpt from dataset):

I want to subtract layer B from layer A by using st_difference. My SQL query looks like this:

CREATE TABLE diff AS
SELECT st_difference(clc.geom, roads.geom) AS geom
FROM clc, roads
WHERE st_intersects(clc.geom, roads.geom);

When looking at the result I was very confused at first:

There was no difference in the upper part! While troubleshooting I changed the layer opacity and I noticed, that the difference was below a polygon with the same shape but without the difference:

Somehow the output of st_difference results in multiple overlapping polygons and some of them are the result of the difference and some are not.

Layer A and B do not have multiple overlapping polygons themselves.

Answer 1

You have done a full cross table join, so you have intersected every road against every land cover polygon in turn, so you have (N-1) polygons that don't intersect for ever polygon that does (where N is the number of roads).

To get the answer you expect you need to union the land cover polygons and then subtract the roads from that polygon (or union the roads and subtract that from each polygon in turn). So something like:

SELECT st_difference(clc.geom, roads.geom) AS geom
FROM clc, (SELECT st_union(geom) from roads) as roads;