Python:离直线最近的点

Python: Closest Point to a line

我有以下问题。我有一个装满坐标和三个点的盒子,它们构成了一条线。现在我想计算所有框坐标到该线的最短距离。我有三种方法可以做到这一点,vtk 和 numpy 版本总是有相同的结果,但不是 shapely 的距离方法。但我需要形状优美的版本,因为我想测量从一个点到整条线的最近距离,而不是到单独的线段。到目前为止,这是一个示例代码。所以问题是“pdist”:

from shapely.geometry import LineString, Point
import vtk, numpy as np
import itertools

import math

from numpy.linalg import norm

x1=np.arange(4,21)
y1=np.arange(4,21)
z1=np.arange(-7,6)

linepoints = np.array([[1,10,0],[10,10,0],[15,15,0]])


for i in itertools.product(x1,y1,z1):
    
    for m in range(len(linepoints)-1):
        
        line3 = LineString([linepoints[m],linepoints[m+1]])
        
        p = Point(i)
        
        d = norm(np.cross(linepoints[m]-linepoints[m+1], linepoints[m]-i))/norm(linepoints[m+1]-linepoints[m])
        
        dist=math.sqrt(vtk.vtkLine().DistanceToLine(i,linepoints[m],linepoints[m+1]))
        
        pdist = p.distance(line3)
        
        print(d,dist,pdist)

问题在于,使用叉积计算的是到由点 linepoints[m]linepoints[m+1] 定义的线段所跨越的直线的正交距离。另一方面,Shapely 计算到线段的距离,即,如果正交投影落在线段的 "outside" 上,它 return 到正交投影或边界点之一的距离。

要获得一致的结果,您可以自己计算正交投影,然后调用 Shapely 距离方法:

import numpy as np
from shapely.geometry import Point, LineString


A = np.array([1,0])
B = np.array([3,0])
C = np.array([0,1])


l = LineString([A, B])
p = Point(C)


d = np.linalg.norm(np.cross(B - A, C - A))/np.linalg.norm(B - A)

n = B - A
v = C - A

z = A + n*(np.dot(v, n)/np.dot(n, n))

print(l.distance(p), d, Point(z).distance(p))
#1.4142135623730951 1.0 1.0

但是请注意,Shapely 有效地忽略了 z 坐标。例如:

import numpy as np
from shapely.geometry import Point, LineString

A = np.array([1,0,1])
B = np.array([0,0,0])

print(Point([1,0,1]).distance(Point([0,0,0])))

return 因为距离只有 1.

编辑: 根据您的评论,这里将是一个计算到段的距离(对于任意数量的维度)的版本:

from shapely.geometry import LineString, Point
import numpy as np
import itertools

import math

from numpy.linalg import norm

x1=np.arange(4,21)
y1=np.arange(4,21)
z1=np.arange(-7,6)

linepoints = np.array([[1,10,0],[10,10,0],[15,15,0]])

def dist(A, B, C):
    """Calculate the distance of point C to line segment spanned by points A, B.

    """

    a = np.asarray(A)
    b = np.asarray(B)
    c = np.asarray(C)

    #project c onto line spanned by a,b but consider the end points
    #should the projection fall "outside" of the segment    
    n, v = b - a, c - a

    #the projection q of c onto the infinite line defined by points a,b
    #can be parametrized as q = a + t*(b - a). In terms of dot-products,
    #the coefficient t is (c - a).(b - a)/( (b-a).(b-a) ). If we want
    #to restrict the "projected" point to belong to the finite segment
    #connecting points a and b, it's sufficient to "clip" it into
    #interval [0,1] - 0 corresponds to a, 1 corresponds to b.

    t = max(0, min(np.dot(v, n)/np.dot(n, n), 1))
    return np.linalg.norm(c - (a + t*n)) #or np.linalg.norm(v - t*n)


for coords in itertools.product(x1,y1,z1):
    for m in range(len(linepoints)-1):

        line3 = LineString([linepoints[m],linepoints[m+1]])
        d = dist(linepoints[m], linepoints[m+1], coords)
        print(coords, d)