Implements CELT_SET_LOSS_PERC
[opus.git] / libcelt / kiss_fft.c
1 /*
2 Copyright (c) 2003-2004, Mark Borgerding
3 Lots of modifications by Jean-Marc Valin
4 Copyright (c) 2005-2007, Xiph.Org Foundation
5 Copyright (c) 2008,      Xiph.Org Foundation, CSIRO
6
7 All rights reserved.
8
9 Redistribution and use in source and binary forms, with or without modification, are permitted provided that the following conditions are met:
10
11     * Redistributions of source code must retain the above copyright notice, this list of conditions and the following disclaimer.
12     * Redistributions in binary form must reproduce the above copyright notice, this list of conditions and the following disclaimer in the documentation and/or other materials provided with the distribution.
13
14 THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
15 */
16
17 #ifndef SKIP_CONFIG_H
18 #  ifdef HAVE_CONFIG_H
19 #    include "config.h"
20 #  endif
21 #endif
22
23 #include "_kiss_fft_guts.h"
24 #include "arch.h"
25 #include "os_support.h"
26 #include "mathops.h"
27 #include "stack_alloc.h"
28
29 /* The guts header contains all the multiplication and addition macros that are defined for
30    complex numbers.  It also delares the kf_ internal functions.
31 */
32
33 static void kf_bfly2(
34                      kiss_fft_cpx * Fout,
35                      const size_t fstride,
36                      const kiss_fft_state *st,
37                      int m,
38                      int N,
39                      int mm
40                     )
41 {
42    kiss_fft_cpx * Fout2;
43    const kiss_twiddle_cpx * tw1;
44    int i,j;
45    kiss_fft_cpx * Fout_beg = Fout;
46    for (i=0;i<N;i++)
47    {
48       Fout = Fout_beg + i*mm;
49       Fout2 = Fout + m;
50       tw1 = st->twiddles;
51       for(j=0;j<m;j++)
52       {
53          kiss_fft_cpx t;
54          Fout->r = SHR(Fout->r, 1);Fout->i = SHR(Fout->i, 1);
55          Fout2->r = SHR(Fout2->r, 1);Fout2->i = SHR(Fout2->i, 1);
56          C_MUL (t,  *Fout2 , *tw1);
57          tw1 += fstride;
58          C_SUB( *Fout2 ,  *Fout , t );
59          C_ADDTO( *Fout ,  t );
60          ++Fout2;
61          ++Fout;
62       }
63    }
64 }
65
66 static void ki_bfly2(
67                      kiss_fft_cpx * Fout,
68                      const size_t fstride,
69                      const kiss_fft_state *st,
70                      int m,
71                      int N,
72                      int mm
73                     )
74 {
75    kiss_fft_cpx * Fout2;
76    const kiss_twiddle_cpx * tw1;
77    kiss_fft_cpx t;
78    int i,j;
79    kiss_fft_cpx * Fout_beg = Fout;
80    for (i=0;i<N;i++)
81    {
82       Fout = Fout_beg + i*mm;
83       Fout2 = Fout + m;
84       tw1 = st->twiddles;
85       for(j=0;j<m;j++)
86       {
87          C_MULC (t,  *Fout2 , *tw1);
88          tw1 += fstride;
89          C_SUB( *Fout2 ,  *Fout , t );
90          C_ADDTO( *Fout ,  t );
91          ++Fout2;
92          ++Fout;
93       }
94    }
95 }
96
97 static void kf_bfly4(
98                      kiss_fft_cpx * Fout,
99                      const size_t fstride,
100                      const kiss_fft_state *st,
101                      int m,
102                      int N,
103                      int mm
104                     )
105 {
106    const kiss_twiddle_cpx *tw1,*tw2,*tw3;
107    kiss_fft_cpx scratch[6];
108    const size_t m2=2*m;
109    const size_t m3=3*m;
110    int i, j;
111
112    kiss_fft_cpx * Fout_beg = Fout;
113    for (i=0;i<N;i++)
114    {
115       Fout = Fout_beg + i*mm;
116       tw3 = tw2 = tw1 = st->twiddles;
117       for (j=0;j<m;j++)
118       {
119          C_MUL4(scratch[0],Fout[m] , *tw1 );
120          C_MUL4(scratch[1],Fout[m2] , *tw2 );
121          C_MUL4(scratch[2],Fout[m3] , *tw3 );
122              
123          Fout->r = PSHR(Fout->r, 2);
124          Fout->i = PSHR(Fout->i, 2);
125          C_SUB( scratch[5] , *Fout, scratch[1] );
126          C_ADDTO(*Fout, scratch[1]);
127          C_ADD( scratch[3] , scratch[0] , scratch[2] );
128          C_SUB( scratch[4] , scratch[0] , scratch[2] );
129          Fout[m2].r = PSHR(Fout[m2].r, 2);
130          Fout[m2].i = PSHR(Fout[m2].i, 2);
131          C_SUB( Fout[m2], *Fout, scratch[3] );
132          tw1 += fstride;
133          tw2 += fstride*2;
134          tw3 += fstride*3;
135          C_ADDTO( *Fout , scratch[3] );
136              
137          Fout[m].r = scratch[5].r + scratch[4].i;
138          Fout[m].i = scratch[5].i - scratch[4].r;
139          Fout[m3].r = scratch[5].r - scratch[4].i;
140          Fout[m3].i = scratch[5].i + scratch[4].r;
141          ++Fout;
142       }
143    }
144 }
145
146 static void ki_bfly4(
147                      kiss_fft_cpx * Fout,
148                      const size_t fstride,
149                      const kiss_fft_state *st,
150                      int m,
151                      int N,
152                      int mm
153                     )
154 {
155    const kiss_twiddle_cpx *tw1,*tw2,*tw3;
156    kiss_fft_cpx scratch[6];
157    const size_t m2=2*m;
158    const size_t m3=3*m;
159    int i, j;
160
161    kiss_fft_cpx * Fout_beg = Fout;
162    for (i=0;i<N;i++)
163    {
164       Fout = Fout_beg + i*mm;
165       tw3 = tw2 = tw1 = st->twiddles;
166       for (j=0;j<m;j++)
167       {
168          C_MULC(scratch[0],Fout[m] , *tw1 );
169          C_MULC(scratch[1],Fout[m2] , *tw2 );
170          C_MULC(scratch[2],Fout[m3] , *tw3 );
171              
172          C_SUB( scratch[5] , *Fout, scratch[1] );
173          C_ADDTO(*Fout, scratch[1]);
174          C_ADD( scratch[3] , scratch[0] , scratch[2] );
175          C_SUB( scratch[4] , scratch[0] , scratch[2] );
176          C_SUB( Fout[m2], *Fout, scratch[3] );
177          tw1 += fstride;
178          tw2 += fstride*2;
179          tw3 += fstride*3;
180          C_ADDTO( *Fout , scratch[3] );
181              
182          Fout[m].r = scratch[5].r - scratch[4].i;
183          Fout[m].i = scratch[5].i + scratch[4].r;
184          Fout[m3].r = scratch[5].r + scratch[4].i;
185          Fout[m3].i = scratch[5].i - scratch[4].r;
186          ++Fout;
187       }
188    }
189 }
190
191 #ifndef RADIX_TWO_ONLY
192
193 static void kf_bfly3(
194                      kiss_fft_cpx * Fout,
195                      const size_t fstride,
196                      const kiss_fft_state *st,
197                      int m,
198                      int N,
199                      int mm
200                     )
201 {
202    int i;
203    size_t k;
204    const size_t m2 = 2*m;
205    const kiss_twiddle_cpx *tw1,*tw2;
206    kiss_fft_cpx scratch[5];
207    kiss_twiddle_cpx epi3;
208
209    kiss_fft_cpx * Fout_beg = Fout;
210    epi3 = st->twiddles[fstride*m];
211    for (i=0;i<N;i++)
212    {
213       Fout = Fout_beg + i*mm;
214       tw1=tw2=st->twiddles;
215       k=m;
216       do {
217          C_FIXDIV(*Fout,3); C_FIXDIV(Fout[m],3); C_FIXDIV(Fout[m2],3);
218
219          C_MUL(scratch[1],Fout[m] , *tw1);
220          C_MUL(scratch[2],Fout[m2] , *tw2);
221
222          C_ADD(scratch[3],scratch[1],scratch[2]);
223          C_SUB(scratch[0],scratch[1],scratch[2]);
224          tw1 += fstride;
225          tw2 += fstride*2;
226
227          Fout[m].r = Fout->r - HALF_OF(scratch[3].r);
228          Fout[m].i = Fout->i - HALF_OF(scratch[3].i);
229
230          C_MULBYSCALAR( scratch[0] , epi3.i );
231
232          C_ADDTO(*Fout,scratch[3]);
233
234          Fout[m2].r = Fout[m].r + scratch[0].i;
235          Fout[m2].i = Fout[m].i - scratch[0].r;
236
237          Fout[m].r -= scratch[0].i;
238          Fout[m].i += scratch[0].r;
239
240          ++Fout;
241       } while(--k);
242    }
243 }
244
245 static void ki_bfly3(
246                      kiss_fft_cpx * Fout,
247                      const size_t fstride,
248                      const kiss_fft_state *st,
249                      size_t m,
250                      int N,
251                      int mm
252                     )
253 {
254    size_t i, k;
255    const size_t m2 = 2*m;
256    const kiss_twiddle_cpx *tw1,*tw2;
257    kiss_fft_cpx scratch[5];
258    kiss_twiddle_cpx epi3;
259
260    kiss_fft_cpx * Fout_beg = Fout;
261    epi3 = st->twiddles[fstride*m];
262    for (i=0;i<N;i++)
263    {
264       Fout = Fout_beg + i*mm;
265       tw1=tw2=st->twiddles;
266       k=m;
267       do{
268
269          C_MULC(scratch[1],Fout[m] , *tw1);
270          C_MULC(scratch[2],Fout[m2] , *tw2);
271
272          C_ADD(scratch[3],scratch[1],scratch[2]);
273          C_SUB(scratch[0],scratch[1],scratch[2]);
274          tw1 += fstride;
275          tw2 += fstride*2;
276
277          Fout[m].r = Fout->r - HALF_OF(scratch[3].r);
278          Fout[m].i = Fout->i - HALF_OF(scratch[3].i);
279
280          C_MULBYSCALAR( scratch[0] , -epi3.i );
281
282          C_ADDTO(*Fout,scratch[3]);
283
284          Fout[m2].r = Fout[m].r + scratch[0].i;
285          Fout[m2].i = Fout[m].i - scratch[0].r;
286
287          Fout[m].r -= scratch[0].i;
288          Fout[m].i += scratch[0].r;
289
290          ++Fout;
291       }while(--k);
292    }
293 }
294
295
296 static void kf_bfly5(
297                      kiss_fft_cpx * Fout,
298                      const size_t fstride,
299                      const kiss_fft_state *st,
300                      int m,
301                      int N,
302                      int mm
303                     )
304 {
305    kiss_fft_cpx *Fout0,*Fout1,*Fout2,*Fout3,*Fout4;
306    int i, u;
307    kiss_fft_cpx scratch[13];
308    const kiss_twiddle_cpx * twiddles = st->twiddles;
309    const kiss_twiddle_cpx *tw;
310    kiss_twiddle_cpx ya,yb;
311    kiss_fft_cpx * Fout_beg = Fout;
312
313    ya = twiddles[fstride*m];
314    yb = twiddles[fstride*2*m];
315    tw=st->twiddles;
316
317    for (i=0;i<N;i++)
318    {
319       Fout = Fout_beg + i*mm;
320       Fout0=Fout;
321       Fout1=Fout0+m;
322       Fout2=Fout0+2*m;
323       Fout3=Fout0+3*m;
324       Fout4=Fout0+4*m;
325
326       for ( u=0; u<m; ++u ) {
327          C_FIXDIV( *Fout0,5); C_FIXDIV( *Fout1,5); C_FIXDIV( *Fout2,5); C_FIXDIV( *Fout3,5); C_FIXDIV( *Fout4,5);
328          scratch[0] = *Fout0;
329
330          C_MUL(scratch[1] ,*Fout1, tw[u*fstride]);
331          C_MUL(scratch[2] ,*Fout2, tw[2*u*fstride]);
332          C_MUL(scratch[3] ,*Fout3, tw[3*u*fstride]);
333          C_MUL(scratch[4] ,*Fout4, tw[4*u*fstride]);
334
335          C_ADD( scratch[7],scratch[1],scratch[4]);
336          C_SUB( scratch[10],scratch[1],scratch[4]);
337          C_ADD( scratch[8],scratch[2],scratch[3]);
338          C_SUB( scratch[9],scratch[2],scratch[3]);
339
340          Fout0->r += scratch[7].r + scratch[8].r;
341          Fout0->i += scratch[7].i + scratch[8].i;
342
343          scratch[5].r = scratch[0].r + S_MUL(scratch[7].r,ya.r) + S_MUL(scratch[8].r,yb.r);
344          scratch[5].i = scratch[0].i + S_MUL(scratch[7].i,ya.r) + S_MUL(scratch[8].i,yb.r);
345
346          scratch[6].r =  S_MUL(scratch[10].i,ya.i) + S_MUL(scratch[9].i,yb.i);
347          scratch[6].i = -S_MUL(scratch[10].r,ya.i) - S_MUL(scratch[9].r,yb.i);
348
349          C_SUB(*Fout1,scratch[5],scratch[6]);
350          C_ADD(*Fout4,scratch[5],scratch[6]);
351
352          scratch[11].r = scratch[0].r + S_MUL(scratch[7].r,yb.r) + S_MUL(scratch[8].r,ya.r);
353          scratch[11].i = scratch[0].i + S_MUL(scratch[7].i,yb.r) + S_MUL(scratch[8].i,ya.r);
354          scratch[12].r = - S_MUL(scratch[10].i,yb.i) + S_MUL(scratch[9].i,ya.i);
355          scratch[12].i = S_MUL(scratch[10].r,yb.i) - S_MUL(scratch[9].r,ya.i);
356
357          C_ADD(*Fout2,scratch[11],scratch[12]);
358          C_SUB(*Fout3,scratch[11],scratch[12]);
359
360          ++Fout0;++Fout1;++Fout2;++Fout3;++Fout4;
361       }
362    }
363 }
364
365 static void ki_bfly5(
366                      kiss_fft_cpx * Fout,
367                      const size_t fstride,
368                      const kiss_fft_state *st,
369                      int m,
370                      int N,
371                      int mm
372                     )
373 {
374    kiss_fft_cpx *Fout0,*Fout1,*Fout2,*Fout3,*Fout4;
375    int i, u;
376    kiss_fft_cpx scratch[13];
377    const kiss_twiddle_cpx * twiddles = st->twiddles;
378    const kiss_twiddle_cpx *tw;
379    kiss_twiddle_cpx ya,yb;
380    kiss_fft_cpx * Fout_beg = Fout;
381
382    ya = twiddles[fstride*m];
383    yb = twiddles[fstride*2*m];
384    tw=st->twiddles;
385
386    for (i=0;i<N;i++)
387    {
388       Fout = Fout_beg + i*mm;
389       Fout0=Fout;
390       Fout1=Fout0+m;
391       Fout2=Fout0+2*m;
392       Fout3=Fout0+3*m;
393       Fout4=Fout0+4*m;
394
395       for ( u=0; u<m; ++u ) {
396          scratch[0] = *Fout0;
397
398          C_MULC(scratch[1] ,*Fout1, tw[u*fstride]);
399          C_MULC(scratch[2] ,*Fout2, tw[2*u*fstride]);
400          C_MULC(scratch[3] ,*Fout3, tw[3*u*fstride]);
401          C_MULC(scratch[4] ,*Fout4, tw[4*u*fstride]);
402
403          C_ADD( scratch[7],scratch[1],scratch[4]);
404          C_SUB( scratch[10],scratch[1],scratch[4]);
405          C_ADD( scratch[8],scratch[2],scratch[3]);
406          C_SUB( scratch[9],scratch[2],scratch[3]);
407
408          Fout0->r += scratch[7].r + scratch[8].r;
409          Fout0->i += scratch[7].i + scratch[8].i;
410
411          scratch[5].r = scratch[0].r + S_MUL(scratch[7].r,ya.r) + S_MUL(scratch[8].r,yb.r);
412          scratch[5].i = scratch[0].i + S_MUL(scratch[7].i,ya.r) + S_MUL(scratch[8].i,yb.r);
413
414          scratch[6].r = -S_MUL(scratch[10].i,ya.i) - S_MUL(scratch[9].i,yb.i);
415          scratch[6].i =  S_MUL(scratch[10].r,ya.i) + S_MUL(scratch[9].r,yb.i);
416
417          C_SUB(*Fout1,scratch[5],scratch[6]);
418          C_ADD(*Fout4,scratch[5],scratch[6]);
419
420          scratch[11].r = scratch[0].r + S_MUL(scratch[7].r,yb.r) + S_MUL(scratch[8].r,ya.r);
421          scratch[11].i = scratch[0].i + S_MUL(scratch[7].i,yb.r) + S_MUL(scratch[8].i,ya.r);
422          scratch[12].r =  S_MUL(scratch[10].i,yb.i) - S_MUL(scratch[9].i,ya.i);
423          scratch[12].i = -S_MUL(scratch[10].r,yb.i) + S_MUL(scratch[9].r,ya.i);
424
425          C_ADD(*Fout2,scratch[11],scratch[12]);
426          C_SUB(*Fout3,scratch[11],scratch[12]);
427
428          ++Fout0;++Fout1;++Fout2;++Fout3;++Fout4;
429       }
430    }
431 }
432
433 #endif
434
435 static void kf_work(
436         kiss_fft_cpx * Fout,
437         const kiss_fft_cpx * f,
438         size_t fstride,
439         int in_stride,
440         const celt_int16 * factors,
441         const kiss_fft_state *st,
442         int N,
443         int m2
444         )
445 {
446     const int p=*factors++; /* the radix  */
447     const int m=*factors++; /* stage's fft length/p */
448     /*printf ("fft %d %d %d %d %d %d %d\n", p*m, m, p, s2, fstride*in_stride, N, m2);*/
449     if (m!=1) 
450         kf_work( Fout , f, fstride*p, in_stride, factors,st, N*p, m);
451
452     /* Compensate for longer twiddles table (when sharing) */
453     if (st->shift>0)
454        fstride <<= st->shift;
455     switch (p) {
456         case 2: kf_bfly2(Fout,fstride,st,m, N, m2); break;
457         case 4: kf_bfly4(Fout,fstride,st,m, N, m2); break;
458 #ifndef RADIX_TWO_ONLY
459         case 3: kf_bfly3(Fout,fstride,st,m, N, m2); break;
460         case 5: kf_bfly5(Fout,fstride,st,m, N, m2); break;
461 #endif
462     }    
463 }
464
465
466 static void ki_work(
467              kiss_fft_cpx * Fout,
468              const kiss_fft_cpx * f,
469              size_t fstride,
470              int in_stride,
471              const celt_int16 * factors,
472              const kiss_fft_state *st,
473              int N,
474              int m2
475             )
476 {
477    const int p=*factors++; /* the radix  */
478    const int m=*factors++; /* stage's fft length/p */
479    /*printf ("fft %d %d %d %d %d %d %d\n", p*m, m, p, s2, fstride*in_stride, N, m2);*/
480    if (m!=1) 
481       ki_work( Fout , f, fstride*p, in_stride, factors,st, N*p, m);
482
483    /* Compensate for longer twiddles table (when sharing) */
484    if (st->shift>0)
485       fstride <<= st->shift;
486    switch (p) {
487       case 2: ki_bfly2(Fout,fstride,st,m, N, m2); break;
488       case 4: ki_bfly4(Fout,fstride,st,m, N, m2); break;
489 #ifndef RADIX_TWO_ONLY
490       case 3: ki_bfly3(Fout,fstride,st,m, N, m2); break;
491       case 5: ki_bfly5(Fout,fstride,st,m, N, m2); break;
492 #endif
493    }    
494 }
495
496
497 #ifdef CUSTOM_MODES
498
499 static
500 void compute_bitrev_table(
501          int Fout,
502          celt_int16 *f,
503          const size_t fstride,
504          int in_stride,
505          celt_int16 * factors,
506          const kiss_fft_state *st
507             )
508 {
509    const int p=*factors++; /* the radix  */
510    const int m=*factors++; /* stage's fft length/p */
511
512     /*printf ("fft %d %d %d %d %d %d\n", p*m, m, p, s2, fstride*in_stride, N);*/
513    if (m==1)
514    {
515       int j;
516       for (j=0;j<p;j++)
517       {
518          *f = Fout+j;
519          f += fstride*in_stride;
520       }
521    } else {
522       int j;
523       for (j=0;j<p;j++)
524       {
525          compute_bitrev_table( Fout , f, fstride*p, in_stride, factors,st);
526          f += fstride*in_stride;
527          Fout += m;
528       }
529    }
530 }
531
532
533 /*  facbuf is populated by p1,m1,p2,m2, ...
534     where 
535     p[i] * m[i] = m[i-1]
536     m0 = n                  */
537 static 
538 int kf_factor(int n,celt_int16 * facbuf)
539 {
540     int p=4;
541
542     /*factor out powers of 4, powers of 2, then any remaining primes */
543     do {
544         while (n % p) {
545             switch (p) {
546                 case 4: p = 2; break;
547                 case 2: p = 3; break;
548                 default: p += 2; break;
549             }
550             if (p>32000 || (celt_int32)p*(celt_int32)p > n)
551                 p = n;          /* no more factors, skip to end */
552         }
553         n /= p;
554 #ifdef RADIX_TWO_ONLY
555         if (p!=2 && p != 4)
556 #else
557         if (p>5)
558 #endif
559         {
560            return 0;
561         }
562         *facbuf++ = p;
563         *facbuf++ = n;
564     } while (n > 1);
565     return 1;
566 }
567
568 static void compute_twiddles(kiss_twiddle_cpx *twiddles, int nfft)
569 {
570    int i;
571 #ifdef FIXED_POINT
572    for (i=0;i<nfft;++i) {
573       celt_word32 phase = -i;
574       kf_cexp2(twiddles+i, DIV32(SHL32(phase,17),nfft));
575    }
576 #else
577    for (i=0;i<nfft;++i) {
578       const double pi=3.14159265358979323846264338327;
579       double phase = ( -2*pi /nfft ) * i;
580       kf_cexp(twiddles+i, phase );
581    }
582 #endif
583 }
584
585
586 /*
587  *
588  * User-callable function to allocate all necessary storage space for the fft.
589  *
590  * The return value is a contiguous block of memory, allocated with malloc.  As such,
591  * It can be freed with free(), rather than a kiss_fft-specific function.
592  * */
593 kiss_fft_state *kiss_fft_alloc_twiddles(int nfft,void * mem,size_t * lenmem,  const kiss_fft_state *base)
594 {
595     kiss_fft_state *st=NULL;
596     size_t memneeded = sizeof(struct kiss_fft_state); /* twiddle factors*/
597
598     if ( lenmem==NULL ) {
599         st = ( kiss_fft_state*)KISS_FFT_MALLOC( memneeded );
600     }else{
601         if (mem != NULL && *lenmem >= memneeded)
602             st = (kiss_fft_state*)mem;
603         *lenmem = memneeded;
604     }
605     if (st) {
606         celt_int16 *bitrev;
607         kiss_twiddle_cpx *twiddles;
608
609         st->nfft=nfft;
610 #ifndef FIXED_POINT
611         st->scale = 1./nfft;
612 #endif
613         if (base != NULL)
614         {
615            st->twiddles = base->twiddles;
616            st->shift = 0;
617            while (nfft<<st->shift != base->nfft && st->shift < 32)
618               st->shift++;
619            /* FIXME: Report error and do proper cleanup */
620            if (st->shift>=32)
621               return NULL;
622         } else {
623            st->twiddles = twiddles = (kiss_twiddle_cpx*)KISS_FFT_MALLOC(sizeof(kiss_twiddle_cpx)*nfft);
624            compute_twiddles(twiddles, nfft);
625            st->shift = -1;
626         }
627         if (!kf_factor(nfft,st->factors))
628         {
629            kiss_fft_free(st);
630            return NULL;
631         }
632
633         /* bitrev */
634         st->bitrev = bitrev = (celt_int16*)KISS_FFT_MALLOC(sizeof(celt_int16)*nfft);
635         compute_bitrev_table(0, bitrev, 1,1, st->factors,st);
636     }
637     return st;
638 }
639
640 kiss_fft_state *kiss_fft_alloc(int nfft,void * mem,size_t * lenmem )
641 {
642    return kiss_fft_alloc_twiddles(nfft, mem, lenmem, NULL);
643 }
644
645 void kiss_fft_free(const kiss_fft_state *cfg)
646 {
647    if (cfg)
648    {
649       celt_free((celt_int16*)cfg->bitrev);
650       if (cfg->shift < 0)
651          celt_free((kiss_twiddle_cpx*)cfg->twiddles);
652       celt_free((kiss_fft_state*)cfg);
653    }
654 }
655
656 #endif /* CUSTOM_MODES */
657
658 static void kiss_fft_stride(const kiss_fft_state *st,const kiss_fft_cpx *fin,kiss_fft_cpx *fout,int in_stride)
659 {
660     int i;
661     celt_assert2 (fin != fout, "In-place FFT not supported");
662     /* Bit-reverse the input */
663     for (i=0;i<st->nfft;i++)
664     {
665        fout[st->bitrev[i]] = fin[i];
666 #ifndef FIXED_POINT
667        fout[st->bitrev[i]].r *= st->scale;
668        fout[st->bitrev[i]].i *= st->scale;
669 #endif
670     }
671     kf_work( fout, fin, 1,in_stride, st->factors,st, 1, 1);
672 }
673
674 void kiss_fft(const kiss_fft_state *cfg,const kiss_fft_cpx *fin,kiss_fft_cpx *fout)
675 {
676     kiss_fft_stride(cfg,fin,fout,1);
677 }
678
679 static void kiss_ifft_stride(const kiss_fft_state *st,const kiss_fft_cpx *fin,kiss_fft_cpx *fout,int in_stride)
680 {
681    int i;
682    celt_assert2 (fin != fout, "In-place FFT not supported");
683    /* Bit-reverse the input */
684    for (i=0;i<st->nfft;i++)
685       fout[st->bitrev[i]] = fin[i];
686    ki_work( fout, fin, 1,in_stride, st->factors,st, 1, 1);
687 }
688
689 void kiss_ifft(const kiss_fft_state *cfg,const kiss_fft_cpx *fin,kiss_fft_cpx *fout)
690 {
691    kiss_ifft_stride(cfg,fin,fout,1);
692 }
693