2017-05-23 88 views
0

我有以下問題。我有一個裝滿座標和三點的盒子,它們組成了一條線。現在我想計算所有框座標到該線的最短距離。我有三種方法來做到這一點,vtk和numpy版本總是有相同的結果,但沒有勻稱的距離方法。但我需要這個造型的版本,因爲我想測量從一個點到hwole線的最近距離,而不是單獨的線段。以下是一個示例代碼。所以,問題是「pdist」:Python:最近指向一行

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) 

回答

1

的問題是,與跨產品您正在計算通過由點linepoints[m]linepoints[m+1]確定的段橫跨的線正交距離。另一方面,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 

但是請注意,身材勻稱有效地忽略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) 
+0

這是個問題!那麼最好的方法是計算所有3D點與多段線的最近距離?因爲在某些情況下,點與線段的無限延伸距離最近,但我需要距離線段/線上的點最近的距離,當然也需要與z方向上的點最接近的距離。我認爲這是你用正交投影描述的東西嗎? – Varlor

+1

@Varlor我已經用我認爲是一個通用的方法更新了答案。在'dist'函數中,將參數't'限制在區間'[0,1]'中,確保它計算到有限區段的距離。沒有這個限制,它會返回到無窮遠線的距離... – ewcz

+0

啊,很好,它的工作原理。你能解釋一下數學上究竟有多少行: n,v = b - a,c - a 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) – Varlor