Раскрасить диаграмму Вороного
Я пытаюсь раскрасить диаграмму Вороного, созданную с помощью scipy.spatial.Voronoi
. Здесь мой код:
import numpy as np
import matplotlib.pyplot as plt
from scipy.spatial import Voronoi, voronoi_plot_2d
# make up data points
points = np.random.rand(15,2)
# compute Voronoi tesselation
vor = Voronoi(points)
# plot
voronoi_plot_2d(vor)
# colorize
for region in vor.regions:
if not -1 in region:
polygon = [vor.vertices[i] for i in region]
plt.fill(*zip(*polygon))
plt.show()
Полученное изображение:
![Voronoi Diagram]()
Как вы можете видеть, некоторые районы Вороного на границе изображения не окрашены. Это связано с тем, что некоторые индексы к вершинам Вороного для этих областей установлены на -1
, т.е. Для тех вершин вне диаграммы Вороного. Согласно документам:
области: (список списка ints, shape (nregions, *)) Индексы верных Вороных, образующих каждую область Вороного. -1 указывает вершину вне диаграммы Вороного.
Чтобы раскрасить эти регионы, я попытался просто удалить эти "внешние" вершины из многоугольника, но это не сработало. Я думаю, мне нужно заполнить некоторые точки на границе области изображения, но я не могу понять, как это сделать разумно.
Может ли кто-нибудь помочь?
Ответы
Ответ 1
Структура данных Voronoi содержит всю необходимую информацию для построения позиций для "точек на бесконечности". Qhull также сообщает их просто как индексы -1
, поэтому Scipy не вычисляет их для вас.
https://gist.github.com/pv/8036995
http://nbviewer.ipython.org/gist/pv/8037100
import numpy as np
import matplotlib.pyplot as plt
from scipy.spatial import Voronoi
def voronoi_finite_polygons_2d(vor, radius=None):
"""
Reconstruct infinite voronoi regions in a 2D diagram to finite
regions.
Parameters
----------
vor : Voronoi
Input diagram
radius : float, optional
Distance to 'points at infinity'.
Returns
-------
regions : list of tuples
Indices of vertices in each revised Voronoi regions.
vertices : list of tuples
Coordinates for revised Voronoi vertices. Same as coordinates
of input vertices, with 'points at infinity' appended to the
end.
"""
if vor.points.shape[1] != 2:
raise ValueError("Requires 2D input")
new_regions = []
new_vertices = vor.vertices.tolist()
center = vor.points.mean(axis=0)
if radius is None:
radius = vor.points.ptp().max()
# Construct a map containing all ridges for a given point
all_ridges = {}
for (p1, p2), (v1, v2) in zip(vor.ridge_points, vor.ridge_vertices):
all_ridges.setdefault(p1, []).append((p2, v1, v2))
all_ridges.setdefault(p2, []).append((p1, v1, v2))
# Reconstruct infinite regions
for p1, region in enumerate(vor.point_region):
vertices = vor.regions[region]
if all(v >= 0 for v in vertices):
# finite region
new_regions.append(vertices)
continue
# reconstruct a non-finite region
ridges = all_ridges[p1]
new_region = [v for v in vertices if v >= 0]
for p2, v1, v2 in ridges:
if v2 < 0:
v1, v2 = v2, v1
if v1 >= 0:
# finite ridge: already in the region
continue
# Compute the missing endpoint of an infinite ridge
t = vor.points[p2] - vor.points[p1] # tangent
t /= np.linalg.norm(t)
n = np.array([-t[1], t[0]]) # normal
midpoint = vor.points[[p1, p2]].mean(axis=0)
direction = np.sign(np.dot(midpoint - center, n)) * n
far_point = vor.vertices[v2] + direction * radius
new_region.append(len(new_vertices))
new_vertices.append(far_point.tolist())
# sort region counterclockwise
vs = np.asarray([new_vertices[v] for v in new_region])
c = vs.mean(axis=0)
angles = np.arctan2(vs[:,1] - c[1], vs[:,0] - c[0])
new_region = np.array(new_region)[np.argsort(angles)]
# finish
new_regions.append(new_region.tolist())
return new_regions, np.asarray(new_vertices)
# make up data points
np.random.seed(1234)
points = np.random.rand(15, 2)
# compute Voronoi tesselation
vor = Voronoi(points)
# plot
regions, vertices = voronoi_finite_polygons_2d(vor)
print "--"
print regions
print "--"
print vertices
# colorize
for region in regions:
polygon = vertices[region]
plt.fill(*zip(*polygon), alpha=0.4)
plt.plot(points[:,0], points[:,1], 'ko')
plt.xlim(vor.min_bound[0] - 0.1, vor.max_bound[0] + 0.1)
plt.ylim(vor.min_bound[1] - 0.1, vor.max_bound[1] + 0.1)
plt.show()
![введите описание изображения здесь]()
Ответ 2
Я не думаю, что из данных, доступных в структуре vor, достаточно информации, чтобы понять это, не выполняя, по крайней мере, некоторые вычисления voronoi. Так как это случай, вот соответствующие части исходной функции voronoi_plot_2d, которую вы должны использовать, чтобы извлечь точки, которые пересекаются с vor.max_bound или vor.min_bound, которые являются нижним левым и верхним правыми углами диаграммы в укажите другие координаты для ваших полигонов.
for simplex in vor.ridge_vertices:
simplex = np.asarray(simplex)
if np.all(simplex >= 0):
ax.plot(vor.vertices[simplex,0], vor.vertices[simplex,1], 'k-')
ptp_bound = vor.points.ptp(axis=0)
center = vor.points.mean(axis=0)
for pointidx, simplex in zip(vor.ridge_points, vor.ridge_vertices):
simplex = np.asarray(simplex)
if np.any(simplex < 0):
i = simplex[simplex >= 0][0] # finite end Voronoi vertex
t = vor.points[pointidx[1]] - vor.points[pointidx[0]] # tangent
t /= np.linalg.norm(t)
n = np.array([-t[1], t[0]]) # normal
midpoint = vor.points[pointidx].mean(axis=0)
direction = np.sign(np.dot(midpoint - center, n)) * n
far_point = vor.vertices[i] + direction * ptp_bound.max()
ax.plot([vor.vertices[i,0], far_point[0]],
[vor.vertices[i,1], far_point[1]], 'k--')
Ответ 3
Просто продолжение решения этой проблемы. Я пытался повторить результат, представленный выше пользователем PV. Результат, который я получаю при запуске его кода: ![enter image description here]()
Мне интересно, было ли обновление qhull или scipy.spatial. Может быть, моя установка не работает, я не знаю.
В любом случае, поскольку у меня та же проблема, и я не могу заставить работать приведенный выше код, у меня есть альтернативное решение. Добавляя кольцо призрачных точек вокруг вашего распределения точек, вы получаете ячейки Вороного, которые выходят за ваши границы, и вы выдвигаете проблему за пределы области построения.
n = 15
points = np.random.rand(n, 2)
# Determine the radius of the circle (in this case the diagonal)
radius = np.sqrt(2.)
# Determine the number of points on the circle
# I use the number of points per dimension for each side
n_circle = 2 * n
for i in range(n_circle):
angle = i * 2*np.pi/n_circle
points = np.concatenate((points, [[radius*np.cos(angle)+0.5,radius*np.sin(angle)+0.5]]))
vor = Voronoi(points)
for r, p in enumerate(vor.point_region):
region = vor.regions[p]
if not -1 in region and r < n:
polygon = [vor.vertices[i] for i in region]
plt.fill(*zip(*polygon))
plt.xlim(0, 1)
plt.ylim(0, 1)
Вот примеры без кругового распределения и с круговым распределением точек:
![enter image description here]()
![enter image description here]()
Ответ 4
У меня есть намного более простое решение этой проблемы, то есть добавить 4 удаленных фиктивных точки в ваш список точек перед вызовом алгоритма Вороного.
Основываясь на ваших кодах, я добавил две строки.
import numpy as np
import matplotlib.pyplot as plt
from scipy.spatial import Voronoi, voronoi_plot_2d
# make up data points
points = np.random.rand(15,2)
# add 4 distant dummy points
points = np.append(points, [[999,999], [-999,999], [999,-999], [-999,-999]], axis = 0)
# compute Voronoi tesselation
vor = Voronoi(points)
# plot
voronoi_plot_2d(vor)
# colorize
for region in vor.regions:
if not -1 in region:
polygon = [vor.vertices[i] for i in region]
plt.fill(*zip(*polygon))
# fix the range of axes
plt.xlim([0,1]), plt.ylim([0,1])
plt.show()
Тогда полученная цифра выглядит следующим образом. ![enter image description here]()