4

This question is related to Seeking PostGIS equivalent to GeoPandas union overlay?

I have two tables in PostGIS, each with two features(rows) and a rectangular geometry (polygon).

enter image description here

I would like to obtain the Symmetrical Difference from those two layers so that the result would be:

enter image description here

the resulting table has the attributes from the first table (df1), or the second (df2)

enter image description here

However when I try to use ST_SymDifference in PostGIS I get a different result:

CREATE TABLE test.symdifference_v04 AS
SELECT 
    test.df1.df1,
    test.df2.df2,
    ST_SymDifference(test.df1.geom,test.df2.geom)
FROM test.df1, test.df2

enter image description here

It seems that the ST_SymDifference is applied on a polygon to polygon basis i.e. df1 row 1 vs df2 row 1
df1 row 1 vs df2 row 2
df1 row 2 vs df2 row 1
df1 row 2 vs df2 row 2

with Multipolygons as a result.

enter image description here enter image description here enter image description here enter image description here

How would I obtain the first mentioned result with the 6 polygons (fancy colors) instead of the 4 multiploygons in PostGIS?

PolyGeo
  • 65,136
  • 29
  • 109
  • 338
RutgerH
  • 3,285
  • 3
  • 27
  • 55

2 Answers2

3

This query should help you solve the problem:

variant 1 -

1) Create wanted parts:

CREATE TABLE table_cut as SELECT DISTINCT a.id, ST_Difference(a.geom, b.geom) AS geom FROM table1 AS a JOIN table2 AS b ON ST_Intersects(a.geom, b.geom) GROUP BY a.id, a.geom, b.geom;

CREATE TABLE table_cut2 as SELECT DISTINCT a.id, ST_Difference(a.geom, b.geom) AS geom FROM table2 AS a JOIN table1 AS b ON ST_Intersects(a.geom, b.geom) GROUP BY a.id, a.geom, b.geom;

CREATE TABLE table_cut1 as SELECT DISTINCT a.id, ST_Intersection(a.geom, b.geom) AS geom FROM table1 AS a JOIN table2 AS b ON ST_Intersects(a.geom, b.geom) GROUP BY a.id, a.geom, b.geom;

2) Create a results tables

CREATE TABLE table3 AS SELECT DISTINCT a.id, (ST_Dump(ST_Difference(a.geom, b.geom))).geom FROM table_cut2 a JOIN table_cut1 b ON ST_Intersects(a.geom, b.geom) GROUP BY a.id, a.geom, b.geom;

CREATE TABLE table4 AS SELECT DISTINCT a.id, (ST_Dump(ST_Difference(a.geom, b.geom))).geom FROM table_cut a JOIN table_cut1 b ON ST_Intersects(a.geom, b.geom) GROUP BY a.id, a.geom, b.geom

3) Remove the two internal parts and symbolize the cells ...

enter image description here

Interesting geometric puzzle: -) ...

option 2 -

WITH
tbla AS (SELECT * FROM table1 UNION SELECT * FROM table2),
tblb AS (SELECT a.id, ST_Intersection(a.geom, b.geom) AS geom FROM table1 AS a JOIN table2 AS b ON ST_Intersects(a.geom, b.geom) GROUP BY a.id, a.geom, b.geom)
(SELECT a.id, (ST_Dump(ST_Difference(a.geom, b.geom))).geom AS geom FROM tbla AS a CROSS JOIN LATERAL (SELECT ST_Collect(geom) AS geom FROM tblb WHERE ST_Intersects(a.geom, geom)) AS b);

...

Cyril Mikhalchenko
  • 4,397
  • 7
  • 14
  • 45
2

You can use st_difference to get your result:

First we create our two dummy tables:

-- FIRST SET OF POLYGON
CREATE TABLE rect1 (name varchar, geom geometry);

INSERT INTO rect1 VALUES
  (1, 'POLYGON((0 0, 1 0, 1 1, 0 1, 0 0))'),
  (2, 'POLYGON((1 1, 2 1, 2 2, 1 2, 1 1))');

SELECT name, ST_AsText(geom) FROM rect1;

-- SECOND SET OF POLYGON
CREATE TABLE rect2 (name varchar, geom geometry);

INSERT INTO rect2 VALUES
  (1, 'POLYGON((0.5 0.5, 1.5 0.5, 1.5 1.5, 0.5 1.5, 0.5 0.5))'),
  (2, 'POLYGON((1.5 1.5, 2.5 1.5, 2.5 2.5, 1.5 2.5, 1.5 1.5))');

SELECT name, ST_AsText(geom) FROM rect2;

Now we can produce the result by using st_difference two time:

  1. The difference between rect1 and rect2
  2. The difference between rect2 and rect1

Noticed that we have to group all the "cutting" geometry with st_union.

WITH T1 AS
(
    SELECT st_difference(r1.geom,st_union(r2.geom)) as geom
    FROM rect1 r1, rect2 r2
    GROUP BY r1.geom

    UNION

    SELECT st_difference(r2.geom,st_union(r1.geom)) as geom
    FROM rect1 r1, rect2 r2
    GROUP BY r2.geom
)
SELECT (st_dump(geom)).geom from T1

We use ST_Dump to split our multipart polygon to singlepart polygon.

As expected we obtain:

enter image description here

obchardon
  • 1,724
  • 2
  • 13
  • 21