给您一个我在
《观音》 中用的FFT。
void fft(int n, double theta, double ar[], double ai[])
{
int m, mh, i, j, k;
double wr, wi, xr, xi;
for (m = n; (mh = m >> 1) >= 1; m = mh)
{
for (i = 0; i < mh; i++)
{
wr = cos(theta * i);
wi = sin(theta * i);
for (j = i; j < n; j += m)
{
k = j + mh;
xr = ar[j] - ar[k];
xi = ai[j] - ai[k];
ar[j] += ar[k];
ai[j] += ai[k];
ar[k] = wr * xr - wi * xi;
ai[k] = wr * xi + wi * xr;
}
}
theta *= 2;
}
i = 0;
for (j = 1; j < n - 1; j++)
{
for (k = n >> 1; k > (i ^= k); k >>= 1);
if (j < i)
{
xr = ar[j];
xi = ai[j];
ar[j] = ar
; ai[j] = ai; ar = xr; ai = xi; } }}使用例:低域通過濾波器 (Low-pass filter, LPF)#define FFT_SIZE 1024void LPF(ar, ai) { double theta = -8 * atan(1.0) / FFT_SIZE; fft(FFT_SIZE, theta, ar, ai); //高周波数部分にゼロを代入する for (i = iFreqHightIndex; i < FFT_SIZE / 2; i ++) { ar = ai = ar[FFT_SIZE - i - 1] = ai[FFT_SIZE - i - 1] = 0.0; } theta = 8 * atan(1.0) / FFT_SIZE; fft(FFT_SIZE, theta, ar, ai); for (i = 0; i < FFT_SIZE; i ++) { ar /= FFT_SIZE; }