简介：面向嵌入式信号处理初学者与需要做频域测量的工程师，这份资料围绕STM32F407平台，用标准库实现FFT并完成正弦波幅值、频率与相位差的测量。内容涉及ADC采样、窗函数处理、FFT变换、频谱分析与结果输出，适合电机控制、音频分析和无线通信等实时信号处理场景的学习与验证。压缩包共283个文件，约9.98MB，包含C/H源文件、Keil工程配置、编译中间文件及可执行输出等，便于直接打开工程、对照源码理解算法流程与工程结构。已有1668人学习下载，可从中获得一套较完整的STM32 FFT测量实现参考，并据此梳理采样、频谱计算和结果输出之间的衔接关系。
1. 一块最小系统板加两路正弦波输入，就能得到这样的读数
一块 STM32 最小系统板，两路同频正弦波分别从 PA0、PA1 输入，串口每隔 100 ms 输出一行结果：A=1.002V F=1000.31Hz dPhi=44.97deg——这是我把 FFT 真正跑在单片机上之后一直关注的东西。很多人对快速傅里叶变换的印象还停留在上位机里那条漂亮的频谱曲线，可一旦搬到 STM32 上，难点根本不在蝶形公式本身，而在采样率怎么定、直流偏置怎么去掉、频率没有正好落在某个频率间隔上时幅值为什么忽大忽小、两个通道不是同一时刻采样时相位差为什么会偏。这套基于标准外设库的 FFT 测参方案做的就是整条链路：定时器触发 ADC 采 1024 点，一次基 2 变换，把频域里的峰值还原成幅值、频率和双通道相位差。它适合手上已经有 F103 或 F407 板子、想把示波器那点测量能力加入自己产品的开发者，也适合刚啃完 DSP 教材、想找个能上板验证的练手项目的人。
2. 信号链与采样前端：把连续正弦波变成 FFT 能处理的数组
采样前端这一层最容易被轻视。算法写得再漂亮，进来的数组本身就是乱的，后面全是白费功夫。这一章先把参数算清楚，再把定时器、ADC、DMA 这套配置搭出来，最后做码值到浮点数组的预处理。
2.1 采样率、点数和分辨率三者互相制约
给定采样率fs和点数N，频率分辨率Δf = fs / N，一帧的时间长度T = N / fs，能无混叠测量的最高频率是fs / 2。测 1 kHz 正弦，取fs = 20 kHz、N = 1024，得到Δf ≈ 19.53 Hz、T = 51.2 ms。这 51.2 ms 意味着一帧内 1 kHz 信号转了 51 个周期，这一点很关键：周期数越多，主瓣越窄，幅值受频谱泄漏影响越小。反过来只采 5 个周期，就算加了窗，幅值读数也会随信号频率抖动，这类精度问题很容易被忽略。
| 参数 | 典型取值 | 决定什么 |
|---|---|---|
FFT_N | 1024 | 频率分辨率与运算量，1024 点基 2 是 5120 次蝶形 |
FS_HZ | 20000 | 可测上限 10 kHz，触发间隔 50 µs |
Δf | 19.53 Hz | 频率的量化单位 |
T | 51.2 ms | 数据刷新率的上限 |
| ADC 时钟 | 12 MHz | F103 上不得超过 14 MHz |
采样率不能只满足奈奎斯特，实际工程里我一般会留出 5～10 倍余量，理由是两个：一是抗混叠滤波器不可能做成理想的矩形，滚降需要过渡带；二是每周期点数太少时，ADC 的孔径抖动和触发抖动会被放大成相位噪声。20 kHz 采 1 kHz 得到每周期 20 点，配上 1024 点一帧，算是精度和运算量之间比较舒服的一档。
抗混叠那一级千万别省。PA0 前面串一个 1 kΩ 电阻、对地挂 33 nF，截止频率约 4.8 kHz，正好压在信号带宽之上、fs/2之下，把 10 kHz 以上会折叠回来的成分先衰减掉。没有这一级，示波器探头夹上去看到的“干干净净的正弦”，在频谱上照样能冒出一堆莫名其妙的谱线。
注意：ADC 采样时间要跟信号源内阻匹配。41.5 个采样周期在 12 MHz 下是 3.46 µs，源阻抗几 kΩ 以内都够用；如果前面加了很大的分压电阻，采样电容充不满，幅值会系统性偏小，而且偏得很“稳定”，特别容易误判成算法问题。
2.2 定时器触发 + ADC 扫描 + DMA 搬运的标准库写法
采样必须由硬件定时器触发，绝不能在主循环里手动启动转换，否则采样间隔抖动直接变成相位误差。标准外设库这套初始化比较长，但每一步都有明确目的，照着改宏定义就能用。
#define FFT_N 1024u #define FS_HZ 20000.0f #define VREF 3.3f #define ADC_FULL 4095.0f uint16_t adc_raw[FFT_N * 2]; /* 两通道交替存放：CH0,CH1,CH0,CH1... */ volatile uint8_t frame_ready = 0; void sample_init(void) { GPIO_InitTypeDef gpio; ADC_InitTypeDef adc; DMA_InitTypeDef dma; TIM_TimeBaseInitTypeDef tim; RCC_APB2PeriphClockCmd(RCC_APB2Periph_GPIOA | RCC_APB2Periph_ADC1, ENABLE); RCC_AHBPeriphClockCmd(RCC_AHBPeriph_DMA1, ENABLE); RCC_APB1PeriphClockCmd(RCC_APB1Periph_TIM2, ENABLE); RCC_ADCCLKConfig(RCC_PCLK2_Div6); /* 72MHz/6 = 12MHz，别超 14MHz */ gpio.GPIO_Pin = GPIO_Pin_0 | GPIO_Pin_1; gpio.GPIO_Mode = GPIO_Mode_AIN; GPIO_Init(GPIOA, &gpio); ADC_DeInit(ADC1); adc.ADC_Mode = ADC_Mode_Independent; adc.ADC_ScanConvMode = ENABLE; /* 两通道轮流转换 */ adc.ADC_ContinuousConvMode = DISABLE; /* 由 TIM2 触发，不能连续 */ adc.ADC_ExternalTrigConv = ADC_ExternalTrigConv_T2_TRGO; adc.ADC_DataAlign = ADC_DataAlign_Right; adc.ADC_NbrOfChannel = 2; ADC_Init(ADC1, &adc); ADC_RegularChannelConfig(ADC1, ADC_Channel_0, 1, ADC_SampleTime_41Cycles5); ADC_RegularChannelConfig(ADC1, ADC_Channel_1, 2, ADC_SampleTime_41Cycles5); ADC_DMACmd(ADC1, ENABLE); ADC_ExternalTrigConvCmd(ADC1, ENABLE); DMA_DeInit(DMA1_Channel1); dma.DMA_PeripheralBaseAddr = (uint32_t)&ADC1->DR; dma.DMA_MemoryBaseAddr = (uint32_t)adc_raw; dma.DMA_DIR = DMA_DIR_PeripheralSRC; dma.DMA_BufferSize = FFT_N * 2; /* 半字为单位 = 2*N */ dma.DMA_PeripheralInc = DMA_PeripheralInc_Disable; dma.DMA_MemoryInc = DMA_MemoryInc_Enable; dma.DMA_PeripheralDataSize = DMA_PeripheralDataSize_HalfWord; dma.DMA_MemoryDataSize = DMA_MemoryDataSize_HalfWord; dma.DMA_Mode = DMA_Mode_Circular; dma.DMA_Priority = DMA_Priority_High; DMA_Init(DMA1_Channel1, &dma); DMA_ITConfig(DMA1_Channel1, DMA_IT_TC, ENABLE); DMA_Cmd(DMA1_Channel1, ENABLE); /* 72MHz / 20000Hz = 3600 个计数，向上计数，更新事件作 TRGO */ tim.TIM_Period = 3600 - 1; tim.TIM_Prescaler = 0; tim.TIM_ClockDivision = TIM_CKD_DIV1; tim.TIM_CounterMode = TIM_CounterMode_Up; TIM_TimeBaseInit(TIM2, &tim); TIM_SelectOutputTrigger(TIM2, TIM_TRGO_Update); TIM_Cmd(TIM2, ENABLE); ADC_Cmd(ADC1, ENABLE); ADC_ResetCalibration(ADC1); while (ADC_GetResetCalibrationStatus(ADC1)); ADC_StartCalibration(ADC1); while (ADC_GetCalibrationStatus(ADC1)); ADC_SoftwareStartConvCmd(ADC1, ENABLE); /* 首次启动，后续靠 TRGO */ }几个参数值得单独说。TIM_Period = 3600 - 1是 72 MHz 除以 20 kHz 的结果，改采样率就改这个数，改完记得同步FS_HZ宏，两处不一致时频率读数会整体缩放，这种错误很典型。ADC_ExternalTrigConv_T2_TRGO配上TIM_SelectOutputTrigger(TIM2, TIM_TRGO_Update)才是完整的一对，只写一半的话 ADC 会一直等不到触发，现象是 DMA 一个数都不搬。DMA_BufferSize的单位是传输元素个数而不是字节，配置成半字的时候就是2 * FFT_N，写成4 * FFT_N会多搬一倍，数组越界却不报错。
还有一点：扫描模式下 CH0 和 CH1 不是同时采的，CH1 比 CH0 晚了“采样时间 + 转换时间”这么久——41.5 + 12.5 = 54 个 ADC 周期，在 12 MHz 下正好 4.5 µs。这个固定偏斜后面算相位差时必须补回来，第 4 章会用到。
2.3 码值到浮点数组：去直流比标定更要紧
转换这一步做两件事：把 12 位码值映射成电压，把直流分量减掉。第二件事如果忘了，0 号频点的能量会大到把自动量程、峰值搜索全部带跑偏。
static void raw_to_float(const uint16_t *raw, float *x0, float *x1) { float sum0 = 0.0f, sum1 = 0.0f; for (uint16_t i = 0; i < FFT_N; i++) { x0[i] = raw[2 * i] * (VREF / ADC_FULL); /* CH0 */ x1[i] = raw[2 * i + 1] * (VREF / ADC_FULL); /* CH1 */ sum0 += x0[i]; sum1 += x1[i]; } float dc0 = sum0 / (float)FFT_N; /* 均值即直流分量 */ float dc1 = sum1 / (float)FFT_N; for (uint16_t i = 0; i < FFT_N; i++) { x0[i] -= dc0; x1[i] -= dc1; } }用均值当直流分量是最省事的做法，前提是一帧里包含整数个信号周期，否则均值本身会带一点残余。测 1 kHz、采 51.2 ms 这种整周期关系时误差可以忽略；信号频率未知或者不是整周期，就该换成硬件上的中点偏置——两路信号先经过 1.65 V 的偏置网络再进 ADC，把 0 V 抬到量程中央，然后减固定常数VREF/2，比减均值稳定得多。
硬件上还得多一句：STM32 的 ADC 输入不能承受负电压，正弦波必须抬到 0～3.3 V 之间。信号源直接输出 ±1 V 就进管脚，轻则波形下半周削平、频谱上出现一堆奇次谐波，重则把管脚打坏。我一般用一只 10 kΩ 电阻把信号拉到 1.65 V 分压点，串 100 nF 隔直，简单可靠。
3. FFT 核心：在 STM32 上跑 1024 点基 2 快速傅里叶变换
前端交出两个干净的浮点数组，接下来是纯算法部分。1024 点基 2 是 10 级、每级 512 个蝶形，一共 5120 次蝶形运算。手写一遍的价值在于你能确切知道每个下标在干什么，出了问题也能定位到具体哪一级。
3.1 位反转置换与旋转因子表
基 2 的时域抽取（DIT）要求输入按位反转顺序排列，输出才是自然顺序。位反转这段代码短，但下标逻辑很容易写错，写成死循环或者只对一部分元素生效都不报错，只是频谱结果完全不对。
static float fft_re[FFT_N], fft_im[FFT_N]; static float tw_re[FFT_N / 2], tw_im[FFT_N / 2]; static void fft_init_twiddle(void) { for (uint16_t i = 0; i < FFT_N / 2; i++) { float ang = -2.0f * 3.14159265358979f * (float)i / (float)FFT_N; tw_re[i] = cosf(ang); tw_im[i] = sinf(ang); } } static void bit_reverse(float *re, float *im, uint16_t n) { uint16_t j = 0; for (uint16_t i = 1; i < n; i++) { uint16_t bit = n >> 1; while (j & bit) { /* 找到 j 的下一个反转值 */ j ^= bit; bit >>= 1; } j ^= bit; if (i < j) { /* 只在前者小于后者时交换，保证每对换一次 */ float t = re[i]; re[i] = re[j]; re[j] = t; t = im[i]; im[i] = im[j]; im[j] = t; } } }if (i < j)这个判断是必须的，去掉它会把每对元素交换两次，等于什么都没做，而数组里只有实部非零、虚部全 0 时，错误又不会崩溃，属于典型的“静默错误”。旋转因子我放在const大数组里烧到 Flash 更省 RAM，1024 点需要 512 个float共 2 KB，F103C8 那点 20 KB RAM 里塞得下但也不宽裕，可以按需只存第一象限再用对称性推出其余三个象限。
3.2 三层循环的蝶形运算
蝶形那段是整个算法最容易写错下标的地方。以 DIT 为例，外层按级数增长块长度len，中间层遍历所有块，内层在块内做一半的蝶形。
void fft_radix2(void) { bit_reverse(fft_re, fft_im, FFT_N); for (uint16_t len = 2; len <= FFT_N; len <<= 1) { uint16_t half = len >> 1; uint16_t step = FFT_N / len; /* 旋转因子表的下标步长 */ for (uint16_t i = 0; i < FFT_N; i += len) { for (uint16_t j = 0; j < half; j++) { uint16_t k = i + j; /* 上半个蝶形 */ uint16_t m = k + half; /* 下半个蝶形 */ uint16_t ti = j * step; /* 取表下标，不是 j */ float wr = tw_re[ti], wi = tw_im[ti]; float mr = fft_re[m], mi = fft_im[m]; float tr = mr * wr - mi * wi; /* 复数乘 W */ float tv = mr * wi + mi * wr; fft_re[m] = fft_re[k] - tr; fft_im[m] = fft_im[k] - tv; fft_re[k] += tr; fft_im[k] += tv; } } } }ti = j * step是最关键的一行。每一级的旋转因子只用到表里的一部分，第一级（len = 2）只用W[0]，最后一级（len = N）要用完整 512 个，step就是用来按级抽取的。写成ti = j的话，第一级就会取到W[1]、W[2]，结果看着像频谱、频率也对得上，但幅值和相位全都慢慢偏掉——这类错误最难排查，因为波形“看起来是对的”。
先做减法再赋值也需要注意顺序：fft_re[k]被更新之前必须先算出tr，上面代码已经先存了中间变量，如果把fft_re[k] += tr写在前面，fft_re[m]那一行用的就是脏数据。
RAM 方面，三个 1024 点浮点数组（实部、虚部，加上输入的两路）就是 12 KB 以上，F103C8 的 20 KB 一下就去掉大半。真要在小容量芯片上跑，输入可以直接原地变换，两路信号分时处理，反正 FFT 是线性变换，两路之间没有耦合。
3.3 浮点、Q15 与 CMSIS-DSP 的三条路
同等正确性下，摆在面前的选择有三条。第一条是上面这种手写单精度浮点：代码透明，改参数方便，缺点是没有 FPU 的 Cortex-M3 上全靠软件浮点库，指令数大约是整数的三到五倍。我在 72 MHz 的 F103 上跑 1024 点单精度，量级落在几十毫秒，做 10 Hz 级的慢速监测完全够，做实时音频就吃力了。带 FPU 的 F407 是另一个数量级，1 ms 以内能出结果，代价是芯片贵一截、功耗高一点。
第二条是 Q15 定点。所有数据缩放到 ±1.0 对应的 ±32767，蝶形每级的增益是 2 倍，所以每做一级就要右移一位防止溢出，1024 点做完 10 级，整体相当于右移了 10 位，小信号的信噪比掉得厉害。想要精度就得用块浮点（每级统一判断是否溢出再决定移几位），复杂度立刻上一个台阶。
第三条是直接用 CMSIS-DSP，这也是我最常用的方案。它是 ARM 的官方 DSP 库，和标准外设库并不冲突，Keil 工程里把arm_cortexM3l_math.lib（M4F 用arm_cortexM4lf_math.lib）加进来就行。
#include "arm_math.h" static arm_rfft_fast_instance_f32 rfft; static float fft_in[FFT_N], fft_out[FFT_N]; void dsp_fft_prepare(void) { arm_rfft_fast_init_f32(&rfft, FFT_N); /* 只需初始化一次 */ } /* 输入实数序列，输出打包格式的频谱 */ void dsp_fft_run(const float *x, float *mag, uint16_t half_n) { memcpy(fft_in, x, sizeof(float) * FFT_N); arm_rfft_fast_f32(&rfft, fft_in, fft_out, 0); /* 0 = 正变换 */ /* 打包格式：out[0]=DC，out[1]=Nyquist， 之后 out[2k]=Re[k]，out[2k+1]=Im[k]，k 从 1 到 N/2-1 */ mag[0] = fabsf(fft_out[0]); mag[half_n] = fabsf(fft_out[1]); for (uint16_t k = 1; k < half_n; k++) { float re = fft_out[2 * k]; float im = fft_out[2 * k + 1]; mag[k] = sqrtf(re * re + im * im); } }arm_rfft_fast_f32用的是实数 FFT，比起把虚部填零再跑复数 FFT，运算量直接减半，这个便宜不占白不占。它输出的打包格式是刚接触时最容易踩坑的地方：不是[Re0, Im0, Re1, Im1, ...]那种规整排列，直流和奈奎斯特两个点被单独拎出来塞在开头。照着这个顺序去取sqrt(re²+im²)，出来的幅值会整体错位，读数完全对不上。
提示：
arm_rfft_fast_f32只支持 32 位单精度。M3 上没有 FPU 的时候，arm_rfft_q15配合arm_cfft_q15会快很多，但要把缩放因子一路管到底，调试期间我一般不这么干。
如果只是教学或者想弄清楚每一步，手写那一版更好；如果项目要交付、要跑得快，直接用 CMSIS-DSP。这两条路我都留着，因为前者用来对照后者，能很快判定问题出在算法还是出在调用。
4. 三参数提取：从频谱到幅值、频率、相位差
频谱算出来只是起点。真正要交给上层的是三个数：幅值多少伏、频率多少赫兹、两路之间相位差多少度。这三件事里，只有频率是自然落在频率间隔上的，另外两个都需要做额外处理。
4.1 峰值搜索与单边谱幅值还原
实数信号的 FFT 是共轭对称的，能量在正负频率各占一半，所以正频率那一侧找到的幅值要乘 2 才是真实幅值。加上窗函数之后还要除掉窗的相干增益，汉宁窗是 0.5。
typedef struct { float amp; /* 幅值，单位 V */ float freq; /* 频率，单位 Hz */ float phase; /* 初相位，单位 rad */ } wave_t; void spectrum_measure(const float *x, wave_t *out) { float win_sum = 0.0f; for (uint16_t i = 0; i < FFT_N; i++) { /* 汉宁窗，N-1 做分母保证首尾对称 */ float w = 0.5f - 0.5f * cosf(2.0f * 3.14159265f * (float)i / (float)(FFT_N - 1)); fft_re[i] = x[i] * w; fft_im[i] = 0.0f; win_sum += w; } fft_radix2(); uint16_t kmax = 1; float mmax = 0.0f; for (uint16_t k = 1; k < FFT_N / 2; k++) { /* 跳过 DC 与奈奎斯特 */ float m = sqrtf(fft_re[k] * fft_re[k] + fft_im[k] * fft_im[k]); if (m > mmax) { mmax = m; kmax = k; } } /* 单边谱 *2，再除以窗相干增益 win_sum/N（汉宁窗约 0.5） */ out->amp = 2.0f * mmax / win_sum; out->freq = (float)kmax * FS_HZ / (float)FFT_N; out->phase = atan2f(fft_im[kmax], fft_re[kmax]); }2.0f * mmax / win_sum这一句把窗增益归一化合并进去了，比硬写mmax / (FFT_N * 0.5f)更通用——换布莱克曼窗或者平顶窗的时候只要改窗函数，不用回头改幅值公式。平顶窗的主瓣很宽，频率分辨能力差，但幅值精度高得多，专门用来做幅值表的场合我会换成它。
峰值搜索的范围同理，从k = 1开始、到N/2 - 1结束。包含 0 号或者N/2那个点的话，只要信号里有一点残余直流，峰值就会被吸到 0 号上去，频率直接读成 0 Hz。
4.2 频率细化：插值法和相位差分法
kmax * fs / N只能分辨到 19.53 Hz 这一档。想再细，得用频谱插值把频率间隔内部的位置估出来。最简单的是对数幅值的抛物线插值：
/* 在峰值周围三点做抛物线拟合，求小数偏移 delta ∈ [-1, 1] */ static float freq_delta(uint16_t k, const float *mag) { float a = mag[k - 1], b = mag[k], c = mag[k + 1]; float denom = a - 2.0f * b + c; if (fabsf(denom) < 1e-12f) return 0.0f; /* 三点等高，无法拟合 */ return 0.5f * (a - c) / denom; } float f_refined = ((float)kmax + freq_delta(kmax, mag)) * FS_HZ / (float)FFT_N;插值法简单、单帧就能算，但精度依赖信噪比，通常能做到Δf的 1/10 到 1/20，也就是 1～2 Hz。要更好的效果就用相位差分法：隔M个采样点取两帧，同一频点k的相位差Δθ和频率偏移δ满足Δθ = 2πδM/N，于是
/* theta1、theta2 是相隔 M 个采样点两帧在同一 bin 上的相位 */ float dtheta = theta2 - theta1; while (dtheta > 3.14159265f) dtheta -= 2.0f * 3.14159265f; while (dtheta < -3.14159265f) dtheta += 2.0f * 3.14159265f; float delta = dtheta * (float)FFT_N / (2.0f * 3.14159265f * (float)M); float f_fine = ((float)kmax + delta) * FS_HZ / (float)FFT_N;代价是要采两帧，刷新率降一半；换来的是在中等信噪比下能把误差压到 0.1 Hz 量级。M越大越灵敏，但2πδM/N超过 π 就会产生相位模糊，M的上限大约是N/(2|δ|max)，实际取N/4比较稳妥。
4.3 双通道相位差：atan2、频率偏移补偿与采样偏斜
两路各自用atan2f(im, re)算出峰值处的相位，相减就是相位差。但这里有两处系统误差必须补。
第一处来自频率没有落在频率间隔正中。设真实频率是(k + δ)·fs/N，那么这个 bin 上测到的相位并不是信号的真实初相，而是被旋转了一个角度。矩形窗和汉宁窗在这一点上关系一致：
/* 把 bin k 上的相位换算回信号真实初相 */ static float phase_correct(float phi_bin, uint16_t k, float delta) { return phi_bin - 3.14159265f * delta * (float)(FFT_N - 1) / (float)FFT_N; }π · δ · (N-1)/N这个量在N = 1024、δ = 0.5时大约是 1.57 rad，将近 90°，不补的话相位差读数会随信号频率漂移，看着像“相位不稳定”，实际是计算漏了一步。
第二处来自两通道不是同时采样。2.2 节算过，扫描模式下 CH1 比 CH0 晚 4.5 µs，这等价于 CH1 的相位被额外超前了 `2πf·
本文还有配套的精品资源，点击获取
