~ppsspp/ppsspp/ffmpeg

« back to all changes in this revision

Viewing changes to libavcodec/fft_template.c

  • Committer: Henrik Rydgård
  • Date: 2014-01-03 10:44:32 UTC
  • Revision ID: git-v1:87c6c126784b1718bfa448ecf2e6a9fef781eb4e
Update our ffmpeg snapshot to a clone of the official repository.

This is because Maxim's at3plus support has been officially merged!

Show diffs side-by-side

added added

removed removed

Lines of Context:
 
1
/*
 
2
 * FFT/IFFT transforms
 
3
 * Copyright (c) 2008 Loren Merritt
 
4
 * Copyright (c) 2002 Fabrice Bellard
 
5
 * Partly based on libdjbfft by D. J. Bernstein
 
6
 *
 
7
 * This file is part of FFmpeg.
 
8
 *
 
9
 * FFmpeg is free software; you can redistribute it and/or
 
10
 * modify it under the terms of the GNU Lesser General Public
 
11
 * License as published by the Free Software Foundation; either
 
12
 * version 2.1 of the License, or (at your option) any later version.
 
13
 *
 
14
 * FFmpeg is distributed in the hope that it will be useful,
 
15
 * but WITHOUT ANY WARRANTY; without even the implied warranty of
 
16
 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
 
17
 * Lesser General Public License for more details.
 
18
 *
 
19
 * You should have received a copy of the GNU Lesser General Public
 
20
 * License along with FFmpeg; if not, write to the Free Software
 
21
 * Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA
 
22
 */
 
23
 
 
24
/**
 
25
 * @file
 
26
 * FFT/IFFT transforms.
 
27
 */
 
28
 
 
29
#include <stdlib.h>
 
30
#include <string.h>
 
31
#include "libavutil/mathematics.h"
 
32
#include "fft.h"
 
33
#include "fft-internal.h"
 
34
 
 
35
#if CONFIG_FFT_FIXED_32
 
36
#include "fft_table.h"
 
37
#else /* CONFIG_FFT_FIXED_32 */
 
38
 
 
39
/* cos(2*pi*x/n) for 0<=x<=n/4, followed by its reverse */
 
40
#if !CONFIG_HARDCODED_TABLES
 
41
COSTABLE(16);
 
42
COSTABLE(32);
 
43
COSTABLE(64);
 
44
COSTABLE(128);
 
45
COSTABLE(256);
 
46
COSTABLE(512);
 
47
COSTABLE(1024);
 
48
COSTABLE(2048);
 
49
COSTABLE(4096);
 
50
COSTABLE(8192);
 
51
COSTABLE(16384);
 
52
COSTABLE(32768);
 
53
COSTABLE(65536);
 
54
#endif
 
55
COSTABLE_CONST FFTSample * const FFT_NAME(ff_cos_tabs)[] = {
 
56
    NULL, NULL, NULL, NULL,
 
57
    FFT_NAME(ff_cos_16),
 
58
    FFT_NAME(ff_cos_32),
 
59
    FFT_NAME(ff_cos_64),
 
60
    FFT_NAME(ff_cos_128),
 
61
    FFT_NAME(ff_cos_256),
 
62
    FFT_NAME(ff_cos_512),
 
63
    FFT_NAME(ff_cos_1024),
 
64
    FFT_NAME(ff_cos_2048),
 
65
    FFT_NAME(ff_cos_4096),
 
66
    FFT_NAME(ff_cos_8192),
 
67
    FFT_NAME(ff_cos_16384),
 
68
    FFT_NAME(ff_cos_32768),
 
69
    FFT_NAME(ff_cos_65536),
 
70
};
 
71
 
 
72
#endif /* CONFIG_FFT_FIXED_32 */
 
73
 
 
74
static void fft_permute_c(FFTContext *s, FFTComplex *z);
 
75
static void fft_calc_c(FFTContext *s, FFTComplex *z);
 
76
 
 
77
static int split_radix_permutation(int i, int n, int inverse)
 
