float.c revision c74663799493f2b1e6123c18def94295d0afab7b
1/* libFLAC - Free Lossless Audio Codec library 2 * Copyright (C) 2004,2005,2006,2007 Josh Coalson 3 * 4 * Redistribution and use in source and binary forms, with or without 5 * modification, are permitted provided that the following conditions 6 * are met: 7 * 8 * - Redistributions of source code must retain the above copyright 9 * notice, this list of conditions and the following disclaimer. 10 * 11 * - Redistributions in binary form must reproduce the above copyright 12 * notice, this list of conditions and the following disclaimer in the 13 * documentation and/or other materials provided with the distribution. 14 * 15 * - Neither the name of the Xiph.org Foundation nor the names of its 16 * contributors may be used to endorse or promote products derived from 17 * this software without specific prior written permission. 18 * 19 * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS 20 * ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT 21 * LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR 22 * A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE FOUNDATION OR 23 * CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, 24 * EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, 25 * PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR 26 * PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF 27 * LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING 28 * NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS 29 * SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. 30 */ 31 32#if HAVE_CONFIG_H 33# include <config.h> 34#endif 35 36#include "FLAC/assert.h" 37 38#include "private/float.h" 39 40#ifdef FLAC__INTEGER_ONLY_LIBRARY 41 42/* adjust for compilers that can't understand using LLU suffix for uint64_t literals */ 43#ifdef _MSC_VER 44#define FLAC__U64L(x) x 45#else 46#define FLAC__U64L(x) x##LLU 47#endif 48 49const FLAC__fixedpoint FLAC__FP_ZERO = 0; 50const FLAC__fixedpoint FLAC__FP_ONE_HALF = 0x00008000; 51const FLAC__fixedpoint FLAC__FP_ONE = 0x00010000; 52const FLAC__fixedpoint FLAC__FP_LN2 = 45426; 53const FLAC__fixedpoint FLAC__FP_E = 178145; 54 55/* Lookup tables for Knuth's logarithm algorithm */ 56#define LOG2_LOOKUP_PRECISION 16 57static const FLAC__uint32 log2_lookup[][LOG2_LOOKUP_PRECISION] = { 58 { 59 /* 60 * 0 fraction bits 61 */ 62 /* undefined */ 0x00000000, 63 /* lg(2/1) = */ 0x00000001, 64 /* lg(4/3) = */ 0x00000000, 65 /* lg(8/7) = */ 0x00000000, 66 /* lg(16/15) = */ 0x00000000, 67 /* lg(32/31) = */ 0x00000000, 68 /* lg(64/63) = */ 0x00000000, 69 /* lg(128/127) = */ 0x00000000, 70 /* lg(256/255) = */ 0x00000000, 71 /* lg(512/511) = */ 0x00000000, 72 /* lg(1024/1023) = */ 0x00000000, 73 /* lg(2048/2047) = */ 0x00000000, 74 /* lg(4096/4095) = */ 0x00000000, 75 /* lg(8192/8191) = */ 0x00000000, 76 /* lg(16384/16383) = */ 0x00000000, 77 /* lg(32768/32767) = */ 0x00000000 78 }, 79 { 80 /* 81 * 4 fraction bits 82 */ 83 /* undefined */ 0x00000000, 84 /* lg(2/1) = */ 0x00000010, 85 /* lg(4/3) = */ 0x00000007, 86 /* lg(8/7) = */ 0x00000003, 87 /* lg(16/15) = */ 0x00000001, 88 /* lg(32/31) = */ 0x00000001, 89 /* lg(64/63) = */ 0x00000000, 90 /* lg(128/127) = */ 0x00000000, 91 /* lg(256/255) = */ 0x00000000, 92 /* lg(512/511) = */ 0x00000000, 93 /* lg(1024/1023) = */ 0x00000000, 94 /* lg(2048/2047) = */ 0x00000000, 95 /* lg(4096/4095) = */ 0x00000000, 96 /* lg(8192/8191) = */ 0x00000000, 97 /* lg(16384/16383) = */ 0x00000000, 98 /* lg(32768/32767) = */ 0x00000000 99 }, 100 { 101 /* 102 * 8 fraction bits 103 */ 104 /* undefined */ 0x00000000, 105 /* lg(2/1) = */ 0x00000100, 106 /* lg(4/3) = */ 0x0000006a, 107 /* lg(8/7) = */ 0x00000031, 108 /* lg(16/15) = */ 0x00000018, 109 /* lg(32/31) = */ 0x0000000c, 110 /* lg(64/63) = */ 0x00000006, 111 /* lg(128/127) = */ 0x00000003, 112 /* lg(256/255) = */ 0x00000001, 113 /* lg(512/511) = */ 0x00000001, 114 /* lg(1024/1023) = */ 0x00000000, 115 /* lg(2048/2047) = */ 0x00000000, 116 /* lg(4096/4095) = */ 0x00000000, 117 /* lg(8192/8191) = */ 0x00000000, 118 /* lg(16384/16383) = */ 0x00000000, 119 /* lg(32768/32767) = */ 0x00000000 120 }, 121 { 122 /* 123 * 12 fraction bits 124 */ 125 /* undefined */ 0x00000000, 126 /* lg(2/1) = */ 0x00001000, 127 /* lg(4/3) = */ 0x000006a4, 128 /* lg(8/7) = */ 0x00000315, 129 /* lg(16/15) = */ 0x0000017d, 130 /* lg(32/31) = */ 0x000000bc, 131 /* lg(64/63) = */ 0x0000005d, 132 /* lg(128/127) = */ 0x0000002e, 133 /* lg(256/255) = */ 0x00000017, 134 /* lg(512/511) = */ 0x0000000c, 135 /* lg(1024/1023) = */ 0x00000006, 136 /* lg(2048/2047) = */ 0x00000003, 137 /* lg(4096/4095) = */ 0x00000001, 138 /* lg(8192/8191) = */ 0x00000001, 139 /* lg(16384/16383) = */ 0x00000000, 140 /* lg(32768/32767) = */ 0x00000000 141 }, 142 { 143 /* 144 * 16 fraction bits 145 */ 146 /* undefined */ 0x00000000, 147 /* lg(2/1) = */ 0x00010000, 148 /* lg(4/3) = */ 0x00006a40, 149 /* lg(8/7) = */ 0x00003151, 150 /* lg(16/15) = */ 0x000017d6, 151 /* lg(32/31) = */ 0x00000bba, 152 /* lg(64/63) = */ 0x000005d1, 153 /* lg(128/127) = */ 0x000002e6, 154 /* lg(256/255) = */ 0x00000172, 155 /* lg(512/511) = */ 0x000000b9, 156 /* lg(1024/1023) = */ 0x0000005c, 157 /* lg(2048/2047) = */ 0x0000002e, 158 /* lg(4096/4095) = */ 0x00000017, 159 /* lg(8192/8191) = */ 0x0000000c, 160 /* lg(16384/16383) = */ 0x00000006, 161 /* lg(32768/32767) = */ 0x00000003 162 }, 163 { 164 /* 165 * 20 fraction bits 166 */ 167 /* undefined */ 0x00000000, 168 /* lg(2/1) = */ 0x00100000, 169 /* lg(4/3) = */ 0x0006a3fe, 170 /* lg(8/7) = */ 0x00031513, 171 /* lg(16/15) = */ 0x00017d60, 172 /* lg(32/31) = */ 0x0000bb9d, 173 /* lg(64/63) = */ 0x00005d10, 174 /* lg(128/127) = */ 0x00002e59, 175 /* lg(256/255) = */ 0x00001721, 176 /* lg(512/511) = */ 0x00000b8e, 177 /* lg(1024/1023) = */ 0x000005c6, 178 /* lg(2048/2047) = */ 0x000002e3, 179 /* lg(4096/4095) = */ 0x00000171, 180 /* lg(8192/8191) = */ 0x000000b9, 181 /* lg(16384/16383) = */ 0x0000005c, 182 /* lg(32768/32767) = */ 0x0000002e 183 }, 184 { 185 /* 186 * 24 fraction bits 187 */ 188 /* undefined */ 0x00000000, 189 /* lg(2/1) = */ 0x01000000, 190 /* lg(4/3) = */ 0x006a3fe6, 191 /* lg(8/7) = */ 0x00315130, 192 /* lg(16/15) = */ 0x0017d605, 193 /* lg(32/31) = */ 0x000bb9ca, 194 /* lg(64/63) = */ 0x0005d0fc, 195 /* lg(128/127) = */ 0x0002e58f, 196 /* lg(256/255) = */ 0x0001720e, 197 /* lg(512/511) = */ 0x0000b8d8, 198 /* lg(1024/1023) = */ 0x00005c61, 199 /* lg(2048/2047) = */ 0x00002e2d, 200 /* lg(4096/4095) = */ 0x00001716, 201 /* lg(8192/8191) = */ 0x00000b8b, 202 /* lg(16384/16383) = */ 0x000005c5, 203 /* lg(32768/32767) = */ 0x000002e3 204 }, 205 { 206 /* 207 * 28 fraction bits 208 */ 209 /* undefined */ 0x00000000, 210 /* lg(2/1) = */ 0x10000000, 211 /* lg(4/3) = */ 0x06a3fe5c, 212 /* lg(8/7) = */ 0x03151301, 213 /* lg(16/15) = */ 0x017d6049, 214 /* lg(32/31) = */ 0x00bb9ca6, 215 /* lg(64/63) = */ 0x005d0fba, 216 /* lg(128/127) = */ 0x002e58f7, 217 /* lg(256/255) = */ 0x001720da, 218 /* lg(512/511) = */ 0x000b8d87, 219 /* lg(1024/1023) = */ 0x0005c60b, 220 /* lg(2048/2047) = */ 0x0002e2d7, 221 /* lg(4096/4095) = */ 0x00017160, 222 /* lg(8192/8191) = */ 0x0000b8ad, 223 /* lg(16384/16383) = */ 0x00005c56, 224 /* lg(32768/32767) = */ 0x00002e2b 225 } 226}; 227 228#if 0 229static const FLAC__uint64 log2_lookup_wide[] = { 230 { 231 /* 232 * 32 fraction bits 233 */ 234 /* undefined */ 0x00000000, 235 /* lg(2/1) = */ FLAC__U64L(0x100000000), 236 /* lg(4/3) = */ FLAC__U64L(0x6a3fe5c6), 237 /* lg(8/7) = */ FLAC__U64L(0x31513015), 238 /* lg(16/15) = */ FLAC__U64L(0x17d60497), 239 /* lg(32/31) = */ FLAC__U64L(0x0bb9ca65), 240 /* lg(64/63) = */ FLAC__U64L(0x05d0fba2), 241 /* lg(128/127) = */ FLAC__U64L(0x02e58f74), 242 /* lg(256/255) = */ FLAC__U64L(0x01720d9c), 243 /* lg(512/511) = */ FLAC__U64L(0x00b8d875), 244 /* lg(1024/1023) = */ FLAC__U64L(0x005c60aa), 245 /* lg(2048/2047) = */ FLAC__U64L(0x002e2d72), 246 /* lg(4096/4095) = */ FLAC__U64L(0x00171600), 247 /* lg(8192/8191) = */ FLAC__U64L(0x000b8ad2), 248 /* lg(16384/16383) = */ FLAC__U64L(0x0005c55d), 249 /* lg(32768/32767) = */ FLAC__U64L(0x0002e2ac) 250 }, 251 { 252 /* 253 * 48 fraction bits 254 */ 255 /* undefined */ 0x00000000, 256 /* lg(2/1) = */ FLAC__U64L(0x1000000000000), 257 /* lg(4/3) = */ FLAC__U64L(0x6a3fe5c60429), 258 /* lg(8/7) = */ FLAC__U64L(0x315130157f7a), 259 /* lg(16/15) = */ FLAC__U64L(0x17d60496cfbb), 260 /* lg(32/31) = */ FLAC__U64L(0xbb9ca64ecac), 261 /* lg(64/63) = */ FLAC__U64L(0x5d0fba187cd), 262 /* lg(128/127) = */ FLAC__U64L(0x2e58f7441ee), 263 /* lg(256/255) = */ FLAC__U64L(0x1720d9c06a8), 264 /* lg(512/511) = */ FLAC__U64L(0xb8d8752173), 265 /* lg(1024/1023) = */ FLAC__U64L(0x5c60aa252e), 266 /* lg(2048/2047) = */ FLAC__U64L(0x2e2d71b0d8), 267 /* lg(4096/4095) = */ FLAC__U64L(0x1716001719), 268 /* lg(8192/8191) = */ FLAC__U64L(0xb8ad1de1b), 269 /* lg(16384/16383) = */ FLAC__U64L(0x5c55d640d), 270 /* lg(32768/32767) = */ FLAC__U64L(0x2e2abcf52) 271 } 272}; 273#endif 274 275FLAC__uint32 FLAC__fixedpoint_log2(FLAC__uint32 x, unsigned fracbits, unsigned precision) 276{ 277 const FLAC__uint32 ONE = (1u << fracbits); 278 const FLAC__uint32 *table = log2_lookup[fracbits >> 2]; 279 280 FLAC__ASSERT(fracbits < 32); 281 FLAC__ASSERT((fracbits & 0x3) == 0); 282 283 if(x < ONE) 284 return 0; 285 286 if(precision > LOG2_LOOKUP_PRECISION) 287 precision = LOG2_LOOKUP_PRECISION; 288 289 /* Knuth's algorithm for computing logarithms, optimized for base-2 with lookup tables */ 290 { 291 FLAC__uint32 y = 0; 292 FLAC__uint32 z = x >> 1, k = 1; 293 while (x > ONE && k < precision) { 294 if (x - z >= ONE) { 295 x -= z; 296 z = x >> k; 297 y += table[k]; 298 } 299 else { 300 z >>= 1; 301 k++; 302 } 303 } 304 return y; 305 } 306} 307 308#endif /* defined FLAC__INTEGER_ONLY_LIBRARY */ 309