pixman/pixman/pixman-radial-gradient.c

   1 /*
   2  *
   3  * Copyright © 2000 Keith Packard, member of The XFree86 Project, Inc.
   4  * Copyright © 2000 SuSE, Inc.
   5  *             2005 Lars Knoll & Zack Rusin, Trolltech
   6  * Copyright © 2007 Red Hat, Inc.
   7  *
   8  *
   9  * Permission to use, copy, modify, distribute, and sell this software and its
  10  * documentation for any purpose is hereby granted without fee, provided that
  11  * the above copyright notice appear in all copies and that both that
  12  * copyright notice and this permission notice appear in supporting
  13  * documentation, and that the name of Keith Packard not be used in
  14  * advertising or publicity pertaining to distribution of the software without
  15  * specific, written prior permission.  Keith Packard makes no
  16  * representations about the suitability of this software for any purpose.  It
  17  * is provided "as is" without express or implied warranty.
  18  *
  19  * THE COPYRIGHT HOLDERS DISCLAIM ALL WARRANTIES WITH REGARD TO THIS
  20  * SOFTWARE, INCLUDING ALL IMPLIED WARRANTIES OF MERCHANTABILITY AND
  21  * FITNESS, IN NO EVENT SHALL THE COPYRIGHT HOLDERS BE LIABLE FOR ANY
  22  * SPECIAL, INDIRECT OR CONSEQUENTIAL DAMAGES OR ANY DAMAGES
  23  * WHATSOEVER RESULTING FROM LOSS OF USE, DATA OR PROFITS, WHETHER IN
  24  * AN ACTION OF CONTRACT, NEGLIGENCE OR OTHER TORTIOUS ACTION, ARISING
  25  * OUT OF OR IN CONNECTION WITH THE USE OR PERFORMANCE OF THIS
  26  * SOFTWARE.
  27  */
  28
  29 #ifdef HAVE_CONFIG_H
  30 #include <config.h>
  31 #endif
  32 #include <stdlib.h>
  33 #include <math.h>
  34 #include "pixman-private.h"
  35
  36 static void
  37 radial_gradient_get_scanline_32 (pixman_image_t *image,
  38                                  int             x,
  39                                  int             y,
  40                                  int             width,
  41                                  uint32_t *      buffer,
  42                                  const uint32_t *mask,
  43                                  uint32_t        mask_bits)
  44 {
  45     /*
  46      * In the radial gradient problem we are given two circles (c₁,r₁) and
  47      * (c₂,r₂) that define the gradient itself. Then, for any point p, we
  48      * must compute the value(s) of t within [0.0, 1.0] representing the
  49      * circle(s) that would color the point.
  50      *
  51      * There are potentially two values of t since the point p can be
  52      * colored by both sides of the circle, (which happens whenever one
  53      * circle is not entirely contained within the other).
  54      *
  55      * If we solve for a value of t that is outside of [0.0, 1.0] then we
  56      * use the extend mode (NONE, REPEAT, REFLECT, or PAD) to map to a
  57      * value within [0.0, 1.0].
  58      *
  59      * Here is an illustration of the problem:
  60      *
  61      *              p₂
  62      *           p  •
  63      *           •   ╲
  64      *        ·       ╲r₂
  65      *  p₁ ·           ╲
  66      *  •              θ╲
  67      *   ╲             ╌╌•
  68      *    ╲r₁        ·   c₂
  69      *    θ╲    ·
  70      *    ╌╌•
  71      *      c₁
  72      *
  73      * Given (c₁,r₁), (c₂,r₂) and p, we must find an angle θ such that two
  74      * points p₁ and p₂ on the two circles are collinear with p. Then, the
  75      * desired value of t is the ratio of the length of p₁p to the length
  76      * of p₁p₂.
  77      *
  78      * So, we have six unknown values: (p₁x, p₁y), (p₂x, p₂y), θ and t.
  79      * We can also write six equations that constrain the problem:
  80      *
  81      * Point p₁ is a distance r₁ from c₁ at an angle of θ:
  82      *
  83      *  1. p₁x = c₁x + r₁·cos θ
  84      *  2. p₁y = c₁y + r₁·sin θ
  85      *
  86      * Point p₂ is a distance r₂ from c₂ at an angle of θ:
  87      *
  88      *  3. p₂x = c₂x + r2·cos θ
  89      *  4. p₂y = c₂y + r2·sin θ
  90      *
  91      * Point p lies at a fraction t along the line segment p₁p₂:
  92      *
  93      *  5. px = t·p₂x + (1-t)·p₁x
  94      *  6. py = t·p₂y + (1-t)·p₁y
  95      *
  96      * To solve, first subtitute 1-4 into 5 and 6:
  97      *
  98      * px = t·(c₂x + r₂·cos θ) + (1-t)·(c₁x + r₁·cos θ)
  99      * py = t·(c₂y + r₂·sin θ) + (1-t)·(c₁y + r₁·sin θ)
 100      *
 101      * Then solve each for cos θ and sin θ expressed as a function of t:
 102      *
 103      * cos θ = (-(c₂x - c₁x)·t + (px - c₁x)) / ((r₂-r₁)·t + r₁)
 104      * sin θ = (-(c₂y - c₁y)·t + (py - c₁y)) / ((r₂-r₁)·t + r₁)
 105      *
 106      * To simplify this a bit, we define new variables for several of the
 107      * common terms as shown below:
 108      *
 109      *              p₂
 110      *           p  •
 111      *           •   ╲
 112      *        ·  ┆    ╲r₂
 113      *  p₁ ·     ┆     ╲
 114      *  •     pdy┆      ╲
 115      *   ╲       ┆       •c₂
 116      *    ╲r₁    ┆   ·   ┆
 117      *     ╲    ·┆       ┆cdy
 118      *      •╌╌╌╌┴╌╌╌╌╌╌╌┘
 119      *    c₁  pdx   cdx
 120      *
 121      * cdx = (c₂x - c₁x)
 122      * cdy = (c₂y - c₁y)
 123      *  dr =  r₂-r₁
 124      * pdx =  px - c₁x
 125      * pdy =  py - c₁y
 126      *
 127      * Note that cdx, cdy, and dr do not depend on point p at all, so can
 128      * be pre-computed for the entire gradient. The simplifed equations
 129      * are now:
 130      *
 131      * cos θ = (-cdx·t + pdx) / (dr·t + r₁)
 132      * sin θ = (-cdy·t + pdy) / (dr·t + r₁)
 133      *
 134      * Finally, to get a single function of t and eliminate the last
 135      * unknown θ, we use the identity sin²θ + cos²θ = 1. First, square
 136      * each equation, (we knew a quadratic was coming since it must be
 137      * possible to obtain two solutions in some cases):
 138      *
 139      * cos²θ = (cdx²t² - 2·cdx·pdx·t + pdx²) / (dr²·t² + 2·r₁·dr·t + r₁²)
 140      * sin²θ = (cdy²t² - 2·cdy·pdy·t + pdy²) / (dr²·t² + 2·r₁·dr·t + r₁²)
 141      *
 142      * Then add both together, set the result equal to 1, and express as a
 143      * standard quadratic equation in t of the form At² + Bt + C = 0
 144      *
 145      * (cdx² + cdy² - dr²)·t² - 2·(cdx·pdx + cdy·pdy + r₁·dr)·t + (pdx² + pdy² - r₁²) = 0
 146      *
 147      * In other words:
 148      *
 149      * A = cdx² + cdy² - dr²
 150      * B = -2·(pdx·cdx + pdy·cdy + r₁·dr)
 151      * C = pdx² + pdy² - r₁²
 152      *
 153      * And again, notice that A does not depend on p, so can be
 154      * precomputed. From here we just use the quadratic formula to solve
 155      * for t:
 156      *
 157      * t = (-2·B ± ⎷(B² - 4·A·C)) / 2·A
 158      */
 159
 160     gradient_t *gradient = (gradient_t *)image;
 161     source_image_t *source = (source_image_t *)image;
 162     radial_gradient_t *radial = (radial_gradient_t *)image;
 163     uint32_t *end = buffer + width;
 164     pixman_gradient_walker_t walker;
 165     pixman_bool_t affine = TRUE;
 166     double cx = 1.;
 167     double cy = 0.;
 168     double cz = 0.;
 169     double rx = x + 0.5;
 170     double ry = y + 0.5;
 171     double rz = 1.;
 172
 173     _pixman_gradient_walker_init (&walker, gradient, source->common.repeat);
 174
 175     if (source->common.transform)
 176     {
 177         pixman_vector_t v;
 178         /* reference point is the center of the pixel */
 179         v.vector[0] = pixman_int_to_fixed (x) + pixman_fixed_1 / 2;
 180         v.vector[1] = pixman_int_to_fixed (y) + pixman_fixed_1 / 2;
 181         v.vector[2] = pixman_fixed_1;
 182
 183         if (!pixman_transform_point_3d (source->common.transform, &v))
 184             return;
 185
 186         cx = source->common.transform->matrix[0][0] / 65536.;
 187         cy = source->common.transform->matrix[1][0] / 65536.;
 188         cz = source->common.transform->matrix[2][0] / 65536.;
 189
 190         rx = v.vector[0] / 65536.;
 191         ry = v.vector[1] / 65536.;
 192         rz = v.vector[2] / 65536.;
 193
 194         affine =
 195             source->common.transform->matrix[2][0] == 0 &&
 196             v.vector[2] == pixman_fixed_1;
 197     }
 198
 199     if (affine)
 200     {
 201         /* When computing t over a scanline, we notice that some expressions
 202          * are constant so we can compute them just once. Given:
 203          *
 204          * t = (-2·B ± ⎷(B² - 4·A·C)) / 2·A
 205          *
 206          * where
 207          *
 208          * A = cdx² + cdy² - dr² [precomputed as radial->A]
 209          * B = -2·(pdx·cdx + pdy·cdy + r₁·dr)
 210          * C = pdx² + pdy² - r₁²
 211          *
 212          * Since we have an affine transformation, we know that (pdx, pdy)
 213          * increase linearly with each pixel,
 214          *
 215          * pdx = pdx₀ + n·cx,
 216          * pdy = pdy₀ + n·cy,
 217          *
 218          * we can then express B in terms of an linear increment along
 219          * the scanline:
 220          *
 221          * B = B₀ + n·cB, with
 222          * B₀ = -2·(pdx₀·cdx + pdy₀·cdy + r₁·dr) and
 223          * cB = -2·(cx·cdx + cy·cdy)
 224          *
 225          * Thus we can replace the full evaluation of B per-pixel (4 multiplies,
 226          * 2 additions) with a single addition.
 227          */
 228         double r1   = radial->c1.radius / 65536.;
 229         double r1sq = r1 * r1;
 230         double pdx  = rx - radial->c1.x / 65536.;
 231         double pdy  = ry - radial->c1.y / 65536.;
 232         double A = radial->A;
 233         double invA = -65536. / (2. * A);
 234         double A4 = -4. * A;
 235         double B  = -2. * (pdx*radial->cdx + pdy*radial->cdy + r1*radial->dr);
 236         double cB = -2. *  (cx*radial->cdx +  cy*radial->cdy);
 237         pixman_bool_t invert = A * radial->dr < 0;
 238
 239         while (buffer < end)
 240         {
 241             if (!mask || *mask++ & mask_bits)
 242             {
 243                 pixman_fixed_48_16_t t;
 244                 double det = B * B + A4 * (pdx * pdx + pdy * pdy - r1sq);
 245                 if (det <= 0.)
 246                     t = (pixman_fixed_48_16_t) (B * invA);
 247                 else if (invert)
 248                     t = (pixman_fixed_48_16_t) ((B + sqrt (det)) * invA);
 249                 else
 250                     t = (pixman_fixed_48_16_t) ((B - sqrt (det)) * invA);
 251
 252                 *buffer = _pixman_gradient_walker_pixel (&walker, t);
 253             }
 254             ++buffer;
 255
 256             pdx += cx;
 257             pdy += cy;
 258             B += cB;
 259         }
 260     }
 261     else
 262     {
 263         /* projective */
 264         while (buffer < end)
 265         {
 266             if (!mask || *mask++ & mask_bits)
 267             {
 268                 double pdx, pdy;
 269                 double B, C;
 270                 double det;
 271                 double c1x = radial->c1.x / 65536.0;
 272                 double c1y = radial->c1.y / 65536.0;
 273                 double r1  = radial->c1.radius / 65536.0;
 274                 pixman_fixed_48_16_t t;
 275                 double x, y;
 276
 277                 if (rz != 0)
 278                 {
 279                     x = rx / rz;
 280                     y = ry / rz;
 281                 }
 282                 else
 283                 {
 284                     x = y = 0.;
 285                 }
 286
 287                 pdx = x - c1x;
 288                 pdy = y - c1y;
 289
 290                 B = -2 * (pdx * radial->cdx +
 291                           pdy * radial->cdy +
 292                           r1 * radial->dr);
 293                 C = (pdx * pdx + pdy * pdy - r1 * r1);
 294
 295                 det = (B * B) - (4 * radial->A * C);
 296                 if (det < 0.0)
 297                     det = 0.0;
 298
 299                 if (radial->A * radial->dr < 0)
 300                     t = (pixman_fixed_48_16_t) ((-B - sqrt (det)) / (2.0 * radial->A) * 65536);
 301                 else
 302                     t = (pixman_fixed_48_16_t) ((-B + sqrt (det)) / (2.0 * radial->A) * 65536);
 303
 304                 *buffer = _pixman_gradient_walker_pixel (&walker, t);
 305             }
 306
 307             ++buffer;
 308
 309             rx += cx;
 310             ry += cy;
 311             rz += cz;
 312         }
 313     }
 314 }
 315
 316 static void
 317 radial_gradient_property_changed (pixman_image_t *image)
 318 {
 319     image->common.get_scanline_32 = radial_gradient_get_scanline_32;
 320     image->common.get_scanline_64 = _pixman_image_get_scanline_generic_64;
 321 }
 322
 323 PIXMAN_EXPORT pixman_image_t *
 324 pixman_image_create_radial_gradient (pixman_point_fixed_t *        inner,
 325                                      pixman_point_fixed_t *        outer,
 326                                      pixman_fixed_t                inner_radius,
 327                                      pixman_fixed_t                outer_radius,
 328                                      const pixman_gradient_stop_t *stops,
 329                                      int                           n_stops)
 330 {
 331     pixman_image_t *image;
 332     radial_gradient_t *radial;
 333
 334     return_val_if_fail (n_stops >= 2, NULL);
 335
 336     image = _pixman_image_allocate ();
 337
 338     if (!image)
 339         return NULL;
 340
 341     radial = &image->radial;
 342
 343     if (!_pixman_init_gradient (&radial->common, stops, n_stops))
 344     {
 345         free (image);
 346         return NULL;
 347     }
 348
 349     image->type = RADIAL;
 350
 351     radial->c1.x = inner->x;
 352     radial->c1.y = inner->y;
 353     radial->c1.radius = inner_radius;
 354     radial->c2.x = outer->x;
 355     radial->c2.y = outer->y;
 356     radial->c2.radius = outer_radius;
 357     radial->cdx = pixman_fixed_to_double (radial->c2.x - radial->c1.x);
 358     radial->cdy = pixman_fixed_to_double (radial->c2.y - radial->c1.y);
 359     radial->dr = pixman_fixed_to_double (radial->c2.radius - radial->c1.radius);
 360     radial->A = (radial->cdx * radial->cdx +
 361                  radial->cdy * radial->cdy -
 362                  radial->dr  * radial->dr);
 363
 364     image->common.property_changed = radial_gradient_property_changed;
 365
 366     return image;
 367 }
 368