78
{
 
79
    int m;
 
80
    if(n <= 2) return i&1;
 
81
    m = n >> 1;
 
82
    if(!(i&m))            return split_radix_permutation(i, m, inverse)*2;
 
83
    m >>= 1;
 
84
    if(inverse == !(i&m)) return split_radix_permutation(i, m, inverse)*4 + 1;
 
85
    else                  return split_radix_permutation(i, m, inverse)*4 - 1;
 
86
}
 
87
 
 
88
av_cold void ff_init_ff_cos_tabs(int index)
 
89
{
 
90
#if (!CONFIG_HARDCODED_TABLES) && (!CONFIG_FFT_FIXED_32)
 
91
    int i;
 
92
    int m = 1<<index;
 
93
    double freq = 2*M_PI/m;
 
94
    FFTSample *tab = FFT_NAME(ff_cos_tabs)[index];
 
95
    for(i=0; i<=m/4; i++)
 
96
        tab[i] = FIX15(cos(i*freq));
 
97
    for(i=1; i<m/4; i++)
 
98
        tab[m/2-i] = tab[i];
 
99
#endif
 
100
}
 
101
 
 
102
static const int avx_tab[] = {
 
103
    0, 4, 1, 5, 8, 12, 9, 13, 2, 6, 3, 7, 10, 14, 11, 15
 
104
};
 
105
 
 
106
static int is_second_half_of_fft32(int i, int n)
 
107
{
 
108
    if (n <= 32)
 
109
        return i >= 16;
 
110
    else if (i < n/2)
 
111
        return is_second_half_of_fft32(i, n/2);
 
112
    else if (i < 3*n/4)
 
113
        return is_second_half_of_fft32(i - n/2, n/4);
 
114
    else
 
115
        return is_second_half_of_fft32(i - 3*n/4, n/4);
 
116
}
 
117
 
 
118
static av_cold void fft_perm_avx(FFTContext *s)
 
119
{
 
120
    int i;
 
121
    int n = 1 << s->nbits;
 
122
 
 
123
    for (i = 0; i < n; i += 16) {
 
124
        int k;
 
125
        if (is_second_half_of_fft32(i, n)) {
 
126
            for (k = 0; k < 16; k++)
 
127
                s->revtab[-split_radix_permutation(i + k, n, s->inverse) & (n - 1)] =
 
128
                    i + avx_tab[k];
 
129
 
 
130
        } else {
 
131
            for (k = 0; k < 16; k++) {
 
132
                int j = i + k;
 
133
                j = (j & ~7) | ((j >> 1) & 3) | ((j << 2) & 4);
 
134
                s->revtab[-split_radix_permutation(i + k, n, s->inverse) & (n - 1)] = j;
 
135
            }
 
136
        }
 
137
    }
 
138
}
 
139
 
 
140
av_cold int ff_fft_init(FFTContext *s, int nbits, int inverse)
 
