【随手笔记】FFT资料整理（一）

LongRunning

已于 2025-03-02 17:14:40 修改

阅读量1.2k

点赞数 30

CC 4.0 BY-SA版权

分类专栏：软件文章标签：笔记单片机

于 2025-01-23 15:48:09 首次发布

本文链接：https://blog.youkuaiyun.com/LongRunning/article/details/145200522

软件专栏收录该内容

9 篇文章

订阅专栏

（一）结果验证

函数波形示例1

#define Fs  44800
#define NPT 256
void InitBufInArray()
{
    int i = 0;
    float fx = 0;
    for(i=0; i<NPT; i++)
    {
        // fx = 1500 * sin(2*PI * i * 350.0 / Fs) +
        //      2700 * sin(2*PI * i * 8400.0 / Fs) +
        //      4000 * sin(2*PI * i * 18725.0 / Fs);

        fx = 1500 * sin(2*PI * i * 350.0 / Fs);//频率350hz幅值峰值为1500
        adc_buff[i] = ((uint16_t)fx) ;
			
		uart_output[3] = adc_buff[i] >> 8;
        uart_output[4] = adc_buff[i] & 0xff;
      	HAL_UART_Transmit(&huart1, (uint8_t *)uart_output, 5, 0xFFFF);
      	HAL_Delay(50);
		//printf("[%d]\r\n", adc_buff[i]);
    }
}

正弦-波峰1500-频率350
在这里插入图片描述

函数波形示例2

void InitBufInArray()
{
    int i = 0;
    float fx = 0;
    for(i=0; i<NPT; i++)
    {
        // fx = 1500 * sin(2*PI * i * 350.0 / Fs) +
        //      2700 * sin(2*PI * i * 8400.0 / Fs) +
        //      4000 * sin(2*PI * i * 18725.0 / Fs);

        //fx = 1500 * sin(2*PI * i * 350.0 / Fs);//频率350hz幅值峰值为1500
        fx = 2700 * sin(2*PI * i * 8400.0 / Fs);//频率8400hz幅值峰值为2700
        adc_buff[i] = ((uint16_t)fx) ;
			
		 uart_output[3] = adc_buff[i] >> 8;
         uart_output[4] = adc_buff[i] & 0xff;
         HAL_UART_Transmit(&huart1, (uint8_t *)uart_output, 5, 0xFFFF);
         HAL_Delay(50);
			//printf("[%d]\r\n", adc_buff[i]);
    }
}

正弦波形-峰值2700-频率8400
在这里插入图片描述

函数波形示例3

#define Fs  44800
#define NPT 256
void InitBufInArray()
{
    int i = 0;
    float fx = 0;
    for(i=0; i<NPT; i++)
    {
        // fx = 1500 * sin(2*PI * i * 350.0 / Fs) +
        //      2700 * sin(2*PI * i * 8400.0 / Fs) +
        //      4000 * sin(2*PI * i * 18725.0 / Fs);

        //fx = 1500 * sin(2*PI * i * 350.0 / Fs);//频率350hz幅值峰值为1500
        //fx = 2700 * sin(2*PI * i * 8400.0 / Fs);//频率8400hz幅值峰值为2700
        fx = 4000 * sin(2*PI * i * 18725.0 / Fs);//频率18725hz幅值峰值为4000
        adc_buff[i] = ((uint16_t)fx) ;
			
		uart_output[3] = adc_buff[i] >> 8;
        uart_output[4] = adc_buff[i] & 0xff;
        HAL_UART_Transmit(&huart1, (uint8_t *)uart_output, 5, 0xFFFF);
        HAL_Delay(50);
		//printf("[%d]\r\n", adc_buff[i]);
    }
}

正弦波形-峰值4000-频率18725
在这里插入图片描述

上述的三种波形整合之后

