我有以下问题。我有一个充满坐标的盒子和三个点,它们组成一条线。现在,我想计算所有框坐标到该线的最短距离。我有三种方法可以做到这一点,而vtk和numpy版本始终具有相同的结果,但是Shape的distance方法却没有。但是我需要匀称的版本,因为我想测量一点到最短距离而不是分开的线段的最近距离。这是到目前为止的示例代码。所以问题出在“ 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计算到分段的距离,即,如果正交投影落在分段的“外部”,它会返回该距离到正交投影或边界点之一。

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

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])))


仅作为距离返回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)

关于python - Python:最接近点的一行,我们在Stack Overflow上找到一个类似的问题:https://stackoverflow.com/questions/44129007/

10-12 17:02