/src/ghostpdl/base/gsmatrix.c
Line | Count | Source (jump to first uncovered line) |
1 | | /* Copyright (C) 2001-2023 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., 39 Mesa Street, Suite 108A, San Francisco, |
13 | | CA 94129, USA, for further information. |
14 | | */ |
15 | | |
16 | | |
17 | | /* Matrix operators for Ghostscript library */ |
18 | | #include "math_.h" |
19 | | #include "memory_.h" |
20 | | #include "gx.h" |
21 | | #include "gserrors.h" |
22 | | #include "gxfarith.h" |
23 | | #include "gxfixed.h" |
24 | | #include "gxmatrix.h" |
25 | | #include "stream.h" |
26 | | |
27 | | /* The identity matrix */ |
28 | | static const gs_matrix gs_identity_matrix = |
29 | | {identity_matrix_body}; |
30 | | |
31 | | /* ------ Matrix creation ------ */ |
32 | | |
33 | | /* Create an identity matrix */ |
34 | | void |
35 | | gs_make_identity(gs_matrix * pmat) |
36 | 100M | { |
37 | 100M | *pmat = gs_identity_matrix; |
38 | 100M | } |
39 | | |
40 | | /* Create a translation matrix */ |
41 | | int |
42 | | gs_make_translation(double dx, double dy, gs_matrix * pmat) |
43 | 136k | { |
44 | 136k | *pmat = gs_identity_matrix; |
45 | 136k | pmat->tx = dx; |
46 | 136k | pmat->ty = dy; |
47 | 136k | return 0; |
48 | 136k | } |
49 | | |
50 | | /* Create a scaling matrix */ |
51 | | int |
52 | | gs_make_scaling(double sx, double sy, gs_matrix * pmat) |
53 | 8.56M | { |
54 | 8.56M | *pmat = gs_identity_matrix; |
55 | 8.56M | pmat->xx = sx; |
56 | 8.56M | pmat->yy = sy; |
57 | 8.56M | return 0; |
58 | 8.56M | } |
59 | | |
60 | | /* Create a rotation matrix. */ |
61 | | /* The angle is in degrees. */ |
62 | | int |
63 | | gs_make_rotation(double ang, gs_matrix * pmat) |
64 | 1.63M | { |
65 | 1.63M | gs_sincos_t sincos; |
66 | | |
67 | 1.63M | gs_sincos_degrees(ang, &sincos); |
68 | 1.63M | pmat->yy = pmat->xx = sincos.cos; |
69 | 1.63M | pmat->xy = sincos.sin; |
70 | 1.63M | pmat->yx = -sincos.sin; |
71 | 1.63M | pmat->tx = pmat->ty = 0.0; |
72 | 1.63M | return 0; |
73 | 1.63M | } |
74 | | |
75 | | /* ------ Matrix arithmetic ------ */ |
76 | | |
77 | | /* Multiply two matrices. We should check for floating exceptions, */ |
78 | | /* but for the moment it's just too awkward. */ |
79 | | /* Since this is used heavily, we check for shortcuts. */ |
80 | | int |
81 | | gs_matrix_multiply(const gs_matrix * pm1, const gs_matrix * pm2, gs_matrix * pmr) |
82 | 58.0M | { |
83 | 58.0M | double xx1 = pm1->xx, yy1 = pm1->yy; |
84 | 58.0M | double tx1 = pm1->tx, ty1 = pm1->ty; |
85 | 58.0M | double xx2 = pm2->xx, yy2 = pm2->yy; |
86 | 58.0M | double xy2 = pm2->xy, yx2 = pm2->yx; |
87 | | |
88 | 58.0M | if (is_xxyy(pm1)) { |
89 | 56.2M | pmr->tx = tx1 * xx2 + pm2->tx; |
90 | 56.2M | pmr->ty = ty1 * yy2 + pm2->ty; |
91 | 56.2M | if (is_fzero(xy2)) |
92 | 55.1M | pmr->xy = 0; |
93 | 1.12M | else |
94 | 1.12M | pmr->xy = xx1 * xy2, |
95 | 1.12M | pmr->ty += tx1 * xy2; |
96 | 56.2M | pmr->xx = xx1 * xx2; |
97 | 56.2M | if (is_fzero(yx2)) |
98 | 55.1M | pmr->yx = 0; |
99 | 1.09M | else |
100 | 1.09M | pmr->yx = yy1 * yx2, |
101 | 1.09M | pmr->tx += ty1 * yx2; |
102 | 56.2M | pmr->yy = yy1 * yy2; |
103 | 56.2M | } else { |
104 | 1.76M | double xy1 = pm1->xy, yx1 = pm1->yx; |
105 | | |
106 | 1.76M | pmr->xx = xx1 * xx2 + xy1 * yx2; |
107 | 1.76M | pmr->xy = xx1 * xy2 + xy1 * yy2; |
108 | 1.76M | pmr->yy = yx1 * xy2 + yy1 * yy2; |
109 | 1.76M | pmr->yx = yx1 * xx2 + yy1 * yx2; |
110 | 1.76M | pmr->tx = tx1 * xx2 + ty1 * yx2 + pm2->tx; |
111 | 1.76M | pmr->ty = tx1 * xy2 + ty1 * yy2 + pm2->ty; |
112 | 1.76M | } |
113 | 58.0M | return 0; |
114 | 58.0M | } |
115 | | int |
116 | | gs_matrix_multiply_double(const gs_matrix_double * pm1, const gs_matrix * pm2, gs_matrix_double * pmr) |
117 | 1.79M | { |
118 | 1.79M | double xx1 = pm1->xx, yy1 = pm1->yy; |
119 | 1.79M | double tx1 = pm1->tx, ty1 = pm1->ty; |
120 | 1.79M | double xx2 = pm2->xx, yy2 = pm2->yy; |
121 | 1.79M | double xy2 = pm2->xy, yx2 = pm2->yx; |
122 | | |
123 | 1.79M | if (is_xxyy(pm1)) { |
124 | 1.61M | pmr->tx = tx1 * xx2 + pm2->tx; |
125 | 1.61M | pmr->ty = ty1 * yy2 + pm2->ty; |
126 | 1.61M | if (is_fzero(xy2)) |
127 | 1.22M | pmr->xy = 0; |
128 | 384k | else |
129 | 384k | pmr->xy = xx1 * xy2, |
130 | 384k | pmr->ty += tx1 * xy2; |
131 | 1.61M | pmr->xx = xx1 * xx2; |
132 | 1.61M | if (is_fzero(yx2)) |
133 | 1.28M | pmr->yx = 0; |
134 | 326k | else |
135 | 326k | pmr->yx = yy1 * yx2, |
136 | 326k | pmr->tx += ty1 * yx2; |
137 | 1.61M | pmr->yy = yy1 * yy2; |
138 | 1.61M | } else { |
139 | 181k | double xy1 = pm1->xy, yx1 = pm1->yx; |
140 | | |
141 | 181k | pmr->xx = xx1 * xx2 + xy1 * yx2; |
142 | 181k | pmr->xy = xx1 * xy2 + xy1 * yy2; |
143 | 181k | pmr->yy = yx1 * xy2 + yy1 * yy2; |
144 | 181k | pmr->yx = yx1 * xx2 + yy1 * yx2; |
145 | 181k | pmr->tx = tx1 * xx2 + ty1 * yx2 + pm2->tx; |
146 | 181k | pmr->ty = tx1 * xy2 + ty1 * yy2 + pm2->ty; |
147 | 181k | } |
148 | 1.79M | return 0; |
149 | 1.79M | } |
150 | | |
151 | | /* Invert a matrix. Return gs_error_undefinedresult if not invertible. */ |
152 | | int |
153 | | gs_matrix_invert(const gs_matrix * pm, gs_matrix * pmr) |
154 | 36.3M | { /* We have to be careful about fetch/store order, */ |
155 | | /* because pm might be the same as pmr. */ |
156 | 36.3M | if (is_xxyy(pm)) { |
157 | 25.2M | if (is_fzero(pm->xx) || is_fzero(pm->yy)) |
158 | 24.6k | return_error(gs_error_undefinedresult); |
159 | 25.2M | pmr->tx = -(pmr->xx = 1.0 / pm->xx) * pm->tx; |
160 | 25.2M | pmr->xy = 0.0; |
161 | 25.2M | pmr->yx = 0.0; |
162 | 25.2M | pmr->ty = -(pmr->yy = 1.0 / pm->yy) * pm->ty; |
163 | 25.2M | } else { |
164 | 11.0M | float mxx = pm->xx, myy = pm->yy, mxy = pm->xy, myx = pm->yx; |
165 | 11.0M | float mtx = pm->tx, mty = pm->ty; |
166 | | /* we declare det as double since on at least some computer (i.e. peeves) |
167 | | declaring it as a float results in different values for pmr depending |
168 | | on whether or not optimization is turned on. I believe this is caused |
169 | | by the compiler keeping the det value in an internal register when |
170 | | optimization is enable. As evidence of this if you add a debugging |
171 | | statement to print out det the optimized code acts the same as the |
172 | | unoptimized code. declearing det as double does not change the CET 10-09.ps |
173 | | output. */ |
174 | 11.0M | double det = (float)(mxx * myy) - (float)(mxy * myx); |
175 | | |
176 | | /* |
177 | | * We are doing the math as floats instead of doubles to reproduce |
178 | | * the results in page 1 of CET 10-09.ps |
179 | | */ |
180 | 11.0M | if (det == 0) |
181 | 9.34k | return_error(gs_error_undefinedresult); |
182 | 11.0M | pmr->xx = myy / det; |
183 | 11.0M | pmr->xy = -mxy / det; |
184 | 11.0M | pmr->yx = -myx / det; |
185 | 11.0M | pmr->yy = mxx / det; |
186 | 11.0M | pmr->tx = (((float)(mty * myx) - (float)(mtx * myy))) / det; |
187 | 11.0M | pmr->ty = (((float)(mtx * mxy) - (float)(mty * mxx))) / det; |
188 | 11.0M | } |
189 | 36.3M | return 0; |
190 | 36.3M | } |
191 | | int |
192 | | gs_matrix_invert_to_double(const gs_matrix * pm, gs_matrix_double * pmr) |
193 | 1.79M | { /* We have to be careful about fetch/store order, */ |
194 | | /* because pm might be the same as pmr. */ |
195 | 1.79M | if (is_xxyy(pm)) { |
196 | 1.61M | if (is_fzero(pm->xx) || is_fzero(pm->yy)) |
197 | 121 | return_error(gs_error_undefinedresult); |
198 | 1.61M | pmr->tx = -(pmr->xx = 1.0 / pm->xx) * pm->tx; |
199 | 1.61M | pmr->xy = 0.0; |
200 | 1.61M | pmr->yx = 0.0; |
201 | 1.61M | pmr->ty = -(pmr->yy = 1.0 / pm->yy) * pm->ty; |
202 | 1.61M | } else { |
203 | 181k | double mxx = pm->xx, myy = pm->yy, mxy = pm->xy, myx = pm->yx; |
204 | 181k | double mtx = pm->tx, mty = pm->ty; |
205 | 181k | double det = (mxx * myy) - (mxy * myx); |
206 | | |
207 | | /* |
208 | | * We are doing the math as floats instead of doubles to reproduce |
209 | | * the results in page 1 of CET 10-09.ps |
210 | | */ |
211 | 181k | if (det == 0) |
212 | 4 | return_error(gs_error_undefinedresult); |
213 | 181k | pmr->xx = myy / det; |
214 | 181k | pmr->xy = -mxy / det; |
215 | 181k | pmr->yx = -myx / det; |
216 | 181k | pmr->yy = mxx / det; |
217 | 181k | pmr->tx = (((mty * myx) - (mtx * myy))) / det; |
218 | 181k | pmr->ty = (((mtx * mxy) - (mty * mxx))) / det; |
219 | 181k | } |
220 | 1.79M | return 0; |
221 | 1.79M | } |
222 | | |
223 | | /* Translate a matrix, possibly in place. */ |
224 | | int |
225 | | gs_matrix_translate(const gs_matrix * pm, double dx, double dy, gs_matrix * pmr) |
226 | 75.1k | { |
227 | 75.1k | gs_point trans; |
228 | 75.1k | int code = gs_distance_transform(dx, dy, pm, &trans); |
229 | | |
230 | 75.1k | if (code < 0) |
231 | 0 | return code; |
232 | 75.1k | if (pmr != pm) |
233 | 75.1k | *pmr = *pm; |
234 | 75.1k | pmr->tx += trans.x; |
235 | 75.1k | pmr->ty += trans.y; |
236 | 75.1k | return 0; |
237 | 75.1k | } |
238 | | |
239 | | /* Scale a matrix, possibly in place. */ |
240 | | int |
241 | | gs_matrix_scale(const gs_matrix * pm, double sx, double sy, gs_matrix * pmr) |
242 | 3.70M | { |
243 | 3.70M | pmr->xx = pm->xx * sx; |
244 | 3.70M | pmr->xy = pm->xy * sx; |
245 | 3.70M | pmr->yx = pm->yx * sy; |
246 | 3.70M | pmr->yy = pm->yy * sy; |
247 | 3.70M | if (pmr != pm) { |
248 | 185k | pmr->tx = pm->tx; |
249 | 185k | pmr->ty = pm->ty; |
250 | 185k | } |
251 | 3.70M | return 0; |
252 | 3.70M | } |
253 | | |
254 | | /* Rotate a matrix, possibly in place. The angle is in degrees. */ |
255 | | int |
256 | | gs_matrix_rotate(const gs_matrix * pm, double ang, gs_matrix * pmr) |
257 | 2.81M | { |
258 | 2.81M | double mxx, mxy; |
259 | 2.81M | gs_sincos_t sincos; |
260 | | |
261 | 2.81M | gs_sincos_degrees(ang, &sincos); |
262 | 2.81M | mxx = pm->xx, mxy = pm->xy; |
263 | 2.81M | pmr->xx = sincos.cos * mxx + sincos.sin * pm->yx; |
264 | 2.81M | pmr->xy = sincos.cos * mxy + sincos.sin * pm->yy; |
265 | 2.81M | pmr->yx = sincos.cos * pm->yx - sincos.sin * mxx; |
266 | 2.81M | pmr->yy = sincos.cos * pm->yy - sincos.sin * mxy; |
267 | 2.81M | if (pmr != pm) { |
268 | 0 | pmr->tx = pm->tx; |
269 | 0 | pmr->ty = pm->ty; |
270 | 0 | } |
271 | 2.81M | return 0; |
272 | 2.81M | } |
273 | | |
274 | | /* ------ Coordinate transformations (floating point) ------ */ |
275 | | |
276 | | /* Note that all the transformation routines take separate */ |
277 | | /* x and y arguments, but return their result in a point. */ |
278 | | |
279 | | /* Transform a point. */ |
280 | | int |
281 | | gs_point_transform(double x, double y, const gs_matrix * pmat, |
282 | | gs_point * ppt) |
283 | 2.60G | { |
284 | | /* |
285 | | * The float casts are there to reproduce results in CET 10-01.ps |
286 | | * page 4. |
287 | | */ |
288 | 2.60G | ppt->x = (float)(x * pmat->xx) + pmat->tx; |
289 | 2.60G | ppt->y = (float)(y * pmat->yy) + pmat->ty; |
290 | 2.60G | if (!is_fzero(pmat->yx)) |
291 | 148M | ppt->x += (float)(y * pmat->yx); |
292 | 2.60G | if (!is_fzero(pmat->xy)) |
293 | 1.97G | ppt->y += (float)(x * pmat->xy); |
294 | 2.60G | return 0; |
295 | 2.60G | } |
296 | | |
297 | | /* Inverse-transform a point. */ |
298 | | /* Return gs_error_undefinedresult if the matrix is not invertible. */ |
299 | | int |
300 | | gs_point_transform_inverse(double x, double y, const gs_matrix * pmat, |
301 | | gs_point * ppt) |
302 | 219M | { |
303 | 219M | if (is_xxyy(pmat)) { |
304 | 211M | if (is_fzero(pmat->xx) || is_fzero(pmat->yy)) |
305 | 6.12k | return_error(gs_error_undefinedresult); |
306 | 211M | ppt->x = (x - pmat->tx) / pmat->xx; |
307 | 211M | ppt->y = (y - pmat->ty) / pmat->yy; |
308 | 211M | return 0; |
309 | 211M | } else if (is_xyyx(pmat)) { |
310 | 1.02M | if (is_fzero(pmat->xy) || is_fzero(pmat->yx)) |
311 | 347 | return_error(gs_error_undefinedresult); |
312 | 1.02M | ppt->x = (y - pmat->ty) / pmat->xy; |
313 | 1.02M | ppt->y = (x - pmat->tx) / pmat->yx; |
314 | 1.02M | return 0; |
315 | 6.95M | } else { /* There are faster ways to do this, */ |
316 | | /* but we won't implement one unless we have to. */ |
317 | 6.95M | gs_matrix imat; |
318 | 6.95M | int code = gs_matrix_invert(pmat, &imat); |
319 | | |
320 | 6.95M | if (code < 0) |
321 | 679 | return code; |
322 | 6.95M | return gs_point_transform(x, y, &imat, ppt); |
323 | 6.95M | } |
324 | 219M | } |
325 | | |
326 | | /* Transform a distance. */ |
327 | | int |
328 | | gs_distance_transform(double dx, double dy, const gs_matrix * pmat, |
329 | | gs_point * pdpt) |
330 | 247M | { |
331 | 247M | pdpt->x = dx * pmat->xx; |
332 | 247M | pdpt->y = dy * pmat->yy; |
333 | 247M | if (!is_fzero(pmat->yx)) |
334 | 58.6M | pdpt->x += dy * pmat->yx; |
335 | 247M | if (!is_fzero(pmat->xy)) |
336 | 58.0M | pdpt->y += dx * pmat->xy; |
337 | 247M | return 0; |
338 | 247M | } |
339 | | |
340 | | /* Inverse-transform a distance. */ |
341 | | /* Return gs_error_undefinedresult if the matrix is not invertible. */ |
342 | | int |
343 | | gs_distance_transform_inverse(double dx, double dy, |
344 | | const gs_matrix * pmat, gs_point * pdpt) |
345 | 77.9M | { |
346 | 77.9M | if (is_xxyy(pmat)) { |
347 | 24.2M | if (is_fzero(pmat->xx) || is_fzero(pmat->yy)) |
348 | 64.2k | return_error(gs_error_undefinedresult); |
349 | 24.2M | pdpt->x = dx / pmat->xx; |
350 | 24.2M | pdpt->y = dy / pmat->yy; |
351 | 53.6M | } else if (is_xyyx(pmat)) { |
352 | 409k | if (is_fzero(pmat->xy) || is_fzero(pmat->yx)) |
353 | 45.9k | return_error(gs_error_undefinedresult); |
354 | 363k | pdpt->x = dy / pmat->xy; |
355 | 363k | pdpt->y = dx / pmat->yx; |
356 | 53.2M | } else { |
357 | 53.2M | double det = pmat->xx * pmat->yy - pmat->xy * pmat->yx; |
358 | | |
359 | 53.2M | if (det == 0) |
360 | 887k | return_error(gs_error_undefinedresult); |
361 | 52.3M | pdpt->x = (dx * pmat->yy - dy * pmat->yx) / det; |
362 | 52.3M | pdpt->y = (dy * pmat->xx - dx * pmat->xy) / det; |
363 | 52.3M | } |
364 | 76.9M | return 0; |
365 | 77.9M | } |
366 | | |
367 | | /* Compute the bounding box of 4 points. */ |
368 | | int |
369 | | gs_points_bbox(const gs_point pts[4], gs_rect * pbox) |
370 | 27.9M | { |
371 | 27.9M | #define assign_min_max(vmin, vmax, v0, v1)\ |
372 | 111M | if ( v0 < v1 ) vmin = v0, vmax = v1; else vmin = v1, vmax = v0 |
373 | 27.9M | #define assign_min_max_4(vmin, vmax, v0, v1, v2, v3)\ |
374 | 55.8M | { double min01, max01, min23, max23;\ |
375 | 55.8M | assign_min_max(min01, max01, v0, v1);\ |
376 | 55.8M | assign_min_max(min23, max23, v2, v3);\ |
377 | 55.8M | vmin = min(min01, min23);\ |
378 | 55.8M | vmax = max(max01, max23);\ |
379 | 55.8M | } |
380 | 27.9M | assign_min_max_4(pbox->p.x, pbox->q.x, |
381 | 27.9M | pts[0].x, pts[1].x, pts[2].x, pts[3].x); |
382 | 27.9M | assign_min_max_4(pbox->p.y, pbox->q.y, |
383 | 27.9M | pts[0].y, pts[1].y, pts[2].y, pts[3].y); |
384 | 27.9M | #undef assign_min_max |
385 | 27.9M | #undef assign_min_max_4 |
386 | 27.9M | return 0; |
387 | 27.9M | } |
388 | | |
389 | | /* Transform or inverse-transform a bounding box. */ |
390 | | /* Return gs_error_undefinedresult if the matrix is not invertible. */ |
391 | | static int |
392 | | bbox_transform_either_only(const gs_rect * pbox_in, const gs_matrix * pmat, |
393 | | gs_point pts[4], |
394 | | int (*point_xform) (double, double, const gs_matrix *, gs_point *)) |
395 | 27.9M | { |
396 | 27.9M | int code; |
397 | | |
398 | 27.9M | if ((code = (*point_xform) (pbox_in->p.x, pbox_in->p.y, pmat, &pts[0])) < 0 || |
399 | 27.9M | (code = (*point_xform) (pbox_in->p.x, pbox_in->q.y, pmat, &pts[1])) < 0 || |
400 | 27.9M | (code = (*point_xform) (pbox_in->q.x, pbox_in->p.y, pmat, &pts[2])) < 0 || |
401 | 27.9M | (code = (*point_xform) (pbox_in->q.x, pbox_in->q.y, pmat, &pts[3])) < 0 |
402 | 27.9M | ) |
403 | 27.9M | DO_NOTHING; |
404 | 27.9M | return code; |
405 | 27.9M | } |
406 | | |
407 | | static int |
408 | | bbox_transform_either(const gs_rect * pbox_in, const gs_matrix * pmat, |
409 | | gs_rect * pbox_out, |
410 | | int (*point_xform) (double, double, const gs_matrix *, gs_point *)) |
411 | 27.8M | { |
412 | 27.8M | int code; |
413 | | |
414 | | /* |
415 | | * In principle, we could transform only one point and two |
416 | | * distance vectors; however, because of rounding, we will only |
417 | | * get fully consistent results if we transform all 4 points. |
418 | | * We must compute the max and min after transforming, |
419 | | * since a rotation may be involved. |
420 | | */ |
421 | 27.8M | gs_point pts[4]; |
422 | | |
423 | 27.8M | if ((code = bbox_transform_either_only(pbox_in, pmat, pts, point_xform)) < 0) |
424 | 923 | return code; |
425 | 27.8M | return gs_points_bbox(pts, pbox_out); |
426 | 27.8M | } |
427 | | int |
428 | | gs_bbox_transform(const gs_rect * pbox_in, const gs_matrix * pmat, |
429 | | gs_rect * pbox_out) |
430 | 18.5M | { |
431 | 18.5M | return bbox_transform_either(pbox_in, pmat, pbox_out, |
432 | 18.5M | gs_point_transform); |
433 | 18.5M | } |
434 | | int |
435 | | gs_bbox_transform_only(const gs_rect * pbox_in, const gs_matrix * pmat, |
436 | | gs_point points[4]) |
437 | 28.8k | { |
438 | 28.8k | return bbox_transform_either_only(pbox_in, pmat, points, |
439 | 28.8k | gs_point_transform); |
440 | 28.8k | } |
441 | | int |
442 | | gs_bbox_transform_inverse(const gs_rect * pbox_in, const gs_matrix * pmat, |
443 | | gs_rect * pbox_out) |
444 | 9.36M | { |
445 | 9.36M | int code = bbox_transform_either(pbox_in, pmat, pbox_out, |
446 | 9.36M | gs_point_transform_inverse); |
447 | | |
448 | 9.36M | return code; |
449 | 9.36M | } |
450 | | |
451 | | /* ------ Coordinate transformations (to fixed point) ------ */ |
452 | | |
453 | 52.3M | #define f_fits_in_fixed(f) f_fits_in_bits(f, fixed_int_bits) |
454 | | |
455 | | /* Make a gs_matrix_fixed from a gs_matrix. */ |
456 | | int |
457 | | gs_matrix_fixed_from_matrix(gs_matrix_fixed *pfmat, const gs_matrix *pmat) |
458 | 25.1M | { |
459 | 25.1M | *(gs_matrix *)pfmat = *pmat; |
460 | 25.1M | if (f_fits_in_fixed(pmat->tx) && f_fits_in_fixed(pmat->ty)) { |
461 | 25.1M | pfmat->tx = fixed2float(pfmat->tx_fixed = float2fixed(pmat->tx)); |
462 | 25.1M | pfmat->ty = fixed2float(pfmat->ty_fixed = float2fixed(pmat->ty)); |
463 | 25.1M | pfmat->txy_fixed_valid = true; |
464 | 25.1M | } else { |
465 | 0 | pfmat->txy_fixed_valid = false; |
466 | 0 | } |
467 | 25.1M | return 0; |
468 | 25.1M | } |
469 | | |
470 | | /* Transform a point with a fixed-point result. */ |
471 | | int |
472 | | gs_point_transform2fixed(const gs_matrix_fixed * pmat, |
473 | | double x, double y, gs_fixed_point * ppt) |
474 | 15.0M | { |
475 | 15.0M | fixed px, py, t; |
476 | 15.0M | double xtemp, ytemp; |
477 | 15.0M | int code; |
478 | | |
479 | 15.0M | if (!pmat->txy_fixed_valid) { /* The translation is out of range. Do the */ |
480 | | /* computation in floating point, and convert to */ |
481 | | /* fixed at the end. */ |
482 | 272k | gs_point fpt; |
483 | | |
484 | 272k | gs_point_transform(x, y, (const gs_matrix *)pmat, &fpt); |
485 | 272k | if (!(f_fits_in_fixed(fpt.x) && f_fits_in_fixed(fpt.y))) |
486 | 347 | return_error(gs_error_limitcheck); |
487 | 271k | ppt->x = float2fixed(fpt.x); |
488 | 271k | ppt->y = float2fixed(fpt.y); |
489 | 271k | return 0; |
490 | 272k | } |
491 | 14.8M | if (!is_fzero(pmat->xy)) { /* Hope for 90 degree rotation */ |
492 | 176k | if ((code = CHECK_DFMUL2FIXED_VARS(px, y, pmat->yx, xtemp)) < 0 || |
493 | 176k | (code = CHECK_DFMUL2FIXED_VARS(py, x, pmat->xy, ytemp)) < 0 |
494 | 176k | ) |
495 | 10 | return code; |
496 | 176k | FINISH_DFMUL2FIXED_VARS(px, xtemp); |
497 | 176k | FINISH_DFMUL2FIXED_VARS(py, ytemp); |
498 | 176k | if (!is_fzero(pmat->xx)) { |
499 | 151k | if ((code = CHECK_DFMUL2FIXED_VARS(t, x, pmat->xx, xtemp)) < 0) |
500 | 0 | return code; |
501 | 151k | FINISH_DFMUL2FIXED_VARS(t, xtemp); |
502 | 151k | if ((code = CHECK_SET_FIXED_SUM(px, px, t)) < 0) |
503 | 0 | return code; |
504 | 151k | } |
505 | 176k | if (!is_fzero(pmat->yy)) { |
506 | 150k | if ((code = CHECK_DFMUL2FIXED_VARS(t, y, pmat->yy, ytemp)) < 0) |
507 | 0 | return code; |
508 | 150k | FINISH_DFMUL2FIXED_VARS(t, ytemp); |
509 | 150k | if ((code = CHECK_SET_FIXED_SUM(py, py, t)) < 0) |
510 | 0 | return code; |
511 | 150k | } |
512 | 14.6M | } else { |
513 | 14.6M | if ((code = CHECK_DFMUL2FIXED_VARS(px, x, pmat->xx, xtemp)) < 0 || |
514 | 14.6M | (code = CHECK_DFMUL2FIXED_VARS(py, y, pmat->yy, ytemp)) < 0 |
515 | 14.6M | ) |
516 | 7.30k | return code; |
517 | 14.6M | FINISH_DFMUL2FIXED_VARS(px, xtemp); |
518 | 14.6M | FINISH_DFMUL2FIXED_VARS(py, ytemp); |
519 | 14.6M | if (!is_fzero(pmat->yx)) { |
520 | 616 | if ((code = CHECK_DFMUL2FIXED_VARS(t, y, pmat->yx, ytemp)) < 0) |
521 | 0 | return code; |
522 | 616 | FINISH_DFMUL2FIXED_VARS(t, ytemp); |
523 | 616 | if ((code = CHECK_SET_FIXED_SUM(px, px, t)) < 0) |
524 | 0 | return code; |
525 | 616 | } |
526 | 14.6M | } |
527 | 14.8M | if (((code = CHECK_SET_FIXED_SUM(ppt->x, px, pmat->tx_fixed)) < 0) || |
528 | 14.8M | ((code = CHECK_SET_FIXED_SUM(ppt->y, py, pmat->ty_fixed)) < 0) ) |
529 | 85 | return code; |
530 | 14.8M | return 0; |
531 | 14.8M | } |
532 | | |
533 | | #if PRECISE_CURRENTPOINT |
534 | | /* Transform a point with a fixed-point result. */ |
535 | | /* Used for the best precision of the current point, |
536 | | see comment in clamp_point_aux. */ |
537 | | int |
538 | | gs_point_transform2fixed_rounding(const gs_matrix_fixed * pmat, |
539 | | double x, double y, gs_fixed_point * ppt) |
540 | 762k | { |
541 | 762k | gs_point fpt; |
542 | | |
543 | 762k | gs_point_transform(x, y, (const gs_matrix *)pmat, &fpt); |
544 | 762k | if (!(f_fits_in_fixed(fpt.x) && f_fits_in_fixed(fpt.y))) |
545 | 94 | return_error(gs_error_limitcheck); |
546 | 762k | ppt->x = float2fixed_rounded(fpt.x); |
547 | 762k | ppt->y = float2fixed_rounded(fpt.y); |
548 | 762k | return 0; |
549 | 762k | } |
550 | | #endif |
551 | | |
552 | | /* Transform a distance with a fixed-point result. */ |
553 | | int |
554 | | gs_distance_transform2fixed(const gs_matrix_fixed * pmat, |
555 | | double dx, double dy, gs_fixed_point * ppt) |
556 | 81.1M | { |
557 | 81.1M | fixed px, py, t; |
558 | 81.1M | double xtemp, ytemp; |
559 | 81.1M | int code; |
560 | | |
561 | 81.1M | if ((code = CHECK_DFMUL2FIXED_VARS(px, dx, pmat->xx, xtemp)) < 0 || |
562 | 81.1M | (code = CHECK_DFMUL2FIXED_VARS(py, dy, pmat->yy, ytemp)) < 0 |
563 | 81.1M | ) |
564 | 436k | return code; |
565 | 80.7M | FINISH_DFMUL2FIXED_VARS(px, xtemp); |
566 | 80.7M | FINISH_DFMUL2FIXED_VARS(py, ytemp); |
567 | 80.7M | if (!is_fzero(pmat->yx)) { |
568 | 1.93M | if ((code = CHECK_DFMUL2FIXED_VARS(t, dy, pmat->yx, ytemp)) < 0) |
569 | 36.2k | return code; |
570 | 1.90M | FINISH_DFMUL2FIXED_VARS(t, ytemp); |
571 | 1.90M | if ((code = CHECK_SET_FIXED_SUM(px, px, t)) < 0) |
572 | 12 | return code; |
573 | 1.90M | } |
574 | 80.7M | if (!is_fzero(pmat->xy)) { |
575 | 1.67M | if ((code = CHECK_DFMUL2FIXED_VARS(t, dx, pmat->xy, xtemp)) < 0) |
576 | 13.8k | return code; |
577 | 1.66M | FINISH_DFMUL2FIXED_VARS(t, xtemp); |
578 | 1.66M | if ((code = CHECK_SET_FIXED_SUM(py, py, t)) < 0) |
579 | 0 | return code; |
580 | 1.66M | } |
581 | 80.7M | ppt->x = px; |
582 | 80.7M | ppt->y = py; |
583 | 80.7M | return 0; |
584 | 80.7M | } |
585 | | |
586 | | /* ------ Serialization ------ */ |
587 | | |
588 | | /* |
589 | | * For maximum conciseness in band lists, we write a matrix as a control |
590 | | * byte followed by 0 to 6 values. The control byte has the format |
591 | | * AABBCD00. AA and BB control (xx,yy) and (xy,yx) as follows: |
592 | | * 00 = values are (0.0, 0.0) |
593 | | * 01 = values are (V, V) [1 value follows] |
594 | | * 10 = values are (V, -V) [1 value follows] |
595 | | * 11 = values are (U, V) [2 values follow] |
596 | | * C and D control tx and ty as follows: |
597 | | * 0 = value is 0.0 |
598 | | * 1 = value follows |
599 | | * The following code is the only place that knows this representation. |
600 | | */ |
601 | | |
602 | | /* Put a matrix on a stream. */ |
603 | | int |
604 | | sput_matrix(stream *s, const gs_matrix *pmat) |
605 | 11.1M | { |
606 | 11.1M | byte buf[1 + 6 * sizeof(float)]; |
607 | 11.1M | byte *cp = buf + 1; |
608 | 11.1M | byte b = 0; |
609 | 11.1M | float coeff[6]; |
610 | 11.1M | int i; |
611 | 11.1M | uint ignore; |
612 | | |
613 | 11.1M | coeff[0] = pmat->xx; |
614 | 11.1M | coeff[1] = pmat->xy; |
615 | 11.1M | coeff[2] = pmat->yx; |
616 | 11.1M | coeff[3] = pmat->yy; |
617 | 11.1M | coeff[4] = pmat->tx; |
618 | 11.1M | coeff[5] = pmat->ty; |
619 | 33.5M | for (i = 0; i < 4; i += 2) { |
620 | 22.3M | float u = coeff[i], v = coeff[i ^ 3]; |
621 | | |
622 | 22.3M | b <<= 2; |
623 | 22.3M | if (u != 0 || v != 0) { |
624 | 13.7M | memcpy(cp, &u, sizeof(float)); |
625 | 13.7M | cp += sizeof(float); |
626 | | |
627 | 13.7M | if (v == u) |
628 | 1.61M | b += 1; |
629 | 12.1M | else if (v == -u) |
630 | 5.18M | b += 2; |
631 | 6.94M | else { |
632 | 6.94M | b += 3; |
633 | 6.94M | memcpy(cp, &v, sizeof(float)); |
634 | 6.94M | cp += sizeof(float); |
635 | 6.94M | } |
636 | 13.7M | } |
637 | 22.3M | } |
638 | 33.5M | for (; i < 6; ++i) { |
639 | 22.3M | float v = coeff[i]; |
640 | | |
641 | 22.3M | b <<= 1; |
642 | 22.3M | if (v != 0) { |
643 | 19.7M | ++b; |
644 | 19.7M | memcpy(cp, &v, sizeof(float)); |
645 | 19.7M | cp += sizeof(float); |
646 | 19.7M | } |
647 | 22.3M | } |
648 | 11.1M | buf[0] = b << 2; |
649 | 11.1M | return sputs(s, buf, cp - buf, &ignore); |
650 | 11.1M | } |
651 | | |
652 | | /* Get a matrix from a stream. */ |
653 | | int |
654 | | sget_matrix(stream *s, gs_matrix *pmat) |
655 | 16.2M | { |
656 | 16.2M | int b = sgetc(s); |
657 | 16.2M | float coeff[6]; |
658 | 16.2M | int i; |
659 | 16.2M | int status; |
660 | 16.2M | uint nread; |
661 | | |
662 | 16.2M | if (b < 0) |
663 | 0 | return b; |
664 | 48.7M | for (i = 0; i < 4; i += 2, b <<= 2) |
665 | 32.5M | if (!(b & 0xc0)) |
666 | 15.0M | coeff[i] = coeff[i ^ 3] = 0.0; |
667 | 17.4M | else { |
668 | 17.4M | float value; |
669 | | |
670 | 17.4M | status = sgets(s, (byte *)&value, sizeof(value), &nread); |
671 | 17.4M | if (status < 0 && status != EOFC) |
672 | 0 | return_error(gs_error_ioerror); |
673 | 17.4M | coeff[i] = value; |
674 | 17.4M | switch ((b >> 6) & 3) { |
675 | 3.44M | case 1: |
676 | 3.44M | coeff[i ^ 3] = value; |
677 | 3.44M | break; |
678 | 8.40M | case 2: |
679 | 8.40M | coeff[i ^ 3] = -value; |
680 | 8.40M | break; |
681 | 5.62M | case 3: |
682 | 5.62M | status = sgets(s, (byte *)&coeff[i ^ 3], |
683 | 5.62M | sizeof(coeff[0]), &nread); |
684 | 5.62M | if (status < 0 && status != EOFC) |
685 | 0 | return_error(gs_error_ioerror); |
686 | 17.4M | } |
687 | 17.4M | } |
688 | 48.7M | for (; i < 6; ++i, b <<= 1) |
689 | 32.5M | if (b & 0x80) { |
690 | 24.5M | status = sgets(s, (byte *)&coeff[i], sizeof(coeff[0]), &nread); |
691 | 24.5M | if (status < 0 && status != EOFC) |
692 | 0 | return_error(gs_error_ioerror); |
693 | 24.5M | } else |
694 | 7.93M | coeff[i] = 0.0; |
695 | 16.2M | pmat->xx = coeff[0]; |
696 | 16.2M | pmat->xy = coeff[1]; |
697 | 16.2M | pmat->yx = coeff[2]; |
698 | 16.2M | pmat->yy = coeff[3]; |
699 | 16.2M | pmat->tx = coeff[4]; |
700 | 16.2M | pmat->ty = coeff[5]; |
701 | 16.2M | return 0; |
702 | 16.2M | } |
703 | | |
704 | | /* Compare two matrices */ |
705 | | int |
706 | 17.1M | gs_matrix_compare(const gs_matrix *pmat1, const gs_matrix *pmat2) { |
707 | 17.1M | if (pmat1->xx != pmat2->xx) |
708 | 242k | return(1); |
709 | 16.9M | if (pmat1->xy != pmat2->xy) |
710 | 1.14k | return(1); |
711 | 16.9M | if (pmat1->yx != pmat2->yx) |
712 | 130 | return(1); |
713 | 16.9M | if (pmat1->yy != pmat2->yy) |
714 | 16.6k | return(1); |
715 | 16.8M | if (pmat1->tx != pmat2->tx) |
716 | 12.6M | return(1); |
717 | 4.28M | if (pmat1->ty != pmat2->ty) |
718 | 67.6k | return(1); |
719 | 4.21M | return(0); |
720 | 4.28M | } |