Fonctions géométriques

Nous allons lister dans un fichier CSV les adresses des Hautes-Alpes se trouvant à 2kms des limites d'un glacier.

Au terme de cet exercice, il nous faut un code complet permettant de générer ce fichier à partir d'un projet QGIS vide, et d'y inclure le nombre d'adresses concernées.

Buffer

En guise d'entrainement, commencez par créer une zone tampon de 5 mètres autour de nos sommets.

vlayer1 = QgsVectorLayer("C:/Users/Georges/Downloads/qgis/peaks/peaks.shp", "Mes sommets", "ogr")
QgsProject.instance().addMapLayer(vlayer1)
 
bufferPath = "C:/Users/Georges/Downloads/qgis/buffers/buffer1.shp"
 
processing.run('native:buffer', {"INPUT": vlayer1, "DISTANCE": 5, \
"SEGMENTS": 10, "END_CAP_STYLE": 0, "DISSOLVE": False, "OUTPUT": bufferPath})
 
vbuffer1 = QgsVectorLayer(bufferPath, "Buffer 1", "ogr")
QgsProject.instance().addMapLayer(vbuffer1)
 
# Passer Buffer en bas
root = QgsProject.instance().layerTreeRoot()
myAlayer = root.findLayer(vbuffer1.id())
myClone = myAlayer.clone()
parent = myAlayer.parent()
parent.insertChildNode(-1, myClone)
parent.removeChildNode(myAlayer)

Les 2 premières lignes appellent la couche peaks, la 4^ème stocke le chemin de notre future zone tampon.

C'est les lignes 6 et 7 qui à elles-seules créent le buffer grâce à l'instruction native/buffer. D'autres paramètres sont disponibles.

Les lignes 9 et 10 appellent ensuite notre couche contenant la zone tampon dans le canvas. Le reste du code est là pour passer les tampons sous les sommets, afin d'y voir clair.

Mais le code ci-dessus fonctionne très mal n'est-ce pas ?

• Essayez de tricher
• Essayez de comprendre pourquoi ça ne marche pas
• Inspectez les unités géométriques utilisées dans votre canvas
• Contrôlez vos logs dans View / Panels / Log Messages

Utilisez plutôt le code suivant pour créer des buffer de 2 kms autour de nos 108 000 adresses des Hautes-Alpes.

vlayer2 = QgsVectorLayer("C:/Users/Georges/Downloads/qgis/data_BDTOPO_V3_Dep05_adresse/data_BDTOPO_V3_Dep05_adresse/ADRESSE_05.shp", "Addresses", "ogr")
QgsProject.instance().addMapLayer(vlayer2)
 
bufferPath = "C:/Users/Georges/Downloads/qgis/buffers/buffer1.shp"
 
processing.run('native:buffer', {"INPUT": vlayer2, "DISTANCE": 2000, \
"SEGMENTS": 10, "END_CAP_STYLE": 0, "DISSOLVE": False, "OUTPUT": bufferPath})
 
vbuffer1 = QgsVectorLayer(bufferPath, "Buffer 1", "ogr")
QgsProject.instance().addMapLayer(vbuffer1)
 
# Passer Buffer en bas
root = QgsProject.instance().layerTreeRoot()
myAlayer = root.findLayer(vbuffer1.id())
myClone = myAlayer.clone()
parent = myAlayer.parent()
parent.insertChildNode(-1, myClone)
parent.removeChildNode(myAlayer)

Glaciers Alpins

Nous allons compter combien d'adresses se trouvent à 2 kms de distance des limites d'un glacier.

Téléchargez les glaciers contenus dans l'emprise de nos adresses Hautes-Alpines en utilisant l'extension QuickOSM (filtre natural=glaciers). Une couche sur ce modèle est disponible dans les pièces jointes de cet article.

Enregistrez-en un shape dans le même système de projection que vos adresses/tampons.

Pour contrôler la projection d'un layer.

{layer}.crs()

Intersection

L'algorithme de calcul d'intersection se présente de façon assez similaire que celle du buffer. D'autres paramètres sont disponibles.

processing.run('qgis:intersection', {\
"INPUT": {LayerInput},\
"OVERLAY": {LayerOverlay},\
"INPUT_FIELDS": {[list]},\
"OVERLAY_FIELDS": {[list]},\
"OVERLAY_FIELDS_PREFIX": {string},\
"OUTPUT": {File}})

Introduisez l'algorithme ci-dessus dans votre code, exécutez-le puis zoomez, ci-dessous, sur le hameau de la Bérarde, au cœur de la zone de glaciers des Écrins, pour contrôler vos résultats.

# Zoom sur la Berarde
myselect = layerAdresse.getFeatures( QgsFeatureRequest().setFilterExpression ( u'"NOM_1" = \'LA BERARDE\'' ) )
layerAdresse.selectByIds( [ f.id() for f in myselect ] )
iface.mapCanvas().zoomToSelected(layerAdresse)
iface.mapCanvas().zoomScale(100000)

Vos intersections doivent exister pour chaque adresse intersectées (présence de doublons opportuns).

Si votre vérification est bonne (même nombre d'entités intersectées dans un groupe de zones tampons proches que d'adresses au sein de ce groupe), alors il ne vous reste plus qu'à compter le nombre total d'intersections.

{Layer}.featureCount()

Export CSV

Vous savez déjà exporter un fichier CSV (voir chapitre Exporter des fichiers plus haut dans ce court), à vous de jouer !

N'oubliez pas d'inclure le nombre d'adresses concernées dans le nom du fichier.

Nettoyage

La couche de nos zone tampons et celle de nos intersections sont à présent inutiles. Nous avons fait le job, nous pouvons donc supprimez ces couches mais aussi fermer QGIS.

from osgeo import ogr
...
# Nettoyage
QgsProject.instance().removeMapLayer(layerBuffer)
QgsProject.instance().removeMapLayer(intersection)
 
driver = ogr.GetDriverByName("ESRI Shapefile")
 
driver.DeleteDataSource(bufferPath)
driver.DeleteDataSource(interstectPath)
 
os._exit(0)

Il semble toute fois que la suppression totale des shapes soit impossible à réaliser à la suite du code complet. Sans doute une mémoire cache à vider dans QGIS ou Python.

Code complet de la génération des buffers, des intersections et du fichier CSV (modifiez les chemins, fonctionne dans l'éditeur Python exclusivement)

• Note personnelle : aller plus loin en se passant de la génération des shapes des buffers (couche virtuelle à la place), des glaciers (API Overpass d'OSM) et des intersections (couche virtuelles). Afin d'obtenir le simple fichier CSV des adresses à 2 kms d'un glacier, et leur nombre.

Liens ou pièces jointes

chemin_de_fer.zip	[ ]	0 Ko
data_BDTOPO_V3_Dep05_adresse.zip	[ ]	3889 Ko
data_IRIS_2019.zip	[ ]	45905 Ko
decathlon_france.zip	[308 magasins Décathlon français depuis OSM le 27 décembre 2020]	11 Ko
eau.zip	[ ]	0 Ko
glaciers.zip	[ ]	231 Ko
iso_iris.zip	[Des zones isochrones à 15 minutes autour de 308 POIs.]	12125 Ko
Koln gml.zip	[ ]	2818 Ko
peaks.zip	[ ]	14 Ko
peaks_selection.zip	[ ]	1 Ko
simple_countries.zip	[ ]	1880 Ko
sol.zip	[ ]	0 Ko
troncons_routes.zip	[ ]	0 Ko