htswanalysis/MACS/lib/gsl/gsl-1.11/ieee-utils/make_rep.c

   1 /* ieee-utils/make_rep.c
   2  *
   3  * Copyright (C) 1996, 1997, 1998, 1999, 2000, 2007 Brian Gough
   4  *
   5  * This program is free software; you can redistribute it and/or modify
   6  * it under the terms of the GNU General Public License as published by
   7  * the Free Software Foundation; either version 3 of the License, or (at
   8  * your option) any later version.
   9  *
  10  * This program is distributed in the hope that it will be useful, but
  11  * WITHOUT ANY WARRANTY; without even the implied warranty of
  12  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
  13  * General Public License for more details.
  14  *
  15  * You should have received a copy of the GNU General Public License
  16  * along with this program; if not, write to the Free Software
  17  * Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA.
  18  */
  19
  20 #include <config.h>
  21 #include <gsl/gsl_ieee_utils.h>
  22
  23 #include "endian.c"
  24 #include "standardize.c"
  25
  26 static void sprint_nybble(int i, char *s) ;
  27 static void sprint_byte(int i, char *s) ;
  28 static int determine_ieee_type (int non_zero, int exponent, int max_exponent);
  29
  30
  31 /* For the IEEE float format the bits are found from the following
  32    masks,
  33
  34    sign      = 0x80000000
  35    exponent  = 0x7f800000
  36    mantisssa = 0x007fffff
  37
  38    For the IEEE double format the masks are,
  39
  40    sign      = 0x8000000000000000
  41    exponent  = 0x7ff0000000000000
  42    mantissa  = 0x000fffffffffffff
  43
  44    */
  45
  46 void
  47 gsl_ieee_float_to_rep (const float * x, gsl_ieee_float_rep * r)
  48 {
  49   int e, non_zero;
  50
  51   union {
  52     float f;
  53     struct  {
  54       unsigned char byte[4] ;
  55     } ieee ;
  56   } u;
  57
  58   u.f = *x ;
  59
  60   if (little_endian_p())
  61     make_float_bigendian(&(u.f)) ;
  62
  63   /* note that r->sign is signed, u.ieee.byte is unsigned */
  64
  65   if (u.ieee.byte[3]>>7)
  66     {
  67       r->sign = 1 ;
  68     }
  69   else
  70     {
  71       r->sign = 0 ;
  72     }
  73
  74   e = (u.ieee.byte[3] & 0x7f) << 1 | (u.ieee.byte[2] & 0x80)>>7 ;
  75
  76   r->exponent = e - 127 ;
  77
  78   sprint_byte((u.ieee.byte[2] & 0x7f) << 1,r->mantissa) ;
  79   sprint_byte(u.ieee.byte[1],r->mantissa + 7) ;
  80   sprint_byte(u.ieee.byte[0],r->mantissa + 15) ;
  81
  82   r->mantissa[23] = '\0' ;
  83
  84   non_zero = u.ieee.byte[0] || u.ieee.byte[1] || (u.ieee.byte[2] & 0x7f);
  85
  86   r->type = determine_ieee_type (non_zero, e, 255) ;
  87 }
  88
  89 void
  90 gsl_ieee_double_to_rep (const double * x, gsl_ieee_double_rep * r)
  91 {
  92
  93   int e, non_zero;
  94
  95   union
  96   {
  97     double d;
  98     struct  {
  99       unsigned char byte[8];
 100     } ieee ;
 101   } u;
 102
 103   u.d= *x ;
 104
 105   if (little_endian_p())
 106     make_double_bigendian(&(u.d)) ;
 107
 108   /* note that r->sign is signed, u.ieee.byte is unsigned */
 109
 110   if (u.ieee.byte[7]>>7)
 111     {
 112       r->sign = 1 ;
 113     }
 114   else
 115     {
 116       r->sign = 0 ;
 117     }
 118
 119
 120   e =(u.ieee.byte[7] & 0x7f)<<4 ^ (u.ieee.byte[6] & 0xf0)>>4 ;
 121
 122   r->exponent = e - 1023 ;
 123
 124   sprint_nybble(u.ieee.byte[6],r->mantissa) ;
 125   sprint_byte(u.ieee.byte[5],r->mantissa + 4) ;
 126   sprint_byte(u.ieee.byte[4],r->mantissa + 12) ;
 127   sprint_byte(u.ieee.byte[3],r->mantissa + 20) ;
 128   sprint_byte(u.ieee.byte[2],r->mantissa + 28) ;
 129   sprint_byte(u.ieee.byte[1],r->mantissa + 36) ;
 130   sprint_byte(u.ieee.byte[0],r->mantissa + 44) ;
 131
 132   r->mantissa[52] = '\0' ;
 133
 134   non_zero = (u.ieee.byte[0] || u.ieee.byte[1] || u.ieee.byte[2]
 135               || u.ieee.byte[3] || u.ieee.byte[4] || u.ieee.byte[5]
 136               || (u.ieee.byte[6] & 0x0f)) ;
 137
 138   r->type = determine_ieee_type (non_zero, e, 2047) ;
 139 }
 140
 141 /* A table of character representations of nybbles */
 142
 143 static char nybble[16][5]={ /* include space for the \0 */
 144   "0000", "0001", "0010", "0011",
 145   "0100", "0101", "0110", "0111",
 146   "1000", "1001", "1010", "1011",
 147   "1100", "1101", "1110", "1111"
 148 }  ;
 149
 150 static void
 151 sprint_nybble(int i, char *s)
 152 {
 153   char *c ;
 154   c=nybble[i & 0x0f ];
 155   *s=c[0] ;  *(s+1)=c[1] ;  *(s+2)=c[2] ;  *(s+3)=c[3] ;
 156 }
 157
 158 static void
 159 sprint_byte(int i, char *s)
 160 {
 161   char *c ;
 162   c=nybble[(i & 0xf0)>>4];
 163   *s=c[0] ;  *(s+1)=c[1] ;  *(s+2)=c[2] ;  *(s+3)=c[3] ;
 164   c=nybble[i & 0x0f];
 165   *(s+4)=c[0] ;  *(s+5)=c[1] ;  *(s+6)=c[2] ;  *(s+7)=c[3] ;
 166 }
 167
 168 static int
 169 determine_ieee_type (int non_zero, int exponent, int max_exponent)
 170 {
 171   if (exponent == max_exponent)
 172     {
 173       if (non_zero)
 174         {
 175           return GSL_IEEE_TYPE_NAN ;
 176         }
 177       else
 178         {
 179           return GSL_IEEE_TYPE_INF ;
 180         }
 181     }
 182   else if (exponent == 0)
 183     {
 184       if (non_zero)
 185         {
 186           return GSL_IEEE_TYPE_DENORMAL ;
 187         }
 188       else
 189         {
 190           return GSL_IEEE_TYPE_ZERO ;
 191         }
 192     }
 193   else
 194     {
 195       return GSL_IEEE_TYPE_NORMAL ;
 196     }
 197 }