ネットワーク分析ライブラリ

改訂 ee19294562 から以降（QGIS >= 1.8）では、新しいネットワーク解析ライブラリがQGISコア解析ライブラリに追加されました。このライブラリは：

  • 地理データから数学的なグラフ（ポリラインベクターレイヤー）を作成します

  • グラフ理論からの基本的なメソッドを実装します（現在はダイクストラ法のみ）

ネットワーク解析ライブラリはRoadGraphコアプラグインから基本機能をエクスポートすることによって作成されました。今はそのメソッドをプラグインで、またはPythonのコンソールから直接使用できます。

一般情報

手短に言えば、一般的なユースケースは、次のように記述できます。

  1. 地理データから地理的なグラフ（たいていはポリラインベクターレイヤー）を作成する

  2. グラフ分析の実行

  3. 分析結果の利用（例えば、これらの可視化）

グラフを構築する

最初にする必要がある事は—入力データを準備することです、つまりベクターレイヤーをグラフに変換することです。これからのすべての操作は、レイヤーではなく、このグラフを使用します。

ソースとしてどんなポリラインベクターレイヤーも使用できます。ポリラインの頂点は、グラフの頂点となり、ポリラインのセグメントは、グラフの辺です。いくつかのノードが同じ座標を持っている場合、それらは同じグラフの頂点です。だから共通のノードを持つ2つの線は接続しています。

さらに、グラフの作成時には、入力ベクターレイヤーに好きな数だけ追加の点を「固定する」（「結びつける」）ことが可能です。追加の点それぞれについて、対応箇所—最も近いグラフの頂点または最も近いグラフの辺、が探し出されます。後者の場合、辺は分割されて新しい頂点が追加されるでしょう。

ベクターレイヤー属性と辺の長さは、辺の性質として使用できます。

Converting from a vector layer to the graph is done using the Builder programming pattern. A graph is constructed using a so-called Director. There is only one Director for now: QgsLineVectorLayerDirector. The director sets the basic settings that will be used to construct a graph from a line vector layer, used by the builder to create the graph. Currently, as in the case with the director, only one builder exists: QgsGraphBuilder, that creates QgsGraph objects. You may want to implement your own builders that will build a graphs compatible with such libraries as BGL or NetworkX.

To calculate edge properties the programming pattern strategy is used. For now only QgsDistanceArcProperter strategy is available, that takes into account the length of the route. You can implement your own strategy that will use all necessary parameters. For example, RoadGraph plugin uses a strategy that computes travel time using edge length and speed value from attributes.

では手順に行きましょう。

まず第一に、このライブラリを使用するために、networkanalysisモジュールをインポートする必要があります

from qgis.networkanalysis import *

それからディレクターを作成するためのいくつかの例

# don't use information about road direction from layer attributes,
# all roads are treated as two-way
director = QgsLineVectorLayerDirector(vLayer, -1, '', '', '', 3)

# 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 = QgsLineVectorLayerDirector(vLayer, 5, 'yes', '1', 'no', 3)

ディレクターを構築するために、ベクターレイヤーを渡さなければなりません。これはグラフ構造および各道路セグメント上で許される移動（一方向または双方向の動き、順または逆の方向）についての情報のソースとして使用されるでしょう。呼び出しは次のようになります

director = QgsLineVectorLayerDirector(vl, directionFieldId,
                                      directDirectionValue,
                                      reverseDirectionValue,
                                      bothDirectionValue,
                                      defaultDirection)

