一维快速傅里叶变换(FFT)
在一维DFT的基础上,推导一维快速傅里叶变换的方法。要利用$e^{-jk\frac{2\pi}{N}n}$的性质,将其记为$W_N^{kn}$主要有以下性质:
周期性:$W_N^{a+N}=W_N^{a}$
对称性:$W_N^{a+\frac{N}{2}}=W_N^{a}$
缩放性:$W_N^{2}=W_{\frac{N}{2}}^{1}$
证明方法就是按旋转因子定义,直接拆开就行,就是代数变换。以上性质意味着值相同的就不用了再多计算一遍了,这就能简化DFT的计算过程。
FFT的递推公式就是基于上述性质对DFT进行变换得到的。
首先,把输入的时域信号x[n] , n = 0 , 1 , . . . ,N−1根据索引分为奇偶两部分
对DFT公式进行化简:
根据旋转因子的缩放性,可以进一步换算:
进一步得到递推公式:
其中,$F_{even}[k]$为偶数项的DFT结果,$F_{odd}[k]$为奇数项DFT的结果。$W_N^{k}$是$e^{-j\frac{2\pi}{N}k}$,是一个复数。有了递推公式就可以递归地实现该算法。然而递归是非常慢的操作,观察公式可以发现。要计算$X[k]$,只需要将输入序列分为偶数项和奇数项,计算偶数项的DFT,可以继续划分,直到最后只剩两项时,计算两个数的DFT。然后向上层传递。
假设有输入序列有8个数。x[0],x[1],x[2]….x[7],划分层次如下:
最底下的数就是输入序列,按这个方式组合计算即可获得更高层次的结果,所有只需要自底向上层层按公式计算即可。在这里,不要被上面的(0,4),(2,6),(1,5),(3,7)的组合放置的顺序误导了,这里的遍历顺序仍然是0,1,2,3。具体的可以参考后面给出的代码实现。这里想要重点提一点的是,上面的$x[n]$虽然输入是一个实数,但其实在代码实现中,可以转成复数计算,这样一来,整个傅里叶变换过程和逆变换过程,都是在复数空间中计算的,计算结果也是复数,这样去看公式,就更统一和好理解。
更进一步地,分析快速傅里叶逆变换(IFFT),这里直接上公式,就不进行推导了,经过了前面的洗礼,推导的过程应该是很简单的。
由此可见,IFFT的公式除了左边的$\frac{1}{N}$之外,还是一个复数乘积的求和公式。但是$W_N^{-kn}$与FFT的相比,反过来了,$X[k]$就是FFT计算得到的序列。而快速计算的推导方法和FFT相同,就不再赘述。下面实现的代码仓库。
FFT和IFFT代码实现的仓库地址
二维快速傅里叶变换(FFT)
二维FFT是在一维FFT基础上实现,实现过程为:
对二维输入数据的每一行进行FFT,变换结果仍然按行存入二维数组中。结果是一行复数序列。
在1的结果基础上,对每一列进行FFT,得到二维FFT结果。
这里仍然拿之前二维DFT的示例程序进行修改。得到的结果与DFT的结果一致。
对于二维快速傅里叶逆变换,上面两点结论改改依然成立:
- 对二维FFT的每一行进行IFFT,结果是一行复数序列,变换结果仍然按行存入二维数组中。
- 在1的结果基础上,对每一列进行IFFT,得到二维IFFT结果。