Loading...
Searching...
No Matches
EST_wave_utils.cc
/*************************************************************************/
/* */
/* Centre for Speech Technology Research */
/* University of Edinburgh, UK */
/* Copyright (c) 1996 */
/* All Rights Reserved. */
/* */
/* Permission is hereby granted, free of charge, to use and distribute */
/* this software and its documentation without restriction, including */
/* without limitation the rights to use, copy, modify, merge, publish, */
/* distribute, sublicense, and/or sell copies of this work, and to */
/* permit persons to whom this work is furnished to do so, subject to */
/* the following conditions: */
/* 1. The code must retain the above copyright notice, this list of */
/* conditions and the following disclaimer. */
/* 2. Any modifications must be clearly marked as such. */
/* 3. Original authors' names are not deleted. */
/* 4. The authors' names are not used to endorse or promote products */
/* derived from this software without specific prior written */
/* permission. */
/* */
/* THE UNIVERSITY OF EDINBURGH AND THE CONTRIBUTORS TO THIS WORK */
/* DISCLAIM ALL WARRANTIES WITH REGARD TO THIS SOFTWARE, INCLUDING */
/* ALL IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS, IN NO EVENT */
/* SHALL THE UNIVERSITY OF EDINBURGH NOR THE CONTRIBUTORS BE LIABLE */
/* FOR ANY SPECIAL, INDIRECT OR CONSEQUENTIAL DAMAGES OR ANY DAMAGES */
/* WHATSOEVER RESULTING FROM LOSS OF USE, DATA OR PROFITS, WHETHER IN */
/* AN ACTION OF CONTRACT, NEGLIGENCE OR OTHER TORTIOUS ACTION, */
/* ARISING OUT OF OR IN CONNECTION WITH THE USE OR PERFORMANCE OF */
/* THIS SOFTWARE. */
/* */
/*************************************************************************/
/* Author : Alan Black and Paul Taylor */
/* Date : June 1996 */
/*-----------------------------------------------------------------------*/
/* Various low-level waveform conversion routines and file format */
/* independent i/o functions */
/* */
/* Acknowledgements */
/* ulaw conversion code provided by */
/* Craig Reese: IDA/Supercomputing Research Center */
/* IEEE extended conversion */
/* Apple Computer, Inc. */
/* */
/*=======================================================================*/
#include <cstdio>
#include <cstdlib>
#include "EST_unix.h"
#include <cstring>
#include <cmath>
#include "EST_wave_utils.h"
#include "EST_wave_aux.h"
#include "EST_error.h"
static short st_ulaw_to_short(unsigned char ulawbyte);
static short st_alaw_to_short(unsigned char alawbyte);
static unsigned char st_short_to_ulaw(short sample);
static unsigned char st_short_to_alaw(short sample);
/*
* This table is
* Copyright 1992 by Jutta Degener and Carsten Bormann, Technische
* Universitaet Berlin. See the accompanying file "COPYRIGHT" for
* details. THERE IS ABSOLUTELY NO WARRANTY FOR THIS SOFTWARE.
* take from speak_freely-7.2/gsm/src/toast_alaw.c
*/
static unsigned short a2s[] = {
5120,60160, 320,65200,20480,44032, 1280,64192,
2560,62848, 64,65456,10240,54784, 640,64864,
7168,58112, 448,65072,28672,35840, 1792,63680,
3584,61824, 192,65328,14336,50688, 896,64608,
4096,61184, 256,65264,16384,48128, 1024,64448,
2048,63360, 0,65520, 8192,56832, 512,64992,
6144,59136, 384,65136,24576,39936, 1536,63936,
3072,62336, 128,65392,12288,52736, 768,64736,
5632,59648, 352,65168,22528,41984, 1408,64064,
2816,62592, 96,65424,11264,53760, 704,64800,
7680,57600, 480,65040,30720,33792, 1920,63552,
3840,61568, 224,65296,15360,49664, 960,64544,
4608,60672, 288,65232,18432,46080, 1152,64320,
2304,63104, 32,65488, 9216,55808, 576,64928,
6656,58624, 416,65104,26624,37888, 1664,63808,
3328,62080, 160,65360,13312,51712, 832,64672,
5376,59904, 336,65184,21504,43008, 1344,64128,
2688,62720, 80,65440,10752,54272, 672,64832,
7424,57856, 464,65056,29696,34816, 1856,63616,
3712,61696, 208,65312,14848,50176, 928,64576,
4352,60928, 272,65248,17408,47104, 1088,64384,
2176,63232, 16,65504, 8704,56320, 544,64960,
6400,58880, 400,65120,25600,38912, 1600,63872,
3200,62208, 144,65376,12800,52224, 800,64704,
5888,59392, 368,65152,23552,40960, 1472,64000,
2944,62464, 112,65408,11776,53248, 736,64768,
7936,57344, 496,65024,31744,32768, 1984,63488,
3968,61440, 240,65280,15872,49152, 992,64512,
4864,60416, 304,65216,19456,45056, 1216,64256,
2432,62976, 48,65472, 9728,55296, 608,64896,
6912,58368, 432,65088,27648,36864, 1728,63744,
3456,61952, 176,65344,13824,51200, 864,64640
};
#define st_alaw_to_short(a) (a2s[(unsigned char)a])
void ulaw_to_short(const unsigned char *ulaw,short *data,int length)
{
/* Convert ulaw to shorts */
int i;
for (i=0; i<length; i++)
data[i] = st_ulaw_to_short(ulaw[i]); /* ulaw convert */
}
void alaw_to_short(const unsigned char *alaw,short *data,int length)
{
/* Convert ulaw to shorts */
int i;
for (i=0; i<length; i++)
{
data[i] = st_alaw_to_short(alaw[i])-32768; /* alaw convert */
}
}
void shorten_to_short(unsigned char *ulaw, short *data, int length)
{
/* Convert Tony Robinson's shorten format to shorts */
(void)ulaw;
(void)data;
(void)length;
}
void uchar_to_short(const unsigned char *chars,short *data,int length)
{
/* Convert 8 bit data to shorts UNSIGNED CHAR */
int i;
for (i=0; i<length; i++)
{
data[i] = (((int)chars[i])-128)*256;
}
}
void schar_to_short(const unsigned char *chars,short *data,int length)
{
/* Convert 8 bit data to shorts SIGNED CHAR */
int i;
for (i=0; i<length; i++)
data[i] = (((unsigned char)chars[i]))*256;
}
void short_to_uchar(const short *data,unsigned char *chars,int length)
{
/* Convert shorts to 8 bit UNSIGNED CHAR */
int i;
for (i=0; i<length; i++)
chars[i] = (data[i]/256)+128;
}
void short_to_schar(const short *data,unsigned char *chars,int length)
{
/* Convert shorts to 8 bit SIGNED CHAR */
int i;
for (i=0; i<length; i++)
chars[i] = (data[i]/256);
}
#if 0
void short_to_adpcm(short *data,signed char *chars,int length)
{
struct adpcm_state state;
state.valprev = 0;
state.index = 0;
adpcm_coder(data,chars,length,&state);
}
void adpcm_to_short(signed char *chars, short *data,int length)
{
struct adpcm_state state;
state.valprev = 0;
state.index = 0;
adpcm_decoder(chars,data,length/2,&state);
}
#endif
void short_to_ulaw(const short *data,unsigned char *ulaw,int length)
{
/* Convert ulaw to shorts */
int i;
for (i=0; i<length; i++)
ulaw[i] = st_short_to_ulaw(data[i]); /* ulaw convert */
}
void short_to_alaw(const short *data,unsigned char *alaw,int length)
{
/* Convert alaw to shorts */
int i;
for (i=0; i<length; i++)
alaw[i] = st_short_to_alaw(data[i]); /* alaw convert */
}
short *convert_raw_data(unsigned char *file_data,int data_length,
enum EST_sample_type_t sample_type,int bo)
{
/* converts data in file_data to native byte order shorts */
/* frees file_data if its not returned */
short *d=NULL;
if (sample_type == st_short)
{
// d = walloc(short,data_length);
if (bo != EST_NATIVE_BO)
swap_bytes_short((short *)file_data,data_length);
return (short *)file_data;
}
else if (sample_type == st_mulaw)
{
d = walloc(short,data_length);
ulaw_to_short(file_data,d,data_length);
wfree(file_data);
return d;
}
else if (sample_type == st_alaw)
{
d = walloc(short,data_length);
alaw_to_short(file_data,d,data_length);
wfree(file_data);
return d;
}
else if (sample_type == st_alaw)
{
d = walloc(short,data_length);
alaw_to_short(file_data,d,data_length);
wfree(file_data);
return d;
}
#if 0
else if (sample_type == st_adpcm)
{ /* Not really checked yet */
d = walloc(short,data_length);
adpcm_to_short((signed char *)file_data,d,data_length);
wfree(file_data);
return d;
}
#endif
else if (sample_type == st_schar)
{
d = walloc(short,data_length);
schar_to_short((unsigned char *)file_data,d,data_length);
wfree(file_data);
return d;
}
else if (sample_type == st_uchar)
{
d = walloc(short,data_length);
uchar_to_short((unsigned char *)file_data,d,data_length);
wfree(file_data);
return d;
}
else
EST_error("Convert raw data: unsupported sample type %s(%d)",
EST_sample_type_map.name(sample_type), sample_type);
/* NOTREACHED */
return NULL;
}
enum EST_write_status save_raw_data(FILE *fp, const short *data, int offset,
int num_samples, int num_channels,
enum EST_sample_type_t sample_type,
int bo)
{
int i;
int n;
if (sample_type == st_mulaw)
{
unsigned char *ulaw = walloc(unsigned char,num_samples*num_channels);
short_to_ulaw(data+(offset*num_channels),
ulaw,num_samples*num_channels);
n = fwrite(ulaw,1,num_channels * num_samples,fp);
wfree(ulaw);
if (n != (num_channels * num_samples))
return misc_write_error;
}
else if (sample_type == st_alaw)
{
unsigned char *alaw = walloc(unsigned char,num_samples*num_channels);
short_to_alaw(data+(offset*num_channels),
alaw,num_samples*num_channels);
n = fwrite(alaw,1,num_channels * num_samples,fp);
wfree(alaw);
if (n != (num_channels * num_samples))
return misc_write_error;
}
else if (sample_type == st_ascii)
{
for (i=offset*num_channels; i < num_samples*num_channels; i++)
fprintf(fp,"%d\n",data[i]);
}
else if (sample_type == st_schar)
{
unsigned char *chars = walloc(unsigned char,num_samples*num_channels);
short_to_schar(data+(offset*num_channels),
chars,num_samples*num_channels);
n = fwrite(chars,1,num_channels * num_samples,fp);
wfree(chars);
if (n != (num_channels * num_samples))
return misc_write_error;
}
else if (sample_type == st_uchar)
{
unsigned char *chars = walloc(unsigned char,num_samples*num_channels);
short_to_uchar(data+(offset*num_channels),
chars,num_samples*num_channels);
n = fwrite(chars,1,num_channels * num_samples,fp);
wfree(chars);
if ( n != (num_channels * num_samples))
return misc_write_error;
}
#if 0
else if (sample_type == st_adpcm)
{
signed char *chars = walloc(signed char,num_samples*num_channels);
short_to_adpcm(data+(offset*num_channels),
chars,num_samples*num_channels);
n = fwrite(chars,1,num_channels * num_samples,fp);
wfree(chars);
if ( n != (num_channels * num_samples))
return misc_write_error;
}
#endif
else if (sample_type == st_short)
{
if (bo != EST_NATIVE_BO)
{
short *xdata = walloc(short,num_channels*num_samples);
memmove(xdata,data+(offset*num_channels),
num_channels*num_samples*sizeof(short));
swap_bytes_short(xdata,num_channels*num_samples);
n = fwrite(xdata, sizeof(short),num_channels * num_samples, fp);
wfree(xdata);
}
else
n = fwrite(&data[offset], sizeof(short),
num_channels * num_samples, fp);
if (n != (num_channels * num_samples))
return misc_write_error;
}
else
{
fprintf(stderr,"save data file: unsupported sample type\n");
return misc_write_error;
}
return write_ok;
}
int get_word_size(enum EST_sample_type_t sample_type)
{
/* Returns word size from type */
int word_size;
switch (sample_type)
{
case st_unknown:
word_size = 2; break;
case st_uchar:
case st_schar:
word_size = 1; break;
case st_mulaw:
case st_alaw:
word_size = 1; break;
#if 0
case st_adpcm: /* maybe I mean 0.5 */
word_size = 1; break;
#endif
case st_short:
word_size = 2; break;
case st_int:
/* Yes I mean 4 not sizeof(int) these are machine independent defs */
word_size = 4; break;
case st_float:
word_size = 4; break;
case st_double:
word_size = 8; break;
default:
fprintf(stderr,"Unknown encoding format error\n");
word_size = 2;
}
return word_size;
}
enum EST_sample_type_t str_to_sample_type(const char *type)
{
/* string to internal value */
if (streq(type,"short"))
return st_short;
if (streq(type,"shorten"))
return st_shorten;
else if ((streq(type,"ulaw")) || (streq(type,"mulaw")))
return st_mulaw;
else if ((streq(type,"char")) || (streq(type,"byte")) ||
(streq(type,"8bit")))
return st_schar;
else if ((streq(type,"unsignedchar")) || (streq(type,"unsignedbyte")) ||
(streq(type,"unsigned8bit")))
return st_uchar;
else if (streq(type,"int"))
return st_int;
#if 0
else if (streq(type,"adpcm"))
return st_adpcm;
#endif
else if ((streq(type,"real")) || (streq(type,"float")) ||
(streq(type,"real4")))
return st_float;
else if ((streq(type,"real8")) || (streq(type,"double")))
return st_double;
else if (streq(type,"alaw"))
return st_alaw;
else if (streq(type,"ascii"))
return st_ascii;
else
{
fprintf(stderr,"Unknown sample type: \"%s\"\n",type);
return st_unknown;
}
}
const char *sample_type_to_str(enum EST_sample_type_t type)
{
switch (type)
{
case st_short: return "short";
case st_shorten: return "shorten";
case st_mulaw: return "ulaw";
case st_alaw: return "alaw";
case st_schar: return "char";
case st_uchar: return "unsignedchar";
case st_int: return "int";
#if 0
case st_adpcm: return "adpcm";
#endif
case st_float: return "float";
case st_double: return "double";
case st_ascii: return "ascii";
case st_unknown: return "unknown";
default:
fprintf(stderr,"Unknown sample_type %d\n",type);
return "very_unknown";
}
}
/*
** This routine converts from linear to ulaw.
**
** Craig Reese: IDA/Supercomputing Research Center
** Joe Campbell: Department of Defense
** 29 September 1989
**
** References:
** 1) CCITT Recommendation G.711 (very difficult to follow)
** 2) "A New Digital Technique for Implementation of Any
** Continuous PCM Companding Law," Villeret, Michel,
** et al. 1973 IEEE Int. Conf. on Communications, Vol 1,
** 1973, pg. 11.12-11.17
** 3) MIL-STD-188-113,"Interoperability and Performance Standards
** for Analog-to_Digital Conversion Techniques,"
** 17 February 1987
**
** Input: Signed 16 bit linear sample
** Output: 8 bit ulaw sample
*/
#define ZEROTRAP /* turn on the trap as per the MIL-STD */
#define BIAS 0x84 /* define the add-in bias for 16 bit samples */
#define CLIP 32635
static unsigned char st_short_to_ulaw(short sample)
{
static int exp_lut[256] = {0,0,1,1,2,2,2,2,3,3,3,3,3,3,3,3,
4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,
5,5,5,5,5,5,5,5,5,5,5,5,5,5,5,5,
5,5,5,5,5,5,5,5,5,5,5,5,5,5,5,5,
6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,
6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,
6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,
6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,
7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,
7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,
7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,
7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,
7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,
7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,
7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,
7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7};
int sign, exponent, mantissa;
unsigned char ulawbyte;
/* Get the sample into sign-magnitude. */
sign = (sample >> 8) & 0x80; /* set aside the sign */
if ( sign != 0 ) sample = -sample; /* get magnitude */
if ( sample > CLIP ) sample = CLIP; /* clip the magnitude */
/* Convert from 16 bit linear to ulaw. */
sample = sample + BIAS;
exponent = exp_lut[( sample >> 7 ) & 0xFF];
mantissa = ( sample >> ( exponent + 3 ) ) & 0x0F;
ulawbyte = ~ ( sign | ( exponent << 4 ) | mantissa );
#ifdef ZEROTRAP
if ( ulawbyte == 0 ) ulawbyte = 0x02; /* optional CCITT trap */
#endif
return ulawbyte;
}
/*
** This routine converts from ulaw to 16 bit linear.
**
** Craig Reese: IDA/Supercomputing Research Center
** 29 September 1989
**
** References:
** 1) CCITT Recommendation G.711 (very difficult to follow)
** 2) MIL-STD-188-113,"Interoperability and Performance Standards
** for Analog-to_Digital Conversion Techniques,"
** 17 February 1987
**
** Input: 8 bit ulaw sample
** Output: signed 16 bit linear sample
*/
static short st_ulaw_to_short( unsigned char ulawbyte )
{
static int exp_lut[8] = { 0, 132, 396, 924, 1980, 4092, 8316, 16764 };
int sign, exponent, mantissa;
short sample;
ulawbyte = ~ ulawbyte;
sign = ( ulawbyte & 0x80 );
exponent = ( ulawbyte >> 4 ) & 0x07;
mantissa = ulawbyte & 0x0F;
sample = exp_lut[exponent] + ( mantissa << ( exponent + 3 ) );
if ( sign != 0 ) sample = -sample;
return sample;
}
/* The following is copied from sox:g711.c */
/*
* This source code is a product of Sun Microsystems, Inc. and is provided
* for unrestricted use. Users may copy or modify this source code without
* charge.
* [ no warranties or liabilities whatsoever; see there for details. ]
*/
/* copy from CCITT G.711 specifications */
static unsigned char _u2a[128] = { /* u- to A-law conversions */
1, 1, 2, 2, 3, 3, 4, 4,
5, 5, 6, 6, 7, 7, 8, 8,
9, 10, 11, 12, 13, 14, 15, 16,
17, 18, 19, 20, 21, 22, 23, 24,
25, 27, 29, 31, 33, 34, 35, 36,
37, 38, 39, 40, 41, 42, 43, 44,
46, 48, 49, 50, 51, 52, 53, 54,
55, 56, 57, 58, 59, 60, 61, 62,
64, 65, 66, 67, 68, 69, 70, 71,
72, 73, 74, 75, 76, 77, 78, 79,
80, 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};
static unsigned char _a2u[128] = { /* A- to u-law conversions */
1, 3, 5, 7, 9, 11, 13, 15,
16, 17, 18, 19, 20, 21, 22, 23,
24, 25, 26, 27, 28, 29, 30, 31,
32, 32, 33, 33, 34, 34, 35, 35,
36, 37, 38, 39, 40, 41, 42, 43,
44, 45, 46, 47, 48, 48, 49, 49,
50, 51, 52, 53, 54, 55, 56, 57,
58, 59, 60, 61, 62, 63, 64, 64,
65, 66, 67, 68, 69, 70, 71, 72,
73, 74, 75, 76, 77, 78, 79, 80,
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};
/* A-law to u-law conversion */
static inline unsigned char st_alaw2ulaw(
unsigned char aval)
{
aval &= 0xff;
return (unsigned char) ((aval & 0x80) ? (0xFF ^ _a2u[aval ^ 0xD5]) :
(0x7F ^ _a2u[aval ^ 0x55]));
}
/* u-law to A-law conversion */
static inline unsigned char st_ulaw2alaw(
unsigned char uval)
{
uval &= 0xff;
return (unsigned char) ((uval & 0x80) ? (0xD5 ^ (_u2a[0xFF ^ uval] - 1)) :
(unsigned char) (0x55 ^ (_u2a[0x7F ^ uval] - 1)));
}
/* end of Sun code */
/* This is somewhat simple-minded, but ... */
static unsigned char st_short_to_alaw(short sample)
{
return st_ulaw2alaw(st_short_to_ulaw(sample));
}
/*
* C O N V E R T T O I E E E E X T E N D E D
*/
/* 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 FloatToUnsigned(f) ((unsigned long)(((long)(f - 2147483648.0)) + 2147483647L) + 1)
void ConvertToIeeeExtended(double num,unsigned 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;
}
/*
* C O N V E R T F R O M I E E E E X T E N D E D
*/
/*
* 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)
/****************************************************************
* Extended precision IEEE floating-point conversion routine.
****************************************************************/
double ConvertFromIeeeExtended(unsigned char *bytes)
{
double f;
int 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((double)(hiMant), expon-=31);
f += ldexp((double)(loMant), expon-=32);
}
}
if (bytes[0] & 0x80)
return -f;
else
return f;
}