#define Fs  44800
#define NPT 256
void InitBufInArray()
{
    int i = 0;
    float fx = 0;
    for(i=0; i<NPT; i++)
    {
        fx = 1500 * sin(2*PI * i * 350.0 / Fs) +
             2700 * sin(2*PI * i * 8400.0 / Fs) +
             4000 * sin(2*PI * i * 18725.0 / Fs);

        //fx = 1500 * sin(2*PI * i * 350.0 / Fs);//频率350hz幅值峰值为1500
        //fx = 2700 * sin(2*PI * i * 8400.0 / Fs);//频率8400hz幅值峰值为2700
        //fx = 4000 * sin(2*PI * i * 18725.0 / Fs);//频率18725hz幅值峰值为4000

        adc_buff[i] = ((uint16_t)fx) ;
			
		uart_output[3] = adc_buff[i] >> 8;
        uart_output[4] = adc_buff[i] & 0xff;
        HAL_UART_Transmit(&huart1, (uint8_t *)uart_output, 5, 0xFFFF);
        HAL_Delay(50);
		//printf("[%d]\r\n", adc_buff[i]);
    }
}

波形的相互叠加和整合
在这里插入图片描述

一种波形的FFT结果

一个基波幅度为15,
波形1: 频率100Hz是幅度为10，
波形2: 频率150Hz是幅度为5.5
的信号

#define FFT_LENGTH  1000
        fx = 15 + 10*arm_sin_f32(2*PI*i*100/FFT_LENGTH) + 
        		  5.5*arm_sin_f32(2*PI*i*150/FFT_LENGTH);

下图可以看到FFT之后的频率图
显示的结果和网上基本一致
波形1和波形2的幅值和频率结果正确
在这里插入图片描述

函数整合波形验证-1

//已经验证和网上波形一致-
//基波幅度为15
//100HZ:幅度是10
//150HZ:幅度是5.5
fx = 15 + 10 * arm_sin_f32(2*PI * i* 100/ FFT_LENGTH) + 
			  5.5* arm_sin_f32(2*PI * i* 150/FFT_LENGTH);

原始波形是这样的
在这里插入图片描述
FFT结果是这样的

分解其中一个波形分析-2

fx = 15 + 10 * arm_sin_f32(2*PI * i* 100/ FFT_LENGTH) ;

原始波形是这样的
在这里插入图片描述
FFT结果是这样的

分解其中一个波形分析-3

fx = 15 +5.5* arm_sin_f32(2*PI * i* 150/FFT_LENGTH);

原始波形
在这里插入图片描述
FFT结果是这样的

函数整合波形验证-4

fx  = 100 +  10*   arm_sin_f32(2*PI* 10*i/ Fs)+
		     30*   arm_sin_f32(2*PI* 40*i/ Fs)+
	         50*   arm_sin_f32(2*PI* 80*i/ Fs);

原始波形是这样的
在这里插入图片描述
FFT结果是这样的

（二）理论参考+标题是链接

音频相关基础知识（采样率、位深度、通道数、PCM、AAC）

人耳听觉范围是20Hz~20kHz
根据香农采样定理(也叫奈奎斯特采样定理)
采集频率是原始信号频率的两倍才会尽可能还原
现一般的专业设备的采样频率为44.1kHz。44.1kHz是专业音频中的最低采样率
位深度：一个采样点用多少位表示
码率：码率=采样率x位深度x声道数
双声道的码率是44.1x16x2=1411.2kbps

音频可视化：采样、频率和傅里叶变换

声音振动的两个关键属性是频率和振幅，
频率是指一秒钟振动多少次，对应于音高，频率越高声音越尖锐刺耳。
振幅则表示最大的位移值，对应于音量，振幅越大声音越响
其中 AAC 和 MP3 格式是有损的，而 WAV 则是无损的
傅里叶变换是一个函数，输入一串数字代表样本值，输出一串复数代表频率分量，每个数字具体代表的频率可以根据样本数量和采样频率计算得知，我们不关心复数的方向，我们只关心复数的模，
傅里叶输出的结果是左右对称的，因此只有一半的信息有价值
随音乐跳动的每个柱子，对应的是一个频率或一组频率，而柱子的高度则是频率的分量大小，这两个信息傅里叶变换都能给到

常用的采样率

8,000 Hz - 电话所用采样率, 对于人的说话已经足够
11,025 Hz-AM调幅广播所用采样率
22,050 Hz和24,000 Hz- FM调频广播所用采样率
32,000 Hz - miniDV 数码视频 camcorder、DAT (LP mode)所用采样率
44,100 Hz - 音频 CD, 也常用于 MPEG-1 音频（VCD, SVCD, MP3）所用采样率

在STM32单片机上使用傅里叶解析信号

