ale_real.h

   1 // Copyright 2002, 2004 David Hilvert <dhilvert@auricle.dyndns.org>,
   2 //                                    <dhilvert@ugcs.caltech.edu>
   3
   4 /*  This file is part of the Anti-Lamenessing Engine.
   5
   6     The Anti-Lamenessing Engine is free software; you can redistribute it and/or modify
   7     it under the terms of the GNU General Public License as published by
   8     the Free Software Foundation; either version 2 of the License, or
   9     (at your option) any later version.
  10
  11     The Anti-Lamenessing Engine is distributed in the hope that it will be useful,
  12     but WITHOUT ANY WARRANTY; without even the implied warranty of
  13     MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
  14     GNU General Public License for more details.
  15
  16     You should have received a copy of the GNU General Public License
  17     along with the Anti-Lamenessing Engine; if not, write to the Free Software
  18     Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
  19 */
  20
  21 #ifndef __ale_real_h__
  22 #define __ale_real_h__
  23
  24 #define SINGLE 1
  25 #define DOUBLE 2
  26 #define HALF 3
  27
  28 /*
  29  * Real-valued type used to represent the range of an image (colors, weights,
  30  * etc.).
  31  *
  32  * ale_real is used in computation.
  33  * ale_sreal is used for storage.
  34  */
  35
  36 #if ALE_COLORS == SINGLE
  37
  38 typedef float ale_real;
  39 typedef float ale_sreal;
  40
  41 #define ALE_REAL_PRECISION_STRING "SINGLE"
  42
  43 #elif ALE_COLORS == DOUBLE
  44
  45 typedef double ale_real;
  46 typedef double ale_sreal;
  47
  48 #define ALE_REAL_PRECISION_STRING "DOUBLE"
  49
  50 #elif ALE_COLORS == HALF
  51
  52 /*
  53  * What follows is one approach to packing a floating point
  54  * number into 16 bits.  This implementation is very slow.
  55  */
  56
  57 #define MANTISSA_BITS (9)
  58 #define EXP_BITS (15 - MANTISSA_BITS)
  59 #define EXP_SPECIAL (1 << (EXP_BITS - 1))
  60 #define EXP_MAX (EXP_SPECIAL - 1)
  61 #define EXP_MIN (-EXP_MAX)
  62
  63 typedef float ale_real;
  64
  65 class ale_sreal {
  66
  67         union {
  68                 uint16_t bits;
  69                 struct {
  70                         uint16_t sign:1;
  71                         uint16_t mant:MANTISSA_BITS;
  72                         int16_t  exp :EXP_BITS;
  73                 } fpr;
  74         } u;
  75 public:
  76         ale_sreal() {
  77                 u.bits = 0;
  78         }
  79         ale_sreal operator=(float v) {
  80                 if (v == 0) {
  81                         u.bits = 0;
  82                 } else if (isnan(v)) {
  83                         u.fpr.exp = EXP_SPECIAL;
  84                         u.fpr.mant = 1;
  85                 } else {
  86
  87                         if (v > 0)
  88                                 u.fpr.sign = 0;
  89                         else if (v < 0) {
  90                                 u.fpr.sign = 1;
  91                                 v = -v;
  92                         } else
  93                                 assert(0);
  94
  95                         /*
  96                          * Get the exponent.
  97                          */
  98
  99                         int log2 = (int) floor (log(v) / log(2));
 100
 101                         /*
 102                          * Test the exponent against the largest expressible
 103                          * exponent for ale_sreal.
 104                          */
 105
 106                         if (log2 > EXP_MAX) {
 107                                 /*
 108                                  * Infinity
 109                                  */
 110
 111                                 u.fpr.exp = EXP_SPECIAL;
 112                                 u.fpr.mant = 0;
 113
 114                                 return *this;
 115                         }
 116
 117                         /*
 118                          * Test the exponent against the smallest expressible
 119                          * exponent for ale_sreal.
 120                          */
 121
 122                         if (log2 < EXP_MIN) {
 123                                 /*
 124                                  * Zero
 125                                  */
 126
 127                                 u.fpr.exp = 0x0;
 128                                 u.fpr.mant = 0;
 129
 130                                 return *this;
 131                         }
 132
 133                         /*
 134                          * The exponent is in range, so use it.
 135                          */
 136
 137                         u.fpr.exp = log2;
 138
 139                         u.fpr.mant = (uint16_t) floor(v / pow(2, log2) * (1 << (MANTISSA_BITS - 1)));
 140                 }
 141
 142                 return *this;
 143         }
 144
 145         operator float() const {
 146                 float result = 3.14159;
 147
 148                 if (((uint16_t) u.fpr.exp == EXP_SPECIAL) && (u.fpr.mant == 1)) {
 149
 150                         /*
 151                          * NaN
 152                          */
 153
 154                         float a = 0;
 155                         float b = 0;
 156
 157                         result = a / b;
 158
 159                 } else if (((uint16_t) u.fpr.exp == EXP_SPECIAL) && (u.fpr.mant == 0)) {
 160
 161                         /*
 162                          * Infinity
 163                          */
 164
 165                         float a = 1;
 166                         float b = 0;
 167
 168                         result = (a / b);
 169
 170                 } else if ((uint16_t) u.fpr.exp != EXP_SPECIAL) {
 171
 172                         /*
 173                          * Value is finite.
 174                          */
 175
 176                         result = u.fpr.mant / ((double) (1 << (MANTISSA_BITS - 1)))
 177                                * pow(2, u.fpr.exp);
 178
 179                 } else
 180                         assert(0);
 181
 182                 if (u.fpr.sign)
 183                         result = -result;
 184
 185                 return result;
 186         }
 187
 188         ale_sreal operator-=(float r) {
 189                 *this = (float) *this - (float) r;
 190                 return *this;
 191         }
 192         ale_sreal operator/=(float r) {
 193                 *this = (float) *this / (float) r;
 194                 return *this;
 195         }
 196         ale_sreal operator*=(float r) {
 197                 *this = (float) *this * (float) r;
 198                 return *this;
 199         }
 200         ale_sreal operator+=(float r) {
 201                 *this = (float) *this + (float) r;
 202                 return *this;
 203         }
 204 };
 205
 206 #undef MANTISSA_BITS
 207 #undef EXP_BITS
 208 #undef EXP_SPECIAL
 209 #undef EXP_MAX
 210 #undef EXP_MIN
 211
 212 #define ALE_REAL_PRECISION_STRING "HALF"
 213
 214 #else
 215
 216 #warning Unknown precision in ale_real.h: Choosing PRECISION=SINGLE.
 217
 218 typedef float ale_real;
 219 typedef float ale_sreal;
 220
 221 #define ALE_REAL_PRECISION_STRING "SINGLE"
 222
 223 #endif
 224
 225 #undef SINGLE
 226 #undef DOUBLE
 227 #undef HALF
 228
 229 #endif