快速傅利葉變換

計算復序列的快速傅利葉變換。

序列\(x(n)(n=0,1,...,n-1)\)的離散傅利葉變換定義為

\[x(k)=\sum_^x(n)w_^, \qquad k=0,1,...,n-1

\]其中\(w_^=e^}\)，將序列\(x(n)\)按序號\(n\)的奇偶分成兩組，即

\[\left.\begin\beginx_(n)=&x(2n)\\ x_(n)=&x(2n+1)\end\end\right\} \qquad n=0,1,...,\frac-1

\]因此，\(x(n)\)的傅利葉變換可寫成

\[\beginx(k) &= \sum_^x(2n)w^_ + \sum_^x(2n+1)w^_\\&= \sum_^x_(n)w^_ + w_^\sum_^x_(n)w^_\end

\]由此可得\(x(k)=x_(k)+w_^x_(k), \qquad k = 0,1,...,\frac\)，式中

\[\beginx_(k)&=\sum_^x(2n)w^_\\x_(k)&=\sum_^x(2n+1)w^_\end

\]他們分別是\(x_1(n)\)和\(x_2(n)\)的\(n/2\)點dft。上面的推導表明：乙個\(n\)點dft被分解為兩個\(n/2\)點dft，這兩個\(n/2\)點dft又可合成乙個\(n\)點dft。但上面給出的公式僅能得到\(x(k)\)的前\(n/2\)點的值，要用\(x_(k)\)和\(x_(k)\)來表達\(x(k)\)的後半部分的值，還必須運用權係數\(w_n\)的週期性與對稱性，即

\[w_^=w_^, \quad w_^=-w_^

\]因此，\(x(k)\)的後\(n/2\)點的值可表示為

\[\beginx(k+\frac)&=x_(k+\frac)+w_^x_(k+\frac)\\&=x_(k)-w_^x_(k), \ k=0,1,...,\frac-1\end

\]通過上面的推導可以看出，\(n\)點的dft可以分解為兩個\(n/2\)點dft，每個\(n/2\)點dft又可以分解為兩個\(n/4\)點dft。依此類推，當\(n\)為2的整數次冪時（\(n=2^m\)），由於每分解一次降低一階冪次，所以通過\(m\)次分解，最後全部成為一系列2點dft運算。以上就是按時間抽取的快速傅利葉變換（fft）演算法。

序列\(x(k)\)的離散傅利葉反變換定義為

\[x(n)=\frac\sum_^x(k)w_^, \qquad n=0,1,...,n-1

\]它與離散傅利葉正變換的區別在於將\(w_n\)改變為\(w_n^\)，並多了乙個除以\(n\)的運算。因為\(w_n\)和\(w_n^\)對於推導按時間抽取的快速傅利葉變換演算法並無實質性區別，因此可將fft和快速傅利葉反變換（ifft）演算法合併在同一程式中。

是用c語言實現快速傅利葉變換（fft）的方法如下：

/************************************
x ---一維陣列，長度為n，開始時存放要變換資料的實部，最後存放變換結果的實部。
y ---一維陣列，長度為n，開始時存放要變換資料的虛部，最後存放變換結果的虛部。
n ---資料長度，必須是2的整數次冪。
sign ---當sign=1時，子函式計算離散傅利葉正變換；當sign=-1時，子函式計算離散傅利葉反變換
************************************/
#include "math.h"
void fft(double *x, double *y, int n, int sign)
n1 = n - 1;
for(j = 0, i = 0; i < n1; i++) 
k = n / 2;
while(k < (j + 1)) 
j = j + k;
} n1 = 1;
for(l = 1; l <= m; l++) 
t = c;
c = c * c1 - s * s1;
s = t * s1 + s * c1;
} }if(sign == -1) 
}}

快速傅利葉變換

傅利葉變換與快速傅利葉變換

快速傅利葉變換

快速傅利葉變換

相關推薦