我相信有很多关于FFT的在线教程，我也看过很多在线教程，都不是很好理解。在对FFT及其编程实现进行了基本研究后。我决定写一个关于FFT的通俗易懂的解释。所以我会在接下来的叙述中尽量用一种很通俗很详细的方式来解释。
我第一次知道傅立叶变换的时候，我沉迷于音乐的频谱，我想做一个音乐频谱显示器。看了很多资料，了解了傅立叶变换和FFT。当然，这一切都发生在之前。经过两周的系统学习和研究，我自制了一个FFT算法。我不敢说它在速度上高人一等，但我个人认为它是一种简单易懂的实现方法。经过实际的VC仿真实验，用STM32运行也很成功。
首先，要了解FFT，必须对DFT有所了解，因为两者本质上是一样的。在此之前，列出离散傅里叶变换对(DFT):
，k=0，1，…N-1
n=0，1…N-1
其中包括：
但FFT之所以叫快速傅立叶变换，是利用了以下性质(最重要的！)
周期性：
对称性：
还原性：
先把这三个公式放在这里，然后你会用到它们。
根据这些特性，一个长的DFT运算可以分解成几个短的DFT运算序列的组合，从而减少计算量。这里，为了便于理解，我使用了时间抽取快速傅立叶变换(DIT-FFT)方法。
首先，把一个序列x(n)一分为二。
对于，k=0，1，n-1，设n是2的整数次幂，即n=2 m。
按奇偶性分组x(n)。
x(2r)=x1(r)
x(2r 1)=x2(r) r=0，1，….(N/2)-1
然后DFT也分两组做预算。
(第一项为偶数，第二项为奇数)
这时，上面提到的三个属性的可还原性将发挥作用：
回顾：
那么这个公式可以简化为：
(公式1-1)
那么x1(k)和X2(k)是长度为N/2的序列X1(k)和x2(k)的N/2个点的离散傅立叶变换。
其中包括：
K=0，1，2…N/2-1
此时，一个N点的DFT被分解成两个N/2的DFT。但是X1(k)，X2(k)只有N/2个点，也就是说，公式(1-1)只是DFT的前半部分。要求DFT的后半部分可以利用其对称性得到，如下：
(公式1-2)
那么表达式(1-1)和(1-2)可以由蝶形信号流图符号来表示。如图所示：
以N=8为例，可以用下图表示：
通过这种分解，每N/2点DFT只需要(n 2)/4次复数乘法，明显节省了计算量。
依此类推，继续把得到的X1(k)和X2(k)按奇偶性分解，方法同上。然后，我们得到了这张可以在百度上找到的图片的大推送。(笑)
对于这张图，我想强调的是，在进行二阶蝶形运算的时候，不应该WN1后面是WN0吗？为什么是W2N？这个问题已经困扰我一段时间了，但是唐别忘了前者的W0N展开式是W0N/2因为此时N已经按宇称分开了，所以是N/2而W2N/2也就是W2N是按可约性得到的。它可以不要在这里混淆。
运营效率就不用提了。
那 FFT算法的理论内容，然后它用C语言实现了这个算法。对于FFT算法的C语言实现，在线方法层出不穷。因为我我很懒(我我懒得看别人程序)，再加上自给自足的原则，我我写了一个我个人认为很容易理解的程序，为了帮助读者理解，我我有意尽可能减少库函数的使用。一些基本函数是自己写的(难免有很多bug)，但是作为FFT算法已经足够了。目前这个程序只能处理2 N个数据。理论上，如果不够，可以按以下顺序加0。我没有为了简单起见，我们不实现这个操作，但是主要的函数都是可用的。代码如下：
#包括
#包括
#定义PI 3.1415926
#定义FFT_LENGTH 8
double input[FFT_LENGTH]={1，1，1，1，1，1，1 }；
Struct complex1{ //定义复杂结构
双实；//实部
双重形象；//虚部
};
//将输入的实数结果存储为复数
struct complex1结果_ dat[8]；
结构复杂1 con_complex(结构复杂1 a，结构复杂1 b){
结构复杂1温度；
temp . real=(a . real * b . real)-(a . image * b . image)；
temp . image=(a . image * b . real)(a . real * b . image)；
返回温度；
}
int mypow(int a，int b){
int i，sum=a；
if(b==0)返回1；
for(I=1；我
sum *=a；
}
返回总和；
}
int log2(int n){
无符号I=1；
int sum=1；
对于(我；i ){
sum *=2；
if(sum=n)break；
}
返回I；
}
void交换(结构复杂1 *a，结构复杂1 *b){
结构复杂1温度；
temp=* a；
* a=* b；
* b=温度；
}
void fft(struct complex1 dat[]，无符号字符N){
struct complex1 dat_buf，dat _ org
/* L是几个蝶形运算，也代表二进制位数。
n是当前蝶泳所需的次数。N最初是N/2，每次蝶形运算后需要/2。
I是反转要用的自增符号，我也用L盘系列。j是计算当前光盘类别的计算次数。
Re_i i_copy是反演中使用的变量。
k是当前圆盘cos(2*pi/N*k)，k也是E的k(-J2 * pi/N)* k
*/
无符号char L，I，j，re_i=0，i_copy=0，k=0，FFT _ flag=1；
//经过观察发现，每个蝶形运算都需要N/2次运算，总共N/2*log2N次运算。
无符号字符FFT _ counter=0；
//这里要补充的2 N，这里必须省略2 n。
//蝴蝶系列(L级)
l=log2(N)；
//计算每一关蝴蝶计算的次数(这里只是一个初始值)然后每次求/2。
//N=N/2；
//反转dat的顺序
for(I=1；ii _ copy=I；
re _ I=0；
for(j=L-1；j0；j - ){
//判断I的副本的最低位数，并将其移动到最高位。
//re_i是交换的号码。每次它的数字不能移动，循环后都要清零。
re _ I |=((I _ copy0x 01)I _ copy=1；
}
swap(dat[i]，dat[re _ I])；
}
//执行fft计算
for(I=0；IFFT _ flag=1；
FFT _ counter=0；
for(j=0；Jif (FFT _ counter==mypower (2，I)){//每隔几次控制，每隔几次计算，
FFT _ flag=0；
} else if(FFT _ counter==0){//rest结束，继续运算。
FFT _ flag=1；
}
//当这个语句没有被判断时，fft_flag保持，这样操作就可以继续。
if(fft_flag){
dat _ buf . real=cos((2 * PI * k)/(N/mypow(2，L-I-1)))；
dat _ buf . image=-sin((2 * PI * k)/(N/mypow(2，L-I-1)))；
dat _ buf=con _ complex(dat[j mypow(2，i)]，dat _ buf)；
//计算当前蝶形运算的奇数项与w的乘积。
dat_org.real=dat[j]。真实；
dat_org.image=dat[j]。形象；//临时存储
dat[j]。real=dat _ org . real dat _ buf . real；
dat[j]。image=dat _ org . image dat _ buf . image；
//实部加实部加虚部加虚部
dat[j mypow(2，i)]real=dat _ org . real-dat _ buf . real；
dat[j mypow(2，i)]image=dat _ org . image-dat _ buf . image；
//实部减实部减虚部
k；
fft _ counter
}否则{
FFT _ counter-；//算几次，就休息几次。
k=0；
}
}
}
}
void main(){
int I；
//首先将输入信号转换成复数。
for(I=0；iresult_dat[i]。image=0；
//输入信号是二维的，暂时没有复数。
result_dat[i]。real=input[I]；
//result_dat[i]。实数=10；
//输入信号都是实数
}
fft(result_dat，FFT _ LENGTH)；
for(I=0；iinput[i]=sqrt(result_dat[i])。real*result_dat[i]。真实结果_dat[i]。image*result_dat[i]。图片)；
//拿模具
printf(% lf \ n ，输入[I])；
}
while(1)；
}
在这个程序中，输入的数组是输入信号，这里只取8个模拟样本，输出数据也是输入的。如果你想看到其他的序列，你可以改变FFT_LENGTH的值和输入的内容。程序输出实部和虚部的模块。如果你只是单纯的想看实部或虚部的幅度，请自行修改程序~这就是这篇关于FFT的简介和FFT算法的基本实现的全部内容。参考书：数字信号处理
标签：FFTIn