compas
compas copied to clipboard
compas-dev
Line closest_point potential bug
I think the line.closest_point method gives bad results for finding the closest parameter on a line.
Example:
Given a line: Line(Point(x=-6.0, y=0.0, z=0.0), Point(x=7.5, y=0.0, z=0.0))
And a point: Point(x=-5.0, y=0.0, z=0.0)
The result is 1.0
But it should be: 0.07407407407407407
After looking at OpenNurbs C++ implementation, this is an equivalent in Python:
def closest_point(line, point):
t = 0.0
v = line.end - line.start
vov = v.length**2
if vov > 0.0:
if (point - line.start).length**2 <= (point - line.end).length**2:
t = v.dot(point - line.start) / vov
else:
t = 1.0 + v.dot(point - line.end) / vov
else:
t = 0.0 # degenerate line
return line.point_at(t), t
this is original one in compas:
a = self.start
vector = point - a
direction = self.direction
t = vector.dot(direction)
c = a + direction * t
if return_parameter:
return c, t
return c
i would double check your test because the result seems correct to me...
from compas.geometry import Point, Line
from compas_view2.app import App
line = Line(Point(-6.0, 0.0, 0.0), Point(7.5, 0.0, 0.0))
point = Point(-5, 0, 0)
closest, t = line.closest_point(point, True)
print(t)
viewer = App()
viewer.add(line)
viewer.add(closest)
viewer.run()
t is expressed wrt to the direction vector, which is a unit vector.
therefore 1.0 seems to make sense, since Point(5, 0, 0) is Point(-6, 0, 0) + Vector(1, 0, 0) * 1.0...
the difference with the C++ implementation lies mostly in the division by the squared length of the vector of the line, but that is implicitly included in the Python implementation because of the use of the unit direction vector.
I see, it is not normalized parameter, but is a scale value of line.direction.
I thought that the t parameter can be directly used in point_on_line = line.point_at(t) method.
This was the confusion.
Which means that instead, it has to be used like this: point_on_line = line.start+line.direction*t
In this case it is not a bug, but intended behavior. I guess I was too much used to Rhino implementation.
no, you're right. it should be possible to combine the two. it was not clear to me from your question/statements that this was the problem...
the domain in also defined as 0.0 => start and 1.0 => end so we should indeed update the implementation.
if don't know why the second "if" is necessary though. t = v.dot(point - line.start) / vov and t = 1.0 + v.dot(point - line.end) / vov seem equivalent to me...
start = self.start
vector = point - start
direction = self.vector
t = vector.dot(direction) / self.length**2
closest = start + direction * t
if return_parameter:
return closest, t
return closest
if we don't care about the lookups, can be simplified to
vector = point - self.start
t = vector.dot(self.vector) / self.length**2
closest = self.start + self.vector * t
if return_parameter:
return closest, t
return closest
if don't know why the second "if" is necessary though. t = v.dot(point - line.start) / vov and t = 1.0 + v.dot(point - line.end) / vov seem equivalent to me...
I checked different scenarios, I agree, there is no need to do that because the result is always the same. The dot product projects one vector to the other, using the same starting point, so it is the same. Also, there is never a need to know to which end the point is closer.
Thank you for changing to the normalized result: t = vector.dot(self.vector) / self.length**2
I think it is more clear when using normalization.
