2 Sep
2020
2 Sep
'20
2:37 p.m.
Message posté par : hadri45 (hadrien.siney(a)me.com)
----------------------------------------
Bonjour,
Je dois réaliser des croisements de couches forestières afin d'enrichir la base de
données des Forêts officielle. Pour cette opération, le mieux est d'utiliser la
fonction Union d'ArcGIS mais j'ai l'obligation de réaliser cela sous
Postgres/PostGIS. Pour se faire, j'ai donc vu qu'il fallait réaliser cela en 2
étapes avec d'abord l'intersection pour les deux couches qui s'intersectent
puis avec INSERT INTO, récupérer celle qui ne s'intersecte pas avec le ST_DIFFERENCE.
Voici le code où je dois réaliser deux croisements tout en gardant les attributs
(code_psg,code_onf et tfv) pour des analyses complémentaires.
Ce code me donne le bon résultat mais vu que je dois réaliser ce traitement pour la France
entière il me faut l'optimiser le plus possible sachant que sur un département cela
prend plus de 2h30 alors que sur ArcGIS le traitement prend 4 minutes.
-----------------
Code :
-- TABLE DU CROISEMENT PRO / PSG
DROP TABLE IF EXISTS surface_unitaire_01.propsg_01;
-- 1. Intersection
CREATE TABLE surface_unitaire_01.propsg_01 AS
SELECT o.code_onf, p.code AS code_psg, ST_SetSRID((ST_Dump(ST_Intersection(o.geom,
p.geom))).geom, 932006) AS geom
FROM surface_unitaire_01.pro_2018_simple_01 o
INNER JOIN surface_unitaire_01.psg_2018_simple_01 p ON (o.geom && p.geom AND
ST_Intersects(o.geom, p.geom));
-- Nettoyage des géométries invalides
DELETE FROM surface_unitaire_01.propsg_01 WHERE ST_GeometryType(geom) !=
'ST_Polygon';
CREATE INDEX sp_idx_propsg_01_geom ON surface_unitaire_01.propsg_01 USING GIST (geom);
ANALYZE surface_unitaire_01.propsg_01;
-- 2. Différences PRO - PSG et PSG - PRO
CREATE TEMP TABLE pro_union AS
SELECT code_onf, ST_Union(geom) AS geom
FROM surface_unitaire_01.propsg_01
GROUP BY code_onf;
CREATE INDEX sp_idx_pro_union_geom ON pro_union USING GIST (geom);
ANALYZE pro_union;
CREATE TEMP TABLE psg_union AS
SELECT code_psg, ST_Union(geom) AS geom
FROM surface_unitaire_01.propsg_01
GROUP BY code_psg;
CREATE INDEX sp_idx_psg_union_geom ON psg_union USING GIST (geom);
ANALYZE psg_union;
DROP INDEX surface_unitaire_01.sp_idx_propsg_01_geom;
INSERT INTO surface_unitaire_01.propsg_01
SELECT o.code_onf, NULL::VARCHAR(3) AS code_psg, (ST_Dump(ST_Difference(o.geom,
ST_SetSRID(COALESCE(t1.geom, 'GEOMETRYCOLLECTION EMPTY'::GEOMETRY),
932006)))).geom AS geom
FROM surface_unitaire_01.pro_2018_simple_01 o
LEFT JOIN pro_union t1 ON o.code_onf = t1.code_onf
UNION
SELECT NULL::VARCHAR(3) AS code_onf, p.code AS code_psg, (ST_Dump(ST_Difference(p.geom,
ST_SetSRID(COALESCE(t2.geom, 'GEOMETRYCOLLECTION EMPTY'::GEOMETRY),
932006)))).geom AS geom
FROM surface_unitaire_01.psg_2018_simple_01 p
LEFT JOIN psg_union t2 ON p.code = t2.code_psg;
-- Nettoyage des géométries invalides
SELECT DISTINCT ST_GeometryType(geom) FROM surface_unitaire_01.propsg_01;
SELECT DISTINCT ST_SRID(geom) FROM surface_unitaire_01.propsg_01;
SELECT COUNT(*) FROM surface_unitaire_01.propsg_01 WHERE NOT ST_IsValid(geom);
DELETE FROM surface_unitaire_01.propsg_01 WHERE ST_GeometryType(geom) !=
'ST_Polygon';
CREATE INDEX sp_idx_propsg_01_geom ON surface_unitaire_01.propsg_01 USING GIST (geom);
ANALYZE surface_unitaire_01.propsg_01;
DROP TABLE pro_union;
DROP TABLE psg_union;
-- contrôles de surfaces
SELECT 'couche onf' AS origine, code_onf, SUM(ST_Area(geom)) AS surface_total
FROM surface_unitaire_01.pro_2018_simple_01
GROUP BY origine, code_onf
UNION
SELECT 'croisement' AS origine, code_onf, SUM(ST_Area(geom)) AS surface_total
FROM surface_unitaire_01.propsg_01
GROUP BY origine, code_onf
ORDER BY code_onf, origine;
SELECT 'couche psg' AS origine, code AS code_psg, SUM(ST_Area(geom)) AS
surface_total
FROM surface_unitaire_01.psg_2018_simple_01
GROUP BY origine, code_psg
UNION
SELECT 'croisement' AS origine, code_psg, SUM(ST_Area(geom)) AS surface_total
FROM surface_unitaire_01.propsg_01
GROUP BY origine, code_psg
ORDER BY code_psg, origine;
SELECT DISTINCT st_geometrytype(geom)
FROM surface_unitaire_01.bd_foret_01_simple
SELECT id_foret, ST_NPoints(geom), ST_MemSize(geom)
FROM surface_unitaire_01.bd_foret_01_simple
ORDER BY ST_NPoints DESC;
SELECT COUNT(*), COUNT(*) FILTER (WHERE ST_MemSize(geom) > 8192)
FROM surface_unitaire_01.bd_foret_01_simple;
CREATE TABLE surface_unitaire_01.bd_foret_01_simple_sub AS
SELECT id_foret, tfv, ST_SubDivide(geom) AS geom
FROM surface_unitaire_01.bd_foret_01_simple;
CREATE INDEX idx_foret_geom_sub_01 ON surface_unitaire_01.bd_foret_01_simple_sub USING
GIST (geom);
CREATE TABLE surface_unitaire_01.propsg_01_sub
AS SELECT ST_SUBDIVIDE(geom) as geom, code_onf,code_psg
FROM surface_unitaire_01.propsg_01
-- TABLE DU CROISEMENT BD_FORET / PRO_PSG
DROP TABLE IF EXISTS surface_unitaire_01.final_cedric_01;
-- 1. Intersection
CREATE TABLE surface_unitaire_01.final_cedric_01 AS
SELECT p.code_onf, p.code_psg, b.tfv,
ST_SetSRID((ST_Dump(ST_Intersection(ST_Buffer(p.geom, 0.0000001), b.geom))).geom, 932006)
AS geom
FROM surface_unitaire_01.propsg_01 p
INNER JOIN surface_unitaire_01.bd_foret_01_simple_sub b ON (p.geom && b.geom
AND ST_Intersects(p.geom, b.geom));
-- Nettoyage des géométries invalides
SELECT DISTINCT ST_GeometryType(geom) FROM surface_unitaire_01.final_cedric_01;
DELETE FROM surface_unitaire_01.final_cedric_01 WHERE ST_GeometryType(geom) !=
'ST_Polygon';
CREATE INDEX sp_idx_final_cedric_01_geom ON surface_unitaire_01.final_cedric_01 USING GIST
(geom);
ANALYZE surface_unitaire_01.final_cedric_01;
-- 2. Différences BD_FORET - PRO_PSG et PRO_PSG - BD_FORET
DROP TABLE IF EXISTS pro_psg;
CREATE TEMP TABLE pro_psg AS
SELECT code_onf, code_psg,geom
FROM surface_unitaire_01.final_cedric_01;
CREATE INDEX sp_idx_pro_psg_geom ON pro_psg USING GIST (geom);
ANALYZE pro_psg;
DROP TABLE IF EXISTS tfv;
CREATE TEMP TABLE tfv AS
SELECT tfv, geom
FROM surface_unitaire_01.final_cedric_01;
CREATE INDEX sp_idx_tfv_geom ON tfv USING GIST (geom);
ANALYZE tfv;
DROP INDEX surface_unitaire_01.sp_idx_final_cedric_01_geom;
INSERT INTO surface_unitaire_01.final_cedric_01
SELECT p.code_onf, p.code_psg, NULL::VARCHAR(10) AS tfv, (ST_Dump(ST_Difference(p.geom,
ST_SetSRID(COALESCE(t1.geom, 'GEOMETRYCOLLECTION EMPTY'::GEOMETRY),
932006)))).geom AS geom
FROM surface_unitaire_01.propsg_01_sub p
LEFT JOIN pro_psg as t1 ON p.code_onf = t1.code_onf AND p.code_psg = t1.code_psg
UNION
SELECT NULL::VARCHAR(3) AS code_onf, NULL::VARCHAR(3) AS code_psg, t.tfv,
(ST_Dump(ST_Difference(t.geom, ST_SetSRID(COALESCE(t2.geom, 'GEOMETRYCOLLECTION
EMPTY'::GEOMETRY), 932006)))).geom AS geom
FROM surface_unitaire_01.bd_foret_01_simple_sub t
LEFT JOIN tfv t2 ON t.tfv = t2.tfv;
CREATE INDEX sp_idx_final_cedric_01_geom ON surface_unitaire_01.final_cedric_01 USING GIST
(geom);
ANALYZE surface_unitaire_01.propsg_01;
DROP TABLE pro_psg_union;
DROP TABLE IF EXISTS tfv_union;
DROP TABLE tfv;
-----------------
J'ai donc essayé de subdiviser ce qui accélère les temps de traitements mais le plus
gros problème reste le dernier ST_DIFFERENCE entre la bd_foret et pro_psg qui est toujours
aussi long malgré cette subdivision et ne permet pas de valider le code. J'ai
également essayé de transformer les géométries en tuile (MapVector Tile) avec ST_ASMVTGeom
mais cela est encore plus long. J'ai également essayé de faire sans le ST_UNION ou
ST_COLLECT mais cela ne change pas beaucoup le temps de traitement.
Ma couche principale est de l'ordre de plus de 25000 entités et les autres couches à
croiser sont autour de 2000.
Pour la configuration, je suis sous PostgreSQL 12 et PostGIS 3.0.
Auriez-vous donc une idée afin d'optimiser le traitement ?
Merci,
----------------------------------------
Le message est situé https://georezo.net/forum/viewtopic.php?pid=335100#p335100
Pour y répondre : geomatique(a)ml.georezo.net ou reply de votre messagerie
Pour vous désabonner connectez-vous sur le forum puis Profil / Abonnement
--
Association GeoRezo - le portail géomatique
https://georezo.net