141
{
 
142
    int i, j, n;
 
143
 
 
144
    if (nbits < 2 || nbits > 16)
 
145
        goto fail;
 
146
    s->nbits = nbits;
 
147
    n = 1 << nbits;
 
148
 
 
149
    s->revtab = av_malloc(n * sizeof(uint16_t));
 
150
    if (!s->revtab)
 
151
        goto fail;
 
152
    s->tmp_buf = av_malloc(n * sizeof(FFTComplex));
 
153
    if (!s->tmp_buf)
 
154
        goto fail;
 
155
    s->inverse = inverse;
 
156
    s->fft_permutation = FF_FFT_PERM_DEFAULT;
 
157
 
 
158
    s->fft_permute = fft_permute_c;
 
159
    s->fft_calc    = fft_calc_c;
 
160
#if CONFIG_MDCT
 
161
    s->imdct_calc  = ff_imdct_calc_c;
 
162
    s->imdct_half  = ff_imdct_half_c;
 
163
    s->mdct_calc   = ff_mdct_calc_c;
 
164
#endif
 
165
 
 
166
#if CONFIG_FFT_FIXED_32
 
167
    {
 
168
        int n=0;
 
169
        ff_fft_lut_init(fft_offsets_lut, 0, 1 << 16, &n);
 
170
    }
 
171
#else /* CONFIG_FFT_FIXED_32 */
 
172
#if CONFIG_FFT_FLOAT
 
173
    if (ARCH_ARM)     ff_fft_init_arm(s);
 
174
    if (ARCH_PPC)     ff_fft_init_ppc(s);
 
175
    if (ARCH_X86)     ff_fft_init_x86(s);
 
176
    if (CONFIG_MDCT)  s->mdct_calcw = s->mdct_calc;
 
177
    if (HAVE_MIPSFPU) ff_fft_init_mips(s);
 
178
#else
 
179
    if (CONFIG_MDCT)  s->mdct_calcw = ff_mdct_calcw_c;
 
180
    if (ARCH_ARM)     ff_fft_fixed_init_arm(s);
 
181
#endif
 
182
    for(j=4; j<=nbits; j++) {
 
183
        ff_init_ff_cos_tabs(j);
 
184
    }
 
185
#endif /* CONFIG_FFT_FIXED_32 */
 
186
 
 
187
 
 
188
    if (s->fft_permutation == FF_FFT_PERM_AVX) {
 
189
        fft_perm_avx(s);
 
190
    } else {
 
191
        for(i=0; i<n; i++) {
 
192
            j = i;
 
193
            if (s->fft_permutation == FF_FFT_PERM_SWAP_LSBS)
 
194
                j = (j&~3) | ((j>>1)&1) | ((j<<1)&2);
 
195
            s->revtab[-split_radix_permutation(i, n, s->inverse) & (n-1)] = j;
 
196
        }
 
197
    }
 
198
 
 
199
    return 0;
 
200
 fail:
 
201
    av_freep(&s->revtab);
 
202
    av_freep(&s->tmp_buf);
 
203
    return -1;
 
204
}
 
205
 
 
206
static void fft_permute_c(FFTContext *s, FFTComplex *z)
 
207
{
 
208
    int j, np;
 
209
    const uint16_t *revtab = s->revtab;
 
210
    np = 1 << s->nbits;
 
211
    /* TODO: handle split-radix permute in a more optimal way, probably in-place */
 
212
    for(j=0;j<np;j++) s->tmp_buf[revtab[j]] = z[j];
 
213
    memcpy(z, s->tmp_buf, np * sizeof(FFTComplex));
 
214
}
 
215
 
 
216
av_cold void ff_fft_end(FFTContext *s)
 
217
{
 
218
    av_freep(&s->revtab);
 
219
    av_freep(&s->tmp_buf);
 
220
}
 
221
 
 
222
#if CONFIG_FFT_FIXED_32
 
