Tutoriel: Créer une carte du niveau d’eaux souterraines à partir de données de forage et d’un MNT

1. Introduction

2. Charger un ensemble de données de forage à partir d’un GeoNode SDI

Nous débuterons avec le chargement d’un ensemble de données de forage à partir d’un GeoNode SDI vers QGIS.

Ce tutoriel vous aidera à utiliser des bases de données de forage de l’aquifère transfrontalier de Stampriet à partir du SIG du bassin versant de Orange-Sengu, qui est un GeoNode SDI

2. Démarrer QGIS avec un nouveau projet

.

4. Choisir l’onglet GeoNode.

5. Cliquer sur nouveau bouton pour créer une nouvelle connexion de service.

6. Dans le dialogue de Créer un nouveau GeoNode taper ORASECOM pour le nom et copier http://gis.orasecom.org.

7. Cliquer sur tester la connexion.

Si le test est concluant vous verrez cette fenêtre apparaitre:

Si la connexion est un échec, vérifier votre connexion internet et l’URL.

8. Cliquer Ok pour fermer la fenêtre.

9. Cliquer Ok dans la boite de dialogue pour la fermer et la nouvelle connexion y est ajoutée.

10. Cliquer sur le bouton de connexion.

11. Vous recevrez une erreur dans une fenêtre que vous pouvez ignorer.

Ve verrez alors les couches sur le GeoNode listées.

Nous allons utiliser les données de l’aquifère de Stampriet.

12. Taper stampriet au filtre

13. Cliquer sur la couche stasboreholes WFS (indiquée dans la figure ci-dessous)

Nous utilisons un service Web WFS, parce que nous voulons accéder aux données et non pas à une figure comme avec des données WMS.

14. Cliquer sur Add pour ajouter la couche à la toile de la carte.

15. Cliquer sur close pour fermer la fenêtre de dialogue

Le chargement depuis internet prendra un moment.

Patienter jusqu’à ce que le chargement soit terminé.

Dans la section suivante nous allons enregistrer cette couche WFS dans on Geopackage pour une utilisation locale ultérieure.

3. Exporter des données WFS dans un GeoPackage

4. Choisir les forages dans une zone donnée

Pour ce tutoriel nous avons uniquement besoin des forages situés dans une zone autour de Stampriet.

Un polygone des limites imaginaires est fourni dans les données du toturiel.

1. Ajouter le shapefile des limites de Stampriet au projet.

2. Le panneau des couches, faire un clic droit sur les limites de Stampriet et choisir zoom sur la couche.

3. styler avec un extérieur rouge.

Nous allons à présent couper les forages qui tombent à l’intérieur des limites et les enregistrer dans notre Geopackage.

4. A partir du menu principal choisir Vecteur/coupper (clip)

5. Dans fenêtre de coupe choisir Stampriet_Data Boreholes comme données d’entrée et Stampriet_boundary comme couche de superposition. Enregistrer le résultat de la coupe dans le Stampriet_Data GeoPackage avec le nom boreholes_clipped (coupé).

6. Cliquer sur exécuter. Cliquer fermer pour fermer la fenêtre de dialogue après le traitement.

Dans la prochaine section nous allons calculer la densité des forages dans cette zone.

5. calculer la densité des forages

Nous avons à présent les forages au sein de notre zone d’étude (ceci aurait aussi pu être un aquifère). Nous pouvons à présent calculer la densité de forages.

1. Ouvrir la boite à outils de traitement : dans le menu principal aller à Processing | Toolbox

3. La fenêtre compter des points choisir Stampriet_boundary comme Polygone et boreholes_clipped comme Points. Changer le nom du champ compter à No of boreholes. Enregistrer la sortie du comptage de Stampriet_Data GeoPackage sous le nom de couche de sortie Stampriet study area.

4. Cliquer exécuter. Cliquer fermer après traitement.

Remarquez que vos entrées ont de multiple polygones (e.x. plusieurs aquifères). Il effectue le count par polygone. Vous pouvez donc appliquer la prochaine étape pour calculer la densité de forage pour tous les polygones.

5. Copier le style de Stampriet_boundary à Stampriet study area : dans le panneau des couches faire clic droit Stampriet_boundary et choisir Styles | Copy Style | All Style Categories. Cliquer ensuite sur Stampriet study area et choisir choose Styles | Paste Style | All Style Categories

6. Retirer Stampriet_boundary du panneau des couches

