1 /* cairo - a vector graphics library with display and print output
3 * Copyright © 2002 University of Southern California
5 * This library is free software; you can redistribute it and/or
6 * modify it either under the terms of the GNU Lesser General Public
7 * License version 2.1 as published by the Free Software Foundation
8 * (the "LGPL") or, at your option, under the terms of the Mozilla
9 * Public License Version 1.1 (the "MPL"). If you do not alter this
10 * notice, a recipient may use your version of this file under either
11 * the MPL or the LGPL.
13 * You should have received a copy of the LGPL along with this library
14 * in the file COPYING-LGPL-2.1; if not, write to the Free Software
15 * Foundation, Inc., 51 Franklin Street, Suite 500, Boston, MA 02110-1335, USA
16 * You should have received a copy of the MPL along with this library
17 * in the file COPYING-MPL-1.1
19 * The contents of this file are subject to the Mozilla Public License
20 * Version 1.1 (the "License"); you may not use this file except in
21 * compliance with the License. You may obtain a copy of the License at
22 * http://www.mozilla.org/MPL/
24 * This software is distributed on an "AS IS" basis, WITHOUT WARRANTY
25 * OF ANY KIND, either express or implied. See the LGPL or the MPL for
26 * the specific language governing rights and limitations.
28 * The Original Code is the cairo graphics library.
30 * The Initial Developer of the Original Code is University of Southern
34 * Carl D. Worth <cworth@cworth.org>
39 #include "cairo-box-inline.h"
40 #include "cairo-slope-private.h"
43 _cairo_spline_intersects (const cairo_point_t *a,
44 const cairo_point_t *b,
45 const cairo_point_t *c,
46 const cairo_point_t *d,
47 const cairo_box_t *box)
51 if (_cairo_box_contains_point (box, a) ||
52 _cairo_box_contains_point (box, b) ||
53 _cairo_box_contains_point (box, c) ||
54 _cairo_box_contains_point (box, d))
59 bounds.p2 = bounds.p1 = *a;
60 _cairo_box_add_point (&bounds, b);
61 _cairo_box_add_point (&bounds, c);
62 _cairo_box_add_point (&bounds, d);
64 if (bounds.p2.x <= box->p1.x || bounds.p1.x >= box->p2.x ||
65 bounds.p2.y <= box->p1.y || bounds.p1.y >= box->p2.y)
70 #if 0 /* worth refining? */
71 bounds.p2 = bounds.p1 = *a;
72 _cairo_box_add_curve_to (&bounds, b, c, d);
73 if (bounds.p2.x <= box->p1.x || bounds.p1.x >= box->p2.x ||
74 bounds.p2.y <= box->p1.y || bounds.p1.y >= box->p2.y)
84 _cairo_spline_init (cairo_spline_t *spline,
85 cairo_spline_add_point_func_t add_point_func,
87 const cairo_point_t *a, const cairo_point_t *b,
88 const cairo_point_t *c, const cairo_point_t *d)
90 /* If both tangents are zero, this is just a straight line */
91 if (a->x == b->x && a->y == b->y && c->x == d->x && c->y == d->y)
94 spline->add_point_func = add_point_func;
95 spline->closure = closure;
100 spline->knots.d = *d;
102 if (a->x != b->x || a->y != b->y)
103 _cairo_slope_init (&spline->initial_slope, &spline->knots.a, &spline->knots.b);
104 else if (a->x != c->x || a->y != c->y)
105 _cairo_slope_init (&spline->initial_slope, &spline->knots.a, &spline->knots.c);
106 else if (a->x != d->x || a->y != d->y)
107 _cairo_slope_init (&spline->initial_slope, &spline->knots.a, &spline->knots.d);
111 if (c->x != d->x || c->y != d->y)
112 _cairo_slope_init (&spline->final_slope, &spline->knots.c, &spline->knots.d);
113 else if (b->x != d->x || b->y != d->y)
114 _cairo_slope_init (&spline->final_slope, &spline->knots.b, &spline->knots.d);
116 return FALSE; /* just treat this as a straight-line from a -> d */
118 /* XXX if the initial, final and vector are all equal, this is just a line */
123 static cairo_status_t
124 _cairo_spline_add_point (cairo_spline_t *spline,
125 const cairo_point_t *point,
126 const cairo_point_t *knot)
131 prev = &spline->last_point;
132 if (prev->x == point->x && prev->y == point->y)
133 return CAIRO_STATUS_SUCCESS;
135 _cairo_slope_init (&slope, point, knot);
137 spline->last_point = *point;
138 return spline->add_point_func (spline->closure, point, &slope);
142 _lerp_half (const cairo_point_t *a, const cairo_point_t *b, cairo_point_t *result)
144 result->x = a->x + ((b->x - a->x) >> 1);
145 result->y = a->y + ((b->y - a->y) >> 1);
149 _de_casteljau (cairo_spline_knots_t *s1, cairo_spline_knots_t *s2)
151 cairo_point_t ab, bc, cd;
152 cairo_point_t abbc, bccd;
155 _lerp_half (&s1->a, &s1->b, &ab);
156 _lerp_half (&s1->b, &s1->c, &bc);
157 _lerp_half (&s1->c, &s1->d, &cd);
158 _lerp_half (&ab, &bc, &abbc);
159 _lerp_half (&bc, &cd, &bccd);
160 _lerp_half (&abbc, &bccd, &final);
172 /* Return an upper bound on the error (squared) that could result from
173 * approximating a spline as a line segment connecting the two endpoints. */
175 _cairo_spline_error_squared (const cairo_spline_knots_t *knots)
177 double bdx, bdy, berr;
178 double cdx, cdy, cerr;
180 /* We are going to compute the distance (squared) between each of the the b
181 * and c control points and the segment a-b. The maximum of these two
182 * distances will be our approximation error. */
184 bdx = _cairo_fixed_to_double (knots->b.x - knots->a.x);
185 bdy = _cairo_fixed_to_double (knots->b.y - knots->a.y);
187 cdx = _cairo_fixed_to_double (knots->c.x - knots->a.x);
188 cdy = _cairo_fixed_to_double (knots->c.y - knots->a.y);
190 if (knots->a.x != knots->d.x || knots->a.y != knots->d.y) {
191 /* Intersection point (px):
192 * px = p1 + u(p2 - p1)
193 * (p - px) ∙ (p2 - p1) = 0
195 * u = ((p - p1) ∙ (p2 - p1)) / ∥p2 - p1∥²;
200 dx = _cairo_fixed_to_double (knots->d.x - knots->a.x);
201 dy = _cairo_fixed_to_double (knots->d.y - knots->a.y);
202 v = dx * dx + dy * dy;
204 u = bdx * dx + bdy * dy;
217 u = cdx * dx + cdy * dy;
231 berr = bdx * bdx + bdy * bdy;
232 cerr = cdx * cdx + cdy * cdy;
239 static cairo_status_t
240 _cairo_spline_decompose_into (cairo_spline_knots_t *s1,
241 double tolerance_squared,
242 cairo_spline_t *result)
244 cairo_spline_knots_t s2;
245 cairo_status_t status;
247 if (_cairo_spline_error_squared (s1) < tolerance_squared)
248 return _cairo_spline_add_point (result, &s1->a, &s1->b);
250 _de_casteljau (s1, &s2);
252 status = _cairo_spline_decompose_into (s1, tolerance_squared, result);
253 if (unlikely (status))
256 return _cairo_spline_decompose_into (&s2, tolerance_squared, result);
260 _cairo_spline_decompose (cairo_spline_t *spline, double tolerance)
262 cairo_spline_knots_t s1;
263 cairo_status_t status;
266 spline->last_point = s1.a;
267 status = _cairo_spline_decompose_into (&s1, tolerance * tolerance, spline);
268 if (unlikely (status))
271 return spline->add_point_func (spline->closure,
272 &spline->knots.d, &spline->final_slope);
275 /* Note: this function is only good for computing bounds in device space. */
277 _cairo_spline_bound (cairo_spline_add_point_func_t add_point_func,
279 const cairo_point_t *p0, const cairo_point_t *p1,
280 const cairo_point_t *p2, const cairo_point_t *p3)
282 double x0, x1, x2, x3;
283 double y0, y1, y2, y3;
287 cairo_status_t status;
289 x0 = _cairo_fixed_to_double (p0->x);
290 y0 = _cairo_fixed_to_double (p0->y);
291 x1 = _cairo_fixed_to_double (p1->x);
292 y1 = _cairo_fixed_to_double (p1->y);
293 x2 = _cairo_fixed_to_double (p2->x);
294 y2 = _cairo_fixed_to_double (p2->y);
295 x3 = _cairo_fixed_to_double (p3->x);
296 y3 = _cairo_fixed_to_double (p3->y);
298 /* The spline can be written as a polynomial of the four points:
300 * (1-t)³p0 + 3t(1-t)²p1 + 3t²(1-t)p2 + t³p3
302 * for 0≤t≤1. Now, the X and Y components of the spline follow the
303 * same polynomial but with x and y replaced for p. To find the
304 * bounds of the spline, we just need to find the X and Y bounds.
305 * To find the bound, we take the derivative and equal it to zero,
306 * and solve to find the t's that give the extreme points.
308 * Here is the derivative of the curve, sorted on t:
310 * 3t²(-p0+3p1-3p2+p3) + 2t(3p0-6p1+3p2) -3p0+3p1
320 * a.t² + 2b.t + c = 0
326 * the extreme points are at -c/2b if a is zero, at (-b±√delta)/a if
327 * delta is positive, and at -b/a if delta is zero.
333 if (0 < _t0 && _t0 < 1) \
337 #define FIND_EXTREMES(a,b,c) \
344 double delta = b2 - a * c; \
346 cairo_bool_t feasible; \
347 double _2ab = 2 * a * b; \
348 /* We are only interested in solutions t that satisfy 0<t<1 \
349 * here. We do some checks to avoid sqrt if the solutions \
350 * are not in that range. The checks can be derived from: \
352 * 0 < (-b±√delta)/a < 1 \
355 feasible = delta > b2 && delta < a*a + b2 + _2ab; \
356 else if (-b / a >= 1) \
357 feasible = delta < b2 && delta > a*a + b2 + _2ab; \
359 feasible = delta < b2 || delta < a*a + b2 + _2ab; \
361 if (unlikely (feasible)) { \
362 double sqrt_delta = sqrt (delta); \
363 ADD ((-b - sqrt_delta) / a); \
364 ADD ((-b + sqrt_delta) / a); \
366 } else if (delta == 0) { \
372 /* Find X extremes */
373 a = -x0 + 3*x1 - 3*x2 + x3;
376 FIND_EXTREMES (a, b, c);
378 /* Find Y extremes */
379 a = -y0 + 3*y1 - 3*y2 + y3;
382 FIND_EXTREMES (a, b, c);
384 status = add_point_func (closure, p0, NULL);
385 if (unlikely (status))
388 for (i = 0; i < t_num; i++) {
393 double t_3_0, t_2_1_3, t_1_2_3, t_0_3;
395 t_1_0 = t[i]; /* t */
396 t_0_1 = 1 - t_1_0; /* (1 - t) */
398 t_2_0 = t_1_0 * t_1_0; /* t * t */
399 t_0_2 = t_0_1 * t_0_1; /* (1 - t) * (1 - t) */
401 t_3_0 = t_2_0 * t_1_0; /* t * t * t */
402 t_2_1_3 = t_2_0 * t_0_1 * 3; /* t * t * (1 - t) * 3 */
403 t_1_2_3 = t_1_0 * t_0_2 * 3; /* t * (1 - t) * (1 - t) * 3 */
404 t_0_3 = t_0_1 * t_0_2; /* (1 - t) * (1 - t) * (1 - t) */
406 /* Bezier polynomial */
416 p.x = _cairo_fixed_from_double (x);
417 p.y = _cairo_fixed_from_double (y);
418 status = add_point_func (closure, &p, NULL);
419 if (unlikely (status))
423 return add_point_func (closure, p3, NULL);