build: updated gschemas.compiled
[adg.git] / src / cpml / cpml-curve.c
blob1c37a4351a1f3ac6553610d22df3404fc0dc8705
1 /* CPML - Cairo Path Manipulation Library
2 * Copyright (C) 2007,2008,2009,2010,2011,2012,2013 Nicola Fontana <ntd at entidi.it>
4 * This library is free software; you can redistribute it and/or
5 * modify it under the terms of the GNU Lesser General Public
6 * License as published by the Free Software Foundation; either
7 * version 2 of the License, or (at your option) any later version.
9 * This library is distributed in the hope that it will be useful,
10 * but WITHOUT ANY WARRANTY; without even the implied warranty of
11 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
12 * Lesser General Public License for more details.
14 * You should have received a copy of the GNU Lesser General Public
15 * License along with this library; if not, write to the
16 * Free Software Foundation, Inc., 51 Franklin Street, Fifth Floor,
17 * Boston, MA 02110-1301, USA.
21 /**
22 * SECTION:cpml-curve
23 * @Section_Id:CpmlCurve
24 * @title: CpmlCurve
25 * @short_description: Bézier cubic curve primitive management
27 * The following functions manipulate %CPML_CURVE #CpmlPrimitive.
28 * No validation is made on the input so use the following methods
29 * only when you are sure the <varname>primitive</varname> argument
30 * is effectively a cubic Bézier curve.
32 * <important>
33 * <title>TODO</title>
34 * <itemizedlist>
35 * <listitem>the get_length() method must be implemented;</listitem>
36 * <listitem>actually the put_extents() method is implemented by computing
37 * the bounding box of the control polygon and this will likely
38 * include some empty space: there is room for improvements;</listitem>
39 * <listitem>the put_pair_at() method must be implemented;</listitem>
40 * <listitem>the put_vector_at() method must be implemented;</listitem>
41 * <listitem>the get_closest_pos() method must be implemented;</listitem>
42 * <listitem>the put_intersections() method must be implemented;</listitem>
43 * <listitem>by default, the offset curve is calculated by using the point
44 * at t=0.5 as reference: use a better candidate;</listitem>
45 * <listitem>in the offset() implementation, when the equations are
46 * inconsistent, the alternative approach performs very bad
47 * if <varname>v0</varname> and <varname>v3</varname> are
48 * opposite or staggered.</listitem>
49 * </itemizedlist>
50 * </important>
51 * </para>
52 * <refsect2 id="offset">
53 * <title>Offseting algorithm</title>
54 * <para>
55 * Given a cubic Bézier primitive, it must be found the approximated
56 * Bézier curve parallel to the original one at a specific distance
57 * (the so called "offset curve"). The four points needed to build
58 * the new curve must be returned.
60 * To solve the offset problem, a custom algorithm is used. First, the
61 * resulting curve MUST have the same slope at the start and end point.
62 * These constraints are not sufficient to resolve the system, so I let
63 * the curve pass thought a given point (pm, known and got from the
64 * original curve) at a given time (m, now hardcoded to 0.5).
66 * Firstly, some useful variables are defined:
68 * |[
69 * v0 = unitvector(p[1] &minus; p[0]) &times; offset;
70 * v3 = unitvector(p[3] &minus; p[2]) &times; offset;
71 * p0 = p[0] + normal v0;
72 * p3 = p[3] + normal v3.
73 * ]|
75 * The resulting curve must have the same slopes than the original
76 * one at the start and end points. Forcing the same slopes means:
78 * |[
79 * p1 = p0 + k0 v0.
80 * ]|
82 * where %k0 is an arbitrary factor. Decomposing for %x and %y:
84 * |[
85 * p1.x = p0.x + k0 v0.x;
86 * p1.y = p0.y + k0 v0.y.
87 * ]|
89 * and doing the same for the end point:
91 * |[
92 * p2.x = p3.x + k3 v3.x;
93 * p2.y = p3.y + k3 v3.y.
94 * ]|
96 * This does not give a resolvable system, though. The curve will be
97 * interpolated by forcing its path to pass throught %pm when
98 * <varname>time</varname> is %m, where <code>0 &le; m &le; 1</code>.
99 * Knowing the function of the cubic Bézier:
101 * |[
102 * C(t) = (1 &minus; t)³p0 + 3t(1 &minus; t)²p1 + 3t²(1 &minus; t)p2 + t³p3.
103 * ]|
105 * and forcing <code>t = m</code> and <code>C(t) = pm</code>:
107 * |[
108 * pm = (1 &minus; m)³p0 + 3m(1 &minus; m)²p1 + 3m²(1 &minus; m)p2 + m³p3.
110 * (1 &minus; m) p1 + m p2 = (pm &minus; (1 &minus; m)³p0 &minus; m³p3) / (3m (1 &minus; m)).
111 * ]|
113 * gives this final system:
115 * |[
116 * p1.x = p0.x + k0 v0.x;
117 * p1.y = p0.y + k0 v0.y;
118 * p2.x = p3.x + k3 v3.x;
119 * p2.y = p3.y + k3 v3.y;
120 * (1 &minus; m) p1.x + m p2.x = (pm.x &minus; (1 &minus; m)³p0.x &minus; m³p3.x) / (3m (1 &minus; m));
121 * (1 &minus; m) p1.y + m p2.y = (pm.y &minus; (1 &minus; m)³p0.y &minus; m³p3.y) / (3m (1 &minus; m)).
122 * ]|
124 * Substituting and resolving for %k0 and %k3:
126 * |[
127 * (1 &minus; m) k0 v0.x + m k3 v3.x = (pm.x &minus; (1 &minus; m)³p0.x &minus; m³p3.x) / (3m (1 &minus; m)) &minus; (1 &minus; m) p0.x &minus; m p3.x;
128 * (1 &minus; m) k0 v0.y + m k3 v3.y = (pm.y &minus; (1 &minus; m)³p0.y &minus; m³p3.y) / (3m (1 &minus; m)) &minus; (1 &minus; m) p0.y &minus; m p3.y.
130 * (1 &minus; m) k0 v0.x + m k3 v3.x = (pm.x &minus; (1 &minus; m)²(1+2m) p0.x &minus; m²(3 &minus; 2m) p3.x) / (3m (1 &minus; m));
131 * (1 &minus; m) k0 v0.y + m k3 v3.y = (pm.y &minus; (1 &minus; m)²(1+2m) p0.y &minus; m²(3 &minus; 2m) p3.y) / (3m (1 &minus; m)).
132 * ]|
134 * Letting:
136 * |[
137 * pk = (pm &minus; (1 &minus; m)²(1+2m) p0 &minus; m²(3 &minus; 2m) p3) / (3m (1 &minus; m)).
138 * ]|
140 * reduces the above to this final equations:
142 * |[
143 * (1 &minus; m) k0 v0.x + m k3 v3.x = pk.x;
144 * (1 &minus; m) k0 v0.y + m k3 v3.y = pk.y.
145 * ]|
147 * If <code>v0.x &ne; 0</code>, the system can be resolved for %k0 and
148 * %k3 calculated accordingly:
150 * |[
151 * k0 = (pk.x &minus; m k3 v3.x) / ((1 &minus; m) v0.x);
152 * (pk.x &minus; m k3 v3.x) v0.y / v0.x + m k3 v3.y = pk.y.
154 * k0 = (pk.x &minus; m k3 v3.x) / ((1 &minus; m) v0.x);
155 * k3 m (v3.y &minus; v3.x v0.y / v0.x) = pk.y &minus; pk.x v0.y / v0.x.
157 * k3 = (pk.y &minus; pk.x v0.y / v0.x) / (m (v3.y &minus; v3.x v0.y / v0.x));
158 * k0 = (pk.x &minus; m k3 v3.x) / ((1 &minus; m) v0.x).
159 * ]|
161 * Otherwise, if <code>v3.x &ne; 0</code>, the system can be solved
162 * for %k3 and %k0 calculated accordingly:
164 * |[
165 * k3 = (pk.x &minus; (1 &minus; m) k0 v0.x) / (m v3.x);
166 * (1 &minus; m) k0 v0.y + (pk.x &minus; (1 &minus; m) k0 v0.x) v3.y / v3.x = pk.y.
168 * k3 = (pk.x &minus; (1 &minus; m) k0 v0.x) / (m v3.x);
169 * k0 (1 &minus; m) (v0.y &minus; k0 v0.x v3.y / v3.x) = pk.y &minus; pk.x v3.y / v3.x.
171 * k0 = (pk.y &minus; pk.x v3.y / v3.x) / ((1 &minus; m) (v0.y &minus; v0.x v3.y / v3.x));
172 * k3 = (pk.x &minus; (1 &minus; m) k0 v0.x) / (m v3.x).
173 * ]|
175 * The whole process must be guarded against division by 0 exceptions.
176 * If either <code>v0.x</code> and <code>v3.x</code> are %0, the first
177 * equation will be inconsistent. More in general, the
178 * <code>v0.x &times; v3.y = v3.x &times; v3.y</code> condition must
179 * be avoided. This is the first situation to avoid, in which case
180 * an alternative approach should be used.
182 * </para>
183 * </refsect2>
184 * <para>
186 * Since: 1.0
190 #include "cpml-internal.h"
191 #include "cpml-extents.h"
192 #include "cpml-segment.h"
193 #include "cpml-primitive.h"
194 #include "cpml-primitive-private.h"
195 #include "cpml-curve.h"
198 static void put_extents (const CpmlPrimitive *curve,
199 CpmlExtents *extents);
200 static void offset (CpmlPrimitive *curve,
201 double offset);
204 const _CpmlPrimitiveClass *
205 _cpml_curve_get_class(void)
207 static _CpmlPrimitiveClass *p_class = NULL;
209 if (p_class == NULL) {
210 static _CpmlPrimitiveClass class_data = {
211 "curve to", 4,
212 NULL,
213 put_extents,
214 NULL,
215 NULL,
216 NULL,
217 NULL,
218 offset,
219 NULL
221 p_class = &class_data;
224 return p_class;
229 * cpml_curve_put_pair_at_time:
230 * @curve: the #CpmlPrimitive curve data
231 * @t: the "time" value
232 * @pair: the destination pair
234 * Given the @curve Bézier cubic, finds the coordinates at time @t
235 * (where 0 is the start and 1 is the end) and stores the result
236 * in @pair. Keep in mind @t is not homogeneous, so 0.5 does not
237 * necessarily means the mid point.
239 * Since: 1.0
241 void
242 cpml_curve_put_pair_at_time(const CpmlPrimitive *curve, double t,
243 CpmlPair *pair)
245 CpmlPair p1, p2, p3, p4;
246 double t_2, t_3, t1, t1_2, t1_3;
248 cpml_primitive_put_point(curve, 0, &p1);
249 cpml_primitive_put_point(curve, 1, &p2);
250 cpml_primitive_put_point(curve, 2, &p3);
251 cpml_primitive_put_point(curve, 3, &p4);
253 t_2 = t * t;
254 t_3 = t_2 * t;
255 t1 = 1 - t;
256 t1_2 = t1 * t1;
257 t1_3 = t1_2 * t1;
259 pair->x = t1_3 * p1.x + 3 * t1_2 * t * p2.x
260 + 3 * t1 * t_2 * p3.x + t_3 * p4.x;
261 pair->y = t1_3 * p1.y + 3 * t1_2 * t * p2.y
262 + 3 * t1 * t_2 * p3.y + t_3 * p4.y;
266 * cpml_curve_put_vector_at_time:
267 * @curve: the #CpmlPrimitive curve data
268 * @t: the "time" value
269 * @vector: the destination vector
271 * Given the @curve Bézier cubic, finds the slope at time @t
272 * (where 0 is the start and 1 is the end) and stores the result
273 * in @vector. Keep in mind @t is not homogeneous, so 0.5
274 * does not necessarily means the mid point.
276 * Since: 1.0
278 void
279 cpml_curve_put_vector_at_time(const CpmlPrimitive *curve,
280 double t, CpmlVector *vector)
282 CpmlPair p1, p2, p3, p4;
283 CpmlPair p21, p32, p43;
284 double t1, t1_2, t_2;
286 cpml_primitive_put_point(curve, 0, &p1);
287 cpml_primitive_put_point(curve, 1, &p2);
288 cpml_primitive_put_point(curve, 2, &p3);
289 cpml_primitive_put_point(curve, 3, &p4);
291 p21.x = p2.x - p1.x;
292 p21.y = p2.y - p1.y;
293 p32.x = p3.x - p2.x;
294 p32.y = p3.y - p2.y;
295 p43.x = p4.x - p3.x;
296 p43.y = p4.y - p3.y;
298 t1 = 1 - t;
299 t1_2 = t1 * t1;
300 t_2 = t * t;
302 vector->x = 3 * t1_2 * p21.x + 6 * t1 * t * p32.x + 3 * t_2 * p43.x;
303 vector->y = 3 * t1_2 * p21.y + 6 * t1 * t * p32.y + 3 * t_2 * p43.y;
307 static void
308 put_extents(const CpmlPrimitive *curve, CpmlExtents *extents)
310 CpmlPair p1, p2, p3, p4;
312 extents->is_defined = 0;
314 cpml_primitive_put_point(curve, 0, &p1);
315 cpml_primitive_put_point(curve, 1, &p2);
316 cpml_primitive_put_point(curve, 2, &p3);
317 cpml_primitive_put_point(curve, 3, &p4);
319 cpml_extents_pair_add(extents, &p1);
320 cpml_extents_pair_add(extents, &p2);
321 cpml_extents_pair_add(extents, &p3);
322 cpml_extents_pair_add(extents, &p4);
325 static void
326 offset(CpmlPrimitive *curve, double offset)
328 double m, mm;
329 CpmlVector v0, v3, vm, vtmp;
330 CpmlPair p0, p1, p2, p3, pm;
332 m = 0.5;
333 mm = 1-m;
335 /* Firstly, convert the curve points from cairo format to cpml format
336 * and store them (temporary) in p0..p3 */
337 cpml_pair_from_cairo(&p0, curve->org);
338 cpml_pair_from_cairo(&p1, &curve->data[1]);
339 cpml_pair_from_cairo(&p2, &curve->data[2]);
340 cpml_pair_from_cairo(&p3, &curve->data[3]);
342 /* v0 = p1-p0 */
343 v0.x = p1.x - p0.x;
344 v0.y = p1.y - p0.y;
346 /* v3 = p3-p2 */
347 v3.x = p3.x - p2.x;
348 v3.y = p3.y - p2.y;
350 /* pm = point in C(m) offseted the requested @offset distance */
351 cpml_curve_put_vector_at_time(curve, m, &vm);
352 cpml_vector_set_length(&vm, offset);
353 cpml_vector_normal(&vm);
354 cpml_curve_put_pair_at_time(curve, m, &pm);
355 pm.x += vm.x;
356 pm.y += vm.y;
358 /* p0 = p0 + normal of v0 of @offset magnitude (exact value) */
359 cpml_pair_copy(&vtmp, &v0);
360 cpml_vector_set_length(&vtmp, offset);
361 cpml_vector_normal(&vtmp);
362 p0.x += vtmp.x;
363 p0.y += vtmp.y;
365 /* p3 = p3 + normal of v3 of @offset magnitude, as done for p0 */
366 cpml_pair_copy(&vtmp, &v3);
367 cpml_vector_set_length(&vtmp, offset);
368 cpml_vector_normal(&vtmp);
369 p3.x += vtmp.x;
370 p3.y += vtmp.y;
372 if (v0.x*v3.y == v3.x*v0.y) {
373 /* Inconsistent equations: use the alternative approach */
374 p1.x = p0.x + v0.x + vm.x * 4/3;
375 p1.y = p0.y + v0.y + vm.y * 4/3;
376 p2.x = p3.x - v3.x + vm.x * 4/3;
377 p2.y = p3.y - v3.y + vm.y * 4/3;
378 } else {
379 CpmlPair pk;
380 double k0, k3;
382 pk.x = (pm.x - mm*mm*(1+m+m)*p0.x - m*m*(1+mm+mm)*p3.x) / (3*m*(1-m));
383 pk.y = (pm.y - mm*mm*(1+m+m)*p0.y - m*m*(1+mm+mm)*p3.y) / (3*m*(1-m));
385 if (v0.x != 0) {
386 k3 = (pk.y - pk.x*v0.y / v0.x) / (m*(v3.y - v3.x*v0.y / v0.x));
387 k0 = (pk.x - m*k3*v3.x) / (mm*v0.x);
388 } else {
389 k0 = (pk.y - pk.x*v3.y / v3.x) / (mm*(v0.y - v0.x*v3.y / v3.x));
390 k3 = (pk.x - mm*k0*v0.x) / (m*v3.x);
393 p1.x = p0.x + k0*v0.x;
394 p1.y = p0.y + k0*v0.y;
395 p2.x = p3.x + k3*v3.x;
396 p2.y = p3.y + k3*v3.y;
399 /* Return the new curve in the original array */
400 cpml_pair_to_cairo(&p0, curve->org);
401 cpml_pair_to_cairo(&p1, &curve->data[1]);
402 cpml_pair_to_cairo(&p2, &curve->data[2]);
403 cpml_pair_to_cairo(&p3, &curve->data[3]);