← Back to branch summary

~noskcaj/ubuntu/saucy/sflphone/merge-1.2.3-2

« back to all changes in this revision

Viewing changes to daemon/libs/pjproject-2.0.1/third_party/speex/libspeex/smallft.c

  • Committer: Jackson Doak
  • Date: 2013-07-10 21:04:46 UTC
  • mfrom: (20.1.3 sid)
  • Revision ID: noskcaj@ubuntu.com-20130710210446-y8f587vza807icr9
https://launchpad.net/bugs/1196978
Properly merged from upstream.

Show diffs side-by-side

added added

removed removed

Lines of Context:
daemon/libs/pjproject-2.0.1/third_party/speex/libspeex/smallft.c
 
1
/********************************************************************
 
2
 *                                                                  *
 
3
 * THIS FILE IS PART OF THE OggVorbis SOFTWARE CODEC SOURCE CODE.   *
 
4
 * USE, DISTRIBUTION AND REPRODUCTION OF THIS LIBRARY SOURCE IS     *
 
5
 * GOVERNED BY A BSD-STYLE SOURCE LICENSE INCLUDED WITH THIS SOURCE *
 
6
 * IN 'COPYING'. PLEASE READ THESE TERMS BEFORE DISTRIBUTING.       *
 
7
 *                                                                  *
 
8
 * THE OggVorbis SOURCE CODE IS (C) COPYRIGHT 1994-2001             *
 
9
 * by the XIPHOPHORUS Company http://www.xiph.org/                  *
 
10
 *                                                                  *
 
11
 ********************************************************************
 
12
 
 
13
 function: *unnormalized* fft transform
 
14
 last mod: $Id: smallft.c 2002 2008-06-10 14:09:37Z nanang $
 
15
 
 
16
 ********************************************************************/
 
17
 
 
18
/* FFT implementation from OggSquish, minus cosine transforms,
 
19
 * minus all but radix 2/4 case.  In Vorbis we only need this
 
20
 * cut-down version.
 
21
 *
 
22
 * To do more than just power-of-two sized vectors, see the full
 
23
 * version I wrote for NetLib.
 
24
 *
 
25
 * Note that the packing is a little strange; rather than the FFT r/i
 
26
 * packing following R_0, I_n, R_1, I_1, R_2, I_2 ... R_n-1, I_n-1,
 
27
 * it follows R_0, R_1, I_1, R_2, I_2 ... R_n-1, I_n-1, I_n like the
 
28
 * FORTRAN version
 
29
 */
 
30
 
 
31
#ifdef HAVE_CONFIG_H
 
32
#include "config.h"
 
33
#endif
 
34
 
 
35
#include <math.h>
 
36
#include "smallft.h"
 
37
#include "arch.h"
 
38
#include "os_support.h"
 
39
 
 
40
static void drfti1(int n, float *wa, int *ifac){
 
41
  static int ntryh[4] = { 4,2,3,5 };
 
42
  static float tpi = 6.28318530717958648f;
 
43
  float arg,argh,argld,fi;
 
44
  int ntry=0,i,j=-1;
 
45
  int k1, l1, l2, ib;
 
46
  int ld, ii, ip, is, nq, nr;
 
47
  int ido, ipm, nfm1;
 
48
  int nl=n;
 
49
  int nf=0;
 
50
 
 
51
 L101:
 
52
  j++;
 
53
  if (j < 4)
 
54
    ntry=ntryh[j];
 
55
  else
 
56
    ntry+=2;
 
57
 
 
58
 L104:
 
59
  nq=nl/ntry;
 
60
  nr=nl-ntry*nq;
 
61
  if (nr!=0) goto L101;
 
62
 
 
63
  nf++;
 
64
  ifac[nf+1]=ntry;
 
65
  nl=nq;
 
66
  if(ntry!=2)goto L107;
 
67
  if(nf==1)goto L107;
 
68
 
 
69
  for (i=1;i<nf;i++){
 
70
    ib=nf-i+1;
 
71
    ifac[ib+1]=ifac[ib];
 
72
  }
 
73
  ifac[2] = 2;
 
74
 
 
75
 L107:
 
76
  if(nl!=1)goto L104;
 
77
  ifac[0]=n;
 
78
  ifac[1]=nf;
 
79
  argh=tpi/n;
 
80
  is=0;
 
81
  nfm1=nf-1;
 
82
  l1=1;
 
83
 
 
84
  if(nfm1==0)return;
 
85
 
 
86
  for (k1=0;k1<nfm1;k1++){
 
87
    ip=ifac[k1+2];
 
88
    ld=0;
 
89
    l2=l1*ip;
 
90
    ido=n/l2;
 
91
    ipm=ip-1;
 
92
 
 
93
    for (j=0;j<ipm;j++){
 
94
      ld+=l1;
 
95
      i=is;
 
96
      argld=(float)ld*argh;
 
97
      fi=0.f;
 
98
      for (ii=2;ii<ido;ii+=2){
 
99
        fi+=1.f;
 
100
        arg=fi*argld;
 
101
        wa[i++]=cos(arg);
 
102
        wa[i++]=sin(arg);
 
103
      }
 
104
      is+=ido;
 
105
    }
 
106
    l1=l2;
 
107
  }
 
108
}
 
109
 
 
110
static void fdrffti(int n, float *wsave, int *ifac){
 
111
 
 
112
  if (n == 1) return;
 
113
  drfti1(n, wsave+n, ifac);
 
114
}
 
115
 
 
116
static void dradf2(int ido,int l1,float *cc,float *ch,float *wa1){
 
117
  int i,k;
 
118
  float ti2,tr2;
 
119
  int t0,t1,t2,t3,t4,t5,t6;
 
120
 
 
121
  t1=0;
 
122
  t0=(t2=l1*ido);
 
123
  t3=ido<<1;
 
124
  for(k=0;k<l1;k++){
 
125
    ch[t1<<1]=cc[t1]+cc[t2];
 
126
    ch[(t1<<1)+t3-1]=cc[t1]-cc[t2];
 
127
    t1+=ido;
 
128
    t2+=ido;
 
129
  }
 
130
 
 
131
  if(ido<2)return;
 
132
  if(ido==2)goto L105;
 
133
 
 
134
  t1=0;
 
135
  t2=t0;
 
136
  for(k=0;k<l1;k++){
 
137
    t3=t2;
 
138
    t4=(t1<<1)+(ido<<1);
 
139
    t5=t1;
 
140
    t6=t1+t1;
 
141
    for(i=2;i<ido;i+=2){
 
142
      t3+=2;
 
143
      t4-=2;
 
144
      t5+=2;
 
145
      t6+=2;
 
146
      tr2=wa1[i-2]*cc[t3-1]+wa1[i-1]*cc[t3];
 
147
      ti2=wa1[i-2]*cc[t3]-wa1[i-1]*cc[t3-1];
 
148
      ch[t6]=cc[t5]+ti2;
 
149
      ch[t4]=ti2-cc[t5];
 
150
      ch[t6-1]=cc[t5-1]+tr2;
 
151
      ch[t4-1]=cc[t5-1]-tr2;
 
152
    }
 
153
    t1+=ido;
 
154
    t2+=ido;
 
155
  }
 
156
 
 
157
  if(ido%2==1)return;
 
158
 
 
159
 L105:
 
160
  t3=(t2=(t1=ido)-1);
 
161
  t2+=t0;
 
162
  for(k=0;k<l1;k++){
 
163
    ch[t1]=-cc[t2];
 
164
    ch[t1-1]=cc[t3];
 
165
    t1+=ido<<1;
 
166
    t2+=ido;
 
167
    t3+=ido;
 
168
  }
 
169
}
 