223
 
 
224
static void fft_calc_c(FFTContext *s, FFTComplex *z) {
 
225
 
 
226
    int nbits, i, n, num_transforms, offset, step;
 
227
    int n4, n2, n34;
 
228
    FFTSample tmp1, tmp2, tmp3, tmp4, tmp5, tmp6, tmp7, tmp8;
 
229
    FFTComplex *tmpz;
 
230
    FFTSample w_re, w_im;
 
231
    FFTSample *w_re_ptr, *w_im_ptr;
 
232
    const int fft_size = (1 << s->nbits);
 
233
    int64_t accu;
 
234
 
 
235
    num_transforms = (0x2aab >> (16 - s->nbits)) | 1;
 
236
 
 
237
    for (n=0; n<num_transforms; n++){
 
238
        offset = fft_offsets_lut[n] << 2;
 
239
        tmpz = z + offset;
 
240
 
 
241
        tmp1 = tmpz[0].re + tmpz[1].re;
 
242
        tmp5 = tmpz[2].re + tmpz[3].re;
 
243
        tmp2 = tmpz[0].im + tmpz[1].im;
 
244
        tmp6 = tmpz[2].im + tmpz[3].im;
 
245
        tmp3 = tmpz[0].re - tmpz[1].re;
 
246
        tmp8 = tmpz[2].im - tmpz[3].im;
 
247
        tmp4 = tmpz[0].im - tmpz[1].im;
 
248
        tmp7 = tmpz[2].re - tmpz[3].re;
 
249
 
 
250
        tmpz[0].re = tmp1 + tmp5;
 
251
        tmpz[2].re = tmp1 - tmp5;
 
252
        tmpz[0].im = tmp2 + tmp6;
 
253
        tmpz[2].im = tmp2 - tmp6;
 
254
        tmpz[1].re = tmp3 + tmp8;
 
255
        tmpz[3].re = tmp3 - tmp8;
 
256
        tmpz[1].im = tmp4 - tmp7;
 
257
        tmpz[3].im = tmp4 + tmp7;
 
258
    }
 
259
 
 
260
    if (fft_size < 8)
 
261
        return;
 
262
 
 
263
    num_transforms = (num_transforms >> 1) | 1;
 
264
 
 
265
    for (n=0; n<num_transforms; n++){
 
266
        offset = fft_offsets_lut[n] << 3;
 
267
        tmpz = z + offset;
 
268
 
 
269
        tmp1 = tmpz[4].re + tmpz[5].re;
 
270
        tmp3 = tmpz[6].re + tmpz[7].re;
 
271
        tmp2 = tmpz[4].im + tmpz[5].im;
 
272
        tmp4 = tmpz[6].im + tmpz[7].im;
 
273
        tmp5 = tmp1 + tmp3;
 
274
        tmp7 = tmp1 - tmp3;
 
275
        tmp6 = tmp2 + tmp4;
 
276
        tmp8 = tmp2 - tmp4;
 
277
 
 
278
        tmp1 = tmpz[4].re - tmpz[5].re;
 
279
        tmp2 = tmpz[4].im - tmpz[5].im;
 
280
        tmp3 = tmpz[6].re - tmpz[7].re;
 
281
        tmp4 = tmpz[6].im - tmpz[7].im;
 
282
 
 
283
        tmpz[4].re = tmpz[0].re - tmp5;
 
284
        tmpz[0].re = tmpz[0].re + tmp5;
 
285
        tmpz[4].im = tmpz[0].im - tmp6;
 
286
        tmpz[0].im = tmpz[0].im + tmp6;
 
287
        tmpz[6].re = tmpz[2].re - tmp8;
 
288
        tmpz[2].re = tmpz[2].re + tmp8;
 
289
        tmpz[6].im = tmpz[2].im + tmp7;
 
290
        tmpz[2].im = tmpz[2].im - tmp7;
 
291
 
 
292
        accu = (int64_t)Q31(M_SQRT1_2)*(tmp1 + tmp2);
 
293
        tmp5 = (int32_t)((accu + 0x40000000) >> 31);
 
294
        accu = (int64_t)Q31(M_SQRT1_2)*(tmp3 - tmp4);
 
295
        tmp7 = (int32_t)((accu + 0x40000000) >> 31);
 
296
        accu = (int64_t)Q31(M_SQRT1_2)*(tmp2 - tmp1);
 
297
        tmp6 = (int32_t)((accu + 0x40000000) >> 31);
 
298
        accu = (int64_t)Q31(M_SQRT1_2)*(tmp3 + tmp4);
 
299
        tmp8 = (int32_t)((accu + 0x40000000) >> 31);
 
300
        tmp1 = tmp5 + tmp7;
 
301
        tmp3 = tmp5 - tmp7;
 
302
        tmp2 = tmp6 + tmp8;
 
303
        tmp4 = tmp6 - tmp8;
 
304
 
 
305
        tmpz[5].re = tmpz[1].re - tmp1;
 
306
        tmpz[1].re = tmpz[1].re + tmp1;
 
307
        tmpz[5].im = tmpz[1].im - tmp2;
 
308
        tmpz[1].im = tmpz[1].im + tmp2;
 
309
        tmpz[7].re = tmpz[3].re - tmp4;
 
310
        tmpz[3].re = tmpz[3].re + tmp4;
 
311
        tmpz[7].im = tmpz[3].im + tmp3;
 
312
        tmpz[3].im = tmpz[3].im - tmp3;
 
313
    }
 
314
 
 
315
    step = 1 << ((MAX_LOG2_NFFT-4) - 4);
 
316
    n4 = 4;
 
317
 
 
318
    for (nbits=4; nbits<=s->nbits; nbits++){
 
319
        n2  = 2*n4;
 
320
        n34 = 3*n4;
 
321
        num_transforms = (num_transforms >> 1) | 1;
 
322
 
 
323
        for (n=0; n<num_transforms; n++){
 
324
            offset = fft_offsets_lut[n] << nbits;
 
325
            tmpz = z + offset;
 
326
 
 
327
            tmp5 = tmpz[ n2].re + tmpz[n34].re;
 
328
            tmp1 = tmpz[ n2].re - tmpz[n34].re;
 
329
            tmp6 = tmpz[ n2].im + tmpz[n34].im;
 
330
            tmp2 = tmpz[ n2].im - tmpz[n34].im;
 
331
 
 
332
            tmpz[ n2].re = tmpz[ 0].re - tmp5;
 
333
            tmpz[  0].re = tmpz[ 0].re + tmp5;
 
334
            tmpz[ n2].im = tmpz[ 0].im - tmp6;
 
335
            tmpz[  0].im = tmpz[ 0].im + tmp6;
 
336
            tmpz[n34].re = tmpz[n4].re - tmp2;
 
337
            tmpz[ n4].re = tmpz[n4].re + tmp2;
 
338
            tmpz[n34].im = tmpz[n4].im + tmp1;
 
339
            tmpz[ n4].im = tmpz[n4].im - tmp1;
 
340
 
 
341
            w_re_ptr = w_tab_sr + step;
 
342
            w_im_ptr = w_tab_sr + MAX_FFT_SIZE/(4*16) - step;
 
343
 
 
344
            for (i=1; i<n4; i++){
 
345
                w_re = w_re_ptr[0];
 
346
                w_im = w_im_ptr[0];
 
347
                accu  = (int64_t)w_re*tmpz[ n2+i].re;
 
348
                accu += (int64_t)w_im*tmpz[ n2+i].im;
 
349
                tmp1 = (int32_t)((accu + 0x40000000) >> 31);
 
350
                accu  = (int64_t)w_re*tmpz[ n2+i].im;
 
351
                accu -= (int64_t)w_im*tmpz[ n2+i].re;
 
352
                tmp2 = (int32_t)((accu + 0x40000000) >> 31);
 
353
                accu  = (int64_t)w_re*tmpz[n34+i].re;
 
354
                accu -= (int64_t)w_im*tmpz[n34+i].im;
 
355
                tmp3 = (int32_t)((accu + 0x40000000) >> 31);
 
356
                accu  = (int64_t)w_re*tmpz[n34+i].im;
 
357
                accu += (int64_t)w_im*tmpz[n34+i].re;
 
358
                tmp4 = (int32_t)((accu + 0x40000000) >> 31);
 
359
 
 
360
                tmp5 = tmp1 + tmp3;
 
361
                tmp1 = tmp1 - tmp3;
 
362
                tmp6 = tmp2 + tmp4;
 
363
                tmp2 = tmp2 - tmp4;
 
364
 
 
365
                tmpz[ n2+i].re = tmpz[   i].re - tmp5;
 
366
                tmpz[    i].re = tmpz[   i].re + tmp5;
 
367
                tmpz[ n2+i].im = tmpz[   i].im - tmp6;
 
368
                tmpz[    i].im = tmpz[   i].im + tmp6;
 
369
                tmpz[n34+i].re = tmpz[n4+i].re - tmp2;
 
370
                tmpz[ n4+i].re = tmpz[n4+i].re + tmp2;
 
371
                tmpz[n34+i].im = tmpz[n4+i].im + tmp1;
 
372
                tmpz[ n4+i].im = tmpz[n4+i].im - tmp1;
 
373
 
 
374
                w_re_ptr += step;
 
375
                w_im_ptr -= step;
 
376
            }
 
377
        }
 
378
        step >>= 1;
 
379
        n4   <<= 1;
 
380
    }
 
381
}
 
