Como fazer um loop recursivo pelas interseções de polígonos pai para obter polígonos menores (filhos) sem sobreposições?

11

Estou lutando com um problema por alguns dias e percebi que muitas pessoas também ficam presas quando o tópico é interseções no PostGIS (v2.5). Por isso, decidi pedir uma pergunta comum mais detalhada e genérica.

Eu tenho a seguinte tabela:

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 com isso:

começar

Quero obter todos os polígonos filhos com base na interseção dos polígonos pai. Para o resultado, seria esperado:

  • Os polígonos filhos sem sobreposição entre eles.
  • Uma coluna contendo a soma do valor de seus polígonos principais,
  • Uma coluna contendo a contagem de polígonos pai de uma categoria
  • Uma coluna contendo a contagem de outra categoria
  • Uma coluna contendo a categoria do polígono filho, com base na seguinte regra: -Se TODOS os polígonos pai são de uma classe, o polígono filho também possui essa classe. Senão, a categoria do polígono filho é uma terceira categoria.

Então parece:

resultado

Então, no fim, a tabela de saída gerado (para este exemplo) terá 7 linhas (todas as 7, polígonos, que não se sobrepõem criança), contendo colunas de category, sum_value, ct_overlap_cat1,ct_overlap_cat2

O código a seguir, iniciado, fornece as interseções individuais, comparando um pai com outro.

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);

Como faço um loop recursivo através do resultado desse código mencionado, que, independentemente do número de polígonos sobrepostos, sempre obtenho seus polígonos 'menores' (filho) (Fig. 2)?

Matt_Geo
fonte

Respostas:

8

Tente o seguinte:

Faça o download dos complementos do PostGIS neste link: https://github.com/pedrogit/postgisaddons

Instale executando o arquivo postgis_addons.sql para obter a função ST_SplitAgg ().

Teste executando o arquivo postgis_addons_test.sql.

Aqui está sua 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
fonte
Eu olhei para o seu addon git repo antes. Coisas inspiradoras.
19418 John Powell
Uau ótima solução. Você fez um trabalho incrível ao criar esses complementos também. Antes de clicar para conceder essa resposta, apenas uma para garantir que algo esteja me incomodando. Ao executar o código que você fornece, o 'polígono 5' (da segunda figura da pergunta) parece não reconhecer a sobreposição com outro polígono ('ct_overlap_cat2 = 1'; 'ct_overlap_cat2 = 0'). Portanto, esse polígono acaba sendo classificado como 'cat1' e 'sum = 2', em vez de 'cat3' e 'sum = 7'. Estou com um pouco de dificuldade para depurar esse pequeno problema. Você poderia me ajudar?
8898 MattelGeo #
11
O único problema com esta solução é que as instruções de caso são codificadas. Em princípio, provavelmente deveria ser capaz de lidar com um número arbitrário de categorias.
19419 John Powell
@Matt_Geo Recebo 6 polígonos na tabela resultante. O polígono do triângulo é dividido em três. Um com soma = 2, um com soma = 7 e um com soma = 8, como na figura desejada.
Pierre Racine
11
E se você substituir ST_Centroid (geom) por ST_Area (geom)?
Pierre Racine
1

Suponho que se você usar o tipo de geometria de polígono em vez de MultiPolygon, tudo se encaixará:

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 );

O resultado são 9 entradas que correspondem às diferentes opções de interseção no exemplo fornecido por você.

Cyril Mikhalchenko
fonte