170
 
 
171
static void dradf4(int ido,int l1,float *cc,float *ch,float *wa1,
 
172
            float *wa2,float *wa3){
 
173
  static float hsqt2 = .70710678118654752f;
 
174
  int i,k,t0,t1,t2,t3,t4,t5,t6;
 
175
  float ci2,ci3,ci4,cr2,cr3,cr4,ti1,ti2,ti3,ti4,tr1,tr2,tr3,tr4;
 
176
  t0=l1*ido;
 
177
 
 
178
  t1=t0;
 
179
  t4=t1<<1;
 
180
  t2=t1+(t1<<1);
 
181
  t3=0;
 
182
 
 
183
  for(k=0;k<l1;k++){
 
184
    tr1=cc[t1]+cc[t2];
 
185
    tr2=cc[t3]+cc[t4];
 
186
 
 
187
    ch[t5=t3<<2]=tr1+tr2;
 
188
    ch[(ido<<2)+t5-1]=tr2-tr1;
 
189
    ch[(t5+=(ido<<1))-1]=cc[t3]-cc[t4];
 
190
    ch[t5]=cc[t2]-cc[t1];
 
191
 
 
192
    t1+=ido;
 
193
    t2+=ido;
 
194
    t3+=ido;
 
195
    t4+=ido;
 
196
  }
 
197
 
 
198
  if(ido<2)return;
 
199
  if(ido==2)goto L105;
 
200
 
 
201
 
 
202
  t1=0;
 
203
  for(k=0;k<l1;k++){
 
204
    t2=t1;
 
205
    t4=t1<<2;
 
206
    t5=(t6=ido<<1)+t4;
 
207
    for(i=2;i<ido;i+=2){
 
208
      t3=(t2+=2);
 
209
      t4+=2;
 
210
      t5-=2;
 
211
 
 
212
      t3+=t0;
 
213
      cr2=wa1[i-2]*cc[t3-1]+wa1[i-1]*cc[t3];
 
214
      ci2=wa1[i-2]*cc[t3]-wa1[i-1]*cc[t3-1];
 
215
      t3+=t0;
 
216
      cr3=wa2[i-2]*cc[t3-1]+wa2[i-1]*cc[t3];
 
217
      ci3=wa2[i-2]*cc[t3]-wa2[i-1]*cc[t3-1];
 
218
      t3+=t0;
 
219
      cr4=wa3[i-2]*cc[t3-1]+wa3[i-1]*cc[t3];
 
220
      ci4=wa3[i-2]*cc[t3]-wa3[i-1]*cc[t3-1];
 
221
 
 
222
      tr1=cr2+cr4;
 
223
      tr4=cr4-cr2;
 
224
      ti1=ci2+ci4;
 
225
      ti4=ci2-ci4;
 
226
 
 
227
      ti2=cc[t2]+ci3;
 
228
      ti3=cc[t2]-ci3;
 
229
      tr2=cc[t2-1]+cr3;
 
230
      tr3=cc[t2-1]-cr3;
 
231
 
 
232
      ch[t4-1]=tr1+tr2;
 
233
      ch[t4]=ti1+ti2;
 
234
 
 
235
      ch[t5-1]=tr3-ti4;
 
236
      ch[t5]=tr4-ti3;
 
237
 
 
238
      ch[t4+t6-1]=ti4+tr3;
 
239
      ch[t4+t6]=tr4+ti3;
 
240
 
 
241
      ch[t5+t6-1]=tr2-tr1;
 
242
      ch[t5+t6]=ti1-ti2;
 
243
    }
 
244
    t1+=ido;
 
245
  }
 
246
  if(ido&1)return;
 
247
 
 
248
 L105:
 
249
 
 
250
  t2=(t1=t0+ido-1)+(t0<<1);
 
251
  t3=ido<<2;
 
252
  t4=ido;
 
253
  t5=ido<<1;
 
254
  t6=ido;
 
255
 
 
256
  for(k=0;k<l1;k++){
 
257
    ti1=-hsqt2*(cc[t1]+cc[t2]);
 
258
    tr1=hsqt2*(cc[t1]-cc[t2]);
 
259
 
 
260
    ch[t4-1]=tr1+cc[t6-1];
 
261
    ch[t4+t5-1]=cc[t6-1]-tr1;
 
262
 
 
263
    ch[t4]=ti1-cc[t1+t0];
 
264
    ch[t4+t5]=ti1+cc[t1+t0];
 
265
 
 
266
    t1+=ido;
 
267
    t2+=ido;
 
268
    t4+=t3;
 
269
    t6+=ido;
 
270
  }
 
271
}
 
