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

Site: OpenCourseWare for GIS
Cours: Formation SIG pour les applications hydrogéologiques
Livre: Tutoriel: Créer une carte du niveau d’eaux souterraines à partir de données de forage et d’un MNT
Imprimé par: Guest user
Date: jeudi 26 décembre 2024, 17:40

Description

1. Introduction

A la fin de ce tutorial, vous serez capable de :

  • Charger un ensemble de données de forages/puits depuis un SDI vers un QGIS
  • Calculer la densité de forages/puits dans un polygone
  • Télécharger des données MNT de SRTM 1-Arc Second
  • Couper et reprojeter un MNT
  • Styler un MNT
  • Echantillonner l’élévation à l’emplacement d’un forage et le comparer à l’attribut de l’élévation
  • Faire de correction et des calculs dans la table des attributs
  • Interpoler les niveaux des forages en raster

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

1. Vérifier les métadonnées de la base de données de forage de l’aquifère transfrontalier de Stampriet. Lire l’onglet d’information et vérifier les attributs qui sont dans la couche.


stasboreholes in GeoNode

2. Démarrer QGIS avec un nouveau projet

3. Cliquer sur le bouton d’ouverture du manager de sources de données

.

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é.

Loading the stasboreholes

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

Notre ensemble de données de forage est toujours au format WFS et dans un Système de Coordonnées Géographiques (SCG) avec les coordonnées en latitude et longitude. Pour nous besoins d’interpolation des élévations nous devons stocker le fichier localement et le reprojeter en UTM Zone 34S/WGS-84.

1. Dans le champ des couches faire un clic droit sur la couche geonode:stasboreholes

2. Choisir Export/ Enregistrer éléments en tant que….

3. Dans la fenêtre de dialogue enregistrer la couche vecteur en tant que …choisir GeoPackage en tant format. Cliquer sur le bouton pour naviguer à votre disque dur où vous avez stoker le the GeoPackage (e.x. Z:\Stampriet). Le renommer Stampriet_Data.gpkg. Pour le nom de la couche taper Borehole.

4. Cliquer sur le bouton pour changer la projection de couche de sortie.

5. La fenêtre de dialogue de séletion du Système de Coordonnés taper 32734 (qui est le code EPSG) dans le filtre pour rechercher la projection de WGS 84 / UTM Zone 34S.

6. Cliquer Ok pour réaliser l’exportation.

7. Après que la couche de données de forage de Stampriet ait été ajoutée à la toile de la carte, retirer la couche geonode:stasboreholes en faisant un clic droit sur ça et choisissant Retirer la couche… Confirmer avec OK.

Vous avez une copie locale de la couche de forage.

8. Changer la projection OTF du projet : dans le panneau des couches faire un clic droit sur la donnée de forage de Stampriet et choisir régler CRS/ Régler le CRS du Projet à partir d’une couche.


9. Enregistrer le projet.

Dans la prochaine section nous allons sélectionner les forages dans une zone spécifique.


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.

Stampriet boundary styled

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.

7. Retirer la couche Stampriet_Data Boreholes du panneau des couches

Selected boreholes Stampriet

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

2. Dans la boite à outils de traitement choisir analyse de vecteur/compter points dans un polygone


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.

Nous utilisons ce no de sortie car cet outil fait des copies du polygone Stampriet_boundary et y ajoute le champ No du boreholes. Nous n’avions pas encore les limites de la zone d’étude dans notre GeoPackage, nous alors l’y ajouter de cette façon.


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

 

8. Dans la table des attributs cliquer le bouton pour basculer activer le mode édition.

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

Utiliser le bouton pour ajouter la division.

12. Cliquer sur OK

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

13. Cliquer pour ouvrir le calculateur de champ.

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)"

Rappelez vous que vous pouvez ajouter des champs en double cliquant les noms de champs sous Champs et valeur dans le panneau du milieu. Utiliser le bouton pour ajouter la division.


15. Cliquer sur Ok

16. Basculer en mode arrêter éditions en cliquant et cliquer Enregistrer dans la fenêtre pour enregistrer les résultats.

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

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

6. Pour la première fois une fenêtre apparaitra pour les informations de connexion. Pour télécharger des feuille SRTM ou avez besoin d’un compte à  https://urs.earthdata.nasa.gov/users/new.

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.

Comme vous pouvez le voir, la feuille de STRM est plus grande que notre zone d’étude. Dans la prochaine section nous allons couper et reprojeter le MNT.

SRTM tile added to project


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.

1. Dans le panneau des couches faire clic droit sur S25E018.hgt et choisir Export | Save As... (enregistrer sous)

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.

2. Renommer la copie du MNT Hillshade.

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

4. Changer le renderer à Singleband pseudocolor.

5. Au niveau de Couleur Ramp cliquer sur la flèche et choisir créer une nouvelle Color Ramp.



6. Dans la fenêtre de type de color ramp choisir Catalog: cpt-city et clicquer sur OK


7. Dans la fenêtre Cpt-city Color Ramp sélectionner une des palettes sous Topography, e.x. elevation, et cliquer sur OK.


8. Dans le panneau Layer styling cliquer sur classify et applique la color ramp.

9. Sous Layer Rendering change le mode mélange à Multiply. Ceci est nécessaire pour la prochaine étape ou nous allons styler Hillshade et voulons la voir à travers le MNT.



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

11. Dans le panneau Layer Styling assurez vous que vous êtes en train de styler la couche Hillshade à présent. Changer le rééchantillonnage (resampling) à Bilinear pour Zoomed in et Average pour Zoomed out. Cela va créer une visualisation lisse. Garder les autres réglages par défaut.

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

DEM styled

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.

Install Point sampling tool

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.

2. Cliquer le boutons Select by attributes (selection basée sur les attributs) .

3. Ecrire l’expression

"elevation" < 0

Select 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.

Scatter plot elevations

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

L’expression devrait ressembler à ceci dans la table des attributs:


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.

Nous allons utiliser deux méthodes : pondération de la distance inverse (IDW) et les Polygones de Thiessen.

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.

4. Sélectionner la couche GW_level_IDW layer dans le panneau des couches et cliquer pour ouvrir le panneau pour styler les couches (Styling panel)


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.

IDW interpolation result

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. 

Thiessen interpolation result

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