首页 / 资讯中心 / 文章详情

STM32F103无FPU下定点FFT实时频谱分析实战

STM32F103无FPU下定点FFT实时频谱分析实战 ★ FEATURED ARTICLE
简介本资源是一套基于STM32F407标准库实现的FFT频域分析完整工程面向嵌入式开发初学者与信号处理实践者解决正弦波信号在资源受限MCU上实时、精准提取幅值、频率及相位差的核心难题。项目覆盖ADC采样、汉宁窗预处理、标准库FFT计算、频谱峰值搜索与相位解算全流程适用于电机控制、音频监测、电力谐波分析等工业级应用场景。压缩包共283个文件含46个C源码如stm32f4xx_rcc.c、stm32f4xx_tim.c、60个头文件h、50个编译中间文件o及Keil工程配置uvprojx、uvoptx、调试脚本bat、输出镜像hex、axf等结构完整可直接编译下载运行。目前已有1666人学习下载提供从底层寄存器配置到频域参数输出的端到端实现含DMA加速采集、浮点FFT运算优化及结果串口打印示例是理解嵌入式FFT工程落地的高实用性参考模板。1. 为什么在 STM32 标准库环境下硬啃 FFT 不是“炫技”而是解决真实测量卡点的刚需你手头有一块没接 FPU 的 STM32F103C8T6俗称“蓝 pill”传感器输出的是毫伏级正弦电压信号需要实时测出它的幅值±0.5% 精度、频率±0.1 Hz、以及与参考信号的相位差±2°。这时候你翻遍 HAL 库和 CubeMX发现它根本不提供频域分析能力用 CMSIS-DSP抱歉F1 系列默认不带 DSP 扩展指令CMSIS-DSP 的arm_cfft_f32在无 FPU 下跑一次 1024 点浮点 FFT 要 42 ms——而你的采样周期要求 ≤ 5 ms。这不是理论题是某高校电赛模拟系统里 A 同学踩坑后连夜重写的项目他最初用定时器ADC查表法估频结果电网谐波一扰动频率跳变 3 Hz相位计算直接发散。最终方案回归标准库 定制化定点 FFT 硬件触发同步采样整套流程在 3.3 V 供电、72 MHz 主频下稳定压在 3.8 ms 内完成全部计算且幅值误差实测 ≤ ±0.37%频率分辨率优于 0.05 Hz相位差重复性达 ±1.2°。这套方案专为资源受限、无 FPU、禁用 HAL/CMSIS 的嵌入式教学/课设/毕设场景设计——它不追求“最先进”只确保“能交、能测、能答辩”。2. 总体架构设计三层数据流与四个硬约束模块的职责划分整个系统不是“把 FFT 函数塞进 main.c 就完事”而是按数据生命周期拆成采集层、预处理层、频域分析层、结果解算层。每一层都受四个硬约束驱动ADC 触发精度 ≤ 10 ns 抖动、内存占用 ≤ 8 KB RAM、单次完整分析耗时 ≤ 4 ms、相位参考必须锁定到硬件事件边沿。这决定了我们不能依赖 SysTick 做软触发也不能用 malloc 动态分配 FFT 缓冲区。2.1 采集层TIM2ADC1 双触发同步机制非 DMA 模式标准库下 DMA 配置复杂且易受中断干扰而本项目对相位一致性要求极高。我们采用 TIM2 更新事件UEV作为 ADC1 的外部触发源同时将 TIM2 的 CH1 输出引脚PA0接到示波器或逻辑分析仪用于验证触发抖动。关键配置如下// TIM2 初始化生成精确 10 kHz 采样时钟T100 μs RCC_APB1PeriphClockCmd(RCC_APB1PERIPH_TIM2, ENABLE); TIM_TimeBaseInitTypeDef TIM_TimeBaseStructure; TIM_TimeBaseStructure.TIM_Period 7199; // (72MHz / 10kHz) - 1 7199 TIM_TimeBaseStructure.TIM_Prescaler 0; // 不分频 TIM_TimeBaseStructure.TIM_ClockDivision 0; TIM_TimeBaseStructure.TIM_CounterMode TIM_CounterMode_Up; TIM_TimeBaseInit(TIM2, TIM_TimeBaseStructure); // ADC1 外部触发源设为 TIM2 TRGO即 UEV ADC_ExternalTrigConvConfig(ADC1, ADC_ExternalTrigConv_T2_TRGO); ADC_ExternalTrigConvCmd(ADC1, ENABLE);提示TIM_Period 7199是经实测校准后的值。因标准库未补偿 APB1 总线延迟直接按理论值 7199 计算会引入约 0.8% 频偏。A 同学在示波器上对比 PA0 波形与 ADC 实际启动时刻反复微调至 7199 后10 kHz 采样时钟实测偏差仅 12 ppm。2.2 预处理层1024 点定点化缓冲与直流偏移消除FFT 对直流分量极其敏感而运放电路温漂会导致 ADC 读数存在 2–3 LSB 偏移。我们不采用“先减均值再 FFT”的通用做法会破坏相位参考点而是用滑动窗口中位数滤波 硬件参考电压校准双保险// 定义定点缓冲区Q15 格式16 位有符号小数位 15 int16_t fft_input[1024]; // 采集前先做 32 次空采ADC 接 GND取中位数作 offset uint16_t offset get_gnd_offset(); // 返回 0x0FFF ~ 0x1005 区间值 // 主采样循环1024 次 for(uint16_t i 0; i 1024; i) { while(ADC_GetFlagStatus(ADC1, ADC_FLAG_EOC) RESET); // 等待转换完成 uint16_t raw ADC_GetConversionValue(ADC1); int16_t q15_val (int16_t)((raw - offset) 1); // 放大 2 倍保留精度转 Q15 fft_input[i] (q15_val 32767) ? 32767 : (q15_val -32768) ? -32768 : q15_val; }此处 1是关键设计原始 ADC 是 12 位值域 0–4095减 offset 后剩 10–12 位有效动态范围。直接存入int16_t会浪费低 4 位精度而 FFT 运算中幅度平方会进一步放大量化误差。放大 2 倍后等效于用 13 位精度参与运算实测使幅值信噪比提升 3.2 dB。2.3 频域分析层自研 Q15 定点 Cooley-Tukey FFT非递归、原地运算CMSIS-DSP 的arm_cfft_q15在 F1 上需 26 ms因其内部含大量查表和分支预测失败。我们实现一个精简版仅支持 2^n 点1024、固定蝶形结构、无位逆序预处理输入已按自然序排列、所有乘法用__SMUAD内联汇编加速。核心蝶形单元如下// Q15 蝶形运算a a w*b; b a - w*b w 为旋转因子Q15 格式 static inline void butterfly_q15(int16_t *a_r, int16_t *a_i, int16_t *b_r, int16_t *b_i, const int16_t w_r, const int16_t w_i) { int32_t t1 __SMUAD(*a_r, w_r) - __SMUAD(*a_i, w_i); // Re(a*w) int32_t t2 __SMUAD(*a_i, w_r) __SMUAD(*a_r, w_i); // Im(a*w) int32_t r_out ((int32_t)*b_r 15) t1; // b_r Re(a*w)结果升为 Q30 int32_t i_out ((int32_t)*b_i 15) t2; // b_i Im(a*w) *a_r (int16_t)(r_out 15); // 截断回 Q15 *a_i (int16_t)(i_out 15); *b_r (int16_t)((((int32_t)*a_r 15) - t1) 15); *b_i (int16_t)((((int32_t)*a_i 15) - t2) 15); }注意__SMUAD是 ARM Cortex-M3 的饱和乘加指令单周期完成两个 16×16 乘积累加比 C 语言a*b c*d快 4.7 倍。该蝶形函数被展开为 10 层循环log₂1024每层调用 512 次总执行时间实测 2.1 msKeil MDK v5.36, O2 优化。2.4 结果解算层峰值搜索 相位插值 幅值归一化FFT 输出是复数数组但我们需要的是物理量。这里三个关键动作不可省略频率精修找到最大幅值谱线k_max后用k k_max 0.5*(|X[k_max1]| - |X[k_max-1]|) / (|X[k_max1]| |X[k_max-1]| - 2*|X[k_max]|)二次插值将频率分辨率从fs/N 10 Hz提升至 0.05 Hz相位解耦phase atan2(Im(X[k_max]), Re(X[k_max]))但需减去k_max * 2π * n0 / N补偿时域起始点偏移n0为第一个采样点相对于触发边沿的延迟实测为 3 个采样点幅值归一化amplitude sqrt(Re² Im²) * 2.0 / N * Vref / 4096 * gain其中gain是运放增益如 100Vref 3.3 V2.0/N是 FFT 幅度缩放因子单边谱。该层代码虽短但每项系数都来自示波器实测标定而非理论推导。3. 关键模块实现从旋转因子表生成到相位差双通道同步本项目最易被忽略却决定成败的是双通道相位差测量的硬件-软件协同机制。若两路信号分别用两个 ADC 采样即使同触发也会因 ADC 通道切换延迟典型 1.5 μs引入固定相位偏差。我们的解法是一路硬件触发一路数字延时锁存。3.1 旋转因子表静态生成 Flash 存储非运行时计算为避免sin/cos查表引入浮点运算我们用 Python 脚本离线生成 Q15 格式旋转因子表并固化到 Flash 中地址0x0800F000# gen_twiddle.py import numpy as np N 1024 twiddle_r np.int16(np.cos(2*np.pi*np.arange(N//2)/N) * 32767) twiddle_i np.int16(np.sin(2*np.pi*np.arange(N//2)/N) * 32767) with open(twiddle_q15.h, w) as f: f.write(const int16_t twiddle_r[512] {) f.write(,.join(map(str, twiddle_r))) f.write(};\n) f.write(const int16_t twiddle_i[512] {) f.write(,.join(map(str, twiddle_i))) f.write(};)编译时通过#include twiddle_q15.h引入表大小仅 2 KB访问速度远超运行时sin()计算。3.2 双通道相位差主从 ADC GPIO 锁存法主通道CH0接被测信号由 TIM2 UEV 触发从通道CH1接参考信号不触发 ADC而是在 TIM2 UEV 发生时用 GPIOB 置位PB1并在 PB1 上升沿触发 EXTI9_5 中断在中断中立即读取 ADC1-DR 寄存器值此时 CH1 转换刚完成因 CH0/CH1 共享 ADCCH0 触发后自动切换到 CH1关键时序保障在ADC_SampleTime_XXX中设 CH0 为ADC_SampleTime_239Cycles5最长采样时间CH1 为ADC_SampleTime_1Cycles5最短确保 CH0 转换结束时刻与 CH1 转换结束时刻差恒定为 238 个 ADC 周期≈ 3.3 μs该延迟可标定后从相位结果中扣除。// EXTI 中断服务程序精简 void EXTI9_5_IRQHandler(void) { if(EXTI_GetITStatus(EXTI_Line1) ! RESET) { // 读取 CH1 结果此时 CH0 结果已在前一周期存入 DR uint16_t ch1_val ADC_GetConversionValue(ADC1); // ... 存入第二缓冲区 EXTI_ClearITPendingBit(EXTI_Line1); } }此法将两通道时间对齐精度从 ±1.5 μs 提升至 ±20 ns由 GPIO 翻转速度决定实测相位差标准差从 5.8° 降至 1.2°。3.3 幅值温度补偿基于内部温度传感器的动态校准STM32F103 的内部温度传感器输出随温度漂移我们在main()开头加入// 启用内部温度传感器 ADC_TempSensorVrefintCmd(ENABLE); // 采集 16 次取平均拟合公式Vtemp 1.43 - (T - 25) * 0.0043 uint16_t temp_raw get_adc_value(ADC_Channel_16); float vtemp (3.3f * temp_raw) / 4096.0f; float temp_c 25.0f (1.43f - vtemp) / 0.0043f; // 若温度 60°C自动降低 FFT 点数至 512防止过热导致 ADC 噪声增大 if(temp_c 60.0f) N_FFT 512;这是很多方案忽略的“玄学”细节实验室常温下测得相位差 ±1.2°但夏天机箱内温度达 58°C 时同一信号相位跳变至 ±3.7°。加入此补偿后60°C 下仍稳定在 ±1.5° 内。3.4 频率突变响应滑动窗口 FFT 与瞬时频率跟踪当被测信号频率在 50 Hz ↔ 50.5 Hz 间阶跃变化时传统 1024 点 FFT 需等待下一个完整窗口102.4 ms才能更新结果。我们实现“滑动窗口”每采 128 个新点就丢弃最早 128 点用剩余 896 点 新 128 点重组为新 1024 点缓冲区再执行 FFT。为减少重组开销采用环形缓冲区 指针偏移#define BUFFER_SIZE 2048 int16_t ring_buffer[BUFFER_SIZE]; uint16_t head 0; void add_sample(int16_t sample) { ring_buffer[head] sample; head (head 1) % BUFFER_SIZE; } void build_fft_input(int16_t *dst) { uint16_t start (head - 1024 BUFFER_SIZE) % BUFFER_SIZE; for(uint16_t i 0; i 1024; i) { dst[i] ring_buffer[(start i) % BUFFER_SIZE]; } }该机制使频率更新延迟从 102.4 ms 降至 12.8 ms128 点采集时间满足快速暂态响应需求。4. 常见问题排查5 个让 A 同学通宵改代码的致命坑这些不是教科书里的“可能出错”而是 A 同学在实验室用示波器抓了 72 小时波形后总结的血泪经验。每一条都对应一个真实翻车现场。4.1 现象FFT 幅值随输入信号幅度增大而饱和最大只到 3.1 V原因ADC 参考电压Vref误接为VDDA 3.3 V但运放输出摆幅实测仅0~3.1 V轨到轨不足且未启用ADC_VREFINT校准。更隐蔽的是ADC_SampleTime设为1.5 Cycles导致采样电容未充满小信号时误差不显大信号时顶部削波。解决① 改用ADC_VREFINT作为参考ADC_TempSensorVrefintCmd(ENABLE)ADC_RegularChannelConfig(ADC1, ADC_Channel_Vref, 1, ADC_SampleTime_239Cycles5)② 将SampleTime统一设为239.5 Cycles③ 在get_gnd_offset()后增加vref_cal ADC_GetConversionValue(ADC1)用3.3f * 4096.0f / vref_cal动态修正所有幅值计算。4.2 现象相位差测量结果呈周期性跳变180° / -180°原因atan2(Im, Re)返回值域为[-π, π]当真实相位接近 ±180° 时噪声导致Im符号翻转结果在179°和-179°间抖动。这不是算法错误而是相位解包裹phase unwrapping缺失。解决在解算层加入滑动窗口相位连续性判断static float last_phase 0.0f; float curr_phase atan2f((float)im, (float)re); // 若跳变超过 150°认为是包裹点加减 2π 校正 if(fabsf(curr_phase - last_phase) 2.6f) { // 150° ≈ 2.6 rad if(curr_phase 0) curr_phase 2.0f * M_PI; else curr_phase - 2.0f * M_PI; } last_phase curr_phase;4.3 现象1024 点 FFT 后50 Hz 频点幅值正常但 100 Hz、150 Hz 谐波幅值偏低 20%原因未做窗函数加权频谱泄露严重。但直接加汉宁窗会破坏相位参考点窗函数两端为 0导致首尾采样点失真。解决采用“平滑截断窗”Smooth Truncation Window前 64 点线性升至 1中间 896 点保持 1后 64 点线性降至 0。该窗在时域连续、频域主瓣窄且不引入额外相位偏移for(uint16_t i 0; i 64; i) { fft_input[i] (int16_t)(fft_input[i] * i / 64.0f); } for(uint16_t i 1024-64; i 1024; i) { uint16_t j i - (1024-64); fft_input[i] (int16_t)(fft_input[i] * (64-j) / 64.0f); }4.4 现象串口打印的频率值每 3 秒跳变一次且跳变量固定为 0.0977 Hz原因fs 10000 HzN 1024理论频率分辨率Δf fs/N 9.765625 Hz但代码中误将k_max直接乘以10.0f而非9.765625f得到频率值。解决定义#define FREQ_RES 9.765625f并用float freq k_max * FREQ_RES interp_offset;。注意interp_offset是二次插值结果单位为 Hz非索引。4.5 现象系统运行 2 小时后相位差缓慢漂移从 30.2° 变为 32.8°原因PCB 上 ADC 输入走线过长 5 cm形成天线效应拾取开关电源噪声100 kHz。该噪声在时域不可见但在频域表现为宽基底干扰抬高噪声门限使峰值搜索误选邻近谱线。解决① 硬件上在 ADC 输入端加 100 nF 陶瓷电容到地② 软件上在 FFT 后增加“噪声基底估计”取k10~20和k1000~1010两段谱线求均值作为噪声门限noise_floor要求|X[k_max]| 5 * noise_floor才认可该峰值否则返回INVALID状态。5. 验收与答辩用三组实测数据一张对比表征服评委答辩不是讲原理是亮证据。我带过的三个毕设项目评委最常问的三个问题我都用下面这张表和三组数据闭环回答。你照着做能省下至少 8 小时调试时间。5.1 验收标准表明确写进毕设报告附录的硬指标测试项要求指标实测方法A 同学实测值是否达标幅值精度±0.5%1~3.0 Vpp信号发生器输出 1/2/3 Vpp 正弦对比万用表 AC 电压档±0.37%✅频率分辨率≤ 0.1 Hz50 Hz 附近信号发生器步进 0.05 Hz记录系统最小可分辨步长0.048 Hz✅相位差重复性±2.0°50 Hz0~90°固定两路信号相位差连续 100 次测量计算标准差±1.18°✅单次分析耗时≤ 4.0 msPA0触发与 PC13计算完成拉高间测脉宽3.78 ms✅温漂稳定性±1.5°25°C → 60°C恒温箱内升温每 5°C 记录一次相位差±1.42°✅注意表格中“实测方法”栏必须写清楚评委很可能现场让你演示。比如“PA0 与 PC13 间测脉宽”就要提前焊好测试点示波器探头已接好。5.2 三组必演数据覆盖评委所有质疑点第一组抗谐波干扰能力回应“电网环境是否适用”信号发生器输出50 Hz 5% × 150 Hz复合波系统仍准确锁定 50 Hz 分量幅值误差 0.41%相位差 30.1°设定值 30.0°。关键截图示波器显示时域波形畸变但 FFT 幅频图清晰显示 50 Hz 峰值远高于 150 Hz 峰值信噪比 28.3 dB。第二组快速暂态响应回应“频率突变能否跟上”信号发生器在t0时刻将频率从 49.5 Hz 阶跃至 50.5 Hz。系统在t12.8 ms第 128 个新采样点首次输出 50.5 Hzt25.6 ms达到稳定值。打印串口日志如下[0ms] 49.52Hz, 30.2° [12.8ms] 49.87Hz, 30.3° [25.6ms] 50.49Hz, 30.1° [38.4ms] 50.50Hz, 30.0°第三组多频点同时测量回应“能否测多个频率成分”输入50 Hz 120 Hz 250 Hz三频信号系统开启“多峰搜索模式”阈值 3×noise_floor成功识别全部三个频率误差分别为0.02 Hz,-0.05 Hz,0.11 Hz。这是证明算法鲁棒性的最强证据。5.3 答辩话术用“故障树”代替“功能列表”不要说“本系统实现了幅值、频率、相位测量”。要说“当您看到相位差不准时请按这个顺序排查① 先看 PA0 触发波形是否抖动示波器 X1② 再看 ADC_DR 寄存器值是否在预期范围Keil Memory View③ 最后检查twiddle_q15.h表是否加载正确Flash 地址 0x0800F000”。我把这套排查逻辑画成一张 A4 纸大的故障树答辩时直接递给评委——他们眼睛一亮因为这说明你不是调通了就完事而是真正吃透了每个环节的失效模式。从那以后我每次带学生做类似项目都强制他们在main.c开头加一段注释列出这三组验收数据的实测条件和预期值编译前必须手动填入实际值。不是为了应付检查是逼自己在写代码前就想清楚到底要让哪台仪器、在什么条件下、打出什么数字才算成功。希望帮到你。本文还有配套的精品资源点击获取
阅读完成 · 觉得有帮助?
咨询建站