softfloat.h 8.42 KB
Newer Older
1 2 3
/*
 * Copyright (c) 2006 Michael Niedermayer <michaelni@gmx.at>
 *
4 5 6
 * This file is part of FFmpeg.
 *
 * FFmpeg is free software; you can redistribute it and/or
7 8
 * modify it under the terms of the GNU Lesser General Public
 * License as published by the Free Software Foundation; either
9
 * version 2.1 of the License, or (at your option) any later version.
10
 *
11
 * FFmpeg is distributed in the hope that it will be useful,
12 13 14 15 16
 * but WITHOUT ANY WARRANTY; without even the implied warranty of
 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
 * Lesser General Public License for more details.
 *
 * You should have received a copy of the GNU Lesser General Public
17
 * License along with FFmpeg; if not, write to the Free Software
18 19
 * Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA
 */
20

21 22
#ifndef AVUTIL_SOFTFLOAT_H
#define AVUTIL_SOFTFLOAT_H
23

24
#include <stdint.h>
25
#include "common.h"
26

27
#include "avassert.h"
28
#include "softfloat_tables.h"
29

30
#define MIN_EXP -149
31 32 33 34 35
#define MAX_EXP  126
#define ONE_BITS 29

typedef struct SoftFloat{
    int32_t mant;
36
    int32_t  exp;
37 38
}SoftFloat;

39 40 41 42 43 44 45
static const SoftFloat FLOAT_0          = {          0,   MIN_EXP};             ///< 0.0
static const SoftFloat FLOAT_05         = { 0x20000000,   0};                   ///< 0.5
static const SoftFloat FLOAT_1          = { 0x20000000,   1};                   ///< 1.0
static const SoftFloat FLOAT_EPSILON    = { 0x29F16B12, -16};                   ///< A small value
static const SoftFloat FLOAT_1584893192 = { 0x32B771ED,   1};                   ///< 1.584893192 (10^.2)
static const SoftFloat FLOAT_100000     = { 0x30D40000,  17};                   ///< 100000
static const SoftFloat FLOAT_0999999    = { 0x3FFFFBCE,   0};                   ///< 0.999999
46
static const SoftFloat FLOAT_MIN        = { 0x20000000,   MIN_EXP};
47

48 49 50 51

/**
 * Convert a SoftFloat to a double precision float.
 */
52 53
static inline av_const double av_sf2double(SoftFloat v) {
    v.exp -= ONE_BITS +1;
54
    return ldexp(v.mant, v.exp);
55 56
}

57
static av_const SoftFloat av_normalize_sf(SoftFloat a){
58 59
    if(a.mant){
#if 1
60
        while((a.mant + 0x1FFFFFFFU)<0x3FFFFFFFU){
61 62 63 64
            a.mant += a.mant;
            a.exp  -= 1;
        }
#else
65
        int s=ONE_BITS - av_log2(FFABS(a.mant));
66 67 68 69 70 71 72 73 74 75 76 77 78
        a.exp   -= s;
        a.mant <<= s;
#endif
        if(a.exp < MIN_EXP){
            a.exp = MIN_EXP;
            a.mant= 0;
        }
    }else{
        a.exp= MIN_EXP;
    }
    return a;
}

79
static inline av_const SoftFloat av_normalize1_sf(SoftFloat a){
80
#if 1
81
    if((int32_t)(a.mant + 0x40000000U) <= 0){
82 83 84
        a.exp++;
        a.mant>>=1;
    }
85
    av_assert2(a.mant < 0x40000000 && a.mant > -0x40000000);
86
    av_assert2(a.exp <= MAX_EXP);
87 88 89
    return a;
#elif 1
    int t= a.mant + 0x40000000 < 0;
90
    return (SoftFloat){ a.mant>>t, a.exp+t};
91
#else
92
    int t= (a.mant + 0x3FFFFFFFU)>>31;
93
    return (SoftFloat){a.mant>>t, a.exp+t};
94 95 96 97
#endif
}

/**
98
 * @return Will not be more denormalized than a*b. So if either input is
99 100
 *         normalized, then the output will not be worse then the other input.
 *         If both are normalized, then the output will be normalized.
101
 */