准备采样数据格式
傅里叶变换
取模
归一化
打印幅值

stm32f1单片机上用FFT测量信号频率（高精度、过程详细）

由理论知识可知，幅值最大值的横坐标对应信号频率，纵坐标对应幅度

【经验分享】如何使用STM32提供的DSP库进行FFT

信号频谱图FFT原理与分析

一个模拟信号，经过ADC采样之后，就变成了数字信号
奈奎斯特采样定理采样频率要大于信号频率的两倍
采样得到的数字信号，就可以做FFT变换了。N个采样点，经过FFT之后，就可以得到N个点的FFT结果。为了方便进行FFT运算，通常N取2的整数次方
假设采样频率为Fs，信号频率F，采样点数为N
那么FFT之后结果就是一个为N点的复数
复数 Z = a + b i
每一个点就对应着一个频率点
模值
这个点的模值，就是该频率值下的幅度特性
假设原始信号的峰值为A
FFT之后第一个点就是直流分量，它的模值就是直流分量的N倍
直流分量 A = 模值 / N
其余点幅值 A = 模值 / (N/2)
频率
第一个点表示 0HZ
最后一个点的下一个点 N+1个点表示频率 Fs
某点的频率 FN
FN = (N-1) * (Fs / N)
频率最小的分辨率 = Fs / N
假如采样频率1024，采样点1024，分辨到1HZ
假如采样频率2048，采样点1024，分辨率2HZ
假如采样频率1024，采样点2048，分辨率0.5HZ
所以要增加频率的精确度则要增加采样点数
FFT之后某点用复数表示 z = a + b i
这个复数的模值 A = 根号 aa + bb
由于FFT结果具有对称性，通常只取前半部分的结果

我所理解的快速傅里叶变换（FFT）

FFT是离散傅立叶变换(DFT)的快速算法
正弦曲线保真度。一个正弦曲线信号输入后，输出的仍是正弦曲线

在这里插入图片描述

产生频谱泄露的主要原因是采样频率和原始信号频率不同步，造成周期的采样信号的相位在始端和终端不连续

如果看了这篇文章你还不懂傅里叶变换，那就过来掐死我吧

在这里插入图片描述

基于STM32的FFT频谱分析+波形识别

FFT官方库

按照FFT官方库的说明，lBufOutArray和lBufInArray都必须是32位的数据类型，其中高16位存储实部，低16位存储虚部

基于STM32F103标准库实现FFT，并实现音乐频谱绘制

从时域到频域的转变，示波器FFT功能及运用

栅栏效应是指对一函数进行采样，即是抽取采样点上的对应的函数值。其效果如同透过栅栏的缝隙观看外景一样，只有落在缝隙前的少数景象被看到
谐波失真测量

谐波

是指对周期性非正弦交流量进行傅里叶级数分解所得到的大于基波频率整数倍的各次分量。谐波分量出现在基频的整数倍处

基波

最长的周期（最小的频率）最大的振幅的正弦波

谐波

频率等于基波频率整数倍的正弦波分量
在这里插入图片描述
在原始振荡波形中，包含基波和谐波
基波：和原始波形一样，最大周期，最小频率，最大振幅的正弦波分量
谐波：频率等于基波频率的整数倍的正弦波分量
例如工频为50HZ的交流电，其二次谐波频率为100HZ，三次谐波为150HZ
在这里插入图片描述

AM-调幅

振幅传播信号-频率保持不变

FM-调频

调制频率

STM32实现FFT，求取幅度频谱

在这里插入图片描述
峰峰值？

STM32上实现FFT算法精准测量正弦波信号的幅值、频率和相位差（标准库）

FFT 变换后，得到的是复数Re+j⋅Im ，包括实部和虚部：

实部 (Re)：与余弦分量相关，决定信号的振幅。
虚部 (Im)：与正弦分量相关，影响信号的相位。
计算幅值和相位：

幅值：表示信号在某频率下的强度，通过以下公式计算：幅值 = sqrt(Re^2 + Im^2)

相位：表示信号相对于参考信号的偏移，通过以下公式计算：相位 = atan2(Im, Re)

STM32F4使用FPU+DSP库进行FFT运算的测试过程一

STM32F4使用FPU+DSP库进行FFT运算的测试过程二

