/* @(#)cvb.c 17.1.1.1 (ES0-DMD) 01/25/02 17:39:08 */ /*=========================================================================== Copyright (C) 1995 European Southern Observatory (ESO) This program is free software; you can redistribute it and/or modify it under the terms of the GNU General Public License as published by the Free Software Foundation; either version 2 of the License, or (at your option) any later version. This program is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License for more details. You should have received a copy of the GNU General Public License along with this program; if not, write to the Free Software Foundation, Inc., 675 Massachusetss Ave, Cambridge, MA 02139, USA. Corresponding concerning ESO-MIDAS should be addressed as follows: Internet e-mail: midas@eso.org Postal address: European Southern Observatory Data Management Division Karl-Schwarzschild-Strasse 2 D 85748 Garching bei Muenchen GERMANY ===========================================================================*/ /*++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ .COPYRIGHT (c) 1995 European Southern Observatory .IDENTIFIC cvb.c .LANGUAGE C .AUTHOR P. Grosbol ESO/IPG .KEYWORDS binary number conversion, IEEE-754, floating point .PURPOSE convert array of binary numbers from standard FITS format to local binary representation .VERSION 1.0 1988-Nov-06 : Creation, PJG .VERSION 1.1 1989-Oct-10 : Byte swap for SP-FP, PJG .VERSION 2.0 1990-Feb-04 : Change calling sequence, PJG .VERSION 2.1 1990-Feb-12 : Change byte-swap algorithm, PJG .VERSION 2.2 1990-Mar-30 : Correct union pointers, PJG .VERSION 2.3 1990-Sep-28 : Correct IEEE -> VAXG, PJG .VERSION 2.4 1991-Feb-12 : Include fp-swap, PJG .VERSION 2.5 1991-Mar-22 : Remove a register variable 'v', PJG .VERSION 2.6 1993-Sep-22 : Correct IEEE NULL check + bswap, PJG .VERSION 2.7 1993-Oct-28 : Change to simple types in calls, PJG .VERSION 2.8 1994-Mar-03 : Include math.h and minor changes, PJG .VERSION 2.9 1994-May-09 : Explicit cast for ANSI-C, PJG .VERSION 3.0 1994-Nov-14 : Check range of VAXFLOAT, PJG .VERSION 3.1 1994-Dec-06 : Initiate ifpx in cvr8, PJG .VERSION 3.2 1995-Feb-02 : Correct ifpx initiation in cvr8, PJG ---------------------------------------------------------------------*/ #include /* Math-library definitions */ #include /* define types of binary data format */ #include /* define data conversion formats */ static int samei2; /* same 2-byte integer format */ static int samei4; /* same 4-byte integer format */ static DFMT efmt; /* external data format definition */ static int ls0,ls1,ls2,ls3; /* 32-bit integer byte swap order */ static int fs0,fs1,fs2,fs3; /* 32-bit float byte swap order */ static int ds0,ds1,ds2,ds3; /* 64-bit double byte swap order */ static int ds4,ds5,ds6,ds7; /* 64-bit double byte swap order */ int cvinit(extfmt) /*++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ .PURPOSE initiate data conversion routines .RETURN error code: ok=0, 1=unknown format .ALGORITM Analyze internal and external data format definitions and setup static varibles to define conversion needed. --------------------------------------------------------------------*/ int extfmt; /* External data format FITS/IHAP */ { int i, lbo[8], xbo[8], sa[8]; long n; switch (extfmt) { case FITS : efmt.ifmt = TWOS_COMP; efmt.fpfmt = IEEEFLOAT; efmt.bos = 12; efmt.bol = 1234; efmt.bof = 1234; efmt.bod = 12345678; break; case IHAP : efmt.ifmt = TWOS_COMP; efmt.fpfmt = HPFLOAT; efmt.bos = 12; efmt.bol = 3412; efmt.bof = 1234; efmt.bod = 12345678; break; default : return 1; } samei2 = ((cpu.ifmt == efmt.ifmt) && (cpu.bos == efmt.bos)); samei4 = ((cpu.ifmt == efmt.ifmt) && (cpu.bol == efmt.bol)); /* Setup byte swap pattern for long int */ n = cpu.bol; for (i=0; i<4; i++) { lbo[3-i] = n%10; n /= 10; } n = efmt.bol; for (i=0; i<4; i++) { xbo[3-i] = n%10; n /= 10; } for (i=0; i<4; i++) for (n=0; n<4; n++) if (lbo[i]==xbo[n]) sa[i] = n; ls0 = sa[0]; ls1 = sa[1]; ls2 = sa[2]; ls3 = sa[3]; /* Setup byte swap pattern for float */ n = cpu.bof; for (i=0; i<4; i++) { lbo[3-i] = n%10; n /= 10; } n = efmt.bof; for (i=0; i<4; i++) { xbo[3-i] = n%10; n /= 10; } for (i=0; i<4; i++) for (n=0; n<4; n++) if (lbo[i]==xbo[n]) sa[i] = n; fs0 = sa[0]; fs1 = sa[1]; fs2 = sa[2]; fs3 = sa[3]; /* Setup byte swap pattern for double */ n = cpu.bod; for (i=0; i<8; i++) { lbo[7-i] = n%10; n /= 10; } n = efmt.bod; for (i=0; i<8; i++) { xbo[7-i] = n%10; n /= 10; } for (i=0; i<8; i++) for (n=0; n<8; n++) if (lbo[i]==xbo[n]) sa[i] = n; ds0 = sa[0]; ds1 = sa[1]; ds2 = sa[2]; ds3 = sa[3]; ds4 = sa[4]; ds5 = sa[5]; ds6 = sa[6]; ds7 = sa[7]; return 0; } typedef union { /* union for conversion */ UINT1 c[2]; /* bytes */ INT2 s; /* 2 byte integer */ } VI2; int cvi2(pbuf,no,to) /*++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ .PURPOSE convert 2 byte integer array between different computers .RETURN error code: ok=0, 1=cannot convert .ALGORITM swap bytes between VAX - nonVAX machines ---------------------------------------------------------------------*/ INT2 *pbuf; /* pointer to data array */ int no; /* no. of values to convert */ int to; /* true if convert to ext.fmt */ { register UINT1 byte; register int n; register VI2 *pv; if (no<1 || samei2) return 0; /* check if conversion needed */ if (cpu.ifmt!=efmt.ifmt) return 1; /* no format conversion */ if (cpu.bos!=efmt.bos) { /* byte swap needed ! */ n = no; pv = (VI2 *) pbuf; while (n--) { /* loop through data array */ byte = pv->c[0]; pv->c[0] = pv->c[1]; pv->c[1] = byte; pv++; } } return 0; } typedef union { /* union for conversion */ UINT1 c[4]; /* bytes */ INT4 i; /* 4 byte integer */ } VI4; int cvi4(pbuf,no,to) /*++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ .PURPOSE convert 4 byte integer array between different computers .RETURN error code: ok=0, 1=cannot convert .ALGORITM swap bytes and words between VAX - nonVAX machines ---------------------------------------------------------------------*/ INT4 *pbuf; /* pointer to data array */ int no; /* no. of values to convert */ int to; /* true if convert to ext.fmt */ { register UINT1 *pb; register VI4 *pv; register int n; VI4 v; if (no<1 || samei4) return 0; /* check if conversion needed */ if (cpu.ifmt!=efmt.ifmt) return 1; /* no format conversion */ if (cpu.bol!=efmt.bol) { /* byte swap needed ! */ n = no; pv = (VI4 *) pbuf; pb = (UINT1 *) pbuf; if (to) while (n--) { /* loop through data array */ v.c[ls0] = *pb++; v.c[ls1] = *pb++; v.c[ls2] = *pb++; v.c[ls3] = *pb++; (pv++)->i = v.i; } else while (n--) { /* loop through data array */ v.i = (pv++)->i; *pb++ = v.c[ls0]; *pb++ = v.c[ls1]; *pb++ = v.c[ls2]; *pb++ = v.c[ls3]; } } return 0; } typedef union { /* union for conversion */ UINT1 c[4]; /* bytes */ INT4 i; /* 4 byte integer */ REAL4 f; /* 4 byte floating point */ } VR4; #ifdef __STDC__ int cvr4(REAL4 * pbuf , int no , int to) #else int cvr4(pbuf, no, to) /*++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ .PURPOSE convert 4 byte real array between different computers .RETURN error code: ok=0, 1=cannot convert .ALGORITM swap bytes between VAX - nonVAX machines in addition some bit change due to different FP def's. ---------------------------------------------------------------------*/ REAL4 *pbuf; /* pointer to data array */ int no; /* no. of values to convert */ int to; /* true if convert to ext.fmt */ #endif { register UINT1 *pb, byte; register VR4 *pv; register int n,i; VR4 v; int ifpx, bswap; double d; if (no<1) return 0; /* check if conversion needed */ bswap = (cpu.bof != efmt.bof); if (to) { /* convert to external format */ switch (cpu.fpfmt) { /* check internal format */ case IEEEFLOAT : switch (efmt.fpfmt) { case IEEEFLOAT : n = no; pv = (VR4 *) pbuf; while (n--) { if (isNULLFLOAT(pv->f)) pv->i = 0xFFFFFFFFL; pv++; } break; case VAXFLOAT : case VAXGFLOAT : default : return 1; } break; case VAXFLOAT : case VAXGFLOAT : switch (efmt.fpfmt) { case IEEEFLOAT : n = no; ifpx = 0; pv = (VR4 *) pbuf; while (n--) { if (isNULLFLOAT(pv->f)) pv->i = 0xFFFFFFFFL; else { d = frexp((double) pv->f, &ifpx); if (ifpx<-125 || !(pv->i & 0x7F80)) pv->f = 0.0; else pv->f *= 0.25; } pv++; } break; case VAXFLOAT : case VAXGFLOAT : break; default : return 1; } break; default : return 1; } if (bswap) { /* check if byte swap machine */ n = no; pv = (VR4 *) pbuf; pb = (UINT1 *) pbuf; while (n--) { /* loop through data array */ v.c[fs0] = *pb++; v.c[fs1] = *pb++; v.c[fs2] = *pb++; v.c[fs3] = *pb++; (pv++)->i = v.i; } } } else { /* copy from external to internal format */ switch (efmt.fpfmt) { case IEEEFLOAT : switch (cpu.fpfmt) { case IEEEFLOAT : n = no; pv = (VR4 *) pbuf; pb = (UINT1 *) pbuf; while (n--) { if ((pv->i & FPXM) == FPXM) { /* check NaN + Inf */ toNULLFLOAT(pv->f); pb += 4; } else if (bswap) { v.i = pv->i; *pb++ = v.c[fs0]; *pb++ = v.c[fs1]; *pb++ = v.c[fs2]; *pb++ = v.c[fs3]; } pv++; } break; case VAXFLOAT : case VAXGFLOAT : n = no; pv = (VR4 *) pbuf; pb = (UINT1 *) pbuf; while (n--) { v.i = pv->i & FPXM; byte = v.c[0]; v.c[0] = v.c[1]; v.c[1] = byte; ifpx = v.i >> 7; v.i = pv->i; *pb++ = v.c[fs0]; *pb++ = v.c[fs1]; *pb++ = v.c[fs2]; *pb++ = v.c[fs3]; if (!ifpx) pv->f = 0.0; else if (254<=ifpx) toNULLFLOAT(pv->f); else if (!ifpx) pv->f = 0.0; else pv->f *= 4.0; pv++; } break; default : return 1; } break; case HPFLOAT : switch (cpu.fpfmt) { case IEEEFLOAT : case VAXFLOAT : case VAXGFLOAT : n = no; pv = (VR4 *) pbuf; pb = (UINT1 *) pbuf; while (n--) { i = *pb++; i = (i<<8) | *pb++; i = (i<<8) | *pb++; i = (i<<8) | *pb++; ifpx = (i & 1) ? (i>>1) | ~0x7F : (i>>1) & 0x7F; pv->f = ldexp((double)(i & (int)0xFFFFFF00),ifpx-31); pv++; } break; default : return 1; } break; default : return 1; } } return 0; } typedef union { /* union for conversion */ UINT1 c[8]; /* bytes */ UINT2 s[4]; /* 2 byte integers */ INT4 i[2]; /* 4 byte integers */ REAL8 d; /* 8 byte floating point */ } VR8; int cvr8(pbuf,no,to) /*++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ .PURPOSE convert 8 byte real array between different computers .RETURN error code: ok=0, 1=cannot convert .ALGORITM swap bytes between VAX - nonVAX machines in addition change some bits due to differences in the internal FP definition. ---------------------------------------------------------------------*/ REAL8 *pbuf; /* pointer to data array */ int no; /* no. of values to convert */ int to; /* true if convert to ext.fmt */ { register UINT1 *pb, byte; register int n; register VR8 *pv; VR8 v; int ifpx, bswap; double d; if (no<1) return 0; /* check if conversion needed */ bswap = (cpu.bod != efmt.bod); if (to) { /* convert to external format */ switch (cpu.fpfmt) { /* check internal format */ case IEEEFLOAT : switch (efmt.fpfmt) { case IEEEFLOAT : n = no; pv = (VR8 *) pbuf; while (n--) { if (isNULLDOUBLE(pv->d)) pv->i[0] = 0xFFFFFFFFL, pv->i[1] = 0xFFFFFFFFL; } break; default : return 1; } break; case VAXFLOAT : switch (efmt.fpfmt) { case VAXFLOAT : break; case IEEEFLOAT : n = no; pv = (VR8 *) pbuf; while (n--) { if (isNULLDOUBLE(pv->d)) pv->i[0] = 0xFFFFFFFFL, pv->i[1] = 0xFFFFFFFFL; else { pv->s[3] >>= 3; pv->c[7] |= (pv->c[4] & 0x07) << 5; pv->s[2] >>= 3; pv->c[5] |= (pv->c[2] & 0x07) << 5; pv->s[1] >>= 3; pv->c[3] |= (pv->c[0] & 0x07) << 5; pv->s[0] >>= 3; if (pv->c[1] & 0x10) pv->s[0] = (pv->s[0] & 0x0FFF) | 0x8000; pv->s[0] += 0x37E0; } pv++; } break; default : return 1; } break; case VAXGFLOAT : switch (efmt.fpfmt) { case IEEEFLOAT : n = no; ifpx = 0; pv = (VR8 *) pbuf; while (n--) { if (isNULLDOUBLE(pv->d)) pv->i[0] = 0xFFFFFFFFL, pv->i[1] = 0xFFFFFFFFL; else { d = frexp(pv->d, &ifpx); if (ifpx<-1021 || !(pv->i[0] & 0x7FF0)) pv->d = 0.0; else pv->d *= 0.25; } pv++; } break; case VAXGFLOAT : break; default : return 1; } break; default : return 1; } if (bswap) { /* check if byte swap machine */ n = no; pv = (VR8 *) pbuf; pb = (UINT1 *) pbuf; while (n--) { v.c[ds0] = *pb++; v.c[ds1] = *pb++; v.c[ds2] = *pb++; v.c[ds3] = *pb++; v.c[ds4] = *pb++; v.c[ds5] = *pb++; v.c[ds6] = *pb++; v.c[ds7] = *pb++; pv->i[0] = v.i[0]; (pv++)->i[1] = v.i[1]; } } } else { /* copy from external to internal format */ switch (efmt.fpfmt) { case IEEEFLOAT : switch (cpu.fpfmt) { case IEEEFLOAT : n = no; pv = (VR8 *) pbuf; pb = (UINT1 *) pbuf; while (n--) { if ((pv->i[0] & DPXM) == DPXM) { /* check NaN + Inf */ toNULLDOUBLE(pv->d); pb += 8; } else if (bswap) { v.i[0] = pv->i[0]; v.i[1] = pv->i[1]; *pb++ = v.c[ds0]; *pb++ = v.c[ds1]; *pb++ = v.c[ds2]; *pb++ = v.c[ds3]; *pb++ = v.c[ds4]; *pb++ = v.c[ds5]; *pb++ = v.c[ds6]; *pb++ = v.c[ds7]; } pv++; } break; case VAXFLOAT : /* not tested - check byte order */ n = no; pv = (VR8 *) pbuf; while (n--) { v.i[0] = pv->i[0] & DPXM; byte = v.c[0]; v.c[0] = v.c[1]; v.c[1] = byte; ifpx = v.i[0] >> 4; byte = pv->c[0]; pv->c[0] = pv->c[1]; pv->c[1] = byte; byte = pv->c[2]; pv->c[2] = pv->c[3]; pv->c[3] = byte; byte = pv->c[4]; pv->c[4] = pv->c[5]; pv->c[5] = byte; byte = pv->c[6]; pv->c[6] = pv->c[7]; pv->c[7] = byte; if (ifpx<=770) pv->d = 0.0; else if (1278<=ifpx) toNULLDOUBLE(pv->d); else { byte = pv->c[1] & 0x80; pv->s[0] = ((pv->s[0] - 0x37E0) << 3) & 0x7FF8; if (byte) pv->s[0] |= 0x8000; pv->s[0] |= (pv->c[3] >> 5) & 0x07; pv->s[1] <<= 3; pv->s[1] |= (pv->c[5] >> 5) & 0x07; pv->s[2] <<= 3; pv->s[2] |= (pv->c[7] >> 5) & 0x07; pv->s[3] <<= 3; } pv++; } break; case VAXGFLOAT : n = no; pv = (VR8 *) pbuf; while (n--) { v.i[0] = pv->i[0] & DPXM; byte = v.c[0]; v.c[0] = v.c[1]; v.c[1] = byte; ifpx = v.i[0] >> 4; byte = pv->c[0]; pv->c[0] = pv->c[1]; pv->c[1] = byte; byte = pv->c[2]; pv->c[2] = pv->c[3]; pv->c[3] = byte; byte = pv->c[4]; pv->c[4] = pv->c[5]; pv->c[5] = byte; byte = pv->c[6]; pv->c[6] = pv->c[7]; pv->c[7] = byte; if (!ifpx) pv->d = 0.0; else if (2046<=ifpx) toNULLDOUBLE(pv->d); else pv->d *= 4.0; pv++; } break; default : return 1; } break; default : return 1; } } return 0; }