LCOV - code coverage report
Current view: top level - gfx/cairo/libpixman/src - pixman-radial-gradient.c (source / functions) Hit Total Coverage
Test: output.info Lines: 0 200 0.0 %
Date: 2017-07-14 16:53:18 Functions: 0 8 0.0 %
Legend: Lines: hit not hit

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

Generated by: LCOV version 1.13