ieee754io.c revision 1.5
1 /* $NetBSD: ieee754io.c,v 1.5 2020/05/25 20:47:25 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 #include "lib_strbuf.h" 56 57 static char * 58 fmt_blong( 59 unsigned long val, 60 int cnt 61 ) 62 { 63 char *buf, *s; 64 int i = cnt; 65 66 val <<= 32 - cnt; 67 LIB_GETBUF(buf); 68 s = buf; 69 70 while (i--) 71 { 72 if (val & 0x80000000) 73 { 74 *s++ = '1'; 75 } 76 else 77 { 78 *s++ = '0'; 79 } 80 val <<= 1; 81 } 82 *s = '\0'; 83 return buf; 84 } 85 86 static char * 87 fmt_flt( 88 unsigned int sign, 89 unsigned long mh, 90 unsigned long ml, 91 unsigned long ch 92 ) 93 { 94 char *buf; 95 96 LIB_GETBUF(buf); 97 snprintf(buf, LIB_BUFLENGTH, "%c %s %s %s", sign ? '-' : '+', 98 fmt_blong(ch, 11), 99 fmt_blong(mh, 20), 100 fmt_blong(ml, 32)); 101 102 return buf; 103 } 104 105 static char * 106 fmt_hex( 107 unsigned char *bufp, 108 int length 109 ) 110 { 111 char * buf; 112 char hex[4]; 113 int i; 114 115 LIB_GETBUF(buf); 116 buf[0] = '\0'; 117 for (i = 0; i < length; i++) { 118 snprintf(hex, sizeof(hex), "%02x", bufp[i]); 119 strlcat(buf, hex, LIB_BUFLENGTH); 120 } 121 122 return buf; 123 } 124 125 #endif 126 127 static unsigned char 128 get_byte( 129 unsigned char *bufp, 130 offsets_t offset, 131 int *fieldindex 132 ) 133 { 134 unsigned char val; 135 136 val = *(bufp + offset[*fieldindex]); 137 #ifdef LIBDEBUG 138 if (debug > 4) 139 printf("fetchieee754: getbyte(0x%08x, %d) = 0x%02x\n", (unsigned int)(bufp)+offset[*fieldindex], *fieldindex, val); 140 #endif 141 (*fieldindex)++; 142 return val; 143 } 144 145 #ifdef __not_yet__ 146 static void 147 put_byte( 148 unsigned char *bufp, 149 offsets_t offsets, 150 int *fieldindex, 151 unsigned char val 152 ) 153 { 154 *(bufp + offsets[*fieldindex]) = val; 155 (*fieldindex)++; 156 } 157 #endif 158 159 /* 160 * make conversions to and from external IEEE754 formats and internal 161 * NTP FP format. 162 */ 163 int 164 fetch_ieee754( 165 unsigned char **buffpp, 166 int size, 167 l_fp *lfpp, 168 offsets_t offsets 169 ) 170 { 171 unsigned char *bufp = *buffpp; 172 unsigned int sign; 173 unsigned int bias; 174 unsigned int maxexp; 175 int mbits; 176 u_long mantissa_low; 177 u_long mantissa_high; 178 u_long characteristic; 179 long exponent; 180 #ifdef LIBDEBUG 181 int length; 182 #endif 183 unsigned char val; 184 int fieldindex = 0; 185 186 switch (size) 187 { 188 case IEEE_DOUBLE: 189 #ifdef LIBDEBUG 190 length = 8; 191 #endif 192 mbits = 52; 193 bias = 1023; 194 maxexp = 2047; 195 break; 196 197 case IEEE_SINGLE: 198 #ifdef LIBDEBUG 199 length = 4; 200 #endif 201 mbits = 23; 202 bias = 127; 203 maxexp = 255; 204 break; 205 206 default: 207 return IEEE_BADCALL; 208 } 209 210 val = get_byte(bufp, offsets, &fieldindex); /* fetch sign byte & first part of characteristic */ 211 212 sign = (val & 0x80) != 0; 213 characteristic = (val & 0x7F); 214 215 val = get_byte(bufp, offsets, &fieldindex); /* fetch rest of characteristic and start of mantissa */ 216 217 switch (size) 218 { 219 case IEEE_SINGLE: 220 characteristic <<= 1; 221 characteristic |= (val & 0x80) != 0; /* grab last characteristic bit */ 222 223 mantissa_high = 0; 224 225 mantissa_low = (val &0x7F) << 16; 226 mantissa_low |= (u_long)get_byte(bufp, offsets, &fieldindex) << 8; 227 mantissa_low |= get_byte(bufp, offsets, &fieldindex); 228 break; 229 230 case IEEE_DOUBLE: 231 characteristic <<= 4; 232 characteristic |= (val & 0xF0) >> 4; /* grab lower characteristic bits */ 233 234 mantissa_high = (val & 0x0F) << 16; 235 mantissa_high |= (u_long)get_byte(bufp, offsets, &fieldindex) << 8; 236 mantissa_high |= get_byte(bufp, offsets, &fieldindex); 237 238 mantissa_low = (u_long)get_byte(bufp, offsets, &fieldindex) << 24; 239 mantissa_low |= (u_long)get_byte(bufp, offsets, &fieldindex) << 16; 240 mantissa_low |= (u_long)get_byte(bufp, offsets, &fieldindex) << 8; 241 mantissa_low |= get_byte(bufp, offsets, &fieldindex); 242 break; 243 244 default: 245 return IEEE_BADCALL; 246 } 247 #ifdef LIBDEBUG 248 if (debug > 4) 249 { 250 double d; 251 float f; 252 253 if (size == IEEE_SINGLE) 254 { 255 int i; 256 257 for (i = 0; i < length; i++) 258 { 259 *((unsigned char *)(&f)+i) = *(*buffpp + offsets[i]); 260 } 261 d = f; 262 } 263 else 264 { 265 int i; 266 267 for (i = 0; i < length; i++) 268 { 269 *((unsigned char *)(&d)+i) = *(*buffpp + offsets[i]); 270 } 271 } 272 273 printf("fetchieee754: FP: %s -> %s -> %e(=%s)\n", fmt_hex(*buffpp, length), 274 fmt_flt(sign, mantissa_high, mantissa_low, characteristic), 275 d, fmt_hex((unsigned char *)&d, length)); 276 } 277 #endif 278 279 *buffpp += fieldindex; 280 281 /* 282 * detect funny numbers 283 */ 284 if (characteristic == maxexp) 285 { 286 /* 287 * NaN or Infinity 288 */ 289 if (mantissa_low || mantissa_high) 290 { 291 /* 292 * NaN 293 */ 294 return IEEE_NAN; 295 } 296 else 297 { 298 /* 299 * +Inf or -Inf 300 */ 301 return sign ? IEEE_NEGINFINITY : IEEE_POSINFINITY; 302 } 303 } 304 else 305 { 306 /* 307 * collect real numbers 308 */ 309 310 L_CLR(lfpp); 311 312 /* 313 * check for overflows 314 */ 315 exponent = characteristic - bias; 316 317 if (exponent > 31) /* sorry - hardcoded */ 318 { 319 /* 320 * overflow only in respect to NTP-FP representation 321 */ 322 return sign ? IEEE_NEGOVERFLOW : IEEE_POSOVERFLOW; 323 } 324 else 325 { 326 int frac_offset; /* where the fraction starts */ 327 328 frac_offset = mbits - exponent; 329 330 if (characteristic == 0) 331 { 332 /* 333 * de-normalized or tiny number - fits only as 0 334 */ 335 return IEEE_OK; 336 } 337 else 338 { 339 /* 340 * adjust for implied 1 341 */ 342 if (mbits > 31) 343 mantissa_high |= 1 << (mbits - 32); 344 else 345 mantissa_low |= 1 << mbits; 346 347 /* 348 * take mantissa apart - if only all machine would support 349 * 64 bit operations 8-( 350 */ 351 if (frac_offset > mbits) 352 { 353 lfpp->l_ui = 0; /* only fractional number */ 354 frac_offset -= mbits + 1; /* will now contain right shift count - 1*/ 355 if (mbits > 31) 356 { 357 lfpp->l_uf = mantissa_high << (63 - mbits); 358 lfpp->l_uf |= mantissa_low >> (mbits - 33); 359 lfpp->l_uf >>= frac_offset; 360 } 361 else 362 { 363 lfpp->l_uf = mantissa_low >> frac_offset; 364 } 365 } 366 else 367 { 368 if (frac_offset > 32) 369 { 370 /* 371 * must split in high word 372 */ 373 lfpp->l_ui = mantissa_high >> (frac_offset - 32); 374 lfpp->l_uf = (mantissa_high & ((1 << (frac_offset - 32)) - 1)) << (64 - frac_offset); 375 lfpp->l_uf |= mantissa_low >> (frac_offset - 32); 376 } 377 else 378 { 379 /* 380 * must split in low word 381 */ 382 lfpp->l_ui = mantissa_high << (32 - frac_offset); 383 lfpp->l_ui |= (mantissa_low >> frac_offset) & ((1 << (32 - frac_offset)) - 1); 384 lfpp->l_uf = (mantissa_low & ((1 << frac_offset) - 1)) << (32 - frac_offset); 385 } 386 } 387 388 /* 389 * adjust for sign 390 */ 391 if (sign) 392 { 393 L_NEG(lfpp); 394 } 395 396 return IEEE_OK; 397 } 398 } 399 } 400 } 401 402 int 403 put_ieee754( 404 unsigned char **bufpp, 405 int size, 406 l_fp *lfpp, 407 offsets_t offsets 408 ) 409 { 410 l_fp outlfp; 411 #ifdef LIBDEBUG 412 unsigned int sign; 413 unsigned int bias; 414 #endif 415 /*unsigned int maxexp;*/ 416 int mbits; 417 int msb; 418 u_long mantissa_low = 0; 419 u_long mantissa_high = 0; 420 #ifdef LIBDEBUG 421 u_long characteristic = 0; 422 long exponent; 423 #endif 424 /*int length;*/ 425 unsigned long mask; 426 427 outlfp = *lfpp; 428 429 switch (size) 430 { 431 case IEEE_DOUBLE: 432 /*length = 8;*/ 433 mbits = 52; 434 #ifdef LIBDEBUG 435 bias = 1023; 436 #endif 437 /*maxexp = 2047;*/ 438 break; 439 440 case IEEE_SINGLE: 441 /*length = 4;*/ 442 mbits = 23; 443 #ifdef LIBDEBUG 444 bias = 127; 445 #endif 446 /*maxexp = 255;*/ 447 break; 448 449 default: 450 return IEEE_BADCALL; 451 } 452 453 /* 454 * find sign 455 */ 456 if (L_ISNEG(&outlfp)) 457 { 458 L_NEG(&outlfp); 459 #ifdef LIBDEBUG 460 sign = 1; 461 #endif 462 } 463 else 464 { 465 #ifdef LIBDEBUG 466 sign = 0; 467 #endif 468 } 469 470 if (L_ISZERO(&outlfp)) 471 { 472 #ifdef LIBDEBUG 473 exponent = mantissa_high = mantissa_low = 0; /* true zero */ 474 #endif 475 } 476 else 477 { 478 /* 479 * find number of significant integer bits 480 */ 481 mask = 0x80000000; 482 if (outlfp.l_ui) 483 { 484 msb = 63; 485 while (mask && ((outlfp.l_ui & mask) == 0)) 486 { 487 mask >>= 1; 488 msb--; 489 } 490 } 491 else 492 { 493 msb = 31; 494 while (mask && ((outlfp.l_uf & mask) == 0)) 495 { 496 mask >>= 1; 497 msb--; 498 } 499 } 500 501 switch (size) 502 { 503 case IEEE_SINGLE: 504 mantissa_high = 0; 505 if (msb >= 32) 506 { 507 mantissa_low = (outlfp.l_ui & ((1 << (msb - 32)) - 1)) << (mbits - (msb - 32)); 508 mantissa_low |= outlfp.l_uf >> (mbits - (msb - 32)); 509 } 510 else 511 { 512 mantissa_low = (outlfp.l_uf << (mbits - msb)) & ((1 << mbits) - 1); 513 } 514 break; 515 516 case IEEE_DOUBLE: 517 if (msb >= 32) 518 { 519 mantissa_high = (outlfp.l_ui << (mbits - msb)) & ((1 << (mbits - 32)) - 1); 520 mantissa_high |= outlfp.l_uf >> (32 - (mbits - msb)); 521 mantissa_low = (outlfp.l_ui & ((1 << (msb - mbits)) - 1)) << (32 - (msb - mbits)); 522 mantissa_low |= outlfp.l_uf >> (msb - mbits); 523 } 524 else 525 { 526 mantissa_high = outlfp.l_uf << (mbits - 32 - msb); 527 mantissa_low = outlfp.l_uf << (mbits - 32); 528 } 529 } 530 531 #ifdef LIBDEBUG 532 exponent = msb - 32; 533 characteristic = exponent + bias; 534 535 if (debug > 4) 536 printf("FP: %s\n", fmt_flt(sign, mantissa_high, mantissa_low, characteristic)); 537 #endif 538 } 539 return IEEE_OK; 540 } 541 542 543 #if defined(DEBUG) && defined(LIBDEBUG) 544 int main( 545 int argc, 546 char **argv 547 ) 548 { 549 static offsets_t native_off = { 0, 1, 2, 3, 4, 5, 6, 7 }; 550 double f = 1.0; 551 double *f_p = &f; 552 l_fp fp; 553 554 if (argc == 2) 555 { 556 if (sscanf(argv[1], "%lf", &f) != 1) 557 { 558 printf("cannot convert %s to a float\n", argv[1]); 559 return 1; 560 } 561 } 562 563 printf("double: %s %s\n", fmt_blong(*(unsigned long *)&f, 32), fmt_blong(*(unsigned long *)((char *)(&f)+4), 32)); 564 printf("fetch from %f = %d\n", f, fetch_ieee754((void *)&f_p, IEEE_DOUBLE, &fp, native_off)); 565 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)); 566 f_p = &f; 567 put_ieee754((void *)&f_p, IEEE_DOUBLE, &fp, native_off); 568 569 return 0; 570 } 571 572 #endif 573 /* 574 * History: 575 * 576 * ieee754io.c,v 577 * Revision 4.12 2005/04/16 17:32:10 kardel 578 * update copyright 579 * 580 * Revision 4.11 2004/11/14 15:29:41 kardel 581 * support PPSAPI, upgrade Copyright to Berkeley style 582 * 583 * Revision 4.8 1999/02/21 12:17:36 kardel 584 * 4.91f reconcilation 585 * 586 * Revision 4.7 1999/02/21 11:26:03 kardel 587 * renamed index to fieldindex to avoid index() name clash 588 * 589 * Revision 4.6 1998/11/15 20:27:52 kardel 590 * Release 4.0.73e13 reconcilation 591 * 592 * Revision 4.5 1998/08/16 19:01:51 kardel 593 * debug information only compile for LIBDEBUG case 594 * 595 * Revision 4.4 1998/08/09 09:39:28 kardel 596 * Release 4.0.73e2 reconcilation 597 * 598 * Revision 4.3 1998/06/13 11:56:19 kardel 599 * disabled putbute() for the time being 600 * 601 * Revision 4.2 1998/06/12 15:16:58 kardel 602 * ansi2knr compatibility 603 * 604 * Revision 4.1 1998/05/24 07:59:56 kardel 605 * conditional debug support 606 * 607 * Revision 4.0 1998/04/10 19:46:29 kardel 608 * Start 4.0 release version numbering 609 * 610 * Revision 1.1 1998/04/10 19:27:46 kardel 611 * initial NTP VERSION 4 integration of PARSE with GPS166 binary support 612 * 613 * Revision 1.1 1997/10/06 21:05:45 kardel 614 * new parse structure 615 * 616 */ 617
Indexes created Mon Apr 13 00:22:24 UTC 2026