arm_cfft_radix4_instance_f32 * S,float32_t * pSrc)，这个函数在STM32F4xx_DSP_StdPeriph_Lib_V1.4.0库说明中，描述为–不要使用该功能，已经被arm_cfft_f32()替代

STM32F4 FFT 音乐频谱不要太easy！

STM32 DSP库的快速添加基于cubemx 调用，使用DSP库

使用CubeMX 添加
通过keil添加

STM32 CubeMX配置ADC+DMA进行FFT（1）

在RT-Thread studio中为STM32系列开启DSP支持

FFT—快速傅里叶变换算法——STM32F1+DSP库实现（2）

ARM DSP库CMSIS-DSP的使用——以STM32F4浮点FFT为例

【STM32】STM32F4调用DSP库实现FFT运算

音频频谱跳动的柱子是什么呢，Y轴是幅值，X轴是频率
时间花费在验证一个波形函数的FFT正确结果是什么这个事情上
知道这个就可以验证程序结果是否正确
复数 Z = a + bi
复数的模记为 |z|
|z| = 开根号（aa + bb）
直流偏置采集的波形信号没有负值比如12位ADC 那采样结果 - 2048

代码整理（三）

方法一

//原始数据
for (int i = 0; i < FFT_LENGTH; i++)
{
    fft_inputbuf[i * 2]     = adcBuff[i] ;//实部赋值
    fft_inputbuf[i * 2 + 1] = 0;//虚部赋值，固定为0.
}
//傅里叶处理
arm_cfft_f32(&arm_cfft_sR_f32_len1024, fft_inputbuf, 0, 1);
//取模
arm_cmplx_mag_f32(fft_inputbuf, fft_outputbuf, FFT_LENGTH); 

//求幅值
fft_outputbuf[0] /= 1024;

for (int i = 1; i < FFT_LENGTH; i++)//输出各次谐波幅值
{
    fft_outputbuf[i] /= 512;
}


  
* arm_cfft_f32(&arm_cfft_sR_f32_len1024, fft_inputbuf, 0, 1);
输入数据进行傅里叶变换。它的输入参数含义如下

arm_cfft_sR_f32_len1024为傅里叶变换结构体，1024是要运算的点数。
我们在进行其他点数的计算时，比如说32个点的FFT，就可以用arm_cfft_sR_f32_len32。
fft_inputbuf为傅里叶变换需要处理的数据的首地址
第三个参数0是正变换，1是反变换
第四个参数一般是1 为1时是按位取反，系统配置为1。

经过傅里叶变换后的结果仍然为复数，虚部和实部的比可以计算出频率点的相位，
这个在这里不进行考虑。我们直接对复数取模。

取模是实部和虚部的平方和取平均来算，不过我们没必要自己去这样写，
因为DSP库为我们提供了取模函数：
* arm_cmplx_mag_f32(fft_inputbuf, fft_outputbuf, FFT_LENGTH);
参数含义如下：

fft_inputbuf源数据，形式为复数
fft_outputbuf取完模后的数据，形式为实数
FFT_LENGTH是取模的点数

方法二


#definen FFT_LENGTH 4096
float fft_inputbuf[FFT_LENGTH*2];	 // FFT输入数组，用于存放复数
float fft_outputbuf[FFT_LENGTH];	 // FFT输出数组，存放幅值
arm_cfft_radix4_instance_f32 scfft;   // FFT实例结构体

arm_cfft_radix4_init_f32(&scfft,FFT_LENGTH,0,1);
//初始化scfft结构体，设定FFT相关参数     

for(idex=0;idex<FFT_LENGTH;idex++) //
{			
	fft_inputbuf[2*idex]=(u16)(sampledata[idex]);    //生成输入信号实部
	fft_inputbuf[2*idex+1]=0;//虚部全部为0
}
//fft运算
arm_cfft_radix4_f32(&scfft,fft_inputbuf);  
//把运算结果复数求模得幅值
arm_cmplx_mag_f32(fft_inputbuf,fft_outputbuf,FFT_LENGTH);	

//求幅值
fft_outputbuf[0] /= 1024;

for (int i = 1; i < FFT_LENGTH; i++)//输出各次谐波幅值
{
    fft_outputbuf[i] /= 512;
}

