aboutsummaryrefslogtreecommitdiff
path: root/pd/src/d_mayer_fft.c
diff options
context:
space:
mode:
Diffstat (limited to 'pd/src/d_mayer_fft.c')
-rw-r--r--pd/src/d_mayer_fft.c268
1 files changed, 134 insertions, 134 deletions
diff --git a/pd/src/d_mayer_fft.c b/pd/src/d_mayer_fft.c
index ea884018..0e2242b7 100644
--- a/pd/src/d_mayer_fft.c
+++ b/pd/src/d_mayer_fft.c
@@ -64,29 +64,29 @@
#if defined(GOOD_TRIG)
#define FHT_SWAP(a,b,t) {(t)=(a);(a)=(b);(b)=(t);}
-#define TRIG_VARS \
+#define TRIG_VARS \
int t_lam=0;
-#define TRIG_INIT(k,c,s) \
- { \
- int i; \
- for (i=2 ; i<=k ; i++) \
- {coswrk[i]=costab[i];sinwrk[i]=sintab[i];} \
- t_lam = 0; \
- c = 1; \
- s = 0; \
+#define TRIG_INIT(k,c,s) \
+ { \
+ int i; \
+ for (i=2 ; i<=k ; i++) \
+ {coswrk[i]=costab[i];sinwrk[i]=sintab[i];} \
+ t_lam = 0; \
+ c = 1; \
+ s = 0; \
}
-#define TRIG_NEXT(k,c,s) \
- { \
- int i,j; \
- (t_lam)++; \
- for (i=0 ; !((1<<i)&t_lam) ; i++); \
- i = k-i; \
- s = sinwrk[i]; \
- c = coswrk[i]; \
- if (i>1) \
- { \
- for (j=k-i+2 ; (1<<j)&t_lam ; j++); \
- j = k - j; \
+#define TRIG_NEXT(k,c,s) \
+ { \
+ int i,j; \
+ (t_lam)++; \
+ for (i=0 ; !((1<<i)&t_lam) ; i++); \
+ i = k-i; \
+ s = sinwrk[i]; \
+ c = coswrk[i]; \
+ if (i>1) \
+ { \
+ for (j=k-i+2 ; (1<<j)&t_lam ; j++); \
+ j = k - j; \
sinwrk[i] = halsec[i] * (sinwrk[i-1] + sinwrk[j]); \
coswrk[i] = halsec[i] * (coswrk[i-1] + coswrk[j]); \
} \
@@ -95,20 +95,20 @@
#endif
#if defined(FAST_TRIG)
-#define TRIG_VARS \
+#define TRIG_VARS \
REAL t_c,t_s;
-#define TRIG_INIT(k,c,s) \
- { \
- t_c = costab[k]; \
- t_s = sintab[k]; \
- c = 1; \
- s = 0; \
+#define TRIG_INIT(k,c,s) \
+ { \
+ t_c = costab[k]; \
+ t_s = sintab[k]; \
+ c = 1; \
+ s = 0; \
}
-#define TRIG_NEXT(k,c,s) \
+#define TRIG_NEXT(k,c,s) \
{ \
REAL t = c; \
- c = t*t_c - s*t_s; \
- s = t*t_s + s*t_c; \
+ c = t*t_c - s*t_s; \
+ s = t*t_s + s*t_c; \
}
#define TRIG_RESET(k,c,s)
#endif
@@ -227,58 +227,58 @@ REAL c1,s1,s2,c2,s3,c3,s4,c4;
REAL aa;
for (k=n>>1; (!((k2^=k)&k)); k>>=1);
if (k1>k2)
- {
- aa=fz[k1];fz[k1]=fz[k2];fz[k2]=aa;
- }
+ {
+ aa=fz[k1];fz[k1]=fz[k2];fz[k2]=aa;
+ }
}
for ( k=0 ; (1<<k)<n ; k++ );
k &= 1;
if (k==0)
{
- for (fi=fz,fn=fz+n;fi<fn;fi+=4)
- {
- REAL f0,f1,f2,f3;
- f1 = fi[0 ]-fi[1 ];
- f0 = fi[0 ]+fi[1 ];
- f3 = fi[2 ]-fi[3 ];
- f2 = fi[2 ]+fi[3 ];
- fi[2 ] = (f0-f2);
- fi[0 ] = (f0+f2);
- fi[3 ] = (f1-f3);
- fi[1 ] = (f1+f3);
- }
+ for (fi=fz,fn=fz+n;fi<fn;fi+=4)
+ {
+ REAL f0,f1,f2,f3;
+ f1 = fi[0 ]-fi[1 ];
+ f0 = fi[0 ]+fi[1 ];
+ f3 = fi[2 ]-fi[3 ];
+ f2 = fi[2 ]+fi[3 ];
+ fi[2 ] = (f0-f2);
+ fi[0 ] = (f0+f2);
+ fi[3 ] = (f1-f3);
+ fi[1 ] = (f1+f3);
+ }
}
else
{
- for (fi=fz,fn=fz+n,gi=fi+1;fi<fn;fi+=8,gi+=8)
- {
- REAL bs1,bc1,bs2,bc2,bs3,bc3,bs4,bc4,
- bg0,bf0,bf1,bg1,bf2,bg2,bf3,bg3;
- bc1 = fi[0 ] - gi[0 ];
- bs1 = fi[0 ] + gi[0 ];
- bc2 = fi[2 ] - gi[2 ];
- bs2 = fi[2 ] + gi[2 ];
- bc3 = fi[4 ] - gi[4 ];
- bs3 = fi[4 ] + gi[4 ];
- bc4 = fi[6 ] - gi[6 ];
- bs4 = fi[6 ] + gi[6 ];
- bf1 = (bs1 - bs2);
- bf0 = (bs1 + bs2);
- bg1 = (bc1 - bc2);
- bg0 = (bc1 + bc2);
- bf3 = (bs3 - bs4);
- bf2 = (bs3 + bs4);
- bg3 = SQRT2*bc4;
- bg2 = SQRT2*bc3;
- fi[4 ] = bf0 - bf2;
- fi[0 ] = bf0 + bf2;
- fi[6 ] = bf1 - bf3;
- fi[2 ] = bf1 + bf3;
- gi[4 ] = bg0 - bg2;
- gi[0 ] = bg0 + bg2;
- gi[6 ] = bg1 - bg3;
- gi[2 ] = bg1 + bg3;
- }
+ for (fi=fz,fn=fz+n,gi=fi+1;fi<fn;fi+=8,gi+=8)
+ {
+ REAL bs1,bc1,bs2,bc2,bs3,bc3,bs4,bc4,
+ bg0,bf0,bf1,bg1,bf2,bg2,bf3,bg3;
+ bc1 = fi[0 ] - gi[0 ];
+ bs1 = fi[0 ] + gi[0 ];
+ bc2 = fi[2 ] - gi[2 ];
+ bs2 = fi[2 ] + gi[2 ];
+ bc3 = fi[4 ] - gi[4 ];
+ bs3 = fi[4 ] + gi[4 ];
+ bc4 = fi[6 ] - gi[6 ];
+ bs4 = fi[6 ] + gi[6 ];
+ bf1 = (bs1 - bs2);
+ bf0 = (bs1 + bs2);
+ bg1 = (bc1 - bc2);
+ bg0 = (bc1 + bc2);
+ bf3 = (bs3 - bs4);
+ bf2 = (bs3 + bs4);
+ bg3 = SQRT2*bc4;
+ bg2 = SQRT2*bc3;
+ fi[4 ] = bf0 - bf2;
+ fi[0 ] = bf0 + bf2;
+ fi[6 ] = bf1 - bf3;
+ fi[2 ] = bf1 + bf3;
+ gi[4 ] = bg0 - bg2;
+ gi[0 ] = bg0 + bg2;
+ gi[6 ] = bg1 - bg3;
+ gi[2 ] = bg1 + bg3;
+ }
}
if (n<16) return;
@@ -292,71 +292,71 @@ REAL c1,s1,s2,c2,s3,c3,s4,c4;
k4 = k2 << 1;
k3 = k2 + k1;
kx = k1 >> 1;
- fi = fz;
- gi = fi + kx;
- fn = fz + n;
- do
- {
- REAL g0,f0,f1,g1,f2,g2,f3,g3;
- f1 = fi[0 ] - fi[k1];
- f0 = fi[0 ] + fi[k1];
- f3 = fi[k2] - fi[k3];
- f2 = fi[k2] + fi[k3];
- fi[k2] = f0 - f2;
- fi[0 ] = f0 + f2;
- fi[k3] = f1 - f3;
- fi[k1] = f1 + f3;
- g1 = gi[0 ] - gi[k1];
- g0 = gi[0 ] + gi[k1];
- g3 = SQRT2 * gi[k3];
- g2 = SQRT2 * gi[k2];
- gi[k2] = g0 - g2;
- gi[0 ] = g0 + g2;
- gi[k3] = g1 - g3;
- gi[k1] = g1 + g3;
- gi += k4;
- fi += k4;
- } while (fi<fn);
+ fi = fz;
+ gi = fi + kx;
+ fn = fz + n;
+ do
+ {
+ REAL g0,f0,f1,g1,f2,g2,f3,g3;
+ f1 = fi[0 ] - fi[k1];
+ f0 = fi[0 ] + fi[k1];
+ f3 = fi[k2] - fi[k3];
+ f2 = fi[k2] + fi[k3];
+ fi[k2] = f0 - f2;
+ fi[0 ] = f0 + f2;
+ fi[k3] = f1 - f3;
+ fi[k1] = f1 + f3;
+ g1 = gi[0 ] - gi[k1];
+ g0 = gi[0 ] + gi[k1];
+ g3 = SQRT2 * gi[k3];
+ g2 = SQRT2 * gi[k2];
+ gi[k2] = g0 - g2;
+ gi[0 ] = g0 + g2;
+ gi[k3] = g1 - g3;
+ gi[k1] = g1 + g3;
+ gi += k4;
+ fi += k4;
+ } while (fi<fn);
TRIG_INIT(k,c1,s1);
for (ii=1;ii<kx;ii++)
{
- REAL c2,s2;
+ REAL c2,s2;
TRIG_NEXT(k,c1,s1);
c2 = c1*c1 - s1*s1;
s2 = 2*(c1*s1);
- fn = fz + n;
- fi = fz +ii;
- gi = fz +k1-ii;
- do
- {
- REAL a,b,g0,f0,f1,g1,f2,g2,f3,g3;
- b = s2*fi[k1] - c2*gi[k1];
- a = c2*fi[k1] + s2*gi[k1];
- f1 = fi[0 ] - a;
- f0 = fi[0 ] + a;
- g1 = gi[0 ] - b;
- g0 = gi[0 ] + b;
- b = s2*fi[k3] - c2*gi[k3];
- a = c2*fi[k3] + s2*gi[k3];
- f3 = fi[k2] - a;
- f2 = fi[k2] + a;
- g3 = gi[k2] - b;
- g2 = gi[k2] + b;
- b = s1*f2 - c1*g3;
- a = c1*f2 + s1*g3;
- fi[k2] = f0 - a;
- fi[0 ] = f0 + a;
- gi[k3] = g1 - b;
- gi[k1] = g1 + b;
- b = c1*g2 - s1*f3;
- a = s1*g2 + c1*f3;
- gi[k2] = g0 - a;
- gi[0 ] = g0 + a;
- fi[k3] = f1 - b;
- fi[k1] = f1 + b;
- gi += k4;
- fi += k4;
- } while (fi<fn);
+ fn = fz + n;
+ fi = fz +ii;
+ gi = fz +k1-ii;
+ do
+ {
+ REAL a,b,g0,f0,f1,g1,f2,g2,f3,g3;
+ b = s2*fi[k1] - c2*gi[k1];
+ a = c2*fi[k1] + s2*gi[k1];
+ f1 = fi[0 ] - a;
+ f0 = fi[0 ] + a;
+ g1 = gi[0 ] - b;
+ g0 = gi[0 ] + b;
+ b = s2*fi[k3] - c2*gi[k3];
+ a = c2*fi[k3] + s2*gi[k3];
+ f3 = fi[k2] - a;
+ f2 = fi[k2] + a;
+ g3 = gi[k2] - b;
+ g2 = gi[k2] + b;
+ b = s1*f2 - c1*g3;
+ a = c1*f2 + s1*g3;
+ fi[k2] = f0 - a;
+ fi[0 ] = f0 + a;
+ gi[k3] = g1 - b;
+ gi[k1] = g1 + b;
+ b = c1*g2 - s1*f3;
+ a = s1*g2 + c1*f3;
+ gi[k2] = g0 - a;
+ gi[0 ] = g0 + a;
+ fi[k3] = f1 - b;
+ fi[k1] = f1 + b;
+ gi += k4;
+ fi += k4;
+ } while (fi<fn);
}
TRIG_RESET(k,c1,s1);
} while (k4<n);