STM32 FFT频谱分析:频率分辨率、幅值还原与嵌入式工程实践
发布时间:2026/10/24 9:36:06 锦皓数字建站

简介一份面向嵌入式初学者的STM32 FFT算法实现说明文档。内容从DFT与FFT基本原理讲起结合STM32官方DSP库以识别21Hz正弦信号为例完整演示原始数据生成、FFT变换调用、频率与幅值计算、频谱打印及两种算法耗时对比官方约9ms、非官方约411ms并给出了关键代码与处理思路便于在振动监测、音频分析等场景中迁移使用。压缩包内含1个docx格式文档大小244KB整体结构完整、目录清晰已有2700人学习浏览。通过阅读可以掌握FFT在STM32上的落地流程包括实部/虚部存储格式、采样点数与频率分辨率的关系、直流分量修正等易错细节适合作为入门学习和项目开发参考。1. 一个反直觉现象FFT 算出来的不是频率表在 STM32 上把 FFT 跑完你拿到的不是一张“频率-幅值”表而是 1024 个格式别扭的复数。21Hz 的正弦信号不会出现在第 21 号元素上而是躺在第 210 号元素里2048 的振幅并不会直接显示得自己乘 2、再除 N还得单独照顾 0 号元素这个直流特例。这套映射规则正是 FFT 在嵌入式系统里能 work 的核心机制。本文用一个识别正弦信号频率的 STM32 工程切口信号按 102.4Hz 采样 1024 点FFT 稳定还原出 21Hz 主频和 2048 幅值再对比官方库 9ms 与通用实现 411ms 的耗时差距。搞懂这些拿到运算结果后你才不至于两眼一黑。2. 采样参数设计为什么 0.1Hz 分辨率必须攒够 10 秒数据STM32 例程里最先值得看的两个数字是 102.4Hz 和 1024。前者是采样频率 Fs后者是数据点数 NPT它们不是随手填的而是围绕“要把 0.1Hz 分辨出来”这个硬指标反推出来的。很多开发者以为采样率越高分辨率越好实际上 FFT 的频率分辨率跟采样率没关系它只由采样总时长决定。这个反直觉结论是理解整条采样链的第一步。2.1 频率分辨率被采样时间锁死不是被采样率锁死DFT 把 N 个采样点映射成 N 个频域格点频率轴从 0Hz 均匀铺到采样率 Fs格点间距是 Fs/N。把 Fs 写成 N/T其中 T 是采样总时长格点间距就变成 1/T。也就是说想分辨 0.1Hz就得攒够 10 秒的数据想分辨 0.01Hz就得采集 100 秒。频率分辨率是采样时间的倒数这是一条逃不掉的物理约束。例程里 Fs102.4Hz、NPT1024采样时长正好是 1024÷102.410 秒分辨率102.4÷10240.1Hz。这三个数满足“分辨率1/采样时间”的等式。如果把采样率改成 1024Hz、点数仍是 1024采样时长缩到 1 秒分辨率立刻退化成 1Hz——虽然数据点一样多但 21Hz 和 21.5Hz 就分不开了。下表可以直观看出这个权衡采样率 Fs点数 NPT采样时长频率分辨率能否分辨 21Hz 与 21.3Hz102.4Hz102410s0.1Hz能102.4Hz2562.5s0.4Hz勉强格点不齐1024Hz10241s1Hz不能所以当你决定“我要 0.1Hz 的分辨率”时10 秒采集时长已经是板上钉钉的事能调整的只是点数与采样率的搭配。点数越多频率格点越密、可分析的最高频率也越高但 RAM 占用和 FFT 耗时同步上涨。这就是为什么官方库把点数做成 64、256、1024 三档——每一档都是基 4 算法的最优规模这点后面细讲。2.2 奈奎斯特那半边频谱只取前一半的机制FFT 输出的 N 个复数中只有前 N/2 个点携带独立信息后半部分是对称镜像。原因是实信号的频谱具有共轭对称性采样率 Fs 对应的最高可分析频率是 Fs/2也就是奈奎斯特频率超过它的成分会被折叠回低频区。例程里 Fs102.4Hz奈奎斯特频率就是 51.2HzFFT 结果只在 0 到 51.2Hz 之间有物理意义。这就是为什么幅值计算的循环只跑到 NPT/2for(i0;iNPT/2;i)后面那 512 个点不过是前 512 个点的镜像重复计算只是浪费 CPU。这里有个新手常犯的误会以为只算一半会丢信息其实完整的 N 点 FFT 输出里后一半存放的是负频率分量对实数信号来说它和正频率部分互为正共轭。取前一半正好把 0Hz 到奈奎斯特频率这段看得清清楚楚。2.3 例程参数为什么这么配21/15/25Hz 全落在整数格点上频率格点的映射公式是 Fn i × Fs / NPTi 从 0 数到 NPT/2-1。把 Fs102.4、NPT1024 代进去每个格点的频率就是 i×0.1Hz。于是三个信号分量的落点非常整齐15Hz 落在第 150 个格点21Hz 落在第 210 个格点25Hz 落在第 250 个格点。这不是巧合。设计例程时特意选了能被 0.1Hz 整除的信号频率目的就是让每个频率分量都精确对准整数格点。只有这样FFT 结果才能在单个格点上形成干净的尖峰能量不会泄漏到相邻格点。如果你的信号频率是 21.3Hz而格点间距是 0.1Hz21.3Hz 会同时泄漏到 21.2 和 21.4 两个格点上主峰幅值随之下降频谱图出现双肩包。这个现象在避坑章节里会展开但采样参数的设计原则现在就该记住先定分辨率需求再定采样时长然后选择让目标频率尽量对准整数格点的 Fs 与 NPT 组合。3. 数据格式与库函数32 位打包和基 4 蝶形如何让 FFT 快到 9ms官方 DSP 库的接口看起来是一个函数调用但前后要做不少手脚输入数组用什么格式、实部虚部怎么塞进 32 位整数、库内部为什么能快到令人咋舌。这些细节直接决定你调用 cr4_fft_1024_stm32 之后拿到的数据能不能用。下面把这条链路逐层拆开。3.1 32 位长整型打包高 16 位实部低 16 位虚部先看例程生成原始信号的代码。这个函数模拟了 ADC 采样结果主信号是 21Hz、幅值 2048 的正弦叠加 15Hz 幅值 1000 和 25Hz 幅值 100 的次分量最后加 0 到 100 的随机噪声。for(i0; iNPT; i) { // FX 直流2048 主信号2048*sin(21Hz) 15Hz分量 25Hz分量 随机噪声 float fx 2048 2048*sin(PI2*i*21/Fs) 1000*sin(PI2*i*15/Fs) 100*sin(PI2*i*25/Fs) rand()%100; lBUFIN[i] ((s16)fx) 16; // 高16位放实部低16位清零当虚部 }关键就在最后一行fx 是 float 类型的采样值先强转成 s16 截断小数部分再左移 16 位塞进 32 位整型的高半段。低 16 位自动补 0表示虚部初始为 0。之所以要预留虚部是因为蝶形运算本质是复数乘法即便输入是实信号运算中间过程的复数结构也省不掉。官方库为了统一接口干脆规定输入输出都用“高 16 位实部、低 16 位虚部”的 Q15 定点格式打包。Q15 是定点数格式里最常见的一种把 [-1, 1] 区间的实数映射到 [-32768, 32767]。例程里 fx 最大约 5295只用了 Q15 量程的六分之一左右剩余的动态范围给噪声和叠加分量留了余量。如果你把满幅 4096 的 ADC 信号直接左移却忘了检查最大值是否超过 32767就会发生定点溢出频谱上所有峰都会变得面目全非。3.2 基 4 蝶形与旋转因子表官方库快的三层原因官方 DSP 库提供了三个汇编函数cr4_fft_64_stm32、cr4_fft_256_stm32、cr4_fft_1024_stm32。注意点数清一色是 4 的幂644³2564⁴10244⁵。这个前缀 cr4 指的就是面向 Cortex-M4 内核优化的基 4 FFT。基 4 的含义是每一级蝶形运算同时处理 4 路输入相比基 2 每级处理 2 路蝶形级数直接减半。1024 点用基 2 需要 10 级用基 4 只需要 5 级。级别越少中间层的读写次数和复数乘加次数越少。再加上三点叠加旋转因子全部提前算好存成查找表运行时不调用任何 sin/cos 函数核心循环用汇编手写寄存器分配完全针对 Cortex-M4 的流水线特性定制定点运算避免浮点开销。三层优化叠加9ms 的 1024 点 FFT 就是这么来的。提示官方库只提供 64、256、1024 三种点数本质原因是基 4 算法要求点数是 4 的幂。如果你的应用需要 512 点或 2048 点要么补零到最近的 4 次幂要么换用非官方实现。3.3 频谱读法第 210 号格点为什么就是 21HzFFT 做完后lBUFOUT 数组里每个元素对应一个频率格点。用前面那个映射公式反算i 是数组下标i×Fs/NPT 就是该格点的物理频率。i210 时210×102.4/102421.0Hz分毫不差。这个例程的信号分量与格点映射关系如下信号分量真实幅值格点索引 i格点频率 Fn直流204800Hz15Hz 正弦100015015.0Hz21Hz 正弦204821021.0Hz25Hz 正弦10025025.0HzFFT 把时域信号拆成这四路分量每路能量集中到自己对应的格点上。21Hz 主分量的模值会在第 210 号格点出现明显尖峰相邻格点接近 0。这正是“找频率”的完整逻辑先算每个格点的物理频率再看哪个格点模值最大就说明信号里存在对应频率的分量。4. 幅值还原为什么模值要乘 2 再除 NFFT 输出的复数模值并不直接等于信号振幅。例程里的核心代码如下注释里的公式是 lBUFMAG[i] sqrt(lX² lY²) × 2 / NPT。为什么是“乘 2 再除 N”这是本章要拆解的机制。4.1 从 DFT 累加看 N/2 倍峰值的来路DFT 的求和过程可以理解为把信号与一组同频率的复指数逐点相乘再累加。一个振幅为 A 的正弦信号用欧拉公式展开后是两个旋转方向相反的复指数各自振幅为 A/2。当 DFT 的扫描频率与信号频率对齐时其中一个复指数在 N 次累加里始终保持同相逐次叠加成 N×(A/2)另一个反向旋转的复指数则均匀转圈累加后相互抵消接近 0。所以对齐频率格点上FFT 输出的模值是 A×N/2。这就是“乘 2 除 N”的数学根源乘 2 是把只保留了一半能量的单边频谱补回来除 N 是把 N 次同相累加的效果归一化回原始振幅。三个运算各司其职少一个都还原不出真实幅值。4.2 直流分量是特例0 号格点不需要乘 2第 0 号格点是直流分量它没有正负频谱之分。直流信号是常数DFT 求和时每个采样点的贡献都是同相的N 个点累加出 N 倍的直流值。还原直流分量只需要除 N不需要乘 2。例程里专门做了这个处理lBUFMAG[0] lBUFMAG[0] / 2; // 直流分量不乘2按2折算后要还原这段代码的逻辑是幅值计算统一按“乘 2 再除 N”流程处理了所有格点但直流格点本身只该除 N所以多乘的那个 2 要再除回去。如果忘记这一步0Hz 处的读数会是真实直流幅值的两倍频谱图上直流尖峰比所有信号峰都高很容易误判。4.3 逐行拆解 dsp_asm_powerMag定点数如何变回真实幅值官方库输出的是 Q15 定点的实部和虚部要变成可打印的幅值需要经过换算。例程的幅值计算函数如下void dsp_asm_powerMag(void) { s16 lX, lY; // lX取低16位是虚部lY取高16位是实部 u32 i; for(i0; iNPT/2; i) // 频谱共轭对称只算前一半格点 { lX (lBUFOUT[i] 16) 16; // 低16位虚部带符号 lY (lBUFOUT[i] 16); // 高16位实部 { float X NPT * ((float)lX) / 32768; // Q15定点还原 float Y NPT * ((float)lY) / 32768; float Mag sqrt(X*X Y*Y) / NPT; // 按DFT累加关系归一化 lBUFMAG[i] (u32)(Mag * 65536); // 放大为定点数方便打印 } } lBUFMAG[0] lBUFMAG[0] / 2; // 直流分量不需要乘2 }这里的 32768 是 Q15 格式的满幅值lX、lY 除以 32768 才还原成真实的浮点比例。乘以 NPT 再开方、最后再除以 NPT对应的正是第 4.1 节推导的累加关系。最后乘以 65536 存成 u32是为了让幅值以“定点放大”的形式方便串口 printf 输出。值得注意的是lBUFMAG[i] 里存的是放大 65536 倍后的定点值打印后如果要还原物理幅值需要除以 65536——但如果你想直接看相对大小比较频谱尖峰这种放大格式反而更直观。4.4 只算一半为什么够共轭对称的工程红利实信号 FFT 后第 k 个格点和第 N-k 个格点的模值相等、相位相反频谱关于奈奎斯特频率镜像对称。因此循环只跑 NPT/2把后一半的运算量整个省掉。这不仅是性能优化也是理解幅值公式里“乘 2”的物理含义的关键第 k 个格点的模值只携带了正频率那部分的能量另一半能量放在镜像格点里所以必须乘 2 才能补全。5. 避坑指南FFT 结果翻车的四个典型现场FFT 的坑不在调用本身而在调用前后那些你以为理所当然的假设。以下四个场景是嵌入式环境里反复出现的翻车现场每一条都是“现象—原因—解决”的完整链路。5.1 主峰幅值明显减半21Hz 旁边多出两个对称包现象信号明明是 21Hz、振幅 2048频谱主峰却只有 1024 左右而且 21Hz 两侧各出现一个小峰像两个肩膀。原因输入信号频率与整数格点错位。FFT 的幅值公式“乘 2 除 N”默认信号频率精确落在某一格点上但真实信号频率往往在格点之间能量被泄漏到相邻的几个格点形成主瓣变宽、幅度下降的栅栏效应。例程特意选了 15/21/25Hz 这些能被 0.1Hz 整除的频率就是为了避开这个坑。解决调整 Fs 和 NPT 的搭配让目标频率尽量对准整数格点。如果信号频率无法精确控制就加汉宁窗或布莱克曼窗抑制旁瓣但要注意加窗本身会改变主峰幅值需要按窗函数的幅值恢复系数做补偿。赶上只差半个格点的时候最简单的办法是加密格点——增加 NPT 即拉长采样时间。5.2 频谱整体形状崩坏0Hz 处顶到天花板现象不管输入什么样的正弦信号频谱都是低频区一团巨大尖峰其他区域几乎平的打印出的幅值完全对不上。原因输入数组的格式填错了。常见的有两种一是把浮点信号直接强转 s32 而没有先转 s16 再左移虚部位置被随机噪声占满复数结构被破坏二是实部虚部放反把信号塞进了低 16 位FFT 算法按“高 16 位实部”解析时拿到的是残缺数据。解决严格按 lBUFIN[i] ((s16)fx) 16 的格式填充确保低 16 位为 0。调完之后喂一个已知频率和幅值的正弦信号验证谱峰位置和幅值是否能对上再接入真实 ADC 数据。这是一条能省下大量调试时间的自测习惯。5.3 高频分量折返25Hz 旁边出现一条幽灵峰现象25Hz 的幅值明明只有 100频谱里却多出一条不属于任何已知分量的峰且位置随采样率变化。原因采样率没有满足奈奎斯特条件。Fs102.4Hz 的上限是 51.2Hz25Hz 本身安全但如果信号里混入了超过 51.2Hz 的高频噪声它不会消失而是会折叠回低频区伪装成一条“幽灵峰”。ADC 前没有抗混叠滤波时这种情况尤其隐蔽。解决硬件上在 ADC 前端加低通滤波器把高于 Fs/2 的成分滤掉软件上先做一次粗略的 FFT 观察高频段是否还有残余能量。工程经验是采样率取最高分析频率的 2.5 倍以上给滤波器的过渡带留出空间。5.4 耗时测试的数字没有说服力9ms 和 411ms 是怎么被污染现象官方库测出 9ms、非官方库测出 411ms数字很漂亮但把 printf 加进测试区间后时间暴涨到一秒以上耗时数据失去意义。原因printf 重定向到串口是阻塞发送波特率 9600 时一个字节差不多 1ms打印 1024 个点就是一秒量级。定时器计时区间如果覆盖了打印操作测的就不是 FFT 而是串口。解决先关定时器计数再执行串口打印。例程用 TIM6 计时时应该把 start 放在 FFT 调用前、stop 放在幅值计算结束后printf 全部挪到计时区间之外。另外确认 TIM6 的中断标志已经清理否则第二次测量会把上次的计数累积进来。6. 9ms 与 411ms官方库和非官方库的设计取舍两者 44 倍的耗时差距不是“官方代码更聪明”一句话能带过的。拆开看差距来自三个层次算法结构、计算策略、实现方式。官方库用基 4 蝶形把 1024 点压缩到 5 级非官方实现大多用基 2 结构总共 10 级蝶形级数翻倍官方库把所有旋转因子提前算成查找表运行时只查表非官方实现每级实时调用 sin/cos 计算旋转因子而三角函数在无 FPU 的内核上本身就是几十微秒量级的开销官方库核心循环是手写汇编非官方实现是 C 循环编译器的寄存器分配再优化也追不上针对流水线手动调度的汇编。这个差距背后是取舍不是优劣。官方库快但点数被锁死在 4 的幂64、256、1024你想做 512 点或 2048 点就得自己补零或移植。非官方实现慢但它可以用任意点数也方便插入窗函数、相位校正等自定义处理。实时性敏感的场合比如电机控制里的电流谐波分析只能选官方库需要灵活的频谱分析工具慢一点但能加窗、能算任意点数非官方实现反而更合适。从那以后我换任何一块新的单片机平台都会先跑一遍同一份 FFT 测时程序确认官方库的实际耗时跟预期一致再确认信号频率落在整数格点上最后才接真实 ADC 数据。这套流程走完FFT 才算真正在这块芯片上落地希望帮到你。本文还有配套的精品资源点击获取
锦
锦皓数字建站
深耕本土企业品牌数字化升级,专注原创端正雅致商务官网,从视觉设计到稳定运维全程保驾护航。