PATCH: Use floats instead of doubles in some spline calculations for improved performance

From 1c77afaa79fc27207c7a5741cb4ceb7c2e8e7db6 Mon Sep 17 00:00:00 2001
From: Marko Kuukkanen <[email protected]>
Date: Tue, 19 Feb 2013 12:54:09 +0200
Subject: [PATCH] Optimized spline calculations

---
 src/cairo-fixed-private.h |    7 ++++
 src/cairo-spline.c        |   85 +++++++++++++++++++++++++++++++++++++-------
 2 files changed, 78 insertions(+), 14 deletions(-)

diff --git a/src/cairo-fixed-private.h b/src/cairo-fixed-private.h
index b6cc6be..d4a2645 100644
--- a/src/cairo-fixed-private.h
+++ b/src/cairo-fixed-private.h
@@ -50,6 +50,7 @@
 
 #define CAIRO_FIXED_ONE        ((cairo_fixed_t)(1 << CAIRO_FIXED_FRAC_BITS))
 #define CAIRO_FIXED_ONE_DOUBLE ((double)(1 << CAIRO_FIXED_FRAC_BITS))
+#define CAIRO_FIXED_ONE_FLOAT  ((float)(1 << CAIRO_FIXED_FRAC_BITS))
 #define CAIRO_FIXED_EPSILON    ((cairo_fixed_t)(1))
 
 #define CAIRO_FIXED_ERROR_DOUBLE (1. / (2 * CAIRO_FIXED_ONE_DOUBLE))
@@ -153,6 +154,12 @@ _cairo_fixed_to_double (cairo_fixed_t f)
     return ((double) f) / CAIRO_FIXED_ONE_DOUBLE;
 }
 
+static inline float
+_cairo_fixed_to_float (cairo_fixed_t f)
+{
+    return ((float) f) / CAIRO_FIXED_ONE_FLOAT;
+}
+
 static inline int
 _cairo_fixed_is_integer (cairo_fixed_t f)
 {
diff --git a/src/cairo-spline.c b/src/cairo-spline.c
index 44634fa..86fb39e 100644
--- a/src/cairo-spline.c
+++ b/src/cairo-spline.c
@@ -169,23 +169,53 @@ _de_casteljau (cairo_spline_knots_t *s1, cairo_spline_knots_t *s2)
     s1->d = final;
 }
 
+/* Should only be used, if we know that there will be no overflow */
+static void
+_lerp_half_fast (const cairo_point_t *a, const cairo_point_t *b, cairo_point_t *result)
+{
+    result->x = (a->x + b->x) >> 1;
+    result->y = (a->y + b->y) >> 1;
+}
+
+static void
+_de_casteljau_fast (cairo_spline_knots_t *s1, cairo_spline_knots_t *s2)
+{
+    cairo_point_t ab, bc, cd;
+    cairo_point_t abbc, bccd;
+    cairo_point_t final;
+
+    _lerp_half_fast (&s1->a, &s1->b, &ab);
+    _lerp_half_fast (&s1->b, &s1->c, &bc);
+    _lerp_half_fast (&s1->c, &s1->d, &cd);
+    _lerp_half_fast (&ab, &bc, &abbc);
+    _lerp_half_fast (&bc, &cd, &bccd);
+    _lerp_half_fast (&abbc, &bccd, &final);
+
+    s2->a = final;
+    s2->b = bccd;
+    s2->c = cd;
+    s2->d = s1->d;
+
+    s1->b = ab;
+    s1->c = abbc;
+    s1->d = final;
+}
+
 /* Return an upper bound on the error (squared) that could result from
  * approximating a spline as a line segment connecting the two endpoints. */