272
 
 
273
static void dradfg(int ido,int ip,int l1,int idl1,float *cc,float *c1,
 
274
                          float *c2,float *ch,float *ch2,float *wa){
 
275
 
 
276
  static float tpi=6.283185307179586f;
 
277
  int idij,ipph,i,j,k,l,ic,ik,is;
 
278
  int t0,t1,t2,t3,t4,t5,t6,t7,t8,t9,t10;
 
279
  float dc2,ai1,ai2,ar1,ar2,ds2;
 
280
  int nbd;
 
281
  float dcp,arg,dsp,ar1h,ar2h;
 
282
  int idp2,ipp2;
 
283
 
 
284
  arg=tpi/(float)ip;
 
285
  dcp=cos(arg);
 
286
  dsp=sin(arg);
 
287
  ipph=(ip+1)>>1;
 
288
  ipp2=ip;
 
289
  idp2=ido;
 
290
  nbd=(ido-1)>>1;
 
291
  t0=l1*ido;
 
292
  t10=ip*ido;
 
293
 
 
294
  if(ido==1)goto L119;
 
295
  for(ik=0;ik<idl1;ik++)ch2[ik]=c2[ik];
 
296
 
 
297
  t1=0;
 
298
  for(j=1;j<ip;j++){
 
299
    t1+=t0;
 
300
    t2=t1;
 
301
    for(k=0;k<l1;k++){
 
302
      ch[t2]=c1[t2];
 
303
      t2+=ido;
 
304
    }
 
305
  }
 
306
 
 
307
  is=-ido;
 
308
  t1=0;
 
309
  if(nbd>l1){
 
310
    for(j=1;j<ip;j++){
 
311
      t1+=t0;
 
312
      is+=ido;
 
313
      t2= -ido+t1;
 
314
      for(k=0;k<l1;k++){
 
315
        idij=is-1;
 
316
        t2+=ido;
 
317
        t3=t2;
 
318
        for(i=2;i<ido;i+=2){
 
319
          idij+=2;
 
320
          t3+=2;
 
321
          ch[t3-1]=wa[idij-1]*c1[t3-1]+wa[idij]*c1[t3];
 
322
          ch[t3]=wa[idij-1]*c1[t3]-wa[idij]*c1[t3-1];
 
323
        }
 
324
      }
 
325
    }
 
326
  }else{
 
327
 
 
328
    for(j=1;j<ip;j++){
 
329
      is+=ido;
 
330
      idij=is-1;
 
331
      t1+=t0;
 
332
      t2=t1;
 
333
      for(i=2;i<ido;i+=2){
 
334
        idij+=2;
 
335
        t2+=2;
 
336
        t3=t2;
 
337
        for(k=0;k<l1;k++){
 
338
          ch[t3-1]=wa[idij-1]*c1[t3-1]+wa[idij]*c1[t3];
 
339
          ch[t3]=wa[idij-1]*c1[t3]-wa[idij]*c1[t3-1];
 
340
          t3+=ido;
 
341
        }
 
342
      }
 
343
    }
 
344
  }
 
345
 
 
346
  t1=0;
 
347
  t2=ipp2*t0;
 
348
  if(nbd<l1){
 
349
    for(j=1;j<ipph;j++){
 
350
      t1+=t0;
 
351
      t2-=t0;
 
352
      t3=t1;
 
353
      t4=t2;
 
354
      for(i=2;i<ido;i+=2){
 
355
        t3+=2;
 
356
        t4+=2;
 
357
        t5=t3-ido;
 
358
        t6=t4-ido;
 
359
        for(k=0;k<l1;k++){
 
360
          t5+=ido;
 
361
          t6+=ido;
 
362
          c1[t5-1]=ch[t5-1]+ch[t6-1];
 
363
          c1[t6-1]=ch[t5]-ch[t6];
 
364
          c1[t5]=ch[t5]+ch[t6];
 
365
          c1[t6]=ch[t6-1]-ch[t5-1];
 
366
        }
 
367
      }
 
368
    }
 
369
  }else{
 
370
    for(j=1;j<ipph;j++){
 
371
      t1+=t0;
 
372
      t2-=t0;
 
373
      t3=t1;
 
374
      t4=t2;
 
375
      for(k=0;k<l1;k++){
 
376
        t5=t3;
 
377
        t6=t4;
 
378
        for(i=2;i<ido;i+=2){
 
379
          t5+=2;
 
380
          t6+=2;
 
381
          c1[t5-1]=ch[t5-1]+ch[t6-1];
 
382
          c1[t6-1]=ch[t5]-ch[t6];
 
383
          c1[t5]=ch[t5]+ch[t6];
 
384
          c1[t6]=ch[t6-1]-ch[t5-1];
 
385
        }
 
386
        t3+=ido;
 
387
        t4+=ido;
 
388
      }
 
389
    }
 
390
  }
 
391
 
 
392
L119:
 
393
  for(ik=0;ik<idl1;ik++)c2[ik]=ch2[ik];
 
394
 
 
395
  t1=0;
 
396
  t2=ipp2*idl1;
 
397
  for(j=1;j<ipph;j++){
 
398
    t1+=t0;
 
399
    t2-=t0;
 
400
    t3=t1-ido;
 
401
    t4=t2-ido;
 
402
    for(k=0;k<l1;k++){
 
403
      t3+=ido;
 
404
      t4+=ido;
 
405
      c1[t3]=ch[t3]+ch[t4];
 
406
      c1[t4]=ch[t4]-ch[t3];
 
407
    }
 
408
  }
 
409
 
 
410
  ar1=1.f;
 
411
  ai1=0.f;
 
412
  t1=0;
 
413
  t2=ipp2*idl1;
 
414
  t3=(ip-1)*idl1;
 
415
  for(l=1;l<ipph;l++){
 
416
    t1+=idl1;
 
417
    t2-=idl1;
 
418
    ar1h=dcp*ar1-dsp*ai1;
 
419
    ai1=dcp*ai1+dsp*ar1;
 
420
    ar1=ar1h;
 
421
    t4=t1;
 
422
    t5=t2;
 
423
    t6=t3;
 
424
    t7=idl1;
 
425
 
 
426
    for(ik=0;ik<idl1;ik++){
 
427
      ch2[t4++]=c2[ik]+ar1*c2[t7++];
 
428
      ch2[t5++]=ai1*c2[t6++];
 
429
    }
 
430
 
 
431
    dc2=ar1;
 
432
    ds2=ai1;
 
433
    ar2=ar1;
 
434
    ai2=ai1;
 
435
 
 
436
    t4=idl1;
 
437
    t5=(ipp2-1)*idl1;
 
438
    for(j=2;j<ipph;j++){
 
439
      t4+=idl1;
 
440
      t5-=idl1;
 
441
 
 
442
      ar2h=dc2*ar2-ds2*ai2;
 
443
      ai2=dc2*ai2+ds2*ar2;
 
444
      ar2=ar2h;
 
445
 
 
446
      t6=t1;
 
447
      t7=t2;
 
448
      t8=t4;
 
449
      t9=t5;
 
450
      for(ik=0;ik<idl1;ik++){
 
451
        ch2[t6++]+=ar2*c2[t8++];
 
452
        ch2[t7++]+=ai2*c2[t9++];
 
453
      }
 
454
    }
 
455
  }
 
456
 
 
457
  t1=0;
 
458
  for(j=1;j<ipph;j++){
 
459
    t1+=idl1;
 
460
    t2=t1;
 
461
    for(ik=0;ik<idl1;ik++)ch2[ik]+=c2[t2++];
 
462
  }
 
463
 
 
464
  if(ido<l1)goto L132;
 
465
 
 
466
  t1=0;
 
467
  t2=0;
 
468
  for(k=0;k<l1;k++){
 
469
    t3=t1;
 
470
    t4=t2;
 
471
    for(i=0;i<ido;i++)cc[t4++]=ch[t3++];
 
472
    t1+=ido;
 
473
    t2+=t10;
 
474
  }
 
475
 
 
476
  goto L135;
 
477
 
 
478
 L132:
 
479
  for(i=0;i<ido;i++){
 
480
    t1=i;
 
481
    t2=i;
 
482
    for(k=0;k<l1;k++){
 
483
      cc[t2]=ch[t1];
 
484
      t1+=ido;
 
485
      t2+=t10;
 
486
    }
 
487
  }
 
488
 
 
489
 L135:
 
490
  t1=0;
 
491
  t2=ido<<1;
 
492
  t3=0;
 
493
  t4=ipp2*t0;
 
494
  for(j=1;j<ipph;j++){
 
495
 
 
496
    t1+=t2;
 
497
    t3+=t0;
 
498
    t4-=t0;
 
499
 
 
500
    t5=t1;
 
501
    t6=t3;
 
502
    t7=t4;
 
503
 
 
504
    for(k=0;k<l1;k++){
 
505
      cc[t5-1]=ch[t6];
 
506
      cc[t5]=ch[t7];
 
507
      t5+=t10;
 
508
      t6+=ido;
 
509
      t7+=ido;
 
510
    }
 
511
  }
 
512
 
 
513
  if(ido==1)return;
 
514
  if(nbd<l1)goto L141;
 
515
 
 
516
  t1=-ido;
 
517
  t3=0;
 
518
  t4=0;
 
519
  t5=ipp2*t0;
 
520
  for(j=1;j<ipph;j++){
 
521
    t1+=t2;
 
522
    t3+=t2;
 
523
    t4+=t0;
 
524
    t5-=t0;
 
525
    t6=t1;
 
526
    t7=t3;
 
527
    t8=t4;
 
528
    t9=t5;
 
529
    for(k=0;k<l1;k++){
 
530
      for(i=2;i<ido;i+=2){
 
531
        ic=idp2-i;
 
532
        cc[i+t7-1]=ch[i+t8-1]+ch[i+t9-1];
 
533
        cc[ic+t6-1]=ch[i+t8-1]-ch[i+t9-1];
 
534
        cc[i+t7]=ch[i+t8]+ch[i+t9];
 
535
        cc[ic+t6]=ch[i+t9]-ch[i+t8];
 
536
      }
 
537
      t6+=t10;
 
538
      t7+=t10;
 
539
      t8+=ido;
 
540
      t9+=ido;
 
541
    }
 
542
  }
 
543
  return;
 
544
 
 
545
 L141:
 
546
 
 
547
  t1=-ido;
 
548
  t3=0;
 
549
  t4=0;
 
550
  t5=ipp2*t0;
 
551
  for(j=1;j<ipph;j++){
 
552
    t1+=t2;
 
553
    t3+=t2;
 
554
    t4+=t0;
 
555
    t5-=t0;
 
556
    for(i=2;i<ido;i+=2){
 
557
      t6=idp2+t1-i;
 
558
      t7=i+t3;
 
559
      t8=i+t4;
 
560
      t9=i+t5;
 
561
      for(k=0;k<l1;k++){
 
562
        cc[t7-1]=ch[t8-1]+ch[t9-1];
 
563
        cc[t6-1]=ch[t8-1]-ch[t9-1];
 
564
        cc[t7]=ch[t8]+ch[t9];
 
565
        cc[t6]=ch[t9]-ch[t8];
 
566
        t6+=t10;
 
567
        t7+=t10;
 
568
        t8+=ido;
 
569
        t9+=ido;
 
570
      }
 
571
    }
 
572
  }
 
573
}
 
