1 /* Copyright (C) 1995, 1996, 1998, 1999 Aladdin Enterprises. All rights reserved.
2
3 This software is provided AS-IS with no warranty, either express or
4 implied.
5
6 This software is distributed under license and may not be copied,
7 modified or distributed except as expressly authorized under the terms
8 of the license contained in the file LICENSE in this distribution.
9
10 For more information about licensing, please refer to
11 http://www.ghostscript.com/licensing/. For information on
12 commercial licensing, go to http://www.artifex.com/licensing/ or
13 contact Artifex Software, Inc., 101 Lucas Valley Road #110,
14 San Rafael, CA 94903, U.S.A., +1(415)492-9861.
15 */
16
17 /* $Id: zdouble.c,v 1.5 2002/06/16 03:43:50 lpd Exp $ */
18 /* Double-precision floating point arithmetic operators */
19 #include "math_.h"
20 #include "memory_.h"
21 #include "string_.h"
22 #include "ctype_.h"
23 #include "ghost.h"
24 #include "gxfarith.h"
25 #include "oper.h"
26 #include "store.h"
27
28 /*
29 * Thanks to Jean-Pierre Demailly of the Institut Fourier of the
30 * Universit\'e de Grenoble I <demailly@fourier.grenet.fr> for proposing
31 * this package and for arranging the funding for its creation.
32 *
33 * These operators work with doubles represented as 8-byte strings. When
34 * applicable, they write their result into a string supplied as an argument.
35 * They also accept ints and reals as arguments.
36 */
37
38 /* Forward references */
39 private int double_params_result(os_ptr, int, double *);
40 private int double_params(os_ptr, int, double *);
41 private int double_result(i_ctx_t *, int, double);
42 private int double_unary(i_ctx_t *, double (*)(double));
43
44 #define dbegin_unary()\
45 os_ptr op = osp;\
46 double num;\
47 int code = double_params_result(op, 1, &num);\
48 \
49 if ( code < 0 )\
50 return code
51
52 #define dbegin_binary()\
53 os_ptr op = osp;\
54 double num[2];\
55 int code = double_params_result(op, 2, num);\
56 \
57 if ( code < 0 )\
58 return code
59
60 /* ------ Arithmetic ------ */
61
62 /* <dnum1> <dnum2> <dresult> .dadd <dresult> */
63 private int
zdadd(i_ctx_t * i_ctx_p)64 zdadd(i_ctx_t *i_ctx_p)
65 {
66 dbegin_binary();
67 return double_result(i_ctx_p, 2, num[0] + num[1]);
68 }
69
70 /* <dnum1> <dnum2> <dresult> .ddiv <dresult> */
71 private int
zddiv(i_ctx_t * i_ctx_p)72 zddiv(i_ctx_t *i_ctx_p)
73 {
74 dbegin_binary();
75 if (num[1] == 0.0)
76 return_error(e_undefinedresult);
77 return double_result(i_ctx_p, 2, num[0] / num[1]);
78 }
79
80 /* <dnum1> <dnum2> <dresult> .dmul <dresult> */
81 private int
zdmul(i_ctx_t * i_ctx_p)82 zdmul(i_ctx_t *i_ctx_p)
83 {
84 dbegin_binary();
85 return double_result(i_ctx_p, 2, num[0] * num[1]);
86 }
87
88 /* <dnum1> <dnum2> <dresult> .dsub <dresult> */
89 private int
zdsub(i_ctx_t * i_ctx_p)90 zdsub(i_ctx_t *i_ctx_p)
91 {
92 dbegin_binary();
93 return double_result(i_ctx_p, 2, num[0] - num[1]);
94 }
95
96 /* ------ Simple functions ------ */
97
98 /* <dnum> <dresult> .dabs <dresult> */
99 private int
zdabs(i_ctx_t * i_ctx_p)100 zdabs(i_ctx_t *i_ctx_p)
101 {
102 return double_unary(i_ctx_p, fabs);
103 }
104
105 /* <dnum> <dresult> .dceiling <dresult> */
106 private int
zdceiling(i_ctx_t * i_ctx_p)107 zdceiling(i_ctx_t *i_ctx_p)
108 {
109 return double_unary(i_ctx_p, ceil);
110 }
111
112 /* <dnum> <dresult> .dfloor <dresult> */
113 private int
zdfloor(i_ctx_t * i_ctx_p)114 zdfloor(i_ctx_t *i_ctx_p)
115 {
116 return double_unary(i_ctx_p, floor);
117 }
118
119 /* <dnum> <dresult> .dneg <dresult> */
120 private int
zdneg(i_ctx_t * i_ctx_p)121 zdneg(i_ctx_t *i_ctx_p)
122 {
123 dbegin_unary();
124 return double_result(i_ctx_p, 1, -num);
125 }
126
127 /* <dnum> <dresult> .dround <dresult> */
128 private int
zdround(i_ctx_t * i_ctx_p)129 zdround(i_ctx_t *i_ctx_p)
130 {
131 dbegin_unary();
132 return double_result(i_ctx_p, 1, floor(num + 0.5));
133 }
134
135 /* <dnum> <dresult> .dsqrt <dresult> */
136 private int
zdsqrt(i_ctx_t * i_ctx_p)137 zdsqrt(i_ctx_t *i_ctx_p)
138 {
139 dbegin_unary();
140 if (num < 0.0)
141 return_error(e_rangecheck);
142 return double_result(i_ctx_p, 1, sqrt(num));
143 }
144
145 /* <dnum> <dresult> .dtruncate <dresult> */
146 private int
zdtruncate(i_ctx_t * i_ctx_p)147 zdtruncate(i_ctx_t *i_ctx_p)
148 {
149 dbegin_unary();
150 return double_result(i_ctx_p, 1, (num < 0 ? ceil(num) : floor(num)));
151 }
152
153 /* ------ Transcendental functions ------ */
154
155 private int
darc(i_ctx_t * i_ctx_p,double (* afunc)(double))156 darc(i_ctx_t *i_ctx_p, double (*afunc)(double))
157 {
158 dbegin_unary();
159 return double_result(i_ctx_p, 1, (*afunc)(num) * radians_to_degrees);
160 }
161 /* <dnum> <dresult> .darccos <dresult> */
162 private int
zdarccos(i_ctx_t * i_ctx_p)163 zdarccos(i_ctx_t *i_ctx_p)
164 {
165 return darc(i_ctx_p, acos);
166 }
167 /* <dnum> <dresult> .darcsin <dresult> */
168 private int
zdarcsin(i_ctx_t * i_ctx_p)169 zdarcsin(i_ctx_t *i_ctx_p)
170 {
171 return darc(i_ctx_p, asin);
172 }
173
174 /* <dnum> <ddenom> <dresult> .datan <dresult> */
175 private int
zdatan(i_ctx_t * i_ctx_p)176 zdatan(i_ctx_t *i_ctx_p)
177 {
178 double result;
179
180 dbegin_binary();
181 if (num[0] == 0) { /* on X-axis, special case */
182 if (num[1] == 0)
183 return_error(e_undefinedresult);
184 result = (num[1] < 0 ? 180 : 0);
185 } else {
186 result = atan2(num[0], num[1]) * radians_to_degrees;
187 if (result < 0)
188 result += 360;
189 }
190 return double_result(i_ctx_p, 2, result);
191 }
192
193 /* <dnum> <dresult> .dcos <dresult> */
194 private int
zdcos(i_ctx_t * i_ctx_p)195 zdcos(i_ctx_t *i_ctx_p)
196 {
197 return double_unary(i_ctx_p, gs_cos_degrees);
198 }
199
200 /* <dbase> <dexponent> <dresult> .dexp <dresult> */
201 private int
zdexp(i_ctx_t * i_ctx_p)202 zdexp(i_ctx_t *i_ctx_p)
203 {
204 double ipart;
205
206 dbegin_binary();
207 if (num[0] == 0.0 && num[1] == 0.0)
208 return_error(e_undefinedresult);
209 if (num[0] < 0.0 && modf(num[1], &ipart) != 0.0)
210 return_error(e_undefinedresult);
211 return double_result(i_ctx_p, 2, pow(num[0], num[1]));
212 }
213
214 private int
dlog(i_ctx_t * i_ctx_p,double (* lfunc)(double))215 dlog(i_ctx_t *i_ctx_p, double (*lfunc)(double))
216 {
217 dbegin_unary();
218 if (num <= 0.0)
219 return_error(e_rangecheck);
220 return double_result(i_ctx_p, 1, (*lfunc)(num));
221 }
222 /* <dposnum> <dresult> .dln <dresult> */
223 private int
zdln(i_ctx_t * i_ctx_p)224 zdln(i_ctx_t *i_ctx_p)
225 {
226 return dlog(i_ctx_p, log);
227 }
228 /* <dposnum> <dresult> .dlog <dresult> */
229 private int
zdlog(i_ctx_t * i_ctx_p)230 zdlog(i_ctx_t *i_ctx_p)
231 {
232 return dlog(i_ctx_p, log10);
233 }
234
235 /* <dnum> <dresult> .dsin <dresult> */
236 private int
zdsin(i_ctx_t * i_ctx_p)237 zdsin(i_ctx_t *i_ctx_p)
238 {
239 return double_unary(i_ctx_p, gs_sin_degrees);
240 }
241
242 /* ------ Comparison ------ */
243
244 private int
dcompare(i_ctx_t * i_ctx_p,int mask)245 dcompare(i_ctx_t *i_ctx_p, int mask)
246 {
247 os_ptr op = osp;
248 double num[2];
249 int code = double_params(op, 2, num);
250
251 if (code < 0)
252 return code;
253 make_bool(op - 1,
254 (mask & (num[0] < num[1] ? 1 : num[0] > num[1] ? 4 : 2))
255 != 0);
256 pop(1);
257 return 0;
258 }
259 /* <dnum1> <dnum2> .deq <bool> */
260 private int
zdeq(i_ctx_t * i_ctx_p)261 zdeq(i_ctx_t *i_ctx_p)
262 {
263 return dcompare(i_ctx_p, 2);
264 }
265 /* <dnum1> <dnum2> .dge <bool> */
266 private int
zdge(i_ctx_t * i_ctx_p)267 zdge(i_ctx_t *i_ctx_p)
268 {
269 return dcompare(i_ctx_p, 6);
270 }
271 /* <dnum1> <dnum2> .dgt <bool> */
272 private int
zdgt(i_ctx_t * i_ctx_p)273 zdgt(i_ctx_t *i_ctx_p)
274 {
275 return dcompare(i_ctx_p, 4);
276 }
277 /* <dnum1> <dnum2> .dle <bool> */
278 private int
zdle(i_ctx_t * i_ctx_p)279 zdle(i_ctx_t *i_ctx_p)
280 {
281 return dcompare(i_ctx_p, 3);
282 }
283 /* <dnum1> <dnum2> .dlt <bool> */
284 private int
zdlt(i_ctx_t * i_ctx_p)285 zdlt(i_ctx_t *i_ctx_p)
286 {
287 return dcompare(i_ctx_p, 1);
288 }
289 /* <dnum1> <dnum2> .dne <bool> */
290 private int
zdne(i_ctx_t * i_ctx_p)291 zdne(i_ctx_t *i_ctx_p)
292 {
293 return dcompare(i_ctx_p, 5);
294 }
295
296 /* ------ Conversion ------ */
297
298 /* Take the easy way out.... */
299 #define MAX_CHARS 50
300
301 /* <dnum> <dresult> .cvd <dresult> */
302 private int
zcvd(i_ctx_t * i_ctx_p)303 zcvd(i_ctx_t *i_ctx_p)
304 {
305 dbegin_unary();
306 return double_result(i_ctx_p, 1, num);
307 }
308
309 /* <string> <dresult> .cvsd <dresult> */
310 private int
zcvsd(i_ctx_t * i_ctx_p)311 zcvsd(i_ctx_t *i_ctx_p)
312 {
313 os_ptr op = osp;
314 int code = double_params_result(op, 0, NULL);
315 double num;
316 char buf[MAX_CHARS + 2];
317 char *str = buf;
318 uint len;
319 char end;
320
321 if (code < 0)
322 return code;
323 check_read_type(op[-1], t_string);
324 len = r_size(op - 1);
325 if (len > MAX_CHARS)
326 return_error(e_limitcheck);
327 memcpy(str, op[-1].value.bytes, len);
328 /*
329 * We check syntax in the following way: we remove whitespace,
330 * verify that the string contains only [0123456789+-.dDeE],
331 * then append a $ and then check that the next character after
332 * the scanned number is a $.
333 */
334 while (len > 0 && isspace(*str))
335 ++str, --len;
336 while (len > 0 && isspace(str[len - 1]))
337 --len;
338 str[len] = 0;
339 if (strspn(str, "0123456789+-.dDeE") != len)
340 return_error(e_syntaxerror);
341 strcat(str, "$");
342 if (sscanf(str, "%lf%c", &num, &end) != 2 || end != '$')
343 return_error(e_syntaxerror);
344 return double_result(i_ctx_p, 1, num);
345 }
346
347 /* <dnum> .dcvi <int> */
348 private int
zdcvi(i_ctx_t * i_ctx_p)349 zdcvi(i_ctx_t *i_ctx_p)
350 {
351 os_ptr op = osp;
352 #define alt_min_long (-1L << (arch_sizeof_long * 8 - 1))
353 #define alt_max_long (~(alt_min_long))
354 static const double min_int_real = (alt_min_long * 1.0 - 1);
355 static const double max_int_real = (alt_max_long * 1.0 + 1);
356 double num;
357 int code = double_params(op, 1, &num);
358
359 if (code < 0)
360 return code;
361
362 if (num < min_int_real || num > max_int_real)
363 return_error(e_rangecheck);
364 make_int(op, (long)num); /* truncates toward 0 */
365 return 0;
366 }
367
368 /* <dnum> .dcvr <real> */
369 private int
zdcvr(i_ctx_t * i_ctx_p)370 zdcvr(i_ctx_t *i_ctx_p)
371 {
372 os_ptr op = osp;
373 #define b30 (0x40000000L * 1.0)
374 #define max_mag (0xffffff * b30 * b30 * b30 * 0x4000)
375 static const float min_real = -max_mag;
376 static const float max_real = max_mag;
377 #undef b30
378 #undef max_mag
379 double num;
380 int code = double_params(op, 1, &num);
381
382 if (code < 0)
383 return code;
384 if (num < min_real || num > max_real)
385 return_error(e_rangecheck);
386 make_real(op, (float)num);
387 return 0;
388 }
389
390 /* <dnum> <string> .dcvs <substring> */
391 private int
zdcvs(i_ctx_t * i_ctx_p)392 zdcvs(i_ctx_t *i_ctx_p)
393 {
394 os_ptr op = osp;
395 double num;
396 int code = double_params(op - 1, 1, &num);
397 char str[MAX_CHARS + 1];
398 int len;
399
400 if (code < 0)
401 return code;
402 check_write_type(*op, t_string);
403 /*
404 * To get fully accurate output results for IEEE double-
405 * precision floats (53 bits of mantissa), the ANSI
406 * %g default of 6 digits is not enough; 16 are needed.
407 * Unfortunately, using %.16g produces unfortunate artifacts such as
408 * 1.2 printing as 1.200000000000005. Therefore, we print using %g,
409 * and if the result isn't accurate enough, print again
410 * using %.16g.
411 */
412 {
413 double scanned;
414
415 sprintf(str, "%g", num);
416 sscanf(str, "%lf", &scanned);
417 if (scanned != num)
418 sprintf(str, "%.16g", num);
419 }
420 len = strlen(str);
421 if (len > r_size(op))
422 return_error(e_rangecheck);
423 memcpy(op->value.bytes, str, len);
424 op[-1] = *op;
425 r_set_size(op - 1, len);
426 pop(1);
427 return 0;
428 }
429
430 /* ------ Initialization table ------ */
431
432 /* We need to split the table because of the 16-element limit. */
433 const op_def zdouble1_op_defs[] = {
434 /* Arithmetic */
435 {"3.dadd", zdadd},
436 {"3.ddiv", zddiv},
437 {"3.dmul", zdmul},
438 {"3.dsub", zdsub},
439 /* Comparison */
440 {"2.deq", zdeq},
441 {"2.dge", zdge},
442 {"2.dgt", zdgt},
443 {"2.dle", zdle},
444 {"2.dlt", zdlt},
445 {"2.dne", zdne},
446 /* Conversion */
447 {"2.cvd", zcvd},
448 {"2.cvsd", zcvsd},
449 {"1.dcvi", zdcvi},
450 {"1.dcvr", zdcvr},
451 {"2.dcvs", zdcvs},
452 op_def_end(0)
453 };
454 const op_def zdouble2_op_defs[] = {
455 /* Simple functions */
456 {"2.dabs", zdabs},
457 {"2.dceiling", zdceiling},
458 {"2.dfloor", zdfloor},
459 {"2.dneg", zdneg},
460 {"2.dround", zdround},
461 {"2.dsqrt", zdsqrt},
462 {"2.dtruncate", zdtruncate},
463 /* Transcendental functions */
464 {"2.darccos", zdarccos},
465 {"2.darcsin", zdarcsin},
466 {"3.datan", zdatan},
467 {"2.dcos", zdcos},
468 {"3.dexp", zdexp},
469 {"2.dln", zdln},
470 {"2.dlog", zdlog},
471 {"2.dsin", zdsin},
472 op_def_end(0)
473 };
474
475 /* ------ Internal procedures ------ */
476
477 /* Get some double arguments. */
478 private int
double_params(os_ptr op,int count,double * pval)479 double_params(os_ptr op, int count, double *pval)
480 {
481 pval += count;
482 while (--count >= 0) {
483 switch (r_type(op)) {
484 case t_real:
485 *--pval = op->value.realval;
486 break;
487 case t_integer:
488 *--pval = op->value.intval;
489 break;
490 case t_string:
491 if (!r_has_attr(op, a_read) ||
492 r_size(op) != sizeof(double)
493 )
494 return_error(e_typecheck);
495 --pval;
496 memcpy(pval, op->value.bytes, sizeof(double));
497 break;
498 case t__invalid:
499 return_error(e_stackunderflow);
500 default:
501 return_error(e_typecheck);
502 }
503 op--;
504 }
505 return 0;
506 }
507
508 /* Get some double arguments, and check for a double result. */
509 private int
double_params_result(os_ptr op,int count,double * pval)510 double_params_result(os_ptr op, int count, double *pval)
511 {
512 check_write_type(*op, t_string);
513 if (r_size(op) != sizeof(double))
514 return_error(e_typecheck);
515 return double_params(op - 1, count, pval);
516 }
517
518 /* Return a double result. */
519 private int
double_result(i_ctx_t * i_ctx_p,int count,double result)520 double_result(i_ctx_t *i_ctx_p, int count, double result)
521 {
522 os_ptr op = osp;
523 os_ptr op1 = op - count;
524
525 ref_assign_inline(op1, op);
526 memcpy(op1->value.bytes, &result, sizeof(double));
527 pop(count);
528 return 0;
529 }
530
531 /* Apply a unary function to a double operand. */
532 private int
double_unary(i_ctx_t * i_ctx_p,double (* func)(double))533 double_unary(i_ctx_t *i_ctx_p, double (*func)(double))
534 {
535 dbegin_unary();
536 return double_result(i_ctx_p, 1, (*func)(num));
537 }
538