A présent jetons un coup d’œil sur la table des attributs de Stampriet study area

7. La panneau des couches, faire un clic droit sur Stampriet study area et choisir ouvrir la table des attributs.

Vous verrez que nous avons 64 forages dans la zone d’étude.

Nous devons à présent calculer la superficie du polygone de la zone d’étude

9. Cliquer pour ouvrir le calculateur de champs.

10. Dans la fenêtre de calcul de champ, créer un nouveau champ avec un nom de champ de sorite area en (Km2) et type de sortie nombre décimal (réel).

11. Sous Expression choisir parmi les fonctions dans le panneau central de géométrie | $area et doucle cliquer sur $area pour l’ajouter à l’expression. $area  calcule la superficie du polygone suivant l’unité de la carte. Dans notre cas en m2. Compléter l’expression comme suit :

$area / 1000000

12. Cliquer sur OK

La prochaine étape sert à calculer la densité des forages.

14. Dans la fenêtre du calculateur de champ créer un nouveau champ avec le nom de sortie Borehole density per km2 et type de champ de sortie nombre décimal (réel).

Sous Expression avec l’équation suivante : "No of boreholes” / “Area (km2)"

15. Cliquer sur Ok

Nous avons calculé que la densité des forages dans la zone d’étude est 0.2 par km².

Dans la prochaine section, nous allons ajouter un MNT au projet.

6. Télécharger le MNT SRTM 1-Arc Second

Dans cette section nous allons ajouter un MNT au projet. Nous allons utiliser le produit de Shuttle Radar Topography Mission (SRTM) à la résolution 1 Arc Second. Ceci est un MNT gratuit qui couvre de 56°S à 60°N et a une résolution spatiale d’environ 30 m au niveau de l’équateur.

Vous pouvez télécharger les feuilles de SRTM o partir du site web de USGS Earth Explorer grâce à un plugin. Nous allons utiliser le plugin de téléchargement de SRTM dans QGIS.

1.Dans le menu principal aller à plugins| Manage and Install Plugins...

2. Rechercher SRTM-Downloader plugin, installer le plugin et retourner au projet.

3. Dans la boite à outils cliquer sur le bouton pour ouvrir la fenêtre du téléchargeur de SRTM.

4. Dans la fenêtre du téléchargeur de SRTM régler l’étendu de la toile (assurez vous que vous êtes toujours zoomés aux limites de la zone d’étude). Pour le Output-path naviguer au dossier ou les données de ce projet sont stockées.

5. Click Download.

5. Cliquer sur télécharger

7. Renseigner votre nom d’utilisateur et mot de passe. Cliquer sur Save credentials et cliquer sur OK

Il commencera à présent à télécharger la feuille SRTM qui couvre la zone d’étude.

8. Cliquer OK sur la fenêtre.

9. Fermer la fenêtre.

7. Couper et reprojeter le MNT

La feuille SRTM est ajoutée à la toile de la carte. Si vous zoomez à la couche S25E018.hgt vous verrez qu’elle couvre une bien plus grande zone. La projection est également Système de Projection Géographique avec les coordonnées latitude/longitude en degré (EPSG :4326).

Nous devons donc couper et reprojeter le MNT et l’ajouter à notre GeoPackage.

2. Dans la fenêtre de dialogue « save raster layer as… » choir GeoPackage comme format de sorie. Naviguer au fichier Stampriet_Data.gpkg et choisir DEM comme nom de la couche. Changer le CRS au CRS du projet: EPSG:32734.

Sous Extension choisir Calculer à partir d’une couche et sélectionner Stampriet study area

Sous Resolution changer la résolution horizontale et verticale à 30m

Vérifier les valeurs No data et ajouter -9999 sous de…. à …… Ceci sera la valeur hors gamme pour les cellules qui seront hors de la zone d’étude après la coupe.

3. Laisser les autres réglages par défaut et cliquer sur OK.

4. Retirer la couche the S25E018.hgt du panneau des couches et zoomer aux limites de la zone d’étude pour voir le résultat.

Dans la prochaine section nous allons styler le MNT.

8. Styler le MNT

Pour la visualisation il est commode de styler le MNT avec une rampe de couleur et des effets ombrage

1. Faire un clic droit sur le MNT et choisir dupliquer couche.

3. Sélectionner DEM dans le panneau des couches et cliquer pour ouvrir le panneau Layer Style.

10. A partir du panneau des couches activer la couche Hillshade.