574
 
 
575
static void drftf1(int n,float *c,float *ch,float *wa,int *ifac){
 
576
  int i,k1,l1,l2;
 
577
  int na,kh,nf;
 
578
  int ip,iw,ido,idl1,ix2,ix3;
 
579
 
 
580
  nf=ifac[1];
 
581
  na=1;
 
582
  l2=n;
 
583
  iw=n;
 
584
 
 
585
  for(k1=0;k1<nf;k1++){
 
586
    kh=nf-k1;
 
587
    ip=ifac[kh+1];
 
588
    l1=l2/ip;
 
589
    ido=n/l2;
 
590
    idl1=ido*l1;
 
591
    iw-=(ip-1)*ido;
 
592
    na=1-na;
 
593
 
 
594
    if(ip!=4)goto L102;
 
595
 
 
596
    ix2=iw+ido;
 
597
    ix3=ix2+ido;
 
598
    if(na!=0)
 
599
      dradf4(ido,l1,ch,c,wa+iw-1,wa+ix2-1,wa+ix3-1);
 
600
    else
 
601
      dradf4(ido,l1,c,ch,wa+iw-1,wa+ix2-1,wa+ix3-1);
 
602
    goto L110;
 
603
 
 
604
 L102:
 
605
    if(ip!=2)goto L104;
 
606
    if(na!=0)goto L103;
 
607
 
 
608
    dradf2(ido,l1,c,ch,wa+iw-1);
 
609
    goto L110;
 
610
 
 
611
  L103:
 
612
    dradf2(ido,l1,ch,c,wa+iw-1);
 
613
    goto L110;
 
614
 
 
615
  L104:
 
616
    if(ido==1)na=1-na;
 
617
    if(na!=0)goto L109;
 
618
 
 
619
    dradfg(ido,ip,l1,idl1,c,c,c,ch,ch,wa+iw-1);
 
620
    na=1;
 
621
    goto L110;
 
622
 
 
623
  L109:
 
624
    dradfg(ido,ip,l1,idl1,ch,ch,ch,c,c,wa+iw-1);
 
625
    na=0;
 
626
 
 
627
  L110:
 
628
    l2=l1;
 
629
  }
 
630
 
 
631
  if(na==1)return;
 
632
 
 
633
  for(i=0;i<n;i++)c[i]=ch[i];
 
634
}
 
