最近,因为项目需要在STM32F103系列处理器上,对采集的音频信号进行FFT运算,然而STM32F103毕竟不是STM32F4系列的处理器,对于一般的FFT运算程序还是比较缓慢的。
幸亏官方提供了针对FFT的官方库,但是去官网找了半天居然找不到那个库的下载,好像官方早就把那个库下架了,估计是为了给带DSP指令集的F4系列让路。然后就只好从别人的项目中把这个官方库给扒出来了……
下载地址:https://pan.baidu.com/s/1GiUJgEkQxDAk79iddXKsaA 提取码: dn2d
FFT的意义
对于很多人来说,采样频率和FFT点数之间的关系可能还是不太清楚。下面就来简单分析一下:
根据采样定理,采样频率必须是被采样信号最高频率的2倍。比如,需要采集音频信号,并且需要被观察到的音频频率的频率范围是20Hz到20KHz,那么使用的采样频率就必须大于40kHz。如果需要观察到的音频频率范围为0Hz到600Hz,那么使用的采样频率只需要大于1200Hz即可。
而FFT点数与采样频率之间有什么关系呢?本质上并没有什么关系,但是FFT点数的大小直接关系到频率分辨率。怎么来说呢?
假设采样频率为Fs,信号频率F,采样点数为N。那么FFT之后结果就是一个为N点的复数。每一个点就对应着一个频率点。这个点的模值,就是该频率值下的幅度特性。同时,FFT后的N个点,开始的那个点表示直流分量(即0Hz),而最后的那个点的再下一个点表示采样频率Fs,这中间被N-1个点平均分成N等份,每个点的频率依次增加。即,某点n所表示的频率为:Fn=(n-1)*Fs/N。这就表示,Fs/N就是频率分辨率。
不太理解的,可以查看博客:FFT后的物理意义:https://blog.csdn.net/iloveyoumj/article/details/53308142。
FFT官方库的使用条件
FFT官方库在使用上并不灵活:属于基4的FFT,即FFT点数必须是4^n。也就是说,如果要做512点或2048点的FFT,那么对不起,没法使用官方库了;
FFT官方库的输入输出是等长的,即256点的FFT输入也必须是256点,如果你的输入小于这个长度,是没有任何性能提升的。
FFT官方库的使用
准备工作
下载得到STM32的DSP库之后,就可以将其添加到自己的工程项目中了。
其中,stm32_dsp.h和table_fft.h两个文件是必须添加的。stm32_dsp.h是STM32的DSP库的头文件。
另外,对于**.s文件可以有选择的添加**(用到那个添加那个即可)。由于本文只用到了256点的FFT,所以这里只添加了cr4_fft_256_stm32.s文件。
FFT函数说明
进行256点的FFT,只需要调用STM32 DSP库函数中的cr4_fft_256_stm32()函数即可。该函数的原型为:
void cr4_fft_256_stm32(void *pssOUT, void *pssIN, uint16_t Nbin);
其中,参数pssOUT表示FFT输出数组指针,参数pssIN表示要进行FFT运算的输入数组指针,参数Nbin表示了点数。至于该函数的具体实现,因为是用汇编语言编写的,我也不懂,这里就不妄谈了。
需要说明的是:按照FFT官方库的说明,pssOUT和pssIN都必须是32位的数据类型,其中高16位存储实部,低16位存储虚部。对于pssIN来说,低16位存储的虚部总是为0。
代码实例
假设ADC采样的声音数据为adc_buf[NPT],FFT运算的输入数组为lBufInArray[NPT]。由于FFT计算出来的数据是对称的,因此通常而言输出数组取一半的数据,为lBufOutArray[NPT/2]。除此之外,还定义各次谐波幅值lBufMagArray[NPT/2]。即:
#define NPT 256 uint32_t adc_buf[NPT]={0}; long lBufInArray[NPT]; long lBufOutArray[NPT/2]; long lBufMagArray[NPT/2];
调用FFT官方库的代码为:
//填充数组 for(i=0;i<NPT;i++) //这里因为单片机的ADC只能测正的电压 所以需要前级加直流偏执 //加入直流偏执后,需要在软件上减去2048即一半,达到负半周期测量的目的(需要根据具体情况来进行配置) lBufInArray[i] = ((signed short)(adc_buf[i]-2048)) << 16; cr4_fft_256_stm32(lBufOutArray, lBufInArray, NPT);
同时,计算各次谐波幅值的函数为:
void GetPowerMag() { signed short lX,lY; float X,Y,Mag; unsigned short i; for(i=0; i < NPT/2; i++) { lX = (lBufOutArray[i] << 16) >> 16; lY = (lBufOutArray[i] >> 16); //除以32768再乘65536是为了符合浮点数计算规律 X = NPT * ((float)lX) / 32768; Y = NPT * ((float)lY) / 32768; Mag = sqrt(X * X + Y * Y) / NPT; if(i == 0) lBufMagArray[i] = (unsigned long)(Mag * 32768); else lBufMagArray[i] = (unsigned long)(Mag * 65536); } }
---------------------
版权声明:本文为CSDN博主「Yngz_Miao」的原创文章,遵循CC 4.0 by-sa版权协议,转载请附上原文出处链接及本声明。
原文链接:https://blog.csdn.net/qq_38410730/article/details/90116695