使用 Scipy 在凸包上找到一个点的投影

Question

从一组点来看，我得到了scipy.spatial带有Delaunay或的凸包ConvexHull（来自 qhull 库）。现在我想将这个凸包外的点投影到船体上（即船体上与外部点距离最小的点）。

这是我到目前为止的代码：

from scipy.spatial import Delaunay, ConvexHull
import numpy as np

hu = np.random.rand(10, 2) ## the set of points to get the hull from
pt = np.array([1.1, 0.5]) ## a point outside
pt2 = np.array([0.4, 0.4]) ## a point inside
hull = ConvexHull(hu) ## get only the convex hull
#hull2 = Delaunay(hu) ## or get the full Delaunay triangulation

import matplotlib.pyplot as plt
plt.plot(hu[:,0], hu[:,1], "ro") ## plot all points
#plt.triplot(hu[:,0], hu[:,1], hull2.simplices.copy()) ## plot the Delaunay triangulation
## Plot the convexhull
for simplex in hull.simplices:
    plt.plot(hu[simplex,0], hu[simplex,1], "ro-")

## Plot the points inside and outside the convex hull 
plt.plot(pt[0], pt[1], "bs")
plt.plot(pt2[0], pt2[1], "bs")
plt.show()

使用图片可能更容易，我想从凸包外的蓝点获得绿色的 x 和 y 坐标。这个例子是二维的，但我也需要在更高的维度上应用它。谢谢您的帮助。

编辑：这里解决了问题，但我在实现它时遇到了麻烦：https : //mathoverflow.net/questions/118088/projection-of-a-point-to-a-convex-hull-in-d-dimensions

Answer 1

我在回答自己。正如 0Tech 指出的那样，ConvexHull.equations为您提供每个平面的平面方程（在 2d 中——因此是一条线），形式为 : [A, B, C]。因此平面定义为

A*x + B*y + C = 0

在平面上投影点 P=(x0, y0) 在这里解释。您希望向量上的一个点平行于平面向量 (A, B) 并通过该点投影 P，这条线由 t 参数化：

P_proj = (x, y) = (x0 + A*t, y0 + B*t)

然后您希望您的点在平面上并使用完整平面方程来做到这一点：

A*(x0 + A*t) + B*(y0 + B*t) + C = 0
=> t=-(C + A*x0 + B*y0)/(A**2+B**2)

在（笨拙的）python 中，它给出了任何维度：

from scipy.spatial import Delaunay, ConvexHull
import numpy as np

hu = np.random.rand(10, 2) ## the set of points to get the hull from
pt = np.array([1.1, 0.5]) ## a point outside
pt2 = np.array([0.4, 0.4]) ## a point inside
hull = ConvexHull(hu) ## get only the convex hull
#hull2 = Delaunay(hu) ## or get the full Delaunay triangulation

import matplotlib.pyplot as plt
plt.plot(hu[:,0], hu[:,1], "ro") ## plot all points
#plt.triplot(hu[:,0], hu[:,1], hull2.simplices.copy()) ## plot the Delaunay triangulation
## Plot the convexhull
for simplex in hull.simplices:
    plt.plot(hu[simplex,0], hu[simplex,1], "ro-")

## Plot the points inside and outside the convex hull 
plt.plot(pt[0], pt[1], "bs")
plt.plot(pt2[0], pt2[1], "bs")

for eq in hull.equations:
    t = -(eq[-1] + np.dot(eq[:-1], pt))/(np.sum(eq[:-1]**2))
    pt_proj = pt + eq[:-1]*t
    plt.plot(pt_proj[0], pt_proj[1], "gD-.")
plt.show()

浏览stackoverflow，让我想到了另一种解决方案，它具有使用段而不是线的优势，因此其中一个段的投影始终位于该段上：

def min_distance(pt1, pt2, p):
    """ return the projection of point p (and the distance) on the closest edge formed by the two points pt1 and pt2"""
    l = np.sum((pt2-pt1)**2) ## compute the squared distance between the 2 vertices
    t = np.max([0., np.min([1., np.dot(p-pt1, pt2-pt1) /l])]) # I let the answer of question 849211 explains this
    proj = pt1 + t*(pt2-pt1) ## project the point
    return proj, np.sum((proj-p)**2) ## return the projection and the point

然后我们可以浏览每个顶点并投影点：

for i in range(len(hull.vertices)):
    pt_proj, d = min_distance(hu[hull.vertices[i]], hu[hull.vertices[(i+1)%len(hull.vertices)]], pt)
    plt.plot([pt[0], pt_proj[0]], [pt[1], pt_proj[1]], "c<:")

和图片，每个平面（线）右侧的蓝点投影为绿色（第一种方法）和青色（第二种方法）：