Importante

La traduzione è uno sforzo comunitario you can join. Questa pagina è attualmente tradotta al 62.50%.

19. Libreria analisi di rete

Suggerimento

I frammenti di codice in questa pagina hanno bisogno delle seguenti importazioni se sei al di fuori della console di pyqgis:

from qgis.core import (
  QgsVectorLayer,
  QgsPointXY,
)

La libreria di analisi di rete può essere utilizzata per:

• creare un grafico matematico a partire da dati geografici (layer vettoriali polilineari)

• implementare i metodi di base della teoria dei grafi (attualmente solo l’algoritmo di Dijkstra)

In breve, un caso d’uso tipico può essere descritto come segue:

1. creare un grafo da geodati (solitamente layer vettoriali polilineari)

2. eseguire l’analisi del grafo

3. utilizzare i risultati dell’analisi (ad esempio, visualizzarli)

19.1. Costruire un grafo

La prima cosa da fare è preparare i dati di input, cioè convertire un layer vettoriale in un grafo. Tutte le azioni successive utilizzeranno questo grafo, non il layer.

As a source we can use any polyline vector layer. Nodes of the polylines become graph vertices, and segments of the polylines are graph edges. If several nodes have the same coordinates then they are the same graph vertex. So two lines that have a common node become connected to each other.

Inoltre, durante la creazione del grafo è possibile «fissare» («legare») al layer del vettore di input un numero qualsiasi di punti aggiuntivi. Per ogni punto aggiuntivo verrà trovata una corrispondenza — con il vertice del grafo più vicino o con il bordo del grafo più vicino. In quest’ultimo caso, il bordo verrà diviso e verrà aggiunto un nuovo vertice.

Gli attributi del layer vettoriale e la lunghezza di un bordo possono essere utilizzati come proprietà di un bordo.

La conversione da un layer vettoriale a grafo viene effettuata utilizzando il modello di programmazione Builder. Un grafo viene costruito usando un cosiddetto Director. Per ora esiste un solo Director: QgsVectorLayerDirector. Il Director imposta le impostazioni di base che verranno utilizzate per costruire un grafico da un layer vettoriale di linee, utilizzato dal costruttore per creare il grafico. Attualmente, come nel caso del director, esiste un solo costruttore: QgsGraphBuilder, che crea oggetti QgsGraph. Puoi implementare i propri costruttori per creare un grafo compatibile con librerie come BGL o NetworkX.

To calculate edge properties the programming pattern strategy is used. For now only QgsNetworkDistanceStrategy strategy (that takes into account the length of the route) and QgsNetworkSpeedStrategy (that also considers the speed) are available. You can implement your own strategy that will use all necessary parameters.

È il momento di addentrarsi nel processo.

1. First of all, to use this library we should import the analysis module:

   from qgis.analysis import *
   

2. Then some examples for creating a director:

   # Don't use information about road direction from layer attributes,
   # all roads are treated as two-way
   director = QgsVectorLayerDirector(
       vectorLayer, -1, "", "", "", QgsVectorLayerDirector.DirectionBoth
   )
   

   # Use field with index 5 as source of information about road direction.
   # one-way roads with direct direction have attribute value "yes",
   # one-way roads with reverse direction have the value "1", and accordingly
   # bidirectional roads have "no". By default roads are treated as two-way.
   # This scheme can be used with OpenStreetMap data
   director = QgsVectorLayerDirector(
       vectorLayer, 5, "yes", "1", "no", QgsVectorLayerDirector.DirectionBoth
   )
   

   To construct a director, we should pass a vector layer that will be used as the source for the graph structure and information about allowed movement on each road segment (one-way or bidirectional movement, direct or reverse direction). The call looks like this (find more details on the parameters at qgis.analysis.QgsVectorLayerDirector):

   director = QgsVectorLayerDirector(
       vectorLayer,
       directionFieldId,
       directDirectionValue,
       reverseDirectionValue,
       bothDirectionValue,
       defaultDirection,
   )
   

3. È quindi necessario creare una strategia per il calcolo delle proprietà dei bordi

   # The index of the field that contains information about the edge speed
   attributeId = 1
   # Default speed value
   defaultValue = 50
   # Conversion from speed to metric units ('1' means no conversion)
   toMetricFactor = 1
   strategy = QgsNetworkSpeedStrategy(attributeId, defaultValue, toMetricFactor)
   

4. E informare il director di questa strategia

   director = QgsVectorLayerDirector(vectorLayer, -1, "", "", "", 3)
   director.addStrategy(strategy)
   

5. Now we can use the builder, which will create the graph, using the QgsGraphBuilder class constructor.

   # only CRS is set, all other values are defaults
   builder = QgsGraphBuilder(vectorLayer.crs())
   

6. Also we can define several points, which will be used in the analysis. For example:

   startPoint = QgsPointXY(1179720.1871, 5419067.3507)
   endPoint = QgsPointXY(1180616.0205, 5419745.7839)
   

