Le portail francophone de la géomatique

meonais · Tue 17 January 2017 20:47

Bonjour,

Je m'excuse par avance de la question qui risque d'être très basique, mais je débute avec postgis même si j'ai quelques notions de sql (autoformation). Je n'ai pas eu le temps d'affiner mon autoformation sur postgis avec le suivi strict de tutoriels (hum...) je l'avoue...

J'ai une couche de multilinestring (cours d'eau a, x, y) que je souhaite couper par une couche de points (confluences .) pour obtenir de nouvelles entités multilinestring (tronçons de cours d'eau : 1, 2, 3, x, y) :

l l
x x
l l
-----a----.-----a-----.-------a----- ==> ---1---.-- -- --2 -- -- --. -----3--------
l l
y y
l l

J'ai déjà ma couche de points, générée avec "intersection de ligne" par QGis.

J'ai démarré le code suivant :

Code:

CREATE TABLE sometable AS 
with 

test as (select t.code_hydro as code, (t.geom) as mline from public."test" t),
confrm as (select (c.geom) as pt from public."conftest" c),

foo as (SELECT te.code as code1, te.mline as mli, co.pt as pts FROM test te, confrm co),
wktcut as (select (ST_Split(f.mli, f.pts)) as spl, code1 from foo f)
    

  SELECT w.code1, w.spl
        FROM wktcut w;

Mais je n'arrive pas à charger la nouvelle table dans QGIS (il n'y a a priori pas de « type spatial » pour la nouvelle couche créée) bien que la table existe avec le champ 'code'.

