一、FFT是什么，能干啥？
FFT(fast Fourier transform)快速傅里叶变换，把时域信号变换到频域，至于换到频域后能做啥，那就多了,你能搜到这篇文章，说明你已经知道他是干啥的了，，，，


二、先用matlab来认识一下FFT
1.上代码
代码如下：

1. fs = 1024;           %采样频率
   
2. NPT = 1024;          %样点个数
   
3. t=0:1/fs:1;          %采样步长,模拟单片机采样频率
   
4. y= 2+200*cos(2*pi*15*t)+40*cos(2*pi*100*t+pi*90/180)+15*cos(2*pi*250*t+pi*90/180);
   
5. plot(t,y);                     %画出未变换前的时域信号波形
   
6. fbl=fs/NPT;          %分辨率
   
7. f=(0:NPT)*fbl;       %每个点所代表的频率
   
8. Y=fft(y);                     %进行FFT变换，matlab中fft函数，输入多少采样数据，就输出多少幅频值(是个复数)
   
9.                                      %所以length(y)=length(Y);
   
10. figure(2)
    
11. plot(f(1:NPT/2),abs(Y(1:NPT/2)));
    
12.                                      %以每个点所代表的频率画出其对应的赋值
    
13. 
    

复制代码

2.代码关键字讲解
思想都是一样的，matab知道了，在32里面也差不多

采样频率： 用单片机ADC采样，结果是离散的数字信号，参考采样定理，采样频率至少大于信号频率的两倍，在这里我们采样频率1024,远大于基波频率2，选取1024也是为了方便自己，可以更直观的观察输出（2的整数倍都可）。
分辨率： 输出时，最小频率间隔(Fs/N)，第n个点所代表的频率即为Fn=(n-1)*Fs/N,在这个程序中，采样频率是1024，采样点数是1024，所以，分辨率为1Hz
输出赋值 fft输出就是复数，取模值即可


3.程序结果
未变换前的时域信号波形图：



转换后输出的幅频函数波形图：



由图可看出输入信号大概有三个频率，15，100，250.与输入信号相符
二、淦! 上32
STM32官方给了三个函数库文件，分别是：
cr4_fft_64_stm32.s
cr4_fft_256_stm32.s
cr4_fft_1024_stm32.s
（这些都是汇编语言写的，好在我们不需要看懂，只需要会调用就好。）
分别对应这采样点数为64，256，1024。库函数的获取，可以在官网上获取我们的目标是计算周期为1ms波形的THD值！，使用采样点数为64位的函数

1、代码移植
新建文件夹，把
cr4_fft_64_stm32.s
cr4_fft_256_stm32.s
cr4_fft_1024_stm32.s
和stm32_dsp.h
文件添加到工程里（就像添加.C文件一样），注：stm32_dsp.h文件里第27行需要根据自己单片机型号来修改




2，ADC采集
ADC采集没有什么好说的，我想让大家注意一点就是采样周期，在这里，为了方便使用，我们设ADC采样频率为64KHz，注意，这里直接调用正点原子的Get_Adc();函数的话，他的执行周期很长，需要改短（自己跳进Get_Adc(）函数里，看函数定义，再改），要不然，在定时器中，ADC采集时间会超过中断事件，导致采样周期不是64KHz,最好在中断中设置一个LED灯反转的程序，检测ADC采样周期，确保ADC采集频率稳定


3、上函数
在上函数之前，我们把函数的每一个输入输出变量的含义说一下
在我的程序里是这样调用的：cr4_fft_64_stm32(output,input,NPT);
她的输入变量是input：是一个int32_t类型的数组，高十六位是实部，第十六位是虚部，显然，我们ADC采集只能采集实数。所以我对采集到的数据进行这样的处理：

1.   ADC_Value[i]=Get_Adc(ADC_Channel_5);
   
2.          input[i] = ((s16)ADC_Value[i])<<16;

复制代码

至于虚部，就不用管啦，移位后默认是0
他的另一个输入是NPT ，采样点数。最好是定义成宏定义，因为要用到他定义的数组太多，这里是64
他的输出是output,也是个复数，高十六位代表实部，低十六位代表虚部，他们的模值就是我们要的幅值
呐，官方给了专门求赋值的函数

1. void dsp_asm_powerMag(int32_t *IBUFOUT)
   
2. {
   
3.         s16 lx,ly;
   
4.         u32 i;
   
5.         for(i=0;i<NPT/2;i++)//由于FFT的频域结果是关于奈奎斯特频率对称的，所以只计算一半
   
6.         {
   
7.                 lx = (IBUFOUT[i] << 16)>>16;      //取低十六位，虚部
   
8.                 ly = (IBUFOUT[i] >> 16);          //取高十六位，实部
   
9.                 float X = NPT*((float)lx)/32768;
   
10.                 float Y = NPT*((float)ly)/32768;
    
11.                 float Mag = sqrt(X*X+Y*Y)/NPT;
    
12.                 OUTPUT_MAG[i] = (u32)(Mag*65536);
    
13.         }//这些就是计算振幅IBUFOUT[i] = sqrt(lx*lx+ly*ly)*2/NPT
    
14.          //之所以先除以32768后乘以65536是为了符合浮点数的计算规律
    
15.         OUTPUT_MAG[0] = OUTPUT_MAG[0]/2;//这个是直流分量，不需要乘以2
    
16. }
    

复制代码

直接用就很OK了

1. dsp_asm_powerMag(output);               
   

复制代码

那么我们调用了函数，也转换成幅值了，所以，我该怎么用这些输出的数呢。他输出的是什么***，接下来就需要对应起来，
我们需要知道的是 每个频率点对应的幅值，现在有幅值，怎么把他对应的点找到，这样我们就大功告成了
还记得前面讲的分辨率嘛，对，我们采样点是64个，采样频率也是64KHz,那么分辨率就是fs/NPT即为1KHz，也就是说输出的output每隔一个对应着一个相差位1Hz的频率，那么output[1]就对应着1KHz，output[2]就对应着2KHz，output[3]就对应着3KHz，，，，我们只取前NPT/2个，不要问为什么，问奈奎斯特吧。
还有就是output[0]对应着直流分量

我只用了前几个意思意思。

1. cr4_fft_64_stm32(output,input,NPT);
   
2.         dsp_asm_powerMag(output);               
   
3. 
   
4.         LCD_ShowxNum(100,0,    OUTPUT_MAG[0]  ,7,16,0);        
   
5.         LCD_ShowxNum(100,16,   OUTPUT_MAG[1]  ,7,16,0);               
   
6.         LCD_ShowxNum(100,32,   OUTPUT_MAG[2]  ,7,16,0);               
   
7.         LCD_ShowxNum(100,48,   OUTPUT_MAG[3]  ,7,16,0);               
   
8.         LCD_ShowxNum(100,64,   OUTPUT_MAG[4]  ,7,16,0);               
   
9.         LCD_ShowxNum(100,80,   OUTPUT_MAG[5]  ,7,16,0);        
   
10.     THD=100*(float)sqrt((Voltage[4]*Voltage[4]+Voltage[2]*Voltage[2]+Voltage[3]*Voltage[3]+Voltage[5]*Voltage[5]))/(float)(Voltage[1]);
    

复制代码

————————————————
版权声明：ddydy111