102
static inline av_const SoftFloat av_mul_sf(SoftFloat a, SoftFloat b){
103
    a.exp += b.exp;
104
    av_assert2((int32_t)((a.mant * (int64_t)b.mant) >> ONE_BITS) == (a.mant * (int64_t)b.mant) >> ONE_BITS);
105
    a.mant = (a.mant * (int64_t)b.mant) >> ONE_BITS;
106 107 108 109
    a = av_normalize1_sf((SoftFloat){a.mant, a.exp - 1});
    if (!a.mant || a.exp < MIN_EXP)
        return FLOAT_0;
    return a;
110 111 112
}

/**
113 114
 * b has to be normalized and not zero.
 * @return Will not be more denormalized than a.
115
 */
116
static inline av_const SoftFloat av_div_sf(SoftFloat a, SoftFloat b){
117 118
    int64_t temp = (int64_t)a.mant * (1<<(ONE_BITS+1));
    temp /= b.mant;
119
    a.exp -= b.exp;
120 121 122 123 124 125
    a.mant = temp;
    while (a.mant != temp) {
        temp /= 2;
        a.exp--;
        a.mant = temp;
    }
126 127 128 129
    a = av_normalize1_sf(a);
    if (!a.mant || a.exp < MIN_EXP)
        return FLOAT_0;
    return a;
130 131
}

132 133 134 135 136 137
/**
 * Compares two SoftFloats.
 * @returns < 0 if the first is less
 *          > 0 if the first is greater
 *            0 if they are equal
 */
138
static inline av_const int av_cmp_sf(SoftFloat a, SoftFloat b){
139
    int t= a.exp - b.exp;
140 141 142 143
    if      (t <-31) return                  -  b.mant      ;
    else if (t <  0) return (a.mant >> (-t)) -  b.mant      ;
    else if (t < 32) return  a.mant          - (b.mant >> t);
    else             return  a.mant                         ;
144 145
}

146 147 148 149
/**
 * Compares two SoftFloats.
 * @returns 1 if a is greater than b, 0 otherwise
 */
150 151 152
static inline av_const int av_gt_sf(SoftFloat a, SoftFloat b)
{
    int t= a.exp - b.exp;
153
    if      (t <-31) return 0                >  b.mant      ;
154 155
    else if (t <  0) return (a.mant >> (-t)) >  b.mant      ;
    else if (t < 32) return  a.mant          > (b.mant >> t);
156
    else             return  a.mant          >  0           ;
157 158
}

159 160 161
/**
 * @returns the sum of 2 SoftFloats.
 */
162
static inline av_const SoftFloat av_add_sf(SoftFloat a, SoftFloat b){
163
    int t= a.exp - b.exp;
164
    if      (t <-31) return b;
165 166
    else if (t <  0) return av_normalize_sf(av_normalize1_sf((SoftFloat){ b.mant + (a.mant >> (-t)), b.exp}));
    else if (t < 32) return av_normalize_sf(av_normalize1_sf((SoftFloat){ a.mant + (b.mant >>   t ), a.exp}));
167
    else             return a;
168 169
}

170 171 172
/**
 * @returns the difference of 2 SoftFloats.
 */
173
static inline av_const SoftFloat av_sub_sf(SoftFloat a, SoftFloat b){
174
    return av_add_sf(a, (SoftFloat){ -b.mant, b.exp});
175 176
}

177
//FIXME log, exp, pow
178

179
/**
180 181 182 183
 * Converts a mantisse and exponent to a SoftFloat.
 * This converts a fixed point value v with frac_bits fractional bits to a
 * SoftFloat.
 * @returns a SoftFloat with value v * 2^-frac_bits
184
 */
185
static inline av_const SoftFloat av_int2sf(int v, int frac_bits){
186
    int exp_offset = 0;
187
    if(v <= INT_MIN + 1){
188 189 190 191
        exp_offset = 1;
        v>>=1;
    }
    return av_normalize_sf(av_normalize1_sf((SoftFloat){v, ONE_BITS + 1 - frac_bits + exp_offset}));
192 193 194
}