La couche cours d'eau multilinestring "test" est composée d'un champ identifiant unique (format text) et d'un champ geom (ainsi que d'autres champs "descriptifs")
La couche points "conftest" est composée d'un identifiant unique (format num), d'un champ geom et des champs id_coursdeau (issus de l'outil intersection de lignes)

Merci d'avance de votre aide pour le complément / correction de ce bout de code afin de bien découper les cours d'eau en plusieurs entités en fonction des zones confluences (avec les champs de la couche initiale si possible) et afin de pouvoir charger et visualiser cette couche dans Qgis.

Je suis sous windows XP, QGis 2.2, Postgresql 9.4 et Postgis 2.1

Dernière modification par meonais (Wed 18 January 2017 12:31)

MathieuR · Tue 24 January 2017 13:44

Bonjour,

Pour information, il est conseillé d'appliquer un st_snap sur les données avant d'utiliser st_split.

Voici le code SQL modifié, mais non testé.

st_dump permettra d'extraire les éléments de la geometrycollection, celle-ci résultant du st_split

Code:

Code:

create table line_split as (
select st_dump(st_split(st_snap(a.geom, pt.geom, 0.1), pt.geom)) geom, 
a.code_hydro, a.classe, a.toponyme, a.candidat
from test as a, 
(select (st_dump(geom)).geom from conftest) as pt)

Dernière modification par MathieuR (Tue 24 January 2017 13:44)

meonais · Wed 28 June 2017 10:41

Bonjour,

Je reviens sur mon post quelques mois après. J'ai travaillé sur plusieurs choses mais je bloque encore sur ce sujet, malgré des recherches complémentaires, par exemple :
-- https://georezo.net/forum/viewtopic.php?id=84911
-- ou https://georezo.net/forum/viewtopic.php?id=102330
-- et https://georezo.net/forum/viewtopic.php?id=106691

Je suis un peu perdue dans les st_...

J'ai besoin de l'ensemble des points amont-aval de chaque ligne et des points le long de la ligne pour effectuer une requête complémentaire, mais seulement après avoir découpé mes lignes entre chaque points (ceux aux extrémités amont-aval et ceux le long).

Donc pour cette découpe je travaille aujourd'hui avec :

-- une couche de points le long des lignes, snappés

Code:

DROP TABLE IF EXISTS ptsnap CASCADE;
CREATE TEMP TABLE ptsnap AS(
SELECT o.cd as id, r.classement as za, r.eu as eu, r.__gid as idesp, ST_Closestpoint(ST_Collect(r.geom), o.geom) AS geom
FROM public."couche_points" o
INNER JOIN public."couche_lignes" r
ON ST_Dwithin(r.geom, o.geom, 50)
GROUP BY o.cd, r.classement, r.eu, r.__gid, o.geom);

Nota : je n'arrive pas à faire fonctionner st_snap correctement, au contraire du code ci-dessus

-- une autre couche de points des extrémités de lignes obtenues par UNION de deux tables :

Code:

DROP TABLE IF EXISTS pt_extrem_ligne CASCADE;
CREATE TEMP  TABLE pt_extrem_ligne AS (
(SELECT  'd'||gid::TEXT as idpt,
__gid as id_esp,
st_startpoint((ST_Dump(couche_lignes.geom)).geom) as pt_geom 
FROM couche_lignes)
UNION
(SELECT  'f'||gid::TEXT as idpt,
__gid as id_esp,
st_endpoint((ST_Dump(couche_lignes.geom)).geom) as pt_geom 
FROM couche_lignes)

-- l'UNION de ces deux couches de points (points le longs et points aux extrémités)

-- la couche de lignes couche_lignes précédemment citée

Quelle que soit l'adaptation de la requête trouvée, j'obtiens :
-- soit la même ligne (même longueur) superposée autant de fois que le nombre de couple de points (=de découpe ?) possible
-- soit une couche geometrycollection que je n'arrive pas à afficher dans QGis.

Merci de votre retour, même juste pour m'aiguiller !

Je suis sous windows 8.1, QGis 2.18.7, Postgresql 9.5 et Postgis 2.3.2

Nicolas Ribot · Wed 28 June 2017 11:36

Bonjour,

Je fais un peu une fixette dessus , mais pour ces pb, jouez un peu avec les fonctions de référencement linéaires plutôt qu'avec du snap ou du découpage par intersection.
Le principe général: voyez vos lignes comme une séquence de points: vous pourrez alors à loisirs ajouter des points, en supprimer, faire plusieurs lignes à partir d'une seule, etc.

Nicolas

meonais · Mon 03 July 2017 13:54

Bonjour,

J'ai résolu mon problème, merci pour l'axe de réflexion, ça m'a aidé à comprendre votre code :-) et à l'appliquer.

Par contre, j'ai créé une couche de lignes dumped avant pour pouvoir travailler. Et j'ai supprimé la couche de points aux extrémités de mes lignes (puisqu'ils sont intégrés dans les locus si j'ai bien compris) sinon la découpe ne se faisait pas au début de chaque ligne.

Donc dans l'ordre :

Dump de la couche multilignes :

Code:

DROP TABLE IF EXISTS lignes_dumped CASCADE;
CREATE  TABLE lignes_dumped AS
SELECT gid ,  __gid ,   name  ,  code ,   (ST_Dump(geom)).geom as geom
FROM lignes;

Création de la couche de points, snappés sur les lignes :

Code:

DROP TABLE IF EXISTS ptsnap CASCADE;
CREATE  TABLE ptsnapAS(
SELECT o.cd as id,  r.__gid as zid, ST_Closestpoint(ST_Collect(r.geom), o.geom) AS geom
FROM public."pt" o
INNER JOIN public."lignes" r
ON ST_Dwithin(r.geom, o.geom, 50)
GROUP BY o.cd, r.__gid, o.geom);

Je travaille sur des cours d'eau, y compris grands cours d'eau. Donc certains points sont géolocalisés au delà de 50m rive droite/gauche de la ligne du cours d'eau. J'ai donc ajouté à cette couche ceux qui intersectent une couche surfacique (hydrographie_surfacique de la bd carthage).

Enfin création d'une nouvelle couche de lignes découpées par la localisation des points le long de la ligne (Merci Nicolas Ribot !). :

Code:

-- premiere table des positions des points par rapport a "leurs" lignes
DROP TABLE IF EXISTS locus CASCADE;
CREATE TEMP TABLE locus as (
        -- on genere la position du debut de la ligne (index 0)  pour le premier segment à découper
    select l.gid, 0 as l
    from lignes_dumped l

    UNION ALL

        -- on genere la position de la fin de la ligne (index 1)  pour le dernier segment à découper
    select l.gid, 1 as l
    from lignes_dumped l

    UNION ALL

        -- calcul des positions des points par rapport aux lignes
    select l.gid, st_linelocatepoint(l.geom, p.geom) as  l
    from lignes_dumped l, ptsnap p
    where st_dwithin(l.geom, p.geom, 0.01)
    
    order by gid, l
);
-- deuxieme table qui ajoute un index croissant par position croissante, pour permettre de faire la jointure plus tard sur elle-meme:
DROP  TABLE IF EXISTS loc_with_idx CASCADE;
CREATE TEMP TABLE
loc_with_idx as (
    select gid, l, rank() over (partition by gid order by l) as idx
    from locus
) ;

-- decoupage final de chaque ligne avec les positions des points sur la ligne, grace a st_line_substring--

DROP TABLE IF EXISTS splitang CASCADE;
CREATE TABLE splitang AS (select l.gid, l.__gid ,  l.name  ,    l.code, st_linesubstring(l.geom, loc1.l, loc2.l) as geom 
from loc_with_idx loc1, loc_with_idx loc2, lignes_dumped l
where l.gid = loc1.gid
and loc1.gid = loc2.gid
and loc2.idx = loc1.idx+1
GROUP BY l.gid, l.__gid ,  l.name   l.code, l.geom, loc1.l, loc2.l );

Encore merci

Dernière modification par meonais (Thu 06 July 2017 16:11)

meonais · Mon 03 July 2017 17:44

Re bonjour !

En fait je suis allée un peu vite...
Globalement la requête fonctionne mais j'ai quelques trous dans la découpe : la couche "splitang" est tronquée de quelques tronçons, pourtant présents dans la couche initiale (dumped) et sans qu'il y ait forcément de point de découpe sur le tracé...

J'ai regardé un peu le résultat des sous-requêtes et je vois l'anomalie sans arriver à la comprendre : certains locus se doublent à 0 ou 1 pour une même ligne (même gid)
Ce sera peut être plus parlant avec des extraits / exemples : le résultat des sous requêtes (locus et id_with_idx) et le résultat final (traits noirs) de la découpe sur le tracé initial (traits verts) dans QGis.

Auriez-vous une idée sur le pourquoi du comment ? :-)

Merci

Dernière modification par meonais (Mon 03 July 2017 17:55)

Fichier(s) joint(s) :
Pour accéder aux fichiers vous devez vous inscrire.

meonais · Mon 03 July 2017 18:58

Bonsoir,

J'ai trouvé les sources d'erreur :

1/ en faisant le dump de la couche de lignes, je gardais les anciens identifiants ==> chaque multiligne (une entité) dumpée devient plusieurs entités, donc avec le même identifiant (gid).
J'ai créé un nouveau champ identifiant SERIAL lors du dump. Ensuite j'applique la découpe sur ce champ

2/ j'ai également inséré un GROUP BY dans le calcul des positions de points par rapport aux lignes, afin de traiter les zones où un point est exactement confondu avec l'extrémité de ligne et/ou deux points se superposent exactement (probablement une résultante du snap).

Du coup la requête totale que j'applique est :

Dump de la couche de multilignes (avec id unique)

Code:

DROP TABLE IF EXISTS lignes_dumped CASCADE;
CREATE  TABLE lignes_dumped AS
SELECT gid ,  __gid ,   name  ,  code ,   (ST_Dump(geom)).geom as geom
FROM lignes;
ALTER TABLE l_gm_ang_l_rm_dumped ADD COLUMN id SERIAL PRIMARY KEY;

Création de la couche de points, snappés sur les lignes :

Code:

DROP TABLE IF EXISTS ptsnap CASCADE;
CREATE  TABLE ptsnapAS(
SELECT o.cd as id,  r.__gid as zid, ST_Closestpoint(ST_Collect(r.geom), o.geom) AS geom
FROM public."pt" o
INNER JOIN public."lignes" r
ON ST_Dwithin(r.geom, o.geom, 50)
GROUP BY o.cd, r.__gid, o.geom);

Enfin création d'une nouvelle couche de lignes découpées par la localisation des points le long de la ligne, avec un group by pour éviter les points (de découpe) superposés entre eux :

Code:

-- premiere table des positions des points par rapport a "leurs" lignes
DROP TABLE IF EXISTS locus CASCADE;
CREATE TEMP TABLE locus as (
        -- on genere la position du debut de la ligne (index 0)  pour le premier segment à découper
    select l.gid, 0 as l
    from lignes_dumped l

    UNION ALL

        -- on genere la position de la fin de la ligne (index 1)  pour le dernier segment à découper
    select l.gid, 1 as l
    from lignes_dumped l

    UNION ALL

        -- calcul des positions des points par rapport aux lignes
    select l.gid, st_linelocatepoint(l.geom, p.geom) as  l
    from lignes_dumped l, ptsnap p
    where st_dwithin(l.geom, p.geom, 0.01)
    GROUP BY l.id, l
    order by gid, l
);
-- deuxieme table qui ajoute un index croissant par position croissante, pour permettre de faire la jointure plus tard sur elle-meme:
DROP  TABLE IF EXISTS loc_with_idx CASCADE;
CREATE TEMP TABLE
loc_with_idx as (
    select gid, l, rank() over (partition by gid order by l) as idx
    from locus
) ;

-- decoupage final de chaque ligne avec les positions des points sur la ligne, grace a st_line_substring--

DROP TABLE IF EXISTS splitang CASCADE;
CREATE TABLE splitang AS (select l.gid, l.__gid ,  l.name  ,    l.code, st_linesubstring(l.geom, loc1.l, loc2.l) as geom 
from loc_with_idx loc1, loc_with_idx loc2, lignes_dumped l
where l.gid = loc1.gid
and loc1.gid = loc2.gid
and loc2.idx = loc1.idx+1
GROUP BY l.gid, l.__gid ,  l.name   l.code, l.geom, loc1.l, loc2.l );

A bientôt

Dernière modification par meonais (Thu 06 July 2017 16:19)

Communauté	Ressources	L'emploi	Le marché	L'association	Fils RSS
Forum Blogs Wiki Planet Sites amis	Agenda Biblio Fiches AMO-CNIG Paris PCGIS Téléchargements	Offres d'emploi Banque de CV Formations Démarche Métiers	Geo-Entreprises Geo-Communiqué Appels d'offres	L'association Qui sommes-nous ? Nous aider Bienfaiteurs Partenaires Historique Statuts	Tout le site Job Géomatique Services GeoBlogs ...	Faire un don Mentions légales

Le portail francophone de la géomatique

Annonce

#1 Tue 17 January 2017 20:47

Découper lignes par points

Code:

#2 Tue 24 January 2017 13:44

Re: Découper lignes par points

Code:

#3 Wed 28 June 2017 10:41

Re: Découper lignes par points

Code:

Code:

#4 Wed 28 June 2017 11:36

Re: Découper lignes par points

#5 Mon 03 July 2017 13:54

Re: Découper lignes par points

Code:

Code:

Code:

#6 Mon 03 July 2017 17:44

Re: Découper lignes par points

#7 Mon 03 July 2017 18:58

Re: Découper lignes par points

Code:

Code:

Code:

Pied de page des forums