7. Now all is in place so we can build the graph and «tie» these points to it:

   tiedPoints = director.makeGraph(builder, [startPoint, endPoint])
   

   Building the graph can take some time (which depends on the number of features in a layer and layer size). tiedPoints is a list with coordinates of «tied» points.

8. When the build operation is finished we can get the graph and use it for the analysis:

   graph = builder.graph()
   

9. With the next code we can get the vertex indexes of our points:

   startId = graph.findVertex(tiedPoints[0])
   endId = graph.findVertex(tiedPoints[1])
   

19.2. Analisi grafo

Networks analysis is used to find answers to two questions: which vertices are connected and how to find a shortest path. To solve these problems the network analysis library provides Dijkstra’s algorithm.

Dijkstra’s algorithm finds the shortest route from one of the vertices of the graph to all the others and the values of the optimization parameters. The results can be represented as a shortest path tree.

L’albero del percorso più breve è un grafo diretto pesato (o più precisamente un albero) con le seguenti proprietà:

• solo un vertice non ha bordi in entrata: la radice dell’albero.

• all other vertices have only one incoming edge

• se il vertice B è raggiungibile dal vertice A, allora il percorso da A a B è il solo percorso disponibile ed è ottimale (il più breve) su questo grafo

Per ottenere l’albero del percorso più breve, utilizzare i metodi shortestTree() e dijkstra() della classe QgsGraphAnalyzer. Si consiglia di utilizzare il metodo dijkstra() perché funziona più velocemente e utilizza la memoria in modo più efficiente.

The shortestTree() method is useful when you want to walk around the shortest path tree. It always creates a new graph object (QgsGraph) and accepts three variables:

• source — grafo in ingresso

• startVertexIdx — indice del punto sull” albero (la radice dell’albero)

• criterionNum — numero di proprietà del bordo da usare (a partire da 0).

tree = QgsGraphAnalyzer.shortestTree(graph, startId, 0)

Il metodo dijkstra() ha gli stessi argomenti, ma restituisce una tupla di array:

• Nel primo array, l’elemento n contiene l’indice del bordo in entrata o -1 se non ci sono bordi in entrata.

• Nel secondo array, l’elemento n contiene la distanza dalla radice dell’albero al vertice n o DOUBLE_MAX se il vertice n non è raggiungibile dalla radice.

(tree, cost) = QgsGraphAnalyzer.dijkstra(graph, startId, 0)

Here is some very simple code to display the shortest path tree using the graph created with the shortestTree() or the dijkstra() method (select linestring layer in Layers panel and replace coordinates with your own).

Avvertimento

Usa questo codice solo come esempio; crea molti oggetti QgsRubberBand e può essere lento su grandi insiemi di dati.

shortestTree()dijkstra()

from qgis.core import *
from qgis.gui import *
from qgis.analysis import *
from qgis.PyQt.QtCore import *
from qgis.PyQt.QtGui import *

vectorLayer = QgsVectorLayer(
    "testdata/network.gpkg|layername=network_lines", "lines"
)
director = QgsVectorLayerDirector(
    vectorLayer, -1, "", "", "", QgsVectorLayerDirector.DirectionBoth
)
strategy = QgsNetworkDistanceStrategy()
director.addStrategy(strategy)
builder = QgsGraphBuilder(vectorLayer.crs())

pStart = QgsPointXY(1179661.925139, 5419188.074362)
tiedPoint = director.makeGraph(builder, [pStart])
pStart = tiedPoint[0]

graph = builder.graph()

idStart = graph.findVertex(pStart)

tree = QgsGraphAnalyzer.shortestTree(graph, idStart, 0)

i = 0
while i < tree.edgeCount():
    rb = QgsRubberBand(iface.mapCanvas())
    rb.setColor(Qt.red)
    rb.addPoint(tree.vertex(tree.edge(i).fromVertex()).point())
    rb.addPoint(tree.vertex(tree.edge(i).toVertex()).point())
    i = i + 1

19.2.1. Trovare i percorsi più brevi

Per trovare il percorso ottimale tra due punti si utilizza il seguente approccio. Entrambi i punti (inizio A e fine B) sono «legati» al grafo quando viene costruito. Quindi, utilizzando il metodo shortestTree() o dijkstra() si costruisce l’albero del percorso più breve con radice nel punto iniziale A. Nello stesso albero si trova anche il punto finale B e si inizia a percorrere l’albero dal punto B al punto A. L’intero algoritmo può essere scritto come:

assign T = B
while T != B
    add point T to path
    get incoming edge for point T
    look for point TT, that is start point of this edge
    assign T = TT
add point A to path

At this point we have the path, in the form of the inverted list of vertices (vertices are listed in reversed order from end point to start point) that will be visited during traveling by this path.

Here is the sample code for QGIS Python Console (you may need to load and select a linestring layer in TOC and replace coordinates in the code with yours) that uses the shortestTree() or dijkstra() method:

