2

I have a table of polygons with a column 'class'. In that column are 25 different categories; however, I am only concerned with three polygon classes: marsh, water, and grass. I would like to extract only those polygons that touch each other among those three classes. I'm looking for continuous areas from the water, marsh, and grass areas only. I have an illustration of a test data set (water = blue, marsh = green, grass = yellow).

The following sql code is a step forward, but it doesn't work for three classes only two. It also doesn't select some of the elements:

DROP TABLE IF EXISTS adj;
    CREATE TABLE adj AS
    WITH m1 AS (
    SELECT a.fid, a.class, (ST_Dump(ST_Union(a.geom, b.geom))).geom
    FROM test as a,
         test as b
    WHERE a.fid < b.fid
    AND ST_Touches(a.geom, b.geom)
    AND a.class = 'MARSH'
    AND b.class = 'WATER'
    GROUP BY a.fid, a.class, a.geom, b.geom
        )
    SELECT * FROM m1;

Image of test data

MJM
  • 975
  • 4
  • 17
  • 4
    Is the condition actually ST_Touches - as in those geometries may not share their interiors - or is at least touches what you are looking for, i.e. ST_Intersects? If the latter, you want to use ST_ClusterDBSCAN with eps := 0. – geozelot Dec 15 '23 at 17:19
  • 2
    https://gis.stackexchange.com/questions/301598/only-union-dissolve-intersecting-or-adjacent-features-to-speed-up-query – BERA Dec 15 '23 at 18:14
  • To make your query work, it would have to be recursive (if A touches B touches C, current code would create AB and also BC), and you would have to include all classes AND a.class in ( 'MARSH', 'water','grass') and b.class in ( 'MARSH', 'water','grass'). But really st_clusterDBScan is the way to go – JGH Dec 15 '23 at 19:00

1 Answers1

2

It sounds like you want to find intersecting clusters of polygons from the specified classes which have cluster count > 1.

Here's a solution using the new ST_ClusterIntersectingWin window function (the data query just provides some sample data):

WITH data(fid, class, geom) AS (
  SELECT ROW_NUMBER() OVER (), 
      CASE x WHEN 5 THEN 1 ELSE x END AS class,
      ST_Buffer(ST_Point(x, 1.5 * y), 0.6, 2) AS geom
    FROM        generate_series(1, 10) AS s1(y)
    CROSS JOIN  generate_series(1, 10) AS s2(x)
),
datacls AS (
  SELECT * FROM data WHERE class IN (1, 2, 3)
),
clust AS (
  SELECT ST_ClusterIntersectingWin(geom) OVER () AS clustid, fid, class, geom FROM datacls
),
clustercnt AS (
  SELECT clustid, COUNT(*) AS cnt FROM clust GROUP BY clustid
)
SELECT fid, class, geom FROM clust c
  JOIN clustercnt cs ON c.clustid = cs.clustid
  WHERE cnt > 1;

Note: one possible issue with this solution is that clusters will include polygons which touch at a point, rather than along an edge. This is a limitation of ST_ClusterIntersectingWin, which might be enhanced in a future version.

dr_jts
  • 5,038
  • 18
  • 15
  • I replaced ST_ClusterIntersectingWin(geom) with ST_ClusterDBSCAN(geom, 0,3). I'm having problems with upgrading my postgis from 3.2.0 to 3.4.0. – MJM Dec 16 '23 at 00:25
  • ST_ClusterDBSCAN is handy since it lets you specify the minimum cluster size, which avoids the clustercnt logic. – dr_jts Dec 17 '23 at 01:27