[QGis 2.x – Spatialite 4.x] Reprojeter des points sur un axe.

Objectif : Suite à la création d’un profil d’après un MNT (voir https://archeomatic.wordpress.com/2013/12/12/qgis-2-0-1-profile-tool-points2one-creer-une-coupe-dapres-un-mnt/) Nous disposons d’un axe de coupe et de mobilier archéologique côté sous forme de points en X,Y et Z. Nous allons « reprojeter » ces points selon l’axe de coupe longitudinal de la fosse. 

 

Reprojection du mobilier archéologique (points) sur un axe de coupe.

La solution proposée ici consiste a reporter les points sur l’axe de coupe, a calculer leur distance au point 0 (ici le point A) puis a les afficher dans un repère en X, Z (X étant la distance au point 0 et Z, l’altitude du point).

 

Nous allons pour cela utiliser les opérateurs spatiaux SQL de Spatialite (voir la présentation du logiciel dans ce précédent article du blog : http://wp.me/p29w2o-4I ) et notamment st_line_locate_point et st_point.

 

1)      Postulat :

 

Nous disposons dans notre projet Qgis :

–          d’un shapefile de points qui représente le mobilier archéologique avec une valeur Z en attribut.

–          D’un shapefile de ligne qui correspond à une coupe longitudinale. (Note : c’est ce même axe qui a servi a créer un profil d’après le MNT sous-jacent dans le tuto précédent).

 

2)      Créer une base de données Spatialite:

 

Ayant toujours des bugs avec Qspatialite je passe toujours par Spatialite lui-même pour créer une base de données spatialite… rien de plus simple il suffit d’ouvrir l’executable « spatialite_gui.exe) téléchargé ici : https://www.gaia-gis.it/fossil/spatialite_gui/index

 

Puis d’appuyer sur l’icône 

Creating a new empty SQlite DB.

 

Le logiciel nous propose alors de créer une nouvelle base de données que nous nommerons REPROJECTION.sqlite

 

 

 

Et voilà c’est fini pour Spatialite, nous pouvons quitter le logiciel et revenir à QGis pour…

 

 

 

3)      Se connecter à la base de données Spatialite et charger nos 2 couches :

 

Dans QGis donc, cliquer sur l’icône  , une fenêtre de dialogue s’ouvre alors :

Il suffit d’appuyer sur le bouton [Nouveau] de parcourir pour trouver notre base spatialite nouvellement créée puis d’appuyer sur le bouton [Connecter].

 

Tout le reste des manipulations se feront dorénavant à partir du gestionnaire de Base De Données integré à Qgis 

(dans la barre d’outils « Base de données, la bien nommée !).

 

La fenêtre de dialogue du gestionnaire de Base de données s’ouvre :

 Pour se connecter à notre base de données, il suffit de dérouler le menu spatialite  1 (appuyer sur le petit + ) puis de cliquer sur l’intitulé REPROJECTION.sqlite. Notre base de données est déjà prête et contient une « foultitude » de table par défaut… mais pas encore de données spatiales.

Pour cela nous allons cliquer sur l’icône « Importer une couche ou un fichier » 2

Il faudra ici :

–          Saisie : choisir notre couche (parmi celle présente dans Qgis) ici Mobilier.shp

–          Table : lui donner un nom : ici SP_Mobilier (le préfixe SP pour Spatialite qui permet de la differencier de l’originale)

–          Clé primaire : cocher et écrire gid

–          Colonne de géométrie : cocher et écrire geom

–          Créer des géométries simples au lieu de multipartities : cocher

–          Créer un index spatial : cocher (ça ne peut pas faire de mal !)

Appuyer sur le bouton [OK]

 

Et voilà une belle couche importée. Pour le vérifier, cliquer sur l’intitulé de la base de données puis rafraichir (icône entouré de rouge ci-dessous). On voit bien que notre couche de point SP_Mobilier a été importée :

On fera de même avec notre coupe Coupe.shp que l’on appellera SP_Coupe.

Une fois les couches importées on peut très simplement visualiser nos couche dans l’arborescence 1 et la table attributaire de nos couches en cliquant sur l’onglet « Table » 2.

On voit par exemple ci-dessus que nos points ont comme champs attributaires : gid (la clé primaire), geom (la colonne de geométrie), x,y et z (no comment), type et code (description du mobilier archéo) et fait (le numéro de fait archéologique qui contient le mobilier).

 

4)      La reprojection proprement dite…

 

… va se passer en 2 étapes …

Note :

– Boouuh avec PostGis on peut le faire en une fois !!

–   oui oui mathias je sais –et merci au passage pour ton aide indispensable- mais là pas besoin d’instalation client-serveur compliquée.. et toc !

et le tout dans la fenêtre SQL à laquelle on accède via l’icône  .

 

