PostGIS: corriger les duplicate rings🔗

Posted by Médéric Ribreux 🗓 In blog/ Qgis/

#gis #postgis #sql

PostGIS est plus à cheval sur la forme et la géométrie des objets géographiques comparé à d'autres SIG. Ainsi, il est possible de récupérer des couches "mal formées" issus de logiciels de SIG moins regardants sur la topologie. La conséquence immédiate est l'impossibilité d'effectuer certains traitements avec PostGIS. Sur ces objets, qu'on arrive parfaitement à afficher par ailleurs, on ne peut effectuer d'interrogations complexes et les opérations de calculs ne fonctionnent pas correctement (ST_Area retourne 0).

Pour vérifier qu'un objet est respectueux du modèle géométrique de PostGIS, on peut utiliser la fonction ST_IsValid(champ_géométrique) qui retourne false en cas de problème. Sous la version PostGIS de Debian Squeeze (1.3.6), lorsqu'un objet est mal formé, une erreur est générée indiquant le type d'erreur et son emplacement géographique.

Du fait de la diversité des problèmes de coordonnées géométriques, il est impossible de trouver un traitement global et automatisé de correction. Toutefois, il est bon de disposer d'un moyen de récupérer des objets corrects même légèrement modifiés.

Le code qui suit est ma tentative balbutiante de création d'une fonction de correction des duplicates rings. Le problème corrigé est très précis. Un "Duplicate Rings" est généré lorsqu'un objet de type multi-polygone contient des doublons dans sa collection de doublons. Ça arrive par exemple si un multi-polygone contient 2 polygones strictement identiques.

Voici le code:

CREATE OR REPLACE FUNCTION corrige_duplicate_rings(geometry)
RETURNS geometry AS $$
DECLARE
  inGeom ALIAS for $1;
  outGeom geometry;
  i integer := 0;
  j integer := 0;
  geom_buffer geometry[];
  doublon boolean := false;
BEGIN
  -- Par défaut, la géométrie en sortie sera équivalente à celle qui doit être traitée
  outGeom := inGeom;
  -- on ne corrige que les polygones multiples et non valides !:
  IF (ST_GeometryType(inGeom) = 'ST_MultiPolygon') THEN
    IF NOT ST_isValid(inGeom) THEN
	  -- Extraction des polygones simples de la géométrie dans un tableau:
      FOR i IN 1..ST_NumGeometries(inGeom) LOOP
        geom_buffer[i] := ST_GeometryN(inGeom,i);
      END LOOP;
      -- Boucle de traitement:
      i:=1;
      WHILE i <= array_upper(geom_buffer,1) LOOP
        -- Recherche de doublons:
        j:=1;
        WHILE j <= array_upper(geom_buffer,1) LOOP
          -- Lorsqu'on a un doublon:
          IF NOT j = i THEN
            IF ST_Equals(geom_buffer[i],geom_buffer[j]) THEN
              doublon := true;
              geom_buffer := geom_buffer[1:j-1] || geom_buffer[j+1:array_upper(geom_buffer,1)];
            END IF;
          END IF;
        j:=j+1;    
        END LOOP;
        i:=i+1;
      END LOOP;

      -- A la seule condition d'une présence de doublon, on reconstruit:
      IF doublon THEN
        outGeom := st_collect_garray(geom_buffer);
      END IF;
    END IF;
  END IF;

  RETURN outGeom;

END;
$$ LANGUAGE plpgsql;

À noter que je débute en PL/PGSQL et qu'il est fort possible que la fonction soit perfectible. Mais ce code peut sans doute en inspirer d'autres plus complets et performants… J'ai testé la fonction sur une trentaine de couches mal formées et elle a bien répondu à chaque fois.

La fonction ne fait que corriger les doublons pour ce type d'objet et seulement si elle en rencontre. Un exemple d'utilisation:

SELECT ST_IsValid(corrige_duplicate_rings(wkb_geometry)) FROM "parcelles";

Pour mettre à jour une table:

UPDATE "parcelles" SET wkb_geometry=corrige_duplicate_rings(wkb_geometry) WHERE ST_IsValid(wkb_geometry) = FALSE;