そして、ここでこれらのパラメーターは何を意味するかの完全なリストは以下のとおりです。

  • vl —グラフを構築するために使用される ベクターレイヤー

  • directionFieldId —道路の方向に関する情報が格納されている属性テーブルのフィールドのインデックス。 -1 、その後、まったくこの情報を使用しない場合。整数。

  • directDirectionValue —順方向（線の最初の点から最後の点へ移動）の道に対するフィールド値。文字列。

  • reverseDirectionValue —逆方向（線の最後の点から最初の点へ移動）の道に対するフィールド値。文字列。

  • bothDirectionValue —双方向道路に対するフィールド値（例えば最初の点から最後まで、また最後から最初まで移動できる道に対する）。文字列。

  • defaultDirection —デフォルトの道路の方向。この値は、フィールド directionFieldId が設定されていない、または上で指定された3つの値のいずれとも異なる値を有する道路に使用されるでしょう。整数です。 1 は順方向、 2 は逆方向、 3 は両方向を示します。

それから、辺性質を計算するための戦略を作成することが必要です

properter = QgsDistanceArcProperter()

そして、ディレクターにこの戦略について教えます

director.addProperter(properter)

これで、グラフを作成するビルダーを使用できます。QgsGraphBuilderクラスのコンストラクタには、いくつかの引数を取ります。

  • crs —使用する空間参照系。必須の引数です。

  • otfEnabled — “その場で” 再投影を使用するかどうか。デフォルトでは True (OTFを使用).

  • トポロジ許容値—トポロジ的な許容値です。デフォルト値は0です。

  • 楕円体ID — 使用する楕円体です。デフォルトは “WGS84”です。

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

分析に使用されるいくつかのポイントを定義することもできます。例えば

startPoint = QgsPoint(82.7112, 55.1672)
endPoint = QgsPoint(83.1879, 54.7079)

これですべてが整いましたので、グラフを構築し、それにこれらの点を「結びつける」ことができます。

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

グラフを構築するには（レイヤー中の地物数とレイヤーのサイズに応じて）いくらか時間がかかることがあります。 tiedPoints は「結びつけられた」点の座標を持つリストです。ビルド操作が完了するとグラフが得られ、それを分析のために使用できます

graph = builder.graph()

次のコードで、点の頂点インデックスを取得できます

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

グラフ分析

ネットワーク分析はこの二つの質問に対する答えを見つけるために使用されます：どの頂点が接続されているか、どのように最短経路を検索するか。これらの問題を解決するため、ネットワーク解析ライブラリではダイクストラのアルゴリズムを提供しています。

ダイクストラ法は、グラフの１つの頂点から他のすべての頂点への最短ルートおよび最適化パラメーターの値を見つけます。結果は、最短経路木として表現できます。

最短経路木は、次の性質を有する有向重み付きグラフ（より正確には、木）です：

  • 流入する辺がない頂点が1つだけあります - 木の根

  • 他のすべての頂点には流入する辺が１つだけあります

  • 頂点Bが頂点Aから到達可能である場合、このグラフ上のAからBへの経路は、単一の利用可能な経路であり、それは最適（最短）です

To get the shortest path tree use the methods shortestTree() and dijkstra() of QgsGraphAnalyzer class. It is recommended to use method dijkstra() because it works faster and uses memory more efficiently.

最短経路木の周りを歩くしたい場合 shortestTree() メソッドが便利です。それは、常に新しいグラフオブジェクト（QgsGraph）を作成し、三つの変数を受けとります。

  • source　— 入力グラフ

  • startVertexIdx — 木の上のポイントのインデックス（木の根）

  • criterionNum — 使用する辺プロパティの数（0から始まる）。

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

dijkstra() メソッドは引数は同じですが２つの配列を返します。第１の配列中の要素iには流入辺のインデックス、または流入辺が存在しない場合は-1が入ります。第２の配列の要素iには木の根から頂点iまでの距離が入ります。頂点iが根から到達不能である場合はDOUBLE_MAXが入ります。

