Importante

La traducción es un esfuerzo comunitario puede unirse. Esta página está actualmente traducida en |progreso de traducción|.

19. Biblioteca de análisis de redes

Consejo

Los fragmentos de código en esta página necesitan las siguientes adiciones si está fuera de la consola de pyqgis:

from qgis.core import (
  QgsVectorLayer,
  QgsPointXY,
)

La biblioteca de análisis de red se puede utilizar para:

  • crear un gráfico matemático a partir de datos geográficos (capas vectoriales de polilínea)

  • implementar métodos básicos de la teoría de grafos (actualmente solo el algoritmo de Dijkstra)

Brevemente, un caso de uso típico se puede describir como:

  1. Crear gráfica de geodatos (normalmente de capa vectorial de polilíneas)

  2. ejecutar análisis gráfico

  3. utilizar resultados de análisis (por ejemplo, visualizarlos)

19.1. Contruir un gráfico

Lo primero que hay que hacer — es preparar la entrada de datos, que es convertir una capa vectorial en un gráfico. Todas las acciones adicionales utilizarán esta gráfica, no la capa.

Como fuente podemos utilizar una capa vectorial de polilínea. Los nodos de las polilíneas se convierten en vértices del gráfico, y los segmentos de la polilínea son bordes de gráfico. Si varios nodos tienen la misma coordenada entonces ellos tienen el mimso vértice gráfico. Por lo que dos líneas que tienen un nodo en común se conectaran entre si.

Además durante la creación del gráfico se puede «arreglar» («atar») a la capa vectorial de entrada cualquier número de puntos adicionales. Para cada punto adicional se encontrará una coincidencia — el vértice gráfico más cercano o el borde gráfico más cercano. En el último caso el borde será dividido y un nuevo vértice se añadirá.

Los atributos de la capa vectorial y la longitud de un borde se puede utilizar como las propiedades de un borde.

La conversión de una capa vectorial al gráfico se realiza mediante el patrón de programación Constructor. Un gráfico se construye utilizando el llamado Director. Solo hay un Director por ahora: QgsVectorLayerDirector. El director establece la configuración básica que se utilizará para construir un gráfico a partir de una capa vectorial de línea, utilizada por el constructor para crear el gráfico. Actualmente, como en el caso del director, solo existe un constructor: QgsGraphBuilder, que crea objetos QgsGraph. Es posible que desee implementar sus propios constructores que crearán un gráfico compatible con bibliotecas como 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.

Es tiempo de sumergirse en el proceso.

  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:

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

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

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

    1director = QgsVectorLayerDirector(
2    vectorLayer,
3    directionFieldId,
4    directDirectionValue,
5    reverseDirectionValue,
6    bothDirectionValue,
7    defaultDirection,
8)

  3. Es necesario entonces crear una estrategia para calcular propiedades de borde

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

  4. Y decirle al director sobre esta estrategia

    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. Análisis gráfico

El análisis de redes es utilizado para encontrar respuestas a dos preguntas: que vértices estan conectados y cómo encontrar la ruta más corta. Para resolver estos problemas la librería de análisis de redes proporciona el algoritmo Dijkstra.

El algoritmo Dijkstra encuentra la ruta más corta de uno de los vértices del grafo a todos los otros y los valores de los parámetros de optimización, El resultado puede ser representado como un árbol de la ruta más corta.

El árbol de ruta más corto es un gráfico ponderado dirigido (o más precisamente un árbol) con las siguientes propiedades:

  • sólo un vértice no tiene bordes entrantes — la raíz del árbol

  • todos los otros vértices sólo tienen un borde entrante

  • Si el vértice B es accesible desde el vértice A, entonces el camino de A a B es la única ruta disponible y es optima (más corta) en este grafo

Para obtener la ruta del árbol más corta, use los métodos shortestTree() y dijkstra() de la clase QgsGraphAnalyzer . Es recomendable usar el método dijkstra() porque funciona más rápido y usa la memoria de manera más eficiente.

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 — gráfica entrante

  • startVertexIdx — índice del punto en el árbol (la raíz del árbol)

  • criterionNum — número de propiedad de borde a usar (comenzando desde 0).

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

