1 /* cairo - a vector graphics library with display and print output
2  *
3  * Copyright © 2002 University of Southern California
4  *
5  * This library is free software; you can redistribute it and/or
6  * modify it either under the terms of the GNU Lesser General Public
7  * License version 2.1 as published by the Free Software Foundation
8  * (the "LGPL") or, at your option, under the terms of the Mozilla
9  * Public License Version 1.1 (the "MPL"). If you do not alter this
10  * notice, a recipient may use your version of this file under either
11  * the MPL or the LGPL.
12  *
13  * You should have received a copy of the LGPL along with this library
14  * in the file COPYING-LGPL-2.1; if not, write to the Free Software
15  * Foundation, Inc., 51 Franklin Street, Suite 500, Boston, MA 02110-1335, USA
16  * You should have received a copy of the MPL along with this library
17  * in the file COPYING-MPL-1.1
18  *
19  * The contents of this file are subject to the Mozilla Public License
20  * Version 1.1 (the "License"); you may not use this file except in
21  * compliance with the License. You may obtain a copy of the License at
22  * http://www.mozilla.org/MPL/
23  *
24  * This software is distributed on an "AS IS" basis, WITHOUT WARRANTY
25  * OF ANY KIND, either express or implied. See the LGPL or the MPL for
26  * the specific language governing rights and limitations.
27  *
28  * The Original Code is the cairo graphics library.
29  *
30  * The Initial Developer of the Original Code is University of Southern
31  * California.
32  *
33  * Contributor(s):
34  *	Carl D. Worth <cworth@cworth.org>
35  */
36 
37 #include "cairoint.h"
38 #include "cairo-error-private.h"
39 #include <float.h>
40 
41 #define PIXMAN_MAX_INT ((pixman_fixed_1 >> 1) - pixman_fixed_e) /* need to ensure deltas also fit */
42 
43 /**
44  * SECTION:cairo-matrix
45  * @Title: cairo_matrix_t
46  * @Short_Description: Generic matrix operations
47  * @See_Also: #cairo_t
48  *
49  * #cairo_matrix_t is used throughout cairo to convert between different
50  * coordinate spaces.  A #cairo_matrix_t holds an affine transformation,
51  * such as a scale, rotation, shear, or a combination of these.
52  * The transformation of a point (<literal>x</literal>,<literal>y</literal>)
53  * is given by:
54  *
55  * <programlisting>
56  * x_new = xx * x + xy * y + x0;
57  * y_new = yx * x + yy * y + y0;
58  * </programlisting>
59  *
60  * The current transformation matrix of a #cairo_t, represented as a
61  * #cairo_matrix_t, defines the transformation from user-space
62  * coordinates to device-space coordinates. See cairo_get_matrix() and
63  * cairo_set_matrix().
64  **/
65 
66 static void
67 _cairo_matrix_scalar_multiply (cairo_matrix_t *matrix, double scalar);
68 
69 static void
70 _cairo_matrix_compute_adjoint (cairo_matrix_t *matrix);
71 
72 /**
73  * cairo_matrix_init_identity:
74  * @matrix: a #cairo_matrix_t
75  *
76  * Modifies @matrix to be an identity transformation.
77  *
78  * Since: 1.0
79  **/
80 void
cairo_matrix_init_identity(cairo_matrix_t * matrix)81 cairo_matrix_init_identity (cairo_matrix_t *matrix)
82 {
83     cairo_matrix_init (matrix,
84 		       1, 0,
85 		       0, 1,
86 		       0, 0);
87 }
88 slim_hidden_def(cairo_matrix_init_identity);
89 
90 /**
91  * cairo_matrix_init:
92  * @matrix: a #cairo_matrix_t
93  * @xx: xx component of the affine transformation
94  * @yx: yx component of the affine transformation
95  * @xy: xy component of the affine transformation
96  * @yy: yy component of the affine transformation
97  * @x0: X translation component of the affine transformation
98  * @y0: Y translation component of the affine transformation
99  *
100  * Sets @matrix to be the affine transformation given by
101  * @xx, @yx, @xy, @yy, @x0, @y0. The transformation is given
102  * by:
103  * <programlisting>
104  *  x_new = xx * x + xy * y + x0;
105  *  y_new = yx * x + yy * y + y0;
106  * </programlisting>
107  *
108  * Since: 1.0
109  **/
110 void
cairo_matrix_init(cairo_matrix_t * matrix,double xx,double yx,double xy,double yy,double x0,double y0)111 cairo_matrix_init (cairo_matrix_t *matrix,
112 		   double xx, double yx,
113 
114 		   double xy, double yy,
115 		   double x0, double y0)
116 {
117     matrix->xx = xx; matrix->yx = yx;
118     matrix->xy = xy; matrix->yy = yy;
119     matrix->x0 = x0; matrix->y0 = y0;
120 }
121 slim_hidden_def(cairo_matrix_init);
122 
123 /**
124  * _cairo_matrix_get_affine:
125  * @matrix: a #cairo_matrix_t
126  * @xx: location to store xx component of matrix
127  * @yx: location to store yx component of matrix
128  * @xy: location to store xy component of matrix
129  * @yy: location to store yy component of matrix
130  * @x0: location to store x0 (X-translation component) of matrix, or %NULL
131  * @y0: location to store y0 (Y-translation component) of matrix, or %NULL
132  *
133  * Gets the matrix values for the affine transformation that @matrix represents.
134  * See cairo_matrix_init().
135  *
136  *
137  * This function is a leftover from the old public API, but is still
138  * mildly useful as an internal means for getting at the matrix
139  * members in a positional way. For example, when reassigning to some
140  * external matrix type, or when renaming members to more meaningful
141  * names (such as a,b,c,d,e,f) for particular manipulations.
142  **/
143 void
_cairo_matrix_get_affine(const cairo_matrix_t * matrix,double * xx,double * yx,double * xy,double * yy,double * x0,double * y0)144 _cairo_matrix_get_affine (const cairo_matrix_t *matrix,
145 			  double *xx, double *yx,
146 			  double *xy, double *yy,
147 			  double *x0, double *y0)
148 {
149     *xx  = matrix->xx;
150     *yx  = matrix->yx;
151 
152     *xy  = matrix->xy;
153     *yy  = matrix->yy;
154 
155     if (x0)
156 	*x0 = matrix->x0;
157     if (y0)
158 	*y0 = matrix->y0;
159 }
160 
161 /**
162  * cairo_matrix_init_translate:
163  * @matrix: a #cairo_matrix_t
164  * @tx: amount to translate in the X direction
165  * @ty: amount to translate in the Y direction
166  *
167  * Initializes @matrix to a transformation that translates by @tx and
168  * @ty in the X and Y dimensions, respectively.
169  *
170  * Since: 1.0
171  **/
172 void
cairo_matrix_init_translate(cairo_matrix_t * matrix,double tx,double ty)173 cairo_matrix_init_translate (cairo_matrix_t *matrix,
174 			     double tx, double ty)
175 {
176     cairo_matrix_init (matrix,
177 		       1, 0,
178 		       0, 1,
179 		       tx, ty);
180 }
181 slim_hidden_def(cairo_matrix_init_translate);
182 
183 /**
184  * cairo_matrix_translate:
185  * @matrix: a #cairo_matrix_t
186  * @tx: amount to translate in the X direction
187  * @ty: amount to translate in the Y direction
188  *
189  * Applies a translation by @tx, @ty to the transformation in
190  * @matrix. The effect of the new transformation is to first translate
191  * the coordinates by @tx and @ty, then apply the original transformation
192  * to the coordinates.
193  *
194  * Since: 1.0
195  **/
196 void
cairo_matrix_translate(cairo_matrix_t * matrix,double tx,double ty)197 cairo_matrix_translate (cairo_matrix_t *matrix, double tx, double ty)
198 {
199     cairo_matrix_t tmp;
200 
201     cairo_matrix_init_translate (&tmp, tx, ty);
202 
203     cairo_matrix_multiply (matrix, &tmp, matrix);
204 }
205 slim_hidden_def (cairo_matrix_translate);
206 
207 /**
208  * cairo_matrix_init_scale:
209  * @matrix: a #cairo_matrix_t
210  * @sx: scale factor in the X direction
211  * @sy: scale factor in the Y direction
212  *
213  * Initializes @matrix to a transformation that scales by @sx and @sy
214  * in the X and Y dimensions, respectively.
215  *
216  * Since: 1.0
217  **/
218 void
cairo_matrix_init_scale(cairo_matrix_t * matrix,double sx,double sy)219 cairo_matrix_init_scale (cairo_matrix_t *matrix,
220 			 double sx, double sy)
221 {
222     cairo_matrix_init (matrix,
223 		       sx,  0,
224 		       0, sy,
225 		       0, 0);
226 }
227 slim_hidden_def(cairo_matrix_init_scale);
228 
229 /**
230  * cairo_matrix_scale:
231  * @matrix: a #cairo_matrix_t
232  * @sx: scale factor in the X direction
233  * @sy: scale factor in the Y direction
234  *
235  * Applies scaling by @sx, @sy to the transformation in @matrix. The
236  * effect of the new transformation is to first scale the coordinates
237  * by @sx and @sy, then apply the original transformation to the coordinates.
238  *
239  * Since: 1.0
240  **/
241 void
cairo_matrix_scale(cairo_matrix_t * matrix,double sx,double sy)242 cairo_matrix_scale (cairo_matrix_t *matrix, double sx, double sy)
243 {
244     cairo_matrix_t tmp;
245 
246     cairo_matrix_init_scale (&tmp, sx, sy);
247 
248     cairo_matrix_multiply (matrix, &tmp, matrix);
249 }
250 slim_hidden_def(cairo_matrix_scale);
251 
252 /**
253  * cairo_matrix_init_rotate:
254  * @matrix: a #cairo_matrix_t
255  * @radians: angle of rotation, in radians. The direction of rotation
256  * is defined such that positive angles rotate in the direction from
257  * the positive X axis toward the positive Y axis. With the default
258  * axis orientation of cairo, positive angles rotate in a clockwise
259  * direction.
260  *
261  * Initialized @matrix to a transformation that rotates by @radians.
262  *
263  * Since: 1.0
264  **/
265 void
cairo_matrix_init_rotate(cairo_matrix_t * matrix,double radians)266 cairo_matrix_init_rotate (cairo_matrix_t *matrix,
267 			  double radians)
268 {
269     double  s;
270     double  c;
271 
272     s = sin (radians);
273     c = cos (radians);
274 
275     cairo_matrix_init (matrix,
276 		       c, s,
277 		       -s, c,
278 		       0, 0);
279 }
280 slim_hidden_def(cairo_matrix_init_rotate);
281 
282 /**
283  * cairo_matrix_rotate:
284  * @matrix: a #cairo_matrix_t
285  * @radians: angle of rotation, in radians. The direction of rotation
286  * is defined such that positive angles rotate in the direction from
287  * the positive X axis toward the positive Y axis. With the default
288  * axis orientation of cairo, positive angles rotate in a clockwise
289  * direction.
290  *
291  * Applies rotation by @radians to the transformation in
292  * @matrix. The effect of the new transformation is to first rotate the
293  * coordinates by @radians, then apply the original transformation
294  * to the coordinates.
295  *
296  * Since: 1.0
297  **/
298 void
cairo_matrix_rotate(cairo_matrix_t * matrix,double radians)299 cairo_matrix_rotate (cairo_matrix_t *matrix, double radians)
300 {
301     cairo_matrix_t tmp;
302 
303     cairo_matrix_init_rotate (&tmp, radians);
304 
305     cairo_matrix_multiply (matrix, &tmp, matrix);
306 }
307 
308 /**
309  * cairo_matrix_multiply:
310  * @result: a #cairo_matrix_t in which to store the result
311  * @a: a #cairo_matrix_t
312  * @b: a #cairo_matrix_t
313  *
314  * Multiplies the affine transformations in @a and @b together
315  * and stores the result in @result. The effect of the resulting
316  * transformation is to first apply the transformation in @a to the
317  * coordinates and then apply the transformation in @b to the
318  * coordinates.
319  *
320  * It is allowable for @result to be identical to either @a or @b.
321  *
322  * Since: 1.0
323  **/
324 /*
325  * XXX: The ordering of the arguments to this function corresponds
326  *      to [row_vector]*A*B. If we want to use column vectors instead,
327  *      then we need to switch the two arguments and fix up all
328  *      uses.
329  */
330 void
cairo_matrix_multiply(cairo_matrix_t * result,const cairo_matrix_t * a,const cairo_matrix_t * b)331 cairo_matrix_multiply (cairo_matrix_t *result, const cairo_matrix_t *a, const cairo_matrix_t *b)
332 {
333     cairo_matrix_t r;
334 
335     r.xx = a->xx * b->xx + a->yx * b->xy;
336     r.yx = a->xx * b->yx + a->yx * b->yy;
337 
338     r.xy = a->xy * b->xx + a->yy * b->xy;
339     r.yy = a->xy * b->yx + a->yy * b->yy;
340 
341     r.x0 = a->x0 * b->xx + a->y0 * b->xy + b->x0;
342     r.y0 = a->x0 * b->yx + a->y0 * b->yy + b->y0;
343 
344     *result = r;
345 }
346 slim_hidden_def(cairo_matrix_multiply);
347 
348 void
_cairo_matrix_multiply(cairo_matrix_t * r,const cairo_matrix_t * a,const cairo_matrix_t * b)349 _cairo_matrix_multiply (cairo_matrix_t *r,
350 			const cairo_matrix_t *a,
351 			const cairo_matrix_t *b)
352 {
353     r->xx = a->xx * b->xx + a->yx * b->xy;
354     r->yx = a->xx * b->yx + a->yx * b->yy;
355 
356     r->xy = a->xy * b->xx + a->yy * b->xy;
357     r->yy = a->xy * b->yx + a->yy * b->yy;
358 
359     r->x0 = a->x0 * b->xx + a->y0 * b->xy + b->x0;
360     r->y0 = a->x0 * b->yx + a->y0 * b->yy + b->y0;
361 }
362 
363 /**
364  * cairo_matrix_transform_distance:
365  * @matrix: a #cairo_matrix_t
366  * @dx: X component of a distance vector. An in/out parameter
367  * @dy: Y component of a distance vector. An in/out parameter
368  *
369  * Transforms the distance vector (@dx,@dy) by @matrix. This is
370  * similar to cairo_matrix_transform_point() except that the translation
371  * components of the transformation are ignored. The calculation of
372  * the returned vector is as follows:
373  *
374  * <programlisting>
375  * dx2 = dx1 * a + dy1 * c;
376  * dy2 = dx1 * b + dy1 * d;
377  * </programlisting>
378  *
379  * Affine transformations are position invariant, so the same vector
380  * always transforms to the same vector. If (@x1,@y1) transforms
381  * to (@x2,@y2) then (@x1+@dx1,@y1+@dy1) will transform to
382  * (@x1+@dx2,@y1+@dy2) for all values of @x1 and @x2.
383  *
384  * Since: 1.0
385  **/
386 void
cairo_matrix_transform_distance(const cairo_matrix_t * matrix,double * dx,double * dy)387 cairo_matrix_transform_distance (const cairo_matrix_t *matrix, double *dx, double *dy)
388 {
389     double new_x, new_y;
390 
391     new_x = (matrix->xx * *dx + matrix->xy * *dy);
392     new_y = (matrix->yx * *dx + matrix->yy * *dy);
393 
394     *dx = new_x;
395     *dy = new_y;
396 }
397 slim_hidden_def(cairo_matrix_transform_distance);
398 
399 /**
400  * cairo_matrix_transform_point:
401  * @matrix: a #cairo_matrix_t
402  * @x: X position. An in/out parameter
403  * @y: Y position. An in/out parameter
404  *
405  * Transforms the point (@x, @y) by @matrix.
406  *
407  * Since: 1.0
408  **/
409 void
cairo_matrix_transform_point(const cairo_matrix_t * matrix,double * x,double * y)410 cairo_matrix_transform_point (const cairo_matrix_t *matrix, double *x, double *y)
411 {
412     cairo_matrix_transform_distance (matrix, x, y);
413 
414     *x += matrix->x0;
415     *y += matrix->y0;
416 }
417 slim_hidden_def(cairo_matrix_transform_point);
418 
419 void
_cairo_matrix_transform_bounding_box(const cairo_matrix_t * matrix,double * x1,double * y1,double * x2,double * y2,cairo_bool_t * is_tight)420 _cairo_matrix_transform_bounding_box (const cairo_matrix_t *matrix,
421 				      double *x1, double *y1,
422 				      double *x2, double *y2,
423 				      cairo_bool_t *is_tight)
424 {
425     int i;
426     double quad_x[4], quad_y[4];
427     double min_x, max_x;
428     double min_y, max_y;
429 
430     if (matrix->xy == 0. && matrix->yx == 0.) {
431 	/* non-rotation/skew matrix, just map the two extreme points */
432 
433 	if (matrix->xx != 1.) {
434 	    quad_x[0] = *x1 * matrix->xx;
435 	    quad_x[1] = *x2 * matrix->xx;
436 	    if (quad_x[0] < quad_x[1]) {
437 		*x1 = quad_x[0];
438 		*x2 = quad_x[1];
439 	    } else {
440 		*x1 = quad_x[1];
441 		*x2 = quad_x[0];
442 	    }
443 	}
444 	if (matrix->x0 != 0.) {
445 	    *x1 += matrix->x0;
446 	    *x2 += matrix->x0;
447 	}
448 
449 	if (matrix->yy != 1.) {
450 	    quad_y[0] = *y1 * matrix->yy;
451 	    quad_y[1] = *y2 * matrix->yy;
452 	    if (quad_y[0] < quad_y[1]) {
453 		*y1 = quad_y[0];
454 		*y2 = quad_y[1];
455 	    } else {
456 		*y1 = quad_y[1];
457 		*y2 = quad_y[0];
458 	    }
459 	}
460 	if (matrix->y0 != 0.) {
461 	    *y1 += matrix->y0;
462 	    *y2 += matrix->y0;
463 	}
464 
465 	if (is_tight)
466 	    *is_tight = TRUE;
467 
468 	return;
469     }
470 
471     /* general matrix */
472     quad_x[0] = *x1;
473     quad_y[0] = *y1;
474     cairo_matrix_transform_point (matrix, &quad_x[0], &quad_y[0]);
475 
476     quad_x[1] = *x2;
477     quad_y[1] = *y1;
478     cairo_matrix_transform_point (matrix, &quad_x[1], &quad_y[1]);
479 
480     quad_x[2] = *x1;
481     quad_y[2] = *y2;
482     cairo_matrix_transform_point (matrix, &quad_x[2], &quad_y[2]);
483 
484     quad_x[3] = *x2;
485     quad_y[3] = *y2;
486     cairo_matrix_transform_point (matrix, &quad_x[3], &quad_y[3]);
487 
488     min_x = max_x = quad_x[0];
489     min_y = max_y = quad_y[0];
490 
491     for (i=1; i < 4; i++) {
492 	if (quad_x[i] < min_x)
493 	    min_x = quad_x[i];
494 	if (quad_x[i] > max_x)
495 	    max_x = quad_x[i];
496 
497 	if (quad_y[i] < min_y)
498 	    min_y = quad_y[i];
499 	if (quad_y[i] > max_y)
500 	    max_y = quad_y[i];
501     }
502 
503     *x1 = min_x;
504     *y1 = min_y;
505     *x2 = max_x;
506     *y2 = max_y;
507 
508     if (is_tight) {
509         /* it's tight if and only if the four corner points form an axis-aligned
510            rectangle.
511            And that's true if and only if we can derive corners 0 and 3 from
512            corners 1 and 2 in one of two straightforward ways...
513            We could use a tolerance here but for now we'll fall back to FALSE in the case
514            of floating point error.
515         */
516         *is_tight =
517             (quad_x[1] == quad_x[0] && quad_y[1] == quad_y[3] &&
518              quad_x[2] == quad_x[3] && quad_y[2] == quad_y[0]) ||
519             (quad_x[1] == quad_x[3] && quad_y[1] == quad_y[0] &&
520              quad_x[2] == quad_x[0] && quad_y[2] == quad_y[3]);
521     }
522 }
523 
524 cairo_private void
_cairo_matrix_transform_bounding_box_fixed(const cairo_matrix_t * matrix,cairo_box_t * bbox,cairo_bool_t * is_tight)525 _cairo_matrix_transform_bounding_box_fixed (const cairo_matrix_t *matrix,
526 					    cairo_box_t          *bbox,
527 					    cairo_bool_t *is_tight)
528 {
529     double x1, y1, x2, y2;
530 
531     _cairo_box_to_doubles (bbox, &x1, &y1, &x2, &y2);
532     _cairo_matrix_transform_bounding_box (matrix, &x1, &y1, &x2, &y2, is_tight);
533     _cairo_box_from_doubles (bbox, &x1, &y1, &x2, &y2);
534 }
535 
536 static void
_cairo_matrix_scalar_multiply(cairo_matrix_t * matrix,double scalar)537 _cairo_matrix_scalar_multiply (cairo_matrix_t *matrix, double scalar)
538 {
539     matrix->xx *= scalar;
540     matrix->yx *= scalar;
541 
542     matrix->xy *= scalar;
543     matrix->yy *= scalar;
544 
545     matrix->x0 *= scalar;
546     matrix->y0 *= scalar;
547 }
548 
549 /* This function isn't a correct adjoint in that the implicit 1 in the
550    homogeneous result should actually be ad-bc instead. But, since this
551    adjoint is only used in the computation of the inverse, which
552    divides by det (A)=ad-bc anyway, everything works out in the end. */
553 static void
_cairo_matrix_compute_adjoint(cairo_matrix_t * matrix)554 _cairo_matrix_compute_adjoint (cairo_matrix_t *matrix)
555 {
556     /* adj (A) = transpose (C:cofactor (A,i,j)) */
557     double a, b, c, d, tx, ty;
558 
559     _cairo_matrix_get_affine (matrix,
560 			      &a,  &b,
561 			      &c,  &d,
562 			      &tx, &ty);
563 
564     cairo_matrix_init (matrix,
565 		       d, -b,
566 		       -c, a,
567 		       c*ty - d*tx, b*tx - a*ty);
568 }
569 
570 /**
571  * cairo_matrix_invert:
572  * @matrix: a #cairo_matrix_t
573  *
574  * Changes @matrix to be the inverse of its original value. Not
575  * all transformation matrices have inverses; if the matrix
576  * collapses points together (it is <firstterm>degenerate</firstterm>),
577  * then it has no inverse and this function will fail.
578  *
579  * Returns: If @matrix has an inverse, modifies @matrix to
580  *  be the inverse matrix and returns %CAIRO_STATUS_SUCCESS. Otherwise,
581  *  returns %CAIRO_STATUS_INVALID_MATRIX.
582  *
583  * Since: 1.0
584  **/
585 cairo_status_t
cairo_matrix_invert(cairo_matrix_t * matrix)586 cairo_matrix_invert (cairo_matrix_t *matrix)
587 {
588     double det;
589 
590     /* Simple scaling|translation matrices are quite common... */
591     if (matrix->xy == 0. && matrix->yx == 0.) {
592 	matrix->x0 = -matrix->x0;
593 	matrix->y0 = -matrix->y0;
594 
595 	if (matrix->xx != 1.) {
596 	    if (matrix->xx == 0.)
597 		return _cairo_error (CAIRO_STATUS_INVALID_MATRIX);
598 
599 	    matrix->xx = 1. / matrix->xx;
600 	    matrix->x0 *= matrix->xx;
601 	}
602 
603 	if (matrix->yy != 1.) {
604 	    if (matrix->yy == 0.)
605 		return _cairo_error (CAIRO_STATUS_INVALID_MATRIX);
606 
607 	    matrix->yy = 1. / matrix->yy;
608 	    matrix->y0 *= matrix->yy;
609 	}
610 
611 	return CAIRO_STATUS_SUCCESS;
612     }
613 
614     /* inv (A) = 1/det (A) * adj (A) */
615     det = _cairo_matrix_compute_determinant (matrix);
616 
617     if (! ISFINITE (det))
618 	return _cairo_error (CAIRO_STATUS_INVALID_MATRIX);
619 
620     if (det == 0)
621 	return _cairo_error (CAIRO_STATUS_INVALID_MATRIX);
622 
623     _cairo_matrix_compute_adjoint (matrix);
624     _cairo_matrix_scalar_multiply (matrix, 1 / det);
625 
626     return CAIRO_STATUS_SUCCESS;
627 }
628 slim_hidden_def(cairo_matrix_invert);
629 
630 cairo_bool_t
_cairo_matrix_is_invertible(const cairo_matrix_t * matrix)631 _cairo_matrix_is_invertible (const cairo_matrix_t *matrix)
632 {
633     double det;
634 
635     det = _cairo_matrix_compute_determinant (matrix);
636 
637     return ISFINITE (det) && det != 0.;
638 }
639 
640 cairo_bool_t
_cairo_matrix_is_scale_0(const cairo_matrix_t * matrix)641 _cairo_matrix_is_scale_0 (const cairo_matrix_t *matrix)
642 {
643     return matrix->xx == 0. &&
644            matrix->xy == 0. &&
645            matrix->yx == 0. &&
646            matrix->yy == 0.;
647 }
648 
649 double
_cairo_matrix_compute_determinant(const cairo_matrix_t * matrix)650 _cairo_matrix_compute_determinant (const cairo_matrix_t *matrix)
651 {
652     double a, b, c, d;
653 
654     a = matrix->xx; b = matrix->yx;
655     c = matrix->xy; d = matrix->yy;
656 
657     return a*d - b*c;
658 }
659 
660 /**
661  * _cairo_matrix_compute_basis_scale_factors:
662  * @matrix: a matrix
663  * @basis_scale: the scale factor in the direction of basis
664  * @normal_scale: the scale factor in the direction normal to the basis
665  * @x_basis: basis to use.  X basis if true, Y basis otherwise.
666  *
667  * Computes |Mv| and det(M)/|Mv| for v=[1,0] if x_basis is true, and v=[0,1]
668  * otherwise, and M is @matrix.
669  *
670  * Return value: the scale factor of @matrix on the height of the font,
671  * or 1.0 if @matrix is %NULL.
672  **/
673 cairo_status_t
_cairo_matrix_compute_basis_scale_factors(const cairo_matrix_t * matrix,double * basis_scale,double * normal_scale,cairo_bool_t x_basis)674 _cairo_matrix_compute_basis_scale_factors (const cairo_matrix_t *matrix,
675 					   double *basis_scale, double *normal_scale,
676 					   cairo_bool_t x_basis)
677 {
678     double det;
679 
680     det = _cairo_matrix_compute_determinant (matrix);
681 
682     if (! ISFINITE (det))
683 	return _cairo_error (CAIRO_STATUS_INVALID_MATRIX);
684 
685     if (det == 0)
686     {
687 	*basis_scale = *normal_scale = 0;
688     }
689     else
690     {
691 	double x = x_basis != 0;
692 	double y = x == 0;
693 	double major, minor;
694 
695 	cairo_matrix_transform_distance (matrix, &x, &y);
696 	major = hypot (x, y);
697 	/*
698 	 * ignore mirroring
699 	 */
700 	if (det < 0)
701 	    det = -det;
702 	if (major)
703 	    minor = det / major;
704 	else
705 	    minor = 0.0;
706 	if (x_basis)
707 	{
708 	    *basis_scale = major;
709 	    *normal_scale = minor;
710 	}
711 	else
712 	{
713 	    *basis_scale = minor;
714 	    *normal_scale = major;
715 	}
716     }
717 
718     return CAIRO_STATUS_SUCCESS;
719 }
720 
721 cairo_bool_t
_cairo_matrix_is_integer_translation(const cairo_matrix_t * matrix,int * itx,int * ity)722 _cairo_matrix_is_integer_translation (const cairo_matrix_t *matrix,
723 				      int *itx, int *ity)
724 {
725     if (_cairo_matrix_is_translation (matrix))
726     {
727         cairo_fixed_t x0_fixed = _cairo_fixed_from_double (matrix->x0);
728         cairo_fixed_t y0_fixed = _cairo_fixed_from_double (matrix->y0);
729 
730         if (_cairo_fixed_is_integer (x0_fixed) &&
731             _cairo_fixed_is_integer (y0_fixed))
732         {
733             if (itx)
734                 *itx = _cairo_fixed_integer_part (x0_fixed);
735             if (ity)
736                 *ity = _cairo_fixed_integer_part (y0_fixed);
737 
738             return TRUE;
739         }
740     }
741 
742     return FALSE;
743 }
744 
745 #define SCALING_EPSILON _cairo_fixed_to_double(1)
746 
747 /* This only returns true if the matrix is 90 degree rotations or
748  * flips. It appears calling code is relying on this. It will return
749  * false for other rotations even if the scale is one. Approximations
750  * are allowed to handle matricies filled in using trig functions
751  * such as sin(M_PI_2).
752  */
753 cairo_bool_t
_cairo_matrix_has_unity_scale(const cairo_matrix_t * matrix)754 _cairo_matrix_has_unity_scale (const cairo_matrix_t *matrix)
755 {
756     /* check that the determinant is near +/-1 */
757     double det = _cairo_matrix_compute_determinant (matrix);
758     if (fabs (det * det - 1.0) < SCALING_EPSILON) {
759 	/* check that one axis is close to zero */
760 	if (fabs (matrix->xy) < SCALING_EPSILON  &&
761 	    fabs (matrix->yx) < SCALING_EPSILON)
762 	    return TRUE;
763 	if (fabs (matrix->xx) < SCALING_EPSILON  &&
764 	    fabs (matrix->yy) < SCALING_EPSILON)
765 	    return TRUE;
766 	/* If rotations are allowed then it must instead test for
767 	 * orthogonality. This is xx*xy+yx*yy ~= 0.
768 	 */
769     }
770     return FALSE;
771 }
772 
773 /* By pixel exact here, we mean a matrix that is composed only of
774  * 90 degree rotations, flips, and integer translations and produces a 1:1
775  * mapping between source and destination pixels. If we transform an image
776  * with a pixel-exact matrix, filtering is not useful.
777  */
778 cairo_bool_t
_cairo_matrix_is_pixel_exact(const cairo_matrix_t * matrix)779 _cairo_matrix_is_pixel_exact (const cairo_matrix_t *matrix)
780 {
781     cairo_fixed_t x0_fixed, y0_fixed;
782 
783     if (! _cairo_matrix_has_unity_scale (matrix))
784 	return FALSE;
785 
786     x0_fixed = _cairo_fixed_from_double (matrix->x0);
787     y0_fixed = _cairo_fixed_from_double (matrix->y0);
788 
789     return _cairo_fixed_is_integer (x0_fixed) && _cairo_fixed_is_integer (y0_fixed);
790 }
791 
792 /*
793   A circle in user space is transformed into an ellipse in device space.
794 
795   The following is a derivation of a formula to calculate the length of the
796   major axis for this ellipse; this is useful for error bounds calculations.
797 
798   Thanks to Walter Brisken <wbrisken@aoc.nrao.edu> for this derivation:
799 
800   1.  First some notation:
801 
802   All capital letters represent vectors in two dimensions.  A prime '
803   represents a transformed coordinate.  Matrices are written in underlined
804   form, ie _R_.  Lowercase letters represent scalar real values.
805 
806   2.  The question has been posed:  What is the maximum expansion factor
807   achieved by the linear transformation
808 
809   X' = X _R_
810 
811   where _R_ is a real-valued 2x2 matrix with entries:
812 
813   _R_ = [a b]
814         [c d]  .
815 
816   In other words, what is the maximum radius, MAX[ |X'| ], reached for any
817   X on the unit circle ( |X| = 1 ) ?
818 
819   3.  Some useful formulae
820 
821   (A) through (C) below are standard double-angle formulae.  (D) is a lesser
822   known result and is derived below:
823 
824   (A)  sin²(θ) = (1 - cos(2*θ))/2
825   (B)  cos²(θ) = (1 + cos(2*θ))/2
826   (C)  sin(θ)*cos(θ) = sin(2*θ)/2
827   (D)  MAX[a*cos(θ) + b*sin(θ)] = sqrt(a² + b²)
828 
829   Proof of (D):
830 
831   find the maximum of the function by setting the derivative to zero:
832 
833        -a*sin(θ)+b*cos(θ) = 0
834 
835   From this it follows that
836 
837        tan(θ) = b/a
838 
839   and hence
840 
841        sin(θ) = b/sqrt(a² + b²)
842 
843   and
844 
845        cos(θ) = a/sqrt(a² + b²)
846 
847   Thus the maximum value is
848 
849        MAX[a*cos(θ) + b*sin(θ)] = (a² + b²)/sqrt(a² + b²)
850                                    = sqrt(a² + b²)
851 
852   4.  Derivation of maximum expansion
853 
854   To find MAX[ |X'| ] we search brute force method using calculus.  The unit
855   circle on which X is constrained is to be parameterized by t:
856 
857        X(θ) = (cos(θ), sin(θ))
858 
859   Thus
860 
861        X'(θ) = X(θ) * _R_ = (cos(θ), sin(θ)) * [a b]
862                                                [c d]
863              = (a*cos(θ) + c*sin(θ), b*cos(θ) + d*sin(θ)).
864 
865   Define
866 
867        r(θ) = |X'(θ)|
868 
869   Thus
870 
871        r²(θ) = (a*cos(θ) + c*sin(θ))² + (b*cos(θ) + d*sin(θ))²
872              = (a² + b²)*cos²(θ) + (c² + d²)*sin²(θ)
873                  + 2*(a*c + b*d)*cos(θ)*sin(θ)
874 
875   Now apply the double angle formulae (A) to (C) from above:
876 
877        r²(θ) = (a² + b² + c² + d²)/2
878 	     + (a² + b² - c² - d²)*cos(2*θ)/2
879   	     + (a*c + b*d)*sin(2*θ)
880              = f + g*cos(φ) + h*sin(φ)
881 
882   Where
883 
884        f = (a² + b² + c² + d²)/2
885        g = (a² + b² - c² - d²)/2
886        h = (a*c + d*d)
887        φ = 2*θ
888 
889   It is clear that MAX[ |X'| ] = sqrt(MAX[ r² ]).  Here we determine MAX[ r² ]
890   using (D) from above:
891 
892        MAX[ r² ] = f + sqrt(g² + h²)
893 
894   And finally
895 
896        MAX[ |X'| ] = sqrt( f + sqrt(g² + h²) )
897 
898   Which is the solution to this problem.
899 
900   Walter Brisken
901   2004/10/08
902 
903   (Note that the minor axis length is at the minimum of the above solution,
904   which is just sqrt ( f - sqrt(g² + h²) ) given the symmetry of (D)).
905 
906 
907   For another derivation of the same result, using Singular Value Decomposition,
908   see doc/tutorial/src/singular.c.
909 */
910 
911 /* determine the length of the major axis of a circle of the given radius
912    after applying the transformation matrix. */
913 double
_cairo_matrix_transformed_circle_major_axis(const cairo_matrix_t * matrix,double radius)914 _cairo_matrix_transformed_circle_major_axis (const cairo_matrix_t *matrix,
915 					     double radius)
916 {
917     double  a, b, c, d, f, g, h, i, j;
918 
919     if (_cairo_matrix_has_unity_scale (matrix))
920 	return radius;
921 
922     _cairo_matrix_get_affine (matrix,
923                               &a, &b,
924                               &c, &d,
925                               NULL, NULL);
926 
927     i = a*a + b*b;
928     j = c*c + d*d;
929 
930     f = 0.5 * (i + j);
931     g = 0.5 * (i - j);
932     h = a*c + b*d;
933 
934     return radius * sqrt (f + hypot (g, h));
935 
936     /*
937      * we don't need the minor axis length, which is
938      * double min = radius * sqrt (f - sqrt (g*g+h*h));
939      */
940 }
941 
942 static const pixman_transform_t pixman_identity_transform = {{
943         {1 << 16,        0,       0},
944         {       0, 1 << 16,       0},
945         {       0,       0, 1 << 16}
946     }};
947 
948 static cairo_status_t
_cairo_matrix_to_pixman_matrix(const cairo_matrix_t * matrix,pixman_transform_t * pixman_transform,double xc,double yc)949 _cairo_matrix_to_pixman_matrix (const cairo_matrix_t	*matrix,
950 				pixman_transform_t	*pixman_transform,
951 				double xc,
952 				double yc)
953 {
954     cairo_matrix_t inv;
955     unsigned max_iterations;
956 
957     pixman_transform->matrix[0][0] = _cairo_fixed_16_16_from_double (matrix->xx);
958     pixman_transform->matrix[0][1] = _cairo_fixed_16_16_from_double (matrix->xy);
959     pixman_transform->matrix[0][2] = _cairo_fixed_16_16_from_double (matrix->x0);
960 
961     pixman_transform->matrix[1][0] = _cairo_fixed_16_16_from_double (matrix->yx);
962     pixman_transform->matrix[1][1] = _cairo_fixed_16_16_from_double (matrix->yy);
963     pixman_transform->matrix[1][2] = _cairo_fixed_16_16_from_double (matrix->y0);
964 
965     pixman_transform->matrix[2][0] = 0;
966     pixman_transform->matrix[2][1] = 0;
967     pixman_transform->matrix[2][2] = 1 << 16;
968 
969     /* The conversion above breaks cairo's translation invariance:
970      * a translation of (a, b) in device space translates to
971      * a translation of (xx * a + xy * b, yx * a + yy * b)
972      * for cairo, while pixman uses rounded versions of xx ... yy.
973      * This error increases as a and b get larger.
974      *
975      * To compensate for this, we fix the point (xc, yc) in pattern
976      * space and adjust pixman's transform to agree with cairo's at
977      * that point.
978      */
979 
980     if (_cairo_matrix_has_unity_scale (matrix))
981 	return CAIRO_STATUS_SUCCESS;
982 
983     if (unlikely (fabs (matrix->xx) > PIXMAN_MAX_INT ||
984 		  fabs (matrix->xy) > PIXMAN_MAX_INT ||
985 		  fabs (matrix->x0) > PIXMAN_MAX_INT ||
986 		  fabs (matrix->yx) > PIXMAN_MAX_INT ||
987 		  fabs (matrix->yy) > PIXMAN_MAX_INT ||
988 		  fabs (matrix->y0) > PIXMAN_MAX_INT))
989     {
990 	return _cairo_error (CAIRO_STATUS_INVALID_MATRIX);
991     }
992 
993     /* Note: If we can't invert the transformation, skip the adjustment. */
994     inv = *matrix;
995     if (cairo_matrix_invert (&inv) != CAIRO_STATUS_SUCCESS)
996 	return CAIRO_STATUS_SUCCESS;
997 
998     /* find the pattern space coordinate that maps to (xc, yc) */
999     max_iterations = 5;
1000     do {
1001 	double x,y;
1002 	pixman_vector_t vector;
1003 	cairo_fixed_16_16_t dx, dy;
1004 
1005 	vector.vector[0] = _cairo_fixed_16_16_from_double (xc);
1006 	vector.vector[1] = _cairo_fixed_16_16_from_double (yc);
1007 	vector.vector[2] = 1 << 16;
1008 
1009 	/* If we can't transform the reference point, skip the adjustment. */
1010 	if (! pixman_transform_point_3d (pixman_transform, &vector))
1011 	    return CAIRO_STATUS_SUCCESS;
1012 
1013 	x = pixman_fixed_to_double (vector.vector[0]);
1014 	y = pixman_fixed_to_double (vector.vector[1]);
1015 	cairo_matrix_transform_point (&inv, &x, &y);
1016 
1017 	/* Ideally, the vector should now be (xc, yc).
1018 	 * We can now compensate for the resulting error.
1019 	 */
1020 	x -= xc;
1021 	y -= yc;
1022 	cairo_matrix_transform_distance (matrix, &x, &y);
1023 	dx = _cairo_fixed_16_16_from_double (x);
1024 	dy = _cairo_fixed_16_16_from_double (y);
1025 	pixman_transform->matrix[0][2] -= dx;
1026 	pixman_transform->matrix[1][2] -= dy;
1027 
1028 	if (dx == 0 && dy == 0)
1029 	    return CAIRO_STATUS_SUCCESS;
1030     } while (--max_iterations);
1031 
1032     /* We didn't find an exact match between cairo and pixman, but
1033      * the matrix should be mostly correct */
1034     return CAIRO_STATUS_SUCCESS;
1035 }
1036 
1037 static inline double
_pixman_nearest_sample(double d)1038 _pixman_nearest_sample (double d)
1039 {
1040     return ceil (d - .5);
1041 }
1042 
1043 /**
1044  * _cairo_matrix_is_pixman_translation:
1045  * @matrix: a matrix
1046  * @filter: the filter to be used on the pattern transformed by @matrix
1047  * @x_offset: the translation in the X direction
1048  * @y_offset: the translation in the Y direction
1049  *
1050  * Checks if @matrix translated by (x_offset, y_offset) can be
1051  * represented using just an offset (within the range pixman can
1052  * accept) and an identity matrix.
1053  *
1054  * Passing a non-zero value in x_offset/y_offset has the same effect
1055  * as applying cairo_matrix_translate(matrix, x_offset, y_offset) and
1056  * setting x_offset and y_offset to 0.
1057  *
1058  * Upon return x_offset and y_offset contain the translation vector if
1059  * the return value is %TRUE. If the return value is %FALSE, they will
1060  * not be modified.
1061  *
1062  * Return value: %TRUE if @matrix can be represented as a pixman
1063  * translation, %FALSE otherwise.
1064  **/
1065 cairo_bool_t
_cairo_matrix_is_pixman_translation(const cairo_matrix_t * matrix,cairo_filter_t filter,int * x_offset,int * y_offset)1066 _cairo_matrix_is_pixman_translation (const cairo_matrix_t     *matrix,
1067 				     cairo_filter_t            filter,
1068 				     int                      *x_offset,
1069 				     int                      *y_offset)
1070 {
1071     double tx, ty;
1072 
1073     if (!_cairo_matrix_is_translation (matrix))
1074 	return FALSE;
1075 
1076     if (matrix->x0 == 0. && matrix->y0 == 0.)
1077 	return TRUE;
1078 
1079     tx = matrix->x0 + *x_offset;
1080     ty = matrix->y0 + *y_offset;
1081 
1082     if (filter == CAIRO_FILTER_FAST || filter == CAIRO_FILTER_NEAREST) {
1083 	tx = _pixman_nearest_sample (tx);
1084 	ty = _pixman_nearest_sample (ty);
1085     } else if (tx != floor (tx) || ty != floor (ty)) {
1086 	return FALSE;
1087     }
1088 
1089     if (fabs (tx) > PIXMAN_MAX_INT || fabs (ty) > PIXMAN_MAX_INT)
1090 	return FALSE;
1091 
1092     *x_offset = _cairo_lround (tx);
1093     *y_offset = _cairo_lround (ty);
1094     return TRUE;
1095 }
1096 
1097 /**
1098  * _cairo_matrix_to_pixman_matrix_offset:
1099  * @matrix: a matrix
1100  * @filter: the filter to be used on the pattern transformed by @matrix
1101  * @xc: the X coordinate of the point to fix in pattern space
1102  * @yc: the Y coordinate of the point to fix in pattern space
1103  * @out_transform: the transformation which best approximates @matrix
1104  * @x_offset: the translation in the X direction
1105  * @y_offset: the translation in the Y direction
1106  *
1107  * This function tries to represent @matrix translated by (x_offset,
1108  * y_offset) as a %pixman_transform_t and an translation.
1109  *
1110  * Passing a non-zero value in x_offset/y_offset has the same effect
1111  * as applying cairo_matrix_translate(matrix, x_offset, y_offset) and
1112  * setting x_offset and y_offset to 0.
1113  *
1114  * If it is possible to represent the matrix with an identity
1115  * %pixman_transform_t and a translation within the valid range for
1116  * pixman, this function will set @out_transform to be the identity,
1117  * @x_offset and @y_offset to be the translation vector and will
1118  * return %CAIRO_INT_STATUS_NOTHING_TO_DO. Otherwise it will try to
1119  * evenly divide the translational component of @matrix between
1120  * @out_transform and (@x_offset, @y_offset).
1121  *
1122  * Upon return x_offset and y_offset contain the translation vector.
1123  *
1124  * Return value: %CAIRO_INT_STATUS_NOTHING_TO_DO if the out_transform
1125  * is the identity, %CAIRO_STATUS_INVALID_MATRIX if it was not
1126  * possible to represent @matrix as a pixman_transform_t without
1127  * overflows, %CAIRO_STATUS_SUCCESS otherwise.
1128  **/
1129 cairo_status_t
_cairo_matrix_to_pixman_matrix_offset(const cairo_matrix_t * matrix,cairo_filter_t filter,double xc,double yc,pixman_transform_t * out_transform,int * x_offset,int * y_offset)1130 _cairo_matrix_to_pixman_matrix_offset (const cairo_matrix_t	*matrix,
1131 				       cairo_filter_t            filter,
1132 				       double                    xc,
1133 				       double                    yc,
1134 				       pixman_transform_t	*out_transform,
1135 				       int                      *x_offset,
1136 				       int                      *y_offset)
1137 {
1138     cairo_bool_t is_pixman_translation;
1139 
1140     is_pixman_translation = _cairo_matrix_is_pixman_translation (matrix,
1141 								 filter,
1142 								 x_offset,
1143 								 y_offset);
1144 
1145     if (is_pixman_translation) {
1146 	*out_transform = pixman_identity_transform;
1147 	return CAIRO_INT_STATUS_NOTHING_TO_DO;
1148     } else {
1149 	cairo_matrix_t m;
1150 
1151 	m = *matrix;
1152 	cairo_matrix_translate (&m, *x_offset, *y_offset);
1153 	if (m.x0 != 0.0 || m.y0 != 0.0) {
1154 	    double tx, ty, norm;
1155 	    int i, j;
1156 
1157 	    /* pixman also limits the [xy]_offset to 16 bits so evenly
1158 	     * spread the bits between the two.
1159 	     *
1160 	     * To do this, find the solutions of:
1161 	     *   |x| = |x*m.xx + y*m.xy + m.x0|
1162 	     *   |y| = |x*m.yx + y*m.yy + m.y0|
1163 	     *
1164 	     * and select the one whose maximum norm is smallest.
1165 	     */
1166 	    tx = m.x0;
1167 	    ty = m.y0;
1168 	    norm = MAX (fabs (tx), fabs (ty));
1169 
1170 	    for (i = -1; i < 2; i+=2) {
1171 		for (j = -1; j < 2; j+=2) {
1172 		    double x, y, den, new_norm;
1173 
1174 		    den = (m.xx + i) * (m.yy + j) - m.xy * m.yx;
1175 		    if (fabs (den) < DBL_EPSILON)
1176 			continue;
1177 
1178 		    x = m.y0 * m.xy - m.x0 * (m.yy + j);
1179 		    y = m.x0 * m.yx - m.y0 * (m.xx + i);
1180 
1181 		    den = 1 / den;
1182 		    x *= den;
1183 		    y *= den;
1184 
1185 		    new_norm = MAX (fabs (x), fabs (y));
1186 		    if (norm > new_norm) {
1187 			norm = new_norm;
1188 			tx = x;
1189 			ty = y;
1190 		    }
1191 		}
1192 	    }
1193 
1194 	    tx = floor (tx);
1195 	    ty = floor (ty);
1196 	    *x_offset = -tx;
1197 	    *y_offset = -ty;
1198 	    cairo_matrix_translate (&m, tx, ty);
1199 	} else {
1200 	    *x_offset = 0;
1201 	    *y_offset = 0;
1202 	}
1203 
1204 	return _cairo_matrix_to_pixman_matrix (&m, out_transform, xc, yc);
1205     }
1206 }
1207