635
 
 
636
static void dradb2(int ido,int l1,float *cc,float *ch,float *wa1){
 
637
  int i,k,t0,t1,t2,t3,t4,t5,t6;
 
638
  float ti2,tr2;
 
639
 
 
640
  t0=l1*ido;
 
641
 
 
642
  t1=0;
 
643
  t2=0;
 
644
  t3=(ido<<1)-1;
 
645
  for(k=0;k<l1;k++){
 
646
    ch[t1]=cc[t2]+cc[t3+t2];
 
647
    ch[t1+t0]=cc[t2]-cc[t3+t2];
 
648
    t2=(t1+=ido)<<1;
 
649
  }
 
650
 
 
651
  if(ido<2)return;
 
652
  if(ido==2)goto L105;
 
653
 
 
654
  t1=0;
 
655
  t2=0;
 
656
  for(k=0;k<l1;k++){
 
657
    t3=t1;
 
658
    t5=(t4=t2)+(ido<<1);
 
659
    t6=t0+t1;
 
660
    for(i=2;i<ido;i+=2){
 
661
      t3+=2;
 
662
      t4+=2;
 
663
      t5-=2;
 
664
      t6+=2;
 
665
      ch[t3-1]=cc[t4-1]+cc[t5-1];
 
666
      tr2=cc[t4-1]-cc[t5-1];
 
667
      ch[t3]=cc[t4]-cc[t5];
 
668
      ti2=cc[t4]+cc[t5];
 
669
      ch[t6-1]=wa1[i-2]*tr2-wa1[i-1]*ti2;
 
670
      ch[t6]=wa1[i-2]*ti2+wa1[i-1]*tr2;
 
671
    }
 
672
    t2=(t1+=ido)<<1;
 
673
  }
 
674
 
 
675
  if(ido%2==1)return;
 
676
 
 
677
L105:
 
678
  t1=ido-1;
 
679
  t2=ido-1;
 
680
  for(k=0;k<l1;k++){
 
681
    ch[t1]=cc[t2]+cc[t2];
 
682
    ch[t1+t0]=-(cc[t2+1]+cc[t2+1]);
 
683
    t1+=ido;
 
684
    t2+=ido<<1;
 
685
  }
 
686
}
 
687
 
 
688
static void dradb3(int ido,int l1,float *cc,float *ch,float *wa1,
 
689
                          float *wa2){
 
690
  static float taur = -.5f;
 
691
  static float taui = .8660254037844386f;
 
692
  int i,k,t0,t1,t2,t3,t4,t5,t6,t7,t8,t9,t10;
 
693
  float ci2,ci3,di2,di3,cr2,cr3,dr2,dr3,ti2,tr2;
 
694
  t0=l1*ido;
 
695
 
 
696
  t1=0;
 
697
  t2=t0<<1;
 
698
  t3=ido<<1;
 
699
  t4=ido+(ido<<1);
 
700
  t5=0;
 
701
  for(k=0;k<l1;k++){
 
702
    tr2=cc[t3-1]+cc[t3-1];
 
703
    cr2=cc[t5]+(taur*tr2);
 
704
    ch[t1]=cc[t5]+tr2;
 
705
    ci3=taui*(cc[t3]+cc[t3]);
 
706
    ch[t1+t0]=cr2-ci3;
 
707
    ch[t1+t2]=cr2+ci3;
 
708
    t1+=ido;
 
709
    t3+=t4;
 
710
    t5+=t4;
 
711
  }
 
712
 
 
713
  if(ido==1)return;
 
714
 
 
715
  t1=0;
 
716
  t3=ido<<1;
 
717
  for(k=0;k<l1;k++){
 
718
    t7=t1+(t1<<1);
 
719
    t6=(t5=t7+t3);
 
720
    t8=t1;
 
721
    t10=(t9=t1+t0)+t0;
 
722
 
 
723
    for(i=2;i<ido;i+=2){
 
724
      t5+=2;
 
725
      t6-=2;
 
726
      t7+=2;
 
727
      t8+=2;
 
728
      t9+=2;
 
729
      t10+=2;
 
730
      tr2=cc[t5-1]+cc[t6-1];
 
731
      cr2=cc[t7-1]+(taur*tr2);
 
732
      ch[t8-1]=cc[t7-1]+tr2;
 
733
      ti2=cc[t5]-cc[t6];
 
734
      ci2=cc[t7]+(taur*ti2);
 
735
      ch[t8]=cc[t7]+ti2;
 
736
      cr3=taui*(cc[t5-1]-cc[t6-1]);
 
737
      ci3=taui*(cc[t5]+cc[t6]);
 
738
      dr2=cr2-ci3;
 
739
      dr3=cr2+ci3;
 
740
      di2=ci2+cr3;
 
741
      di3=ci2-cr3;
 
742
      ch[t9-1]=wa1[i-2]*dr2-wa1[i-1]*di2;
 
743
      ch[t9]=wa1[i-2]*di2+wa1[i-1]*dr2;
 
744
      ch[t10-1]=wa2[i-2]*dr3-wa2[i-1]*di3;
 
745
      ch[t10]=wa2[i-2]*di3+wa2[i-1]*dr3;
 
746
    }
 
747
    t1+=ido;
 
748
  }
 
749
}
 