El método dijkstra() tiene los mismos argumentos, pero devuelve una tupla de matrices:

  • En la primera matriz, el elemento n contiene el índice del borde entrante o -1 si no hay bordes entrantes.

  • En la segunda matriz, el elemento n contiene la distancia desde la raíz del árbol hasta el vértice n o DOUBLE_MAX si el vértice n es inaccesible desde la raíz.

(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).

Advertencia

Use este código solo como ejemplo, crea muchos objetos QgsRubberBand y puede ser lento en conjuntos de datos extensos.

 1from qgis.core import *
 2from qgis.gui import *
 3from qgis.analysis import *
 4from qgis.PyQt.QtCore import *
 5from qgis.PyQt.QtGui import *
 6
 7vectorLayer = QgsVectorLayer(
 8    "testdata/network.gpkg|layername=network_lines", "lines"
 9)
10director = QgsVectorLayerDirector(
11    vectorLayer, -1, "", "", "", QgsVectorLayerDirector.DirectionBoth
12)
13strategy = QgsNetworkDistanceStrategy()
14director.addStrategy(strategy)
15builder = QgsGraphBuilder(vectorLayer.crs())
16
17pStart = QgsPointXY(1179661.925139, 5419188.074362)
18tiedPoint = director.makeGraph(builder, [pStart])
19pStart = tiedPoint[0]
20
21graph = builder.graph()
22
23idStart = graph.findVertex(pStart)
24
25tree = QgsGraphAnalyzer.shortestTree(graph, idStart, 0)
26
27i = 0
28while i < tree.edgeCount():
29    rb = QgsRubberBand(iface.mapCanvas())
30    rb.setColor(Qt.red)
31    rb.addPoint(tree.vertex(tree.edge(i).fromVertex()).point())
32    rb.addPoint(tree.vertex(tree.edge(i).toVertex()).point())
33    i = i + 1

19.2.1. Encontrar la ruta más corta

Para encontrar la ruta óptima entre dos puntos, se utiliza el siguiente enfoque. Ambos puntos (inicio A y final B) están «vinculados» al gráfico cuando se construye. Luego, usando el método shortestTree () o dijkstra() construimos el árbol de ruta más corto con la raiz en el punto inicial A. En el mismo árbol también encontramos el punto final B y comenzamos a caminar a través del árbol desde el punto B al punto A. Todo el algoritmo se puede escribir como:

1assign T = B
2while T != B
3    add point T to path
4    get incoming edge for point T
5    look for point TT, that is start point of this edge
6    assign T = TT
7add point A to path

En este punto tenemos la ruta, en el formulario de la lista invertida de vértices (los vértices están listados en orden invertida del punto final al punto inicial) que serán visitados durante el viaje por este camino.

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:

 1from qgis.core import *
 2from qgis.gui import *
 3from qgis.analysis import *
 4
 5from qgis.PyQt.QtCore import *
 6from qgis.PyQt.QtGui import *
 7
 8vectorLayer = QgsVectorLayer(
 9    "testdata/network.gpkg|layername=network_lines", "lines"
10)
11director = QgsVectorLayerDirector(
12    vectorLayer, -1, "", "", "", QgsVectorLayerDirector.DirectionBoth
13)
14strategy = QgsNetworkDistanceStrategy()
15director.addStrategy(strategy)
16
17builder = QgsGraphBuilder(vectorLayer.sourceCrs())
18
19startPoint = QgsPointXY(1179661.925139, 5419188.074362)
20endPoint = QgsPointXY(1180942.970617, 5420040.097560)
21
22tiedPoints = director.makeGraph(builder, [startPoint, endPoint])
23tStart, tStop = tiedPoints
24
25graph = builder.graph()
26idxStart = graph.findVertex(tStart)
27
28tree = QgsGraphAnalyzer.shortestTree(graph, idxStart, 0)
29
30idxStart = tree.findVertex(tStart)
31idxEnd = tree.findVertex(tStop)
32
33if idxEnd == -1:
34    raise Exception("No route!")
35
36# Add last point
37route = [tree.vertex(idxEnd).point()]
38
39# Iterate the graph
40while idxEnd != idxStart:
41    edgeIds = tree.vertex(idxEnd).incomingEdges()
42    if len(edgeIds) == 0:
43        break
44    edge = tree.edge(edgeIds[0])
45    route.insert(0, tree.vertex(edge.fromVertex()).point())
46    idxEnd = edge.fromVertex()
47
48# Display
49rb = QgsRubberBand(iface.mapCanvas())
50rb.setColor(Qt.green)
51
52# This may require coordinate transformation if project's CRS
53# is different from layer's CRS
54for p in route:
55    rb.addPoint(p)

