home *** CD-ROM | disk | FTP | other *** search
/ AmigActive 13 / AACD13.ISO / AACD / Sound / LAME / src / ieeefloat.c < prev    next >
C/C++ Source or Header  |  2000-06-07  |  29KB  |  964 lines
  1. /* Copyright (C) 1988-1991 Apple Computer, Inc.
  2.  * All Rights Reserved.
  3.  *
  4.  * Warranty Information
  5.  * Even though Apple has reviewed this software, Apple makes no warranty
  6.  * or representation, either express or implied, with respect to this
  7.  * software, its quality, accuracy, merchantability, or fitness for a
  8.  * particular purpose.  As a result, this software is provided "as is,"
  9.  * and you, its user, are assuming the entire risk as to its quality
  10.  * and accuracy.
  11.  *
  12.  * This code may be used and freely distributed as long as it includes
  13.  * this copyright notice and the warranty information.
  14.  *
  15.  * Machine-independent I/O routines for IEEE floating-point numbers.
  16.  *
  17.  * NaN's and infinities are converted to HUGE_VAL or HUGE, which
  18.  * happens to be infinity on IEEE machines.  Unfortunately, it is
  19.  * impossible to preserve NaN's in a machine-independent way.
  20.  * Infinities are, however, preserved on IEEE machines.
  21.  *
  22.  * These routines have been tested on the following machines:
  23.  *    Apple Macintosh, MPW 3.1 C compiler
  24.  *    Apple Macintosh, THINK C compiler
  25.  *    Silicon Graphics IRIS, MIPS compiler
  26.  *    Cray X/MP and Y/MP
  27.  *    Digital Equipment VAX
  28.  *    Sequent Balance (Multiprocesor 386)
  29.  *    NeXT
  30.  *
  31.  *
  32.  * Implemented by Malcolm Slaney and Ken Turkowski.
  33.  *
  34.  * Malcolm Slaney contributions during 1988-1990 include big- and little-
  35.  * endian file I/O, conversion to and from Motorola's extended 80-bit
  36.  * floating-point format, and conversions to and from IEEE single-
  37.  * precision floating-point format.
  38.  *
  39.  * In 1991, Ken Turkowski implemented the conversions to and from
  40.  * IEEE double-precision format, added more precision to the extended
  41.  * conversions, and accommodated conversions involving +/- infinity,
  42.  * NaN's, and denormalized numbers.
  43.  *
  44.  * $Id: ieeefloat.c,v 1.4 2000/06/07 22:56:02 sbellon Exp $
  45.  *
  46.  * $Log: ieeefloat.c,v $
  47.  * Revision 1.4  2000/06/07 22:56:02  sbellon
  48.  * added support for FPA10 hardware (RISC OS only)
  49.  *
  50.  * Revision 1.3  2000/02/21 23:05:05  markt
  51.  * some 64bit DEC Alpha patches
  52.  *
  53.  * Revision 1.2  2000/02/19 13:32:30  afaber
  54.  * Fixed many warning messages when compiling with MSVC
  55.  *
  56.  * Revision 1.1.1.1  1999/11/24 08:42:58  markt
  57.  * initial checkin of LAME
  58.  * Starting with LAME 3.57beta with some modifications
  59.  *
  60.  * Revision 1.1  1993/06/11  17:45:46  malcolm
  61.  * Initial revision
  62.  *
  63.  */
  64.  
  65. #include        <limits.h>
  66. #include    <stdio.h>
  67. #if defined(__riscos__) && defined(FPA10)
  68. #include    "ymath.h"
  69. #else
  70. #include    <math.h>
  71. #endif
  72. #include    "ieeefloat.h"
  73.  
  74.  
  75. /****************************************************************
  76.  * The following two routines make up for deficiencies in many
  77.  * compilers to convert properly between unsigned integers and
  78.  * floating-point.  Some compilers which have this bug are the
  79.  * THINK_C compiler for the Macintosh and the C compiler for the
  80.  * Silicon Graphics MIPS-based Iris.
  81.  ****************************************************************/
  82.  
  83. #ifdef applec    /* The Apple C compiler works */
  84. # define FloatToUnsigned(f)    ((unsigned long)(f))
  85. # define UnsignedToFloat(u)    ((defdouble)(u))
  86. #else /* applec */
  87. # define FloatToUnsigned(f)    ((unsigned long)(((long)((f) - 2147483648.0)) + 2147483647L + 1))
  88. # define UnsignedToFloat(u)    (((defdouble)((long)((u) - 2147483647L - 1))) + 2147483648.0)
  89. #endif /* applec */
  90.  
  91.  
  92. /****************************************************************
  93.  * Single precision IEEE floating-point conversion routines
  94.  ****************************************************************/
  95.  
  96. #define SEXP_MAX        255
  97. #define SEXP_OFFSET        127
  98. #define SEXP_SIZE        8
  99. #define SEXP_POSITION    (32-SEXP_SIZE-1)
  100.  
  101.  
  102. defdouble
  103. ConvertFromIeeeSingle(char* bytes)
  104. {
  105.     defdouble    f;
  106.     long    mantissa, expon;
  107.     long    bits;
  108.  
  109.     bits =    ((unsigned long)(bytes[0] & 0xFF) << 24)
  110.         |    ((unsigned long)(bytes[1] & 0xFF) << 16)
  111.         |    ((unsigned long)(bytes[2] & 0xFF) << 8)
  112.         |     (unsigned long)(bytes[3] & 0xFF);        /* Assemble bytes into a long */
  113.  
  114.     if ((bits & 0x7FFFFFFF) == 0) {
  115.         f = 0;
  116.     }
  117.  
  118.     else {
  119.         expon = (bits & 0x7F800000) >> SEXP_POSITION;
  120.         if (expon == SEXP_MAX) {        /* Infinity or NaN */
  121.             f = HUGE_VAL;        /* Map NaN's to infinity */
  122.         }
  123.         else {
  124.             if (expon == 0) {    /* Denormalized number */
  125.                 mantissa = (bits & 0x7fffff);
  126.                 f = ldexp((defdouble) mantissa, (int) (expon - SEXP_OFFSET - SEXP_POSITION + 1));
  127.             }
  128.             else {                /* Normalized number */
  129.                 mantissa = (bits & 0x7fffff) + 0x800000;    /* Insert hidden bit */
  130.                 f = ldexp((defdouble) mantissa, (int) (expon - SEXP_OFFSET - SEXP_POSITION));
  131.             }
  132.         }
  133.     }
  134.  
  135.     if (bits & LONG_MIN)
  136.         return -f;
  137.     else
  138.         return f;
  139. }
  140.  
  141.  
  142. /****************************************************************/
  143.  
  144.  
  145. void
  146. ConvertToIeeeSingle(defdouble num, char* bytes)
  147. {
  148.     long    sign;
  149.     register long bits;
  150.  
  151.     if (num < 0) {    /* Can't distinguish a negative zero */
  152.         sign = LONG_MIN;
  153.         num *= -1;
  154.     } else {
  155.         sign = 0;
  156.     }
  157.  
  158.     if (num == 0) {
  159.         bits = 0;
  160.     }
  161.  
  162.     else {
  163.         defdouble fMant;
  164.         int expon;
  165.  
  166.         fMant = frexp(num, &expon);
  167.  
  168.         if ((expon > (SEXP_MAX-SEXP_OFFSET+1)) || !(fMant < 1)) {
  169.             /* NaN's and infinities fail second test */
  170.             bits = sign | 0x7F800000;        /* +/- infinity */
  171.         }
  172.  
  173.         else {
  174.             long mantissa;
  175.  
  176.             if (expon < -(SEXP_OFFSET-2)) {    /* Smaller than normalized */
  177.                 int shift = (SEXP_POSITION+1) + (SEXP_OFFSET-2) + expon;
  178.                 if (shift < 0) {    /* Way too small: flush to zero */
  179.                     bits = sign;
  180.                 }
  181.                 else {            /* Nonzero denormalized number */
  182.                     mantissa = (long)(fMant * (1L << shift));
  183.                     bits = sign | mantissa;
  184.                 }
  185.             }
  186.  
  187.             else {                /* Normalized number */
  188.                 mantissa = (long)floor(fMant * (1L << (SEXP_POSITION+1)));
  189.                 mantissa -= (1L << SEXP_POSITION);            /* Hide MSB */
  190.                 bits = sign | ((long)((expon + SEXP_OFFSET - 1)) << SEXP_POSITION) | mantissa;
  191.             }
  192.         }
  193.     }
  194.  
  195.     bytes[0] = (char)(bits >> 24);    /* Copy to byte string */
  196.     bytes[1] = (char)(bits >> 16);
  197.     bytes[2] = (char)(bits >> 8);
  198.     bytes[3] = (char)(bits);
  199. }
  200.  
  201.  
  202. /****************************************************************
  203.  * Double precision IEEE floating-point conversion routines
  204.  ****************************************************************/
  205.  
  206. #define DEXP_MAX        2047
  207. #define DEXP_OFFSET        1023
  208. #define DEXP_SIZE        11
  209. #define DEXP_POSITION    (32-DEXP_SIZE-1)
  210.  
  211.  
  212. defdouble
  213. ConvertFromIeeeDouble(char* bytes)
  214. {
  215.     defdouble    f;
  216.     long    mantissa, expon;
  217.     unsigned long first, second;
  218.  
  219.     first = ((unsigned long)(bytes[0] & 0xFF) << 24)
  220.         |    ((unsigned long)(bytes[1] & 0xFF) << 16)
  221.         |    ((unsigned long)(bytes[2] & 0xFF) << 8)
  222.         |     (unsigned long)(bytes[3] & 0xFF);
  223.     second= ((unsigned long)(bytes[4] & 0xFF) << 24)
  224.         |    ((unsigned long)(bytes[5] & 0xFF) << 16)
  225.         |    ((unsigned long)(bytes[6] & 0xFF) << 8)
  226.         |     (unsigned long)(bytes[7] & 0xFF);
  227.  
  228.     if (first == 0 && second == 0) {
  229.         f = 0;
  230.     }
  231.  
  232.     else {
  233.         expon = (first & 0x7FF00000) >> DEXP_POSITION;
  234.         if (expon == DEXP_MAX) {        /* Infinity or NaN */
  235.             f = HUGE_VAL;        /* Map NaN's to infinity */
  236.         }
  237.         else {
  238.             if (expon == 0) {    /* Denormalized number */
  239.                 mantissa = (first & 0x000FFFFF);
  240.                 f = ldexp((defdouble) mantissa, (int) (expon - DEXP_OFFSET - DEXP_POSITION + 1));
  241.                 f += ldexp(UnsignedToFloat(second), (int) (expon - DEXP_OFFSET - DEXP_POSITION + 1 - 32));
  242.             }
  243.             else {                /* Normalized number */
  244.                 mantissa = (first & 0x000FFFFF) + 0x00100000;    /* Insert hidden bit */
  245.                 f = ldexp((defdouble) mantissa, (int) (expon - DEXP_OFFSET - DEXP_POSITION));
  246.                 f += ldexp(UnsignedToFloat(second), (int) (expon - DEXP_OFFSET - DEXP_POSITION - 32));
  247.             }
  248.         }
  249.     }
  250.  
  251.     if (first & 0x80000000)
  252.         return -f;
  253.     else
  254.         return f;
  255. }
  256.  
  257.  
  258. /****************************************************************/
  259.  
  260.  
  261. void
  262. ConvertToIeeeDouble(defdouble num, char *bytes)
  263. {
  264.     long    sign;
  265.     long    first, second;
  266.  
  267.     if (num < 0) {    /* Can't distinguish a negative zero */
  268.         sign = LONG_MIN;
  269.         num *= -1;
  270.     } else {
  271.         sign = 0;
  272.     }
  273.  
  274.     if (num == 0) {
  275.         first = 0;
  276.         second = 0;
  277.     }
  278.  
  279.     else {
  280.         defdouble fMant, fsMant;
  281.         int expon;
  282.  
  283.         fMant = frexp(num, &expon);
  284.  
  285.         if ((expon > (DEXP_MAX-DEXP_OFFSET+1)) || !(fMant < 1)) {
  286.             /* NaN's and infinities fail second test */
  287.             first = sign | 0x7FF00000;        /* +/- infinity */
  288.             second = 0;
  289.         }
  290.  
  291.         else {
  292.             long mantissa;
  293.  
  294.             if (expon < -(DEXP_OFFSET-2)) {    /* Smaller than normalized */
  295.                 int shift = (DEXP_POSITION+1) + (DEXP_OFFSET-2) + expon;
  296.                 if (shift < 0) {    /* Too small for something in the MS word */
  297.                     first = sign;
  298.                     shift += 32;
  299.                     if (shift < 0) {    /* Way too small: flush to zero */
  300.                         second = 0;
  301.                     }
  302.                     else {            /* Pretty small demorn */
  303.                         second = FloatToUnsigned(floor(ldexp(fMant, shift)));
  304.                     }
  305.                 }
  306.                 else {            /* Nonzero denormalized number */
  307.                     fsMant = ldexp(fMant, shift);
  308.                     mantissa = (long)floor(fsMant);
  309.                     first = sign | mantissa;
  310.                     second = FloatToUnsigned(floor(ldexp(fsMant - mantissa, 32)));
  311.                 }
  312.             }
  313.  
  314.             else {                /* Normalized number */
  315.                 fsMant = ldexp(fMant, DEXP_POSITION+1);
  316.                 mantissa = (long)floor(fsMant);
  317.                 mantissa -= (1L << DEXP_POSITION);            /* Hide MSB */
  318.                 fsMant -= (1L << DEXP_POSITION);
  319.                 first = sign | ((long)((expon + DEXP_OFFSET - 1)) << DEXP_POSITION) | mantissa;
  320.                 second = FloatToUnsigned(floor(ldexp(fsMant - mantissa, 32)));
  321.             }
  322.         }
  323.     }
  324.  
  325.     bytes[0] = (char)(first >> 24);
  326.     bytes[1] = (char)(first >> 16);
  327.     bytes[2] = (char)(first >> 8);
  328.     bytes[3] = (char)(first);
  329.     bytes[4] = (char)(second >> 24);
  330.     bytes[5] = (char)(second >> 16);
  331.     bytes[6] = (char)(second >> 8);
  332.     bytes[7] = (char)(second);
  333. }
  334.  
  335.  
  336. /****************************************************************
  337.  * Extended precision IEEE floating-point conversion routines
  338.  ****************************************************************/
  339.  
  340. defdouble
  341. ConvertFromIeeeExtended(char* bytes)
  342. {
  343.     defdouble    f;
  344.     long    expon;
  345.     unsigned long hiMant, loMant;
  346.  
  347. #ifdef    TEST
  348. printf("ConvertFromIEEEExtended(%lx,%lx,%lx,%lx,%lx,%lx,%lx,%lx,%lx,%lx\r",
  349.     (long)bytes[0], (long)bytes[1], (long)bytes[2], (long)bytes[3],
  350.     (long)bytes[4], (long)bytes[5], (long)bytes[6],
  351.     (long)bytes[7], (long)bytes[8], (long)bytes[9]);
  352. #endif
  353.  
  354.     expon = ((bytes[0] & 0x7F) << 8) | (bytes[1] & 0xFF);
  355.     hiMant    =    ((unsigned long)(bytes[2] & 0xFF) << 24)
  356.             |    ((unsigned long)(bytes[3] & 0xFF) << 16)
  357.             |    ((unsigned long)(bytes[4] & 0xFF) << 8)
  358.             |    ((unsigned long)(bytes[5] & 0xFF));
  359.     loMant    =    ((unsigned long)(bytes[6] & 0xFF) << 24)
  360.             |    ((unsigned long)(bytes[7] & 0xFF) << 16)
  361.             |    ((unsigned long)(bytes[8] & 0xFF) << 8)
  362.             |    ((unsigned long)(bytes[9] & 0xFF));
  363.  
  364.     if (expon == 0 && hiMant == 0 && loMant == 0) {
  365.         f = 0;
  366.     }
  367.     else {
  368.         if (expon == 0x7FFF) {    /* Infinity or NaN */
  369.             f = HUGE_VAL;
  370.         }
  371.         else {
  372.             expon -= 16383;
  373.             f  = ldexp(UnsignedToFloat(hiMant), (int) (expon -= 31));
  374.             f += ldexp(UnsignedToFloat(loMant), (int) (expon -= 32));
  375.         }
  376.     }
  377.  
  378.     if (bytes[0] & 0x80)
  379.         return -f;
  380.     else
  381.         return f;
  382. }
  383.  
  384.  
  385. /****************************************************************/
  386.  
  387.  
  388. void
  389. ConvertToIeeeExtended(defdouble num, char *bytes)
  390. {
  391.     int    sign;
  392.     int expon;
  393.     defdouble fMant, fsMant;
  394.     unsigned long hiMant, loMant;
  395.  
  396.     if (num < 0) {
  397.         sign = 0x8000;
  398.         num *= -1;
  399.     } else {
  400.         sign = 0;
  401.     }
  402.  
  403.     if (num == 0) {
  404.         expon = 0; hiMant = 0; loMant = 0;
  405.     }
  406.     else {
  407.         fMant = frexp(num, &expon);
  408.         if ((expon > 16384) || !(fMant < 1)) {    /* Infinity or NaN */
  409.             expon = sign|0x7FFF; hiMant = 0; loMant = 0; /* infinity */
  410.         }
  411.         else {    /* Finite */
  412.             expon += 16382;
  413.             if (expon < 0) {    /* denormalized */
  414.                 fMant = ldexp(fMant, expon);
  415.                 expon = 0;
  416.             }
  417.             expon |= sign;
  418.             fMant = ldexp(fMant, 32);          fsMant = floor(fMant); hiMant = FloatToUnsigned(fsMant);
  419.             fMant = ldexp(fMant - fsMant, 32); fsMant = floor(fMant); loMant = FloatToUnsigned(fsMant);
  420.         }
  421.     }
  422.  
  423.     bytes[0] = expon >> 8;
  424.     bytes[1] = expon;
  425.     bytes[2] = (char)(hiMant >> 24);
  426.     bytes[3] = (char)(hiMant >> 16);
  427.     bytes[4] = (char)(hiMant >> 8);
  428.     bytes[5] = (char)(hiMant);
  429.     bytes[6] = (char)(loMant >> 24);
  430.     bytes[7] = (char)(loMant >> 16);
  431.     bytes[8] = (char)(loMant >> 8);
  432.     bytes[9] = (char)(loMant);
  433. }
  434.  
  435. /****************************************************************
  436.  * Testing routines for the floating-point conversions.
  437.  ****************************************************************/
  438.  
  439. #ifdef METROWERKS
  440. #define IEEE
  441. #endif
  442. #ifdef applec
  443. # define IEEE
  444. #endif /* applec */
  445. #ifdef THINK_C
  446. # define IEEE
  447. #endif /* THINK_C */
  448. #ifdef sgi
  449. # define IEEE
  450. #endif /* sgi */
  451. #ifdef sequent
  452. # define IEEE
  453. # define LITTLE_ENDIAN
  454. #endif /* sequent */
  455. #ifdef sun
  456. # define IEEE
  457. #endif /* sun */
  458. #ifdef NeXT
  459. # define IEEE
  460. #endif /* NeXT */
  461.  
  462. #ifdef MAIN
  463.  
  464. union SParts {
  465.     Single s;
  466.     long i;
  467. };
  468. union DParts {
  469.     Double d;
  470.     long i[2];
  471. };
  472. union EParts {
  473.     defdouble e;
  474.     short i[6];
  475. };
  476.  
  477.  
  478. int
  479. GetHexValue(register int x)
  480. {
  481.     x &= 0x7F;
  482.  
  483.     if ('0' <= x && x <= '9')
  484.         x -= '0';
  485.     else if ('a' <= x && x <= 'f')
  486.         x = x - 'a' + 0xA;
  487.     else if ('A' <= x && x <= 'F')
  488.         x = x - 'A' + 0xA;
  489.     else
  490.         x = 0;
  491.  
  492.     return(x);
  493. }
  494.  
  495.  
  496. void
  497. Hex2Bytes(register char *hex, register char *bytes)
  498. {
  499.     for ( ; *hex; hex += 2) {
  500.         *bytes++ = (GetHexValue(hex[0]) << 4) | GetHexValue(hex[1]);
  501.         if (hex[1] == 0)
  502.             break;    /* Guard against odd bytes */
  503.     }
  504. }
  505.  
  506.  
  507. int
  508. GetHexSymbol(register int x)
  509. {
  510.     x &= 0xF;
  511.     if (x <= 9)
  512.         x += '0';
  513.     else
  514.         x += 'A' - 0xA;
  515.     return(x);
  516. }
  517.  
  518.  
  519. void
  520. Bytes2Hex(register char *bytes, register char *hex, register int nBytes)
  521. {
  522.     for ( ; nBytes--; bytes++) {
  523.         *hex++ = GetHexSymbol(*bytes >> 4);
  524.         *hex++ = GetHexSymbol(*bytes);
  525.     }
  526.     *hex = 0;
  527. }
  528.  
  529.  
  530. void
  531. MaybeSwapBytes(char* bytes, int nBytes)
  532. {
  533. #ifdef LITTLE_ENDIAN
  534.     register char *p, *q, t;
  535.     for (p = bytes, q = bytes+nBytes-1; p < q; p++, q--) {
  536.         t = *p;
  537.         *p = *q;
  538.         *q = t;
  539.     }
  540. #else
  541.     if (bytes, nBytes);        /* Just so it's used */
  542. #endif /* LITTLE_ENDIAN */
  543.  
  544. }
  545.  
  546.  
  547. float
  548. MachineIEEESingle(char* bytes)
  549. {
  550.     float t;
  551.     MaybeSwapBytes(bytes, 4);
  552.     t = *((float*)(bytes));
  553.     MaybeSwapBytes(bytes, 4);
  554.     return (t);
  555. }
  556.  
  557.  
  558. Double
  559. MachineIEEEDouble(char* bytes)
  560. {
  561.     Double t;
  562.     MaybeSwapBytes(bytes, 8);
  563.     t = *((Double*)(bytes));
  564.     MaybeSwapBytes(bytes, 8);
  565.     return (t);
  566. }
  567.  
  568.  
  569. void
  570. TestFromIeeeSingle(char *hex)
  571. {
  572.     defdouble f;
  573.     union SParts p;
  574.     char bytes[4];
  575.  
  576.     Hex2Bytes(hex, bytes);
  577.     f = ConvertFromIeeeSingle(bytes);
  578.     p.s = f;
  579.  
  580. #ifdef IEEE
  581.     fprintf(stderr, "IEEE(%g) [%s] --> float(%g) [%08lX]\n",
  582.     MachineIEEESingle(bytes),
  583.     hex, f, p.i);
  584. #else /* IEEE */
  585.     fprintf(stderr, "IEEE[%s] --> float(%g) [%08lX]\n", hex, f, p.i);
  586. #endif /* IEEE */
  587. }
  588.  
  589.  
  590. void
  591. TestToIeeeSingle(defdouble f)
  592. {
  593.     union SParts p;
  594.     char bytes[4];
  595.     char hex[8+1];
  596.  
  597.     p.s = f;
  598.  
  599.     ConvertToIeeeSingle(f, bytes);
  600.     Bytes2Hex(bytes, hex, 4);
  601. #ifdef IEEE
  602.     fprintf(stderr, "float(%g) [%08lX] --> IEEE(%g) [%s]\n",
  603.         f, p.i,
  604.         MachineIEEESingle(bytes),
  605.         hex
  606.     );
  607. #else /* IEEE */
  608.     fprintf(stderr, "float(%g) [%08lX] --> IEEE[%s]\n", f, p.i, hex);
  609. #endif /* IEEE */
  610. }
  611.  
  612.  
  613. void
  614. TestFromIeeeDouble(char *hex)
  615. {
  616.     defdouble f;
  617.     union DParts p;
  618.     char bytes[8];
  619.  
  620.     Hex2Bytes(hex, bytes);
  621.     f = ConvertFromIeeeDouble(bytes);
  622.     p.d = f;
  623.  
  624. #ifdef IEEE
  625.     fprintf(stderr, "IEEE(%g) [%.8s %.8s] --> double(%g) [%08lX %08lX]\n",
  626.     MachineIEEEDouble(bytes),
  627.     hex, hex+8, f, p.i[0], p.i[1]);
  628. #else /* IEEE */
  629.     fprintf(stderr, "IEEE[%.8s %.8s] --> double(%g) [%08lX %08lX]\n",
  630.         hex, hex+8, f, p.i[0], p.i[1]);
  631. #endif /* IEEE */
  632.  
  633. }
  634.  
  635. void
  636. TestToIeeeDouble(defdouble f)
  637. {
  638.     union DParts p;
  639.     char bytes[8];
  640.     char hex[16+1];
  641.  
  642.     p.d = f;
  643.  
  644.     ConvertToIeeeDouble(f, bytes);
  645.     Bytes2Hex(bytes, hex, 8);
  646. #ifdef IEEE
  647.     fprintf(stderr, "double(%g) [%08lX %08lX] --> IEEE(%g) [%.8s %.8s]\n",
  648.         f, p.i[0], p.i[1],
  649.         MachineIEEEDouble(bytes),
  650.         hex, hex+8
  651.     );
  652. #else /* IEEE */
  653.     fprintf(stderr, "double(%g) [%08lX %08lX] --> IEEE[%.8s %.8s]\n",
  654.         f, p.i[0], p.i[1], hex, hex+8
  655.     );
  656. #endif /* IEEE */
  657.  
  658. }
  659.  
  660.  
  661. void
  662. TestFromIeeeExtended(char *hex)
  663. {
  664.     defdouble f;
  665.     union EParts p;
  666.     char bytes[12];
  667.  
  668.     Hex2Bytes(hex, bytes);
  669.     f = ConvertFromIeeeExtended(bytes);
  670.     p.e = f;
  671.  
  672.     bytes[11] = bytes[9];
  673.     bytes[10] = bytes[8];
  674.     bytes[9] = bytes[7];
  675.     bytes[8] = bytes[6];
  676.     bytes[7] = bytes[5];
  677.     bytes[6] = bytes[4];
  678.     bytes[5] = bytes[3];
  679.     bytes[4] = bytes[2];
  680.     bytes[3] = 0;
  681.     bytes[2] = 0;
  682.  
  683. #if defined(applec) || defined(THINK_C) || defined(METROWERKS)
  684.     fprintf(stderr, "IEEE(%g) [%.4s %.8s %.8s] --> extended(%g) [%04X %04X%04X %04X%04X]\n",
  685.         *((defdouble*)(bytes)),
  686.         hex, hex+4, hex+12, f,
  687.         p.i[0]&0xFFFF, p.i[2]&0xFFFF, p.i[3]&0xFFFF, p.i[4]&0xFFFF, p.i[5]&0xFFFF
  688.     );
  689. #else /* !Macintosh */
  690.     fprintf(stderr, "IEEE[%.4s %.8s %.8s] --> extended(%g) [%04X %04X%04X %04X%04X]\n",
  691.         hex, hex+4, hex+12, f,
  692.         p.i[0]&0xFFFF, p.i[2]&0xFFFF, p.i[3]&0xFFFF, p.i[4]&0xFFFF, p.i[5]&0xFFFF
  693.     );
  694. #endif /* Macintosh */
  695. }
  696.  
  697.  
  698. void
  699. TestToIeeeExtended(defdouble f)
  700. {
  701.     char bytes[12];
  702.     char hex[24+1];
  703.  
  704.     ConvertToIeeeExtended(f, bytes);
  705.     Bytes2Hex(bytes, hex, 10);
  706.  
  707.     bytes[11] = bytes[9];
  708.     bytes[10] = bytes[8];
  709.     bytes[9] = bytes[7];
  710.     bytes[8] = bytes[6];
  711.     bytes[7] = bytes[5];
  712.     bytes[6] = bytes[4];
  713.     bytes[5] = bytes[3];
  714.     bytes[4] = bytes[2];
  715.     bytes[3] = 0;
  716.     bytes[2] = 0;
  717.  
  718. #if defined(applec) || defined(THINK_C) || defined(METROWERKS)
  719.     fprintf(stderr, "extended(%g) --> IEEE(%g) [%.4s %.8s %.8s]\n",
  720.         f, *((defdouble*)(bytes)),
  721.         hex, hex+4, hex+12
  722.     );
  723. #else /* !Macintosh */
  724.     fprintf(stderr, "extended(%g) --> IEEE[%.4s %.8s %.8s]\n",
  725.         f,
  726.         hex, hex+4, hex+12
  727.     );
  728. #endif /* Macintosh */
  729. }
  730.  
  731. #include    <signal.h>
  732.  
  733. void SignalFPE(int i, void (*j)())
  734. {
  735.     printf("[Floating Point Interrupt Caught.]\n", i, j);
  736.     signal(SIGFPE, SignalFPE);
  737. }
  738.  
  739. void
  740. main(void)
  741. {
  742.     long d[3];
  743.     char bytes[12];
  744.  
  745.     signal(SIGFPE, SignalFPE);
  746.  
  747.     TestFromIeeeSingle("00000000");
  748.     TestFromIeeeSingle("80000000");
  749.     TestFromIeeeSingle("3F800000");
  750.     TestFromIeeeSingle("BF800000");
  751.     TestFromIeeeSingle("40000000");
  752.     TestFromIeeeSingle("C0000000");
  753.     TestFromIeeeSingle("7F800000");
  754.     TestFromIeeeSingle("FF800000");
  755.     TestFromIeeeSingle("00800000");
  756.     TestFromIeeeSingle("00400000");
  757.     TestFromIeeeSingle("00000001");
  758.     TestFromIeeeSingle("80000001");
  759.     TestFromIeeeSingle("3F8FEDCB");
  760.     TestFromIeeeSingle("7FC00100");    /* Quiet NaN(1) */
  761.     TestFromIeeeSingle("7F800100");    /* Signalling NaN(1) */
  762.  
  763.     TestToIeeeSingle(0.0);
  764.     TestToIeeeSingle(-0.0);
  765.     TestToIeeeSingle(1.0);
  766.     TestToIeeeSingle(-1.0);
  767.     TestToIeeeSingle(2.0);
  768.     TestToIeeeSingle(-2.0);
  769.     TestToIeeeSingle(3.0);
  770.     TestToIeeeSingle(-3.0);
  771. #if !(defined(sgi) || defined(NeXT))
  772.     TestToIeeeSingle(HUGE_VAL);
  773.     TestToIeeeSingle(-HUGE_VAL);
  774. #endif
  775.  
  776. #ifdef IEEE
  777.     /* These only work on big-endian IEEE machines */
  778.     d[0] = 0x00800000L; MaybeSwapBytes(d,4); TestToIeeeSingle(*((float*)(&d[0])));        /* Smallest normalized */
  779.     d[0] = 0x00400000L; MaybeSwapBytes(d,4); TestToIeeeSingle(*((float*)(&d[0])));        /* Almost largest denormalized */
  780.     d[0] = 0x00000001L; MaybeSwapBytes(d,4); TestToIeeeSingle(*((float*)(&d[0])));        /* Smallest denormalized */
  781.     d[0] = 0x00000001L; MaybeSwapBytes(d,4); TestToIeeeSingle(*((float*)(&d[0])) * 0.5);    /* Smaller than smallest denorm */
  782.     d[0] = 0x3F8FEDCBL; MaybeSwapBytes(d,4); TestToIeeeSingle(*((float*)(&d[0])));
  783. #if !(defined(sgi) || defined(NeXT))
  784.     d[0] = 0x7FC00100L; MaybeSwapBytes(d,4); TestToIeeeSingle(*((float*)(&d[0])));        /* Quiet NaN(1) */
  785.     d[0] = 0x7F800100L; MaybeSwapBytes(d,4); TestToIeeeSingle(*((float*)(&d[0])));        /* Signalling NaN(1) */
  786. #endif /* sgi */
  787. #endif /* IEEE */
  788.  
  789.  
  790.  
  791.     TestFromIeeeDouble("0000000000000000");
  792.     TestFromIeeeDouble("8000000000000000");
  793.     TestFromIeeeDouble("3FF0000000000000");
  794.     TestFromIeeeDouble("BFF0000000000000");
  795.     TestFromIeeeDouble("4000000000000000");
  796.     TestFromIeeeDouble("C000000000000000");
  797.     TestFromIeeeDouble("7FF0000000000000");
  798.     TestFromIeeeDouble("FFF0000000000000");
  799.     TestFromIeeeDouble("0010000000000000");
  800.     TestFromIeeeDouble("0008000000000000");
  801.     TestFromIeeeDouble("0000000000000001");
  802.     TestFromIeeeDouble("8000000000000001");
  803.     TestFromIeeeDouble("3FFFEDCBA9876543");
  804.     TestFromIeeeDouble("7FF8002000000000");    /* Quiet NaN(1) */
  805.     TestFromIeeeDouble("7FF0002000000000");    /* Signalling NaN(1) */
  806.  
  807.     TestToIeeeDouble(0.0);
  808.     TestToIeeeDouble(-0.0);
  809.     TestToIeeeDouble(1.0);
  810.     TestToIeeeDouble(-1.0);
  811.     TestToIeeeDouble(2.0);
  812.     TestToIeeeDouble(-2.0);
  813.     TestToIeeeDouble(3.0);
  814.     TestToIeeeDouble(-3.0);
  815. #if !(defined(sgi) || defined(NeXT))
  816.     TestToIeeeDouble(HUGE_VAL);
  817.     TestToIeeeDouble(-HUGE_VAL);
  818. #endif
  819.  
  820. #ifdef IEEE
  821.     /* These only work on big-endian IEEE machines */
  822.     Hex2Bytes("0010000000000000", bytes); MaybeSwapBytes(d,8); TestToIeeeDouble(*((Double*)(bytes)));    /* Smallest normalized */
  823.     Hex2Bytes("0010000080000000", bytes); MaybeSwapBytes(d,8); TestToIeeeDouble(*((Double*)(bytes)));    /* Normalized, problem with unsigned */
  824.     Hex2Bytes("0008000000000000", bytes); MaybeSwapBytes(d,8); TestToIeeeDouble(*((Double*)(bytes)));    /* Almost largest denormalized */
  825.     Hex2Bytes("0000000080000000", bytes); MaybeSwapBytes(d,8); TestToIeeeDouble(*((Double*)(bytes)));    /* Denorm problem with unsigned */
  826.     Hex2Bytes("0000000000000001", bytes); MaybeSwapBytes(d,8); TestToIeeeDouble(*((Double*)(bytes)));    /* Smallest denormalized */
  827.     Hex2Bytes("0000000000000001", bytes); MaybeSwapBytes(d,8); TestToIeeeDouble(*((Double*)(bytes)) * 0.5);    /* Smaller than smallest denorm */
  828.     Hex2Bytes("3FFFEDCBA9876543", bytes); MaybeSwapBytes(d,8); TestToIeeeDouble(*((Double*)(bytes)));    /* accuracy test */
  829. #if !(defined(sgi) || defined(NeXT))
  830.     Hex2Bytes("7FF8002000000000", bytes); MaybeSwapBytes(d,8); TestToIeeeDouble(*((Double*)(bytes)));    /* Quiet NaN(1) */
  831.     Hex2Bytes("7FF0002000000000", bytes); MaybeSwapBytes(d,8); TestToIeeeDouble(*((Double*)(bytes)));    /* Signalling NaN(1) */
  832. #endif /* sgi */
  833. #endif /* IEEE */
  834.  
  835.     TestFromIeeeExtended("00000000000000000000");    /* +0 */
  836.     TestFromIeeeExtended("80000000000000000000");    /* -0 */
  837.     TestFromIeeeExtended("3FFF8000000000000000");    /* +1 */
  838.     TestFromIeeeExtended("BFFF8000000000000000");    /* -1 */
  839.     TestFromIeeeExtended("40008000000000000000");    /* +2 */
  840.     TestFromIeeeExtended("C0008000000000000000");    /* -2 */
  841.     TestFromIeeeExtended("7FFF0000000000000000");    /* +infinity */
  842.     TestFromIeeeExtended("FFFF0000000000000000");    /* -infinity */
  843.     TestFromIeeeExtended("7FFF8001000000000000");    /* Quiet NaN(1) */
  844.     TestFromIeeeExtended("7FFF0001000000000000");    /* Signalling NaN(1) */
  845.     TestFromIeeeExtended("3FFFFEDCBA9876543210");    /* accuracy test */
  846.  
  847.     TestToIeeeExtended(0.0);
  848.     TestToIeeeExtended(-0.0);
  849.     TestToIeeeExtended(1.0);
  850.     TestToIeeeExtended(-1.0);
  851.     TestToIeeeExtended(2.0);
  852.     TestToIeeeExtended(-2.0);
  853. #if !(defined(sgi) || defined(NeXT))
  854.     TestToIeeeExtended(HUGE_VAL);
  855.     TestToIeeeExtended(-HUGE_VAL);
  856. #endif /* sgi */
  857.  
  858. #if defined(applec) || defined(THINK_C) || defined(METROWERKS)
  859.     Hex2Bytes("7FFF00008001000000000000", bytes); TestToIeeeExtended(*((long double*)(bytes)));    /* Quiet NaN(1) */
  860.     Hex2Bytes("7FFF00000001000000000000", bytes); TestToIeeeExtended(*((long double*)(bytes)));    /* Signalling NaN(1) */
  861.     Hex2Bytes("7FFE00008000000000000000", bytes); TestToIeeeExtended(*((long double*)(bytes)));
  862.     Hex2Bytes("000000008000000000000000", bytes); TestToIeeeExtended(*((long double*)(bytes)));
  863.     Hex2Bytes("000000000000000000000001", bytes); TestToIeeeExtended(*((long double*)(bytes)));
  864.     Hex2Bytes("3FFF0000FEDCBA9876543210", bytes); TestToIeeeExtended(*((long double*)(bytes)));
  865. #endif
  866. }
  867.  
  868.  
  869. /* This is the output of the test program on an IEEE machine:
  870. IEEE(0) [00000000] --> float(0) [00000000]
  871. IEEE(-0) [80000000] --> float(-0) [80000000]
  872. IEEE(1) [3F800000] --> float(1) [3F800000]
  873. IEEE(-1) [BF800000] --> float(-1) [BF800000]
  874. IEEE(2) [40000000] --> float(2) [40000000]
  875. IEEE(-2) [C0000000] --> float(-2) [C0000000]
  876. IEEE(INF) [7F800000] --> float(INF) [7F800000]
  877. IEEE(-INF) [FF800000] --> float(-INF) [FF800000]
  878. IEEE(1.17549e-38) [00800000] --> float(1.17549e-38) [00800000]
  879. IEEE(5.87747e-39) [00400000] --> float(5.87747e-39) [00400000]
  880. IEEE(1.4013e-45) [00000001] --> float(1.4013e-45) [00000001]
  881. IEEE(-1.4013e-45) [80000001] --> float(-1.4013e-45) [80000001]
  882. IEEE(1.12444) [3F8FEDCB] --> float(1.12444) [3F8FEDCB]
  883. IEEE(NAN(001)) [7FC00100] --> float(INF) [7F800000]
  884. IEEE(NAN(001)) [7F800100] --> float(INF) [7F800000]
  885. float(0) [00000000] --> IEEE(0) [00000000]
  886. float(-0) [80000000] --> IEEE(0) [00000000]
  887. float(1) [3F800000] --> IEEE(1) [3F800000]
  888. float(-1) [BF800000] --> IEEE(-1) [BF800000]
  889. float(2) [40000000] --> IEEE(2) [40000000]
  890. float(-2) [C0000000] --> IEEE(-2) [C0000000]
  891. float(3) [40400000] --> IEEE(3) [40400000]
  892. float(-3) [C0400000] --> IEEE(-3) [C0400000]
  893. float(INF) [7F800000] --> IEEE(INF) [7F800000]
  894. float(-INF) [FF800000] --> IEEE(-INF) [FF800000]
  895. float(1.17549e-38) [00800000] --> IEEE(1.17549e-38) [00800000]
  896. float(5.87747e-39) [00400000] --> IEEE(5.87747e-39) [00400000]
  897. float(1.4013e-45) [00000001] --> IEEE(1.4013e-45) [00000001]
  898. float(7.00649e-46) [00000000] --> IEEE(0) [00000000]
  899. float(1.12444) [3F8FEDCB] --> IEEE(1.12444) [3F8FEDCB]
  900. float(NAN(001)) [7FC00100] --> IEEE(INF) [7F800000]
  901. float(NAN(001)) [7FC00100] --> IEEE(INF) [7F800000]
  902. IEEE(0) [00000000 00000000] --> double(0) [00000000 00000000]
  903. IEEE(-0) [80000000 00000000] --> double(-0) [80000000 00000000]
  904. IEEE(1) [3FF00000 00000000] --> double(1) [3FF00000 00000000]
  905. IEEE(-1) [BFF00000 00000000] --> double(-1) [BFF00000 00000000]
  906. IEEE(2) [40000000 00000000] --> double(2) [40000000 00000000]
  907. IEEE(-2) [C0000000 00000000] --> double(-2) [C0000000 00000000]
  908. IEEE(INF) [7FF00000 00000000] --> double(INF) [7FF00000 00000000]
  909. IEEE(-INF) [FFF00000 00000000] --> double(-INF) [FFF00000 00000000]
  910. IEEE(2.22507e-308) [00100000 00000000] --> double(2.22507e-308) [00100000 00000000]
  911. IEEE(1.11254e-308) [00080000 00000000] --> double(1.11254e-308) [00080000 00000000]
  912. IEEE(4.94066e-324) [00000000 00000001] --> double(4.94066e-324) [00000000 00000001]
  913. IEEE(-4.94066e-324) [80000000 00000001] --> double(-4.94066e-324) [80000000 00000001]
  914. IEEE(1.99556) [3FFFEDCB A9876543] --> double(1.99556) [3FFFEDCB A9876543]
  915. IEEE(NAN(001)) [7FF80020 00000000] --> double(INF) [7FF00000 00000000]
  916. IEEE(NAN(001)) [7FF00020 00000000] --> double(INF) [7FF00000 00000000]
  917. double(0) [00000000 00000000] --> IEEE(0) [00000000 00000000]
  918. double(-0) [80000000 00000000] --> IEEE(0) [00000000 00000000]
  919. double(1) [3FF00000 00000000] --> IEEE(1) [3FF00000 00000000]
  920. double(-1) [BFF00000 00000000] --> IEEE(-1) [BFF00000 00000000]
  921. double(2) [40000000 00000000] --> IEEE(2) [40000000 00000000]
  922. double(-2) [C0000000 00000000] --> IEEE(-2) [C0000000 00000000]
  923. double(3) [40080000 00000000] --> IEEE(3) [40080000 00000000]
  924. double(-3) [C0080000 00000000] --> IEEE(-3) [C0080000 00000000]
  925. double(INF) [7FF00000 00000000] --> IEEE(INF) [7FF00000 00000000]
  926. double(-INF) [FFF00000 00000000] --> IEEE(-INF) [FFF00000 00000000]
  927. double(2.22507e-308) [00100000 00000000] --> IEEE(2.22507e-308) [00100000 00000000]
  928. double(2.22507e-308) [00100000 80000000] --> IEEE(2.22507e-308) [00100000 80000000]
  929. double(1.11254e-308) [00080000 00000000] --> IEEE(1.11254e-308) [00080000 00000000]
  930. double(1.061e-314) [00000000 80000000] --> IEEE(1.061e-314) [00000000 80000000]
  931. double(4.94066e-324) [00000000 00000001] --> IEEE(4.94066e-324) [00000000 00000001]
  932. double(4.94066e-324) [00000000 00000001] --> IEEE(4.94066e-324) [00000000 00000001]
  933. double(1.99556) [3FFFEDCB A9876543] --> IEEE(1.99556) [3FFFEDCB A9876543]
  934. double(NAN(001)) [7FF80020 00000000] --> IEEE(INF) [7FF00000 00000000]
  935. double(NAN(001)) [7FF80020 00000000] --> IEEE(INF) [7FF00000 00000000]
  936. IEEE(0) [0000 00000000 00000000] --> extended(0) [0000 00000000 00000000]
  937. IEEE(-0) [8000 00000000 00000000] --> extended(-0) [8000 00000000 00000000]
  938. IEEE(1) [3FFF 80000000 00000000] --> extended(1) [3FFF 80000000 00000000]
  939. IEEE(-1) [BFFF 80000000 00000000] --> extended(-1) [BFFF 80000000 00000000]
  940. IEEE(2) [4000 80000000 00000000] --> extended(2) [4000 80000000 00000000]
  941. IEEE(-2) [C000 80000000 00000000] --> extended(-2) [C000 80000000 00000000]
  942. IEEE(INF) [7FFF 00000000 00000000] --> extended(INF) [7FFF 00000000 00000000]
  943. IEEE(-INF) [FFFF 00000000 00000000] --> extended(-INF) [FFFF 00000000 00000000]
  944. IEEE(NAN(001)) [7FFF 80010000 00000000] --> extended(INF) [7FFF 00000000 00000000]
  945. IEEE(NAN(001)) [7FFF 00010000 00000000] --> extended(INF) [7FFF 00000000 00000000]
  946. IEEE(1.99111) [3FFF FEDCBA98 76543210] --> extended(1.99111) [3FFF FEDCBA98 76543210]
  947. extended(0) --> IEEE(0) [0000 00000000 00000000]
  948. extended(-0) --> IEEE(0) [0000 00000000 00000000]
  949. extended(1) --> IEEE(1) [3FFF 80000000 00000000]
  950. extended(-1) --> IEEE(-1) [BFFF 80000000 00000000]
  951. extended(2) --> IEEE(2) [4000 80000000 00000000]
  952. extended(-2) --> IEEE(-2) [C000 80000000 00000000]
  953. extended(INF) --> IEEE(INF) [7FFF 00000000 00000000]
  954. extended(-INF) --> IEEE(-INF) [FFFF 00000000 00000000]
  955. extended(NAN(001)) --> IEEE(INF) [7FFF 00000000 00000000]
  956. extended(NAN(001)) --> IEEE(INF) [7FFF 00000000 00000000]
  957. extended(5.94866e+4931) --> IEEE(5.94866e+4931) [7FFE 80000000 00000000]
  958. extended(1e-4927) --> IEEE(1e-4927) [0000 80000000 00000000]
  959. extended(1e-4927) --> IEEE(1e-4927) [0000 00000000 00000001]
  960. extended(1.99111) --> IEEE(1.99111) [3FFF FEDCBA98 76543210]
  961. */
  962.  
  963. #endif /* TEST_FP */
  964.