ieeexdbl.c
4.44 KB
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
/***************************************************************************
* SGI Audio File Library internal subroutines
*
* _AFConvertFromIeeeExtended
* _AFConvertToIeeeExtended
*
* convert to/from ieee 80-bit extended floating point
* from apple computer, june 91
***************************************************************************/
#ident "$Revision: 1.1.1.1 $"
#include <math.h>
/*
* Copyright (C) 1988-1991 Apple Computer, Inc.
* All rights reserved.
*
* Machine-independent I/O routines for IEEE floating-point numbers.
*
* NaN's and infinities are converted to HUGE_VAL or HUGE, which
* happens to be infinity on IEEE machines. Unfortunately, it is
* impossible to preserve NaN's in a machine-independent way.
* Infinities are, however, preserved on IEEE machines.
*
* These routines have been tested on the following machines:
* Apple Macintosh, MPW 3.1 C compiler
* Apple Macintosh, THINK C compiler
* Silicon Graphics IRIS, MIPS compiler
* Cray X/MP and Y/MP
* Digital Equipment VAX
*
*
* Implemented by Malcolm Slaney and Ken Turkowski.
*
* Malcolm Slaney contributions during 1988-1990 include big- and little-
* endian file I/O, conversion to and from Motorola's extended 80-bit
* floating-point format, and conversions to and from IEEE single-
* precision floating-point format.
*
* In 1991, Ken Turkowski implemented the conversions to and from
* IEEE double-precision format, added more precision to the extended
* conversions, and accommodated conversions involving +/- infinity,
* NaN's, and denormalized numbers.
*/
#ifndef HUGE_VAL
# define HUGE_VAL HUGE
#endif /*HUGE_VAL*/
# define UnsignedToFloat(u) \
(((double)((long)(u - 2147483647L - 1))) + 2147483648.0)
# define FloatToUnsigned(f) \
((unsigned long)(((long)(f - 2147483648.0)) + 2147483647L + 1))
/****************************************************************
* Extended precision IEEE floating-point conversion routines.
****************************************************************/
/*
* C O N V E R T F R O M I E E E E X T E N D E D
*/
double
_AFConvertFromIeeeExtended(char *bytes)
{
double f;
long expon;
unsigned long hiMant, loMant;
expon = ((bytes[0] & 0x7F) << 8) | (bytes[1] & 0xFF);
hiMant = ((unsigned long)(bytes[2] & 0xFF) << 24)
| ((unsigned long)(bytes[3] & 0xFF) << 16)
| ((unsigned long)(bytes[4] & 0xFF) << 8)
| ((unsigned long)(bytes[5] & 0xFF));
loMant = ((unsigned long)(bytes[6] & 0xFF) << 24)
| ((unsigned long)(bytes[7] & 0xFF) << 16)
| ((unsigned long)(bytes[8] & 0xFF) << 8)
| ((unsigned long)(bytes[9] & 0xFF));
if (expon == 0 && hiMant == 0 && loMant == 0) {
f = 0;
}
else {
if (expon == 0x7FFF) { /* Infinity or NaN */
f = HUGE_VAL;
}
else {
expon -= 16383;
f = ldexp(UnsignedToFloat(hiMant), expon-=31);
f += ldexp(UnsignedToFloat(loMant), expon-=32);
}
}
if (bytes[0] & 0x80)
return -f;
else
return f;
}
/*
* C O N V E R T T O I E E E E X T E N D E D
*/
void
_AFConvertToIeeeExtended(double num, char *bytes)
{
int sign;
int expon;
double fMant, fsMant;
unsigned long hiMant, loMant;
if (num < 0) {
sign = 0x8000;
num *= -1;
} else {
sign = 0;
}
if (num == 0) {
expon = 0; hiMant = 0; loMant = 0;
}
else {
fMant = frexp(num, &expon);
if ((expon > 16384) || !(fMant < 1)) { /* Infinity or NaN */
expon = sign|0x7FFF; hiMant = 0; loMant = 0; /* infinity */
}
else { /* Finite */
expon += 16382;
if (expon < 0) { /* denormalized */
fMant = ldexp(fMant, expon);
expon = 0;
}
expon |= sign;
fMant = ldexp(fMant, 32);
fsMant = floor(fMant);
hiMant = FloatToUnsigned(fsMant);
fMant = ldexp(fMant - fsMant, 32);
fsMant = floor(fMant);
loMant = FloatToUnsigned(fsMant);
}
}
bytes[0] = expon >> 8;
bytes[1] = expon;
bytes[2] = hiMant >> 24;
bytes[3] = hiMant >> 16;
bytes[4] = hiMant >> 8;
bytes[5] = hiMant;
bytes[6] = loMant >> 24;
bytes[7] = loMant >> 16;
bytes[8] = loMant >> 8;
bytes[9] = loMant;
}