1 /* $NetBSD: ieee754io.c,v 1.6 2024/08/18 20:47:17 christos Exp $ */ 2 3 /* 4 * /src/NTP/ntp4-dev/libntp/ieee754io.c,v 4.12 2005/04/16 17:32:10 kardel RELEASE_20050508_A 5 * 6 * ieee754io.c,v 4.12 2005/04/16 17:32:10 kardel RELEASE_20050508_A 7 * 8 * $Created: Sun Jul 13 09:12:02 1997 $ 9 * 10 * Copyright (c) 1997-2005 by Frank Kardel <kardel <AT> ntp.org> 11 * 12 * Redistribution and use in source and binary forms, with or without 13 * modification, are permitted provided that the following conditions 14 * are met: 15 * 1. Redistributions of source code must retain the above copyright 16 * notice, this list of conditions and the following disclaimer. 17 * 2. Redistributions in binary form must reproduce the above copyright 18 * notice, this list of conditions and the following disclaimer in the 19 * documentation and/or other materials provided with the distribution. 20 * 3. Neither the name of the author nor the names of its contributors 21 * may be used to endorse or promote products derived from this software 22 * without specific prior written permission. 23 * 24 * THIS SOFTWARE IS PROVIDED BY THE AUTHOR AND CONTRIBUTORS ``AS IS'' AND 25 * ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE 26 * IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE 27 * ARE DISCLAIMED. IN NO EVENT SHALL THE AUTHOR OR CONTRIBUTORS BE LIABLE 28 * FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL 29 * DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS 30 * OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) 31 * HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT 32 * LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY 33 * OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF 34 * SUCH DAMAGE. 35 * 36 */ 37 38 #ifdef HAVE_CONFIG_H 39 #include "config.h" 40 #endif 41 42 #include <stdio.h> 43 #include "l_stdlib.h" 44 #include "ntp_stdlib.h" 45 #include "ntp_fp.h" 46 #include "ieee754io.h" 47 48 static unsigned char get_byte (unsigned char *, offsets_t, int *); 49 #ifdef __not_yet__ 50 static void put_byte (unsigned char *, offsets_t, int *, unsigned char); 51 #endif 52 53 #ifdef LIBDEBUG 54 55 static char * 56 fmt_blong( 57 unsigned long val, 58 int cnt 59 ) 60 { 61 char *buf, *s; 62 int i = cnt; 63 64 val <<= 32 - cnt; 65 LIB_GETBUF(buf); 66 s = buf; 67 68 while (i--) 69 { 70 if (val & 0x80000000) 71 { 72 *s++ = '1'; 73 } 74 else 75 { 76 *s++ = '0'; 77 } 78 val <<= 1; 79 } 80 *s = '\0'; 81 return buf; 82 } 83 84 static char * 85 fmt_flt( 86 unsigned int sign, 87 unsigned long mh, 88 unsigned long ml, 89 unsigned long ch 90 ) 91 { 92 char *buf; 93 94 LIB_GETBUF(buf); 95 snprintf(buf, LIB_BUFLENGTH, "%c %s %s %s", sign ? '-' : '+', 96 fmt_blong(ch, 11), 97 fmt_blong(mh, 20), 98 fmt_blong(ml, 32)); 99 100 return buf; 101 } 102 103 static char * 104 fmt_hex( 105 unsigned char *bufp, 106 int length 107 ) 108 { 109 char * buf; 110 char hex[4]; 111 int i; 112 113 LIB_GETBUF(buf); 114 buf[0] = '\0'; 115 for (i = 0; i < length; i++) { 116 snprintf(hex, sizeof(hex), "%02x", bufp[i]); 117 strlcat(buf, hex, LIB_BUFLENGTH); 118 } 119 120 return buf; 121 } 122 123 #endif 124 125 static unsigned char 126 get_byte( 127 unsigned char *bufp, 128 offsets_t offset, 129 int *fieldindex 130 ) 131 { 132 unsigned char val; 133 134 val = *(bufp + offset[*fieldindex]); 135 #ifdef LIBDEBUG 136 if (debug > 4) 137 printf("fetchieee754: getbyte(0x%08x, %d) = 0x%02x\n", (unsigned int)(bufp)+offset[*fieldindex], *fieldindex, val); 138 #endif 139 (*fieldindex)++; 140 return val; 141 } 142 143 #ifdef __not_yet__ 144 static void 145 put_byte( 146 unsigned char *bufp, 147 offsets_t offsets, 148 int *fieldindex, 149 unsigned char val 150 ) 151 { 152 *(bufp + offsets[*fieldindex]) = val; 153 (*fieldindex)++; 154 } 155 #endif 156 157 /* 158 * make conversions to and from external IEEE754 formats and internal 159 * NTP FP format. 160 */ 161 int 162 fetch_ieee754( 163 unsigned char **buffpp, 164 int size, 165 l_fp *lfpp, 166 offsets_t offsets 167 ) 168 { 169 unsigned char *bufp = *buffpp; 170 unsigned int sign; 171 unsigned int bias; 172 unsigned int maxexp; 173 int mbits; 174 u_long mantissa_low; 175 u_long mantissa_high; 176 u_long characteristic; 177 long exponent; 178 #ifdef LIBDEBUG 179 int length; 180 #endif 181 unsigned char val; 182 int fieldindex = 0; 183 184 switch (size) 185 { 186 case IEEE_DOUBLE: 187 #ifdef LIBDEBUG 188 length = 8; 189 #endif 190 mbits = 52; 191 bias = 1023; 192 maxexp = 2047; 193 break; 194 195 case IEEE_SINGLE: 196 #ifdef LIBDEBUG 197 length = 4; 198 #endif 199 mbits = 23; 200 bias = 127; 201 maxexp = 255; 202 break; 203 204 default: 205 return IEEE_BADCALL; 206 } 207 208 val = get_byte(bufp, offsets, &fieldindex); /* fetch sign byte & first part of characteristic */ 209 210 sign = (val & 0x80) != 0; 211 characteristic = (val & 0x7F); 212 213 val = get_byte(bufp, offsets, &fieldindex); /* fetch rest of characteristic and start of mantissa */ 214 215 switch (size) 216 { 217 case IEEE_SINGLE: 218 characteristic <<= 1; 219 characteristic |= (val & 0x80) != 0; /* grab last characteristic bit */ 220 221 mantissa_high = 0; 222 223 mantissa_low = (val &0x7F) << 16; 224 mantissa_low |= (u_long)get_byte(bufp, offsets, &fieldindex) << 8; 225 mantissa_low |= get_byte(bufp, offsets, &fieldindex); 226 break; 227 228 case IEEE_DOUBLE: 229 characteristic <<= 4; 230 characteristic |= (val & 0xF0) >> 4; /* grab lower characteristic bits */ 231 232 mantissa_high = (val & 0x0F) << 16; 233 mantissa_high |= (u_long)get_byte(bufp, offsets, &fieldindex) << 8; 234 mantissa_high |= get_byte(bufp, offsets, &fieldindex); 235 236 mantissa_low = (u_long)get_byte(bufp, offsets, &fieldindex) << 24; 237 mantissa_low |= (u_long)get_byte(bufp, offsets, &fieldindex) << 16; 238 mantissa_low |= (u_long)get_byte(bufp, offsets, &fieldindex) << 8; 239 mantissa_low |= get_byte(bufp, offsets, &fieldindex); 240 break; 241 242 default: 243 return IEEE_BADCALL; 244 } 245 #ifdef LIBDEBUG 246 if (debug > 4) 247 { 248 double d; 249 float f; 250 251 if (size == IEEE_SINGLE) 252 { 253 int i; 254 255 for (i = 0; i < length; i++) 256 { 257 *((unsigned char *)(&f)+i) = *(*buffpp + offsets[i]); 258 } 259 d = f; 260 } 261 else 262 { 263 int i; 264 265 for (i = 0; i < length; i++) 266 { 267 *((unsigned char *)(&d)+i) = *(*buffpp + offsets[i]); 268 } 269 } 270 271 printf("fetchieee754: FP: %s -> %s -> %e(=%s)\n", fmt_hex(*buffpp, length), 272 fmt_flt(sign, mantissa_high, mantissa_low, characteristic), 273 d, fmt_hex((unsigned char *)&d, length)); 274 } 275 #endif 276 277 *buffpp += fieldindex; 278 279 /* 280 * detect funny numbers 281 */ 282 if (characteristic == maxexp) 283 { 284 /* 285 * NaN or Infinity 286 */ 287 if (mantissa_low || mantissa_high) 288 { 289 /* 290 * NaN 291 */ 292 return IEEE_NAN; 293 } 294 else 295 { 296 /* 297 * +Inf or -Inf 298 */ 299 return sign ? IEEE_NEGINFINITY : IEEE_POSINFINITY; 300 } 301 } 302 else 303 { 304 /* 305 * collect real numbers 306 */ 307 308 L_CLR(lfpp); 309 310 /* 311 * check for overflows 312 */ 313 exponent = characteristic - bias; 314 315 if (exponent > 31) /* sorry - hardcoded */ 316 { 317 /* 318 * overflow only in respect to NTP-FP representation 319 */ 320 return sign ? IEEE_NEGOVERFLOW : IEEE_POSOVERFLOW; 321 } 322 else 323 { 324 int frac_offset; /* where the fraction starts */ 325 326 frac_offset = mbits - exponent; 327 328 if (characteristic == 0) 329 { 330 /* 331 * de-normalized or tiny number - fits only as 0 332 */ 333 return IEEE_OK; 334 } 335 else 336 { 337 /* 338 * adjust for implied 1 339 */ 340 if (mbits > 31) 341 mantissa_high |= 1 << (mbits - 32); 342 else 343 mantissa_low |= 1 << mbits; 344 345 /* 346 * take mantissa apart - if only all machine would support 347 * 64 bit operations 8-( 348 */ 349 if (frac_offset > mbits) 350 { 351 lfpp->l_ui = 0; /* only fractional number */ 352 frac_offset -= mbits + 1; /* will now contain right shift count - 1*/ 353 if (mbits > 31) 354 { 355 lfpp->l_uf = mantissa_high << (63 - mbits); 356 lfpp->l_uf |= mantissa_low >> (mbits - 33); 357 lfpp->l_uf >>= frac_offset; 358 } 359 else 360 { 361 lfpp->l_uf = mantissa_low >> frac_offset; 362 } 363 } 364 else 365 { 366 if (frac_offset > 32) 367 { 368 /* 369 * must split in high word 370 */ 371 lfpp->l_ui = mantissa_high >> (frac_offset - 32); 372 lfpp->l_uf = (mantissa_high & ((1 << (frac_offset - 32)) - 1)) << (64 - frac_offset); 373 lfpp->l_uf |= mantissa_low >> (frac_offset - 32); 374 } 375 else 376 { 377 /* 378 * must split in low word 379 */ 380 lfpp->l_ui = mantissa_high << (32 - frac_offset); 381 lfpp->l_ui |= (mantissa_low >> frac_offset) & ((1 << (32 - frac_offset)) - 1); 382 lfpp->l_uf = (mantissa_low & ((1 << frac_offset) - 1)) << (32 - frac_offset); 383 } 384 } 385 386 /* 387 * adjust for sign 388 */ 389 if (sign) 390 { 391 L_NEG(lfpp); 392 } 393 394 return IEEE_OK; 395 } 396 } 397 } 398 } 399 400 /* 401 * DLH: This function is currently unused in ntpd. If you think about 402 * using it, be sure it does what you intend. I notice the bufpp arg 403 * is never referenced, and the calculated mantissa_high & mantissa_low 404 * are only referenced in debug output. It seems they're supposed to 405 * be composed into an ieee754-format float and stored at *bufpp or 406 * possibly **bufpp. Brought to my attention by this: 407 * 408 * ieee754io.c:414:10: warning: variable 'mantissa_low' set but not used 409 * [-Wunused-but-set-variable] 410 * 411 * To quiet it I'm #ifdef'ing the function away for now, here and below 412 * the call to it in main(). 413 */ 414 #ifdef PUT_IEEE754_UNUSED_FUNC 415 int 416 put_ieee754( 417 unsigned char **bufpp, 418 int size, 419 l_fp *lfpp, 420 offsets_t offsets 421 ) 422 { 423 l_fp outlfp; 424 #ifdef LIBDEBUG 425 unsigned int sign; 426 unsigned int bias; 427 #endif 428 /*unsigned int maxexp;*/ 429 int mbits; 430 int msb; 431 u_long mantissa_low = 0; 432 u_long mantissa_high = 0; 433 #ifdef LIBDEBUG 434 u_long characteristic = 0; 435 long exponent; 436 #endif 437 /*int length;*/ 438 unsigned long mask; 439 440 outlfp = *lfpp; 441 442 switch (size) 443 { 444 case IEEE_DOUBLE: 445 /*length = 8;*/ 446 mbits = 52; 447 #ifdef LIBDEBUG 448 bias = 1023; 449 #endif 450 /*maxexp = 2047;*/ 451 break; 452 453 case IEEE_SINGLE: 454 /*length = 4;*/ 455 mbits = 23; 456 #ifdef LIBDEBUG 457 bias = 127; 458 #endif 459 /*maxexp = 255;*/ 460 break; 461 462 default: 463 return IEEE_BADCALL; 464 } 465 466 /* 467 * find sign 468 */ 469 if (L_ISNEG(&outlfp)) 470 { 471 L_NEG(&outlfp); 472 #ifdef LIBDEBUG 473 sign = 1; 474 #endif 475 } 476 else 477 { 478 #ifdef LIBDEBUG 479 sign = 0; 480 #endif 481 } 482 483 if (L_ISZERO(&outlfp)) 484 { 485 #ifdef LIBDEBUG 486 exponent = mantissa_high = mantissa_low = 0; /* true zero */ 487 #endif 488 } 489 else 490 { 491 /* 492 * find number of significant integer bits 493 */ 494 mask = 0x80000000; 495 if (outlfp.l_ui) 496 { 497 msb = 63; 498 while (mask && ((outlfp.l_ui & mask) == 0)) 499 { 500 mask >>= 1; 501 msb--; 502 } 503 } 504 else 505 { 506 msb = 31; 507 while (mask && ((outlfp.l_uf & mask) == 0)) 508 { 509 mask >>= 1; 510 msb--; 511 } 512 } 513 514 switch (size) 515 { 516 case IEEE_SINGLE: 517 mantissa_high = 0; 518 if (msb >= 32) 519 { 520 mantissa_low = (outlfp.l_ui & ((1 << (msb - 32)) - 1)) << (mbits - (msb - 32)); 521 mantissa_low |= outlfp.l_uf >> (mbits - (msb - 32)); 522 } 523 else 524 { 525 mantissa_low = (outlfp.l_uf << (mbits - msb)) & ((1 << mbits) - 1); 526 } 527 break; 528 529 case IEEE_DOUBLE: 530 if (msb >= 32) 531 { 532 mantissa_high = (outlfp.l_ui << (mbits - msb)) & ((1 << (mbits - 32)) - 1); 533 mantissa_high |= outlfp.l_uf >> (32 - (mbits - msb)); 534 mantissa_low = (outlfp.l_ui & ((1 << (msb - mbits)) - 1)) << (32 - (msb - mbits)); 535 mantissa_low |= outlfp.l_uf >> (msb - mbits); 536 } 537 else 538 { 539 mantissa_high = outlfp.l_uf << (mbits - 32 - msb); 540 mantissa_low = outlfp.l_uf << (mbits - 32); 541 } 542 } 543 544 #ifdef LIBDEBUG 545 exponent = msb - 32; 546 characteristic = exponent + bias; 547 548 if (debug > 4) 549 printf("FP: %s\n", fmt_flt(sign, mantissa_high, mantissa_low, characteristic)); 550 #endif 551 } 552 return IEEE_OK; 553 } 554 #endif /* PUT_IEEE754_UNUSED_FUNC */ 555 556 557 #if defined(DEBUG) && defined(LIBDEBUG) 558 int main( 559 int argc, 560 char **argv 561 ) 562 { 563 static offsets_t native_off = { 0, 1, 2, 3, 4, 5, 6, 7 }; 564 double f = 1.0; 565 double *f_p = &f; 566 l_fp fp; 567 568 if (argc == 2) 569 { 570 if (sscanf(argv[1], "%lf", &f) != 1) 571 { 572 printf("cannot convert %s to a float\n", argv[1]); 573 return 1; 574 } 575 } 576 577 printf("double: %s %s\n", fmt_blong(*(unsigned long *)&f, 32), fmt_blong(*(unsigned long *)((char *)(&f)+4), 32)); 578 printf("fetch from %f = %d\n", f, fetch_ieee754((void *)&f_p, IEEE_DOUBLE, &fp, native_off)); 579 printf("fp [%s %s] = %s\n", fmt_blong(fp.l_ui, 32), fmt_blong(fp.l_uf, 32), mfptoa(fp.l_ui, fp.l_uf, 15)); 580 f_p = &f; 581 #ifdef PUT_IEEE754_UNUSED_FUNC 582 put_ieee754((void *)&f_p, IEEE_DOUBLE, &fp, native_off); 583 /* there should be a check on *f_p (f) having the expected result here */ 584 #endif /* PUT_IEEE754_UNUSED_FUNC */ 585 586 return 0; 587 } 588 589 #endif 590 /* 591 * History: 592 * 593 * ieee754io.c,v 594 * Revision 4.12 2005/04/16 17:32:10 kardel 595 * update copyright 596 * 597 * Revision 4.11 2004/11/14 15:29:41 kardel 598 * support PPSAPI, upgrade Copyright to Berkeley style 599 * 600 * Revision 4.8 1999/02/21 12:17:36 kardel 601 * 4.91f reconcilation 602 * 603 * Revision 4.7 1999/02/21 11:26:03 kardel 604 * renamed index to fieldindex to avoid index() name clash 605 * 606 * Revision 4.6 1998/11/15 20:27:52 kardel 607 * Release 4.0.73e13 reconcilation 608 * 609 * Revision 4.5 1998/08/16 19:01:51 kardel 610 * debug information only compile for LIBDEBUG case 611 * 612 * Revision 4.4 1998/08/09 09:39:28 kardel 613 * Release 4.0.73e2 reconcilation 614 * 615 * Revision 4.3 1998/06/13 11:56:19 kardel 616 * disabled putbute() for the time being 617 * 618 * Revision 4.2 1998/06/12 15:16:58 kardel 619 * ansi2knr compatibility 620 * 621 * Revision 4.1 1998/05/24 07:59:56 kardel 622 * conditional debug support 623 * 624 * Revision 4.0 1998/04/10 19:46:29 kardel 625 * Start 4.0 release version numbering 626 * 627 * Revision 1.1 1998/04/10 19:27:46 kardel 628 * initial NTP VERSION 4 integration of PARSE with GPS166 binary support 629 * 630 * Revision 1.1 1997/10/06 21:05:45 kardel 631 * new parse structure 632 * 633 */ 634