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