PyQgis et l'analyse réseau

QGIS 3.x dispose d’un module d’analyse réseau qui peut être appelé au moyen de scripts python. Voici quelques balises pour bien comprendre les outils mis à disposition et démarrer vos projets.

Un premier réseau dans PyQGIS

La première étape de l’analyse réseau est de créer un graphe à partir de données sources. Dans QGIS, ce graphe peut être créé à partir de données vecteurs de type LineString en 2D ou en 3D (et, donc, de LineStringZ, …).

Le nom des différentes interfaces a été modifié depuis QGIS 2.8. Au moment de la publication, la documentation n’était pas remise à jour, mais les concepts de base sont restés identiques.

Le graphe créé sera constité de faces (edge - objets QgsGraphEdge), et de sommets (vertex - objets QgsGraphVertex. Chaque vertex et edge comportera un identifiant unique.

Il est également possible de forcer l’ajout de vertex sur le réseau. Dans ce cas, PyQgis projettera chacun des points sur l’edge le plus proche. L’utilisateur pourra récupérer les points ainsi projetés (tied point) dans l’espace Python et les ré-utiliser par la suite.

Par défaut, on pourrait décrire le traitement effectué par PyQGIS pour constituer un réseau:

• Un vertex est créé au niveau de chaque “noeud” (ou sommet) des lignes du réseau (donc à chaque fois qu’une ligne change de direction) ;
• Les sommets sont reliés par un edge (une face). Si deux lignes partagent le même sommet, alors elles s’entrecroisent sur le graphe. Pour un réseau routier, il s’agira par exemple des carrefours qui relient deux rues.

Quelques illustrations:

Dans le graphe, chaque edge se voit attribuer un coût de parcours.

Le graphe est orienté (les edges ne peuvent être parcourus que dans un seul sens). Il est possible de configurer les coûts et de restreindre les directions possibles depuis la couche ligne de départ. Dans ce cas, PyQGIS ne va créer qu’un seul edge entre deux vertex, et non deux edges orientés dans chaque direction. Nous laissons cependant cette possibilité de côté pour le moment.

Interconnexion des lignes

Les lignes qui s’entrecroisent ne s’interconnecteront pas forcément l’une à l’autre. Pour qu’elle partagent un vertex commun (et, donc, qu’elle soient interconnectées dans le graphe), elles doivent partager un noeud commun.

Ainsi, donc, si vous disposez de lignes de réseau routier, il est nécessaire de s’assurer qu’un noeud commun relie chaque carrefour. Sans quoi, les rues ne seront pas connectées les unes aux autres.

Astuce: la barre d’outil “accrochage” de QGIS peut être utilisée pour faciliter l’accrochage des nouveaux noeuds avec les existants.

Construction du graphe en python

Pour construire le graphe en PyQGIS, plusieurs objets sont utilisés:

• un objet QgsGraphBuilder, qui construira le graphe ;
• un objet QgsVectorLayerDirector, qui “dirigera” la construction du graphe à partir de la couche vecteur ;
• enfin, un objet permettra de déterminer le coût du parcours de chaque face. Ici, l’on a choisi la stratégie de distance (la distance de la face détermine son coût): QgsDistanceNetworkStrategy. Il est également possible d’utiliser la vitesse de parcours comme coût.

from qgis.analysis import *
from qgis.core import *


# on appelle la couche active
# **ATTENTION** si la couche sélectionnée n'est pas le réseau, le script échouera silencieusement
layer = iface.activeLayer()
crs = layer.sourceCrs()

# l'objet graphBuilder est instancié
builder = QgsGraphBuilder(crs, False)

# l'objet "director" permet de diriger la construction à partir, ici, d'une couche vecteur
# La signature complète est disponible ici: https://qgis.org/pyqgis/3.8/analysis/QgsVectorLayerDirector.html#qgis.analysis.QgsVectorLayerDirector
director = QgsVectorLayerDirector(layer, -1, '', '', '', QgsVectorLayerDirector.DirectionBoth)
# une distance à chaque face sera ajoutée. Elle est créée à partir de la longueur du edge, ici:
# sans properter, toutes les distances sont égales à zéro
director.addStrategy(QgsNetworkDistanceStrategy())

# on crée le graphe, en lui ajoutant les POI
start = QgsPointXY(-1.265, -3.631) # point "home"
end   = QgsPointXY(5.058, -3.177)  # point "school"
poi_list = director.makeGraph(builder, [start, end])
graph = builder.graph()

L’objet poi_list contient une liste des points projetés (start et end) sur le graphe. Les index des points sont respectés: donc le point start ayant l’index 0, son point projeté sur le graphe sera poi_list[0].

>>> print(poi_list)
[<QgsPointXY: POINT(-1.26481915406766543 -3.63082468926068369)>, <QgsPointXY: POINT(4.90604443387641265 -3.19642270696956432)>]

Calcul du chemin le plus court

Une fois le graphe créé, il est possible de calculer des itinéraires le long de celui-ci.

Pour cela, nous allons utiliser l’algorithme Dijkstra. De manière synthétique, il fonctionne de la manière suivante:

• Le graphe est parcouru depuis le point de départ ;
• Tous les noeuds qui sont connectés à ce premier point sont parcourus. Le coût de la distance est associée à chacun de ces noeuds ;
• Une fois tous les noeuds parcourus, l’algorithme “se déplace” vers le point le plus proche, et recommence à parcourir tous les points connectés. Pour chacun des points, il additionne son propre coût et celui à parcourir pour arriver à ce point, et l’associe au point qu’il atteint. Cependant, si le point comporte déjà un coût moins élevé, il ne va pas modifier le coût associé à ce point, et le laisser tel quel ;
• Puis l’algorithme se déplace sur le point le plus proche, et exécute le même calcul. Chaque point n’est visité qu’une seule fois.

L’algorithme s’arrête quand tous les points sont visités, ou quand le point de destination est atteint.

Il est illustré ici: .

Application de l’algorithme de Dijkstra en PyQGIS

L’algorithme peut être appelé avec l’objet QgsGraphAnalyzer.

L’algorithme est appliqué à partir de l’id du point d’arrivée, et un coût est calculé pour chaque vertex du graph:

# On récupère l'id du vertex du point de départ
start_id = graph.findVertex(points_list[0])
# On récupère l'id du vertex du point d'arrivée
end_id  = graph.findVertex(points_list[1])

# l'algo dijkstra est exécuté
(tree, cost) = QgsGraphAnalyzer.dijkstra(graph, start_id, 0)

A ce stade, l’algorithme renvoie deux objets:

• tree: contient une liste constitué d’entiers. L’indice de chaque liste renvoie à l’indice du vertex dans le graphe. La valeur est l’identifiant de l’edge à emprunter pour retourner au noeud d’origine en empruntant le plus court chemin.

  Donc, si l’on considère le vertex avec l’indice 40, tree[40] contient l’identifiant du edge à emprunter pour retourner au point de départ.

  Si le vertex ne peut pas être atteint via le réseau, la valeur est égale à -1. Cela arrive lorsque, par exemple, le noeud n’est pas connecté au graphe.

• cost contient également une liste, cette fois composée de réels. L’indice correspond toujours à un identifiant de edge, et la valeur correspond au coût pour rejoindre le noeud de départ depuis ce edge.

Reconstituer le parcours

Il est donc possible de reconstituer le chemin à parcourir depuis n’importer quel noeud pour retourner vers le noeud de départ. Pour chaque noeud, on peut retrouver l’identifiant de l’edge à parcourir à partir de tree. Le edge correspondant peut être récupéré dans le graphe, il est une instance de QgsGraphEdge. Il renseigne un noeud de départ et d’arrivée. On peut donc ré-interroger l’objet tree pour obtenir l’id du edge suivant, et ainsi de suite.

Exemple:

>>> # on part du noeud avec l'id 0:
>>> tree[0]
1
>>> # le edge suivant est celui qui a l'id 1
>>> # on récupère le edge avec l'id 1
>>> e = graph.edge(1)
>>> # vertex d'arrivée de ce edge:
>>> e.toVertex()
0
>>> # vertex de départ de ce edge:
>>> e.fromVertex()
1
>>> # on peut donc parcourir l'objet tree à partir du noeud suivant,
>>> # ayant l'id 1, et ainsi de suite
>>> tree[1]
2

Graphiquement, cela donne ceci:

La boucle suivante permet de reconstituer le chemin, et de le dessiner dans le canevas. Une ligne rouge va matérialiser le chemin, et les edge seront marqués par une croix:

# Vérification: peut-on bien obtenir le point d'arrivée ?
if tree[end_id] == -1:
    # note: si la couche sélectionnée n'est pas celle du réseau, vous aurez une
    # erreur ici.
    raise Exception("could not reach end point (is network selected ?)")

# list des points que l'on va récupéré
p = list()

# on se positionne sur le point d'arrivée
cur_pos = end_id

# boucle pour remonter au point de départ
while cur_pos != start_id: 
    # on récupère l'edge précédent
    previous_edge = graph.edge(tree[cur_pos])
    # on récupère le noeud de départ
    starting_vertex_id = previous_edge.toVertex() 
    # on ajoute ce noeud à notre liste
    p.append(graph.vertex(starting_vertex_id).point())
    # on se "déplace" vers le point suivant
    cur_pos = previous_edge.fromVertex()

# on ajoute le point de départ pour fermer le chemin
p.append(points_list[0])

# création d'un Rubber band sur la carte
rb = QgsRubberBand(iface.mapCanvas())
rb.setColor(Qt.red)

import random

for pnt in p:
    # ajout d'un marker à chaque vertex rencontré
    m = QgsVertexMarker(iface.mapCanvas())
    m.setCenter(pnt)
    m.setColor(QColor(random.randint(0,250), random.randint(0, 250), random.randint(0, 250)))
    m.setIconType(QgsVertexMarker.ICON_CROSS)
    m.setPenWidth(1)
    # on étend la bande rouge d'un point
    rb.addPoint(pnt)

Obtenir les coûts de parcours

Le coût du parcours peut également être obtenu à partir de l’objet cost:

>>> # coût du parcours depuis le vertex 0:
>>> cost[0]
808405.6919667209
>>> # coût du parcours de l'edge à emprunter depuis le vertex 0:
>>> graph.edge(tree[0]).cost(0)
52325.11077089761
>>> # coût au vertex suivant (le vertex avec l'id 1)
>>> cost[1]
756080.5811958233
>>> # ce cout est bien la différence entre le coût au vertex 0 et au vertex 1
>>> cost[0] - graph.edge(tree[0]).cost(0) == cost[1]
True

Il est possible d’ajouter plusieurs stratégies de calcul de coût, comme signalé plus haut. Dès lors, pour récupérer le coût de parcours d’un edge, il est nécessaire d’indiquer l’indice de la stratégie indiquée dans l’objet QgsVectorLayerDirector. La fonction stratégie renvoie la liste complète des coûts:

>>> graph.edge(1).strategies()
[52325.11077089761]

Optimisation de l’analyse de réseau

La fonction dijkstra est donc appliquée sur tout le graphe, et renvoie un coût global pour chaque vertex du réseau. Elle n’est pas arrêtée lorsque le point d’arrivée est atteint.

Dès lors, si vous prévoyez d’analyser les calculs vers différents vertex à partir du même point de départ, il suffira d’exécuter une seule fois la fonction dijkstra. Le résultat pourra être ré-utilisé.