4a) Créer une table de points avec un « xrelatif » :

Nous donnons donc ici un nom 1 à notre première requête « Créer_Table_Point » puis après l’avoir tapée 3 le bouton [Stocker] 2 permettra de la réutiliser.

Avant d’appuyer sur le bouton [Executer] 4 et d’avoir vérifier que la requête SQL avait été traitée convenablement 5 attardons nous sur la requête :

 

create table T301 as

on crée une table T301 (T pour table et 301 car c’est le nom de la fosse) avec

Select p.gid, p. »z », p. »type », p. »code », p. »fait »,

les champs gid, z, type, code et fait…(notez les noms de champs avec les guillemets)

st_line_locate_point(c.geom, p.geom)*st_length(c.geom) as xrelatif

… et un champ appellé xrelatif qui correspond à la distance au point 0 de la coup, de chaque point préalablement reporté sur l’axe de coupe

from « SP_Coupe » c, « SP_Mobilier » p

et ce, à partir des couches SP_Coupe (c de son petit nom) et SP_Mobilier (p pour point de son petit nom)

where p. »fait » like ‘301’ and  c. »gid » = 1

avec comme condition que les points appartiennent au fait 301 et que la clé primaire de la coupe (gid) soit égale à 1.(notez la valeur texte encadrée d’apostrophes)

order by xrelatif ;

et trions tout par ordre croissant à partir du champs xrelatif !

 

Vous avez tout compris ? Vous êtes donc autorisé a [Executer] votre requête via le bouton idoine …choisir le SCR… (en passant si la requête a fonctionné appuyer sur le bouton [Stocker] pour la garder en mémoire !)…

…et a regarder la jolie table que vous venez de créer… après avoir actualisé votre base de donnée et cliqué sur l’onglet « Table » bien sûr !

 

Nous avons donc bien une table avec les champs voulus et un nouveau champ : « xrelatif » qui correspond pour chaque point reporté sur l’axe a sa distance par rapport au départ de l’axe (en mètre bien sûr, mon projet QGis et mes couches sont en RGF93/ CC47)

Ne reste plus qu’à…

 4b) Créer une couche de point en Xrelatif, Z :

Rien de plus simple (mais si vous allez voir…)

 

1 Donner un nom à notre requête, ici « Créer_Points_Projetés »

2 Ne pas appuyer sur le bouton [Stocker]… attendre d’avoir écrit la requête et vérifié qu »elle fonctionne proprement 😉

3 Ecrire la requête :

 

Par laquelle on demande

select st_point(xrelatif, z)

de créer les points avec comme abscisse les valeurs du champ xrelatif et en ordonnée les valeurs du champ z (attention pas de create table ici, juste un select)

as geom, gid, « z », « type », « code » , « fait », « xrelatif »

et d’y ajouter comme attributs, les champs geom, gid, z, type, code, fait et xrelatif

from T301

et ce, à partir de la table T301, que l’on vient de créer avec la requête précédente

order by xrelatif

et puis, tant qu’à faire de trier à partir du champ xrelatif

 

4 On [Exécute]

5 On regarde la table résultante…

6 On coche « Charger en tant que nouvelle couche »

7 choisir la clé primaire : gid

8 identifier la colonne qui contient la geométrie : geom

9 donner un nom à la couche, ici « Reprojection_Points »…

10 et on [Charge !]

Conclusion :

Alors, c’est-ti-pas-beau ? avec le profil en plus…

Cette méthode embarquée (pas de système client serveur) permet simplement de faire ce dont je rêvais depuis pas mal de temps… reprojeter des points sur un axe afin d’avoir une vision « en coupe » des altitudes de différents objets archéologiques. On peut évidemment choisir l’axe de reprojection à l’envie mais surtout cette méthode totalement intégrée au logiciel d’information géographique permet évidemment de garder tous les attributs des points d’origine et donc de les représenter avec la même symbologie en coupe et en plan.

Personnellement j’y vois déjà de multiples usages à des échelles différentes et à des fins différentes.. oui oui géomorphos de tous poils et responsables d’opérations des profondes vallées, c’est à vous que je parle !

En tout cas si ce tuto vous a servi, dites-m’en plus !

 

Epilogue :

Vous avez bien compris la combine qui se rapproche de celle utilisée dans le tuto précédent pour créer un profil selon un axe d’après un MNT, on triche puisque l’on crée des points dans un système Xrelatif, Z au lieu de X,Y.

Un des plus gros problème qui en résulte c’est que le composeur de Qgis refuse d’afficher les objets trop proches de 0,0 (Bug ?), pour palier a ce détail technique il suffit de changer la 3eme ligne de la première requête ainsi : (oui bah vous avez qu’à suivre un peu !)

st_line_locate_point(c.geom, p.geom)*st_length(c.geom)+1000000 as xrelatif