750
 
 
751
static void dradb4(int ido,int l1,float *cc,float *ch,float *wa1,
 
752
                          float *wa2,float *wa3){
 
753
  static float sqrt2=1.414213562373095f;
 
754
  int i,k,t0,t1,t2,t3,t4,t5,t6,t7,t8;
 
755
  float ci2,ci3,ci4,cr2,cr3,cr4,ti1,ti2,ti3,ti4,tr1,tr2,tr3,tr4;
 
756
  t0=l1*ido;
 
757
 
 
758
  t1=0;
 
759
  t2=ido<<2;
 
760
  t3=0;
 
761
  t6=ido<<1;
 
762
  for(k=0;k<l1;k++){
 
763
    t4=t3+t6;
 
764
    t5=t1;
 
765
    tr3=cc[t4-1]+cc[t4-1];
 
766
    tr4=cc[t4]+cc[t4];
 
767
    tr1=cc[t3]-cc[(t4+=t6)-1];
 
768
    tr2=cc[t3]+cc[t4-1];
 
769
    ch[t5]=tr2+tr3;
 
770
    ch[t5+=t0]=tr1-tr4;
 
771
    ch[t5+=t0]=tr2-tr3;
 
772
    ch[t5+=t0]=tr1+tr4;
 
773
    t1+=ido;
 
774
    t3+=t2;
 
775
  }
 
776
 
 
777
  if(ido<2)return;
 
778
  if(ido==2)goto L105;
 
779
 
 
780
  t1=0;
 
781
  for(k=0;k<l1;k++){
 
782
    t5=(t4=(t3=(t2=t1<<2)+t6))+t6;
 
783
    t7=t1;
 
784
    for(i=2;i<ido;i+=2){
 
785
      t2+=2;
 
786
      t3+=2;
 
787
      t4-=2;
 
788
      t5-=2;
 
789
      t7+=2;
 
790
      ti1=cc[t2]+cc[t5];
 
791
      ti2=cc[t2]-cc[t5];
 
792
      ti3=cc[t3]-cc[t4];
 
793
      tr4=cc[t3]+cc[t4];
 
794
      tr1=cc[t2-1]-cc[t5-1];
 
795
      tr2=cc[t2-1]+cc[t5-1];
 
796
      ti4=cc[t3-1]-cc[t4-1];
 
797
      tr3=cc[t3-1]+cc[t4-1];
 
798
      ch[t7-1]=tr2+tr3;
 
799
      cr3=tr2-tr3;
 
800
      ch[t7]=ti2+ti3;
 
801
      ci3=ti2-ti3;
 
802
      cr2=tr1-tr4;
 
803
      cr4=tr1+tr4;
 
804
      ci2=ti1+ti4;
 
805
      ci4=ti1-ti4;
 
806
 
 
807
      ch[(t8=t7+t0)-1]=wa1[i-2]*cr2-wa1[i-1]*ci2;
 
808
      ch[t8]=wa1[i-2]*ci2+wa1[i-1]*cr2;
 
809
      ch[(t8+=t0)-1]=wa2[i-2]*cr3-wa2[i-1]*ci3;
 
810
      ch[t8]=wa2[i-2]*ci3+wa2[i-1]*cr3;
 
811
      ch[(t8+=t0)-1]=wa3[i-2]*cr4-wa3[i-1]*ci4;
 
812
      ch[t8]=wa3[i-2]*ci4+wa3[i-1]*cr4;
 
813
    }
 
814
    t1+=ido;
 
815
  }
 
816
 
 
817
  if(ido%2 == 1)return;
 
818
 
 
819
 L105:
 
820
 
 
821
  t1=ido;
 
822
  t2=ido<<2;
 
823
  t3=ido-1;
 
824
  t4=ido+(ido<<1);
 
825
  for(k=0;k<l1;k++){
 
826
    t5=t3;
 
827
    ti1=cc[t1]+cc[t4];
 
828
    ti2=cc[t4]-cc[t1];
 
829
    tr1=cc[t1-1]-cc[t4-1];
 
830
    tr2=cc[t1-1]+cc[t4-1];
 
831
    ch[t5]=tr2+tr2;
 
832
    ch[t5+=t0]=sqrt2*(tr1-ti1);
 
833
    ch[t5+=t0]=ti2+ti2;
 
834
    ch[t5+=t0]=-sqrt2*(tr1+ti1);
 
835
 
 
836
    t3+=ido;
 
837
    t1+=t2;
 
838
    t4+=t2;
 
839
  }
 
840
}
 
