Curve Fitting: minor change to re-fitting method
authorCampbell Barton <ideasman42@gmail.com>
Sat, 29 Apr 2017 14:01:16 +0000 (00:01 +1000)
committerCampbell Barton <ideasman42@gmail.com>
Sat, 29 Apr 2017 14:01:16 +0000 (00:01 +1000)
Avoid calculating a new split-index when re-fitting.

While checking if a knot can be removed, the index with the highest error
can be used as a candidate to replace the knot
(in the case it can't be removed).

extern/curve_fit_nd/curve_fit_nd.h
extern/curve_fit_nd/intern/curve_fit_cubic.c
extern/curve_fit_nd/intern/curve_fit_cubic_refit.c
source/blender/editors/curve/editcurve.c

index 7232f802e288a3d6eec9d02ecaabc70c4c7ce830..18244799b0f485feaa8546fad778e35942a34404 100644 (file)
@@ -36,7 +36,7 @@
 /* curve_fit_cubic.c */
 
 /**
 /* curve_fit_cubic.c */
 
 /**
- * Takes a flat array of points and evalues that to calculate a bezier spline.
+ * Takes a flat array of points and evaluates that to calculate a bezier spline.
  *
  * \param points, points_len: The array of points to calculate a cubics from.
  * \param dims: The number of dimensions for for each element in \a points.
  *
  * \param points, points_len: The array of points to calculate a cubics from.
  * \param dims: The number of dimensions for for each element in \a points.
@@ -82,7 +82,7 @@ int curve_fit_cubic_to_points_fl(
         unsigned int **r_corners_index_array, unsigned int *r_corners_index_len);
 
 /**
         unsigned int **r_corners_index_array, unsigned int *r_corners_index_len);
 
 /**
- * Takes a flat array of points and evalues that to calculate handle lengths.
+ * Takes a flat array of points and evaluates that to calculate handle lengths.
  *
  * \param points, points_len: The array of points to calculate a cubics from.
  * \param dims: The number of dimensions for for each element in \a points.
  *
  * \param points, points_len: The array of points to calculate a cubics from.
  * \param dims: The number of dimensions for for each element in \a points.
@@ -107,7 +107,8 @@ int curve_fit_cubic_to_points_single_db(
 
         double  r_handle_l[],
         double  r_handle_r[],
 
         double  r_handle_l[],
         double  r_handle_r[],
-        double *r_error_sq);
+        double *r_error_sq,
+        unsigned int *r_error_index);
 
 int curve_fit_cubic_to_points_single_fl(
         const float       *points,
 
 int curve_fit_cubic_to_points_single_fl(
         const float       *points,
@@ -120,7 +121,8 @@ int curve_fit_cubic_to_points_single_fl(
 
         float   r_handle_l[],
         float   r_handle_r[],
 
         float   r_handle_l[],
         float   r_handle_r[],
-        float  *r_error_sq);
+        float  *r_error_sq,
+        unsigned int *r_error_index);
 
 enum {
        CURVE_FIT_CALC_HIGH_QUALIY          = (1 << 0),
 
 enum {
        CURVE_FIT_CALC_HIGH_QUALIY          = (1 << 0),
index 0a32f1e796afc63863cd610316750604986a1d8b..ed855d34b0886b3b66d2c3e292aec438af369a2b 100644 (file)
@@ -554,8 +554,8 @@ static void cubic_from_points_fallback(
        r_cubic->orig_span = (points_offset_len - 1);
 #endif
 
        r_cubic->orig_span = (points_offset_len - 1);
 #endif
 
-       /* p1 = p0 - (tan_l * alpha_l);
-        * p2 = p3 + (tan_r * alpha_r);
+       /* p1 = p0 - (tan_l * alpha);
+        * p2 = p3 + (tan_r * alpha);
         */
        msub_vn_vnvn_fl(p1, p0, tan_l, alpha, dims);
        madd_vn_vnvn_fl(p2, p3, tan_r, alpha, dims);
         */
        msub_vn_vnvn_fl(p1, p0, tan_l, alpha, dims);
        madd_vn_vnvn_fl(p2, p3, tan_r, alpha, dims);
@@ -1436,12 +1436,11 @@ int curve_fit_cubic_to_points_single_db(
 
         double  r_handle_l[],
         double  r_handle_r[],
 
         double  r_handle_l[],
         double  r_handle_r[],
-        double  *r_error_max_sq)
+        double *r_error_max_sq,
+        uint   *r_error_index)
 {
        Cubic *cubic = alloca(cubic_alloc_size(dims));
 
 {
        Cubic *cubic = alloca(cubic_alloc_size(dims));
 
-       uint split_index;
-
        /* in this instance theres no advantage in using length cache,
         * since we're not recursively calculating values. */
 #ifdef USE_LENGTH_CACHE
        /* in this instance theres no advantage in using length cache,
         * since we're not recursively calculating values. */
 #ifdef USE_LENGTH_CACHE
@@ -1462,7 +1461,7 @@ int curve_fit_cubic_to_points_single_db(
 #endif
                tan_l, tan_r, error_threshold, dims,
 
 #endif
                tan_l, tan_r, error_threshold, dims,
 
-               cubic, r_error_max_sq, &split_index);
+               cubic, r_error_max_sq, r_error_index);
 
 #ifdef USE_LENGTH_CACHE
        if (points_length_cache_alloc) {
 
 #ifdef USE_LENGTH_CACHE
        if (points_length_cache_alloc) {
@@ -1487,7 +1486,8 @@ int curve_fit_cubic_to_points_single_fl(
 
         float   r_handle_l[],
         float   r_handle_r[],
 
         float   r_handle_l[],
         float   r_handle_r[],
-        float  *r_error_sq)
+        float  *r_error_sq,
+        uint   *r_error_index)
 {
        const uint points_flat_len = points_len * dims;
        double *points_db = malloc(sizeof(double) * points_flat_len);
 {
        const uint points_flat_len = points_len * dims;
        double *points_db = malloc(sizeof(double) * points_flat_len);
@@ -1521,7 +1521,8 @@ int curve_fit_cubic_to_points_single_fl(
                (double)error_threshold,
                tan_l_db, tan_r_db,
                r_handle_l_db, r_handle_r_db,
                (double)error_threshold,
                tan_l_db, tan_r_db,
                r_handle_l_db, r_handle_r_db,
-               &r_error_sq_db);
+               &r_error_sq_db,
+               r_error_index);
 
        free(points_db);
 
 
        free(points_db);
 
index bf1ab99995f68bf4b089d4a3dfc7b8c7b1a06daf..96ec9a33270237093d905e59949b3cf5f479ade8 100644 (file)
@@ -207,7 +207,7 @@ struct KnotCornerState {
 
 /* Utility functions */
 
 
 /* Utility functions */
 
-#ifdef USE_KNOT_REFIT
+#if defined(USE_KNOT_REFIT) && !defined(USE_KNOT_REFIT_REMOVE)
 /**
  * Find the most distant point between the 2 knots.
  */
 /**
  * Find the most distant point between the 2 knots.
  */
@@ -269,7 +269,7 @@ static uint knot_find_split_point(
 
        return split_point;
 }
 
        return split_point;
 }
-#endif  /* USE_KNOT_REFIT */
+#endif  /* USE_KNOT_REFIT && !USE_KNOT_REFIT_REMOVE */
 
 
 #ifdef USE_CORNER_DETECT
 
 
 #ifdef USE_CORNER_DETECT
@@ -322,7 +322,7 @@ static double knot_remove_error_value(
         const double *points_offset_length_cache,
         const uint dims,
         /* Avoid having to re-calculate again */
         const double *points_offset_length_cache,
         const uint dims,
         /* Avoid having to re-calculate again */
-        double r_handle_factors[2])
+        double r_handle_factors[2], uint *r_error_index)
 {
        double error_sq = FLT_MAX;
 
 {
        double error_sq = FLT_MAX;
 
@@ -338,7 +338,7 @@ static double knot_remove_error_value(
                points_offset, points_offset_len, points_offset_length_cache, dims, 0.0,
                tan_l, tan_r,
                handle_factor_l, handle_factor_r,
                points_offset, points_offset_len, points_offset_length_cache, dims, 0.0,
                tan_l, tan_r,
                handle_factor_l, handle_factor_r,
-               &error_sq);
+               &error_sq, r_error_index);
 
        assert(error_sq != FLT_MAX);
 
 
        assert(error_sq != FLT_MAX);
 
@@ -363,6 +363,7 @@ static double knot_calc_curve_error_value(
                ((knot_r->index + pd->points_len) - knot_l->index)) + 1;
 
        if (points_offset_len != 2) {
                ((knot_r->index + pd->points_len) - knot_l->index)) + 1;
 
        if (points_offset_len != 2) {
+               uint error_index_dummy;
                return knot_remove_error_value(
                        tan_l, tan_r,
                        &pd->points[knot_l->index * dims], points_offset_len,
                return knot_remove_error_value(
                        tan_l, tan_r,
                        &pd->points[knot_l->index * dims], points_offset_len,
@@ -372,7 +373,55 @@ static double knot_calc_curve_error_value(
                        NULL,
 #endif
                        dims,
                        NULL,
 #endif
                        dims,
-                       r_handle_factors);
+                       r_handle_factors, &error_index_dummy);
+       }
+       else {
+               /* No points between, use 1/3 handle length with no error as a fallback. */
+               assert(points_offset_len == 2);
+#ifdef USE_LENGTH_CACHE
+               r_handle_factors[0] = r_handle_factors[1] = pd->points_length_cache[knot_l->index] / 3.0;
+#else
+               r_handle_factors[0] = r_handle_factors[1] = len_vnvn(
+                       &pd->points[(knot_l->index + 0) * dims],
+                       &pd->points[(knot_l->index + 1) * dims], dims) / 3.0;
+#endif
+               return 0.0;
+       }
+}
+
+#ifdef USE_KNOT_REFIT_REMOVE
+
+static double knot_calc_curve_error_value_and_index(
+        const struct PointData *pd,
+        const struct Knot *knot_l, const struct Knot *knot_r,
+        const double *tan_l, const double *tan_r,
+        const uint dims,
+        double r_handle_factors[2],
+        uint *r_error_index)
+{
+       const uint points_offset_len = ((knot_l->index < knot_r->index) ?
+               (knot_r->index - knot_l->index) :
+               ((knot_r->index + pd->points_len) - knot_l->index)) + 1;
+
+       if (points_offset_len != 2) {
+               const double error_sq = knot_remove_error_value(
+                       tan_l, tan_r,
+                       &pd->points[knot_l->index * dims], points_offset_len,
+#ifdef USE_LENGTH_CACHE
+                       &pd->points_length_cache[knot_l->index],
+#else
+                       NULL,
+#endif
+                       dims,
+                       r_handle_factors, r_error_index);
+
+               /* Adjust the offset index to the global index & wrap if needed. */
+               *r_error_index += knot_l->index;
+               if (*r_error_index >= pd->points_len) {
+                       *r_error_index -= pd->points_len;
+               }
+
+               return error_sq;
        }
        else {
                /* No points between, use 1/3 handle length with no error as a fallback. */
        }
        else {
                /* No points between, use 1/3 handle length with no error as a fallback. */
@@ -384,9 +433,11 @@ static double knot_calc_curve_error_value(
                        &pd->points[(knot_l->index + 0) * dims],
                        &pd->points[(knot_l->index + 1) * dims], dims) / 3.0;
 #endif
                        &pd->points[(knot_l->index + 0) * dims],
                        &pd->points[(knot_l->index + 1) * dims], dims) / 3.0;
 #endif
+               *r_error_index = 0;
                return 0.0;
        }
 }
                return 0.0;
        }
 }
