gfx/cairo/libpixman/src/pixman-radial-gradient.c

Wed, 31 Dec 2014 06:09:35 +0100

author
Michael Schloh von Bennewitz <michael@schloh.com>
date
Wed, 31 Dec 2014 06:09:35 +0100
changeset 0
6474c204b198
permissions
-rw-r--r--

Cloned upstream origin tor-browser at tor-browser-31.3.0esr-4.5-1-build1
revision ID fc1c9ff7c1b2defdbc039f12214767608f46423f for hacking purpose.

michael@0 1 /* -*- Mode: c; c-basic-offset: 4; tab-width: 8; indent-tabs-mode: t; -*- */
michael@0 2 /*
michael@0 3 *
michael@0 4 * Copyright © 2000 Keith Packard, member of The XFree86 Project, Inc.
michael@0 5 * Copyright © 2000 SuSE, Inc.
michael@0 6 * 2005 Lars Knoll & Zack Rusin, Trolltech
michael@0 7 * Copyright © 2007 Red Hat, Inc.
michael@0 8 *
michael@0 9 *
michael@0 10 * Permission to use, copy, modify, distribute, and sell this software and its
michael@0 11 * documentation for any purpose is hereby granted without fee, provided that
michael@0 12 * the above copyright notice appear in all copies and that both that
michael@0 13 * copyright notice and this permission notice appear in supporting
michael@0 14 * documentation, and that the name of Keith Packard not be used in
michael@0 15 * advertising or publicity pertaining to distribution of the software without
michael@0 16 * specific, written prior permission. Keith Packard makes no
michael@0 17 * representations about the suitability of this software for any purpose. It
michael@0 18 * is provided "as is" without express or implied warranty.
michael@0 19 *
michael@0 20 * THE COPYRIGHT HOLDERS DISCLAIM ALL WARRANTIES WITH REGARD TO THIS
michael@0 21 * SOFTWARE, INCLUDING ALL IMPLIED WARRANTIES OF MERCHANTABILITY AND
michael@0 22 * FITNESS, IN NO EVENT SHALL THE COPYRIGHT HOLDERS BE LIABLE FOR ANY
michael@0 23 * SPECIAL, INDIRECT OR CONSEQUENTIAL DAMAGES OR ANY DAMAGES
michael@0 24 * WHATSOEVER RESULTING FROM LOSS OF USE, DATA OR PROFITS, WHETHER IN
michael@0 25 * AN ACTION OF CONTRACT, NEGLIGENCE OR OTHER TORTIOUS ACTION, ARISING
michael@0 26 * OUT OF OR IN CONNECTION WITH THE USE OR PERFORMANCE OF THIS
michael@0 27 * SOFTWARE.
michael@0 28 */
michael@0 29
michael@0 30 #ifdef HAVE_CONFIG_H
michael@0 31 #include <config.h>
michael@0 32 #endif
michael@0 33 #include <stdlib.h>
michael@0 34 #include <math.h>
michael@0 35 #include "pixman-private.h"
michael@0 36
michael@0 37 #include "pixman-dither.h"
michael@0 38
michael@0 39 static inline pixman_fixed_32_32_t
michael@0 40 dot (pixman_fixed_48_16_t x1,
michael@0 41 pixman_fixed_48_16_t y1,
michael@0 42 pixman_fixed_48_16_t z1,
michael@0 43 pixman_fixed_48_16_t x2,
michael@0 44 pixman_fixed_48_16_t y2,
michael@0 45 pixman_fixed_48_16_t z2)
michael@0 46 {
michael@0 47 /*
michael@0 48 * Exact computation, assuming that the input values can
michael@0 49 * be represented as pixman_fixed_16_16_t
michael@0 50 */
michael@0 51 return x1 * x2 + y1 * y2 + z1 * z2;
michael@0 52 }
michael@0 53
michael@0 54 static inline double
michael@0 55 fdot (double x1,
michael@0 56 double y1,
michael@0 57 double z1,
michael@0 58 double x2,
michael@0 59 double y2,
michael@0 60 double z2)
michael@0 61 {
michael@0 62 /*
michael@0 63 * Error can be unbound in some special cases.
michael@0 64 * Using clever dot product algorithms (for example compensated
michael@0 65 * dot product) would improve this but make the code much less
michael@0 66 * obvious
michael@0 67 */
michael@0 68 return x1 * x2 + y1 * y2 + z1 * z2;
michael@0 69 }
michael@0 70
michael@0 71 static uint32_t
michael@0 72 radial_compute_color (double a,
michael@0 73 double b,
michael@0 74 double c,
michael@0 75 double inva,
michael@0 76 double dr,
michael@0 77 double mindr,
michael@0 78 pixman_gradient_walker_t *walker,
michael@0 79 pixman_repeat_t repeat)
michael@0 80 {
michael@0 81 /*
michael@0 82 * In this function error propagation can lead to bad results:
michael@0 83 * - discr can have an unbound error (if b*b-a*c is very small),
michael@0 84 * potentially making it the opposite sign of what it should have been
michael@0 85 * (thus clearing a pixel that would have been colored or vice-versa)
michael@0 86 * or propagating the error to sqrtdiscr;
michael@0 87 * if discr has the wrong sign or b is very small, this can lead to bad
michael@0 88 * results
michael@0 89 *
michael@0 90 * - the algorithm used to compute the solutions of the quadratic
michael@0 91 * equation is not numerically stable (but saves one division compared
michael@0 92 * to the numerically stable one);
michael@0 93 * this can be a problem if a*c is much smaller than b*b
michael@0 94 *
michael@0 95 * - the above problems are worse if a is small (as inva becomes bigger)
michael@0 96 */
michael@0 97 double discr;
michael@0 98
michael@0 99 if (a == 0)
michael@0 100 {
michael@0 101 double t;
michael@0 102
michael@0 103 if (b == 0)
michael@0 104 return 0;
michael@0 105
michael@0 106 t = pixman_fixed_1 / 2 * c / b;
michael@0 107 if (repeat == PIXMAN_REPEAT_NONE)
michael@0 108 {
michael@0 109 if (0 <= t && t <= pixman_fixed_1)
michael@0 110 return _pixman_gradient_walker_pixel (walker, t);
michael@0 111 }
michael@0 112 else
michael@0 113 {
michael@0 114 if (t * dr >= mindr)
michael@0 115 return _pixman_gradient_walker_pixel (walker, t);
michael@0 116 }
michael@0 117
michael@0 118 return 0;
michael@0 119 }
michael@0 120
michael@0 121 discr = fdot (b, a, 0, b, -c, 0);
michael@0 122 if (discr >= 0)
michael@0 123 {
michael@0 124 double sqrtdiscr, t0, t1;
michael@0 125
michael@0 126 sqrtdiscr = sqrt (discr);
michael@0 127 t0 = (b + sqrtdiscr) * inva;
michael@0 128 t1 = (b - sqrtdiscr) * inva;
michael@0 129
michael@0 130 /*
michael@0 131 * The root that must be used is the biggest one that belongs
michael@0 132 * to the valid range ([0,1] for PIXMAN_REPEAT_NONE, any
michael@0 133 * solution that results in a positive radius otherwise).
michael@0 134 *
michael@0 135 * If a > 0, t0 is the biggest solution, so if it is valid, it
michael@0 136 * is the correct result.
michael@0 137 *
michael@0 138 * If a < 0, only one of the solutions can be valid, so the
michael@0 139 * order in which they are tested is not important.
michael@0 140 */
michael@0 141 if (repeat == PIXMAN_REPEAT_NONE)
michael@0 142 {
michael@0 143 if (0 <= t0 && t0 <= pixman_fixed_1)
michael@0 144 return _pixman_gradient_walker_pixel (walker, t0);
michael@0 145 else if (0 <= t1 && t1 <= pixman_fixed_1)
michael@0 146 return _pixman_gradient_walker_pixel (walker, t1);
michael@0 147 }
michael@0 148 else
michael@0 149 {
michael@0 150 if (t0 * dr >= mindr)
michael@0 151 return _pixman_gradient_walker_pixel (walker, t0);
michael@0 152 else if (t1 * dr >= mindr)
michael@0 153 return _pixman_gradient_walker_pixel (walker, t1);
michael@0 154 }
michael@0 155 }
michael@0 156
michael@0 157 return 0;
michael@0 158 }
michael@0 159
michael@0 160 static uint32_t *
michael@0 161 radial_get_scanline_narrow (pixman_iter_t *iter, const uint32_t *mask)
michael@0 162 {
michael@0 163 /*
michael@0 164 * Implementation of radial gradients following the PDF specification.
michael@0 165 * See section 8.7.4.5.4 Type 3 (Radial) Shadings of the PDF Reference
michael@0 166 * Manual (PDF 32000-1:2008 at the time of this writing).
michael@0 167 *
michael@0 168 * In the radial gradient problem we are given two circles (c₁,r₁) and
michael@0 169 * (c₂,r₂) that define the gradient itself.
michael@0 170 *
michael@0 171 * Mathematically the gradient can be defined as the family of circles
michael@0 172 *
michael@0 173 * ((1-t)·c₁ + t·(c₂), (1-t)·r₁ + t·r₂)
michael@0 174 *
michael@0 175 * excluding those circles whose radius would be < 0. When a point
michael@0 176 * belongs to more than one circle, the one with a bigger t is the only
michael@0 177 * one that contributes to its color. When a point does not belong
michael@0 178 * to any of the circles, it is transparent black, i.e. RGBA (0, 0, 0, 0).
michael@0 179 * Further limitations on the range of values for t are imposed when
michael@0 180 * the gradient is not repeated, namely t must belong to [0,1].
michael@0 181 *
michael@0 182 * The graphical result is the same as drawing the valid (radius > 0)
michael@0 183 * circles with increasing t in [-inf, +inf] (or in [0,1] if the gradient
michael@0 184 * is not repeated) using SOURCE operator composition.
michael@0 185 *
michael@0 186 * It looks like a cone pointing towards the viewer if the ending circle
michael@0 187 * is smaller than the starting one, a cone pointing inside the page if
michael@0 188 * the starting circle is the smaller one and like a cylinder if they
michael@0 189 * have the same radius.
michael@0 190 *
michael@0 191 * What we actually do is, given the point whose color we are interested
michael@0 192 * in, compute the t values for that point, solving for t in:
michael@0 193 *
michael@0 194 * length((1-t)·c₁ + t·(c₂) - p) = (1-t)·r₁ + t·r₂
michael@0 195 *
michael@0 196 * Let's rewrite it in a simpler way, by defining some auxiliary
michael@0 197 * variables:
michael@0 198 *
michael@0 199 * cd = c₂ - c₁
michael@0 200 * pd = p - c₁
michael@0 201 * dr = r₂ - r₁
michael@0 202 * length(t·cd - pd) = r₁ + t·dr
michael@0 203 *
michael@0 204 * which actually means
michael@0 205 *
michael@0 206 * hypot(t·cdx - pdx, t·cdy - pdy) = r₁ + t·dr
michael@0 207 *
michael@0 208 * or
michael@0 209 *
michael@0 210 * ⎷((t·cdx - pdx)² + (t·cdy - pdy)²) = r₁ + t·dr.
michael@0 211 *
michael@0 212 * If we impose (as stated earlier) that r₁ + t·dr >= 0, it becomes:
michael@0 213 *
michael@0 214 * (t·cdx - pdx)² + (t·cdy - pdy)² = (r₁ + t·dr)²
michael@0 215 *
michael@0 216 * where we can actually expand the squares and solve for t:
michael@0 217 *
michael@0 218 * t²cdx² - 2t·cdx·pdx + pdx² + t²cdy² - 2t·cdy·pdy + pdy² =
michael@0 219 * = r₁² + 2·r₁·t·dr + t²·dr²
michael@0 220 *
michael@0 221 * (cdx² + cdy² - dr²)t² - 2(cdx·pdx + cdy·pdy + r₁·dr)t +
michael@0 222 * (pdx² + pdy² - r₁²) = 0
michael@0 223 *
michael@0 224 * A = cdx² + cdy² - dr²
michael@0 225 * B = pdx·cdx + pdy·cdy + r₁·dr
michael@0 226 * C = pdx² + pdy² - r₁²
michael@0 227 * At² - 2Bt + C = 0
michael@0 228 *
michael@0 229 * The solutions (unless the equation degenerates because of A = 0) are:
michael@0 230 *
michael@0 231 * t = (B ± ⎷(B² - A·C)) / A
michael@0 232 *
michael@0 233 * The solution we are going to prefer is the bigger one, unless the
michael@0 234 * radius associated to it is negative (or it falls outside the valid t
michael@0 235 * range).
michael@0 236 *
michael@0 237 * Additional observations (useful for optimizations):
michael@0 238 * A does not depend on p
michael@0 239 *
michael@0 240 * A < 0 <=> one of the two circles completely contains the other one
michael@0 241 * <=> for every p, the radiuses associated with the two t solutions
michael@0 242 * have opposite sign
michael@0 243 */
michael@0 244 pixman_image_t *image = iter->image;
michael@0 245 int x = iter->x;
michael@0 246 int y = iter->y;
michael@0 247 int width = iter->width;
michael@0 248 uint32_t *buffer = iter->buffer;
michael@0 249
michael@0 250 gradient_t *gradient = (gradient_t *)image;
michael@0 251 radial_gradient_t *radial = (radial_gradient_t *)image;
michael@0 252 uint32_t *end = buffer + width;
michael@0 253 pixman_gradient_walker_t walker;
michael@0 254 pixman_vector_t v, unit;
michael@0 255
michael@0 256 /* reference point is the center of the pixel */
michael@0 257 v.vector[0] = pixman_int_to_fixed (x) + pixman_fixed_1 / 2;
michael@0 258 v.vector[1] = pixman_int_to_fixed (y) + pixman_fixed_1 / 2;
michael@0 259 v.vector[2] = pixman_fixed_1;
michael@0 260
michael@0 261 _pixman_gradient_walker_init (&walker, gradient, image->common.repeat);
michael@0 262
michael@0 263 if (image->common.transform)
michael@0 264 {
michael@0 265 if (!pixman_transform_point_3d (image->common.transform, &v))
michael@0 266 return iter->buffer;
michael@0 267
michael@0 268 unit.vector[0] = image->common.transform->matrix[0][0];
michael@0 269 unit.vector[1] = image->common.transform->matrix[1][0];
michael@0 270 unit.vector[2] = image->common.transform->matrix[2][0];
michael@0 271 }
michael@0 272 else
michael@0 273 {
michael@0 274 unit.vector[0] = pixman_fixed_1;
michael@0 275 unit.vector[1] = 0;
michael@0 276 unit.vector[2] = 0;
michael@0 277 }
michael@0 278
michael@0 279 if (unit.vector[2] == 0 && v.vector[2] == pixman_fixed_1)
michael@0 280 {
michael@0 281 /*
michael@0 282 * Given:
michael@0 283 *
michael@0 284 * t = (B ± ⎷(B² - A·C)) / A
michael@0 285 *
michael@0 286 * where
michael@0 287 *
michael@0 288 * A = cdx² + cdy² - dr²
michael@0 289 * B = pdx·cdx + pdy·cdy + r₁·dr
michael@0 290 * C = pdx² + pdy² - r₁²
michael@0 291 * det = B² - A·C
michael@0 292 *
michael@0 293 * Since we have an affine transformation, we know that (pdx, pdy)
michael@0 294 * increase linearly with each pixel,
michael@0 295 *
michael@0 296 * pdx = pdx₀ + n·ux,
michael@0 297 * pdy = pdy₀ + n·uy,
michael@0 298 *
michael@0 299 * we can then express B, C and det through multiple differentiation.
michael@0 300 */
michael@0 301 pixman_fixed_32_32_t b, db, c, dc, ddc;
michael@0 302
michael@0 303 /* warning: this computation may overflow */
michael@0 304 v.vector[0] -= radial->c1.x;
michael@0 305 v.vector[1] -= radial->c1.y;
michael@0 306
michael@0 307 /*
michael@0 308 * B and C are computed and updated exactly.
michael@0 309 * If fdot was used instead of dot, in the worst case it would
michael@0 310 * lose 11 bits of precision in each of the multiplication and
michael@0 311 * summing up would zero out all the bit that were preserved,
michael@0 312 * thus making the result 0 instead of the correct one.
michael@0 313 * This would mean a worst case of unbound relative error or
michael@0 314 * about 2^10 absolute error
michael@0 315 */
michael@0 316 b = dot (v.vector[0], v.vector[1], radial->c1.radius,
michael@0 317 radial->delta.x, radial->delta.y, radial->delta.radius);
michael@0 318 db = dot (unit.vector[0], unit.vector[1], 0,
michael@0 319 radial->delta.x, radial->delta.y, 0);
michael@0 320
michael@0 321 c = dot (v.vector[0], v.vector[1],
michael@0 322 -((pixman_fixed_48_16_t) radial->c1.radius),
michael@0 323 v.vector[0], v.vector[1], radial->c1.radius);
michael@0 324 dc = dot (2 * (pixman_fixed_48_16_t) v.vector[0] + unit.vector[0],
michael@0 325 2 * (pixman_fixed_48_16_t) v.vector[1] + unit.vector[1],
michael@0 326 0,
michael@0 327 unit.vector[0], unit.vector[1], 0);
michael@0 328 ddc = 2 * dot (unit.vector[0], unit.vector[1], 0,
michael@0 329 unit.vector[0], unit.vector[1], 0);
michael@0 330
michael@0 331 while (buffer < end)
michael@0 332 {
michael@0 333 if (!mask || *mask++)
michael@0 334 {
michael@0 335 *buffer = radial_compute_color (radial->a, b, c,
michael@0 336 radial->inva,
michael@0 337 radial->delta.radius,
michael@0 338 radial->mindr,
michael@0 339 &walker,
michael@0 340 image->common.repeat);
michael@0 341 }
michael@0 342
michael@0 343 b += db;
michael@0 344 c += dc;
michael@0 345 dc += ddc;
michael@0 346 ++buffer;
michael@0 347 }
michael@0 348 }
michael@0 349 else
michael@0 350 {
michael@0 351 /* projective */
michael@0 352 /* Warning:
michael@0 353 * error propagation guarantees are much looser than in the affine case
michael@0 354 */
michael@0 355 while (buffer < end)
michael@0 356 {
michael@0 357 if (!mask || *mask++)
michael@0 358 {
michael@0 359 if (v.vector[2] != 0)
michael@0 360 {
michael@0 361 double pdx, pdy, invv2, b, c;
michael@0 362
michael@0 363 invv2 = 1. * pixman_fixed_1 / v.vector[2];
michael@0 364
michael@0 365 pdx = v.vector[0] * invv2 - radial->c1.x;
michael@0 366 /* / pixman_fixed_1 */
michael@0 367
michael@0 368 pdy = v.vector[1] * invv2 - radial->c1.y;
michael@0 369 /* / pixman_fixed_1 */
michael@0 370
michael@0 371 b = fdot (pdx, pdy, radial->c1.radius,
michael@0 372 radial->delta.x, radial->delta.y,
michael@0 373 radial->delta.radius);
michael@0 374 /* / pixman_fixed_1 / pixman_fixed_1 */
michael@0 375
michael@0 376 c = fdot (pdx, pdy, -radial->c1.radius,
michael@0 377 pdx, pdy, radial->c1.radius);
michael@0 378 /* / pixman_fixed_1 / pixman_fixed_1 */
michael@0 379
michael@0 380 *buffer = radial_compute_color (radial->a, b, c,
michael@0 381 radial->inva,
michael@0 382 radial->delta.radius,
michael@0 383 radial->mindr,
michael@0 384 &walker,
michael@0 385 image->common.repeat);
michael@0 386 }
michael@0 387 else
michael@0 388 {
michael@0 389 *buffer = 0;
michael@0 390 }
michael@0 391 }
michael@0 392
michael@0 393 ++buffer;
michael@0 394
michael@0 395 v.vector[0] += unit.vector[0];
michael@0 396 v.vector[1] += unit.vector[1];
michael@0 397 v.vector[2] += unit.vector[2];
michael@0 398 }
michael@0 399 }
michael@0 400
michael@0 401 iter->y++;
michael@0 402 return iter->buffer;
michael@0 403 }
michael@0 404
michael@0 405 static uint32_t *
michael@0 406 radial_get_scanline_16 (pixman_iter_t *iter, const uint32_t *mask)
michael@0 407 {
michael@0 408 /*
michael@0 409 * Implementation of radial gradients following the PDF specification.
michael@0 410 * See section 8.7.4.5.4 Type 3 (Radial) Shadings of the PDF Reference
michael@0 411 * Manual (PDF 32000-1:2008 at the time of this writing).
michael@0 412 *
michael@0 413 * In the radial gradient problem we are given two circles (c₁,r₁) and
michael@0 414 * (c₂,r₂) that define the gradient itself.
michael@0 415 *
michael@0 416 * Mathematically the gradient can be defined as the family of circles
michael@0 417 *
michael@0 418 * ((1-t)·c₁ + t·(c₂), (1-t)·r₁ + t·r₂)
michael@0 419 *
michael@0 420 * excluding those circles whose radius would be < 0. When a point
michael@0 421 * belongs to more than one circle, the one with a bigger t is the only
michael@0 422 * one that contributes to its color. When a point does not belong
michael@0 423 * to any of the circles, it is transparent black, i.e. RGBA (0, 0, 0, 0).
michael@0 424 * Further limitations on the range of values for t are imposed when
michael@0 425 * the gradient is not repeated, namely t must belong to [0,1].
michael@0 426 *
michael@0 427 * The graphical result is the same as drawing the valid (radius > 0)
michael@0 428 * circles with increasing t in [-inf, +inf] (or in [0,1] if the gradient
michael@0 429 * is not repeated) using SOURCE operator composition.
michael@0 430 *
michael@0 431 * It looks like a cone pointing towards the viewer if the ending circle
michael@0 432 * is smaller than the starting one, a cone pointing inside the page if
michael@0 433 * the starting circle is the smaller one and like a cylinder if they
michael@0 434 * have the same radius.
michael@0 435 *
michael@0 436 * What we actually do is, given the point whose color we are interested
michael@0 437 * in, compute the t values for that point, solving for t in:
michael@0 438 *
michael@0 439 * length((1-t)·c₁ + t·(c₂) - p) = (1-t)·r₁ + t·r₂
michael@0 440 *
michael@0 441 * Let's rewrite it in a simpler way, by defining some auxiliary
michael@0 442 * variables:
michael@0 443 *
michael@0 444 * cd = c₂ - c₁
michael@0 445 * pd = p - c₁
michael@0 446 * dr = r₂ - r₁
michael@0 447 * length(t·cd - pd) = r₁ + t·dr
michael@0 448 *
michael@0 449 * which actually means
michael@0 450 *
michael@0 451 * hypot(t·cdx - pdx, t·cdy - pdy) = r₁ + t·dr
michael@0 452 *
michael@0 453 * or
michael@0 454 *
michael@0 455 * ⎷((t·cdx - pdx)² + (t·cdy - pdy)²) = r₁ + t·dr.
michael@0 456 *
michael@0 457 * If we impose (as stated earlier) that r₁ + t·dr >= 0, it becomes:
michael@0 458 *
michael@0 459 * (t·cdx - pdx)² + (t·cdy - pdy)² = (r₁ + t·dr)²
michael@0 460 *
michael@0 461 * where we can actually expand the squares and solve for t:
michael@0 462 *
michael@0 463 * t²cdx² - 2t·cdx·pdx + pdx² + t²cdy² - 2t·cdy·pdy + pdy² =
michael@0 464 * = r₁² + 2·r₁·t·dr + t²·dr²
michael@0 465 *
michael@0 466 * (cdx² + cdy² - dr²)t² - 2(cdx·pdx + cdy·pdy + r₁·dr)t +
michael@0 467 * (pdx² + pdy² - r₁²) = 0
michael@0 468 *
michael@0 469 * A = cdx² + cdy² - dr²
michael@0 470 * B = pdx·cdx + pdy·cdy + r₁·dr
michael@0 471 * C = pdx² + pdy² - r₁²
michael@0 472 * At² - 2Bt + C = 0
michael@0 473 *
michael@0 474 * The solutions (unless the equation degenerates because of A = 0) are:
michael@0 475 *
michael@0 476 * t = (B ± ⎷(B² - A·C)) / A
michael@0 477 *
michael@0 478 * The solution we are going to prefer is the bigger one, unless the
michael@0 479 * radius associated to it is negative (or it falls outside the valid t
michael@0 480 * range).
michael@0 481 *
michael@0 482 * Additional observations (useful for optimizations):
michael@0 483 * A does not depend on p
michael@0 484 *
michael@0 485 * A < 0 <=> one of the two circles completely contains the other one
michael@0 486 * <=> for every p, the radiuses associated with the two t solutions
michael@0 487 * have opposite sign
michael@0 488 */
michael@0 489 pixman_image_t *image = iter->image;
michael@0 490 int x = iter->x;
michael@0 491 int y = iter->y;
michael@0 492 int width = iter->width;
michael@0 493 uint16_t *buffer = iter->buffer;
michael@0 494 pixman_bool_t toggle = ((x ^ y) & 1);
michael@0 495
michael@0 496 gradient_t *gradient = (gradient_t *)image;
michael@0 497 radial_gradient_t *radial = (radial_gradient_t *)image;
michael@0 498 uint16_t *end = buffer + width;
michael@0 499 pixman_gradient_walker_t walker;
michael@0 500 pixman_vector_t v, unit;
michael@0 501
michael@0 502 /* reference point is the center of the pixel */
michael@0 503 v.vector[0] = pixman_int_to_fixed (x) + pixman_fixed_1 / 2;
michael@0 504 v.vector[1] = pixman_int_to_fixed (y) + pixman_fixed_1 / 2;
michael@0 505 v.vector[2] = pixman_fixed_1;
michael@0 506
michael@0 507 _pixman_gradient_walker_init (&walker, gradient, image->common.repeat);
michael@0 508
michael@0 509 if (image->common.transform)
michael@0 510 {
michael@0 511 if (!pixman_transform_point_3d (image->common.transform, &v))
michael@0 512 return iter->buffer;
michael@0 513
michael@0 514 unit.vector[0] = image->common.transform->matrix[0][0];
michael@0 515 unit.vector[1] = image->common.transform->matrix[1][0];
michael@0 516 unit.vector[2] = image->common.transform->matrix[2][0];
michael@0 517 }
michael@0 518 else
michael@0 519 {
michael@0 520 unit.vector[0] = pixman_fixed_1;
michael@0 521 unit.vector[1] = 0;
michael@0 522 unit.vector[2] = 0;
michael@0 523 }
michael@0 524
michael@0 525 if (unit.vector[2] == 0 && v.vector[2] == pixman_fixed_1)
michael@0 526 {
michael@0 527 /*
michael@0 528 * Given:
michael@0 529 *
michael@0 530 * t = (B ± ⎷(B² - A·C)) / A
michael@0 531 *
michael@0 532 * where
michael@0 533 *
michael@0 534 * A = cdx² + cdy² - dr²
michael@0 535 * B = pdx·cdx + pdy·cdy + r₁·dr
michael@0 536 * C = pdx² + pdy² - r₁²
michael@0 537 * det = B² - A·C
michael@0 538 *
michael@0 539 * Since we have an affine transformation, we know that (pdx, pdy)
michael@0 540 * increase linearly with each pixel,
michael@0 541 *
michael@0 542 * pdx = pdx₀ + n·ux,
michael@0 543 * pdy = pdy₀ + n·uy,
michael@0 544 *
michael@0 545 * we can then express B, C and det through multiple differentiation.
michael@0 546 */
michael@0 547 pixman_fixed_32_32_t b, db, c, dc, ddc;
michael@0 548
michael@0 549 /* warning: this computation may overflow */
michael@0 550 v.vector[0] -= radial->c1.x;
michael@0 551 v.vector[1] -= radial->c1.y;
michael@0 552
michael@0 553 /*
michael@0 554 * B and C are computed and updated exactly.
michael@0 555 * If fdot was used instead of dot, in the worst case it would
michael@0 556 * lose 11 bits of precision in each of the multiplication and
michael@0 557 * summing up would zero out all the bit that were preserved,
michael@0 558 * thus making the result 0 instead of the correct one.
michael@0 559 * This would mean a worst case of unbound relative error or
michael@0 560 * about 2^10 absolute error
michael@0 561 */
michael@0 562 b = dot (v.vector[0], v.vector[1], radial->c1.radius,
michael@0 563 radial->delta.x, radial->delta.y, radial->delta.radius);
michael@0 564 db = dot (unit.vector[0], unit.vector[1], 0,
michael@0 565 radial->delta.x, radial->delta.y, 0);
michael@0 566
michael@0 567 c = dot (v.vector[0], v.vector[1],
michael@0 568 -((pixman_fixed_48_16_t) radial->c1.radius),
michael@0 569 v.vector[0], v.vector[1], radial->c1.radius);
michael@0 570 dc = dot (2 * (pixman_fixed_48_16_t) v.vector[0] + unit.vector[0],
michael@0 571 2 * (pixman_fixed_48_16_t) v.vector[1] + unit.vector[1],
michael@0 572 0,
michael@0 573 unit.vector[0], unit.vector[1], 0);
michael@0 574 ddc = 2 * dot (unit.vector[0], unit.vector[1], 0,
michael@0 575 unit.vector[0], unit.vector[1], 0);
michael@0 576
michael@0 577 while (buffer < end)
michael@0 578 {
michael@0 579 if (!mask || *mask++)
michael@0 580 {
michael@0 581 *buffer = dither_8888_to_0565(
michael@0 582 radial_compute_color (radial->a, b, c,
michael@0 583 radial->inva,
michael@0 584 radial->delta.radius,
michael@0 585 radial->mindr,
michael@0 586 &walker,
michael@0 587 image->common.repeat),
michael@0 588 toggle);
michael@0 589 }
michael@0 590
michael@0 591 toggle ^= 1;
michael@0 592 b += db;
michael@0 593 c += dc;
michael@0 594 dc += ddc;
michael@0 595 ++buffer;
michael@0 596 }
michael@0 597 }
michael@0 598 else
michael@0 599 {
michael@0 600 /* projective */
michael@0 601 /* Warning:
michael@0 602 * error propagation guarantees are much looser than in the affine case
michael@0 603 */
michael@0 604 while (buffer < end)
michael@0 605 {
michael@0 606 if (!mask || *mask++)
michael@0 607 {
michael@0 608 if (v.vector[2] != 0)
michael@0 609 {
michael@0 610 double pdx, pdy, invv2, b, c;
michael@0 611
michael@0 612 invv2 = 1. * pixman_fixed_1 / v.vector[2];
michael@0 613
michael@0 614 pdx = v.vector[0] * invv2 - radial->c1.x;
michael@0 615 /* / pixman_fixed_1 */
michael@0 616
michael@0 617 pdy = v.vector[1] * invv2 - radial->c1.y;
michael@0 618 /* / pixman_fixed_1 */
michael@0 619
michael@0 620 b = fdot (pdx, pdy, radial->c1.radius,
michael@0 621 radial->delta.x, radial->delta.y,
michael@0 622 radial->delta.radius);
michael@0 623 /* / pixman_fixed_1 / pixman_fixed_1 */
michael@0 624
michael@0 625 c = fdot (pdx, pdy, -radial->c1.radius,
michael@0 626 pdx, pdy, radial->c1.radius);
michael@0 627 /* / pixman_fixed_1 / pixman_fixed_1 */
michael@0 628
michael@0 629 *buffer = dither_8888_to_0565 (
michael@0 630 radial_compute_color (radial->a, b, c,
michael@0 631 radial->inva,
michael@0 632 radial->delta.radius,
michael@0 633 radial->mindr,
michael@0 634 &walker,
michael@0 635 image->common.repeat),
michael@0 636 toggle);
michael@0 637 }
michael@0 638 else
michael@0 639 {
michael@0 640 *buffer = 0;
michael@0 641 }
michael@0 642 }
michael@0 643
michael@0 644 ++buffer;
michael@0 645 toggle ^= 1;
michael@0 646
michael@0 647 v.vector[0] += unit.vector[0];
michael@0 648 v.vector[1] += unit.vector[1];
michael@0 649 v.vector[2] += unit.vector[2];
michael@0 650 }
michael@0 651 }
michael@0 652
michael@0 653 iter->y++;
michael@0 654 return iter->buffer;
michael@0 655 }
michael@0 656 static uint32_t *
michael@0 657 radial_get_scanline_wide (pixman_iter_t *iter, const uint32_t *mask)
michael@0 658 {
michael@0 659 uint32_t *buffer = radial_get_scanline_narrow (iter, NULL);
michael@0 660
michael@0 661 pixman_expand_to_float (
michael@0 662 (argb_t *)buffer, buffer, PIXMAN_a8r8g8b8, iter->width);
michael@0 663
michael@0 664 return buffer;
michael@0 665 }
michael@0 666
michael@0 667 void
michael@0 668 _pixman_radial_gradient_iter_init (pixman_image_t *image, pixman_iter_t *iter)
michael@0 669 {
michael@0 670 if (iter->iter_flags & ITER_16)
michael@0 671 iter->get_scanline = radial_get_scanline_16;
michael@0 672 else if (iter->iter_flags & ITER_NARROW)
michael@0 673 iter->get_scanline = radial_get_scanline_narrow;
michael@0 674 else
michael@0 675 iter->get_scanline = radial_get_scanline_wide;
michael@0 676 }
michael@0 677
michael@0 678
michael@0 679 PIXMAN_EXPORT pixman_image_t *
michael@0 680 pixman_image_create_radial_gradient (const pixman_point_fixed_t * inner,
michael@0 681 const pixman_point_fixed_t * outer,
michael@0 682 pixman_fixed_t inner_radius,
michael@0 683 pixman_fixed_t outer_radius,
michael@0 684 const pixman_gradient_stop_t *stops,
michael@0 685 int n_stops)
michael@0 686 {
michael@0 687 pixman_image_t *image;
michael@0 688 radial_gradient_t *radial;
michael@0 689
michael@0 690 image = _pixman_image_allocate ();
michael@0 691
michael@0 692 if (!image)
michael@0 693 return NULL;
michael@0 694
michael@0 695 radial = &image->radial;
michael@0 696
michael@0 697 if (!_pixman_init_gradient (&radial->common, stops, n_stops))
michael@0 698 {
michael@0 699 free (image);
michael@0 700 return NULL;
michael@0 701 }
michael@0 702
michael@0 703 image->type = RADIAL;
michael@0 704
michael@0 705 radial->c1.x = inner->x;
michael@0 706 radial->c1.y = inner->y;
michael@0 707 radial->c1.radius = inner_radius;
michael@0 708 radial->c2.x = outer->x;
michael@0 709 radial->c2.y = outer->y;
michael@0 710 radial->c2.radius = outer_radius;
michael@0 711
michael@0 712 /* warning: this computations may overflow */
michael@0 713 radial->delta.x = radial->c2.x - radial->c1.x;
michael@0 714 radial->delta.y = radial->c2.y - radial->c1.y;
michael@0 715 radial->delta.radius = radial->c2.radius - radial->c1.radius;
michael@0 716
michael@0 717 /* computed exactly, then cast to double -> every bit of the double
michael@0 718 representation is correct (53 bits) */
michael@0 719 radial->a = dot (radial->delta.x, radial->delta.y, -radial->delta.radius,
michael@0 720 radial->delta.x, radial->delta.y, radial->delta.radius);
michael@0 721 if (radial->a != 0)
michael@0 722 radial->inva = 1. * pixman_fixed_1 / radial->a;
michael@0 723
michael@0 724 radial->mindr = -1. * pixman_fixed_1 * radial->c1.radius;
michael@0 725
michael@0 726 return image;
michael@0 727 }

mercurial