Python, SciPy, Matplotlibでドロネー図・ボロノイ図をプロット
PythonではSciPyとMatplotlibを使ってドロネー図とボロノイ図を簡単にプロットできる。それぞれドロネー三角形やボロノイ領域の座標などを取得することも可能。
SciPyのドロネー図、ボロノイ図の算出にはQhullというライブラリが使用されている。
ここでは以下の内容について説明する。
- ドロネー図とボロノイ図
- サンプルコードの下準備
- ドロネー図をプロット:
scipy.spatial.delaunay_plot_2d - ドロネー三角形の座標リストを取得:
scipy.spatial.Delaunay - ボロノイ図をプロット:
scipy.spatial.voronoi_plot_2d - ボロノイ領域の座標リストを取得:
scipy.spatial.Voronoi - ボロノイ図を塗りつぶし
- ドロネー図とボロノイ図を同時にプロット
サンプルコードでの各ライブラリのバージョンは以下の通り。
- SciPy: 1.1.0
- Matplotlib: 2.2.2
OpenCVを使ってドロネー図およびボロノイ図を描画することも可能。
ドロネー図とボロノイ図
ドロネー図とボロノイ図についての説明はWikipediaなどを参照。英語版のほうが詳しい。
サンプルコードで最後に作成する図を先に示す。オレンジがドロネー図、黒がボロノイ図。
ドロネー図とは
平面においては、複数個の点に対して、各点を結ぶ三角形の最小角度が最大になるようにして分割したものがドロネー図（ドロネー三角形分割）。
ボロノイ図とは
平面においては、複数個の点に対して、それ以外の他の点がどの点に近いかによって領域分けされた図がボロノイ図。領域の境界線は各々の点の二等分線の一部となる。
サンプルコードの下準備
以下のように各ライブラリをインポートし、対象領域の幅と高さを指定、その領域内のランダムな点を用意する。この点をもとにドロネー図とボロノイ図をプロットする。
n個の点のx, y座標をn行2列のnumpy.ndarrayとする。
from scipy.spatial import Delaunay, delaunay_plot_2d, Voronoi, voronoi_plot_2d
import matplotlib.pyplot as plt
import numpy as np
w = h = 360
n = 6
np.random.seed(0)
pts = np.random.randint(0, w, (n, 2))
print(pts)
# [[172 47]
# [117 192]
# [323 251]
# [195 359]
# [ 9 211]
# [277 242]]
print(type(pts))
# <class 'numpy.ndarray'>
print(pts.shape)
# (6, 2)
ドロネー図をプロット: scipy.spatial.delaunay_plot_2d
まずDelaunay()の引数に複数点の座標を示すnumpy.ndarrayを指定しオブジェクトを生成。そのオブジェクトをdelaunay_plot_2d()の引数に指定するとドロネー図がプロットできる。
delaunay_plot_2d()が返すのはmatplotlib.figure.Figure。savefig()メソッドで画像として保存できる。
tri = Delaunay(pts)
print(type(tri))
# <class 'scipy.spatial.qhull.Delaunay'>
fig = delaunay_plot_2d(tri)
fig.savefig('data/dst/scipy_matplotlib_delaunay.png')
ドロネー三角形の座標リストを取得: scipy.spatial.Delaunay
Delaunay()で生成したオブジェクトからドロネー三角形の座標リストなどを取得できる。
ドロネー点はpoints属性で取得する。これはコンストラクタに与えた座標と等価。なお、pointsの型はfloat。
print(tri.points)
# [[172. 47.]
# [117. 192.]
# [323. 251.]
# [195. 359.]
# [ 9. 211.]
# [277. 242.]]
print(tri.points == pts)
# [[ True True]
# [ True True]
# [ True True]
# [ True True]
# [ True True]
# [ True True]]
どの点を結んで三角形としているかはsimplices属性で取得する。
print(tri.simplices)
# [[0 1 4]
# [1 3 4]
# [5 0 2]
# [5 1 0]
# [3 5 2]
# [5 3 1]]
simplicesの値はpointsの何番目の点かを示している。例えば[0 1 4]はpointsの0番目、1番目、4番目の点で三角形を構成していることを表している。座標リストを取得したい場合は、以下のように元のndarrayかpointsのインデックスに指定すればよい。
print(pts[tri.simplices])
# [[[172 47]
# [117 192]
# [ 9 211]]
#
# [[117 192]
# [195 359]
# [ 9 211]]
#
# [[277 242]
# [172 47]
# [323 251]]
#
# [[277 242]
# [117 192]
# [172 47]]
#
# [[195 359]
# [277 242]
# [323 251]]
#
# [[277 242]
# [195 359]
# [117 192]]]
ボロノイ図をプロット: scipy.spatial.voronoi_plot_2d
ボロノイ図もドロネー図と同じ流れ。
まずVoronoi()の引数に複数点の座標を示すnumpy.ndarrayを指定しオブジェクトを生成。そのオブジェクトをvoronoi_plot_2d()の引数に指定するとボロノイ図がプロットできる。
voronoi_plot_2d()が返すのはmatplotlib.figure.Figure。savefig()メソッドで画像として保存できる。
vor = Voronoi(pts)
print(type(vor))
# <class 'scipy.spatial.qhull.Voronoi'>
fig = voronoi_plot_2d(vor)
fig.savefig('data/dst/scipy_matplotlib_voronoi.png')
ボロノイ領域の座標リストを取得: scipy.spatial.Voronoi
Voronoi()で生成したオブジェクトから、ボロノイ点の座標やボロノイ領域を構成する座標リストなどを取得できる。
ボロノイ点（上図のオレンジ点）の座標はvertices属性で取得する。
print(vor.vertices)
# [[ 41.71518987 80.51265823]
# [ 82.09150528 310.02013526]
# [331.19719626 87.04766355]
# [218.67630058 147.63583815]
# [282.99117647 333.43398693]
# [182.6014461 263.07537248]]
ボロノイ領域を構成するボロノイ点の組み合わせはregions属性で取得する。
print(vor.regions)
# [[-1, 0, 1], [], [5, 3, 2, 4], [3, 0, -1, 2], [4, -1, 2], [5, 1, 0, 3], [5, 1, -1, 4]]
regionsの値はverticesの何番目の点かを示している。ここで、-1は無限遠点（領域外の点）。
領域内の点のみの組み合わせは以下のように取得可能。-1が含まれておらず空ではないリストをリスト内包表記で抽出する。
- 関連記事: Pythonリスト内包表記の使い方
print([r for r in vor.regions if -1 not in r and r])
# [[5, 3, 2, 4], [5, 1, 0, 3]]
これをverticesのインデックスに指定すればボロノイ領域を構成する点の座標が確認できる。
for region in [r for r in vor.regions if -1 not in r and r]:
print(vor.vertices[region])
# [[182.6014461 263.07537248]
# [218.67630058 147.63583815]
# [331.19719626 87.04766355]
# [282.99117647 333.43398693]]
# [[182.6014461 263.07537248]
# [ 82.09150528 310.02013526]
# [ 41.71518987 80.51265823]
# [218.67630058 147.63583815]]
ボロノイ図を塗りつぶし
ボロノイ領域を塗りつぶしたい場合、Matplotlibのfill()メソッドを使う。
voronoi_plot_2d()の第二引数に描画するmatplotlib.axes.Axesを指定できるので、同じAxesにfill()メソッドで塗りつぶしを行えばよい。
fill()の第一引数にx座標のリスト、第二引数にy座標のリスト、第三引数に塗りつぶし色を指定する。
fig, ax = plt.subplots()
voronoi_plot_2d(vor, ax)
for region, c in zip([r for r in vor.regions if -1 not in r and r], ['yellow', 'pink']):
ax.fill(vor.vertices[region][:, 0],
vor.vertices[region][:, 1],
color=c)
fig.savefig('data/dst/scipy_matplotlib_voronoi_fill.png')
外側の領域も塗りつぶしたい場合はOpenCVを使ったほうが簡単。
Matplotlibを使いたい場合は領域外の点の座標を算出する必要がある。voronoi_plot_2d()のソースコードなどが参考になる。
- [Spatial data structures and algorithms (scipy.spatial) — SciPy v1.1.0 Reference Guide](https://docs.scipy.org/doc/scipy/reference/tutorial/spatial.html
- scipy/_plotutils.py at master · scipy/scipy
ドロネー図とボロノイ図を同時にプロット
delaunay_plot_2d()の第二引数にも描画するmatplotlib.axes.Axesを指定可能。
delaunay_plot_2d()とvoronoi_plot_2d()に同じAxesを指定すればドロネー図とボロノイ図を重ねてプロットできる。
その他のグラフの設定も適宜行える。
fig, ax = plt.subplots(figsize=(4, 4))
delaunay_plot_2d(tri, ax)
voronoi_plot_2d(vor, ax, show_vertices=False)
ax.set_xlim(0, w)
ax.set_ylim(0, h)
ax.grid(linestyle='--')
fig.savefig('data/dst/scipy_matplotlib_delaunay_voronoi.png')