* arm_cfft_radix4_init_f32(&scfft,FFT_LENGTH,0,1);
fftLen用于确定FFT变换的长度，
其默认有5种不同的设置，分别是（16/64/256/1024/4096）。
本系统采用4096点傅里叶变换，故将此设置为4096。
ifftFlag这个参数为0或者是1，为0时是傅里叶变换，为1时为傅里叶反变换，系统采用傅里叶变换，固将其配置为0；
bitReverseFlag的取值为0或者1，为1时是按位取反，系统配置为1。
最后这些参数存储在一个结构体指针S里面。
* arm_cfft_radix4_f32(&scfft,fft_inputbuf);  
本系统在计算基-4浮点傅里叶变换运算时就是利用这个函数，
通过指针pSrc传入采集到的输入信号数据（实部和虚部形式）。
这里要特别注意，输入信号应以复数的形式传入，实部放在奇数位，虚部放在偶数位。
进行4096点FFT运算，则pSrc指向4096个复数，实际应用中，
系统通过A/D采样，这样偶数位全部为0，奇数位为A/D转换器采集到的电压值。
同时FFT变换后的数据也存放在pSrc里面，pSrc必须不小于fftLen的长度。
另外，S结构体指针参数是由先该函数设置好，然后传入该函数的。
* void arm_cmplx_mag_f32(float32_t * pSrc,float32_t * pDst,uint32_t numSamples)

系统通过这个函数来计算复数模值，
可以对第二个函数* pSrc指向的输出结果进行傅里叶变换进行取模操作。
pSrc指向第二个函数傅里叶变换的结果，为复数输入数组指针；
pDst为输出数组指针，存储取模后的值；
numSamples是进行取模的点数，
系统进行4096点的傅里叶变换,
在此配置为4096即可。

方法三
FFT官方库

	for(i=0;i<NPT;i++)
  {
    lBufInArray[i]=ADC_Value[i]<<16;
  }
  cr4_fft_1024_stm32(lBufOutArray, lBufInArray, NPT);   
  
  GetPowerMag();
  
void cr4_fft_256_stm32(void *pssOUT, void *pssIN, u16 Nbin);

        参数1：*pssOUT一个用于存放结果的地址，通常为一个长度为Nbin/2的数组。

        参数2：*pssIN一个用于存放待变换数据的地址，通常为一个长度为Nbin的数组。

        参数3：指定待变换的点数。尽管参数在函数定义中被明确指定为 256 点，但这个参数可能用于函数的一般化设计，以便将来可以轻松地修改或扩展为其他点数的 FFT 变换。

　　调用该函数之后，在pssOUT数组中就存放了进行FFT运算之后的结果数据。该数组中每个元素的数据格式为：高16位存储虚部，低16位存储实部
/******************************************************************
函数名称:GetPowerMag()
函数功能:计算各次谐波幅值  [short 的范围，是-32767 到 32767 。也就是 -(2^15 - 1)到(2^15 - 1)。]
参数说明:
备　　注:先将lBufOutArray分解成实部(X)和虚部(Y)，然后计算幅值(sqrt(X*X+Y*Y)
*******************************************************************/
void GetPowerMag(void)
{
    signed short lX,lY;                                                  //算频率的话Fn=i*Fs/NPT  //由于此处i是从0开始的，所以不需要再减1
    float X,Y,Mag;                                                      
    unsigned short i;
    for(i=0; i<NPT/2; i++)                                                  //经过FFT后，每个频率点处的真实幅值  A0=lBufOutArray[0]/NPT
    {                                                                       //                                 Ai=lBufOutArray[i]*2/NPT
        lX  = (lBufOutArray[i] << 16) >> 16;  //lX  = lBufOutArray[i];
        lY  = (lBufOutArray[i] >> 16);
                                    
        X = NPT * ((float)lX) / 32768;//除以32768再乘65536是为了符合浮点数计算规律,不管他
        Y = NPT * ((float)lY) / 32768;
        Mag = sqrt(X * X + Y * Y) / NPT;
        if(i == 0)
            lBufMagArray[i] = (unsigned long)(Mag * 32768);   //0Hz是直流分量，直流分量不需要乘以2
        else
            lBufMagArray[i] = (unsigned long)(Mag * 65536);
    }
}