Re-enabled intra-frame prediction, which seems to have exposed a few issues
[opus.git] / libcelt / fftwrap.c
1 /* Copyright (C) 2005 Jean-Marc Valin 
2    File: fftwrap.c
3
4    Wrapper for various FFTs 
5
6    Redistribution and use in source and binary forms, with or without
7    modification, are permitted provided that the following conditions
8    are met:
9    
10    - Redistributions of source code must retain the above copyright
11    notice, this list of conditions and the following disclaimer.
12    
13    - Redistributions in binary form must reproduce the above copyright
14    notice, this list of conditions and the following disclaimer in the
15    documentation and/or other materials provided with the distribution.
16    
17    - Neither the name of the Xiph.org Foundation nor the names of its
18    contributors may be used to endorse or promote products derived from
19    this software without specific prior written permission.
20    
21    THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
22    ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
23    LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR
24    A PARTICULAR PURPOSE ARE DISCLAIMED.  IN NO EVENT SHALL THE FOUNDATION OR
25    CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
26    EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
27    PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
28    PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
29    LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
30    NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
31    SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
32
33 */
34
35 #ifdef HAVE_CONFIG_H
36 #include "config.h"
37 #endif
38
39 #define USE_SMALLFT
40 /*#define USE_KISS_FFT*/
41
42
43 #include "arch.h"
44 #include "os_support.h"
45
46 #ifdef FIXED_POINT
47 static int maximize_range(spx_word16_t *in, spx_word16_t *out, spx_word16_t bound, int len)
48 {
49    int i, shift;
50    spx_word16_t max_val = 0;
51    for (i=0;i<len;i++)
52    {
53       if (in[i]>max_val)
54          max_val = in[i];
55       if (-in[i]>max_val)
56          max_val = -in[i];
57    }
58    shift=0;
59    while (max_val <= (bound>>1) && max_val != 0)
60    {
61       max_val <<= 1;
62       shift++;
63    }
64    for (i=0;i<len;i++)
65    {
66       out[i] = in[i] << shift;
67    }   
68    return shift;
69 }
70
71 static void renorm_range(spx_word16_t *in, spx_word16_t *out, int shift, int len)
72 {
73    int i;
74    for (i=0;i<len;i++)
75    {
76       out[i] = (in[i] + (1<<(shift-1))) >> shift;
77    }
78 }
79 #endif
80
81 #ifdef USE_SMALLFT
82
83 #include "smallft.h"
84 #include <math.h>
85
86 void *spx_fft_init(int size)
87 {
88    struct drft_lookup *table;
89    table = celt_alloc(sizeof(struct drft_lookup));
90    spx_drft_init((struct drft_lookup *)table, size);
91    return (void*)table;
92 }
93
94 void spx_fft_destroy(void *table)
95 {
96    spx_drft_clear(table);
97    celt_free(table);
98 }
99
100 void spx_fft(void *table, float *in, float *out)
101 {
102    if (in==out)
103    {
104       int i;
105       celt_warning("FFT should not be done in-place");
106       float scale = 1./((struct drft_lookup *)table)->n;
107       for (i=0;i<((struct drft_lookup *)table)->n;i++)
108          out[i] = scale*in[i];
109    } else {
110       int i;
111       float scale = 1./((struct drft_lookup *)table)->n;
112       for (i=0;i<((struct drft_lookup *)table)->n;i++)
113          out[i] = scale*in[i];
114    }
115    spx_drft_forward((struct drft_lookup *)table, out);
116 }
117
118 void spx_ifft(void *table, float *in, float *out)
119 {
120    if (in==out)
121    {
122       celt_warning("FFT should not be done in-place");
123    } else {
124       int i;
125       for (i=0;i<((struct drft_lookup *)table)->n;i++)
126          out[i] = in[i];
127    }
128    spx_drft_backward((struct drft_lookup *)table, out);
129 }
130
131 #elif defined(USE_KISS_FFT)
132
133 #include "kiss_fftr.h"
134 #include "kiss_fft.h"
135
136 struct kiss_config {
137    kiss_fftr_cfg forward;
138    kiss_fftr_cfg backward;
139    kiss_fft_cpx *freq_data;
140    int N;
141 };
142
143 void *spx_fft_init(int size)
144 {
145    struct kiss_config *table;
146    table = speex_alloc(sizeof(struct kiss_config));
147    table->freq_data = speex_alloc(sizeof(kiss_fft_cpx)*((size>>1)+1));
148    table->forward = kiss_fftr_alloc(size,0,NULL,NULL);
149    table->backward = kiss_fftr_alloc(size,1,NULL,NULL);
150    table->N = size;
151    return table;
152 }
153
154 void spx_fft_destroy(void *table)
155 {
156    struct kiss_config *t = (struct kiss_config *)table;
157    kiss_fftr_free(t->forward);
158    kiss_fftr_free(t->backward);
159    speex_free(t->freq_data);
160    speex_free(table);
161 }
162
163 #ifdef FIXED_POINT
164
165 void spx_fft(void *table, spx_word16_t *in, spx_word16_t *out)
166 {
167    int i;
168    int shift;
169    struct kiss_config *t = (struct kiss_config *)table;
170    shift = maximize_range(in, in, 32000, t->N);
171    kiss_fftr(t->forward, in, t->freq_data);
172    out[0] = t->freq_data[0].r;
173    for (i=1;i<t->N>>1;i++)
174    {
175       out[(i<<1)-1] = t->freq_data[i].r;
176       out[(i<<1)] = t->freq_data[i].i;
177    }
178    out[(i<<1)-1] = t->freq_data[i].r;
179    renorm_range(in, in, shift, t->N);
180    renorm_range(out, out, shift, t->N);
181 }
182
183 #else
184
185 void spx_fft(void *table, spx_word16_t *in, spx_word16_t *out)
186 {
187    int i;
188    float scale;
189    struct kiss_config *t = (struct kiss_config *)table;
190    scale = 1./t->N;
191    kiss_fftr(t->forward, in, t->freq_data);
192    out[0] = scale*t->freq_data[0].r;
193    for (i=1;i<t->N>>1;i++)
194    {
195       out[(i<<1)-1] = scale*t->freq_data[i].r;
196       out[(i<<1)] = scale*t->freq_data[i].i;
197    }
198    out[(i<<1)-1] = scale*t->freq_data[i].r;
199 }
200 #endif
201
202 void spx_ifft(void *table, spx_word16_t *in, spx_word16_t *out)
203 {
204    int i;
205    struct kiss_config *t = (struct kiss_config *)table;
206    t->freq_data[0].r = in[0];
207    t->freq_data[0].i = 0;
208    for (i=1;i<t->N>>1;i++)
209    {
210       t->freq_data[i].r = in[(i<<1)-1];
211       t->freq_data[i].i = in[(i<<1)];
212    }
213    t->freq_data[i].r = in[(i<<1)-1];
214    t->freq_data[i].i = 0;
215
216    kiss_fftri(t->backward, t->freq_data, out);
217 }
218
219
220 #else
221
222 #error No other FFT implemented
223
224 #endif
225
226
227 int fixed_point = 1;
228 #ifdef FIXED_POINT
229 #include "smallft.h"
230
231
232 void spx_fft_float(void *table, float *in, float *out)
233 {
234    int i;
235 #ifdef USE_SMALLFT
236    int N = ((struct drft_lookup *)table)->n;
237 #elif defined(USE_KISS_FFT)
238    int N = ((struct kiss_config *)table)->N;
239 #else
240 #endif
241    spx_word16_t _in[N];
242    spx_word16_t _out[N];
243    for (i=0;i<N;i++)
244       _in[i] = (int)floor(.5+in[i]);
245    spx_fft(table, _in, _out);
246    for (i=0;i<N;i++)
247       out[i] = _out[i];
248    if (!fixed_point)
249    {
250       struct drft_lookup t;
251       spx_drft_init(&t, ((struct kiss_config *)table)->N);
252       float scale = 1./((struct kiss_config *)table)->N;
253       for (i=0;i<((struct kiss_config *)table)->N;i++)
254          out[i] = scale*in[i];
255       spx_drft_forward(&t, out);
256       spx_drft_clear(&t);
257    }
258 }
259
260 void spx_ifft_float(void *table, float *in, float *out)
261 {
262    int i;
263 #ifdef USE_SMALLFT
264    int N = ((struct drft_lookup *)table)->n;
265 #elif defined(USE_KISS_FFT)
266    int N = ((struct kiss_config *)table)->N;
267 #else
268 #endif
269    spx_word16_t _in[N];
270    spx_word16_t _out[N];
271    for (i=0;i<N;i++)
272       _in[i] = (int)floor(.5+in[i]);
273    spx_ifft(table, _in, _out);
274    for (i=0;i<N;i++)
275       out[i] = _out[i];
276    if (!fixed_point)
277    {
278       int i;
279       struct drft_lookup t;
280       spx_drft_init(&t, ((struct kiss_config *)table)->N);
281       for (i=0;i<((struct kiss_config *)table)->N;i++)
282          out[i] = in[i];
283       spx_drft_backward(&t, out);
284       spx_drft_clear(&t);
285    }
286 }
287
288 #else
289
290 void spx_fft_float(void *table, float *in, float *out)
291 {
292    spx_fft(table, in, out);
293 }
294 void spx_ifft_float(void *table, float *in, float *out)
295 {
296    spx_ifft(table, in, out);
297 }
298
299 #endif