包含一组给定点的边界
Boundary enclosing a given set of points
我目前使用的算法有点问题。我想让它成为一个边界。
这是当前行为的示例:
这是一个 MSPaint 所需行为的示例:
C#中凸包的当前代码:https://hastebin.com/dudejesuja.cs
所以这是我的问题:
1) 这可能吗?
R: 是
2) 这还叫凸包吗? (我不这么认为)
R: 不,它叫做边界,link: https://www.mathworks.com/help/matlab/ref/boundary.html
3) 这会不会比传统的凸包更不友好?
R:嗯,据我研究应该是一样的性能
4) 该算法的伪代码或类似代码示例?
R:还没有回答或者我还没有找到解决办法
我会使用不同的方法来解决这个问题。由于我们使用的是一组二维点,因此可以直接计算点区域的边界矩形。然后我会通过水平线和垂直线将这个矩形分成“单元格”,对于每个单元格,只需计算位于其边界内的像素数。由于每个单元格只能有 4 个相邻单元格(通过单元格边相邻),因此边界单元格将是至少有一个空相邻单元格或单元格边位于边界矩形边界的单元格。然后将沿着边界单元格边构建边界。边界看起来像一个“楼梯”,但选择较小的像元大小会改善结果。事实上,细胞大小应该通过实验确定;它不能太小,否则区域内可能会出现空单元格。点之间的平均距离可以用作像元大小的下边界。
考虑使用 Alpha 形状,有时称为凹壳。 https://en.wikipedia.org/wiki/Alpha_shape
它可以从 Delaunay 三角剖分中构建,时间为 O(N log N)。
这里是一些 Python 代码,用于计算 alpha 形状(凹包)并仅保留外边界。这大概就是matlab的boundary在里面做的事情。
from scipy.spatial import Delaunay
import numpy as np
def alpha_shape(points, alpha, only_outer=True):
"""
Compute the alpha shape (concave hull) of a set of points.
:param points: np.array of shape (n,2) points.
:param alpha: alpha value.
:param only_outer: boolean value to specify if we keep only the outer border
or also inner edges.
:return: set of (i,j) pairs representing edges of the alpha-shape. (i,j) are
the indices in the points array.
"""
assert points.shape[0] > 3, "Need at least four points"
def add_edge(edges, i, j):
"""
Add an edge between the i-th and j-th points,
if not in the list already
"""
if (i, j) in edges or (j, i) in edges:
# already added
assert (j, i) in edges, "Can't go twice over same directed edge right?"
if only_outer:
# if both neighboring triangles are in shape, it's not a boundary edge
edges.remove((j, i))
return
edges.add((i, j))
tri = Delaunay(points)
edges = set()
# Loop over triangles:
# ia, ib, ic = indices of corner points of the triangle
for ia, ib, ic in tri.vertices:
pa = points[ia]
pb = points[ib]
pc = points[ic]
# Computing radius of triangle circumcircle
# www.mathalino.com/reviewer/derivation-of-formulas/derivation-of-formula-for-radius-of-circumcircle
a = np.sqrt((pa[0] - pb[0]) ** 2 + (pa[1] - pb[1]) ** 2)
b = np.sqrt((pb[0] - pc[0]) ** 2 + (pb[1] - pc[1]) ** 2)
c = np.sqrt((pc[0] - pa[0]) ** 2 + (pc[1] - pa[1]) ** 2)
s = (a + b + c) / 2.0
area = np.sqrt(s * (s - a) * (s - b) * (s - c))
circum_r = a * b * c / (4.0 * area)
if circum_r < alpha:
add_edge(edges, ia, ib)
add_edge(edges, ib, ic)
add_edge(edges, ic, ia)
return edges
如果你 运行 它与下面的测试代码你会得到这个数字,它看起来像你需要的:
from matplotlib.pyplot import *
# Constructing the input point data
np.random.seed(0)
x = 3.0 * np.random.rand(2000)
y = 2.0 * np.random.rand(2000) - 1.0
inside = ((x ** 2 + y ** 2 > 1.0) & ((x - 3) ** 2 + y ** 2 > 1.0)
points = np.vstack([x[inside], y[inside]]).T
# Computing the alpha shape
edges = alpha_shape(points, alpha=0.25, only_outer=True)
# Plotting the output
figure()
axis('equal')
plot(points[:, 0], points[:, 1], '.')
for i, j in edges:
plot(points[[i, j], 0], points[[i, j], 1])
show()
编辑:根据评论中的请求,这里有一些代码 "stitches" 将输出边设置为连续边序列。
def find_edges_with(i, edge_set):
i_first = [j for (x,j) in edge_set if x==i]
i_second = [j for (j,x) in edge_set if x==i]
return i_first,i_second
def stitch_boundaries(edges):
edge_set = edges.copy()
boundary_lst = []
while len(edge_set) > 0:
boundary = []
edge0 = edge_set.pop()
boundary.append(edge0)
last_edge = edge0
while len(edge_set) > 0:
i,j = last_edge
j_first, j_second = find_edges_with(j, edge_set)
if j_first:
edge_set.remove((j, j_first[0]))
edge_with_j = (j, j_first[0])
boundary.append(edge_with_j)
last_edge = edge_with_j
elif j_second:
edge_set.remove((j_second[0], j))
edge_with_j = (j, j_second[0]) # flip edge rep
boundary.append(edge_with_j)
last_edge = edge_with_j
if edge0[0] == last_edge[1]:
break
boundary_lst.append(boundary)
return boundary_lst
然后您可以遍历边界列表的列表并在每条边中追加与第一个索引对应的点以获得边界多边形。
这里是 JavaScript 构建凹包的代码:https://github.com/AndriiHeonia/hull 也许你可以将它移植到 C#。
正如之前大多数专家所指出的,这可能不是凸包而是凹包,或者说是Alpha Shape。 Iddo 提供了一个干净的 Python 代码来获取这个形状。然而,你也可以直接利用一些现有的包来实现这一点,如果你正在处理大量的点云,也许速度更快,计算内存更少。
[1] Alpha Shape Toolbox:生成n维alpha形状的工具箱。
[2] Plotly:它可以生成一个 Mesh3d 对象,根据键值可以是该集合的凸包、Delaunay 三角剖分或 alpha 集合。
我目前使用的算法有点问题。我想让它成为一个边界。
这是当前行为的示例:
这是一个 MSPaint 所需行为的示例:
C#中凸包的当前代码:https://hastebin.com/dudejesuja.cs
所以这是我的问题:
1) 这可能吗?
R: 是
2) 这还叫凸包吗? (我不这么认为)
R: 不,它叫做边界,link: https://www.mathworks.com/help/matlab/ref/boundary.html
3) 这会不会比传统的凸包更不友好?
R:嗯,据我研究应该是一样的性能
4) 该算法的伪代码或类似代码示例?
R:还没有回答或者我还没有找到解决办法
我会使用不同的方法来解决这个问题。由于我们使用的是一组二维点,因此可以直接计算点区域的边界矩形。然后我会通过水平线和垂直线将这个矩形分成“单元格”,对于每个单元格,只需计算位于其边界内的像素数。由于每个单元格只能有 4 个相邻单元格(通过单元格边相邻),因此边界单元格将是至少有一个空相邻单元格或单元格边位于边界矩形边界的单元格。然后将沿着边界单元格边构建边界。边界看起来像一个“楼梯”,但选择较小的像元大小会改善结果。事实上,细胞大小应该通过实验确定;它不能太小,否则区域内可能会出现空单元格。点之间的平均距离可以用作像元大小的下边界。
考虑使用 Alpha 形状,有时称为凹壳。 https://en.wikipedia.org/wiki/Alpha_shape
它可以从 Delaunay 三角剖分中构建,时间为 O(N log N)。
这里是一些 Python 代码,用于计算 alpha 形状(凹包)并仅保留外边界。这大概就是matlab的boundary在里面做的事情。
from scipy.spatial import Delaunay
import numpy as np
def alpha_shape(points, alpha, only_outer=True):
"""
Compute the alpha shape (concave hull) of a set of points.
:param points: np.array of shape (n,2) points.
:param alpha: alpha value.
:param only_outer: boolean value to specify if we keep only the outer border
or also inner edges.
:return: set of (i,j) pairs representing edges of the alpha-shape. (i,j) are
the indices in the points array.
"""
assert points.shape[0] > 3, "Need at least four points"
def add_edge(edges, i, j):
"""
Add an edge between the i-th and j-th points,
if not in the list already
"""
if (i, j) in edges or (j, i) in edges:
# already added
assert (j, i) in edges, "Can't go twice over same directed edge right?"
if only_outer:
# if both neighboring triangles are in shape, it's not a boundary edge
edges.remove((j, i))
return
edges.add((i, j))
tri = Delaunay(points)
edges = set()
# Loop over triangles:
# ia, ib, ic = indices of corner points of the triangle
for ia, ib, ic in tri.vertices:
pa = points[ia]
pb = points[ib]
pc = points[ic]
# Computing radius of triangle circumcircle
# www.mathalino.com/reviewer/derivation-of-formulas/derivation-of-formula-for-radius-of-circumcircle
a = np.sqrt((pa[0] - pb[0]) ** 2 + (pa[1] - pb[1]) ** 2)
b = np.sqrt((pb[0] - pc[0]) ** 2 + (pb[1] - pc[1]) ** 2)
c = np.sqrt((pc[0] - pa[0]) ** 2 + (pc[1] - pa[1]) ** 2)
s = (a + b + c) / 2.0
area = np.sqrt(s * (s - a) * (s - b) * (s - c))
circum_r = a * b * c / (4.0 * area)
if circum_r < alpha:
add_edge(edges, ia, ib)
add_edge(edges, ib, ic)
add_edge(edges, ic, ia)
return edges
如果你 运行 它与下面的测试代码你会得到这个数字,它看起来像你需要的:
from matplotlib.pyplot import *
# Constructing the input point data
np.random.seed(0)
x = 3.0 * np.random.rand(2000)
y = 2.0 * np.random.rand(2000) - 1.0
inside = ((x ** 2 + y ** 2 > 1.0) & ((x - 3) ** 2 + y ** 2 > 1.0)
points = np.vstack([x[inside], y[inside]]).T
# Computing the alpha shape
edges = alpha_shape(points, alpha=0.25, only_outer=True)
# Plotting the output
figure()
axis('equal')
plot(points[:, 0], points[:, 1], '.')
for i, j in edges:
plot(points[[i, j], 0], points[[i, j], 1])
show()
编辑:根据评论中的请求,这里有一些代码 "stitches" 将输出边设置为连续边序列。
def find_edges_with(i, edge_set):
i_first = [j for (x,j) in edge_set if x==i]
i_second = [j for (j,x) in edge_set if x==i]
return i_first,i_second
def stitch_boundaries(edges):
edge_set = edges.copy()
boundary_lst = []
while len(edge_set) > 0:
boundary = []
edge0 = edge_set.pop()
boundary.append(edge0)
last_edge = edge0
while len(edge_set) > 0:
i,j = last_edge
j_first, j_second = find_edges_with(j, edge_set)
if j_first:
edge_set.remove((j, j_first[0]))
edge_with_j = (j, j_first[0])
boundary.append(edge_with_j)
last_edge = edge_with_j
elif j_second:
edge_set.remove((j_second[0], j))
edge_with_j = (j, j_second[0]) # flip edge rep
boundary.append(edge_with_j)
last_edge = edge_with_j
if edge0[0] == last_edge[1]:
break
boundary_lst.append(boundary)
return boundary_lst
然后您可以遍历边界列表的列表并在每条边中追加与第一个索引对应的点以获得边界多边形。
这里是 JavaScript 构建凹包的代码:https://github.com/AndriiHeonia/hull 也许你可以将它移植到 C#。
正如之前大多数专家所指出的,这可能不是凸包而是凹包,或者说是Alpha Shape。 Iddo 提供了一个干净的 Python 代码来获取这个形状。然而,你也可以直接利用一些现有的包来实现这一点,如果你正在处理大量的点云,也许速度更快,计算内存更少。
[1] Alpha Shape Toolbox:生成n维alpha形状的工具箱。
[2] Plotly:它可以生成一个 Mesh3d 对象,根据键值可以是该集合的凸包、Delaunay 三角剖分或 alpha 集合。