841
 
 
842
static void dradbg(int ido,int ip,int l1,int idl1,float *cc,float *c1,
 
843
            float *c2,float *ch,float *ch2,float *wa){
 
844
  static float tpi=6.283185307179586f;
 
845
  int idij,ipph,i,j,k,l,ik,is,t0,t1,t2,t3,t4,t5,t6,t7,t8,t9,t10,
 
846
      t11,t12;
 
847
  float dc2,ai1,ai2,ar1,ar2,ds2;
 
848
  int nbd;
 
849
  float dcp,arg,dsp,ar1h,ar2h;
 
850
  int ipp2;
 
851
 
 
852
  t10=ip*ido;
 
853
  t0=l1*ido;
 
854
  arg=tpi/(float)ip;
 
855
  dcp=cos(arg);
 
856
  dsp=sin(arg);
 
857
  nbd=(ido-1)>>1;
 
858
  ipp2=ip;
 
859
  ipph=(ip+1)>>1;
 
860
  if(ido<l1)goto L103;
 
861
 
 
862
  t1=0;
 
863
  t2=0;
 
864
  for(k=0;k<l1;k++){
 
865
    t3=t1;
 
866
    t4=t2;
 
867
    for(i=0;i<ido;i++){
 
868
      ch[t3]=cc[t4];
 
869
      t3++;
 
870
      t4++;
 
871
    }
 
872
    t1+=ido;
 
873
    t2+=t10;
 
874
  }
 
875
  goto L106;
 
876
 
 
877
 L103:
 
878
  t1=0;
 
879
  for(i=0;i<ido;i++){
 
880
    t2=t1;
 
881
    t3=t1;
 
882
    for(k=0;k<l1;k++){
 
883
      ch[t2]=cc[t3];
 
884
      t2+=ido;
 
885
      t3+=t10;
 
886
    }
 
887
    t1++;
 
888
  }
 
889
 
 
890
 L106:
 
891
  t1=0;
 
892
  t2=ipp2*t0;
 
893
  t7=(t5=ido<<1);
 
894
  for(j=1;j<ipph;j++){
 
895
    t1+=t0;
 
896
    t2-=t0;
 
897
    t3=t1;
 
898
    t4=t2;
 
899
    t6=t5;
 
900
    for(k=0;k<l1;k++){
 
901
      ch[t3]=cc[t6-1]+cc[t6-1];
 
902
      ch[t4]=cc[t6]+cc[t6];
 
903
      t3+=ido;
 
904
      t4+=ido;
 
905
      t6+=t10;
 
906
    }
 
907
    t5+=t7;
 
908
  }
 
909
 
 
910
  if (ido == 1)goto L116;
 
911
  if(nbd<l1)goto L112;
 
912
 
 
913
  t1=0;
 
914
  t2=ipp2*t0;
 
915
  t7=0;
 
916
  for(j=1;j<ipph;j++){
 
917
    t1+=t0;
 
918
    t2-=t0;
 
919
    t3=t1;
 
920
    t4=t2;
 
921
 
 
922
    t7+=(ido<<1);
 
923
    t8=t7;
 
924
    for(k=0;k<l1;k++){
 
925
      t5=t3;
 
926
      t6=t4;
 
927
      t9=t8;
 
928
      t11=t8;
 
929
      for(i=2;i<ido;i+=2){
 
930
        t5+=2;
 
931
        t6+=2;
 
932
        t9+=2;
 
933
        t11-=2;
 
934
        ch[t5-1]=cc[t9-1]+cc[t11-1];
 
935
        ch[t6-1]=cc[t9-1]-cc[t11-1];
 
936
        ch[t5]=cc[t9]-cc[t11];
 
937
        ch[t6]=cc[t9]+cc[t11];
 
938
      }
 
939
      t3+=ido;
 
940
      t4+=ido;
 
941
      t8+=t10;
 
942
    }
 
943
  }
 
944
  goto L116;
 
945
 
 
946
 L112:
 
947
  t1=0;
 
948
  t2=ipp2*t0;
 
949
  t7=0;
 
950
  for(j=1;j<ipph;j++){
 
951
    t1+=t0;
 
952
    t2-=t0;
 
953
    t3=t1;
 
954
    t4=t2;
 
955
    t7+=(ido<<1);
 
956
    t8=t7;
 
957
    t9=t7;
 
958
    for(i=2;i<ido;i+=2){
 
959
      t3+=2;
 
960
      t4+=2;
 
961
      t8+=2;
 
962
      t9-=2;
 
963
      t5=t3;
 
964
      t6=t4;
 
965
      t11=t8;
 
966
      t12=t9;
 
967
      for(k=0;k<l1;k++){
 
968
        ch[t5-1]=cc[t11-1]+cc[t12-1];
 
969
        ch[t6-1]=cc[t11-1]-cc[t12-1];
 
970
        ch[t5]=cc[t11]-cc[t12];
 
971
        ch[t6]=cc[t11]+cc[t12];
 
972
        t5+=ido;
 
973
        t6+=ido;
 
974
        t11+=t10;
 
975
        t12+=t10;
 
976
      }
 
977
    }
 
978
  }
 
979
 
 
980
L116:
 
981
  ar1=1.f;
 
982
  ai1=0.f;
 
983
  t1=0;
 
984
  t9=(t2=ipp2*idl1);
 
985
  t3=(ip-1)*idl1;
 
986
  for(l=1;l<ipph;l++){
 
987
    t1+=idl1;
 
988
    t2-=idl1;
 
989
 
 
990
    ar1h=dcp*ar1-dsp*ai1;
 
991
    ai1=dcp*ai1+dsp*ar1;
 
992
    ar1=ar1h;
 
993
    t4=t1;
 
994
    t5=t2;
 
995
    t6=0;
 
996
    t7=idl1;
 
997
    t8=t3;
 
998
    for(ik=0;ik<idl1;ik++){
 
999
      c2[t4++]=ch2[t6++]+ar1*ch2[t7++];
 
1000
      c2[t5++]=ai1*ch2[t8++];
 
1001
    }
 
1002
    dc2=ar1;
 
1003
    ds2=ai1;
 
1004
    ar2=ar1;
 
1005
    ai2=ai1;
 
1006
 
 
1007
    t6=idl1;
 
1008
    t7=t9-idl1;
 
1009
    for(j=2;j<ipph;j++){
 
1010
      t6+=idl1;
 
1011
      t7-=idl1;
 
1012
      ar2h=dc2*ar2-ds2*ai2;
 
1013
      ai2=dc2*ai2+ds2*ar2;
 
1014
      ar2=ar2h;
 
1015
      t4=t1;
 
1016
      t5=t2;
 
1017
      t11=t6;
 
1018
      t12=t7;
 
1019
      for(ik=0;ik<idl1;ik++){
 
1020
        c2[t4++]+=ar2*ch2[t11++];
 
1021
        c2[t5++]+=ai2*ch2[t12++];
 
1022
      }
 
1023
    }
 
1024
  }
 
1025
 
 
1026
  t1=0;
 
1027
  for(j=1;j<ipph;j++){
 
1028
    t1+=idl1;
 
1029
    t2=t1;
 
1030
    for(ik=0;ik<idl1;ik++)ch2[ik]+=ch2[t2++];
 
1031
  }
 
1032
 
 
1033
  t1=0;
 
1034
  t2=ipp2*t0;
 
1035
  for(j=1;j<ipph;j++){
 
1036
    t1+=t0;
 
1037
    t2-=t0;
 
1038
    t3=t1;
 
1039
    t4=t2;
 
1040
    for(k=0;k<l1;k++){
 
1041
      ch[t3]=c1[t3]-c1[t4];
 
1042
      ch[t4]=c1[t3]+c1[t4];
 
1043
      t3+=ido;
 
1044
      t4+=ido;
 
1045
    }
 
1046
  }
 
1047
 
 
1048
  if(ido==1)goto L132;
 
1049
  if(nbd<l1)goto L128;
 
1050
 
 
1051
  t1=0;
 
1052
  t2=ipp2*t0;
 
1053
  for(j=1;j<ipph;j++){
 
1054
    t1+=t0;
 
1055
    t2-=t0;
 
1056
    t3=t1;
 
1057
    t4=t2;
 
1058
    for(k=0;k<l1;k++){
 
1059
      t5=t3;
 
1060
      t6=t4;
 
1061
      for(i=2;i<ido;i+=2){
 
1062
        t5+=2;
 
1063
        t6+=2;
 
1064
        ch[t5-1]=c1[t5-1]-c1[t6];
 
1065
        ch[t6-1]=c1[t5-1]+c1[t6];
 
1066
        ch[t5]=c1[t5]+c1[t6-1];
 
1067
        ch[t6]=c1[t5]-c1[t6-1];
 
1068
      }
 
1069
      t3+=ido;
 
1070
      t4+=ido;
 
1071
    }
 
1072
  }
 
1073
  goto L132;
 
1074
 
 
1075
 L128:
 
1076
  t1=0;
 
1077
  t2=ipp2*t0;
 
1078
  for(j=1;j<ipph;j++){
 
1079
    t1+=t0;
 
1080
    t2-=t0;
 
1081
    t3=t1;
 
1082
    t4=t2;
 
1083
    for(i=2;i<ido;i+=2){
 
1084
      t3+=2;
 
1085
      t4+=2;
 
1086
      t5=t3;
 
1087
      t6=t4;
 
1088
      for(k=0;k<l1;k++){
 
1089
        ch[t5-1]=c1[t5-1]-c1[t6];
 
1090
        ch[t6-1]=c1[t5-1]+c1[t6];
 
1091
        ch[t5]=c1[t5]+c1[t6-1];
 
1092
        ch[t6]=c1[t5]-c1[t6-1];
 
1093
        t5+=ido;
 
1094
        t6+=ido;
 
1095
      }
 
1096
    }
 
1097
  }
 
1098
 
 
1099
L132:
 
1100
  if(ido==1)return;
 
1101
 
 
1102
  for(ik=0;ik<idl1;ik++)c2[ik]=ch2[ik];
 
1103
 
 
1104
  t1=0;
 
1105
  for(j=1;j<ip;j++){
 
1106
    t2=(t1+=t0);
 
1107
    for(k=0;k<l1;k++){
 
1108
      c1[t2]=ch[t2];
 
1109
      t2+=ido;
 
1110
    }
 
1111
  }
 
1112
 
 
1113
  if(nbd>l1)goto L139;
 
1114
 
 
1115
  is= -ido-1;
 
1116
  t1=0;
 
1117
  for(j=1;j<ip;j++){
 
1118
    is+=ido;
 
1119
    t1+=t0;
 
1120
    idij=is;
 
1121
    t2=t1;
 
1122
    for(i=2;i<ido;i+=2){
 
1123
      t2+=2;
 
1124
      idij+=2;
 
1125
      t3=t2;
 
1126
      for(k=0;k<l1;k++){
 
1127
        c1[t3-1]=wa[idij-1]*ch[t3-1]-wa[idij]*ch[t3];
 
1128
        c1[t3]=wa[idij-1]*ch[t3]+wa[idij]*ch[t3-1];
 
1129
        t3+=ido;
 
1130
      }
 
1131
    }
 
1132
  }
 
1133
  return;
 
1134
 
 
1135
 L139:
 
1136
  is= -ido-1;
 
1137
  t1=0;
 
1138
  for(j=1;j<ip;j++){
 
1139
    is+=ido;
 
1140
    t1+=t0;
 
1141
    t2=t1;
 
1142
    for(k=0;k<l1;k++){
 
1143
      idij=is;
 
1144
      t3=t2;
 
1145
      for(i=2;i<ido;i+=2){
 
1146
        idij+=2;
 
1147
        t3+=2;
 
1148
        c1[t3-1]=wa[idij-1]*ch[t3-1]-wa[idij]*ch[t3];
 
1149
        c1[t3]=wa[idij-1]*ch[t3]+wa[idij]*ch[t3-1];
 
1150
      }
 
1151
      t2+=ido;
 
1152
    }
 
1153
  }
 
1154
}
 