382
 
 
383
#else /* CONFIG_FFT_FIXED_32 */
 
384
 
 
385
#define BUTTERFLIES(a0,a1,a2,a3) {\
 
386
    BF(t3, t5, t5, t1);\
 
387
    BF(a2.re, a0.re, a0.re, t5);\
 
388
    BF(a3.im, a1.im, a1.im, t3);\
 
389
    BF(t4, t6, t2, t6);\
 
390
    BF(a3.re, a1.re, a1.re, t4);\
 
391
    BF(a2.im, a0.im, a0.im, t6);\
 
392
}
 
393
 
 
394
// force loading all the inputs before storing any.
 
395
// this is slightly slower for small data, but avoids store->load aliasing
 
396
// for addresses separated by large powers of 2.
 
397
#define BUTTERFLIES_BIG(a0,a1,a2,a3) {\
 
398
    FFTSample r0=a0.re, i0=a0.im, r1=a1.re, i1=a1.im;\
 
399
    BF(t3, t5, t5, t1);\
 
400
    BF(a2.re, a0.re, r0, t5);\
 
401
    BF(a3.im, a1.im, i1, t3);\
 
402
    BF(t4, t6, t2, t6);\
 
403
    BF(a3.re, a1.re, r1, t4);\
 
404
    BF(a2.im, a0.im, i0, t6);\
 
