递归式FFT

递归式FFT

N点有限长序列x(n)

 

递归式FFT变换方法：

 

FFT函数：

快速傅立叶变换 

//数组a为输入，数组y为输出，2的power次方为数组的长度 

void ImageProcess::fft(const complex<double> a[], complex<double> y[], int power) 

{ 

if(0 == power) 

{ 

y[0] = a[0]; 

return; 

} 

int n = 1 << power; 

double angle = 2 * PI / n; 

complex<double> wn(cos(angle), sin(angle)); 

complex<double> w(1, 0); 

complex<double> *a0 = new complex<double>[n / 2]; 

complex<double> *a1 = new complex<double>[n / 2]; 

complex<double> *y0 = new complex<double>[n / 2]; 

complex<double> *y1 = new complex<double>[n / 2]; 

for(int i = 0; i < n / 2; i ++) 

{ 

a0[i] = a[2 * i]; 

a1[i] = a[2 * i + 1]; 

} 

//分开成两个子fft过程 

fft(a0, y0, power - 1); 

fft(a1, y1, power - 1); 

complex<double> u; 

for(int k = 0; k < n / 2; k++) //蝶形算法 

{ 

u = w * y1[k]; 

y[k] = y0[k] + u; 

y[k + n / 2] = y0[k] - u; 

w = w * wn; 

} 

delete[] a0; 

delete[] a1; 

delete[] y0; 

delete[] y1; 

} 

快速傅立叶反变换 

//y为输入，a为输出，2的power次方为数组的长度 

void ImageProcess::ifft(const complex<double> y[], complex<double> a[], int power) 

{ 

int count = 1 << power; 

complex<double> *x = new complex<double>[count]; 

memcpy(x, y, sizeof(complex<double>) * count); 

int i; 

for(i = 0; i < count; i++) 

{ 

x[i] = complex<double>(x[i].real(), -x[i].imag()); //共轭复数 

} 

fft(x, a, power); //调用快速傅立叶变换算法 

for(i = 0; i < count; i++) 

{ 

a[i] = complex<double>(a[i].real() / count, -a[i].imag() / count); //共轭复数 

} 

delete[] x; 

}