Fusionar todos los polígonos adyacentes

Question

Me gustaría hacer pruebas de adyacencia en una capa de parcelas (polígonos) y fusionarlas si se ajustan a ciertos criterios (podría ser el tamaño). Según la imagen siguiente, me gustaría fusionar los polígonos 1,2,3 y 4, pero no 5.

Tengo dos problemas:

1. ST_TOUCHES devuelve TRUE si sólo se tocan las esquinas y no un segmento de línea. Creo que necesito ST_RELATE para comprobar si hay segmentos de línea compartidos.
2. Idealmente, me gustaría fusionar TODOS los polígonos adyacentes en uno, pero no estoy seguro de cómo escalar más allá de dos - es decir, fusionar 1,2,3 y 4 (y posiblemente más en los datos reales) en una ronda.

La estructura que tengo ahora se basa en una autounión en ST_TOUCHES .

Datos del juguete

CREATE TABLE testpoly AS 
SELECT 
1 AS id, ST_PolyFromText('POLYGON ((0 0, 10 0, 10 20, 00 20, 0 0 ))') AS geom UNION SELECT
2 AS id, ST_PolyFromText('POLYGON ((10 0, 20 0, 20 20, 10 20, 10 0 ))') AS geom UNION SELECT
3 AS id, ST_PolyFromText('POLYGON ((10 -20, 20 -20, 20 0, 10 0, 10 -20 ))') AS geom UNION SELECT
4 AS id, ST_PolyFromText('POLYGON ((20 -20, 30 -20, 30 0, 20 0, 20 -20 ))') AS geom  UNION SELECT 
5 AS id, ST_PolyFromText('POLYGON ((30 0, 40 0, 40 20, 30 20, 30 0 ))') AS geom ;

Selección

SELECT 
    gid, adj_gid,
    st_AStext(st_union(l2.g1,l2.g2)) AS geo_combo
from (
    --level 2
    SELECT
      t1.id AS gid,
      t1.geom AS g1,
      t2.id AS adj_gid,
      t2.geom AS g2
     from
      testpoly  t1,
      testpoly  t2
     where
      ST_Touches( t1.geom, t2.geom ) 
      AND t1.geom && t2.geom 
) 
l2

Este es el resultado:

+-----+---------+-------------------------------------------------------------------------------+
| gid | adj_gid | geo_combo                                                                     |
+-----+---------+-------------------------------------------------------------------------------+
| 1   | 2       | POLYGON((10 0,0 0,0 20,10 20,20 20,20 0,10 0))                                |
+-----+---------+-------------------------------------------------------------------------------+
| 1   | 3       | MULTIPOLYGON(((10 0,0 0,0 20,10 20,10 0)),((10 0,20 0,20 -20,10 -20,10 0)))   |
+-----+---------+-------------------------------------------------------------------------------+
| 2   | 1       | POLYGON((10 20,20 20,20 0,10 0,0 0,0 20,10 20))                               |
+-----+---------+-------------------------------------------------------------------------------+
| 2   | 3       | POLYGON((10 0,10 20,20 20,20 0,20 -20,10 -20,10 0))                           |
+-----+---------+-------------------------------------------------------------------------------+
| 2   | 4       | MULTIPOLYGON(((20 0,10 0,10 20,20 20,20 0)),((20 0,30 0,30 -20,20 -20,20 0))) |
+-----+---------+-------------------------------------------------------------------------------+
| 3   | 1       | MULTIPOLYGON(((10 0,20 0,20 -20,10 -20,10 0)),((10 0,0 0,0 20,10 20,10 0)))   |
+-----+---------+-------------------------------------------------------------------------------+
| 3   | 2       | POLYGON((20 0,20 -20,10 -20,10 0,10 20,20 20,20 0))                           |
+-----+---------+-------------------------------------------------------------------------------+
| 3   | 4       | POLYGON((20 -20,10 -20,10 0,20 0,30 0,30 -20,20 -20))                         |
+-----+---------+-------------------------------------------------------------------------------+
| 4   | 2       | MULTIPOLYGON(((20 0,30 0,30 -20,20 -20,20 0)),((20 0,10 0,10 20,20 20,20 0))) |
+-----+---------+-------------------------------------------------------------------------------+
| 4   | 3       | POLYGON((20 0,30 0,30 -20,20 -20,10 -20,10 0,20 0))                           |
+-----+---------+-------------------------------------------------------------------------------+
| 4   | 5       | MULTIPOLYGON(((30 0,30 -20,20 -20,20 0,30 0)),((30 0,30 20,40 20,40 0,30 0))) |
+-----+---------+-------------------------------------------------------------------------------+
| 5   | 4       | MULTIPOLYGON(((30 0,30 20,40 20,40 0,30 0)),((30 0,30 -20,20 -20,20 0,30 0))) |
+-----+---------+-------------------------------------------------------------------------------+

Observe que el polígono id=3 comparte un punto con id=1 y, por tanto, se devuelve como resultado positivo. Si cambio la cláusula WHERE por ST_Touches( t1.geom, t2.geom ) AND t1.geom && t2.geom AND ST_Relate(t1.geom, t2.geom ,'T*T***T**'); No obtengo ningún registro.

1. Así que primero ¿Cómo puedo especificar ST_Relate para asegurarme de que sólo se tienen en cuenta las parcelas que comparten un segmento de línea?

2. Y entonces, ¿cómo fusionaría los polígonos 1,2,3,4 en una ronda, colapsando los resultados de la llamada anterior, todo ello reconociendo que la adyacencia 1 a 2 es la misma que la inversa?

Actualización

Si añado esto al where Obviamente sólo obtengo polígonos y no multipolígonos, eliminando así los falsos positivos para mis propósitos: los toques en las esquinas serán ignorados.

