1 /* -*- Mode: c; c-basic-offset: 4; tab-width: 8; indent-tabs-mode: t; -*- */
4 * Copyright © 2000 Keith Packard, member of The XFree86 Project, Inc.
5 * Copyright © 2000 SuSE, Inc.
6 * 2005 Lars Knoll & Zack Rusin, Trolltech
7 * Copyright © 2007 Red Hat, Inc.
10 * Permission to use, copy, modify, distribute, and sell this software and its
11 * documentation for any purpose is hereby granted without fee, provided that
12 * the above copyright notice appear in all copies and that both that
13 * copyright notice and this permission notice appear in supporting
14 * documentation, and that the name of Keith Packard not be used in
15 * advertising or publicity pertaining to distribution of the software without
16 * specific, written prior permission. Keith Packard makes no
17 * representations about the suitability of this software for any purpose. It
18 * is provided "as is" without express or implied warranty.
20 * THE COPYRIGHT HOLDERS DISCLAIM ALL WARRANTIES WITH REGARD TO THIS
21 * SOFTWARE, INCLUDING ALL IMPLIED WARRANTIES OF MERCHANTABILITY AND
22 * FITNESS, IN NO EVENT SHALL THE COPYRIGHT HOLDERS BE LIABLE FOR ANY
23 * SPECIAL, INDIRECT OR CONSEQUENTIAL DAMAGES OR ANY DAMAGES
24 * WHATSOEVER RESULTING FROM LOSS OF USE, DATA OR PROFITS, WHETHER IN
25 * AN ACTION OF CONTRACT, NEGLIGENCE OR OTHER TORTIOUS ACTION, ARISING
26 * OUT OF OR IN CONNECTION WITH THE USE OR PERFORMANCE OF THIS
35 #include "pixman-private.h"
37 static inline pixman_fixed_32_32_t
38 dot (pixman_fixed_48_16_t x1
,
39 pixman_fixed_48_16_t y1
,
40 pixman_fixed_48_16_t z1
,
41 pixman_fixed_48_16_t x2
,
42 pixman_fixed_48_16_t y2
,
43 pixman_fixed_48_16_t z2
)
46 * Exact computation, assuming that the input values can
47 * be represented as pixman_fixed_16_16_t
49 return x1
* x2
+ y1
* y2
+ z1
* z2
;
61 * Error can be unbound in some special cases.
62 * Using clever dot product algorithms (for example compensated
63 * dot product) would improve this but make the code much less
66 return x1
* x2
+ y1
* y2
+ z1
* z2
;
70 radial_compute_color (double a
,
76 pixman_gradient_walker_t
*walker
,
77 pixman_repeat_t repeat
)
80 * In this function error propagation can lead to bad results:
81 * - discr can have an unbound error (if b*b-a*c is very small),
82 * potentially making it the opposite sign of what it should have been
83 * (thus clearing a pixel that would have been colored or vice-versa)
84 * or propagating the error to sqrtdiscr;
85 * if discr has the wrong sign or b is very small, this can lead to bad
88 * - the algorithm used to compute the solutions of the quadratic
89 * equation is not numerically stable (but saves one division compared
90 * to the numerically stable one);
91 * this can be a problem if a*c is much smaller than b*b
93 * - the above problems are worse if a is small (as inva becomes bigger)
104 t
= pixman_fixed_1
/ 2 * c
/ b
;
105 if (repeat
== PIXMAN_REPEAT_NONE
)
107 if (0 <= t
&& t
<= pixman_fixed_1
)
108 return _pixman_gradient_walker_pixel (walker
, t
);
113 return _pixman_gradient_walker_pixel (walker
, t
);
119 discr
= fdot (b
, a
, 0, b
, -c
, 0);
122 double sqrtdiscr
, t0
, t1
;
124 sqrtdiscr
= sqrt (discr
);
125 t0
= (b
+ sqrtdiscr
) * inva
;
126 t1
= (b
- sqrtdiscr
) * inva
;
129 * The root that must be used is the biggest one that belongs
130 * to the valid range ([0,1] for PIXMAN_REPEAT_NONE, any
131 * solution that results in a positive radius otherwise).
133 * If a > 0, t0 is the biggest solution, so if it is valid, it
134 * is the correct result.
136 * If a < 0, only one of the solutions can be valid, so the
137 * order in which they are tested is not important.
139 if (repeat
== PIXMAN_REPEAT_NONE
)
141 if (0 <= t0
&& t0
<= pixman_fixed_1
)
142 return _pixman_gradient_walker_pixel (walker
, t0
);
143 else if (0 <= t1
&& t1
<= pixman_fixed_1
)
144 return _pixman_gradient_walker_pixel (walker
, t1
);
148 if (t0
* dr
>= mindr
)
149 return _pixman_gradient_walker_pixel (walker
, t0
);
150 else if (t1
* dr
>= mindr
)
151 return _pixman_gradient_walker_pixel (walker
, t1
);
159 radial_get_scanline_narrow (pixman_iter_t
*iter
, const uint32_t *mask
)
162 * Implementation of radial gradients following the PDF specification.
163 * See section 8.7.4.5.4 Type 3 (Radial) Shadings of the PDF Reference
164 * Manual (PDF 32000-1:2008 at the time of this writing).
166 * In the radial gradient problem we are given two circles (c₁,r₁) and
167 * (c₂,r₂) that define the gradient itself.
169 * Mathematically the gradient can be defined as the family of circles
171 * ((1-t)·c₁ + t·(c₂), (1-t)·r₁ + t·r₂)
173 * excluding those circles whose radius would be < 0. When a point
174 * belongs to more than one circle, the one with a bigger t is the only
175 * one that contributes to its color. When a point does not belong
176 * to any of the circles, it is transparent black, i.e. RGBA (0, 0, 0, 0).
177 * Further limitations on the range of values for t are imposed when
178 * the gradient is not repeated, namely t must belong to [0,1].
180 * The graphical result is the same as drawing the valid (radius > 0)
181 * circles with increasing t in [-inf, +inf] (or in [0,1] if the gradient
182 * is not repeated) using SOURCE operator composition.
184 * It looks like a cone pointing towards the viewer if the ending circle
185 * is smaller than the starting one, a cone pointing inside the page if
186 * the starting circle is the smaller one and like a cylinder if they
187 * have the same radius.
189 * What we actually do is, given the point whose color we are interested
190 * in, compute the t values for that point, solving for t in:
192 * length((1-t)·c₁ + t·(c₂) - p) = (1-t)·r₁ + t·r₂
194 * Let's rewrite it in a simpler way, by defining some auxiliary
200 * length(t·cd - pd) = r₁ + t·dr
202 * which actually means
204 * hypot(t·cdx - pdx, t·cdy - pdy) = r₁ + t·dr
208 * ⎷((t·cdx - pdx)² + (t·cdy - pdy)²) = r₁ + t·dr.
210 * If we impose (as stated earlier) that r₁ + t·dr >= 0, it becomes:
212 * (t·cdx - pdx)² + (t·cdy - pdy)² = (r₁ + t·dr)²
214 * where we can actually expand the squares and solve for t:
216 * t²cdx² - 2t·cdx·pdx + pdx² + t²cdy² - 2t·cdy·pdy + pdy² =
217 * = r₁² + 2·r₁·t·dr + t²·dr²
219 * (cdx² + cdy² - dr²)t² - 2(cdx·pdx + cdy·pdy + r₁·dr)t +
220 * (pdx² + pdy² - r₁²) = 0
222 * A = cdx² + cdy² - dr²
223 * B = pdx·cdx + pdy·cdy + r₁·dr
224 * C = pdx² + pdy² - r₁²
227 * The solutions (unless the equation degenerates because of A = 0) are:
229 * t = (B ± ⎷(B² - A·C)) / A
231 * The solution we are going to prefer is the bigger one, unless the
232 * radius associated to it is negative (or it falls outside the valid t
235 * Additional observations (useful for optimizations):
236 * A does not depend on p
238 * A < 0 <=> one of the two circles completely contains the other one
239 * <=> for every p, the radiuses associated with the two t solutions
242 pixman_image_t
*image
= iter
->image
;
245 int width
= iter
->width
;
246 uint32_t *buffer
= iter
->buffer
;
248 gradient_t
*gradient
= (gradient_t
*)image
;
249 radial_gradient_t
*radial
= (radial_gradient_t
*)image
;
250 uint32_t *end
= buffer
+ width
;
251 pixman_gradient_walker_t walker
;
252 pixman_vector_t v
, unit
;
254 /* reference point is the center of the pixel */
255 v
.vector
[0] = pixman_int_to_fixed (x
) + pixman_fixed_1
/ 2;
256 v
.vector
[1] = pixman_int_to_fixed (y
) + pixman_fixed_1
/ 2;
257 v
.vector
[2] = pixman_fixed_1
;
259 _pixman_gradient_walker_init (&walker
, gradient
, image
->common
.repeat
);
261 if (image
->common
.transform
)
263 if (!pixman_transform_point_3d (image
->common
.transform
, &v
))
266 unit
.vector
[0] = image
->common
.transform
->matrix
[0][0];
267 unit
.vector
[1] = image
->common
.transform
->matrix
[1][0];
268 unit
.vector
[2] = image
->common
.transform
->matrix
[2][0];
272 unit
.vector
[0] = pixman_fixed_1
;
277 if (unit
.vector
[2] == 0 && v
.vector
[2] == pixman_fixed_1
)
282 * t = (B ± ⎷(B² - A·C)) / A
286 * A = cdx² + cdy² - dr²
287 * B = pdx·cdx + pdy·cdy + r₁·dr
288 * C = pdx² + pdy² - r₁²
291 * Since we have an affine transformation, we know that (pdx, pdy)
292 * increase linearly with each pixel,
297 * we can then express B, C and det through multiple differentiation.
299 pixman_fixed_32_32_t b
, db
, c
, dc
, ddc
;
301 /* warning: this computation may overflow */
302 v
.vector
[0] -= radial
->c1
.x
;
303 v
.vector
[1] -= radial
->c1
.y
;
306 * B and C are computed and updated exactly.
307 * If fdot was used instead of dot, in the worst case it would
308 * lose 11 bits of precision in each of the multiplication and
309 * summing up would zero out all the bit that were preserved,
310 * thus making the result 0 instead of the correct one.
311 * This would mean a worst case of unbound relative error or
312 * about 2^10 absolute error
314 b
= dot (v
.vector
[0], v
.vector
[1], radial
->c1
.radius
,
315 radial
->delta
.x
, radial
->delta
.y
, radial
->delta
.radius
);
316 db
= dot (unit
.vector
[0], unit
.vector
[1], 0,
317 radial
->delta
.x
, radial
->delta
.y
, 0);
319 c
= dot (v
.vector
[0], v
.vector
[1],
320 -((pixman_fixed_48_16_t
) radial
->c1
.radius
),
321 v
.vector
[0], v
.vector
[1], radial
->c1
.radius
);
322 dc
= dot (2 * (pixman_fixed_48_16_t
) v
.vector
[0] + unit
.vector
[0],
323 2 * (pixman_fixed_48_16_t
) v
.vector
[1] + unit
.vector
[1],
325 unit
.vector
[0], unit
.vector
[1], 0);
326 ddc
= 2 * dot (unit
.vector
[0], unit
.vector
[1], 0,
327 unit
.vector
[0], unit
.vector
[1], 0);
331 if (!mask
|| *mask
++)
333 *buffer
= radial_compute_color (radial
->a
, b
, c
,
335 radial
->delta
.radius
,
338 image
->common
.repeat
);
351 * error propagation guarantees are much looser than in the affine case
355 if (!mask
|| *mask
++)
357 if (v
.vector
[2] != 0)
359 double pdx
, pdy
, invv2
, b
, c
;
361 invv2
= 1. * pixman_fixed_1
/ v
.vector
[2];
363 pdx
= v
.vector
[0] * invv2
- radial
->c1
.x
;
364 /* / pixman_fixed_1 */
366 pdy
= v
.vector
[1] * invv2
- radial
->c1
.y
;
367 /* / pixman_fixed_1 */
369 b
= fdot (pdx
, pdy
, radial
->c1
.radius
,
370 radial
->delta
.x
, radial
->delta
.y
,
371 radial
->delta
.radius
);
372 /* / pixman_fixed_1 / pixman_fixed_1 */
374 c
= fdot (pdx
, pdy
, -radial
->c1
.radius
,
375 pdx
, pdy
, radial
->c1
.radius
);
376 /* / pixman_fixed_1 / pixman_fixed_1 */
378 *buffer
= radial_compute_color (radial
->a
, b
, c
,
380 radial
->delta
.radius
,
383 image
->common
.repeat
);
393 v
.vector
[0] += unit
.vector
[0];
394 v
.vector
[1] += unit
.vector
[1];
395 v
.vector
[2] += unit
.vector
[2];
404 radial_get_scanline_wide (pixman_iter_t
*iter
, const uint32_t *mask
)
406 uint32_t *buffer
= radial_get_scanline_narrow (iter
, NULL
);
408 pixman_expand_to_float (
409 (argb_t
*)buffer
, buffer
, PIXMAN_a8r8g8b8
, iter
->width
);
415 _pixman_radial_gradient_iter_init (pixman_image_t
*image
, pixman_iter_t
*iter
)
417 if (iter
->iter_flags
& ITER_NARROW
)
418 iter
->get_scanline
= radial_get_scanline_narrow
;
420 iter
->get_scanline
= radial_get_scanline_wide
;
423 PIXMAN_EXPORT pixman_image_t
*
424 pixman_image_create_radial_gradient (const pixman_point_fixed_t
* inner
,
425 const pixman_point_fixed_t
* outer
,
426 pixman_fixed_t inner_radius
,
427 pixman_fixed_t outer_radius
,
428 const pixman_gradient_stop_t
*stops
,
431 pixman_image_t
*image
;
432 radial_gradient_t
*radial
;
434 image
= _pixman_image_allocate ();
439 radial
= &image
->radial
;
441 if (!_pixman_init_gradient (&radial
->common
, stops
, n_stops
))
447 image
->type
= RADIAL
;
449 radial
->c1
.x
= inner
->x
;
450 radial
->c1
.y
= inner
->y
;
451 radial
->c1
.radius
= inner_radius
;
452 radial
->c2
.x
= outer
->x
;
453 radial
->c2
.y
= outer
->y
;
454 radial
->c2
.radius
= outer_radius
;
456 /* warning: this computations may overflow */
457 radial
->delta
.x
= radial
->c2
.x
- radial
->c1
.x
;
458 radial
->delta
.y
= radial
->c2
.y
- radial
->c1
.y
;
459 radial
->delta
.radius
= radial
->c2
.radius
- radial
->c1
.radius
;
461 /* computed exactly, then cast to double -> every bit of the double
462 representation is correct (53 bits) */
463 radial
->a
= dot (radial
->delta
.x
, radial
->delta
.y
, -radial
->delta
.radius
,
464 radial
->delta
.x
, radial
->delta
.y
, radial
->delta
.radius
);
466 radial
->inva
= 1. * pixman_fixed_1
/ radial
->a
;
468 radial
->mindr
= -1. * pixman_fixed_1
* radial
->c1
.radius
;