19.2.2. Áreas de disponibilidad

El área de la disponibilidad para el vértice A es el subconjunto de vértices del grafo que son accesibles desde el vértice A y el costo de los caminos de la A a estos vértices son no es mayor que cierto valor.

Más claramente esto se puede demostrar con el siguiente ejemplo: «Hay una estación de bomberos ¿Qué partes de la ciudad puede un camión de bomberos alcanzar en 5 minutos? 10 minutos? 15 minutos?». Las respuestas a estas preguntas son las zonas de la estación de bomberos de la disponibilidad.

Para encontrar las áreas de disponibilidad podemos usar el método dijkstra() de la clase QgsGraphAnalyzer. Es suficiente comparar los elementos de la matriz de costos con un valor predefinido. Si el costo [i] es menor o igual a un valor predefinido, entonces el vértice i está dentro del área de disponibilidad, de lo contrario, está fuera.

Un problema más difícil es conseguir los límites de la zona de disponibilidad. El borde inferior es el conjunto de vértices que son todavía accesibles, y el borde superior es el conjunto de vértices que no son accesibles. De hecho esto es simple: es la frontera disponibilidad basado en los bordes del árbol de ruta más corta para los que el vértice origen del contorno es más accesible y el vértice destino del borde no lo es.

Here is an example:

 1director = QgsVectorLayerDirector(
 2  vectorLayer, -1, "", "", "", QgsVectorLayerDirector.DirectionBoth
 3)
 4strategy = QgsNetworkDistanceStrategy()
 5director.addStrategy(strategy)
 6builder = QgsGraphBuilder(vectorLayer.crs())
 7
 8
 9pStart = QgsPointXY(1179661.925139, 5419188.074362)
10delta = iface.mapCanvas().getCoordinateTransform().mapUnitsPerPixel() * 1
11
12rb = QgsRubberBand(iface.mapCanvas())
13rb.setColor(Qt.green)
14rb.addPoint(QgsPointXY(pStart.x() - delta, pStart.y() - delta))
15rb.addPoint(QgsPointXY(pStart.x() + delta, pStart.y() - delta))
16rb.addPoint(QgsPointXY(pStart.x() + delta, pStart.y() + delta))
17rb.addPoint(QgsPointXY(pStart.x() - delta, pStart.y() + delta))
18
19tiedPoints = director.makeGraph(builder, [pStart])
20graph = builder.graph()
21tStart = tiedPoints[0]
22
23idStart = graph.findVertex(tStart)
24
25(tree, cost) = QgsGraphAnalyzer.dijkstra(graph, idStart, 0)
26
27upperBound = []
28r = 1500.0
29i = 0
30tree.reverse()
31
32while i < len(cost):
33    if cost[i] > r and tree[i] != -1:
34        outVertexId = graph.edge(tree [i]).toVertex()
35        if cost[outVertexId] < r:
36            upperBound.append(i)
37    i = i + 1
38
39for i in upperBound:
40    centerPoint = graph.vertex(i).point()
41    rb = QgsRubberBand(iface.mapCanvas())
42    rb.setColor(Qt.red)
43    rb.addPoint(QgsPointXY(centerPoint.x() - delta, centerPoint.y() - delta))
44    rb.addPoint(QgsPointXY(centerPoint.x() + delta, centerPoint.y() - delta))
45    rb.addPoint(QgsPointXY(centerPoint.x() + delta, centerPoint.y() + delta))
46    rb.addPoint(QgsPointXY(centerPoint.x() - delta, centerPoint.y() + delta))