(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() method (select linestring layer in TOC and replace coordinates with your own). Warning: use this code only as an example, it creates a lots of QgsRubberBand objects and may be slow on large data-sets.

from PyQt4.QtCore import *
from PyQt4.QtGui import *

from qgis.core import *
from qgis.gui import *
from qgis.networkanalysis import *

vl = qgis.utils.iface.mapCanvas().currentLayer()
director = QgsLineVectorLayerDirector(vl, -1, '', '', '', 3)
properter = QgsDistanceArcProperter()
director.addProperter(properter)
crs = qgis.utils.iface.mapCanvas().mapRenderer().destinationCrs()
builder = QgsGraphBuilder(crs)

pStart = QgsPoint(-0.743804, 0.22954)
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.arcCount()):
  rb = QgsRubberBand(qgis.utils.iface.mapCanvas())
  rb.setColor (Qt.red)
  rb.addPoint (tree.vertex(tree.arc(i).inVertex()).point())
  rb.addPoint (tree.vertex(tree.arc(i).outVertex()).point())
  i = i + 1

同じですが dijkstra() メソッドを使用して

from PyQt4.QtCore import *
from PyQt4.QtGui import *

from qgis.core import *
from qgis.gui import *
from qgis.networkanalysis import *

vl = qgis.utils.iface.mapCanvas().currentLayer()
director = QgsLineVectorLayerDirector(vl, -1, '', '', '', 3)
properter = QgsDistanceArcProperter()
director.addProperter(properter)
crs = qgis.utils.iface.mapCanvas().mapRenderer().destinationCrs()
builder = QgsGraphBuilder(crs)

pStart = QgsPoint(-1.37144, 0.543836)
tiedPoint = director.makeGraph(builder, [pStart])
pStart = tiedPoint[0]

graph = builder.graph()

idStart = graph.findVertex(pStart)

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

for edgeId in tree:
  if edgeId == -1:
    continue
  rb = QgsRubberBand(qgis.utils.iface.mapCanvas())
  rb.setColor (Qt.red)
  rb.addPoint (graph.vertex(graph.arc(edgeId).inVertex()).point())
  rb.addPoint (graph.vertex(graph.arc(edgeId).outVertex()).point())

最短経路を見つける

2点間の最適経路を見つけるためには以下のアプローチが使用されます。両方の点（始点Aおよび終点B）は構築されたときにグラフに「結ばれます」。それから、メソッド shortestTree() または dijkstra() を使用して、開始点Aに根を持つ最短経路木を構築します。また、同じ木の中に終点Bを見つけ、木を点Aから点Bまで歩き始めます。全体のアルゴリズムは次のように書けます

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

この時点において、この経路で走行中に訪問される頂点の反転リストの形（頂点は逆順で終点から始点へと列挙されている）で、経路が得られます。

これはメソッド shortestTree() を使用するQGIS Pythonコンソールのためのサンプルコードです（TOC中のラインストリングレイヤーを選択して、コードの座標を自分のものに置き換える必要があります）

from PyQt4.QtCore import *
from PyQt4.QtGui import *

from qgis.core import *
from qgis.gui import *
from qgis.networkanalysis import *

vl = qgis.utils.iface.mapCanvas().currentLayer()
director = QgsLineVectorLayerDirector(vl, -1, '', '', '', 3)
properter = QgsDistanceArcProperter()
director.addProperter(properter)
crs = qgis.utils.iface.mapCanvas().mapRenderer().destinationCrs()
builder = QgsGraphBuilder(crs)

pStart = QgsPoint(-0.835953, 0.15679)
pStop = QgsPoint(-1.1027, 0.699986)

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

tStart = tiedPoints[0]
tStop = tiedPoints[1]

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

idStart = tree.findVertex(tStart)
idStop = tree.findVertex(tStop)

if idStop == -1:
  print "Path not found"
else:
  p = []
  while (idStart != idStop):
    l = tree.vertex(idStop).inArc()
    if len(l) == 0:
      break
    e = tree.arc(l[0])
    p.insert(0, tree.vertex(e.inVertex()).point())
    idStop = e.outVertex()

  p.insert(0, tStart)
  rb = QgsRubberBand(qgis.utils.iface.mapCanvas())
  rb.setColor(Qt.red)

  for pnt in p:
    rb.addPoint(pnt)

