diff options
Diffstat (limited to 'pixman/pixman/pixman-radial-gradient.c')
-rw-r--r-- | pixman/pixman/pixman-radial-gradient.c | 923 |
1 files changed, 460 insertions, 463 deletions
diff --git a/pixman/pixman/pixman-radial-gradient.c b/pixman/pixman/pixman-radial-gradient.c index 6523b8259..63c712cc2 100644 --- a/pixman/pixman/pixman-radial-gradient.c +++ b/pixman/pixman/pixman-radial-gradient.c @@ -1,463 +1,460 @@ -/* -*- Mode: c; c-basic-offset: 4; tab-width: 8; indent-tabs-mode: t; -*- */ -/* - * - * Copyright © 2000 Keith Packard, member of The XFree86 Project, Inc. - * Copyright © 2000 SuSE, Inc. - * 2005 Lars Knoll & Zack Rusin, Trolltech - * Copyright © 2007 Red Hat, Inc. - * - * - * Permission to use, copy, modify, distribute, and sell this software and its - * documentation for any purpose is hereby granted without fee, provided that - * the above copyright notice appear in all copies and that both that - * copyright notice and this permission notice appear in supporting - * documentation, and that the name of Keith Packard not be used in - * advertising or publicity pertaining to distribution of the software without - * specific, written prior permission. Keith Packard makes no - * representations about the suitability of this software for any purpose. It - * is provided "as is" without express or implied warranty. - * - * THE COPYRIGHT HOLDERS DISCLAIM ALL WARRANTIES WITH REGARD TO THIS - * SOFTWARE, INCLUDING ALL IMPLIED WARRANTIES OF MERCHANTABILITY AND - * FITNESS, IN NO EVENT SHALL THE COPYRIGHT HOLDERS BE LIABLE FOR ANY - * SPECIAL, INDIRECT OR CONSEQUENTIAL DAMAGES OR ANY DAMAGES - * WHATSOEVER RESULTING FROM LOSS OF USE, DATA OR PROFITS, WHETHER IN - * AN ACTION OF CONTRACT, NEGLIGENCE OR OTHER TORTIOUS ACTION, ARISING - * OUT OF OR IN CONNECTION WITH THE USE OR PERFORMANCE OF THIS - * SOFTWARE. - */ - -#ifdef HAVE_CONFIG_H -#include <config.h> -#endif -#include <stdlib.h> -#include <math.h> -#include "pixman-private.h" - -static inline pixman_fixed_32_32_t -dot (pixman_fixed_48_16_t x1, - pixman_fixed_48_16_t y1, - pixman_fixed_48_16_t z1, - pixman_fixed_48_16_t x2, - pixman_fixed_48_16_t y2, - pixman_fixed_48_16_t z2) -{ - /* - * Exact computation, assuming that the input values can - * be represented as pixman_fixed_16_16_t - */ - return x1 * x2 + y1 * y2 + z1 * z2; -} - -static inline double -fdot (double x1, - double y1, - double z1, - double x2, - double y2, - double z2) -{ - /* - * Error can be unbound in some special cases. - * Using clever dot product algorithms (for example compensated - * dot product) would improve this but make the code much less - * obvious - */ - return x1 * x2 + y1 * y2 + z1 * z2; -} - -static uint32_t -radial_compute_color (double a, - double b, - double c, - double inva, - double dr, - double mindr, - pixman_gradient_walker_t *walker, - pixman_repeat_t repeat) -{ - /* - * In this function error propagation can lead to bad results: - * - det can have an unbound error (if b*b-a*c is very small), - * potentially making it the opposite sign of what it should have been - * (thus clearing a pixel that would have been colored or vice-versa) - * or propagating the error to sqrtdet; - * if det has the wrong sign or b is very small, this can lead to bad - * results - * - * - the algorithm used to compute the solutions of the quadratic - * equation is not numerically stable (but saves one division compared - * to the numerically stable one); - * this can be a problem if a*c is much smaller than b*b - * - * - the above problems are worse if a is small (as inva becomes bigger) - */ - double det; - - if (a == 0) - { - double t; - - if (b == 0) - return 0; - - t = pixman_fixed_1 / 2 * c / b; - if (repeat == PIXMAN_REPEAT_NONE) - { - if (0 <= t && t <= pixman_fixed_1) - return _pixman_gradient_walker_pixel (walker, t); - } - else - { - if (t * dr > mindr) - return _pixman_gradient_walker_pixel (walker, t); - } - - return 0; - } - - det = fdot (b, a, 0, b, -c, 0); - if (det >= 0) - { - double sqrtdet, t0, t1; - - sqrtdet = sqrt (det); - t0 = (b + sqrtdet) * inva; - t1 = (b - sqrtdet) * inva; - - if (repeat == PIXMAN_REPEAT_NONE) - { - if (0 <= t0 && t0 <= pixman_fixed_1) - return _pixman_gradient_walker_pixel (walker, t0); - else if (0 <= t1 && t1 <= pixman_fixed_1) - return _pixman_gradient_walker_pixel (walker, t1); - } - else - { - if (t0 * dr > mindr) - return _pixman_gradient_walker_pixel (walker, t0); - else if (t1 * dr > mindr) - return _pixman_gradient_walker_pixel (walker, t1); - } - } - - return 0; -} - -static uint32_t * -radial_get_scanline_narrow (pixman_iter_t *iter, const uint32_t *mask) -{ - /* - * Implementation of radial gradients following the PDF specification. - * See section 8.7.4.5.4 Type 3 (Radial) Shadings of the PDF Reference - * Manual (PDF 32000-1:2008 at the time of this writing). - * - * In the radial gradient problem we are given two circles (c₁,r₁) and - * (c₂,r₂) that define the gradient itself. - * - * Mathematically the gradient can be defined as the family of circles - * - * ((1-t)·c₁ + t·(c₂), (1-t)·r₁ + t·r₂) - * - * excluding those circles whose radius would be < 0. When a point - * belongs to more than one circle, the one with a bigger t is the only - * one that contributes to its color. When a point does not belong - * to any of the circles, it is transparent black, i.e. RGBA (0, 0, 0, 0). - * Further limitations on the range of values for t are imposed when - * the gradient is not repeated, namely t must belong to [0,1]. - * - * The graphical result is the same as drawing the valid (radius > 0) - * circles with increasing t in [-inf, +inf] (or in [0,1] if the gradient - * is not repeated) using SOURCE operatior composition. - * - * It looks like a cone pointing towards the viewer if the ending circle - * is smaller than the starting one, a cone pointing inside the page if - * the starting circle is the smaller one and like a cylinder if they - * have the same radius. - * - * What we actually do is, given the point whose color we are interested - * in, compute the t values for that point, solving for t in: - * - * length((1-t)·c₁ + t·(c₂) - p) = (1-t)·r₁ + t·r₂ - * - * Let's rewrite it in a simpler way, by defining some auxiliary - * variables: - * - * cd = c₂ - c₁ - * pd = p - c₁ - * dr = r₂ - r₁ - * lenght(t·cd - pd) = r₁ + t·dr - * - * which actually means - * - * hypot(t·cdx - pdx, t·cdy - pdy) = r₁ + t·dr - * - * or - * - * ⎷((t·cdx - pdx)² + (t·cdy - pdy)²) = r₁ + t·dr. - * - * If we impose (as stated earlier) that r₁ + t·dr >= 0, it becomes: - * - * (t·cdx - pdx)² + (t·cdy - pdy)² = (r₁ + t·dr)² - * - * where we can actually expand the squares and solve for t: - * - * t²cdx² - 2t·cdx·pdx + pdx² + t²cdy² - 2t·cdy·pdy + pdy² = - * = r₁² + 2·r₁·t·dr + t²·dr² - * - * (cdx² + cdy² - dr²)t² - 2(cdx·pdx + cdy·pdy + r₁·dr)t + - * (pdx² + pdy² - r₁²) = 0 - * - * A = cdx² + cdy² - dr² - * B = pdx·cdx + pdy·cdy + r₁·dr - * C = pdx² + pdy² - r₁² - * At² - 2Bt + C = 0 - * - * The solutions (unless the equation degenerates because of A = 0) are: - * - * t = (B ± ⎷(B² - A·C)) / A - * - * The solution we are going to prefer is the bigger one, unless the - * radius associated to it is negative (or it falls outside the valid t - * range). - * - * Additional observations (useful for optimizations): - * A does not depend on p - * - * A < 0 <=> one of the two circles completely contains the other one - * <=> for every p, the radiuses associated with the two t solutions - * have opposite sign - */ - pixman_image_t *image = iter->image; - int x = iter->x; - int y = iter->y; - int width = iter->width; - uint32_t *buffer = iter->buffer; - - gradient_t *gradient = (gradient_t *)image; - radial_gradient_t *radial = (radial_gradient_t *)image; - uint32_t *end = buffer + width; - pixman_gradient_walker_t walker; - pixman_vector_t v, unit; - - /* reference point is the center of the pixel */ - v.vector[0] = pixman_int_to_fixed (x) + pixman_fixed_1 / 2; - v.vector[1] = pixman_int_to_fixed (y) + pixman_fixed_1 / 2; - v.vector[2] = pixman_fixed_1; - - _pixman_gradient_walker_init (&walker, gradient, image->common.repeat); - - if (image->common.transform) - { - if (!pixman_transform_point_3d (image->common.transform, &v)) - return iter->buffer; - - unit.vector[0] = image->common.transform->matrix[0][0]; - unit.vector[1] = image->common.transform->matrix[1][0]; - unit.vector[2] = image->common.transform->matrix[2][0]; - } - else - { - unit.vector[0] = pixman_fixed_1; - unit.vector[1] = 0; - unit.vector[2] = 0; - } - - if (unit.vector[2] == 0 && v.vector[2] == pixman_fixed_1) - { - /* - * Given: - * - * t = (B ± ⎷(B² - A·C)) / A - * - * where - * - * A = cdx² + cdy² - dr² - * B = pdx·cdx + pdy·cdy + r₁·dr - * C = pdx² + pdy² - r₁² - * det = B² - A·C - * - * Since we have an affine transformation, we know that (pdx, pdy) - * increase linearly with each pixel, - * - * pdx = pdx₀ + n·ux, - * pdy = pdy₀ + n·uy, - * - * we can then express B, C and det through multiple differentiation. - */ - pixman_fixed_32_32_t b, db, c, dc, ddc; - - /* warning: this computation may overflow */ - v.vector[0] -= radial->c1.x; - v.vector[1] -= radial->c1.y; - - /* - * B and C are computed and updated exactly. - * If fdot was used instead of dot, in the worst case it would - * lose 11 bits of precision in each of the multiplication and - * summing up would zero out all the bit that were preserved, - * thus making the result 0 instead of the correct one. - * This would mean a worst case of unbound relative error or - * about 2^10 absolute error - */ - b = dot (v.vector[0], v.vector[1], radial->c1.radius, - radial->delta.x, radial->delta.y, radial->delta.radius); - db = dot (unit.vector[0], unit.vector[1], 0, - radial->delta.x, radial->delta.y, 0); - - c = dot (v.vector[0], v.vector[1], - -((pixman_fixed_48_16_t) radial->c1.radius), - v.vector[0], v.vector[1], radial->c1.radius); - dc = dot (2 * (pixman_fixed_48_16_t) v.vector[0] + unit.vector[0], - 2 * (pixman_fixed_48_16_t) v.vector[1] + unit.vector[1], - 0, - unit.vector[0], unit.vector[1], 0); - ddc = 2 * dot (unit.vector[0], unit.vector[1], 0, - unit.vector[0], unit.vector[1], 0); - - while (buffer < end) - { - if (!mask || *mask++) - { - *buffer = radial_compute_color (radial->a, b, c, - radial->inva, - radial->delta.radius, - radial->mindr, - &walker, - image->common.repeat); - } - - b += db; - c += dc; - dc += ddc; - ++buffer; - } - } - else - { - /* projective */ - /* Warning: - * error propagation guarantees are much looser than in the affine case - */ - while (buffer < end) - { - if (!mask || *mask++) - { - if (v.vector[2] != 0) - { - double pdx, pdy, invv2, b, c; - - invv2 = 1. * pixman_fixed_1 / v.vector[2]; - - pdx = v.vector[0] * invv2 - radial->c1.x; - /* / pixman_fixed_1 */ - - pdy = v.vector[1] * invv2 - radial->c1.y; - /* / pixman_fixed_1 */ - - b = fdot (pdx, pdy, radial->c1.radius, - radial->delta.x, radial->delta.y, - radial->delta.radius); - /* / pixman_fixed_1 / pixman_fixed_1 */ - - c = fdot (pdx, pdy, -radial->c1.radius, - pdx, pdy, radial->c1.radius); - /* / pixman_fixed_1 / pixman_fixed_1 */ - - *buffer = radial_compute_color (radial->a, b, c, - radial->inva, - radial->delta.radius, - radial->mindr, - &walker, - image->common.repeat); - } - else - { - *buffer = 0; - } - } - - ++buffer; - - v.vector[0] += unit.vector[0]; - v.vector[1] += unit.vector[1]; - v.vector[2] += unit.vector[2]; - } - } - - iter->y++; - return iter->buffer; -} - -static uint32_t * -radial_get_scanline_wide (pixman_iter_t *iter, const uint32_t *mask) -{ - uint32_t *buffer = radial_get_scanline_narrow (iter, NULL); - - pixman_expand ((uint64_t *)buffer, buffer, PIXMAN_a8r8g8b8, iter->width); - - return buffer; -} - -void -_pixman_radial_gradient_iter_init (pixman_image_t *image, - pixman_iter_t *iter, - int x, int y, int width, int height, - uint8_t *buffer, iter_flags_t flags) -{ - if (flags & ITER_NARROW) - iter->get_scanline = radial_get_scanline_narrow; - else - iter->get_scanline = radial_get_scanline_wide; -} - -PIXMAN_EXPORT pixman_image_t * -pixman_image_create_radial_gradient (pixman_point_fixed_t * inner, - pixman_point_fixed_t * outer, - pixman_fixed_t inner_radius, - pixman_fixed_t outer_radius, - const pixman_gradient_stop_t *stops, - int n_stops) -{ - pixman_image_t *image; - radial_gradient_t *radial; - - image = _pixman_image_allocate (); - - if (!image) - return NULL; - - radial = &image->radial; - - if (!_pixman_init_gradient (&radial->common, stops, n_stops)) - { - free (image); - return NULL; - } - - image->type = RADIAL; - - radial->c1.x = inner->x; - radial->c1.y = inner->y; - radial->c1.radius = inner_radius; - radial->c2.x = outer->x; - radial->c2.y = outer->y; - radial->c2.radius = outer_radius; - - /* warning: this computations may overflow */ - radial->delta.x = radial->c2.x - radial->c1.x; - radial->delta.y = radial->c2.y - radial->c1.y; - radial->delta.radius = radial->c2.radius - radial->c1.radius; - - /* computed exactly, then cast to double -> every bit of the double - representation is correct (53 bits) */ - radial->a = dot (radial->delta.x, radial->delta.y, -radial->delta.radius, - radial->delta.x, radial->delta.y, radial->delta.radius); - if (radial->a != 0) - radial->inva = 1. * pixman_fixed_1 / radial->a; - - radial->mindr = -1. * pixman_fixed_1 * radial->c1.radius; - - return image; -} - +/* -*- Mode: c; c-basic-offset: 4; tab-width: 8; indent-tabs-mode: t; -*- */
+/*
+ *
+ * Copyright © 2000 Keith Packard, member of The XFree86 Project, Inc.
+ * Copyright © 2000 SuSE, Inc.
+ * 2005 Lars Knoll & Zack Rusin, Trolltech
+ * Copyright © 2007 Red Hat, Inc.
+ *
+ *
+ * Permission to use, copy, modify, distribute, and sell this software and its
+ * documentation for any purpose is hereby granted without fee, provided that
+ * the above copyright notice appear in all copies and that both that
+ * copyright notice and this permission notice appear in supporting
+ * documentation, and that the name of Keith Packard not be used in
+ * advertising or publicity pertaining to distribution of the software without
+ * specific, written prior permission. Keith Packard makes no
+ * representations about the suitability of this software for any purpose. It
+ * is provided "as is" without express or implied warranty.
+ *
+ * THE COPYRIGHT HOLDERS DISCLAIM ALL WARRANTIES WITH REGARD TO THIS
+ * SOFTWARE, INCLUDING ALL IMPLIED WARRANTIES OF MERCHANTABILITY AND
+ * FITNESS, IN NO EVENT SHALL THE COPYRIGHT HOLDERS BE LIABLE FOR ANY
+ * SPECIAL, INDIRECT OR CONSEQUENTIAL DAMAGES OR ANY DAMAGES
+ * WHATSOEVER RESULTING FROM LOSS OF USE, DATA OR PROFITS, WHETHER IN
+ * AN ACTION OF CONTRACT, NEGLIGENCE OR OTHER TORTIOUS ACTION, ARISING
+ * OUT OF OR IN CONNECTION WITH THE USE OR PERFORMANCE OF THIS
+ * SOFTWARE.
+ */
+
+#ifdef HAVE_CONFIG_H
+#include <config.h>
+#endif
+#include <stdlib.h>
+#include <math.h>
+#include "pixman-private.h"
+
+static inline pixman_fixed_32_32_t
+dot (pixman_fixed_48_16_t x1,
+ pixman_fixed_48_16_t y1,
+ pixman_fixed_48_16_t z1,
+ pixman_fixed_48_16_t x2,
+ pixman_fixed_48_16_t y2,
+ pixman_fixed_48_16_t z2)
+{
+ /*
+ * Exact computation, assuming that the input values can
+ * be represented as pixman_fixed_16_16_t
+ */
+ return x1 * x2 + y1 * y2 + z1 * z2;
+}
+
+static inline double
+fdot (double x1,
+ double y1,
+ double z1,
+ double x2,
+ double y2,
+ double z2)
+{
+ /*
+ * Error can be unbound in some special cases.
+ * Using clever dot product algorithms (for example compensated
+ * dot product) would improve this but make the code much less
+ * obvious
+ */
+ return x1 * x2 + y1 * y2 + z1 * z2;
+}
+
+static uint32_t
+radial_compute_color (double a,
+ double b,
+ double c,
+ double inva,
+ double dr,
+ double mindr,
+ pixman_gradient_walker_t *walker,
+ pixman_repeat_t repeat)
+{
+ /*
+ * In this function error propagation can lead to bad results:
+ * - det can have an unbound error (if b*b-a*c is very small),
+ * potentially making it the opposite sign of what it should have been
+ * (thus clearing a pixel that would have been colored or vice-versa)
+ * or propagating the error to sqrtdet;
+ * if det has the wrong sign or b is very small, this can lead to bad
+ * results
+ *
+ * - the algorithm used to compute the solutions of the quadratic
+ * equation is not numerically stable (but saves one division compared
+ * to the numerically stable one);
+ * this can be a problem if a*c is much smaller than b*b
+ *
+ * - the above problems are worse if a is small (as inva becomes bigger)
+ */
+ double det;
+
+ if (a == 0)
+ {
+ double t;
+
+ if (b == 0)
+ return 0;
+
+ t = pixman_fixed_1 / 2 * c / b;
+ if (repeat == PIXMAN_REPEAT_NONE)
+ {
+ if (0 <= t && t <= pixman_fixed_1)
+ return _pixman_gradient_walker_pixel (walker, t);
+ }
+ else
+ {
+ if (t * dr > mindr)
+ return _pixman_gradient_walker_pixel (walker, t);
+ }
+
+ return 0;
+ }
+
+ det = fdot (b, a, 0, b, -c, 0);
+ if (det >= 0)
+ {
+ double sqrtdet, t0, t1;
+
+ sqrtdet = sqrt (det);
+ t0 = (b + sqrtdet) * inva;
+ t1 = (b - sqrtdet) * inva;
+
+ if (repeat == PIXMAN_REPEAT_NONE)
+ {
+ if (0 <= t0 && t0 <= pixman_fixed_1)
+ return _pixman_gradient_walker_pixel (walker, t0);
+ else if (0 <= t1 && t1 <= pixman_fixed_1)
+ return _pixman_gradient_walker_pixel (walker, t1);
+ }
+ else
+ {
+ if (t0 * dr > mindr)
+ return _pixman_gradient_walker_pixel (walker, t0);
+ else if (t1 * dr > mindr)
+ return _pixman_gradient_walker_pixel (walker, t1);
+ }
+ }
+
+ return 0;
+}
+
+static uint32_t *
+radial_get_scanline_narrow (pixman_iter_t *iter, const uint32_t *mask)
+{
+ /*
+ * Implementation of radial gradients following the PDF specification.
+ * See section 8.7.4.5.4 Type 3 (Radial) Shadings of the PDF Reference
+ * Manual (PDF 32000-1:2008 at the time of this writing).
+ *
+ * In the radial gradient problem we are given two circles (c₁,r₁) and
+ * (c₂,r₂) that define the gradient itself.
+ *
+ * Mathematically the gradient can be defined as the family of circles
+ *
+ * ((1-t)·c₁ + t·(c₂), (1-t)·r₁ + t·r₂)
+ *
+ * excluding those circles whose radius would be < 0. When a point
+ * belongs to more than one circle, the one with a bigger t is the only
+ * one that contributes to its color. When a point does not belong
+ * to any of the circles, it is transparent black, i.e. RGBA (0, 0, 0, 0).
+ * Further limitations on the range of values for t are imposed when
+ * the gradient is not repeated, namely t must belong to [0,1].
+ *
+ * The graphical result is the same as drawing the valid (radius > 0)
+ * circles with increasing t in [-inf, +inf] (or in [0,1] if the gradient
+ * is not repeated) using SOURCE operatior composition.
+ *
+ * It looks like a cone pointing towards the viewer if the ending circle
+ * is smaller than the starting one, a cone pointing inside the page if
+ * the starting circle is the smaller one and like a cylinder if they
+ * have the same radius.
+ *
+ * What we actually do is, given the point whose color we are interested
+ * in, compute the t values for that point, solving for t in:
+ *
+ * length((1-t)·c₁ + t·(c₂) - p) = (1-t)·r₁ + t·r₂
+ *
+ * Let's rewrite it in a simpler way, by defining some auxiliary
+ * variables:
+ *
+ * cd = c₂ - c₁
+ * pd = p - c₁
+ * dr = r₂ - r₁
+ * lenght(t·cd - pd) = r₁ + t·dr
+ *
+ * which actually means
+ *
+ * hypot(t·cdx - pdx, t·cdy - pdy) = r₁ + t·dr
+ *
+ * or
+ *
+ * ⎷((t·cdx - pdx)² + (t·cdy - pdy)²) = r₁ + t·dr.
+ *
+ * If we impose (as stated earlier) that r₁ + t·dr >= 0, it becomes:
+ *
+ * (t·cdx - pdx)² + (t·cdy - pdy)² = (r₁ + t·dr)²
+ *
+ * where we can actually expand the squares and solve for t:
+ *
+ * t²cdx² - 2t·cdx·pdx + pdx² + t²cdy² - 2t·cdy·pdy + pdy² =
+ * = r₁² + 2·r₁·t·dr + t²·dr²
+ *
+ * (cdx² + cdy² - dr²)t² - 2(cdx·pdx + cdy·pdy + r₁·dr)t +
+ * (pdx² + pdy² - r₁²) = 0
+ *
+ * A = cdx² + cdy² - dr²
+ * B = pdx·cdx + pdy·cdy + r₁·dr
+ * C = pdx² + pdy² - r₁²
+ * At² - 2Bt + C = 0
+ *
+ * The solutions (unless the equation degenerates because of A = 0) are:
+ *
+ * t = (B ± ⎷(B² - A·C)) / A
+ *
+ * The solution we are going to prefer is the bigger one, unless the
+ * radius associated to it is negative (or it falls outside the valid t
+ * range).
+ *
+ * Additional observations (useful for optimizations):
+ * A does not depend on p
+ *
+ * A < 0 <=> one of the two circles completely contains the other one
+ * <=> for every p, the radiuses associated with the two t solutions
+ * have opposite sign
+ */
+ pixman_image_t *image = iter->image;
+ int x = iter->x;
+ int y = iter->y;
+ int width = iter->width;
+ uint32_t *buffer = iter->buffer;
+
+ gradient_t *gradient = (gradient_t *)image;
+ radial_gradient_t *radial = (radial_gradient_t *)image;
+ uint32_t *end = buffer + width;
+ pixman_gradient_walker_t walker;
+ pixman_vector_t v, unit;
+
+ /* reference point is the center of the pixel */
+ v.vector[0] = pixman_int_to_fixed (x) + pixman_fixed_1 / 2;
+ v.vector[1] = pixman_int_to_fixed (y) + pixman_fixed_1 / 2;
+ v.vector[2] = pixman_fixed_1;
+
+ _pixman_gradient_walker_init (&walker, gradient, image->common.repeat);
+
+ if (image->common.transform)
+ {
+ if (!pixman_transform_point_3d (image->common.transform, &v))
+ return iter->buffer;
+
+ unit.vector[0] = image->common.transform->matrix[0][0];
+ unit.vector[1] = image->common.transform->matrix[1][0];
+ unit.vector[2] = image->common.transform->matrix[2][0];
+ }
+ else
+ {
+ unit.vector[0] = pixman_fixed_1;
+ unit.vector[1] = 0;
+ unit.vector[2] = 0;
+ }
+
+ if (unit.vector[2] == 0 && v.vector[2] == pixman_fixed_1)
+ {
+ /*
+ * Given:
+ *
+ * t = (B ± ⎷(B² - A·C)) / A
+ *
+ * where
+ *
+ * A = cdx² + cdy² - dr²
+ * B = pdx·cdx + pdy·cdy + r₁·dr
+ * C = pdx² + pdy² - r₁²
+ * det = B² - A·C
+ *
+ * Since we have an affine transformation, we know that (pdx, pdy)
+ * increase linearly with each pixel,
+ *
+ * pdx = pdx₀ + n·ux,
+ * pdy = pdy₀ + n·uy,
+ *
+ * we can then express B, C and det through multiple differentiation.
+ */
+ pixman_fixed_32_32_t b, db, c, dc, ddc;
+
+ /* warning: this computation may overflow */
+ v.vector[0] -= radial->c1.x;
+ v.vector[1] -= radial->c1.y;
+
+ /*
+ * B and C are computed and updated exactly.
+ * If fdot was used instead of dot, in the worst case it would
+ * lose 11 bits of precision in each of the multiplication and
+ * summing up would zero out all the bit that were preserved,
+ * thus making the result 0 instead of the correct one.
+ * This would mean a worst case of unbound relative error or
+ * about 2^10 absolute error
+ */
+ b = dot (v.vector[0], v.vector[1], radial->c1.radius,
+ radial->delta.x, radial->delta.y, radial->delta.radius);
+ db = dot (unit.vector[0], unit.vector[1], 0,
+ radial->delta.x, radial->delta.y, 0);
+
+ c = dot (v.vector[0], v.vector[1],
+ -((pixman_fixed_48_16_t) radial->c1.radius),
+ v.vector[0], v.vector[1], radial->c1.radius);
+ dc = dot (2 * (pixman_fixed_48_16_t) v.vector[0] + unit.vector[0],
+ 2 * (pixman_fixed_48_16_t) v.vector[1] + unit.vector[1],
+ 0,
+ unit.vector[0], unit.vector[1], 0);
+ ddc = 2 * dot (unit.vector[0], unit.vector[1], 0,
+ unit.vector[0], unit.vector[1], 0);
+
+ while (buffer < end)
+ {
+ if (!mask || *mask++)
+ {
+ *buffer = radial_compute_color (radial->a, b, c,
+ radial->inva,
+ radial->delta.radius,
+ radial->mindr,
+ &walker,
+ image->common.repeat);
+ }
+
+ b += db;
+ c += dc;
+ dc += ddc;
+ ++buffer;
+ }
+ }
+ else
+ {
+ /* projective */
+ /* Warning:
+ * error propagation guarantees are much looser than in the affine case
+ */
+ while (buffer < end)
+ {
+ if (!mask || *mask++)
+ {
+ if (v.vector[2] != 0)
+ {
+ double pdx, pdy, invv2, b, c;
+
+ invv2 = 1. * pixman_fixed_1 / v.vector[2];
+
+ pdx = v.vector[0] * invv2 - radial->c1.x;
+ /* / pixman_fixed_1 */
+
+ pdy = v.vector[1] * invv2 - radial->c1.y;
+ /* / pixman_fixed_1 */
+
+ b = fdot (pdx, pdy, radial->c1.radius,
+ radial->delta.x, radial->delta.y,
+ radial->delta.radius);
+ /* / pixman_fixed_1 / pixman_fixed_1 */
+
+ c = fdot (pdx, pdy, -radial->c1.radius,
+ pdx, pdy, radial->c1.radius);
+ /* / pixman_fixed_1 / pixman_fixed_1 */
+
+ *buffer = radial_compute_color (radial->a, b, c,
+ radial->inva,
+ radial->delta.radius,
+ radial->mindr,
+ &walker,
+ image->common.repeat);
+ }
+ else
+ {
+ *buffer = 0;
+ }
+ }
+
+ ++buffer;
+
+ v.vector[0] += unit.vector[0];
+ v.vector[1] += unit.vector[1];
+ v.vector[2] += unit.vector[2];
+ }
+ }
+
+ iter->y++;
+ return iter->buffer;
+}
+
+static uint32_t *
+radial_get_scanline_wide (pixman_iter_t *iter, const uint32_t *mask)
+{
+ uint32_t *buffer = radial_get_scanline_narrow (iter, NULL);
+
+ pixman_expand ((uint64_t *)buffer, buffer, PIXMAN_a8r8g8b8, iter->width);
+
+ return buffer;
+}
+
+void
+_pixman_radial_gradient_iter_init (pixman_image_t *image, pixman_iter_t *iter)
+{
+ if (iter->flags & ITER_NARROW)
+ iter->get_scanline = radial_get_scanline_narrow;
+ else
+ iter->get_scanline = radial_get_scanline_wide;
+}
+
+PIXMAN_EXPORT pixman_image_t *
+pixman_image_create_radial_gradient (pixman_point_fixed_t * inner,
+ pixman_point_fixed_t * outer,
+ pixman_fixed_t inner_radius,
+ pixman_fixed_t outer_radius,
+ const pixman_gradient_stop_t *stops,
+ int n_stops)
+{
+ pixman_image_t *image;
+ radial_gradient_t *radial;
+
+ image = _pixman_image_allocate ();
+
+ if (!image)
+ return NULL;
+
+ radial = &image->radial;
+
+ if (!_pixman_init_gradient (&radial->common, stops, n_stops))
+ {
+ free (image);
+ return NULL;
+ }
+
+ image->type = RADIAL;
+
+ radial->c1.x = inner->x;
+ radial->c1.y = inner->y;
+ radial->c1.radius = inner_radius;
+ radial->c2.x = outer->x;
+ radial->c2.y = outer->y;
+ radial->c2.radius = outer_radius;
+
+ /* warning: this computations may overflow */
+ radial->delta.x = radial->c2.x - radial->c1.x;
+ radial->delta.y = radial->c2.y - radial->c1.y;
+ radial->delta.radius = radial->c2.radius - radial->c1.radius;
+
+ /* computed exactly, then cast to double -> every bit of the double
+ representation is correct (53 bits) */
+ radial->a = dot (radial->delta.x, radial->delta.y, -radial->delta.radius,
+ radial->delta.x, radial->delta.y, radial->delta.radius);
+ if (radial->a != 0)
+ radial->inva = 1. * pixman_fixed_1 / radial->a;
+
+ radial->mindr = -1. * pixman_fixed_1 * radial->c1.radius;
+
+ return image;
+}
+
|