Fusionner tous les polygones adjacents

Je voudrais faire des tests d'adjacence sur une couche parcellaire (polygones) et les fusionner s'ils correspondent à certains critères (peut être la taille). Par l'image ci-dessous, je voudrais fusionner les polygones 1, 2, 3 et 4, mais pas 5.

J'ai deux problèmes:

1. ST_TOUCHESrenvoie VRAI si seulement les coins se touchent et non un segment de ligne. Je pense que j'ai besoin de ST_RELATE pour vérifier les segments de ligne partagés.
2. Idéalement, je voudrais fusionner TOUS les polygones adjacents en un seul, mais je ne sais pas comment évoluer au-delà de deux - comme dans, fusionner 1,2,3 et 4 (et peut-être plus sur les données réelles) en un seul tour.

La structure que j'ai maintenant est basée sur une auto-jointure ST_TOUCHES.

Données sur les jouets

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 ;

Sélection

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

Voici la sortie:

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

Notez que le polygone id = 3 partage un point avec id = 1 et est donc retourné comme résultat positif. Si je change la clause WHERE en ST_Touches( t1.geom, t2.geom ) AND t1.geom && t2.geom AND ST_Relate(t1.geom, t2.geom ,'T*T***T**');je n'obtiens aucun enregistrement.

1. Donc, tout d' abord , comment spécifier ST_Relate pour vous assurer que seules les parcelles partageant un segment de ligne sont prises en compte.

2. Et puis, comment fusionnerais-je les polygones 1, 2, 3, 4 en un tour, en regroupant les résultats de l'appel ci-dessus, tout en reconnaissant que la contiguïté 1 à 2 est identique à l'inverse?

Mise à jour

Si j'ajoute ceci à la whereclause, je n'obtiens évidemment que des polygones et non des multipolygones, éliminant ainsi les faux positifs à mes fins - les touches de coin seront ignorées.

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

Bien que ce ne soit pas idéal (je préfère utiliser les vérifications de topologie avec ST_RELATEcomme solution plus générale), c'est une voie à suivre. Reste alors la question de leur déduplication et de leur union. Peut-être, si je pouvais générer une séquence pour que seuls les polygones se touchent, je pourrais l'union là-dessus.

Mise à jour II

Celui-ci semble fonctionner pour sélectionner des polygones partageant des lignes (mais pas des coins) et est donc une solution plus générale que le MULTIPOLYGONtest ci-dessus . Ma clause where ressemble maintenant à ceci:

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

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

Maintenant, ce qui reste est de savoir comment faire la fusion pour plus qu'une paire de polygones, mais pour un nombre arbitraire correspondant aux critères, en une seule fois.

Je ne pouvais pas m'empêcher de penser que votre exemple est en fait un raster et bien que vous ayez mentionné que vous aimeriez fusionner en fonction de "certains critères (pourrait être la taille)", je voudrais lui donner un coup avec une conversion raster.

Pour votre exemple spécifique, cela fonctionnerait:

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)

Ce qui se passe, c'est que puisque vos polygones sont des cellules parfaitement alignées, elles se convertiront bien en un raster (taille de cellule 10x20). Le dumpaspolygons vous aide ici en fusionnant toutes les cellules adjacentes en une seule et en comparant avec les polygones d'origine, vous pourrez même récupérer l'identifiant pour les poly non fusionnés.

Après avoir expliqué cela, je suis très curieux de savoir comment cela évoluerait et quelle serait la taille de votre ensemble de données: D

Voici un exemple de la façon de procéder dans un style procédural avec plusieurs passes sous le capot.

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

Vous devriez pouvoir transporter plus de colonnes et appliquer des critères supplémentaires pour rejoindre en modifiant le fonctionnement de la LIMIT 1sélection ci-dessous:

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;

Exécutez la chose:

SELECT reduce_joined_testpoly();

Unions appropriées, pas de multipolygones:

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

Voici une autre stratégie (ne fonctionne pas) pour référence (que je n'ai pas pu exclure le cas du point de contact unique). Il devrait être plus rapide que mon autre réponse car il ne prend qu'un seul passage.

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;

(n'hésitez pas à modifier et à poster une autre réponse si quelqu'un peut obtenir la géométrie id = 5 dans son propre groupe)

Pour récupérer la liste des identifiants, etc., vous devez utiliser st_containspour rejoindre la table de testpoly comme détaillé dans la réponse suivante: /programming//a/37486732/6691 mais je n'ai pas pu faire fonctionner cela pour les polygones pour une raison quelconque.

Voici un aperçu rapide en utilisant votre requête d'origine un peu ajustée:

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;

Références: https://postgis.net/docs/using_postgis_dbmanagement.html#DE-9IM