Python 着色 Voronoi 图
声明：本页面是StackOverFlow热门问题的中英对照翻译，遵循CC BY-SA 4.0协议，如果您需要使用它，必须同样遵循CC BY-SA许可，注明原文地址和作者信息，同时你必须将它归于原作者（不是我）：StackOverFlow
原文地址: http://stackoverflow.com/questions/20515554/
Warning: these are provided under cc-by-sa 4.0 license. You are free to use/share it, But you must attribute it to the original authors (not me):
StackOverFlow
Colorize Voronoi Diagram
提问by moooeeeep
I'm trying to colorize a Voronoi Diagram created using scipy.spatial.Voronoi. Here's my code:
我正在尝试为使用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()
The resulting image:
结果图像：
As you can see some of the Voronoi regions at the border of the image are not colored. That is because some indices to the Voronoi vertices for these regions are set to -1, i.e., for those vertices outside the Voronoi diagram. According to the docs:
正如您所看到的，图像边界的一些 Voronoi 区域没有着色。那是因为这些区域的 Voronoi 顶点的一些索引被设置为-1，即 Voronoi 图之外的那些顶点。根据文档：
regions:(list of list of ints, shape (nregions, *)) Indices of the Voronoi vertices forming each Voronoi region. -1 indicates vertex outside the Voronoi diagram.
区域：（整数列表，形状（nregions，*））形成每个 Voronoi 区域的 Voronoi 顶点的索引。-1 表示 Voronoi 图之外的顶点。
In order to colorize these regions as well, I've tried to just remove these "outside" vertices from the polygon, but that didn't work. I think, I need to fill in some points at the border of the image region, but I can't seem to figure out how to achieve this reasonably.
为了也为这些区域着色，我试图从多边形中删除这些“外部”顶点，但这不起作用。我想，我需要在图像区域的边界填充一些点，但我似乎无法弄清楚如何合理地实现这一点。
Can anyone help?
任何人都可以帮忙吗？
采纳答案by pv.
The Voronoi data structure contains all the necessary information to construct positions for the "points at infinity". Qhull also reports them simply as -1indices, so Scipy doesn't compute them for you.
Voronoi 数据结构包含构建“无穷远点”位置的所有必要信息。Qhull 还简单地将它们报告为-1索引，因此 Scipy 不会为您计算它们。
https://gist.github.com/pv/8036995
https://gist.github.com/pv/8036995
http://nbviewer.ipython.org/gist/pv/8037100
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()
回答by Troy Rockwood
I don't think there is enough information from the data available in the vor structure to figure this out without doing at least some of the voronoi computation again. Since that's the case, here are the relevant portions of the original voronoi_plot_2d function that you should be able to use to extract the points that intersect with the vor.max_bound or vor.min_bound which are the bottom left and top right corners of the diagram in order figure out the other coordinates for your polygons.
我认为 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--')
回答by Arrows
I have a much simpler solution to this problem, that is to add 4 distant dummy points to your point list before calling the Voronoi algorithm.
对于这个问题，我有一个更简单的解决方案，即在调用 Voronoi 算法之前将 4 个遥远的虚拟点添加到您的点列表中。
Based on your codes, I added two lines.
根据您的代码，我添加了两行。
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()
Then the resulting figure just looks like the following.
然后生成的图形如下所示。