GeometryType(st_union(t1.geom,t2.geom)) != 'MULTIPOLYGON'

Aunque esto no es lo ideal (preferiría utilizar comprobaciones de topología con ST_RELATE como solución más general), es un camino a seguir. Entonces queda la cuestión de de-duping y union'ing estos. Posiblemente, si pudiera generar una secuencia sólo para los polígonos que se tocan, podría unirlos.

Actualización II

Esta parece funcionar para seleccionar polígonos que comparten líneas (pero no esquinas) y, por tanto, es una solución más general que la anterior MULTIPOLYGON test. Mi cláusula where tiene ahora este aspecto:

WHERE
              ST_Touches( t1.geom, t2.geom ) 
              AND t1.geom && t2.geom 

              -- 'overlap' relation
              AND ST_Relate(t1.geom, t2.geom)='FF2F11212') t2 

Ahora lo que queda es cómo hacer la fusión para más de un par de polígonos, sino para un número arbitrario que se ajuste a los criterios, de una sola vez.

Answer 1

No he podido evitar pensar que tu ejemplo es en realidad un raster y aunque has mencionado que te gustaría fusionar basándote en "ciertos criterios (podría ser el tamaño)" Me gustaría intentarlo con una conversión raster.

Para tu ejemplo concreto, esto funcionaría:

WITH rast AS (
  SELECT 
  ST_UNION(ST_AsRaster(geom,10, 20, '2BUI')) r
  FROM testpoly 
)
,p AS (
    SELECT (ST_DumpAsPolygons(r)).geom FROM rast
)
SELECT t.id,p.* 
FROM p
LEFT JOIN testpoly  t ON ST_Equals(p.geom, t.geom)

Lo que ocurre es que como tus polígonos son celdas perfectamente alineadas, se convertirán perfectamente en una trama (tamaño de celda 10x20). El dumpaspolygons le ayuda aquí mediante la fusión de todas las células adyacentes en una sola y mediante la comparación con los polígonos originales que incluso será capaz de obtener la identificación de nuevo por poly's no fusionados.

Una vez explicado esto, tengo mucha curiosidad por saber cómo se escalaría y qué tamaño tiene tu conjunto de datos :D

Answer 2

He aquí un ejemplo de cómo hacer esto en estilo procedimental con múltiples pasadas bajo el capó.

CREATE TABLE joined_testpoly AS SELECT array[id] ids, geom FROM testpoly; 

Debería poder transportar más columnas y aplicar criterios adicionales para la unión modificando la forma en que se utiliza la función LIMIT 1 seleccione a continuación las obras:

CREATE OR REPLACE FUNCTION reduce_joined_testpoly()
RETURNS void
AS $$
DECLARE
  joined_row joined_testpoly%ROWTYPE;
BEGIN
  LOOP
     SELECT array_cat(a.ids, b.ids), st_union(a.geom, b.geom)
         INTO joined_row 
     FROM joined_testpoly a INNER JOIN joined_testpoly b
           on a.ids != b.ids
              and ST_Touches(a.geom, b.geom) and a.geom && b.geom 
              and ST_Relate(a.geom, b.geom)='FF2F11212'
         LIMIT 1;
     IF NOT FOUND THEN
           EXIT;
     END IF;
     INSERT INTO joined_testpoly VALUES (joined_row.ids, joined_row.geom);
     DELETE FROM joined_testpoly
         WHERE joined_testpoly.ids <@ joined_row.ids 
           AND joined_testpoly.ids != joined_row.ids;
  END LOOP;
  RETURN;
END;
$$ LANGUAGE plpgsql;

Corre la cosa:

SELECT reduce_joined_testpoly();

Uniones adecuadas, sin multipolígonos:

SELECT ids, st_geometrytype(geom), st_area(geom), st_numgeometries(geom) 
FROM joined_testpoly;
    ids    | st_geometrytype | st_area | st_numgeometries 
-----------+-----------------+---------+------------------
 {5}       | ST_Polygon      |     200 |                1
 {1,2,3,4} | ST_Polygon      |     800 |                1

Answer 3

Aquí hay otro (no funciona) estrategia de referencia (que no pude llegar a excluir el caso de un solo punto de contacto). Debería ser más rápido que mi otra respuesta, ya que sólo toma 'una pasada'.

SELECT st_numgeometries(g), (SELECT st_union(x.geom) FROM st_dump(g) x GROUP BY g)
FROM (
    SELECT unnest(st_clusterintersecting(geom)) g, id < 100 as other_arbitrary_grouping 
    FROM testpoly
    GROUP BY other_arbitrary_grouping) c;

(siéntase libre de modificar y publicar otra respuesta si alguien puede conseguir la geometría id=5 en su propio grupo)

Para recuperar la lista de ids etc. tendrías que utilizar st_contains para volver a unir en la tabla testpoly como se detalla en la siguiente respuesta: https://stackoverflow.com/a/37486732/6691 pero por alguna razón no pude conseguir que funcionara para los polígonos.

Answer 4

He aquí un intento rápido utilizando tu consulta original ajustada un poco:

with gr as (SELECT 
    gid, adj_gid,
    st_AStext(st_union(l2.g1,l2.g2)) AS geo_combo
from (
    --level 2
    SELECT
      t1.id AS gid,
      t1.geom AS g1,
      t2.id AS adj_gid,
      t2.geom AS g2
     from
      testpoly  t1,
      testpoly  t2
     where
      ST_Touches( t1.geom, t2.geom ) 
      AND ST_Relate(t1.geom,t2.geom, '****1****')
      AND t1.geom && t2.geom 
) 
l2) select ST_AsText(st_union(gr.geo_combo)) from gr;

Referencias: https://postgis.net/docs/using_postgis_dbmanagement.html#DE-9IM