我有以下问题。我有一个充满坐标的盒子和三个点,它们组成一条线。现在,我想计算所有框坐标到该线的最短距离。我有三种方法可以做到这一点,而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/