Votre projet doit ressembler maintenant à l’image ci-dessous.

Dans la prochaine section, nous allons ajouter les élévations du MNT à la table des attributs des forages.

9. Echantillonner les élévations du MNT aux forages

Nous allons à présent échantillonner les élévations MNT au niveau des locations des forages. Pour cela nous avons besoin d’installer le plugin pour un échantillonnage de point. Remarquez que l’on a besoin d’une connexion internet pour installer le plugin.

1. Dans le menu principal aller à plugins/Manage and Intall Plugins…

2. Dans la fenêtre de dialogue du Plugin rechercher l’outil Point sampling et cliquer sur installer.

3. cliquer sur le bouton dans la barre à outils pour ouvrir la fenêtre de l’outil Point Sampling.

4. Dans la fenêtre de dialogue de l’outil Point Sampling dans l’onglet General sélectionner boreholes_clipped comme couche contenant les points d’échantillonnage.  Choisir tous les champs de boreholes_clipped et la couche MNT. Remarquez que ces couches doivent être sélectionner dans le panneau des couches pour être visible dans la liste. Enregistrer la sortie vecteur points en tant que shapefile ESRI, parce que cet outil ne peut l’ajouter à un Geopackage. Le renommer Boreholes_Z.shp.

5. Maintenant cliquer sur l’onglet des champs. Là-bas vous pouvez changer le nom des champs de sortie si nécessaire. Ici nous changeons uniquement le champ du MNT en Z.

6.  Cliquer sur Ok. Cliquer sur fermer pour fermer la fenêtre de dialogue après le traitement.

7. Vérifier la table des attributs de la couche Boreholes_Z.

Vous pouvez voir que le champ Z a été ajouté.

Dans la prochaine section, nous allons ajouter un diagramme de dispersion des élévation rapportées et les valeurs de Z du MNT.

10. Comparer les valeurs rapportées aux valeurs du MNT

Dans la table des attributs de Boreholes_Z nous avons les valeurs d’élévations rapportées et de valeurs du MNT. Comparons les deux champs en faisant un digramme de dispersion.

Le champ elevation contient plusieurs données manquantes (no data). Ignorons les dans le diagramme de dispersion

1. Ouvrir la table des attributs de Boreholes_Z.

3. Ecrire l’expression

"elevation" < 0

4. Cliquer sur Select features (sélectionner les éléments) et sur Close pour fermer la fenêtre de dialogue.

10 éléments ont été sélectionnés. Ont pourrait les retirer de la table des attributs en basculant en mode éditer et cliquant sur la poubelle. Pour le digramme de dispersion cependant, nous pouvons inverser la sélection et indiquer que nous voulons utiliser uniquement les éléments sélectionnés.

5. Cliquer pour inverser la selection.

Pour faire le diagramme de dispersion, nous devons installer le plugin « Dataploty ».

6. Installer le le plugin « Data ploty » à partir du manager de Plugins.

7. Cliquer sur l’icon dans la barre des outils pour ouvrir le panneau Data Plotly.

8. Dans le panneau du Data ploty choisir Scatter plot (diagramme de dispersion) comme type de diagramme (plot). Choisir Boreholes_Z.shp comme couche. Vérifier la boxe pour n’utiliser que les élément choisis. Choisir Elevation pour l champ X et Z pour le champ Y. Garder le style par défaut.

9. Cliquer pour aller à l’onglet de configuration de la légende.

10. Décocher voir la boxe de légende. Changer le titre du diagramme en Borehole vs DEM elevation. Retirer le titre de la légende. Ecrire Borehole elevation (m) pour le nom de X et DEM elevation (m) pour le nom de Y. Garder le style par défaut.

11. Cliquer sur Créer diagramme (create plot) au bas du panneau.

Dans la prochaine section nous allons remplacer les valeurs manquantes avec les élévations du MNT.

11. Remplacer les nodata avec les élévations du MNT

Comme le diagramme de dispersion parrait raisonnable, nous allons remplacer les nodata (données manquntes) du champ elevation de la couche Boreholes_Z.shp à partir des valeurs du MNT.

1. Ouvrir la table des attributs à partir de la couche Boreholes_Z.shp.

2. Cliquer sur tout désélectionner.

3. Cliquer pour basculer en mode édition.

4. Changer la sortie en elevation et cliquer le bouton pour écrire l’expression.

5. Dans la fenêtre de dialogue écrire l’expression suivante :

if( "elevation" < 0,  "Z" ,  "elevation" )

