Correction automatique des erreurs de géométrie avec PostGIS

Ces opérations sont réalisées pour des tables comportant un champs de type géométrie, stockées dans un SGBD PostGreSQL avec l’extension PostGIS installée.

1. Comptage des erreurs GEOS, des géométries nulles et des collections

SELECT 'invalid' AS nb, count(*) FROM my_table WHERE NOT ST_IsValid(the_geom)
UNION
SELECT 'geonul' AS nb, count(*) FROM my_table WHERE the_geom is null
UNION
SELECT 'collection' AS nb, count(*) FROM ma_table WHERE not ST_IsValid(the_geom)
AND ST_GeometryType(ST_MakeValid(the_geom))='ST_GeometryCollection';

Afficher les causes des erreurs

SELECT id, ST_IsValidReason(the_geom) FROM my_table WHERE NOT ST_IsValid(the_geom);

2. Correction des géométries invalides, des collections et des nœuds doubles (sur des couches de polygones)

2.1 Requête globale

UPDATE my_table
SET the_geom =
ST_Multi(ST_Simplify(ST_Multi(ST_CollectionExtract(ST_ForceCollection(ST_MakeVa
lid(the_geom)),3)),0))
WHERE ST_GeometryType(the_geom) = 'ST_MultiPolygon';

2.2 Correction seule de la géométrie de type GEOS

UPDATE my_table SET the_geom = ST_MakeValid(the_geom) WHERE NOT ST_IsValid(the_geom);

2.3 Suppression seule des géométries nulles

DELETE FROM my_table WHERE the_geom IS NULL;

A n’effectuer que si le nombre n’est pas nul à l’étape 1.

2.4 Correction seule des nœuds doubles

UPDATE my_table SET the_geom = ST_Multi(ST_Simplify(the_geom,0));

A n’effectuer que si la liste des erreurs indique la présence de nœuds double et que l’étape 2 n’est pas réalisée.

2.5 Vérifier de nouveau le nombre d’erreurs de type GEOS

Publicités

Installer PostGIS sous Windows

Installer PostGIS :

Télécharger le bundle correspondant à la version installée de PostGreSQL :

Lancer l’exécutable en vérifiant que la version de PostGreSQL indiquée est la bonne.

Activer PostGIS :

Une fois installé, PostGIS doit être activé dans une base de données sous forme d’extension. Il ne faut pas l’activer pour la base postgres :

— Enable PostGIS (includes raster)
CREATE EXTENSION postgis;
— Enable Topology
CREATE EXTENSION postgis_topology;

Connaître la version de PostGIS :

Pour vérifier quelle version de PostGIS est installée, on peut utiliser la commande suivante dans pgAdmin :

SELECT PostGIS_full_version()

La requête doit normalement renvoyer un message de ce type si PostGIS est installé :

POSTGIS= »2.1.3 r12547″ GEOS= »3.4.2-CAPI-1.8.2 r3924″ PROJ= »Rel. 4.8.0, 6 March 2012″ GDAL= »GDAL 1.10.0, released 2013/04/24″ LIBXML= »2.7.8″ LIBJSON= »UNKNOWN » TOPOLOGY RASTER

ou afficher l’erreur suivante dans le cas contraire :

ERREUR: la fonction postgis_full_version() n’existe pas
LINE 1: SELECT PostGIS_full_version();

Sources :

Supprimer les scories à l’intérieur d’un polygon avec PostGIS

L’union de polygones contigus peut parfois générés des scories (lignes résiduelles, trous…) dans le polygone ou le multipolygone assemblé. L’utilisation des fonctions de PostGIS permet de résoudre ces problèmes en quelques étapes.

Dans le cas de polygones simples, les étapes 2 et 3 sont suffisantes.

L’utilisation des fonctions est décomposée en plusieurs étapes mais il est bien sûr possible de les combiner. De même, des tables ou des vues sont créées à chaque étape afin de vérifier visuellement les résultats dans un outil SIG bureautique.

1. Exploser un multipolygone

Si la géométrie à traiter est un multipolygone, il est nécessaire de l’exploser au préalable en plusieurs polygones à l’aide de la fonction ST_Dump :

CREATE TABLE tablename_st_dump AS
SELECT fieldname, (ST_Dump(the_geom)).geom AS the_geom 
FROM tablename;

2. Récupérer la polyligne extérieure d’un polygone

On utilise ensuite la fonction ST_ExteriorRing pour récupérer la polyligne extérieure du polygone ce qui permet de supprimer toutes les scories internes au polygone (c’est-à-dire n’intersectant pas cette ligne) :

CREATE TABLE tablename_st_exterior_ring AS 
SELECT fieldname, ST_ExteriorRing(the_geom) AS the_geom
FROM tablename_st_dump;

3. Transformer la polyligne en polygone

On utilise la fonction ST_MakePolygon pour transformer la polyligne récupérer en polygone :

CREATE TABLE tablename_st_make_polygon AS 
SELECT fieldname, ST_MakePolygon(the_geom) AS the_geom
FROM tablename_st_exterior_ring;

4. Assembler les polygones dans un multipolygone

Si la géométrie initiale était un multipolygone, on peut assembler les différents polygones  avec la fonction ST_Union :

CREATE TABLE tablename_st_make_polygon AS 
SELECT fieldname, ST_Union(the_geom) AS the_geom
FROM tablename_st_make_polygon
GROUP BY fieldname;

Sources :

Calculer le pourcentage d’intersection entre deux polygones

Dans le cas du calcul de l’intersection de deux couches de polygones commune et zonage pour lesquels le contour des polygones zonage peut être constitué d’un assemblage de polygones commune, la fonction ST_Intersects donne également dans le résultat tous les polygones commune touchant la limite d’un polygone zonage sans le couper.

Pour affiner le résultat, il est possible d’utiliser la requête suivante :

CREATE OR REPLACE VIEW my_schema.my_view AS
SELECT *
FROM my_schema.zonage zonage, my_schema.commune commune
WHERE st_intersects(zonage.the_geom, commune.the_geom)
AND (st_area(st_intersection(zonage.the_geom, commune.the_geom)) / st_area(commune.the_geom)) > 0.1::double precision;

Une vue est créée pour faciliter la réutilisation de la requête et notamment le paramétrage du pourcentage de surface minimum pour que le polygone commune soit prise en compte dans le résultat.

Source :