1/*
2 * Copyright © 2008 Keith Packard
3 *
4 * Permission to use, copy, modify, distribute, and sell this software and its
5 * documentation for any purpose is hereby granted without fee, provided that
6 * the above copyright notice appear in all copies and that both that copyright
7 * notice and this permission notice appear in supporting documentation, and
8 * that the name of the copyright holders not be used in advertising or
9 * publicity pertaining to distribution of the software without specific,
10 * written prior permission.  The copyright holders make no representations
11 * about the suitability of this software for any purpose.  It is provided "as
12 * is" without express or implied warranty.
13 *
14 * THE COPYRIGHT HOLDERS DISCLAIM ALL WARRANTIES WITH REGARD TO THIS SOFTWARE,
15 * INCLUDING ALL IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS, IN NO
16 * EVENT SHALL THE COPYRIGHT HOLDERS BE LIABLE FOR ANY SPECIAL, INDIRECT OR
17 * CONSEQUENTIAL DAMAGES OR ANY DAMAGES WHATSOEVER RESULTING FROM LOSS OF USE,
18 * DATA OR PROFITS, WHETHER IN AN ACTION OF CONTRACT, NEGLIGENCE OR OTHER
19 * TORTIOUS ACTION, ARISING OUT OF OR IN CONNECTION WITH THE USE OR PERFORMANCE
20 * OF THIS SOFTWARE.
21 */
22
23/*
24 * Matrix interfaces
25 */
26
27#ifdef HAVE_CONFIG_H
28#include <config.h>
29#endif
30
31#include <math.h>
32#include <string.h>
33#include "pixman-private.h"
34
35#define F(x)    pixman_int_to_fixed (x)
36
37static force_inline int
38count_leading_zeros (uint32_t x)
39{
40#ifdef __GNUC__
41    return __builtin_clz (x);
42#else
43    int n = 0;
44    while (x)
45    {
46        n++;
47        x >>= 1;
48    }
49    return 32 - n;
50#endif
51}
52
53/*
54 * Large signed/unsigned integer division with rounding for the platforms with
55 * only 64-bit integer data type supported (no 128-bit data type).
56 *
57 * Arguments:
58 *     hi, lo - high and low 64-bit parts of the dividend
59 *     div    - 48-bit divisor
60 *
61 * Returns: lowest 64 bits of the result as a return value and highest 64
62 *          bits of the result to "result_hi" pointer
63 */
64
65/* grade-school unsigned division (128-bit by 48-bit) with rounding to nearest */
66static force_inline uint64_t
67rounded_udiv_128_by_48 (uint64_t  hi,
68                        uint64_t  lo,
69                        uint64_t  div,
70                        uint64_t *result_hi)
71{
72    uint64_t tmp, remainder, result_lo;
73    assert(div < ((uint64_t)1 << 48));
74
75    remainder = hi % div;
76    *result_hi = hi / div;
77
78    tmp = (remainder << 16) + (lo >> 48);
79    result_lo = tmp / div;
80    remainder = tmp % div;
81
82    tmp = (remainder << 16) + ((lo >> 32) & 0xFFFF);
83    result_lo = (result_lo << 16) + (tmp / div);
84    remainder = tmp % div;
85
86    tmp = (remainder << 16) + ((lo >> 16) & 0xFFFF);
87    result_lo = (result_lo << 16) + (tmp / div);
88    remainder = tmp % div;
89
90    tmp = (remainder << 16) + (lo & 0xFFFF);
91    result_lo = (result_lo << 16) + (tmp / div);
92    remainder = tmp % div;
93
94    /* round to nearest */
95    if (remainder * 2 >= div && ++result_lo == 0)
96        *result_hi += 1;
97
98    return result_lo;
99}
100
101/* signed division (128-bit by 49-bit) with rounding to nearest */
102static inline int64_t
103rounded_sdiv_128_by_49 (int64_t   hi,
104                        uint64_t  lo,
105                        int64_t   div,
106                        int64_t  *signed_result_hi)
107{
108    uint64_t result_lo, result_hi;
109    int sign = 0;
110    if (div < 0)
111    {
112        div = -div;
113        sign ^= 1;
114    }
115    if (hi < 0)
116    {
117        if (lo != 0)
118            hi++;
119        hi = -hi;
120        lo = -lo;
121        sign ^= 1;
122    }
123    result_lo = rounded_udiv_128_by_48 (hi, lo, div, &result_hi);
124    if (sign)
125    {
126        if (result_lo != 0)
127            result_hi++;
128        result_hi = -result_hi;
129        result_lo = -result_lo;
130    }
131    if (signed_result_hi)
132    {
133        *signed_result_hi = result_hi;
134    }
135    return result_lo;
136}
137
138/*
139 * Multiply 64.16 fixed point value by (2^scalebits) and convert
140 * to 128-bit integer.
141 */
142static force_inline void
143fixed_64_16_to_int128 (int64_t  hi,
144                       int64_t  lo,
145                       int64_t *rhi,
146                       int64_t *rlo,
147                       int      scalebits)
148{
149    /* separate integer and fractional parts */
150    hi += lo >> 16;
151    lo &= 0xFFFF;
152
153    if (scalebits <= 0)
154    {
155        *rlo = hi >> (-scalebits);
156        *rhi = *rlo >> 63;
157    }
158    else
159    {
160        *rhi = hi >> (64 - scalebits);
161        *rlo = (uint64_t)hi << scalebits;
162        if (scalebits < 16)
163            *rlo += lo >> (16 - scalebits);
164        else
165            *rlo += lo << (scalebits - 16);
166    }
167}
168
169/*
170 * Convert 112.16 fixed point value to 48.16 with clamping for the out
171 * of range values.
172 */
173static force_inline pixman_fixed_48_16_t
174fixed_112_16_to_fixed_48_16 (int64_t hi, int64_t lo, pixman_bool_t *clampflag)
175{
176    if ((lo >> 63) != hi)
177    {
178        *clampflag = TRUE;
179        return hi >= 0 ? INT64_MAX : INT64_MIN;
180    }
181    else
182    {
183        return lo;
184    }
185}
186
187/*
188 * Transform a point with 31.16 fixed point coordinates from the destination
189 * space to a point with 48.16 fixed point coordinates in the source space.
190 * No overflows are possible for affine transformations and the results are
191 * accurate including the least significant bit. Projective transformations
192 * may overflow, in this case the results are just clamped to return maximum
193 * or minimum 48.16 values (so that the caller can at least handle the NONE
194 * and PAD repeats correctly) and the return value is FALSE to indicate that
195 * such clamping has happened.
196 */
197PIXMAN_EXPORT pixman_bool_t
198pixman_transform_point_31_16 (const pixman_transform_t    *t,
199                              const pixman_vector_48_16_t *v,
200                              pixman_vector_48_16_t       *result)
201{
202    pixman_bool_t clampflag = FALSE;
203    int i;
204    int64_t tmp[3][2], divint;
205    uint16_t divfrac;
206
207    /* input vector values must have no more than 31 bits (including sign)
208     * in the integer part */
209    assert (v->v[0] <   ((pixman_fixed_48_16_t)1 << (30 + 16)));
210    assert (v->v[0] >= -((pixman_fixed_48_16_t)1 << (30 + 16)));
211    assert (v->v[1] <   ((pixman_fixed_48_16_t)1 << (30 + 16)));
212    assert (v->v[1] >= -((pixman_fixed_48_16_t)1 << (30 + 16)));
213    assert (v->v[2] <   ((pixman_fixed_48_16_t)1 << (30 + 16)));
214    assert (v->v[2] >= -((pixman_fixed_48_16_t)1 << (30 + 16)));
215
216    for (i = 0; i < 3; i++)
217    {
218        tmp[i][0] = (int64_t)t->matrix[i][0] * (v->v[0] >> 16);
219        tmp[i][1] = (int64_t)t->matrix[i][0] * (v->v[0] & 0xFFFF);
220        tmp[i][0] += (int64_t)t->matrix[i][1] * (v->v[1] >> 16);
221        tmp[i][1] += (int64_t)t->matrix[i][1] * (v->v[1] & 0xFFFF);
222        tmp[i][0] += (int64_t)t->matrix[i][2] * (v->v[2] >> 16);
223        tmp[i][1] += (int64_t)t->matrix[i][2] * (v->v[2] & 0xFFFF);
224    }
225
226    /*
227     * separate 64-bit integer and 16-bit fractional parts for the divisor,
228     * which is also scaled by 65536 after fixed point multiplication.
229     */
230    divint  = tmp[2][0] + (tmp[2][1] >> 16);
231    divfrac = tmp[2][1] & 0xFFFF;
232
233    if (divint == pixman_fixed_1 && divfrac == 0)
234    {
235        /*
236         * this is a simple affine transformation
237         */
238        result->v[0] = tmp[0][0] + ((tmp[0][1] + 0x8000) >> 16);
239        result->v[1] = tmp[1][0] + ((tmp[1][1] + 0x8000) >> 16);
240        result->v[2] = pixman_fixed_1;
241    }
242    else if (divint == 0 && divfrac == 0)
243    {
244        /*
245         * handle zero divisor (if the values are non-zero, set the
246         * results to maximum positive or minimum negative)
247         */
248        clampflag = TRUE;
249
250        result->v[0] = tmp[0][0] + ((tmp[0][1] + 0x8000) >> 16);
251        result->v[1] = tmp[1][0] + ((tmp[1][1] + 0x8000) >> 16);
252
253        if (result->v[0] > 0)
254            result->v[0] = INT64_MAX;
255        else if (result->v[0] < 0)
256            result->v[0] = INT64_MIN;
257
258        if (result->v[1] > 0)
259            result->v[1] = INT64_MAX;
260        else if (result->v[1] < 0)
261            result->v[1] = INT64_MIN;
262    }
263    else
264    {
265        /*
266         * projective transformation, analyze the top 32 bits of the divisor
267         */
268        int32_t hi32divbits = divint >> 32;
269        if (hi32divbits < 0)
270            hi32divbits = ~hi32divbits;
271
272        if (hi32divbits == 0)
273        {
274            /* the divisor is small, we can actually keep all the bits */
275            int64_t hi, rhi, lo, rlo;
276            int64_t div = (divint << 16) + divfrac;
277
278            fixed_64_16_to_int128 (tmp[0][0], tmp[0][1], &hi, &lo, 32);
279            rlo = rounded_sdiv_128_by_49 (hi, lo, div, &rhi);
280            result->v[0] = fixed_112_16_to_fixed_48_16 (rhi, rlo, &clampflag);
281
282            fixed_64_16_to_int128 (tmp[1][0], tmp[1][1], &hi, &lo, 32);
283            rlo = rounded_sdiv_128_by_49 (hi, lo, div, &rhi);
284            result->v[1] = fixed_112_16_to_fixed_48_16 (rhi, rlo, &clampflag);
285        }
286        else
287        {
288            /* the divisor needs to be reduced to 48 bits */
289            int64_t hi, rhi, lo, rlo, div;
290            int shift = 32 - count_leading_zeros (hi32divbits);
291            fixed_64_16_to_int128 (divint, divfrac, &hi, &div, 16 - shift);
292
293            fixed_64_16_to_int128 (tmp[0][0], tmp[0][1], &hi, &lo, 32 - shift);
294            rlo = rounded_sdiv_128_by_49 (hi, lo, div, &rhi);
295            result->v[0] = fixed_112_16_to_fixed_48_16 (rhi, rlo, &clampflag);
296
297            fixed_64_16_to_int128 (tmp[1][0], tmp[1][1], &hi, &lo, 32 - shift);
298            rlo = rounded_sdiv_128_by_49 (hi, lo, div, &rhi);
299            result->v[1] = fixed_112_16_to_fixed_48_16 (rhi, rlo, &clampflag);
300        }
301    }
302    result->v[2] = pixman_fixed_1;
303    return !clampflag;
304}
305
306PIXMAN_EXPORT void
307pixman_transform_point_31_16_affine (const pixman_transform_t    *t,
308                                     const pixman_vector_48_16_t *v,
309                                     pixman_vector_48_16_t       *result)
310{
311    int64_t hi0, lo0, hi1, lo1;
312
313    /* input vector values must have no more than 31 bits (including sign)
314     * in the integer part */
315    assert (v->v[0] <   ((pixman_fixed_48_16_t)1 << (30 + 16)));
316    assert (v->v[0] >= -((pixman_fixed_48_16_t)1 << (30 + 16)));
317    assert (v->v[1] <   ((pixman_fixed_48_16_t)1 << (30 + 16)));
318    assert (v->v[1] >= -((pixman_fixed_48_16_t)1 << (30 + 16)));
319
320    hi0  = (int64_t)t->matrix[0][0] * (v->v[0] >> 16);
321    lo0  = (int64_t)t->matrix[0][0] * (v->v[0] & 0xFFFF);
322    hi0 += (int64_t)t->matrix[0][1] * (v->v[1] >> 16);
323    lo0 += (int64_t)t->matrix[0][1] * (v->v[1] & 0xFFFF);
324    hi0 += (int64_t)t->matrix[0][2];
325
326    hi1  = (int64_t)t->matrix[1][0] * (v->v[0] >> 16);
327    lo1  = (int64_t)t->matrix[1][0] * (v->v[0] & 0xFFFF);
328    hi1 += (int64_t)t->matrix[1][1] * (v->v[1] >> 16);
329    lo1 += (int64_t)t->matrix[1][1] * (v->v[1] & 0xFFFF);
330    hi1 += (int64_t)t->matrix[1][2];
331
332    result->v[0] = hi0 + ((lo0 + 0x8000) >> 16);
333    result->v[1] = hi1 + ((lo1 + 0x8000) >> 16);
334    result->v[2] = pixman_fixed_1;
335}
336
337PIXMAN_EXPORT void
338pixman_transform_point_31_16_3d (const pixman_transform_t    *t,
339                                 const pixman_vector_48_16_t *v,
340                                 pixman_vector_48_16_t       *result)
341{
342    int i;
343    int64_t tmp[3][2];
344
345    /* input vector values must have no more than 31 bits (including sign)
346     * in the integer part */
347    assert (v->v[0] <   ((pixman_fixed_48_16_t)1 << (30 + 16)));
348    assert (v->v[0] >= -((pixman_fixed_48_16_t)1 << (30 + 16)));
349    assert (v->v[1] <   ((pixman_fixed_48_16_t)1 << (30 + 16)));
350    assert (v->v[1] >= -((pixman_fixed_48_16_t)1 << (30 + 16)));
351    assert (v->v[2] <   ((pixman_fixed_48_16_t)1 << (30 + 16)));
352    assert (v->v[2] >= -((pixman_fixed_48_16_t)1 << (30 + 16)));
353
354    for (i = 0; i < 3; i++)
355    {
356        tmp[i][0] = (int64_t)t->matrix[i][0] * (v->v[0] >> 16);
357        tmp[i][1] = (int64_t)t->matrix[i][0] * (v->v[0] & 0xFFFF);
358        tmp[i][0] += (int64_t)t->matrix[i][1] * (v->v[1] >> 16);
359        tmp[i][1] += (int64_t)t->matrix[i][1] * (v->v[1] & 0xFFFF);
360        tmp[i][0] += (int64_t)t->matrix[i][2] * (v->v[2] >> 16);
361        tmp[i][1] += (int64_t)t->matrix[i][2] * (v->v[2] & 0xFFFF);
362    }
363
364    result->v[0] = tmp[0][0] + ((tmp[0][1] + 0x8000) >> 16);
365    result->v[1] = tmp[1][0] + ((tmp[1][1] + 0x8000) >> 16);
366    result->v[2] = tmp[2][0] + ((tmp[2][1] + 0x8000) >> 16);
367}
368
369PIXMAN_EXPORT void
370pixman_transform_init_identity (struct pixman_transform *matrix)
371{
372    int i;
373
374    memset (matrix, '\0', sizeof (struct pixman_transform));
375    for (i = 0; i < 3; i++)
376	matrix->matrix[i][i] = F (1);
377}
378
379typedef pixman_fixed_32_32_t pixman_fixed_34_30_t;
380
381PIXMAN_EXPORT pixman_bool_t
382pixman_transform_point_3d (const struct pixman_transform *transform,
383                           struct pixman_vector *         vector)
384{
385    pixman_vector_48_16_t tmp;
386    tmp.v[0] = vector->vector[0];
387    tmp.v[1] = vector->vector[1];
388    tmp.v[2] = vector->vector[2];
389
390    pixman_transform_point_31_16_3d (transform, &tmp, &tmp);
391
392    vector->vector[0] = tmp.v[0];
393    vector->vector[1] = tmp.v[1];
394    vector->vector[2] = tmp.v[2];
395
396    return vector->vector[0] == tmp.v[0] &&
397           vector->vector[1] == tmp.v[1] &&
398           vector->vector[2] == tmp.v[2];
399}
400
401PIXMAN_EXPORT pixman_bool_t
402pixman_transform_point (const struct pixman_transform *transform,
403                        struct pixman_vector *         vector)
404{
405    pixman_vector_48_16_t tmp;
406    tmp.v[0] = vector->vector[0];
407    tmp.v[1] = vector->vector[1];
408    tmp.v[2] = vector->vector[2];
409
410    if (!pixman_transform_point_31_16 (transform, &tmp, &tmp))
411        return FALSE;
412
413    vector->vector[0] = tmp.v[0];
414    vector->vector[1] = tmp.v[1];
415    vector->vector[2] = tmp.v[2];
416
417    return vector->vector[0] == tmp.v[0] &&
418           vector->vector[1] == tmp.v[1] &&
419           vector->vector[2] == tmp.v[2];
420}
421
422PIXMAN_EXPORT pixman_bool_t
423pixman_transform_multiply (struct pixman_transform *      dst,
424                           const struct pixman_transform *l,
425                           const struct pixman_transform *r)
426{
427    struct pixman_transform d;
428    int dx, dy;
429    int o;
430
431    for (dy = 0; dy < 3; dy++)
432    {
433	for (dx = 0; dx < 3; dx++)
434	{
435	    pixman_fixed_48_16_t v;
436	    pixman_fixed_32_32_t partial;
437
438	    v = 0;
439	    for (o = 0; o < 3; o++)
440	    {
441		partial =
442		    (pixman_fixed_32_32_t) l->matrix[dy][o] *
443		    (pixman_fixed_32_32_t) r->matrix[o][dx];
444
445		v += (partial + 0x8000) >> 16;
446	    }
447
448	    if (v > pixman_max_fixed_48_16 || v < pixman_min_fixed_48_16)
449		return FALSE;
450
451	    d.matrix[dy][dx] = (pixman_fixed_t) v;
452	}
453    }
454
455    *dst = d;
456    return TRUE;
457}
458
459PIXMAN_EXPORT void
460pixman_transform_init_scale (struct pixman_transform *t,
461                             pixman_fixed_t           sx,
462                             pixman_fixed_t           sy)
463{
464    memset (t, '\0', sizeof (struct pixman_transform));
465
466    t->matrix[0][0] = sx;
467    t->matrix[1][1] = sy;
468    t->matrix[2][2] = F (1);
469}
470
471static pixman_fixed_t
472fixed_inverse (pixman_fixed_t x)
473{
474    return (pixman_fixed_t) ((((pixman_fixed_48_16_t) F (1)) * F (1)) / x);
475}
476
477PIXMAN_EXPORT pixman_bool_t
478pixman_transform_scale (struct pixman_transform *forward,
479                        struct pixman_transform *reverse,
480                        pixman_fixed_t           sx,
481                        pixman_fixed_t           sy)
482{
483    struct pixman_transform t;
484
485    if (sx == 0 || sy == 0)
486	return FALSE;
487
488    if (forward)
489    {
490	pixman_transform_init_scale (&t, sx, sy);
491	if (!pixman_transform_multiply (forward, &t, forward))
492	    return FALSE;
493    }
494
495    if (reverse)
496    {
497	pixman_transform_init_scale (&t, fixed_inverse (sx),
498	                             fixed_inverse (sy));
499	if (!pixman_transform_multiply (reverse, reverse, &t))
500	    return FALSE;
501    }
502
503    return TRUE;
504}
505
506PIXMAN_EXPORT void
507pixman_transform_init_rotate (struct pixman_transform *t,
508                              pixman_fixed_t           c,
509                              pixman_fixed_t           s)
510{
511    memset (t, '\0', sizeof (struct pixman_transform));
512
513    t->matrix[0][0] = c;
514    t->matrix[0][1] = -s;
515    t->matrix[1][0] = s;
516    t->matrix[1][1] = c;
517    t->matrix[2][2] = F (1);
518}
519
520PIXMAN_EXPORT pixman_bool_t
521pixman_transform_rotate (struct pixman_transform *forward,
522                         struct pixman_transform *reverse,
523                         pixman_fixed_t           c,
524                         pixman_fixed_t           s)
525{
526    struct pixman_transform t;
527
528    if (forward)
529    {
530	pixman_transform_init_rotate (&t, c, s);
531	if (!pixman_transform_multiply (forward, &t, forward))
532	    return FALSE;
533    }
534
535    if (reverse)
536    {
537	pixman_transform_init_rotate (&t, c, -s);
538	if (!pixman_transform_multiply (reverse, reverse, &t))
539	    return FALSE;
540    }
541
542    return TRUE;
543}
544
545PIXMAN_EXPORT void
546pixman_transform_init_translate (struct pixman_transform *t,
547                                 pixman_fixed_t           tx,
548                                 pixman_fixed_t           ty)
549{
550    memset (t, '\0', sizeof (struct pixman_transform));
551
552    t->matrix[0][0] = F (1);
553    t->matrix[0][2] = tx;
554    t->matrix[1][1] = F (1);
555    t->matrix[1][2] = ty;
556    t->matrix[2][2] = F (1);
557}
558
559PIXMAN_EXPORT pixman_bool_t
560pixman_transform_translate (struct pixman_transform *forward,
561                            struct pixman_transform *reverse,
562                            pixman_fixed_t           tx,
563                            pixman_fixed_t           ty)
564{
565    struct pixman_transform t;
566
567    if (forward)
568    {
569	pixman_transform_init_translate (&t, tx, ty);
570
571	if (!pixman_transform_multiply (forward, &t, forward))
572	    return FALSE;
573    }
574
575    if (reverse)
576    {
577	pixman_transform_init_translate (&t, -tx, -ty);
578
579	if (!pixman_transform_multiply (reverse, reverse, &t))
580	    return FALSE;
581    }
582    return TRUE;
583}
584
585PIXMAN_EXPORT pixman_bool_t
586pixman_transform_bounds (const struct pixman_transform *matrix,
587                         struct pixman_box16 *          b)
588
589{
590    struct pixman_vector v[4];
591    int i;
592    int x1, y1, x2, y2;
593
594    v[0].vector[0] = F (b->x1);
595    v[0].vector[1] = F (b->y1);
596    v[0].vector[2] = F (1);
597
598    v[1].vector[0] = F (b->x2);
599    v[1].vector[1] = F (b->y1);
600    v[1].vector[2] = F (1);
601
602    v[2].vector[0] = F (b->x2);
603    v[2].vector[1] = F (b->y2);
604    v[2].vector[2] = F (1);
605
606    v[3].vector[0] = F (b->x1);
607    v[3].vector[1] = F (b->y2);
608    v[3].vector[2] = F (1);
609
610    for (i = 0; i < 4; i++)
611    {
612	if (!pixman_transform_point (matrix, &v[i]))
613	    return FALSE;
614
615	x1 = pixman_fixed_to_int (v[i].vector[0]);
616	y1 = pixman_fixed_to_int (v[i].vector[1]);
617	x2 = pixman_fixed_to_int (pixman_fixed_ceil (v[i].vector[0]));
618	y2 = pixman_fixed_to_int (pixman_fixed_ceil (v[i].vector[1]));
619
620	if (i == 0)
621	{
622	    b->x1 = x1;
623	    b->y1 = y1;
624	    b->x2 = x2;
625	    b->y2 = y2;
626	}
627	else
628	{
629	    if (x1 < b->x1) b->x1 = x1;
630	    if (y1 < b->y1) b->y1 = y1;
631	    if (x2 > b->x2) b->x2 = x2;
632	    if (y2 > b->y2) b->y2 = y2;
633	}
634    }
635
636    return TRUE;
637}
638
639PIXMAN_EXPORT pixman_bool_t
640pixman_transform_invert (struct pixman_transform *      dst,
641                         const struct pixman_transform *src)
642{
643    struct pixman_f_transform m;
644
645    pixman_f_transform_from_pixman_transform (&m, src);
646
647    if (!pixman_f_transform_invert (&m, &m))
648	return FALSE;
649
650    if (!pixman_transform_from_pixman_f_transform (dst, &m))
651	return FALSE;
652
653    return TRUE;
654}
655
656static pixman_bool_t
657within_epsilon (pixman_fixed_t a,
658                pixman_fixed_t b,
659                pixman_fixed_t epsilon)
660{
661    pixman_fixed_t t = a - b;
662
663    if (t < 0)
664	t = -t;
665
666    return t <= epsilon;
667}
668
669#define EPSILON (pixman_fixed_t) (2)
670
671#define IS_SAME(a, b) (within_epsilon (a, b, EPSILON))
672#define IS_ZERO(a)    (within_epsilon (a, 0, EPSILON))
673#define IS_ONE(a)     (within_epsilon (a, F (1), EPSILON))
674#define IS_UNIT(a)			    \
675    (within_epsilon (a, F (1), EPSILON) ||  \
676     within_epsilon (a, F (-1), EPSILON) || \
677     IS_ZERO (a))
678#define IS_INT(a)    (IS_ZERO (pixman_fixed_frac (a)))
679
680PIXMAN_EXPORT pixman_bool_t
681pixman_transform_is_identity (const struct pixman_transform *t)
682{
683    return (IS_SAME (t->matrix[0][0], t->matrix[1][1]) &&
684	    IS_SAME (t->matrix[0][0], t->matrix[2][2]) &&
685	    !IS_ZERO (t->matrix[0][0]) &&
686	    IS_ZERO (t->matrix[0][1]) &&
687	    IS_ZERO (t->matrix[0][2]) &&
688	    IS_ZERO (t->matrix[1][0]) &&
689	    IS_ZERO (t->matrix[1][2]) &&
690	    IS_ZERO (t->matrix[2][0]) &&
691	    IS_ZERO (t->matrix[2][1]));
692}
693
694PIXMAN_EXPORT pixman_bool_t
695pixman_transform_is_scale (const struct pixman_transform *t)
696{
697    return (!IS_ZERO (t->matrix[0][0]) &&
698            IS_ZERO (t->matrix[0][1]) &&
699            IS_ZERO (t->matrix[0][2]) &&
700
701            IS_ZERO (t->matrix[1][0]) &&
702            !IS_ZERO (t->matrix[1][1]) &&
703            IS_ZERO (t->matrix[1][2]) &&
704
705            IS_ZERO (t->matrix[2][0]) &&
706            IS_ZERO (t->matrix[2][1]) &&
707            !IS_ZERO (t->matrix[2][2]));
708}
709
710PIXMAN_EXPORT pixman_bool_t
711pixman_transform_is_int_translate (const struct pixman_transform *t)
712{
713    return (IS_ONE (t->matrix[0][0]) &&
714            IS_ZERO (t->matrix[0][1]) &&
715            IS_INT (t->matrix[0][2]) &&
716
717            IS_ZERO (t->matrix[1][0]) &&
718            IS_ONE (t->matrix[1][1]) &&
719            IS_INT (t->matrix[1][2]) &&
720
721            IS_ZERO (t->matrix[2][0]) &&
722            IS_ZERO (t->matrix[2][1]) &&
723            IS_ONE (t->matrix[2][2]));
724}
725
726PIXMAN_EXPORT pixman_bool_t
727pixman_transform_is_inverse (const struct pixman_transform *a,
728                             const struct pixman_transform *b)
729{
730    struct pixman_transform t;
731
732    if (!pixman_transform_multiply (&t, a, b))
733	return FALSE;
734
735    return pixman_transform_is_identity (&t);
736}
737
738PIXMAN_EXPORT void
739pixman_f_transform_from_pixman_transform (struct pixman_f_transform *    ft,
740                                          const struct pixman_transform *t)
741{
742    int i, j;
743
744    for (j = 0; j < 3; j++)
745    {
746	for (i = 0; i < 3; i++)
747	    ft->m[j][i] = pixman_fixed_to_double (t->matrix[j][i]);
748    }
749}
750
751PIXMAN_EXPORT pixman_bool_t
752pixman_transform_from_pixman_f_transform (struct pixman_transform *        t,
753                                          const struct pixman_f_transform *ft)
754{
755    int i, j;
756
757    for (j = 0; j < 3; j++)
758    {
759	for (i = 0; i < 3; i++)
760	{
761	    double d = ft->m[j][i];
762	    if (d < -32767.0 || d > 32767.0)
763		return FALSE;
764	    d = d * 65536.0 + 0.5;
765	    t->matrix[j][i] = (pixman_fixed_t) floor (d);
766	}
767    }
768
769    return TRUE;
770}
771
772PIXMAN_EXPORT pixman_bool_t
773pixman_f_transform_invert (struct pixman_f_transform *      dst,
774                           const struct pixman_f_transform *src)
775{
776    static const int a[3] = { 2, 2, 1 };
777    static const int b[3] = { 1, 0, 0 };
778    pixman_f_transform_t d;
779    double det;
780    int i, j;
781
782    det = 0;
783    for (i = 0; i < 3; i++)
784    {
785	double p;
786	int ai = a[i];
787	int bi = b[i];
788	p = src->m[i][0] * (src->m[ai][2] * src->m[bi][1] -
789	                    src->m[ai][1] * src->m[bi][2]);
790	if (i == 1)
791	    p = -p;
792	det += p;
793    }
794
795    if (det == 0)
796	return FALSE;
797
798    det = 1 / det;
799    for (j = 0; j < 3; j++)
800    {
801	for (i = 0; i < 3; i++)
802	{
803	    double p;
804	    int ai = a[i];
805	    int aj = a[j];
806	    int bi = b[i];
807	    int bj = b[j];
808
809	    p = (src->m[ai][aj] * src->m[bi][bj] -
810	         src->m[ai][bj] * src->m[bi][aj]);
811
812	    if (((i + j) & 1) != 0)
813		p = -p;
814
815	    d.m[j][i] = det * p;
816	}
817    }
818
819    *dst = d;
820
821    return TRUE;
822}
823
824PIXMAN_EXPORT pixman_bool_t
825pixman_f_transform_point (const struct pixman_f_transform *t,
826                          struct pixman_f_vector *         v)
827{
828    struct pixman_f_vector result;
829    int i, j;
830    double a;
831
832    for (j = 0; j < 3; j++)
833    {
834	a = 0;
835	for (i = 0; i < 3; i++)
836	    a += t->m[j][i] * v->v[i];
837	result.v[j] = a;
838    }
839
840    if (!result.v[2])
841	return FALSE;
842
843    for (j = 0; j < 2; j++)
844	v->v[j] = result.v[j] / result.v[2];
845
846    v->v[2] = 1;
847
848    return TRUE;
849}
850
851PIXMAN_EXPORT void
852pixman_f_transform_point_3d (const struct pixman_f_transform *t,
853                             struct pixman_f_vector *         v)
854{
855    struct pixman_f_vector result;
856    int i, j;
857    double a;
858
859    for (j = 0; j < 3; j++)
860    {
861	a = 0;
862	for (i = 0; i < 3; i++)
863	    a += t->m[j][i] * v->v[i];
864	result.v[j] = a;
865    }
866
867    *v = result;
868}
869
870PIXMAN_EXPORT void
871pixman_f_transform_multiply (struct pixman_f_transform *      dst,
872                             const struct pixman_f_transform *l,
873                             const struct pixman_f_transform *r)
874{
875    struct pixman_f_transform d;
876    int dx, dy;
877    int o;
878
879    for (dy = 0; dy < 3; dy++)
880    {
881	for (dx = 0; dx < 3; dx++)
882	{
883	    double v = 0;
884	    for (o = 0; o < 3; o++)
885		v += l->m[dy][o] * r->m[o][dx];
886	    d.m[dy][dx] = v;
887	}
888    }
889
890    *dst = d;
891}
892
893PIXMAN_EXPORT void
894pixman_f_transform_init_scale (struct pixman_f_transform *t,
895                               double                     sx,
896                               double                     sy)
897{
898    t->m[0][0] = sx;
899    t->m[0][1] = 0;
900    t->m[0][2] = 0;
901    t->m[1][0] = 0;
902    t->m[1][1] = sy;
903    t->m[1][2] = 0;
904    t->m[2][0] = 0;
905    t->m[2][1] = 0;
906    t->m[2][2] = 1;
907}
908
909PIXMAN_EXPORT pixman_bool_t
910pixman_f_transform_scale (struct pixman_f_transform *forward,
911                          struct pixman_f_transform *reverse,
912                          double                     sx,
913                          double                     sy)
914{
915    struct pixman_f_transform t;
916
917    if (sx == 0 || sy == 0)
918	return FALSE;
919
920    if (forward)
921    {
922	pixman_f_transform_init_scale (&t, sx, sy);
923	pixman_f_transform_multiply (forward, &t, forward);
924    }
925
926    if (reverse)
927    {
928	pixman_f_transform_init_scale (&t, 1 / sx, 1 / sy);
929	pixman_f_transform_multiply (reverse, reverse, &t);
930    }
931
932    return TRUE;
933}
934
935PIXMAN_EXPORT void
936pixman_f_transform_init_rotate (struct pixman_f_transform *t,
937                                double                     c,
938                                double                     s)
939{
940    t->m[0][0] = c;
941    t->m[0][1] = -s;
942    t->m[0][2] = 0;
943    t->m[1][0] = s;
944    t->m[1][1] = c;
945    t->m[1][2] = 0;
946    t->m[2][0] = 0;
947    t->m[2][1] = 0;
948    t->m[2][2] = 1;
949}
950
951PIXMAN_EXPORT pixman_bool_t
952pixman_f_transform_rotate (struct pixman_f_transform *forward,
953                           struct pixman_f_transform *reverse,
954                           double                     c,
955                           double                     s)
956{
957    struct pixman_f_transform t;
958
959    if (forward)
960    {
961	pixman_f_transform_init_rotate (&t, c, s);
962	pixman_f_transform_multiply (forward, &t, forward);
963    }
964
965    if (reverse)
966    {
967	pixman_f_transform_init_rotate (&t, c, -s);
968	pixman_f_transform_multiply (reverse, reverse, &t);
969    }
970
971    return TRUE;
972}
973
974PIXMAN_EXPORT void
975pixman_f_transform_init_translate (struct pixman_f_transform *t,
976                                   double                     tx,
977                                   double                     ty)
978{
979    t->m[0][0] = 1;
980    t->m[0][1] = 0;
981    t->m[0][2] = tx;
982    t->m[1][0] = 0;
983    t->m[1][1] = 1;
984    t->m[1][2] = ty;
985    t->m[2][0] = 0;
986    t->m[2][1] = 0;
987    t->m[2][2] = 1;
988}
989
990PIXMAN_EXPORT pixman_bool_t
991pixman_f_transform_translate (struct pixman_f_transform *forward,
992                              struct pixman_f_transform *reverse,
993                              double                     tx,
994                              double                     ty)
995{
996    struct pixman_f_transform t;
997
998    if (forward)
999    {
1000	pixman_f_transform_init_translate (&t, tx, ty);
1001	pixman_f_transform_multiply (forward, &t, forward);
1002    }
1003
1004    if (reverse)
1005    {
1006	pixman_f_transform_init_translate (&t, -tx, -ty);
1007	pixman_f_transform_multiply (reverse, reverse, &t);
1008    }
1009
1010    return TRUE;
1011}
1012
1013PIXMAN_EXPORT pixman_bool_t
1014pixman_f_transform_bounds (const struct pixman_f_transform *t,
1015                           struct pixman_box16 *            b)
1016{
1017    struct pixman_f_vector v[4];
1018    int i;
1019    int x1, y1, x2, y2;
1020
1021    v[0].v[0] = b->x1;
1022    v[0].v[1] = b->y1;
1023    v[0].v[2] = 1;
1024    v[1].v[0] = b->x2;
1025    v[1].v[1] = b->y1;
1026    v[1].v[2] = 1;
1027    v[2].v[0] = b->x2;
1028    v[2].v[1] = b->y2;
1029    v[2].v[2] = 1;
1030    v[3].v[0] = b->x1;
1031    v[3].v[1] = b->y2;
1032    v[3].v[2] = 1;
1033
1034    for (i = 0; i < 4; i++)
1035    {
1036	if (!pixman_f_transform_point (t, &v[i]))
1037	    return FALSE;
1038
1039	x1 = floor (v[i].v[0]);
1040	y1 = floor (v[i].v[1]);
1041	x2 = ceil (v[i].v[0]);
1042	y2 = ceil (v[i].v[1]);
1043
1044	if (i == 0)
1045	{
1046	    b->x1 = x1;
1047	    b->y1 = y1;
1048	    b->x2 = x2;
1049	    b->y2 = y2;
1050	}
1051	else
1052	{
1053	    if (x1 < b->x1) b->x1 = x1;
1054	    if (y1 < b->y1) b->y1 = y1;
1055	    if (x2 > b->x2) b->x2 = x2;
1056	    if (y2 > b->y2) b->y2 = y2;
1057	}
1058    }
1059
1060    return TRUE;
1061}
1062
1063PIXMAN_EXPORT void
1064pixman_f_transform_init_identity (struct pixman_f_transform *t)
1065{
1066    int i, j;
1067
1068    for (j = 0; j < 3; j++)
1069    {
1070	for (i = 0; i < 3; i++)
1071	    t->m[j][i] = i == j ? 1 : 0;
1072    }
1073}
1074