Cela signifie: si elevation à des valeurs inférieures à zero, alors utliser le champ Z, sinon utiliser les valeurs du champ  elevation.  

6. Cliquer sur OK

7. Cliquer sur tout mettre à jour.

A présent toutes les valeurs manquantes ont été remplacées par l’élévation du MNT.

8. Cliquer pour basculer en arrêt du mode édition et cliquer save pour enregistrer les résultats.

Dans la prochaine section, nous allons calculer le niveau absolu de la nappe dans les forages.

12. Calculer le niveau de la nappe dans les forages

Dans la table des attributs des forages nous avons la profondeur de la nappe.

Dans cette section nous allons calculer l’élévation absolue de la nappe dans les forages.

1. Ouvrir la table des attributs de Boreholes_Z.shp.

2. Cliquer pour basculer en mode édition.

3. Cliquer pour ouvrir Select features (sélectionner éléments) en utilisant une fenêtre de dialogue.

4. Ecrire l’expression suivante :

"depth_m" < 0

5. Cliquer sur Select features. Cliquer Close pour fermer la fenêtre de dialogue.

6. Cliquer pour se débarrasser des éléments sélectionnés.

7. Cliquer pour enregistrer les éditions jusqu’à présent.

8. Cliquer pour ouvrir la fenêtre de dialogue du calculateur de champ.

9. Dans la fenêtre de dialogue du calculateur de champ créer un champ avec nom champ de sorite name GW_level_m et type de champ de sortie Nombre décimal (réel). Garder les autres réglages par défaut et taper l’expression suivante :

"elevation" - "depth_m"

10. Cliquer sur OK.

11. Vérifier les résultats et cliquer pour basculer en mode arrêt d’édition. Cliquer sur save pour enregistrer les résultats.

Dans la prochaine section nous allons interpoler les niveaux de la nappe en raster.

13. Interpolation des niveaux de la nappe en raster

Nous avons maintenant les niveaux de la nappe dans les forages. Nous pouvons interpoler les points en raster.

Commençons par l’interpolation par IDW

1. Dans la boite à outils de traitement, choisir Interpolation | IDW interpolation.

2. Dans la fenêtre de dialogue de IDW Interpolation choisir la couche Boreholes_Z.shp couche vecteur. Utiliser GW_level_m comme attribut de l’interpolation. Cliquer pour ajouter ces réglages à la liste. Garder le coefficient de distance P à 2 pour une fonction de décroissance exponentielle. Sélectionner l’étendu (extent) à partir de la couche Stampriet study area. Changer la taille des pixels à 30 m. Enregistrer le raster interpolé dans le dossier du projet en tant que GW_level_IDW.tif.

3. Cliquer sur Exécuter (run). Cliquer sur Close à la fin du traitement.

Stylons le résultat avec un rampe de couleur.

5. Dans le panneau pour styler les couches changer le renderer à Singleband pseudocolor et choisir All Color Ramps | RdBu.

Le panneau de stylisation de la couche devrait ressembler à ceci.

Et la couche interpolée devrait ressembler à la figure ci-dessous. Vous pour essayer plusieurs rampes de couleur et classes pour améliorer la visualisation. Vous pouvez voir que l’interpolateur IDW créé une sortie lisse, mais comme le nombre de forages est faible, il y a de nombreux biais en direction des forage individuels et les aberrations.

Appliquons à présent les polygones de Thiessen (au connu comme voisin le plus proche ou tesselation de Voronoi).

6. Dans le menu principal choisir Raster | Analysis | Grid (Nearest Neighbor).

Dans la grille de dialogue (plus proche voisin) choisir la couche Boreholes_Z.shp comme couche de points. Sous Advanced (avancée) choisir GW_level_m comme champ de valeur Z. Enregistrer la sortie interpolée en tant que GW_level_Thiessen.tif. Garder le reste par défaut.

8. Cliquer sur Exécuter (run). Cliquer sur close pour fermer la fenêtre de dialogue.

9. Styler le résultat.

Cet algorithme ne couvre pas toute la zone d’étude mais se limite à qu’on appelle enveloppe convexe de la couche d’entrée des points. Il donne un résultat moins lisse comparé à l’interpolation IDW. Cependant, avec des données rares, les plus proches valeurs peuvent être la meilleure approximation. 

Nous pouvons terminer le projet en stylant les puits et exportant toutes les couches que nous voulons garder dans le GeoPackage.