1155
 
 
1156
static void drftb1(int n, float *c, float *ch, float *wa, int *ifac){
 
1157
  int i,k1,l1,l2;
 
1158
  int na;
 
1159
  int nf,ip,iw,ix2,ix3,ido,idl1;
 
1160
 
 
1161
  nf=ifac[1];
 
1162
  na=0;
 
1163
  l1=1;
 
1164
  iw=1;
 
1165
 
 
1166
  for(k1=0;k1<nf;k1++){
 
1167
    ip=ifac[k1 + 2];
 
1168
    l2=ip*l1;
 
1169
    ido=n/l2;
 
1170
    idl1=ido*l1;
 
1171
    if(ip!=4)goto L103;
 
1172
    ix2=iw+ido;
 
1173
    ix3=ix2+ido;
 
1174
 
 
1175
    if(na!=0)
 
1176
      dradb4(ido,l1,ch,c,wa+iw-1,wa+ix2-1,wa+ix3-1);
 
1177
    else
 
1178
      dradb4(ido,l1,c,ch,wa+iw-1,wa+ix2-1,wa+ix3-1);
 
1179
    na=1-na;
 
1180
    goto L115;
 
1181
 
 
1182
  L103:
 
1183
    if(ip!=2)goto L106;
 
1184
 
 
1185
    if(na!=0)
 
1186
      dradb2(ido,l1,ch,c,wa+iw-1);
 
1187
    else
 
1188
      dradb2(ido,l1,c,ch,wa+iw-1);
 
1189
    na=1-na;
 
1190
    goto L115;
 
1191
 
 
1192
  L106:
 
1193
    if(ip!=3)goto L109;
 
1194
 
 
1195
    ix2=iw+ido;
 
1196
    if(na!=0)
 
1197
      dradb3(ido,l1,ch,c,wa+iw-1,wa+ix2-1);
 
1198
    else
 
1199
      dradb3(ido,l1,c,ch,wa+iw-1,wa+ix2-1);
 
1200
    na=1-na;
 
1201
    goto L115;
 
1202
 
 
1203
  L109:
 
1204
/*    The radix five case can be translated later..... */
 
1205
/*    if(ip!=5)goto L112;
 
1206
 
 
1207
    ix2=iw+ido;
 
1208
    ix3=ix2+ido;
 
1209
    ix4=ix3+ido;
 
1210
    if(na!=0)
 
1211
      dradb5(ido,l1,ch,c,wa+iw-1,wa+ix2-1,wa+ix3-1,wa+ix4-1);
 
1212
    else
 
1213
      dradb5(ido,l1,c,ch,wa+iw-1,wa+ix2-1,wa+ix3-1,wa+ix4-1);
 
1214
    na=1-na;
 
1215
    goto L115;
 
1216
 
 
1217
  L112:*/
 
1218
    if(na!=0)
 
1219
      dradbg(ido,ip,l1,idl1,ch,ch,ch,c,c,wa+iw-1);
 
1220
    else
 
1221
      dradbg(ido,ip,l1,idl1,c,c,c,ch,ch,wa+iw-1);
 
1222
    if(ido==1)na=1-na;
 
1223
 
 
1224
  L115:
 
1225
    l1=l2;
 
1226
    iw+=(ip-1)*ido;
 
1227
  }
 
1228
 
 
1229
  if(na==0)return;
 
1230
 
 
1231
  for(i=0;i<n;i++)c[i]=ch[i];
 
1232
}
 
1233
 
 
1234
void spx_drft_forward(struct drft_lookup *l,float *data){
 
1235
  if(l->n==1)return;
 
1236
  drftf1(l->n,data,l->trigcache,l->trigcache+l->n,l->splitcache);
 
1237
}
 
1238
 
 
1239
void spx_drft_backward(struct drft_lookup *l,float *data){
 
1240
  if (l->n==1)return;
 
1241
  drftb1(l->n,data,l->trigcache,l->trigcache+l->n,l->splitcache);
 
1242
}
 
1243
 
 
1244
void spx_drft_init(struct drft_lookup *l,int n)
 
1245
{
 
1246
  l->n=n;
 
1247
  l->trigcache=(float*)speex_alloc(3*n*sizeof(*l->trigcache));
 
1248
  l->splitcache=(int*)speex_alloc(32*sizeof(*l->splitcache));
 
1249
  fdrffti(n, l->trigcache, l->splitcache);
 
1250
}
 
1251
 
 
1252
void spx_drft_clear(struct drft_lookup *l)
 
1253
{
 
1254
  if(l)
 
1255
  {
 
1256
    if(l->trigcache)
 
1257
      speex_free(l->trigcache);
 
1258
    if(l->splitcache)
 
1259
      speex_free(l->splitcache);
 
1260
  }
 
1261
}