+#endif  /* USE_KNOT_REFIT_REMOVE */
 
 struct KnotRemove_Params {
        Heap *heap;
 
 struct KnotRemove_Params {
        Heap *heap;
@@ -556,15 +607,18 @@ static void knot_refit_error_recalculate(
        assert(k->can_remove);
 
 #ifdef USE_KNOT_REFIT_REMOVE
        assert(k->can_remove);
 
 #ifdef USE_KNOT_REFIT_REMOVE
+       (void)knots_len;
+
+       uint refit_index = SPLIT_POINT_INVALID;
        {
                double handles[2];
 
                /* First check if we can remove, this allows to refit and remove as we go. */
        {
                double handles[2];
 
                /* First check if we can remove, this allows to refit and remove as we go. */
-               const double cost_sq = knot_calc_curve_error_value(
+               const double cost_sq = knot_calc_curve_error_value_and_index(
                        p->pd, k->prev, k->next,
                        k->prev->tan[1], k->next->tan[0],
                        dims,
                        p->pd, k->prev, k->next,
                        k->prev->tan[1], k->next->tan[0],
                        dims,
-                       handles);
+                       handles, &refit_index);
 
                if (cost_sq < error_sq_max) {
                        struct KnotRefitState *r;
 
                if (cost_sq < error_sq_max) {
                        struct KnotRefitState *r;
@@ -598,13 +652,14 @@ static void knot_refit_error_recalculate(
        }
 #else
        (void)error_sq_max;
        }
 #else
        (void)error_sq_max;
-#endif  /* USE_KNOT_REFIT_REMOVE */
 
        const uint refit_index = knot_find_split_point(
                 p->pd, k->prev, k->next,
                 knots_len,
                 dims);
 
 
        const uint refit_index = knot_find_split_point(
                 p->pd, k->prev, k->next,
                 knots_len,
                 dims);
 
+#endif  /* USE_KNOT_REFIT_REMOVE */
+
        if ((refit_index == SPLIT_POINT_INVALID) ||
            (refit_index == k->index))
        {
        if ((refit_index == SPLIT_POINT_INVALID) ||
            (refit_index == k->index))
        {
index 648fe93003096920b8305503ccda04229accb749..a69264cd0126a93738d85ff46d2ebd1cae4142b2 100644 (file)
@@ -5850,6 +5850,7 @@ static int curve_dissolve_exec(bContext *C, wmOperator *UNUSED(op))
                                        BLI_assert(points_stride + dims == points + (points_len * dims));
 
                                        float tan_l[3], tan_r[3], error_sq_dummy;
                                        BLI_assert(points_stride + dims == points + (points_len * dims));
 
                                        float tan_l[3], tan_r[3], error_sq_dummy;
+                                       unsigned int error_index_dummy;
 
                                        sub_v3_v3v3(tan_l, bezt_prev->vec[1], bezt_prev->vec[2]);
                                        normalize_v3(tan_l);
 
                                        sub_v3_v3v3(tan_l, bezt_prev->vec[1], bezt_prev->vec[2]);
                                        normalize_v3(tan_l);
@@ -5860,7 +5861,7 @@ static int curve_dissolve_exec(bContext *C, wmOperator *UNUSED(op))
                                                points, points_len, NULL, dims, FLT_EPSILON,
                                                tan_l, tan_r,
                                                bezt_prev->vec[2], bezt_next->vec[0],
                                                points, points_len, NULL, dims, FLT_EPSILON,
                                                tan_l, tan_r,
                                                bezt_prev->vec[2], bezt_next->vec[0],
-                                               &error_sq_dummy);
+                                               &error_sq_dummy, &error_index_dummy);
 
                                        if (!ELEM(bezt_prev->h2, HD_FREE, HD_ALIGN)) {
                                                bezt_prev->h2 = (bezt_prev->h2 == HD_VECT) ? HD_FREE : HD_ALIGN;
 
                                        if (!ELEM(bezt_prev->h2, HD_FREE, HD_ALIGN)) {
                                                bezt_prev->h2 = (bezt_prev->h2 == HD_VECT) ? HD_FREE : HD_ALIGN;