そしてこれは同じサンプルですが dijkstra() メソッドを使用しています

from PyQt4.QtCore import *
from PyQt4.QtGui import *

from qgis.core import *
from qgis.gui import *
from qgis.networkanalysis import *

vl = qgis.utils.iface.mapCanvas().currentLayer()
director = QgsLineVectorLayerDirector(vl, -1, '', '', '', 3)
properter = QgsDistanceArcProperter()
director.addProperter(properter)
crs = qgis.utils.iface.mapCanvas().mapRenderer().destinationCrs()
builder = QgsGraphBuilder(crs)

pStart = QgsPoint(-0.835953, 0.15679)
pStop = QgsPoint(-1.1027, 0.699986)

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

tStart = tiedPoints[0]
tStop = tiedPoints[1]

idStart = graph.findVertex(tStart)
idStop = graph.findVertex(tStop)

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

if tree[idStop] == -1:
  print "Path not found"
else:
  p = []
  curPos = idStop
  while curPos != idStart:
    p.append(graph.vertex(graph.arc(tree[curPos]).inVertex()).point())
    curPos = graph.arc(tree[curPos]).outVertex();

  p.append(tStart)

  rb = QgsRubberBand(qgis.utils.iface.mapCanvas())
  rb.setColor(Qt.red)

  for pnt in p:
    rb.addPoint(pnt)

利用可能領域

頂点Aに対する利用可能領域とは、頂点Aから到達可能であり、Aからこれらの頂点までの経路のコストがある値以下になるような、グラフの頂点の部分集合です。

より明確に、これは次の例で示すことができます。「消防署があります。消防車が5分、10分、15分で到達できるのは市内のどの部分ですか？」。これらの質問への回答が消防署の利用可能領域です。

利用可能領域を見つけるためには、 QgsGraphAnalyzer クラスのメソッド dijkstra() を使用できます。事前に定義された値とコスト配列の要素を比較すれば十分です。コスト[i]が所定値以下の場合は、頂点iは利用可能領域内に、そうでない場合は利用可能領域外にあります。

より難しい問題は、利用可能領域の境界を取得することです。下限はまだ到達できる頂点の集合であり、上限は到達できない頂点の集合です。実際にはこれは簡単です：それは、辺の元頂点が到達可能であり辺の先頂点が到達可能ではないような、最短経路木の辺に基づく利用可能境界です。

例です

from PyQt4.QtCore import *
from PyQt4.QtGui import *

from qgis.core import *
from qgis.gui import *
from qgis.networkanalysis import *

vl = qgis.utils.iface.mapCanvas().currentLayer()
director = QgsLineVectorLayerDirector(vl, -1, '', '', '', 3)
properter = QgsDistanceArcProperter()
director.addProperter(properter)
crs = qgis.utils.iface.mapCanvas().mapRenderer().destinationCrs()
builder = QgsGraphBuilder(crs)

pStart = QgsPoint(65.5462, 57.1509)
delta = qgis.utils.iface.mapCanvas().getCoordinateTransform().mapUnitsPerPixel() * 1

rb = QgsRubberBand(qgis.utils.iface.mapCanvas(), True)
rb.setColor(Qt.green)
rb.addPoint(QgsPoint(pStart.x() - delta, pStart.y() - delta))
rb.addPoint(QgsPoint(pStart.x() + delta, pStart.y() - delta))
rb.addPoint(QgsPoint(pStart.x() + delta, pStart.y() + delta))
rb.addPoint(QgsPoint(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 = 2000.0
i = 0
while i < len(cost):
  if cost[i] > r and tree[i] != -1:
    outVertexId = graph.arc(tree [i]).outVertex()
    if cost[outVertexId] < r:
      upperBound.append(i)
  i = i + 1

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