L’objectif de cet exercice est double : manipuler les géométries relatives à chaque entité et éditer une couche pour l’enrichir.
À partir de la couche vectorielle contenant les pays du monde, nous cherchons à trouver tous les polygones limitrophes à chaque polygone (tous les voisins de chaque pays).
La couche d’entrée va être enrichie de deux nouveaux champs, l’un contenant le noms des voisins et l’autre contenant la population totale de ces voisins.
Les données proviennent de Natural Earth [1] :
ne_10m_admin_0_countries.shp
(Polygones - Pays du monde) ➜ TéléchargementLa méthode à utiliser peut être décrite de la manière suivante :
Pour chaque pays P:
trouver les polygones qui intersectent sa bounding box
parmis ces polygones :
trouver ceux qui sont réelement voisins du pays P
enregistrer le nom des voisins de P dans le champ
NEIGHBORS
enregistrer la somme de la population de ces voisins dans le champ
SUM
Cet exerice mobilise plusieurs notions nouvelles :
édition d’une couche
# Au début de l'édition : layer.startEditing() # À la fin de l'édition pour valider : layer.commitChanges() # Pour abandonner les changements en cours : layer.rollBack()
création de nouveaux champs attributaires
from PyQt5.QtCore import QVariant # La couche `layer` doit être en cours d'édition ! # Création d'un champ de chacun des 3 principaux types # qui seront utilisés : layer.dataProvider().addAttributes( [QgsField("Champ_1", QVariant.String), QgsField("Champ_2", QVariant.Int), QgsField("Champ_3", QVariant.Double)]) layer.updateFields()
mise à jour d’une entité après l’avoir modifiée
# La couche `layer` doit être en cours d'édition ! # `ft` est une variable de type QgsFeature layer.updateFeature(ft)
utilisation de prédicats spatiaux (intersects, disjoint, etc.)
# Obtenir la géométrie d'une entité : geom1 = feature1.geometry() geom2 = feature2.geometry() # Prédicats spatiaux : geom1.intersects(geom2) # -> True/False geom1.disjoint(geom2) # -> True/False geom1.touches(geom2) # -> True/False
utilisation d’un index spatial
# Création de l'index : index = QgsSpatialIndex() for f in layer.getFeatures(): index.insertFeature(f) # Interrogation de l'index # (l'argument de la méthode intersects de l'index # spatial doit être de type `QgsRectangle`) intersecting_ids = index.intersects(geom.boundingBox()) # Une des façon d'obtenir les entités qui correspondent à ces ids : request = QgsFeatureRequest().setFilterFids(intersecting_ids) for feature in vector_layer.getFeatures(request): # ...
À vous de jouer ! Essayer d’écrire le script permettant d’atteindre l’objectif de cet exercice
La méthode à utiliser peut être décrite de la manière suivante :
Création de deux champs
NEIGHBORS
etSUM
(attention au type !)Création de l’index spatial
Parcours des entités de la couche, pour chaque pays P:
trouver les polygones qui intersectent sa bounding box (grace à l’index spatial)
parmis ces polygones :
trouver ceux qui sont réelement voisins du pays P
enregistrer le nom des voisins de P dans le champ
NEIGHBORS
enregistrer la somme de la population de ces voisins dans le champ
SUM
mettre à jour l’entité en question
Validation des changements effectués sur la couche.
➜ Quelle est la population des pays voisins de la Hongrie ?
➜ Combien de pays n’ont aucun voisin ?
Footnotes