Tried, but failed, to write a better smooth bezier function

manimlib/mobject/types/vectorized_mobject.py

@@ -9,7 +9,8 @@ from manimlib.constants import *
 from manimlib.mobject.mobject import Mobject
 from manimlib.mobject.mobject import Point
 from manimlib.utils.bezier import bezier
-from manimlib.utils.bezier import get_smooth_handle_points
+from manimlib.utils.bezier import get_smooth_cubic_bezier_handle_points
+from manimlib.utils.bezier import get_smooth_quadratic_bezier_handle_points
 from manimlib.utils.bezier import get_quadratic_approximation_of_cubic
 from manimlib.utils.bezier import interpolate
 from manimlib.utils.bezier import integer_interpolate
@@ -477,7 +478,7 @@ class VMobject(Mobject):
             for subpath in subpaths:
                 anchors = np.vstack([subpath[::nppc], subpath[-1:]])
                 if mode == "smooth":
-                    h1, h2 = get_smooth_handle_points(anchors)
+                    h1, h2 = get_smooth_cubic_bezier_handle_points(anchors)
                     new_subpath = get_quadratic_approximation_of_cubic(
                         anchors[:-1], h1, h2, anchors[1:]
                     )
@@ -523,10 +524,7 @@ class VMobject(Mobject):
 
     #
     def consider_points_equals(self, p0, p1):
-        return np.allclose(
+        return get_norm(p1 - p0) < self.tolerance_for_point_equality
-            p0, p1,
-            atol=self.tolerance_for_point_equality
-        )
 
     # Information about the curve
     def get_bezier_tuples_from_points(self, points):
@@ -543,10 +541,14 @@ class VMobject(Mobject):
 
     def get_subpaths_from_points(self, points):
         nppc = self.n_points_per_curve
-        split_indices = filter(
+        diffs = points[nppc - 1:-1:nppc] - points[nppc::nppc]
-            lambda n: not self.consider_points_equals(points[n - 1], points[n]),
+        splits = (diffs * diffs).sum(1) > self.tolerance_for_point_equality
-            range(nppc, len(points), nppc)
+        split_indices = np.arange(nppc, len(points), nppc, dtype=int)[splits]
-        )
+
+        # split_indices = filter(
+        #     lambda n: not self.consider_points_equals(points[n - 1], points[n]),
+        #     range(nppc, len(points), nppc)
+        # )
         split_indices = [0, *split_indices, len(points)]
         return [
             points[i1:i2]

manimlib/utils/bezier.py

@@ -1,8 +1,11 @@
 from scipy import linalg
 import numpy as np
 
+from manimlib.constants import PI
 from manimlib.utils.simple_functions import choose
+from manimlib.utils.space_ops import rotate_vector
 from manimlib.utils.space_ops import find_intersection
+from manimlib.utils.space_ops import cross
 from manimlib.utils.space_ops import cross2d
 
 CLOSED_THRESHOLD = 0.001
@@ -83,9 +86,36 @@ def match_interpolate(new_start, new_end, old_start, old_end, old_value):
 
 
 # Figuring out which bezier curves most smoothly connect a sequence of points
+def get_smooth_quadratic_bezier_handle_points(points):
+    # Alas, this function does not actually work very well.
+    # 
+    # For each point P_i, where 1 <= i <= n, draw a line through
+    # P_i parallel to the line through (P_{i-1}, P_{i+1}).  The
+    # intersection of these lines form most of the handles.
+    #
+    # What remains are those near the end points.  For that, we want
+    # the handle between P_0 and P_1 to be closest to (P_0 + P_1) / 2,
+    # which will minimize the second derivative of that curve.  Likewise
+    # for the last handle point.
+    t01 = points[1] - points[0]
+    t12 = points[2] - points[1]
+    tm2 = points[-2] - points[-3]
+    tm1 = points[-1] - points[-2]
+    tangents = np.vstack([
+        rotate_vector(t01, PI / 2, cross(t01, t12)),
+        points[2:] - points[:-2],
+        rotate_vector(tm1, PI / 2, cross(tm1, tm2))
+    ])
+    alt_points = np.array(points)
+    alt_points[0] = points[:2].mean(0)
+    alt_points[-1] = points[-2:].mean(0)
+    return find_intersection(
+        alt_points[:-1], tangents[:-1],
+        alt_points[1:], tangents[1:],
+    )
 
 
-def get_smooth_handle_points(points):
+def get_smooth_cubic_bezier_handle_points(points):
     points = np.array(points)
     num_handles = len(points) - 1
     dim = points.shape[1]
@@ -234,7 +264,8 @@ def get_quadratic_approximation_of_cubic(a0, h0, h1, a1):
 
 
 def get_smooth_quadratic_bezier_path_through(points):
-    h0, h1 = get_smooth_handle_points(points)
+    # TODO
+    h0, h1 = get_smooth_cubic_bezier_handle_points(points)
     a0 = points[:-1]
     a1 = points[1:]
     return get_quadratic_approximation_of_cubic(a0, h0, h1, a1)