/**
195
 * Converts a SoftFloat to an integer.
196
 * Rounding is to -inf.
197
 */
198
static inline av_const int av_sf2int(SoftFloat v, int frac_bits){
199
    v.exp += frac_bits - (ONE_BITS + 1);
200 201 202
    if(v.exp >= 0) return v.mant <<  v.exp ;
    else           return v.mant >>(-v.exp);
}
203

204 205 206 207 208 209 210 211
/**
 * Rounding-to-nearest used.
 */
static av_always_inline SoftFloat av_sqrt_sf(SoftFloat val)
{
    int tabIndex, rem;

    if (val.mant == 0)
212
        val.exp = MIN_EXP;
213
    else if (val.mant < 0)
214
        abort();
215 216 217 218 219 220 221 222 223 224 225 226 227 228 229 230 231 232 233 234 235 236 237 238 239
    else
    {
        tabIndex = (val.mant - 0x20000000) >> 20;

        rem = val.mant & 0xFFFFF;
        val.mant  = (int)(((int64_t)av_sqrttbl_sf[tabIndex] * (0x100000 - rem) +
                           (int64_t)av_sqrttbl_sf[tabIndex + 1] * rem +
                           0x80000) >> 20);
        val.mant = (int)(((int64_t)av_sqr_exp_multbl_sf[val.exp & 1] * val.mant +
                          0x10000000) >> 29);

        if (val.mant < 0x40000000)
            val.exp -= 2;
        else
            val.mant >>= 1;

        val.exp = (val.exp >> 1) + 1;
    }

    return val;
}

/**
 * Rounding-to-nearest used.
 */
240 241 242 243 244 245 246
static av_unused void av_sincos_sf(int a, int *s, int *c)
{
    int idx, sign;
    int sv, cv;
    int st, ct;

    idx = a >> 26;
247
    sign = (int32_t)((unsigned)idx << 27) >> 31;
248 249 250 251
    cv = av_costbl_1_sf[idx & 0xf];
    cv = (cv ^ sign) - sign;

    idx -= 8;
252
    sign = (int32_t)((unsigned)idx << 27) >> 31;
253 254 255 256 257 258 259 260 261 262 263 264 265 266 267 268 269 270 271 272 273 274 275 276 277 278 279 280 281 282 283 284 285 286 287
    sv = av_costbl_1_sf[idx & 0xf];
    sv = (sv ^ sign) - sign;

    idx = a >> 21;
    ct = av_costbl_2_sf[idx & 0x1f];
    st = av_sintbl_2_sf[idx & 0x1f];

    idx = (int)(((int64_t)cv * ct - (int64_t)sv * st + 0x20000000) >> 30);

    sv = (int)(((int64_t)cv * st + (int64_t)sv * ct + 0x20000000) >> 30);

    cv = idx;

    idx = a >> 16;
    ct = av_costbl_3_sf[idx & 0x1f];
    st = av_sintbl_3_sf[idx & 0x1f];

    idx = (int)(((int64_t)cv * ct - (int64_t)sv * st + 0x20000000) >> 30);

    sv = (int)(((int64_t)cv * st + (int64_t)sv * ct + 0x20000000) >> 30);
    cv = idx;

    idx = a >> 11;

    ct = (int)(((int64_t)av_costbl_4_sf[idx & 0x1f] * (0x800 - (a & 0x7ff)) +
                (int64_t)av_costbl_4_sf[(idx & 0x1f)+1]*(a & 0x7ff) +
                0x400) >> 11);
    st = (int)(((int64_t)av_sintbl_4_sf[idx & 0x1f] * (0x800 - (a & 0x7ff)) +
                (int64_t)av_sintbl_4_sf[(idx & 0x1f) + 1] * (a & 0x7ff) +
                0x400) >> 11);

    *c = (int)(((int64_t)cv * ct + (int64_t)sv * st + 0x20000000) >> 30);

    *s = (int)(((int64_t)cv * st + (int64_t)sv * ct + 0x20000000) >> 30);
}
288

289
#endif /* AVUTIL_SOFTFLOAT_H */