From 0507e7a5b8dd558e12c906dd398a6e0bd7f76911 Mon Sep 17 00:00:00 2001 From: Campbell Barton Date: Sun, 10 Jul 2011 17:27:57 +0000 Subject: added some more methods of evaluating the curve. --- modules/curve_utils.py | 356 ++++++++++++++++++++++++++++++++++++++----------- 1 file changed, 281 insertions(+), 75 deletions(-) diff --git a/modules/curve_utils.py b/modules/curve_utils.py index fba6489e..529ddd0c 100644 --- a/modules/curve_utils.py +++ b/modules/curve_utils.py @@ -92,10 +92,7 @@ def treat_points(points, def solve_curvature(p1, p2, n1, n2, fac, fallback): """ Add a nice circular curvature on """ - from mathutils import Vector - from mathutils.geometry import (barycentric_transform, - intersect_line_line, - intersect_point_line, + from mathutils.geometry import (intersect_line_line, ) p1_a = p1 + n1 @@ -126,7 +123,7 @@ def solve_curvature(p1, p2, n1, n2, fac, fallback): def points_to_bezier(points_orig, double_limit=0.0001, kink_tolerance=0.25, - bezier_tolerance=0.02, # error distance, scale dependant + bezier_tolerance=0.05, # error distance, scale dependant subdiv=8, angle_span=0.95, # 1.0 tries to evaluate splines of 180d ): @@ -360,28 +357,75 @@ def points_to_bezier(points_orig, self.points[0].is_joint, self.points[-1].is_joint = joint self.calc_all() + # raise Exception("END") - def bezier_solve_(self): - """ Calculate bezier handles, - assume the splines have been broken up. - - http://polymathprogrammer.com/ + def intersect_line(self, l1, l2, reverse=False): + """ Spectial kind of intersection, works in 3d on the plane + defimed by the points normal and the line. """ - p1 = self.points[0] - p2 = self.points[-1] + from mathutils.geometry import (intersect_point_line, + ) + + if reverse: + p_first = self.points[-1] + no = -self.points[-1].no + point_iter = reversed(self.points[:-1]) + else: + p_first = self.points[0] + no = self.points[0].no + point_iter = self.points[1:] - line_ix_p1 = self.points[len(self.points) // 3] - line_ix_p2 = self.points[int((len(self.points) / 3) * 2)] + # calculate the line right angles to the line + bi_no = (no - no.project(l2 - l1)).normalized() - u = 1 / 3 - v = 2 / 3 + bi_l1 = p_first.co + bi_l2 = p_first.co + bi_no - p0x, p0y, p0z = p1.co - p1x, p1y, p1z = line_ix_p1.co - p2x, p2y, p2z = line_ix_p2.co - p3x, p3y, p3z = p2.co - + for p_apex in point_iter: + ix, fac = intersect_point_line(p_apex.co, bi_l1, bi_l2) + + if fac < 0.0001: + + if reverse: + p_apex_other = p_apex.next + else: + p_apex_other = p_apex.prev + + # find the exact point on the line between the apex and + # the middle + p_test_1 = intersect_point_line(p_apex.co, + l1, + l2)[0] + p_test_2 = intersect_point_line(p_apex_other.co, + l1, + l2)[0] + + w1 = (p_test_1 - p_apex.co).length + w2 = (p_test_2 - p_apex_other.co).length + + #assert(w1 + w2 != 0) + try: + fac = w1 / (w1 + w2) + except ZeroDivisionError: + fac = 0.5 + assert(fac >= 0.0 and fac <= 1.0) + + p_apex_co = p_apex.co.lerp(p_apex_other.co, fac) + p_apex_no = p_apex.no.lerp(p_apex_other.no, fac) + p_apex_no.normalize() + + # visualize_line(p_mid.to_3d(), corner.to_3d()) + # visualize_line(p_apex.co.to_3d(), p_apex_co.to_3d()) + + return p_apex_co, p_apex_no, p_apex + + # intersection not found + return None, None, None + + + @staticmethod + def bez_solve(p0, p1, p2, p3, u, v): ui = 1.0 - u vi = 1.0 - v u_p3 = u * u * u @@ -389,14 +433,6 @@ def points_to_bezier(points_orig, ui_p3 = ui * ui * ui vi_p3 = vi * vi * vi - - # --- snip - - q1 = Vector() - q2 = Vector() - - pos = Vector(), Vector(), Vector(), Vector() - a = 3.0 * ui * ui * u b = 3.0 * ui * u * u c = 3.0 * vi * vi * v @@ -407,26 +443,117 @@ def points_to_bezier(points_orig, assert(0) return 0 - q1.x = p1x - (ui_p3 * p0x + u_p3 * p3x) - q1.y = p1y - (ui_p3 * p0y + u_p3 * p3y) - q1.z = p1z - (ui_p3 * p0z + u_p3 * p3z) + q1 = p1 - (ui_p3 * p0 + u_p3 * p3) + q2 = p2 - (vi_p3 * p0 + v_p3 * p3) - q2.x = p2x - (vi_p3 * p0x + v_p3 * p3x) - q2.y = p2y - (vi_p3 * p0y + v_p3 * p3y) - q2.z = p2z - (vi_p3 * p0z + v_p3 * p3z) + return ((d * q1 - b * q2) / det, + (-c * q1 + a * q2) / det + ) - pos[1].x = (d * q1.x - b * q2.x) / det - pos[1].y = (d * q1.y - b * q2.y) / det - pos[1].z = (d * q1.z - b * q2.z) / det + def bezier_solve__math1(self): + """ Calculate bezier handles, + assume the splines have been broken up. - pos[2].x = ((-c) * q1.x + a * q2.x) / det - pos[2].y = ((-c) * q1.y + a * q2.y) / det - pos[2].z = ((-c) * q1.z + a * q2.z) / det + http://polymathprogrammer.com/ + """ - self.handle_left = pos[1] - self.handle_right = pos[2] + def get(f, min=0.0, max=1.0): + f = (f * (max - min) + min) + return self.points[int((len(self.points) - 1) * f)].co + + + p1 = get(0.0) + p2 = get(1.0) + i1 = get(1/3) + i2 = get(2/3) + + pos = __class__.bez_solve(p1, i1, i2, p2, 1.0 / 3.0, 2.0 / 3.0) + self.handle_left = self.points[0].co + (pos[0] - self.points[0].co) + self.handle_right = self.points[-1].co + (pos[1] - self.points[-1].co) + + def bezier_solve__math2(self): + + def get(f, min=0.0, max=1.0): + f = (f * (max - min) + min) + return self.points[int((len(self.points) - 1) * f)].co + + p1 = get(0.0, min=0.0, max=0.5) + p2 = get(1.0, min=0.0, max=0.5) + i1 = get(1/3, min=0.0, max=0.5) + i2 = get(2/3, min=0.0, max=0.5) + + pos_a = __class__.bez_solve(p1, i1, i2, p2, 1.0 / 3.0, 2.0 / 3.0) + + p1 = get(0.0, min=0.5, max=1.0) + p2 = get(1.0, min=0.5, max=1.0) + i1 = get(1/3, min=0.5, max=1.0) + i2 = get(2/3, min=0.5, max=1.0) + + pos_b = __class__.bez_solve(p1, i1, i2, p2, 1.0 / 3.0, 2.0 / 3.0) + + self.handle_left = self.points[0].co + (pos_a[0] - self.points[0].co) * 2 + self.handle_right = self.points[-1].co + (pos_b[1] - self.points[-1].co) * 2 + + def bezier_solve__inkscape(self): + + # static void + # estimate_bi(Point bezier[4], unsigned const ei, + # Point const data[], double const u[], unsigned const len) + def estimate_bi(bezier, ei, data, u): + + def B0(u): return ( ( 1.0 - u ) * ( 1.0 - u ) * ( 1.0 - u ) ) + def B1(u): return ( 3 * u * ( 1.0 - u ) * ( 1.0 - u ) ) + def B2(u): return ( 3 * u * u * ( 1.0 - u ) ) + def B3(u): return ( u * u * u ) + + # assert( not (1 <= ei and ei <= 2)) + oi = 3 - ei + num = [0.0, 0.0, 0.0] + den = 0.0 + + for i in range(len(data)): + ui = u[i]; + b = [ + B0(ui), + B1(ui), + B2(ui), + B3(ui) + ] + + for d in range(3): + num[d] += (b[ei] * (b[0] * bezier[0][d] + + b[oi] * bezier[oi][d] + + b[3] * bezier[3][d] + + - data[i][d])) + + den -= b[ei] * b[ei] + + if den: + for d in range(3): + bezier[ei][d] = num[d] / den + else: + bezier[ei] = (oi * bezier[0] + ei * bezier[3]) / 3.0 + bezier = [ + self.points[0].co, + self.points[0].co.lerp(self.points[-1].co, 1/3), + self.points[0].co.lerp(self.points[-1].co, 2/3), + self.points[-1].co, + ] + data = [p.co for p in self.points] + u = [i / len(self.points) for i in range(len(self.points))] + estimate_bi(bezier, 1, data, u) + estimate_bi(bezier, 2, data, u) + estimate_bi(bezier, 1, data, u) + estimate_bi(bezier, 2, data, u) + estimate_bi(bezier, 1, data, u) + estimate_bi(bezier, 2, data, u) + estimate_bi(bezier, 1, data, u) + estimate_bi(bezier, 2, data, u) + + self.handle_left = bezier[1] + self.handle_right = bezier[2] - def bezier_solve(self): + def bezier_solve_ideasman42(self): from mathutils.geometry import (intersect_point_line, intersect_line_line, ) @@ -449,41 +576,92 @@ def points_to_bezier(points_orig, # visualize_line(p1.co, l1_co) # visualize_line(p2.co, l2_co) - # picking 1/2 and 2/3'rds works best - line_ix_p1 = self.points[int(len(self.points) * (1.0 / 3.0))] - line_ix_p1_co, line_ix_p1_no = line_ix_p1.co, line_ix_p1.no - line_ix_p2 = self.points[int(len(self.points) * (2.0 / 3.0))] - line_ix_p2_co, line_ix_p2_no = line_ix_p2.co, line_ix_p2.no - - # used to seek for the upper most point but this gives mostly - # as good results - p1_apex_co = self.points[int(len(self.points) * (1.0 / 3.0) * 0.75)].co - p2_apex_co = self.points[int(len(self.points) * (1.0 - (1.0 / 3.0) * 0.75))].co + line_ix_p1_co, line_ix_p1_no, line_ix_p1 = \ + self.intersect_line(p1.co, + l1_co, + ) + line_ix_p2_co, line_ix_p2_no, line_ix_p2 = \ + self.intersect_line(p2.co, + l2_co, + reverse=True, + ) + if line_ix_p1_co is None: + line_ix_p1_co, line_ix_p1_no, line_ix_p1 = \ + p1.next.co, p1.next.no, p1.next + if line_ix_p2_co is None: + line_ix_p2_co, line_ix_p2_no, line_ix_p2 = \ + p2.prev.co, p2.prev.no, p2.prev + + # vis_circle_object(line_ix_p1_co) + # vis_circle_object(line_ix_p2_co) + + l1_max = 0.0 + p1_apex_co = None + p = self.points[1] + while p and (not p.is_joint) and p != line_ix_p1: + ix = intersect_point_line(p.co, p1.co, l1_co)[0] + length = (ix - p.co).length + if length > l1_max: + l1_max = length + p1_apex_co = p.co + p = p.next + + l2_max = 0.0 + p2_apex_co = None + p = self.points[-2] + while p and (not p.is_joint) and p != line_ix_p2: + ix = intersect_point_line(p.co, p2.co, l2_co)[0] + length = (ix - p.co).length + if length > l2_max: + l2_max = length + p2_apex_co = p.co + p = p.prev + + if p1_apex_co is None: + p1_apex_co = p1.next.co + if p2_apex_co is None: + p2_apex_co = p2.prev.co l1_tan = (p1.no - p1.no.project(l1_no)).normalized() l2_tan = -(p2.no - p2.no.project(l2_no)).normalized() + # values are good! + visualize_line(p1.co, p1.co + l1_tan) + visualize_line(p2.co, p2.co + l2_tan) + + visualize_line(p1.co, p1.co + l1_no) + visualize_line(p2.co, p2.co + l2_no) + # calculate bias based on the position of the other point allong # the tangent. # first need to reflect the second normal for angle comparison # first fist need the reflection normal + + # angle between - 0 - 1 + from math import pi no_ref = p_vec.cross(p2.no).cross(p_vec).normalized() l2_no_ref = p2.no.reflect(no_ref).normalized() + no_angle = p1.no.angle(l2_no_ref) / pi del no_ref - from math import pi # This could be tweaked but seems to work well - fac_fac = (p1.no.angle(l2_no_ref) / pi) - fac_1 = p1.no.angle(line_ix_p1_co - p1.co) / pi - fac_2 = (-p2.no).angle(line_ix_p2_co - p2.co) / pi + # fac_fac = 1.0 + + print("angle", "%.6f" % no_angle) - # fac_1 = fac_2 = 0.0 - print(fac_1, fac_2) - # why * 3 ? - it just gives best results - h1_fac = ((p1.co - p1_apex_co).length / 0.75) * (1.0 + fac_1 * fac_fac * 3.0) - h2_fac = ((p2.co - p2_apex_co).length / 0.75) * (1.0 + fac_2 * fac_fac * 3.0) + fac_1 = intersect_point_line(p2_apex_co, + p1.co, + p1.co + l1_tan * (p1.co - p1_apex_co).length, + )[1] * (1.0 + no_angle) + fac_2 = intersect_point_line(p1_apex_co, + p2.co, + p2.co + l2_tan * (p2.co - p2_apex_co).length, + )[1] * (1.0 + no_angle) + + h1_fac = abs(fac_1) + h2_fac = abs(fac_2) h1 = p1.co + (p1.no * h1_fac) h2 = p2.co - (p2.no * h2_fac) @@ -491,14 +669,17 @@ def points_to_bezier(points_orig, self.handle_left = h1 self.handle_right = h2 - """ + ''' visualize_line(p1.co, p1_apex_co) visualize_line(p1_apex_co, p2_apex_co) visualize_line(p2.co, p2_apex_co) visualize_line(p1.co, p2.co) - """ + ''' + + def bezier_solve(self): + return self.bezier_solve__inkscape() - def bezier_error(self, error_max=-1.0, test_count=16): + def bezier_error(self, error_max=-1.0, test_count=8): from mathutils.geometry import interpolate_bezier test_points = interpolate_bezier(self.points[0].co, @@ -515,7 +696,6 @@ def points_to_bezier(points_orig, # this is a rough method measuring the error but should be ok # TODO. dont test against every single point. for co in test_points: - co = co # initial values co_best = self.points[0].co @@ -752,11 +932,25 @@ def points_to_bezier(points_orig, for s in curve.splines: s.bezier_solve() ''' + + ''' + def angle_point(s): + a = 0.0 + a_best = len(s.points) // 2 + i = 1 + for p in s.points[2:-2]: + if p.angle > a: + a = p.angle + a_best = i + i += 1 + return a_best + ''' + # or recursively subdivide... curve.split_func_spline(lambda s: - len(s.points) // 2 + len(s.points) // 2 # angle_point(s) if ((s.bezier_solve(), - s.bezier_error(bezier_tolerance))[1] + s.bezier_error(bezier_tolerance))[1] and (len(s.points))) else -1, recursive=True, @@ -774,14 +968,26 @@ def points_to_bezier(points_orig, if __name__ == "__main__": - bpy.ops.wm.open_mainfile(filepath="/root/curve_test2.blend") + if 0: + bpy.ops.wm.open_mainfile(filepath="/root/curve_test3.blend") + + for c in "Curve Curve.001 Curve.002 Curve.003 Curve.004 Curve.005".split(): + print("---", c) + ob = bpy.data.objects[c] + points = [p.co.xyz for s in ob.data.splines for p in s.points] + + print("points_to_bezier 1") + points_to_bezier(points) + print("points_to_bezier 2") + else: + bpy.ops.wm.open_mainfile(filepath="/root/curve_test2.blend") - ob = bpy.data.objects["Curve"] - points = [p.co.xyz for s in ob.data.splines for p in s.points] + ob = bpy.data.objects['Curve'] + points = [p.co.xyz for s in ob.data.splines for p in s.points] - print("points_to_bezier 1") - points_to_bezier(points) - print("points_to_bezier 2") + print("points_to_bezier 1") + points_to_bezier(points) + print("points_to_bezier 2") bpy.ops.wm.save_as_mainfile(filepath="/root/curve_test_edit.blend", copy=True) -- cgit v1.2.3