1 /* Copyright (C) 2001-2012 Artifex Software, Inc.
2 All Rights Reserved.
3
4 This software is provided AS-IS with no warranty, either express or
5 implied.
6
7 This software is distributed under license and may not be copied,
8 modified or distributed except as expressly authorized under the terms
9 of the license contained in the file LICENSE in this distribution.
10
11 Refer to licensing information at http://www.artifex.com or contact
12 Artifex Software, Inc., 7 Mt. Lassen Drive - Suite A-134, San Rafael,
13 CA 94903, U.S.A., +1(415)492-9861, for further information.
14 */
15
16
17 /* Rendering for non-mesh shadings */
18 #include "math_.h"
19 #include "memory_.h"
20 #include "gx.h"
21 #include "gserrors.h"
22 #include "gsmatrix.h" /* for gscoord.h */
23 #include "gscoord.h"
24 #include "gspath.h"
25 #include "gsptype2.h"
26 #include "gxcspace.h"
27 #include "gxdcolor.h"
28 #include "gxfarith.h"
29 #include "gxfixed.h"
30 #include "gxistate.h"
31 #include "gxpath.h"
32 #include "gxshade.h"
33 #include "gxdevcli.h"
34 #include "gxshade4.h"
35 #include "vdtrace.h"
36 #include "gsicc_cache.h"
37
38 #define VD_TRACE_AXIAL_PATCH 1
39 #define VD_TRACE_RADIAL_PATCH 1
40 #define VD_TRACE_FUNCTIONAL_PATCH 1
41
42 /* ---------------- Function-based shading ---------------- */
43
44 typedef struct Fb_frame_s { /* A rudiment of old code. */
45 gs_rect region;
46 gs_client_color cc[4]; /* colors at 4 corners */
47 int state;
48 } Fb_frame_t;
49
50 typedef struct Fb_fill_state_s {
51 shading_fill_state_common;
52 const gs_shading_Fb_t *psh;
53 gs_matrix_fixed ptm; /* parameter space -> device space */
54 Fb_frame_t frame;
55 } Fb_fill_state_t;
56 /****** NEED GC DESCRIPTOR ******/
57
58 static inline void
make_other_poles(patch_curve_t curve[4])59 make_other_poles(patch_curve_t curve[4])
60 {
61 int i, j;
62
63 for (i = 0; i < 4; i++) {
64 j = (i + 1) % 4;
65 curve[i].control[0].x = (curve[i].vertex.p.x * 2 + curve[j].vertex.p.x) / 3;
66 curve[i].control[0].y = (curve[i].vertex.p.y * 2 + curve[j].vertex.p.y) / 3;
67 curve[i].control[1].x = (curve[i].vertex.p.x + curve[j].vertex.p.x * 2) / 3;
68 curve[i].control[1].y = (curve[i].vertex.p.y + curve[j].vertex.p.y * 2) / 3;
69 curve[i].straight = true;
70 }
71 }
72
73 /* Transform a point with a fixed-point result. */
74 static void
gs_point_transform2fixed_clamped(const gs_matrix_fixed * pmat,floatp x,floatp y,gs_fixed_point * ppt)75 gs_point_transform2fixed_clamped(const gs_matrix_fixed * pmat,
76 floatp x, floatp y, gs_fixed_point * ppt)
77 {
78 gs_point fpt;
79
80 gs_point_transform(x, y, (const gs_matrix *)pmat, &fpt);
81 ppt->x = clamp_coord(fpt.x);
82 ppt->y = clamp_coord(fpt.y);
83 }
84
85 static int
Fb_fill_region(Fb_fill_state_t * pfs,const gs_fixed_rect * rect)86 Fb_fill_region(Fb_fill_state_t * pfs, const gs_fixed_rect *rect)
87 {
88 patch_fill_state_t pfs1;
89 patch_curve_t curve[4];
90 Fb_frame_t * fp = &pfs->frame;
91 int code;
92
93 if (VD_TRACE_FUNCTIONAL_PATCH && vd_allowed('s')) {
94 vd_get_dc('s');
95 vd_set_shift(0, 0);
96 vd_set_scale(0.01);
97 vd_set_origin(0, 0);
98 }
99 memcpy(&pfs1, (shading_fill_state_t *)pfs, sizeof(shading_fill_state_t));
100 pfs1.Function = pfs->psh->params.Function;
101 code = init_patch_fill_state(&pfs1);
102 if (code < 0)
103 return code;
104 pfs1.maybe_self_intersecting = false;
105 pfs1.n_color_args = 2;
106 pfs1.rect = *rect;
107 gs_point_transform2fixed(&pfs->ptm, fp->region.p.x, fp->region.p.y, &curve[0].vertex.p);
108 gs_point_transform2fixed(&pfs->ptm, fp->region.q.x, fp->region.p.y, &curve[1].vertex.p);
109 gs_point_transform2fixed(&pfs->ptm, fp->region.q.x, fp->region.q.y, &curve[2].vertex.p);
110 gs_point_transform2fixed(&pfs->ptm, fp->region.p.x, fp->region.q.y, &curve[3].vertex.p);
111 make_other_poles(curve);
112 curve[0].vertex.cc[0] = fp->region.p.x; curve[0].vertex.cc[1] = fp->region.p.y;
113 curve[1].vertex.cc[0] = fp->region.q.x; curve[1].vertex.cc[1] = fp->region.p.y;
114 curve[2].vertex.cc[0] = fp->region.q.x; curve[2].vertex.cc[1] = fp->region.q.y;
115 curve[3].vertex.cc[0] = fp->region.p.x; curve[3].vertex.cc[1] = fp->region.q.y;
116 code = patch_fill(&pfs1, curve, NULL, NULL);
117 if (term_patch_fill_state(&pfs1))
118 return_error(gs_error_unregistered); /* Must not happen. */
119 if (VD_TRACE_FUNCTIONAL_PATCH && vd_allowed('s'))
120 vd_release_dc;
121 return code;
122 }
123
124 int
gs_shading_Fb_fill_rectangle(const gs_shading_t * psh0,const gs_rect * rect,const gs_fixed_rect * rect_clip,gx_device * dev,gs_imager_state * pis)125 gs_shading_Fb_fill_rectangle(const gs_shading_t * psh0, const gs_rect * rect,
126 const gs_fixed_rect * rect_clip,
127 gx_device * dev, gs_imager_state * pis)
128 {
129 const gs_shading_Fb_t * const psh = (const gs_shading_Fb_t *)psh0;
130 gs_matrix save_ctm;
131 int xi, yi, code;
132 float x[2], y[2];
133 Fb_fill_state_t state;
134
135 code = shade_init_fill_state((shading_fill_state_t *) & state, psh0, dev, pis);
136 if (code < 0)
137 return code;
138 state.psh = psh;
139 /****** HACK FOR FIXED-POINT MATRIX MULTIPLY ******/
140 gs_currentmatrix((gs_state *) pis, &save_ctm);
141 gs_concat((gs_state *) pis, &psh->params.Matrix);
142 state.ptm = pis->ctm;
143 gs_setmatrix((gs_state *) pis, &save_ctm);
144 /* Compute the parameter X and Y ranges. */
145 {
146 gs_rect pbox;
147
148 gs_bbox_transform_inverse(rect, &psh->params.Matrix, &pbox);
149 x[0] = max(pbox.p.x, psh->params.Domain[0]);
150 x[1] = min(pbox.q.x, psh->params.Domain[1]);
151 y[0] = max(pbox.p.y, psh->params.Domain[2]);
152 y[1] = min(pbox.q.y, psh->params.Domain[3]);
153 }
154 if (x[0] > x[1] || y[0] > y[1]) {
155 /* The region is outside the shading area. */
156 if (state.icclink != NULL) gsicc_release_link(state.icclink);
157 return 0;
158 }
159 for (xi = 0; xi < 2; ++xi)
160 for (yi = 0; yi < 2; ++yi) {
161 float v[2];
162
163 v[0] = x[xi], v[1] = y[yi];
164 gs_function_evaluate(psh->params.Function, v,
165 state.frame.cc[yi * 2 + xi].paint.values);
166 }
167 state.frame.region.p.x = x[0];
168 state.frame.region.p.y = y[0];
169 state.frame.region.q.x = x[1];
170 state.frame.region.q.y = y[1];
171 code = Fb_fill_region(&state, rect_clip);
172 if (state.icclink != NULL) gsicc_release_link(state.icclink);
173 return code;
174 }
175
176 /* ---------------- Axial shading ---------------- */
177
178 typedef struct A_fill_state_s {
179 const gs_shading_A_t *psh;
180 gs_point delta;
181 double length;
182 double t0, t1;
183 double v0, v1, u0, u1;
184 } A_fill_state_t;
185 /****** NEED GC DESCRIPTOR ******/
186
187 /* Note t0 and t1 vary over [0..1], not the Domain. */
188
189 static int
A_fill_region(A_fill_state_t * pfs,patch_fill_state_t * pfs1)190 A_fill_region(A_fill_state_t * pfs, patch_fill_state_t *pfs1)
191 {
192 const gs_shading_A_t * const psh = pfs->psh;
193 double x0 = psh->params.Coords[0] + pfs->delta.x * pfs->v0;
194 double y0 = psh->params.Coords[1] + pfs->delta.y * pfs->v0;
195 double x1 = psh->params.Coords[0] + pfs->delta.x * pfs->v1;
196 double y1 = psh->params.Coords[1] + pfs->delta.y * pfs->v1;
197 double h0 = pfs->u0, h1 = pfs->u1;
198 patch_curve_t curve[4];
199
200 gs_point_transform2fixed(&pfs1->pis->ctm, x0 + pfs->delta.y * h0, y0 - pfs->delta.x * h0, &curve[0].vertex.p);
201 gs_point_transform2fixed(&pfs1->pis->ctm, x1 + pfs->delta.y * h0, y1 - pfs->delta.x * h0, &curve[1].vertex.p);
202 gs_point_transform2fixed(&pfs1->pis->ctm, x1 + pfs->delta.y * h1, y1 - pfs->delta.x * h1, &curve[2].vertex.p);
203 gs_point_transform2fixed(&pfs1->pis->ctm, x0 + pfs->delta.y * h1, y0 - pfs->delta.x * h1, &curve[3].vertex.p);
204 curve[0].vertex.cc[0] = pfs->t0; /* The element cc[1] is set to a dummy value against */
205 curve[1].vertex.cc[0] = pfs->t1; /* interrupts while an idle priocessing in gxshade.6.c . */
206 curve[2].vertex.cc[0] = pfs->t1;
207 curve[3].vertex.cc[0] = pfs->t0;
208 curve[0].vertex.cc[1] = 0; /* The element cc[1] is set to a dummy value against */
209 curve[1].vertex.cc[1] = 0; /* interrupts while an idle priocessing in gxshade.6.c . */
210 curve[2].vertex.cc[1] = 0;
211 curve[3].vertex.cc[1] = 0;
212 make_other_poles(curve);
213 return patch_fill(pfs1, curve, NULL, NULL);
214 }
215
216 static inline int
gs_shading_A_fill_rectangle_aux(const gs_shading_t * psh0,const gs_rect * rect,const gs_fixed_rect * clip_rect,gx_device * dev,gs_imager_state * pis)217 gs_shading_A_fill_rectangle_aux(const gs_shading_t * psh0, const gs_rect * rect,
218 const gs_fixed_rect *clip_rect,
219 gx_device * dev, gs_imager_state * pis)
220 {
221 const gs_shading_A_t *const psh = (const gs_shading_A_t *)psh0;
222 gs_function_t * const pfn = psh->params.Function;
223 gs_matrix cmat;
224 gs_rect t_rect;
225 A_fill_state_t state;
226 patch_fill_state_t pfs1;
227 float d0 = psh->params.Domain[0], d1 = psh->params.Domain[1];
228 float dd = d1 - d0;
229 double t0, t1;
230 gs_point dist;
231 int code;
232
233 state.psh = psh;
234 code = shade_init_fill_state((shading_fill_state_t *)&pfs1, psh0, dev, pis);
235 if (code < 0)
236 return code;
237 pfs1.Function = pfn;
238 pfs1.rect = *clip_rect;
239 code = init_patch_fill_state(&pfs1);
240 if (code < 0) {
241 if (pfs1.icclink != NULL) gsicc_release_link(pfs1.icclink);
242 return code;
243 }
244 pfs1.maybe_self_intersecting = false;
245 pfs1.function_arg_shift = 1;
246 /*
247 * Compute the parameter range. We construct a matrix in which
248 * (0,0) corresponds to t = 0 and (0,1) corresponds to t = 1,
249 * and use it to inverse-map the rectangle to be filled.
250 */
251 cmat.tx = psh->params.Coords[0];
252 cmat.ty = psh->params.Coords[1];
253 state.delta.x = psh->params.Coords[2] - psh->params.Coords[0];
254 state.delta.y = psh->params.Coords[3] - psh->params.Coords[1];
255 cmat.yx = state.delta.x;
256 cmat.yy = state.delta.y;
257 cmat.xx = cmat.yy;
258 cmat.xy = -cmat.yx;
259 gs_bbox_transform_inverse(rect, &cmat, &t_rect);
260 t0 = min(max(t_rect.p.y, 0), 1);
261 t1 = max(min(t_rect.q.y, 1), 0);
262 state.v0 = t0;
263 state.v1 = t1;
264 state.u0 = t_rect.p.x;
265 state.u1 = t_rect.q.x;
266 state.t0 = t0 * dd + d0;
267 state.t1 = t1 * dd + d0;
268 gs_distance_transform(state.delta.x, state.delta.y, &ctm_only(pis),
269 &dist);
270 state.length = hypot(dist.x, dist.y); /* device space line length */
271 code = A_fill_region(&state, &pfs1);
272 if (psh->params.Extend[0] && t0 > t_rect.p.y) {
273 if (code < 0) {
274 if (pfs1.icclink != NULL) gsicc_release_link(pfs1.icclink);
275 return code;
276 }
277 /* Use the general algorithm, because we need the trapping. */
278 state.v0 = t_rect.p.y;
279 state.v1 = t0;
280 state.t0 = state.t1 = t0 * dd + d0;
281 code = A_fill_region(&state, &pfs1);
282 }
283 if (psh->params.Extend[1] && t1 < t_rect.q.y) {
284 if (code < 0) {
285 if (pfs1.icclink != NULL) gsicc_release_link(pfs1.icclink);
286 return code;
287 }
288 /* Use the general algorithm, because we need the trapping. */
289 state.v0 = t1;
290 state.v1 = t_rect.q.y;
291 state.t0 = state.t1 = t1 * dd + d0;
292 code = A_fill_region(&state, &pfs1);
293 }
294 if (pfs1.icclink != NULL) gsicc_release_link(pfs1.icclink);
295 if (term_patch_fill_state(&pfs1))
296 return_error(gs_error_unregistered); /* Must not happen. */
297 return code;
298 }
299
300 int
gs_shading_A_fill_rectangle(const gs_shading_t * psh0,const gs_rect * rect,const gs_fixed_rect * rect_clip,gx_device * dev,gs_imager_state * pis)301 gs_shading_A_fill_rectangle(const gs_shading_t * psh0, const gs_rect * rect,
302 const gs_fixed_rect * rect_clip,
303 gx_device * dev, gs_imager_state * pis)
304 {
305 int code;
306
307 if (VD_TRACE_AXIAL_PATCH && vd_allowed('s')) {
308 vd_get_dc('s');
309 vd_set_shift(0, 0);
310 vd_set_scale(0.01);
311 vd_set_origin(0, 0);
312 }
313 code = gs_shading_A_fill_rectangle_aux(psh0, rect, rect_clip, dev, pis);
314 if (VD_TRACE_AXIAL_PATCH && vd_allowed('s'))
315 vd_release_dc;
316 return code;
317 }
318
319 /* ---------------- Radial shading ---------------- */
320
321 /* Some notes on what I have struggled to understand about the following
322 * function. This function renders the 'tube' given by interpolating one
323 * circle to another.
324 *
325 * The first circle is at (x0, y0) with radius r0, and has 'color' t0.
326 * The other circle is at (x1, y1) with radius r1, and has 'color' t1.
327 *
328 * We perform this rendering by approximating each quadrant of the 'tube'
329 * by a tensor patch. The tensor patch is formed by taking a curve along
330 * 1/4 of the circumference of the first circle, a straight line to the
331 * equivalent point on the circumference of the second circle, a curve
332 * back along the circumference of the second circle, and then a straight
333 * line back to where we started.
334 *
335 * There is additional logic in this function that forms the directions of
336 * the curves differently for different quadrants. This is done to ensure
337 * that we always paint 'around' the tube from the back towards the front,
338 * so we don't get unexpected regions showing though. This is explained more
339 * below.
340 *
341 * The original code here examined the position change between the two
342 * circles dx and dy. Based upon this vector it would pick which quadrant/
343 * tensor patch to draw first. It would draw the quadrants/tensor patches
344 * in anticlockwise order. Presumably this was intended to be done so that
345 * the 'top' quadrant would be drawn last.
346 *
347 * Unfortunately this did not always work; see bug 692513. If the quadrants
348 * were rendered in the order 0,1,2,3, the rendering of 1 was leaving traces
349 * on top of 0, which was unexpected.
350 *
351 * I have therefore altered the code slightly; rather than picking a start
352 * quadrant and moving anticlockwise, we now draw the 'undermost' quadrant,
353 * then the two adjacent quadrants, then the topmost quadrant.
354 *
355 * For the purposes of explaination, we shall label the octants as below:
356 *
357 * \2|1/ and Quadrants as: |
358 * 3\|/0 Q1 | Q0
359 * ---+--- ----+----
360 * 4/|\7 Q2 | Q3
361 * /5|6\ |
362 *
363 * We find (dx,dy), the difference between the centres of the circles.
364 * We look to see which octant this falls in. Firstly, this tells us which
365 * quadrant of the circle we need to draw first (Octant n, starts with
366 * Quadrant floor(n/2)). Secondly, it tells us which direction to form the
367 * tensor patch in; we always want to draw from the side 'closest' to
368 * dx/dy to the side further away. This ensures that we don't overwrite
369 * pixels in the incorrect order as the patch decomposes.
370 */
371 static int
R_tensor_annulus(patch_fill_state_t * pfs,double x0,double y0,double r0,double t0,double x1,double y1,double r1,double t1)372 R_tensor_annulus(patch_fill_state_t *pfs,
373 double x0, double y0, double r0, double t0,
374 double x1, double y1, double r1, double t1)
375 {
376 double dx = x1 - x0, dy = y1 - y0;
377 double d = hypot(dx, dy);
378 gs_point p0, p1, pc0, pc1;
379 int k, j, code, dirn;
380 bool inside = 0;
381
382 /* pc0 and pc1 are the centres of the respective circles. */
383 pc0.x = x0, pc0.y = y0;
384 pc1.x = x1, pc1.y = y1;
385 /* Set p0 up so it's a unit vector giving the direction of 90 degrees
386 * to the right of the major axis as we move from p0c to p1c. */
387 if (r0 + d <= r1 || r1 + d <= r0) {
388 /* One circle is inside another one.
389 Use any subdivision,
390 but don't depend on dx, dy, which may be too small. */
391 p0.x = 0, p0.y = -1, dirn = 0;
392 /* Align stripes along radii for faster triangulation : */
393 inside = 1;
394 } else {
395 /* Must generate canonic quadrangle arcs,
396 because we approximate them with curves. */
397 if(dx >= 0) {
398 if (dy >= 0)
399 p0.x = 1, p0.y = 0, dirn = (dx >= dy ? 1 : 0);
400 else
401 p0.x = 0, p0.y = -1, dirn = (dx >= -dy ? 0 : 1);
402 } else {
403 if (dy >= 0)
404 p0.x = 0, p0.y = 1, dirn = (-dx >= dy ? 1 : 0);
405 else
406 p0.x = -1, p0.y = 0, dirn = (-dx >= -dy ? 0 : 1);
407 }
408 }
409 /* fixme: wish: cut invisible parts off.
410 Note : when r0 != r1 the invisible part is not a half circle. */
411 for (k = 0; k < 4; k++) {
412 gs_point p[12];
413 patch_curve_t curve[4];
414
415 /* Set p1 to be 90 degrees anticlockwise from p0 */
416 p1.x = -p0.y; p1.y = p0.x;
417 if (dirn == 0) { /* Clockwise */
418 make_quadrant_arc(p + 0, &pc0, &p1, &p0, r0);
419 make_quadrant_arc(p + 6, &pc1, &p0, &p1, r1);
420 } else { /* Anticlockwise */
421 make_quadrant_arc(p + 0, &pc0, &p0, &p1, r0);
422 make_quadrant_arc(p + 6, &pc1, &p1, &p0, r1);
423 }
424 p[4].x = (p[3].x * 2 + p[6].x) / 3;
425 p[4].y = (p[3].y * 2 + p[6].y) / 3;
426 p[5].x = (p[3].x + p[6].x * 2) / 3;
427 p[5].y = (p[3].y + p[6].y * 2) / 3;
428 p[10].x = (p[9].x * 2 + p[0].x) / 3;
429 p[10].y = (p[9].y * 2 + p[0].y) / 3;
430 p[11].x = (p[9].x + p[0].x * 2) / 3;
431 p[11].y = (p[9].y + p[0].y * 2) / 3;
432 for (j = 0; j < 4; j++) {
433 int jj = (j + inside) % 4;
434
435 if (gs_point_transform2fixed(&pfs->pis->ctm, p[j*3 + 0].x, p[j*3 + 0].y, &curve[jj].vertex.p) < 0)
436 gs_point_transform2fixed_clamped(&pfs->pis->ctm, p[j*3 + 0].x, p[j*3 + 0].y, &curve[jj].vertex.p);
437
438 if (gs_point_transform2fixed(&pfs->pis->ctm, p[j*3 + 1].x, p[j*3 + 1].y, &curve[jj].control[0]) < 0)
439 gs_point_transform2fixed_clamped(&pfs->pis->ctm, p[j*3 + 1].x, p[j*3 + 1].y, &curve[jj].control[0]);
440
441 if (gs_point_transform2fixed(&pfs->pis->ctm, p[j*3 + 2].x, p[j*3 + 2].y, &curve[jj].control[1]) < 0)
442 gs_point_transform2fixed_clamped(&pfs->pis->ctm, p[j*3 + 2].x, p[j*3 + 2].y, &curve[jj].control[1]);
443 curve[j].straight = (((j + inside) & 1) != 0);
444 }
445 curve[(0 + inside) % 4].vertex.cc[0] = t0;
446 curve[(1 + inside) % 4].vertex.cc[0] = t0;
447 curve[(2 + inside) % 4].vertex.cc[0] = t1;
448 curve[(3 + inside) % 4].vertex.cc[0] = t1;
449 curve[0].vertex.cc[1] = curve[1].vertex.cc[1] = 0; /* Initialize against FPE. */
450 curve[2].vertex.cc[1] = curve[3].vertex.cc[1] = 0; /* Initialize against FPE. */
451 code = patch_fill(pfs, curve, NULL, NULL);
452 if (code < 0)
453 return code;
454 /* Move p0 to be ready for the next position */
455 if (k == 0) {
456 /* p0 moves clockwise */
457 p1 = p0;
458 p0.x = p1.y; p0.y = -p1.x;
459 dirn = 0;
460 } else if (k == 1) {
461 /* p0 flips sides */
462 p0.x = -p0.x; p0.y = -p0.y;
463 dirn = 1;
464 } else if (k == 2) {
465 /* p0 moves anti-clockwise */
466 p1 = p0;
467 p0.x = -p1.y; p0.y = p1.x;
468 dirn = 0;
469 }
470 }
471 return 0;
472 }
473
474 static int
R_outer_circle(patch_fill_state_t * pfs,const gs_rect * rect,double x0,double y0,double r0,double x1,double y1,double r1,double * x2,double * y2,double * r2)475 R_outer_circle(patch_fill_state_t *pfs, const gs_rect *rect,
476 double x0, double y0, double r0,
477 double x1, double y1, double r1,
478 double *x2, double *y2, double *r2)
479 {
480 double dx = x1 - x0, dy = y1 - y0;
481 double sp, sq, s;
482
483 /* Compute a cone circle, which contacts the rect externally. */
484 /* Don't bother with all 4 sides of the rect,
485 just do with the X or Y span only,
486 so it's not an exact contact, sorry. */
487 if (any_abs(dx) > any_abs(dy)) {
488 /* Solving :
489 x0 + (x1 - x0) * sq + r0 + (r1 - r0) * sq == bbox_px
490 (x1 - x0) * sp + (r1 - r0) * sp == bbox_px - x0 - r0
491 sp = (bbox_px - x0 - r0) / (x1 - x0 + r1 - r0)
492
493 x0 + (x1 - x0) * sq - r0 - (r1 - r0) * sq == bbox_qx
494 (x1 - x0) * sq - (r1 - r0) * sq == bbox_x - x0 + r0
495 sq = (bbox_x - x0 + r0) / (x1 - x0 - r1 + r0)
496 */
497 if (x1 - x0 + r1 - r0 == 0) /* We checked for obtuse cone. */
498 return_error(gs_error_unregistered); /* Must not happen. */
499 if (x1 - x0 - r1 + r0 == 0) /* We checked for obtuse cone. */
500 return_error(gs_error_unregistered); /* Must not happen. */
501 sp = (rect->p.x - x0 - r0) / (x1 - x0 + r1 - r0);
502 sq = (rect->q.x - x0 + r0) / (x1 - x0 - r1 + r0);
503 } else {
504 /* Same by Y. */
505 if (y1 - y0 + r1 - r0 == 0) /* We checked for obtuse cone. */
506 return_error(gs_error_unregistered); /* Must not happen. */
507 if (y1 - y0 - r1 + r0 == 0) /* We checked for obtuse cone. */
508 return_error(gs_error_unregistered); /* Must not happen. */
509 sp = (rect->p.y - y0 - r0) / (y1 - y0 + r1 - r0);
510 sq = (rect->q.y - y0 + r0) / (y1 - y0 - r1 + r0);
511 }
512 if (sp >= 1 && sq >= 1)
513 s = max(sp, sq);
514 else if(sp >= 1)
515 s = sp;
516 else if (sq >= 1)
517 s = sq;
518 else {
519 /* The circle 1 is outside the rect, use it. */
520 s = 1;
521 }
522 if (r0 + (r1 - r0) * s < 0) {
523 /* Passed the cone apex, use the apex. */
524 s = r0 / (r0 - r1);
525 *r2 = 0;
526 } else
527 *r2 = r0 + (r1 - r0) * s;
528 *x2 = x0 + (x1 - x0) * s;
529 *y2 = y0 + (y1 - y0) * s;
530 return 0;
531 }
532
533 static double
R_rect_radius(const gs_rect * rect,double x0,double y0)534 R_rect_radius(const gs_rect *rect, double x0, double y0)
535 {
536 double d, dd;
537
538 dd = hypot(rect->p.x - x0, rect->p.y - y0);
539 d = hypot(rect->p.x - x0, rect->q.y - y0);
540 dd = max(dd, d);
541 d = hypot(rect->q.x - x0, rect->q.y - y0);
542 dd = max(dd, d);
543 d = hypot(rect->q.x - x0, rect->p.y - y0);
544 dd = max(dd, d);
545 return dd;
546 }
547
548 static int
R_fill_triangle_new(patch_fill_state_t * pfs,const gs_rect * rect,double x0,double y0,double x1,double y1,double x2,double y2,double t)549 R_fill_triangle_new(patch_fill_state_t *pfs, const gs_rect *rect,
550 double x0, double y0, double x1, double y1, double x2, double y2, double t)
551 {
552 shading_vertex_t p0, p1, p2;
553 patch_color_t *c;
554 int code;
555 reserve_colors(pfs, &c, 1); /* Can't fail */
556
557 p0.c = c;
558 p1.c = c;
559 p2.c = c;
560 code = gs_point_transform2fixed(&pfs->pis->ctm, x0, y0, &p0.p);
561 if (code >= 0)
562 code = gs_point_transform2fixed(&pfs->pis->ctm, x1, y1, &p1.p);
563 if (code >= 0)
564 code = gs_point_transform2fixed(&pfs->pis->ctm, x2, y2, &p2.p);
565 if (code >= 0) {
566 c->t[0] = c->t[1] = t;
567 patch_resolve_color(c, pfs);
568 code = mesh_triangle(pfs, &p0, &p1, &p2);
569 }
570 release_colors(pfs, pfs->color_stack, 1);
571 return code;
572 }
573
574 static int
R_obtuse_cone(patch_fill_state_t * pfs,const gs_rect * rect,double x0,double y0,double r0,double x1,double y1,double r1,double t0,double r_rect,bool inwards)575 R_obtuse_cone(patch_fill_state_t *pfs, const gs_rect *rect,
576 double x0, double y0, double r0,
577 double x1, double y1, double r1, double t0, double r_rect,
578 bool inwards)
579 {
580 double dx = x1 - x0, dy = y1 - y0, dr = any_abs(r1 - r0);
581 double d = hypot(dx, dy);
582 /* Assuming dr > d / 3 && d > dr + 1e-7 * (d + dr), see the caller. */
583 double r = r_rect * 1.4143; /* A few bigger than sqrt(2). */
584 double ax, ay, as; /* Cone apex. */
585 double g0; /* The distance from apex to the tangent point of the 0th circle. */
586 int code;
587
588 as = r0 / (r0 - r1);
589 ax = x0 + (x1 - x0) * as;
590 ay = y0 + (y1 - y0) * as;
591 g0 = sqrt(dx * dx + dy * dy - dr * dr) * as;
592 if (g0 < 1e-7 * r0) {
593 /* Nearly degenerate, replace with half-plane. */
594 /* Restrict the half plane with triangle, which covers the rect. */
595 gs_point p0, p1, p2; /* Right tangent limit, apex limit, left tangent linit,
596 (right, left == when looking from the apex). */
597
598 p0.x = ax - dy * r / d;
599 p0.y = ay + dx * r / d;
600 p1.x = ax - dx * r / d;
601 p1.y = ay - dy * r / d;
602 p2.x = ax + dy * r / d;
603 p2.y = ay - dx * r / d;
604 /* Split into 2 triangles at the apex,
605 so that the apex is preciselly covered.
606 Especially important when it is not exactly degenerate. */
607 code = R_fill_triangle_new(pfs, rect, ax, ay, p0.x, p0.y, p1.x, p1.y, t0);
608 if (code < 0)
609 return code;
610 return R_fill_triangle_new(pfs, rect, ax, ay, p1.x, p1.y, p2.x, p2.y, t0);
611 } else {
612 /* Compute the "limit" circle so that its
613 tangent points are outside the rect. */
614 /* Note: this branch is executed when the condition above is false :
615 g0 >= 1e-7 * r0 .
616 We believe that computing this branch with doubles
617 provides enough precision after converting coordinates into 'fixed',
618 and that the limit circle radius is not dramatically big.
619 */
620 double es, er; /* The limit circle parameter, radius. */
621 double ex, ey; /* The limit circle centrum. */
622
623 es = as - as * r / g0; /* Always negative. */
624 er = r * r0 / g0 ;
625 ex = x0 + dx * es;
626 ey = y0 + dy * es;
627 /* Fill the annulus: */
628 code = R_tensor_annulus(pfs, x0, y0, r0, t0, ex, ey, er, t0);
629 if (code < 0)
630 return code;
631 /* Fill entire ending circle to ensure entire rect is covered, but
632 * only if we are filling "inwards" (as otherwise we will overwrite
633 * all the hard work we have done to this point) */
634 if (inwards)
635 code = R_tensor_annulus(pfs, ex, ey, er, t0, ex, ey, 0, t0);
636 return code;
637 }
638 }
639
640 static int
R_tensor_cone_apex(patch_fill_state_t * pfs,const gs_rect * rect,double x0,double y0,double r0,double x1,double y1,double r1,double t)641 R_tensor_cone_apex(patch_fill_state_t *pfs, const gs_rect *rect,
642 double x0, double y0, double r0,
643 double x1, double y1, double r1, double t)
644 {
645 double as = r0 / (r0 - r1);
646 double ax = x0 + (x1 - x0) * as;
647 double ay = y0 + (y1 - y0) * as;
648
649 return R_tensor_annulus(pfs, x1, y1, r1, t, ax, ay, 0, t);
650 }
651
652 /*
653 * A map of this code:
654 *
655 * R_extensions
656 * |-> (R_rect_radius)
657 * |-> (R_outer_circle)
658 * |-> R_obtuse_cone
659 * | |-> R_fill_triangle_new
660 * | | '-> mesh_triangle
661 * | | '-> mesh_triangle_rec <--.
662 * | | |--------------------'
663 * | | |-> small_mesh_triangle
664 * | | | '-> fill_triangle
665 * | | | '-> triangle_by_4 <--.
666 * | | | |----------------'
667 * | | | |-> constant_color_triangle
668 * | | | |-> make_wedge_median (etc)
669 * | | '-----------+--------------------.
670 * | '-------------------. |
671 * |-> R_tensor_cone_apex | |
672 * | '-------------------+ |
673 * '-> R_tensor_annulus <--' \|/
674 * |-> (make_quadrant_arc) |
675 * '-> patch_fill |
676 * |-> fill_patch <--. |
677 * | |-------------' |
678 * | |------------------------------------+
679 * | '-> fill_stripe |
680 * | |-----------------------. |
681 * | \|/ | |
682 * |-> fill_wedges | |
683 * '-> fill_wedges_aux <--. | |
684 * |------------------' \|/ |
685 * |----------------> mesh_padding '
686 * | '----------------------------------.
687 * '-> wedge_by_triangles <--. . |
688 * |---------------------' | |
689 * '-> fill_triangle_wedge <----' |
690 * '-> fill_triangle_wedge_aux |
691 * '-> fill_wedge_trap |
692 * '-> wedge_trap_decompose |
693 * '-> linear_color_trapezoid |
694 * '-> decompose_linear_color <--|
695 * |-------------------------'
696 * '-> constant_color_trapezoid
697 */
698 static int
R_extensions(patch_fill_state_t * pfs,const gs_shading_R_t * psh,const gs_rect * rect,double t0,double t1,bool Extend0,bool Extend1)699 R_extensions(patch_fill_state_t *pfs, const gs_shading_R_t *psh, const gs_rect *rect,
700 double t0, double t1, bool Extend0, bool Extend1)
701 {
702 float x0 = psh->params.Coords[0], y0 = psh->params.Coords[1];
703 floatp r0 = psh->params.Coords[2];
704 float x1 = psh->params.Coords[3], y1 = psh->params.Coords[4];
705 floatp r1 = psh->params.Coords[5];
706 double dx = x1 - x0, dy = y1 - y0, dr = any_abs(r1 - r0);
707 double d = hypot(dx, dy), r;
708 int code;
709
710 if (dr >= d - 1e-7 * (d + dr)) {
711 /* Nested circles, or degenerate. */
712 if (r0 > r1) {
713 if (Extend0) {
714 r = R_rect_radius(rect, x0, y0);
715 if (r > r0) {
716 code = R_tensor_annulus(pfs, x0, y0, r, t0, x0, y0, r0, t0);
717 if (code < 0)
718 return code;
719 }
720 }
721 if (Extend1 && r1 > 0)
722 return R_tensor_annulus(pfs, x1, y1, r1, t1, x1, y1, 0, t1);
723 } else {
724 if (Extend1) {
725 r = R_rect_radius(rect, x1, y1);
726 if (r > r1) {
727 code = R_tensor_annulus(pfs, x1, y1, r, t1, x1, y1, r1, t1);
728 if (code < 0)
729 return code;
730 }
731 }
732 if (Extend0 && r0 > 0)
733 return R_tensor_annulus(pfs, x0, y0, r0, t0, x0, y0, 0, t0);
734 }
735 } else if (dr > d / 3) {
736 /* Obtuse cone. */
737 if (r0 > r1) {
738 if (Extend0) {
739 r = R_rect_radius(rect, x0, y0);
740 code = R_obtuse_cone(pfs, rect, x0, y0, r0, x1, y1, r1, t0, r, true);
741 if (code < 0)
742 return code;
743 }
744 if (Extend1 && r1 != 0)
745 return R_tensor_cone_apex(pfs, rect, x0, y0, r0, x1, y1, r1, t1);
746 return 0;
747 } else {
748 if (Extend1) {
749 r = R_rect_radius(rect, x1, y1);
750 code = R_obtuse_cone(pfs, rect, x1, y1, r1, x0, y0, r0, t1, r, false);
751 if (code < 0)
752 return code;
753 }
754 if (Extend0 && r0 != 0)
755 return R_tensor_cone_apex(pfs, rect, x1, y1, r1, x0, y0, r0, t0);
756 }
757 } else {
758 /* Acute cone or cylinder. */
759 double x2, y2, r2, x3, y3, r3;
760
761 if (Extend0) {
762 code = R_outer_circle(pfs, rect, x1, y1, r1, x0, y0, r0, &x3, &y3, &r3);
763 if (code < 0)
764 return code;
765 if (x3 != x1 || y3 != y1) {
766 code = R_tensor_annulus(pfs, x0, y0, r0, t0, x3, y3, r3, t0);
767 if (code < 0)
768 return code;
769 }
770 }
771 if (Extend1) {
772 code = R_outer_circle(pfs, rect, x0, y0, r0, x1, y1, r1, &x2, &y2, &r2);
773 if (code < 0)
774 return code;
775 if (x2 != x0 || y2 != y0) {
776 code = R_tensor_annulus(pfs, x1, y1, r1, t1, x2, y2, r2, t1);
777 if (code < 0)
778 return code;
779 }
780 }
781 }
782 return 0;
783 }
784
785 static int
R_fill_rect_with_const_color(patch_fill_state_t * pfs,const gs_fixed_rect * clip_rect,float t)786 R_fill_rect_with_const_color(patch_fill_state_t *pfs, const gs_fixed_rect *clip_rect, float t)
787 {
788 #if 0 /* Disabled because the clist writer device doesn't pass
789 the clipping path with fill_recatangle. */
790 patch_color_t pc;
791 const gs_color_space *pcs = pfs->direct_space;
792 gx_device_color dc;
793 int code;
794
795 code = gs_function_evaluate(pfs->Function, &t, pc.cc.paint.values);
796 if (code < 0)
797 return code;
798 pcs->type->restrict_color(&pc.cc, pcs);
799 code = patch_color_to_device_color(pfs, &pc, &dc);
800 if (code < 0)
801 return code;
802 return gx_fill_rectangle_device_rop(fixed2int_pixround(clip_rect->p.x), fixed2int_pixround(clip_rect->p.y),
803 fixed2int_pixround(clip_rect->q.x) - fixed2int_pixround(clip_rect->p.x),
804 fixed2int_pixround(clip_rect->q.y) - fixed2int_pixround(clip_rect->p.y),
805 &dc, pfs->dev, pfs->pis->log_op);
806 #else
807 /* Can't apply fill_rectangle, because the clist writer device doesn't pass
808 the clipping path with fill_recatangle. Convert into trapezoids instead.
809 */
810 quadrangle_patch p;
811 shading_vertex_t pp[2][2];
812 const gs_color_space *pcs = pfs->direct_space;
813 patch_color_t pc;
814 int code;
815
816 code = gs_function_evaluate(pfs->Function, &t, pc.cc.paint.values);
817 if (code < 0)
818 return code;
819 pcs->type->restrict_color(&pc.cc, pcs);
820 pc.t[0] = pc.t[1] = t;
821 pp[0][0].p = clip_rect->p;
822 pp[0][1].p.x = clip_rect->q.x;
823 pp[0][1].p.y = clip_rect->p.y;
824 pp[1][0].p.x = clip_rect->p.x;
825 pp[1][0].p.y = clip_rect->q.y;
826 pp[1][1].p = clip_rect->q;
827 pp[0][0].c = pp[0][1].c = pp[1][0].c = pp[1][1].c = &pc;
828 p.p[0][0] = &pp[0][0];
829 p.p[0][1] = &pp[0][1];
830 p.p[1][0] = &pp[1][0];
831 p.p[1][1] = &pp[1][1];
832 return constant_color_quadrangle(pfs, &p, false);
833 #endif
834 }
835
836 typedef struct radial_shading_attrs_s {
837 double x0, y0;
838 double x1, y1;
839 double span[2][2];
840 double apex;
841 bool have_apex;
842 bool have_root[2]; /* ongoing contact, outgoing contact. */
843 bool outer_contact[2];
844 gs_point p[6]; /* 4 corners of the rectangle, p[4] = p[0], p[5] = p[1] */
845 } radial_shading_attrs_t;
846
847 #define Pw2(a) ((a)*(a))
848
849 static void
radial_shading_external_contact(radial_shading_attrs_t * rsa,int point_index,double t,double r0,double r1,bool at_corner,int root_index)850 radial_shading_external_contact(radial_shading_attrs_t *rsa, int point_index, double t, double r0, double r1, bool at_corner, int root_index)
851 {
852 double cx = rsa->x0 + (rsa->x1 - rsa->x0) * t;
853 double cy = rsa->y0 + (rsa->y1 - rsa->y0) * t;
854 double rx = rsa->p[point_index].x - cx;
855 double ry = rsa->p[point_index].y - cy;
856 double dx = rsa->p[point_index - 1].x - rsa->p[point_index].x;
857 double dy = rsa->p[point_index - 1].y - rsa->p[point_index].y;
858
859 if (at_corner) {
860 double Dx = rsa->p[point_index + 1].x - rsa->p[point_index].x;
861 double Dy = rsa->p[point_index + 1].y - rsa->p[point_index].y;
862 bool b1 = (dx * rx + dy * ry >= 0);
863 bool b2 = (Dx * rx + Dy * ry >= 0);
864
865 if (b1 & b2)
866 rsa->outer_contact[root_index] = true;
867 } else {
868 if (rx * dy - ry * dx < 0)
869 rsa->outer_contact[root_index] = true;
870 }
871 }
872
873 static void
store_roots(radial_shading_attrs_t * rsa,const bool have_root[2],const double t[2],double r0,double r1,int point_index,bool at_corner)874 store_roots(radial_shading_attrs_t *rsa, const bool have_root[2], const double t[2], double r0, double r1, int point_index, bool at_corner)
875 {
876 int i;
877
878 for (i = 0; i < 2; i++) {
879 bool good_root;
880
881 if (!have_root[i])
882 continue;
883 good_root = (!rsa->have_apex || (rsa->apex <= 0 || r0 == 0 ? t[i] >= rsa->apex : t[i] <= rsa->apex));
884 if (good_root) {
885 radial_shading_external_contact(rsa, point_index, t[i], r0, r1, at_corner, i);
886 if (!rsa->have_root[i]) {
887 rsa->span[i][0] = rsa->span[i][1] = t[i];
888 rsa->have_root[i] = true;
889 } else {
890 if (rsa->span[i][0] > t[i])
891 rsa->span[i][0] = t[i];
892 if (rsa->span[i][1] < t[i])
893 rsa->span[i][1] = t[i];
894 }
895 }
896 }
897 }
898
899 static void
compute_radial_shading_span_extended_side(radial_shading_attrs_t * rsa,double r0,double r1,int point_index)900 compute_radial_shading_span_extended_side(radial_shading_attrs_t *rsa, double r0, double r1, int point_index)
901 {
902 double cc, c;
903 bool have_root[2] = {false, false};
904 double t[2];
905 bool by_x = (rsa->p[point_index].x != rsa->p[point_index + 1].x);
906 int i;
907
908 /* As t moves from 0 to 1, the circles move from r0 to r1, and from
909 * from position p0 to py. For simplicity, adjust so that p0 is at
910 * the origin. Consider the projection of the circle drawn at any given
911 * time onto the x axis. The range of points would be:
912 * p1x*t +/- (r0+(r1-r0)*t). We are interested in the first (and last)
913 * moments when the range includes a point c on the x axis. So solve for:
914 * p1x*t +/- (r0+(r1-r0)*t) = c. Let cc = p1x.
915 * So p1x*t0 + (r1-r0)*t0 = c - r0 => t0 = (c - r0)/(p1x + r1 - r0)
916 * p1x*t1 - (r1-r0)*t1 = c + r0 => t1 = (c + r0)/(p1x - r1 + r0)
917 */
918 if (by_x) {
919 c = rsa->p[point_index].x - rsa->x0;
920 cc = rsa->x1 - rsa->x0;
921 } else {
922 c = rsa->p[point_index].y - rsa->y0;
923 cc = rsa->y1 - rsa->y0;
924 }
925 t[0] = (c - r0) / (cc + r1 - r0);
926 t[1] = (c + r0) / (cc - r1 + r0);
927 if (t[0] > t[1]) {
928 c = t[0];
929 t[0] = t[1];
930 t[1] = c;
931 }
932 for (i = 0; i < 2; i++) {
933 double d, d0, d1;
934
935 if (by_x) {
936 d = rsa->y1 - rsa->y0 + r0 + (r1 - r0) * t[i];
937 d0 = rsa->p[point_index].y;
938 d1 = rsa->p[point_index + 1].y;
939 } else {
940 d = rsa->x1 - rsa->x0 + r0 + (r1 - r0) * t[i];
941 d0 = rsa->p[point_index].x;
942 d1 = rsa->p[point_index + 1].x;
943 }
944 if (d1 > d0 ? d0 <= d && d <= d1 : d1 <= d && d <= d0)
945 have_root[i] = true;
946 }
947 store_roots(rsa, have_root, t, r0, r1, point_index, false);
948 }
949
950 static int
compute_radial_shading_span_extended_point(radial_shading_attrs_t * rsa,double r0,double r1,int point_index)951 compute_radial_shading_span_extended_point(radial_shading_attrs_t *rsa, double r0, double r1, int point_index)
952 {
953 /* As t moves from 0 to 1, the circles move from r0 to r1, and from
954 * from position p0 to py. At any given time t, therefore, we
955 * paint the points that are distance r0+(r1-r0)*t from point
956 * (p0x+(p1x-p0x)*t,p0y+(p1y-p0y)*t) = P(t).
957 *
958 * To simplify our algebra, adjust so that (p0x, p0y) is at the origin.
959 * To find the time(s) t at which the a point q is painted, we therefore
960 * solve for t in:
961 *
962 * |q-P(t)| = r0+(r1-r0)*t
963 *
964 * (qx-p1x*t)^2 + (qy-p1y*t)^2 - (r0+(r1-r0)*t)^2 = 0
965 * = qx^2 - 2qx.p1x.t + p1x^2.t^2 + qy^2 - 2qy.p1y.t + p1y^2.t^2 -
966 * (r0^2 + 2r0(r1-r0)t + (r1-r0)^2.t^2)
967 * = qx^2 + qy^2 - r0^2
968 * + -2(qx.p1x + qy.p1y + r0(r1-r0)).t
969 * + (p1x^2 + p1y^2 - (r1-r0)^2).t^2
970 *
971 * So solve using the usual t = (-b +/- SQRT(b^2 - 4ac)) where
972 * a = p1x^2 + p1y^2 - (r1-r0)^2
973 * b = -2(qx.p1x + qy.p1y + r0(r1-r0))
974 * c = qx^2 + qy^2 - r0^2
975 */
976 double p1x = rsa->x1 - rsa->x0;
977 double p1y = rsa->y1 - rsa->y0;
978 double qx = rsa->p[point_index].x - rsa->x0;
979 double qy = rsa->p[point_index].y - rsa->y0;
980 double a = (Pw2(p1x) + Pw2(p1y) - Pw2(r0 - r1));
981 bool have_root[2] = {false, false};
982 double t[2];
983
984 if (fabs(a) < 1e-8) {
985 /* Linear equation. */
986 /* This case is always the ongoing ellipse contact. */
987 double cx = rsa->x0 - (rsa->x1 - rsa->x0) * r0 / (r1 - r0);
988 double cy = rsa->y0 - (rsa->y1 - rsa->y0) * r0 / (r1 - r0);
989
990 t[0] = (Pw2(qx) + Pw2(qy))/(cx*qx + cy*qy) / 2;
991 have_root[0] = true;
992 } else {
993 /* Square equation. No solution if b^2 - 4ac = 0. Equivalently if
994 * (b^2)/4 -a.c = 0 === (b/2)^2 - a.c = 0 === (-b/2)^2 - a.c = 0 */
995 double minushalfb = r0*(r1-r0) + p1x*qx + p1y*qy;
996 double c = Pw2(qx) + Pw2(qy) - Pw2(r0);
997 double desc2 = Pw2(minushalfb) - a*c; /* desc2 = 1/4 (b^2-4ac) */
998
999 if (desc2 < 0) {
1000 return -1; /* The point is outside the shading coverage.
1001 Do not shorten, because we didn't observe it in practice. */
1002 } else {
1003 double desc1 = sqrt(desc2); /* desc1 = 1/2 SQRT(b^2-4ac) */
1004
1005 if (a > 0) {
1006 t[0] = (minushalfb - desc1) / a;
1007 t[1] = (minushalfb + desc1) / a;
1008 } else {
1009 t[0] = (minushalfb + desc1) / a;
1010 t[1] = (minushalfb - desc1) / a;
1011 }
1012 have_root[0] = have_root[1] = true;
1013 }
1014 }
1015 store_roots(rsa, have_root, t, r0, r1, point_index, true);
1016 if (have_root[0] && have_root[1])
1017 return 15;
1018 if (have_root[0])
1019 return 15 - 4;
1020 if (have_root[1])
1021 return 15 - 2;
1022 return -1;
1023 }
1024
1025 #undef Pw2
1026
1027 static int
compute_radial_shading_span_extended(radial_shading_attrs_t * rsa,double r0,double r1)1028 compute_radial_shading_span_extended(radial_shading_attrs_t *rsa, double r0, double r1)
1029 {
1030 int span_type0, span_type1;
1031
1032 span_type0 = compute_radial_shading_span_extended_point(rsa, r0, r1, 1);
1033 if (span_type0 == -1)
1034 return -1;
1035 span_type1 = compute_radial_shading_span_extended_point(rsa, r0, r1, 2);
1036 if (span_type0 != span_type1)
1037 return -1;
1038 span_type1 = compute_radial_shading_span_extended_point(rsa, r0, r1, 3);
1039 if (span_type0 != span_type1)
1040 return -1;
1041 span_type1 = compute_radial_shading_span_extended_point(rsa, r0, r1, 4);
1042 if (span_type0 != span_type1)
1043 return -1;
1044 compute_radial_shading_span_extended_side(rsa, r0, r1, 1);
1045 compute_radial_shading_span_extended_side(rsa, r0, r1, 2);
1046 compute_radial_shading_span_extended_side(rsa, r0, r1, 3);
1047 compute_radial_shading_span_extended_side(rsa, r0, r1, 4);
1048 return span_type0;
1049 }
1050
1051 static int
compute_radial_shading_span(radial_shading_attrs_t * rsa,float x0,float y0,floatp r0,float x1,float y1,floatp r1,const gs_rect * rect)1052 compute_radial_shading_span(radial_shading_attrs_t *rsa, float x0, float y0, floatp r0, float x1, float y1, floatp r1, const gs_rect * rect)
1053 {
1054 /* If the shading area is much larger than the path bbox,
1055 we want to shorten the shading for a faster rendering.
1056 If any point of the path bbox falls outside the shading area,
1057 our math is not applicable, and we render entire shading.
1058 If the path bbox is inside the shading area,
1059 we compute 1 or 2 'spans' - the shading parameter intervals,
1060 which covers the bbox. For doing that we need to resolve
1061 a square eqation by the shading parameter
1062 for each corner of the bounding box,
1063 and for each side of the shading bbox.
1064 Note the equation to be solved in the user space.
1065 Since each equation gives 2 roots (because the points are
1066 strongly inside the shading area), we will get 2 parameter intervals -
1067 the 'lower' one corresponds to the first (ongoing) contact of
1068 the running circle, and the second one corresponds to the last (outgoing) contact
1069 (like in a sun eclipse; well our sun is rectangular).
1070
1071 Here are few exceptions.
1072
1073 First, the equation degenerates when the distance sqrt((x1-x0)^2 + (y1-y0)^2)
1074 appears equal to r0-r1. In this case the base circles do contact,
1075 and the running circle does contact at the same point.
1076 The equation degenerates to a linear one.
1077 Since we don't want float precision noize to affect the result,
1078 we compute this condition in 'fixed' coordinates.
1079
1080 Second, Postscript approximates any circle with 3d order beziers.
1081 This approximation may give a 2% error.
1082 Therefore using the precise roots may cause a dropout.
1083 To prevetn them, we slightly modify the base radii.
1084 However the sign of modification smartly depends
1085 on the relative sizes of the base circles,
1086 and on the contact number. Currently we don't want to
1087 define and debug the smart optimal logic for that,
1088 so we simply try all 4 variants for each source equation,
1089 and use the union of intervals.
1090
1091 Third, we could compute which quarter of the circle
1092 really covers the path bbox. Using it we could skip
1093 rendering of uncovering quarters. Currently we do not
1094 implement this optimization. The general tensor patch algorithm
1095 will skip uncovering parts.
1096
1097 Fourth, when one base circle is (almost) inside the other,
1098 the parameter interval must include the shading apex.
1099 To know that, we determine whether the contacting circle
1100 is outside the rectangle (the "outer" contact),
1101 or it is (partially) inside the rectangle.
1102
1103 At last, a small shortening of a shading won't give a
1104 sensible speedup, but it may replace a symmetric function domain
1105 with an assymmetric one, so that the rendering
1106 would be asymmetyric for a symmetric shading.
1107 Therefore we do not perform a small sortening.
1108 Instead we shorten only if the shading span
1109 is much smaller that the shading domain.
1110 */
1111 const double extent = 1.02;
1112 int span_type0, span_type1, span_type;
1113
1114 memset(rsa, 0, sizeof(*rsa));
1115 rsa->x0 = x0;
1116 rsa->y0 = y0;
1117 rsa->x1 = x1;
1118 rsa->y1 = y1;
1119 rsa->p[0] = rsa->p[4] = rect->p;
1120 rsa->p[1].x = rsa->p[5].x = rect->p.x;
1121 rsa->p[1].y = rsa->p[5].y = rect->q.y;
1122 rsa->p[2] = rect->q;
1123 rsa->p[3].x = rect->q.x;
1124 rsa->p[3].y = rect->p.y;
1125 rsa->have_apex = any_abs(r1 - r0) > 1e-7 * any_abs(r1 + r0);
1126 rsa->apex = (rsa->have_apex ? -r0 / (r1 - r0) : 0);
1127 span_type0 = compute_radial_shading_span_extended(rsa, r0 / extent, r1 * extent);
1128 if (span_type0 == -1)
1129 return -1;
1130 span_type1 = compute_radial_shading_span_extended(rsa, r0 / extent, r1 / extent);
1131 if (span_type0 != span_type1)
1132 return -1;
1133 span_type1 = compute_radial_shading_span_extended(rsa, r0 * extent, r1 * extent);
1134 if (span_type0 != span_type1)
1135 return -1;
1136 span_type1 = compute_radial_shading_span_extended(rsa, r0 * extent, r1 / extent);
1137 if (span_type1 == -1)
1138 return -1;
1139 if (r0 < r1) {
1140 if (rsa->have_root[0] && !rsa->outer_contact[0])
1141 rsa->span[0][0] = rsa->apex; /* Likely never happens. Remove ? */
1142 if (rsa->have_root[1] && !rsa->outer_contact[1])
1143 rsa->span[1][0] = rsa->apex;
1144 } else if (r0 > r1) {
1145 if (rsa->have_root[0] && !rsa->outer_contact[0])
1146 rsa->span[0][1] = rsa->apex;
1147 if (rsa->have_root[1] && !rsa->outer_contact[1])
1148 rsa->span[1][1] = rsa->apex; /* Likely never happens. Remove ? */
1149 }
1150 span_type = 0;
1151 if (rsa->have_root[0] && rsa->span[0][0] < 0)
1152 span_type |= 1;
1153 if (rsa->have_root[1] && rsa->span[1][0] < 0)
1154 span_type |= 1;
1155 if (rsa->have_root[0] && rsa->span[0][1] > 0 && rsa->span[0][0] < 1)
1156 span_type |= 2;
1157 if (rsa->have_root[1] && rsa->span[1][1] > 0 && rsa->span[1][0] < 1)
1158 span_type |= 4;
1159 if (rsa->have_root[0] && rsa->span[0][1] > 1)
1160 span_type |= 8;
1161 if (rsa->have_root[1] && rsa->span[1][1] > 1)
1162 span_type |= 8;
1163 return span_type;
1164 }
1165
1166 static bool
shorten_radial_shading(float * x0,float * y0,floatp * r0,float * d0,float * x1,float * y1,floatp * r1,float * d1,double span_[2])1167 shorten_radial_shading(float *x0, float *y0, floatp *r0, float *d0, float *x1, float *y1, floatp *r1, float *d1, double span_[2])
1168 {
1169 double s0 = span_[0], s1 = span_[1], w;
1170
1171 if (s0 < 0)
1172 s0 = 0;
1173 if (s1 < 0)
1174 s1 = 0;
1175 if (s0 > 1)
1176 s0 = 1;
1177 if (s1 > 1)
1178 s1 = 1;
1179 w = s1 - s0;
1180 if (w == 0)
1181 return false; /* Don't pass a degenerate shading. */
1182 if (w > 0.3)
1183 return false; /* The span is big, don't shorten it. */
1184 { /* Do shorten. */
1185 double R0 = *r0, X0 = *x0, Y0 = *y0, D0 = *d0;
1186 double R1 = *r1, X1 = *x1, Y1 = *y1, D1 = *d1;
1187
1188 *r0 = R0 + (R1 - R0) * s0;
1189 *x0 = X0 + (X1 - X0) * s0;
1190 *y0 = Y0 + (Y1 - Y0) * s0;
1191 *d0 = D0 + (D1 - D0) * s0;
1192 *r1 = R0 + (R1 - R0) * s1;
1193 *x1 = X0 + (X1 - X0) * s1;
1194 *y1 = Y0 + (Y1 - Y0) * s1;
1195 *d1 = D0 + (D1 - D0) * s1;
1196 }
1197 return true;
1198 }
1199
1200 static bool inline
is_radial_shading_large(double x0,double y0,double r0,double x1,double y1,double r1,const gs_rect * rect)1201 is_radial_shading_large(double x0, double y0, double r0, double x1, double y1, double r1, const gs_rect * rect)
1202 {
1203 const double d = hypot(x1 - x0, y1 - y0);
1204 const double area0 = M_PI * r0 * r0 / 2;
1205 const double area1 = M_PI * r1 * r1 / 2;
1206 const double area2 = (r0 + r1) / 2 * d;
1207 const double arbitrary = 8;
1208 double areaX, areaY;
1209
1210 /* The shading area is not equal to area0 + area1 + area2
1211 when one circle is (almost) inside the other.
1212 We believe that the 'arbitrary' coefficient recovers that
1213 when it is set greater than 2. */
1214 /* If one dimension is large enough, the shading parameter span is wide. */
1215 areaX = (rect->q.x - rect->p.x) * (rect->q.x - rect->p.x);
1216 if (areaX * arbitrary < area0 + area1 + area2)
1217 return true;
1218 areaY = (rect->q.y - rect->p.y) * (rect->q.y - rect->p.y);
1219 if (areaY * arbitrary < area0 + area1 + area2)
1220 return true;
1221 return false;
1222 }
1223
1224 static int
gs_shading_R_fill_rectangle_aux(const gs_shading_t * psh0,const gs_rect * rect,const gs_fixed_rect * clip_rect,gx_device * dev,gs_imager_state * pis)1225 gs_shading_R_fill_rectangle_aux(const gs_shading_t * psh0, const gs_rect * rect,
1226 const gs_fixed_rect *clip_rect,
1227 gx_device * dev, gs_imager_state * pis)
1228 {
1229 const gs_shading_R_t *const psh = (const gs_shading_R_t *)psh0;
1230 float d0 = psh->params.Domain[0], d1 = psh->params.Domain[1];
1231 float x0 = psh->params.Coords[0], y0 = psh->params.Coords[1];
1232 floatp r0 = psh->params.Coords[2];
1233 float x1 = psh->params.Coords[3], y1 = psh->params.Coords[4];
1234 floatp r1 = psh->params.Coords[5];
1235 radial_shading_attrs_t rsa;
1236 int span_type; /* <0 - don't shorten, 1 - extent0, 2 - first contact, 4 - last contact, 8 - extent1. */
1237 int code;
1238 patch_fill_state_t pfs1;
1239
1240 if (r0 == 0 && r1 == 0)
1241 return 0; /* PLRM requires to paint nothing. */
1242 code = shade_init_fill_state((shading_fill_state_t *)&pfs1, psh0, dev, pis);
1243 if (code < 0)
1244 return code;
1245 pfs1.Function = psh->params.Function;
1246 code = init_patch_fill_state(&pfs1);
1247 if (code < 0) {
1248 if (pfs1.icclink != NULL) gsicc_release_link(pfs1.icclink);
1249 return code;
1250 }
1251 pfs1.function_arg_shift = 1;
1252 pfs1.rect = *clip_rect;
1253 pfs1.maybe_self_intersecting = false;
1254 if (is_radial_shading_large(x0, y0, r0, x1, y1, r1, rect))
1255 span_type = compute_radial_shading_span(&rsa, x0, y0, r0, x1, y1, r1, rect);
1256 else
1257 span_type = -1;
1258 if (span_type < 0) {
1259 code = R_extensions(&pfs1, psh, rect, d0, d1, psh->params.Extend[0], false);
1260 if (code >= 0)
1261 code = R_tensor_annulus(&pfs1, x0, y0, r0, d0, x1, y1, r1, d1);
1262 if (code >= 0)
1263 code = R_extensions(&pfs1, psh, rect, d0, d1, false, psh->params.Extend[1]);
1264 } else if (span_type == 1) {
1265 code = R_fill_rect_with_const_color(&pfs1, clip_rect, d0);
1266 } else if (span_type == 8) {
1267 code = R_fill_rect_with_const_color(&pfs1, clip_rect, d1);
1268 } else {
1269 bool second_interval = true;
1270
1271 code = 0;
1272 if (span_type & 1)
1273 code = R_extensions(&pfs1, psh, rect, d0, d1, psh->params.Extend[0], false);
1274 if ((code >= 0) && (span_type & 2)) {
1275 float X0 = x0, Y0 = y0, D0 = d0, X1 = x1, Y1 = y1, D1 = d1;
1276 floatp R0 = r0, R1 = r1;
1277
1278 if ((span_type & 4) && rsa.span[0][1] >= rsa.span[1][0]) {
1279 double united[2];
1280
1281 united[0] = rsa.span[0][0];
1282 united[1] = rsa.span[1][1];
1283 shorten_radial_shading(&X0, &Y0, &R0, &D0, &X1, &Y1, &R1, &D1, united);
1284 second_interval = false;
1285 } else {
1286 second_interval = shorten_radial_shading(&X0, &Y0, &R0, &D0, &X1, &Y1, &R1, &D1, rsa.span[0]);
1287 }
1288 code = R_tensor_annulus(&pfs1, X0, Y0, R0, D0, X1, Y1, R1, D1);
1289 }
1290 if (code >= 0 && second_interval) {
1291 if (span_type & 4) {
1292 float X0 = x0, Y0 = y0, D0 = d0, X1 = x1, Y1 = y1, D1 = d1;
1293 floatp R0 = r0, R1 = r1;
1294
1295 shorten_radial_shading(&X0, &Y0, &R0, &D0, &X1, &Y1, &R1, &D1, rsa.span[1]);
1296 code = R_tensor_annulus(&pfs1, X0, Y0, R0, D0, X1, Y1, R1, D1);
1297 }
1298 }
1299 if (code >= 0 && (span_type & 8))
1300 code = R_extensions(&pfs1, psh, rect, d0, d1, false, psh->params.Extend[1]);
1301 }
1302 if (pfs1.icclink != NULL) gsicc_release_link(pfs1.icclink);
1303 if (term_patch_fill_state(&pfs1))
1304 return_error(gs_error_unregistered); /* Must not happen. */
1305 return code;
1306 }
1307
1308 int
gs_shading_R_fill_rectangle(const gs_shading_t * psh0,const gs_rect * rect,const gs_fixed_rect * rect_clip,gx_device * dev,gs_imager_state * pis)1309 gs_shading_R_fill_rectangle(const gs_shading_t * psh0, const gs_rect * rect,
1310 const gs_fixed_rect * rect_clip,
1311 gx_device * dev, gs_imager_state * pis)
1312 {
1313 int code;
1314
1315 if (VD_TRACE_RADIAL_PATCH && vd_allowed('s')) {
1316 vd_get_dc('s');
1317 vd_set_shift(0, 0);
1318 vd_set_scale(0.01);
1319 vd_set_origin(0, 0);
1320 }
1321 code = gs_shading_R_fill_rectangle_aux(psh0, rect, rect_clip, dev, pis);
1322 if (VD_TRACE_FUNCTIONAL_PATCH && vd_allowed('s'))
1323 vd_release_dc;
1324 return code;
1325 }
1326