-static double
+static float
 _cairo_spline_error_squared (const cairo_spline_knots_t *knots)
 {
-    double bdx, bdy, berr;
-    double cdx, cdy, cerr;
-
+    float bdx, bdy, berr;
+    float cdx, cdy, cerr;
     /* We are going to compute the distance (squared) between each of the the b
      * and c control points and the segment a-b. The maximum of these two
      * distances will be our approximation error. */
+    bdx = _cairo_fixed_to_float (knots->b.x - knots->a.x);
+    bdy = _cairo_fixed_to_float (knots->b.y - knots->a.y);
 
-    bdx = _cairo_fixed_to_double (knots->b.x - knots->a.x);
-    bdy = _cairo_fixed_to_double (knots->b.y - knots->a.y);
-
-    cdx = _cairo_fixed_to_double (knots->c.x - knots->a.x);
-    cdy = _cairo_fixed_to_double (knots->c.y - knots->a.y);
+    cdx = _cairo_fixed_to_float (knots->c.x - knots->a.x);
+    cdy = _cairo_fixed_to_float (knots->c.y - knots->a.y);
 
     if (knots->a.x != knots->d.x || knots->a.y != knots->d.y) {
 	/* Intersection point (px):
@@ -195,10 +225,10 @@ _cairo_spline_error_squared (const cairo_spline_knots_t *knots)
 	 *     u = ((p - p1) ∙ (p2 - p1)) / ∥p2 - p1∥²;
 	 */
 
-	double dx, dy, u, v;
+	float dx, dy, u, v;
 
-	dx = _cairo_fixed_to_double (knots->d.x - knots->a.x);
-	dy = _cairo_fixed_to_double (knots->d.y - knots->a.y);
+	dx = _cairo_fixed_to_float (knots->d.x - knots->a.x);
+	dy = _cairo_fixed_to_float (knots->d.y - knots->a.y);
 	 v = dx * dx + dy * dy;
 
 	u = bdx * dx + bdy * dy;
@@ -238,7 +268,7 @@ _cairo_spline_error_squared (const cairo_spline_knots_t *knots)
 
 static cairo_status_t
 _cairo_spline_decompose_into (cairo_spline_knots_t *s1,
-			      double tolerance_squared,
+			      float tolerance_squared,
 			      cairo_spline_t *result)
 {
     cairo_spline_knots_t s2;
@@ -256,6 +286,26 @@ _cairo_spline_decompose_into (cairo_spline_knots_t *s1,
     return _cairo_spline_decompose_into (&s2, tolerance_squared, result);
 }
 
+static cairo_status_t
+_cairo_spline_decompose_fast_into (cairo_spline_knots_t *s1,
+			           float tolerance_squared,
+			           cairo_spline_t *result)
+{
+    cairo_spline_knots_t s2;
+    cairo_status_t status;
+
+    if (_cairo_spline_error_squared (s1) < tolerance_squared)
+	return _cairo_spline_add_point (result, &s1->a, &s1->b);
+
+    _de_casteljau_fast (s1, &s2);
+
+    status = _cairo_spline_decompose_fast_into (s1, tolerance_squared, result);
+    if (unlikely (status))
+	return status;
+
+    return _cairo_spline_decompose_fast_into (&s2, tolerance_squared, result);
+}
+
 cairo_status_t
 _cairo_spline_decompose (cairo_spline_t *spline, double tolerance)
 {
@@ -264,7 +314,14 @@ _cairo_spline_decompose (cairo_spline_t *spline, double tolerance)
 
     s1 = spline->knots;
     spline->last_point = s1.a;
-    status = _cairo_spline_decompose_into (&s1, tolerance * tolerance, spline);
+    if (abs(s1.a.x) < INT_MAX && abs(s1.a.y) < INT_MAX &&
+        abs(s1.b.x) < INT_MAX && abs(s1.b.y) < INT_MAX &&
+        abs(s1.c.x) < INT_MAX && abs(s1.c.y) < INT_MAX &&
+        abs(s1.d.x) < INT_MAX && abs(s1.d.y) < INT_MAX)
+        status = _cairo_spline_decompose_fast_into (&s1, (float)(tolerance * tolerance), spline);
+    else
+        status = _cairo_spline_decompose_into (&s1, (float)(tolerance * tolerance), spline);
+
     if (unlikely (status))
 	return status;
 
-- 
1.7.2.5