405
}
 
406
 
 
407
#define TRANSFORM(a0,a1,a2,a3,wre,wim) {\
 
408
    CMUL(t1, t2, a2.re, a2.im, wre, -wim);\
 
409
    CMUL(t5, t6, a3.re, a3.im, wre,  wim);\
 
410
    BUTTERFLIES(a0,a1,a2,a3)\
 
411
}
 
412
 
 
413
#define TRANSFORM_ZERO(a0,a1,a2,a3) {\
 
414
    t1 = a2.re;\
 
415
    t2 = a2.im;\
 
416
    t5 = a3.re;\
 
417
    t6 = a3.im;\
 
418
    BUTTERFLIES(a0,a1,a2,a3)\
 
419
}
 
420
 
 
421
/* z[0...8n-1], w[1...2n-1] */
 
422
#define PASS(name)\
 
423
static void name(FFTComplex *z, const FFTSample *wre, unsigned int n)\
 
424
{\
 
425
    FFTDouble t1, t2, t3, t4, t5, t6;\
 
426
    int o1 = 2*n;\
 
427
    int o2 = 4*n;\
 
428
    int o3 = 6*n;\
 
429
    const FFTSample *wim = wre+o1;\
 
430
    n--;\
 
431
\
 
432
    TRANSFORM_ZERO(z[0],z[o1],z[o2],z[o3]);\
 
433
    TRANSFORM(z[1],z[o1+1],z[o2+1],z[o3+1],wre[1],wim[-1]);\
 
434
    do {\
 
435
        z += 2;\
 
436
        wre += 2;\
 
437
        wim -= 2;\
 
438
        TRANSFORM(z[0],z[o1],z[o2],z[o3],wre[0],wim[0]);\
 
439
        TRANSFORM(z[1],z[o1+1],z[o2+1],z[o3+1],wre[1],wim[-1]);\
 
440
    } while(--n);\
 
441
}
 
442
 
 
443
PASS(pass)
 
444
#undef BUTTERFLIES
 
445
#define BUTTERFLIES BUTTERFLIES_BIG
 
446
PASS(pass_big)
 
447
 
 
448
#define DECL_FFT(n,n2,n4)\
 
449
static void fft##n(FFTComplex *z)\
 