shortestTree()dijkstra()

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

from qgis.PyQt.QtCore import *
from qgis.PyQt.QtGui import *

vectorLayer = QgsVectorLayer(
    "testdata/network.gpkg|layername=network_lines", "lines"
)
director = QgsVectorLayerDirector(
    vectorLayer, -1, "", "", "", QgsVectorLayerDirector.DirectionBoth
)
strategy = QgsNetworkDistanceStrategy()
director.addStrategy(strategy)

builder = QgsGraphBuilder(vectorLayer.sourceCrs())

startPoint = QgsPointXY(1179661.925139, 5419188.074362)
endPoint = QgsPointXY(1180942.970617, 5420040.097560)

tiedPoints = director.makeGraph(builder, [startPoint, endPoint])
tStart, tStop = tiedPoints

graph = builder.graph()
idxStart = graph.findVertex(tStart)

tree = QgsGraphAnalyzer.shortestTree(graph, idxStart, 0)

idxStart = tree.findVertex(tStart)
idxEnd = tree.findVertex(tStop)

if idxEnd == -1:
    raise Exception("No route!")

# Add last point
route = [tree.vertex(idxEnd).point()]

# Iterate the graph
while idxEnd != idxStart:
    edgeIds = tree.vertex(idxEnd).incomingEdges()
    if len(edgeIds) == 0:
        break
    edge = tree.edge(edgeIds[0])
    route.insert(0, tree.vertex(edge.fromVertex()).point())
    idxEnd = edge.fromVertex()

# Display
rb = QgsRubberBand(iface.mapCanvas())
rb.setColor(Qt.green)

# This may require coordinate transformation if project's CRS
# is different from layer's CRS
for p in route:
    rb.addPoint(p)

19.2.2. Zone di disponibilità

The area of availability for vertex A is the subset of graph vertices that are accessible from vertex A and the cost of the paths from A to these vertices are not greater that some value.

Questo può essere mostrato più chiaramente con il seguente esempio: «C’è una stazione dei vigili del fuoco. Quali zone della città può raggiungere un camion dei pompieri in 5 minuti? 10 minuti? 15 minuti?». Le risposte a queste domande sono le zone di disponibilità della caserma dei pompieri.

Per trovare le zone di disponibilità possiamo usare il metodo dijkstra() della classe QgsGraphAnalyzer. È sufficiente confrontare gli elementi dell’array cost con un valore predefinito. Se cost[i] è minore o uguale a un valore predefinito, allora il vertice i è all’interno della zona di disponibilità, altrimenti è fuori.

A more difficult problem is to get the borders of the area of availability. The bottom border is the set of vertices that are still accessible, and the top border is the set of vertices that are not accessible. In fact this is simple: it is the availability border based on the edges of the shortest path tree for which the source vertex of the edge is accessible and the target vertex of the edge is not.

Here is an example:

director = QgsVectorLayerDirector(
  vectorLayer, -1, "", "", "", QgsVectorLayerDirector.DirectionBoth
)
strategy = QgsNetworkDistanceStrategy()
director.addStrategy(strategy)
builder = QgsGraphBuilder(vectorLayer.crs())


pStart = QgsPointXY(1179661.925139, 5419188.074362)
delta = iface.mapCanvas().getCoordinateTransform().mapUnitsPerPixel() * 1

rb = QgsRubberBand(iface.mapCanvas())
rb.setColor(Qt.green)
rb.addPoint(QgsPointXY(pStart.x() - delta, pStart.y() - delta))
rb.addPoint(QgsPointXY(pStart.x() + delta, pStart.y() - delta))
rb.addPoint(QgsPointXY(pStart.x() + delta, pStart.y() + delta))
rb.addPoint(QgsPointXY(pStart.x() - delta, pStart.y() + delta))

tiedPoints = director.makeGraph(builder, [pStart])
graph = builder.graph()
tStart = tiedPoints[0]

idStart = graph.findVertex(tStart)

(tree, cost) = QgsGraphAnalyzer.dijkstra(graph, idStart, 0)

upperBound = []
r = 1500.0
i = 0
tree.reverse()

while i < len(cost):
    if cost[i] > r and tree[i] != -1:
        outVertexId = graph.edge(tree [i]).toVertex()
        if cost[outVertexId] < r:
            upperBound.append(i)
    i = i + 1

for i in upperBound:
    centerPoint = graph.vertex(i).point()
    rb = QgsRubberBand(iface.mapCanvas())
    rb.setColor(Qt.red)
    rb.addPoint(QgsPointXY(centerPoint.x() - delta, centerPoint.y() - delta))
    rb.addPoint(QgsPointXY(centerPoint.x() + delta, centerPoint.y() - delta))
    rb.addPoint(QgsPointXY(centerPoint.x() + delta, centerPoint.y() + delta))
    rb.addPoint(QgsPointXY(centerPoint.x() - delta, centerPoint.y() + delta))