¿Cómo recorrer recursivamente las intersecciones de polígonos primarios para obtener los polígonos más pequeños (secundarios) sin superposiciones?

11

Estoy luchando con un problema durante un par de días y me di cuenta de que muchas personas también se estancan cuando el tema es las intersecciones en PostGIS (v2.5). Es por eso que decidí pedir una pregunta común más detallada y genérica.

Tengo la siguiente tabla:

DROP TABLE IF EXISTS tbl_foo;
CREATE TABLE tbl_foo (
    id bigint NOT NULL,
    geom public.geometry(MultiPolygon, 4326),
    att_category character varying(15),
    att_value integer
);
INSERT INTO tbl_foo (id, geom, att_category, att_value) VALUES 
    (1, ST_SetSRID('MULTIPOLYGON (((0 6, 0 12, 8 9, 0 6)))'::geometry,4326) , 'cat1', 2 );
INSERT INTO tbl_foo (id, geom, att_category, att_value) VALUES 
    (2, ST_SetSRID('MULTIPOLYGON (((5 0, 5 12, 9 12, 9 0, 5 0)))'::geometry,4326), 'cat1', 1 );
INSERT INTO tbl_foo (id, geom, att_category, att_value) VALUES 
    (3, ST_SetSRID('MULTIPOLYGON (((4 4, 3 8, 4 12, 7 14,10 12, 11 8, 10 4, 4 4)))'::geometry,4326) , 'cat2', 5 );

Se parece a esto:

comienzo

Quiero obtener todos los polígonos secundarios basados ​​en la intersección de los polígonos principales. Para el resultado, se esperaría:

  • Los polígonos secundarios sin superposición entre ellos.
  • Una columna que contiene la suma del valor de sus polígonos principales,
  • Una columna que contiene el recuento de polígonos primarios de una categoría.
  • Una columna que contiene el recuento de otra categoría.
  • Una columna que contiene la categoría del polígono secundario, basada en la siguiente regla: -Si TODOS los polígonos principales son de una clase, el polígono secundario también tiene esta clase. De lo contrario, la categoría del polígono hijo es una tercera categoría.

Así se ve así:

salida

Por lo tanto, al final, la tabla de salida generados (para este ejemplo) tendrá 7 filas (todos los 7,, niño polígonos no se solapan), que contienen columnas de category, sum_value, ct_overlap_cat1,ct_overlap_cat2

El siguiente código que comencé, me da las intersecciones individuales, comparando un padre con otro.

SELECT
(ST_Dump(
    ST_SymDifference(a.geom, b.geom) 
)).geom
FROM tbl_foo a, tbl_foo b
WHERE a.ID < b.ID AND ST_INTERSECTS(a.geom, b.geom)
UNION ALL
SELECT
ST_Intersection(a.geom, b.geom) as geom
FROM tbl_foo a, tbl_foo b
WHERE a.ID < b.ID AND ST_INTERSECTS(a.geom, b.geom);

¿Cómo recorro recursivamente el resultado de este código mencionado, que, independientemente del número de polígonos superpuestos, siempre obtengo sus polígonos 'más pequeños' (secundarios) (Fig. 2)?

Matt_Geo
fuente

Respuestas:

8

Prueba esto:

Descargue los complementos PostGIS desde este enlace: https://github.com/pedrogit/postgisaddons

Instale ejecutando el archivo postgis_addons.sql para obtener la función ST_SplitAgg ().

Pruebe ejecutando el archivo postgis_addons_test.sql.

Aquí está su consulta:

WITH  result_table AS (
    WITH  parts AS (
      SELECT a.att_value val,
             CASE WHEN a.att_category = 'cat1' THEN 1 ELSE 0 END cat1,
             CASE WHEN a.att_category = 'cat2' THEN 1 ELSE 0 END cat2,
             unnest(ST_SplitAgg(a.geom, b.geom, 0.00001)) geom
      FROM tbl_foo a,
           tbl_foo b
      WHERE ST_Equals(a.geom, b.geom) OR
            ST_Contains(a.geom, b.geom) OR
            ST_Contains(b.geom, a.geom) OR
            ST_Overlaps(a.geom, b.geom)
      GROUP BY a.id, a.att_category , ST_AsEWKB(a.geom), val
    )
    SELECT CASE WHEN sum(cat2) = 0 THEN 'cat1'
                WHEN sum(cat1) = 0 THEN 'cat2'
                ELSE 'cat3'
           END category, 
           sum(val*1.0) sum_value, 
           sum(cat1) ct_overlap_cat1, 
           sum(cat2) ct_overlap_cat2, 
           ST_Union(geom) geom
    FROM parts
    GROUP BY ST_Area(geom)
)
SELECT category, sum_value, ct_overlap_cat1, ct_overlap_cat2,
(ST_Dump(result_table.geom)).geom as geom
FROM result_table
Pierre Racine
fuente
He visto tus complementos git repo antes. Cosas inspiradoras.
John Powell el
Wow gran solución. También hiciste un trabajo increíble para crear estos complementos. Antes de hacer clic para otorgar esta respuesta, solo uno para asegurarme de una cosa que me está molestando. Al ejecutar el código que proporciona, el 'polígono 5' (de la segunda figura de la pregunta) no parece reconocer la superposición con otro polígono ('ct_overlap_cat2 = 1'; 'ct_overlap_cat2 = 0'). Por lo tanto, este polígono termina clasificado como 'cat1' y 'sum = 2', en lugar de 'cat3' y 'sum = 7'. Me enfrento a un poco de dificultad para depurar este pequeño problema. ¿Usted me podría ayudar?
Matt_Geo
1
El único problema con esta solución es que las declaraciones de casos están codificadas. En principio, probablemente debería poder manejar un número arbitrario de categorías.
John Powell el
@Matt_Geo Obtengo 6 polígonos en la tabla resultante. El polígono triangular se divide en tres. Uno con suma = 2, uno con suma = 7 y uno con suma = 8 como en su figura deseada.
Pierre Racine
1
¿Qué sucede si reemplaza ST_Centroid (geom) con ST_Area (geom)?
Pierre Racine
1

Supongo que si usa el tipo de geometría de polígono en lugar de MultiPolygon, todo encajará:

DROP TABLE IF EXISTS tbl_foo;
CREATE TABLE tbl_foo (
    id bigint NOT NULL,
    geom public.geometry(Polygon, 4326),
    att_category character varying(15),
    att_value integer
);

INSERT INTO tbl_foo (id, geom, att_category, att_value) VALUES 
    (1, ST_SetSRID('POLYGON ((0 6, 0 12, 8 9, 0 6))'::geometry,4326) , 'cat1', 2 );
INSERT INTO tbl_foo (id, geom, att_category, att_value) VALUES 
    (2, ST_SetSRID('POLYGON ((5 0, 5 12, 9 12, 9 0, 5 0))'::geometry,4326), 'cat1', 1 );
INSERT INTO tbl_foo (id, geom, att_category, att_value) VALUES 
    (3, ST_SetSRID('POLYGON ((4 4, 3 8, 4 12, 7 14,10 12, 11 8, 10 4, 4 4))'::geometry,4326) , 'cat2', 5 );

El resultado son 9 entradas que corresponden a las diferentes opciones de intersección en el ejemplo dado por usted.

Cyril Mikhalchenko
fuente