450
{\
 
451
    fft##n2(z);\
 
452
    fft##n4(z+n4*2);\
 
453
    fft##n4(z+n4*3);\
 
454
    pass(z,FFT_NAME(ff_cos_##n),n4/2);\
 
455
}
 
456
 
 
457
static void fft4(FFTComplex *z)
 
458
{
 
459
    FFTDouble t1, t2, t3, t4, t5, t6, t7, t8;
 
460
 
 
461
    BF(t3, t1, z[0].re, z[1].re);
 
462
    BF(t8, t6, z[3].re, z[2].re);
 
463
    BF(z[2].re, z[0].re, t1, t6);
 
464
    BF(t4, t2, z[0].im, z[1].im);
 
465
    BF(t7, t5, z[2].im, z[3].im);
 
466
    BF(z[3].im, z[1].im, t4, t8);
 
467
    BF(z[3].re, z[1].re, t3, t7);
 
468
    BF(z[2].im, z[0].im, t2, t5);
 
469
}
 
470
 
 
471
static void fft8(FFTComplex *z)
 
472
{
 
473
    FFTDouble t1, t2, t3, t4, t5, t6;
 
474
 
 
475
    fft4(z);
 
476
 
 
477
    BF(t1, z[5].re, z[4].re, -z[5].re);
 
478
    BF(t2, z[5].im, z[4].im, -z[5].im);
 
479
    BF(t5, z[7].re, z[6].re, -z[7].re);
 
480
    BF(t6, z[7].im, z[6].im, -z[7].im);
 
481
 
 
482
    BUTTERFLIES(z[0],z[2],z[4],z[6]);
 
483
    TRANSFORM(z[1],z[3],z[5],z[7],sqrthalf,sqrthalf);
 
484
}
 
485
 
 
486
#if !CONFIG_SMALL
 
487
static void fft16(FFTComplex *z)
 
488
{
 
489
    FFTDouble t1, t2, t3, t4, t5, t6;
 
490
    FFTSample cos_16_1 = FFT_NAME(ff_cos_16)[1];
 
491
    FFTSample cos_16_3 = FFT_NAME(ff_cos_16)[3];
 
492
 
 
493
    fft8(z);
 
494
    fft4(z+8);
 
495
    fft4(z+12);
 
496
 
 
497
    TRANSFORM_ZERO(z[0],z[4],z[8],z[12]);
 
498
    TRANSFORM(z[2],z[6],z[10],z[14],sqrthalf,sqrthalf);
 
499
    TRANSFORM(z[1],z[5],z[9],z[13],cos_16_1,cos_16_3);
 
500
    TRANSFORM(z[3],z[7],z[11],z[15],cos_16_3,cos_16_1);
 
501
}
 
502
#else
 
503
DECL_FFT(16,8,4)
 
504
#endif
 
505
DECL_FFT(32,16,8)
 
506
DECL_FFT(64,32,16)
 
507
DECL_FFT(128,64,32)
 
508
DECL_FFT(256,128,64)
 
509
DECL_FFT(512,256,128)
 
510
#if !CONFIG_SMALL
 
511
#define pass pass_big
 
512
#endif
 
513
DECL_FFT(1024,512,256)
 
514
DECL_FFT(2048,1024,512)
 
515
DECL_FFT(4096,2048,1024)
 
516
DECL_FFT(8192,4096,2048)
 
517
DECL_FFT(16384,8192,4096)
 
518
DECL_FFT(32768,16384,8192)
 
519
DECL_FFT(65536,32768,16384)
 
520
 
 
521
static void (* const fft_dispatch[])(FFTComplex*) = {
 
522
    fft4, fft8, fft16, fft32, fft64, fft128, fft256, fft512, fft1024,
 
523
    fft2048, fft4096, fft8192, fft16384, fft32768, fft65536,
 
524
};
 
525
 
 
526
static void fft_calc_c(FFTContext *s, FFTComplex *z)
 
527
{
 
528
    fft_dispatch[s->nbits-2](z);
 
529
}
 
530
#endif /* CONFIG_FFT_FIXED_32 */