Newer
Older
return (3 * ((y3 - y0) * (x1 + x2) - (x3 - x0) * (y1 + y2)
+ y1 * (x0 - x2) - x1 * (y0 - y2)
+ y3 * (x2 + x0 / 3) - x3 * (y2 + y0 / 3)) / 20)
1005
1006
1007
1008
1009
1010
1011
1012
1013
1014
1015
1016
1017
1018
1019
1020
1021
1022
1023
1024
1025
1026
1027
1028
1029
##############################################
def closest_point(self, point):
# n = P3 - 3*P2 + 3*P1 - P0
# r = 3*(P2 - 2*P1 + P0
# s = 3*(P1 - P0)
# v = P0
# Q(t) = n*t**3 + r*t**2 + s*t + v
# Q'(t) = 3*n*t**2 + 2*r*t + s
# n, r, s, v = symbols('n r s v')
# Q = n*t**3 + r*t**2 + s*t + v
# Qp = simplify(Q.diff(t))
# collect(expand((P*Qp - Q*Qp)), t)
# -3*n**2 * t**5
# -5*n*r * t**4
# -2*(2*n*s + r**2) * t**3
# 3*(P*n - n*v - r*s) * t**2
# (2*P*r - 2*r*v - s**2) * t
# P*s - s*v
n = self._p3 - self._p2*3 + self._p1*3 - self._p0
r = (self._p2 - self._p1*2 + self._p0)*3
s = (self._p1 - self._p0)*3
v = self._p0
roots = fifth_root(
-3 * n.magnitude_square,
-5 * n.dot(r),
-2 * (2*n.dot(s) + r.magnitude_square),
3 * (point.dot(n) - n.dot(v) - r.dot(s)),
2*point.dot(r) - 2*r.dot(v) - s.magnitude_square,
point.dot(s) - s.dot(v),
)
# Fixme: to func
t = [root for root in roots if 0 <= root <= 1]
if not t:
return None
elif len(t) > 1:
raise NameError("Found more than one root: {}".format(t))