企业官网建设流程全解析

简介：Savitzky-Golay滤波器实现包面向STM32单片机平台，用于嵌入式系统在数据采集中对信号进行平滑除噪。该算法由Savitzky与Golay于1964年提出，核心为局域多项式最小二乘拟合，能够在抑制噪声的同时保持信号形状与宽度不变，适用于传感器数据、仪表信号等实时处理场景。资源压缩包共4个文件，包含2个C源文件与2个头文件，整体大小仅4KB，C文件完成滤波主流程与矩阵运算，头文件提供接口声明，结构精简，便于快速集成到现有工程。已有1996人学习/下载，适合中高级嵌入式开发者直接使用或作为算法移植参考。代码设计上将矩阵运算与SG滤波器主体拆分为独立模块，读者既能获得一个可直接调用的滤波函数，也能通过阅读实现理解窗宽、多项式阶数等关键参数对滤波效果的影响，便于依据具体项目需求进行二次开发与调优。

1. 为什么在STM32单片机信号链中值得用Savitzky-Golay滤波

STM32单片机的ADC采样数据里做滤波，我先试过滑动平均，也试过一阶低通，最后换到Savitzky-Golay滤波才把波形保住了。滑动平均窗口拉长一点，阶跃沿就被抹成斜坡；一阶低通则会在波形上叠加可观的相位延迟。Savitzky-Golay滤波器最初由Savitzky和Golay于1964年发表在Analytical Chemistry杂志，它基于局域多项式最小二乘拟合，在时域做滑动窗内的数据平滑，最大特点是滤除噪声的同时保持信号的形状、宽度不变。这个项目用纯C语言实现，包含SavitzkyGolayFilter.c、SavitzkyGolayFilter.h和matrix.c、matrix.h，不依赖第三方库，很适合STM32这类资源有限的单片机。接下来我按数学原理、参数配置、ADC实战和验证技巧逐步拆开，每个部分都是我实际编译调试时最在意的点。

2. 局域多项式最小二乘：Savitzky-Golay的数学原理与代码映射

2.1 卷积思想与最小二乘解

Savitzky-Golay滤波器的本质并不复杂。窗口长度为奇数2m+1，在某个输出点n处，取x[n-m]到x[n+m]共2m+1个采样点，用一个K阶多项式去拟合这段局部数据，然后取拟合多项式在窗口中心点n处的值作为y[n]。窗口逐点平移，每个输出点都重复这一过程。如果直接对每个点做最小二乘拟合，计算量偏高；但Savitzky和Golay发现这个拟合过程可以离线算成固定的卷积核，在线滤波就退化为一次乘加运算：

y[n] = Σ_{k=-m}^{m} h[k] · x[n-k]

h[k]就是S-G卷积系数。为了推导h[k]，先定义窗口内相对索引i从-m到m，构造基矩阵A，其中第i行第j列是i的j次方，j从0到K。多项式拟合系数c由最小二乘条件‖A·c - x_window‖²最小确定，正规方程为(AᵀA)c = Aᵀx_window。若输入x_window取单位脉冲向量e，即窗口中心为1、其余为0，则解出的c就是卷积核的中心权重。对每个延迟位置重复该操作，就能得到完整的h数组。

需要说明的是，正规方程里的AᵀA是范德蒙德矩阵的Gram矩阵，当阶数K升高时，矩阵条件数会快速恶化。在STM32的单精度浮点环境下，这个问题会被放大，所以在这个库里用matrix.c做矩阵求解时，要特别关注返回值，一旦出现秩亏缺，就得降低阶数或改用double做初始化阶段的系数计算。

2.2 资源包里的matrix.c和SavitzkyGolayFilter模块结构

这个包的设计思路比较典型：matrix.c负责矩阵创建、乘法、转置、求逆等基础操作，SavitzkyGolayFilter.c负责把窗口参数和矩阵结果封装成滤波器对外接口。我在移植时最关心的是对外函数怎么调用，读懂头文件后基本就清楚了。常用接口可以整理成下表：

调用方式很直接，先用Init指定窗口和阶数，再对数组执行Process：

#include "SavitzkyGolayFilter.h" SGFilter_Typedef filter; float input[256], output[256]; /* 7点窗口，3阶多项式，适合传感器平滑 */ SGFilter_Init(&filter, 7, 3); SGFilter_Process(&filter, input, output, 256);

参数说明：window必须是正奇数，order必须小于window-1，否则内部矩阵奇异，Init会返回错误。实际使用时，7点3阶是综合性能较好的起点。如果order设为0，就退化为滑动平均；如果window设为1，则滤波器直通，输出等于输入，这两个边界情况在调试时容易让人迷惑。

2.3 求解卷积核时的C代码逻辑

matrix.c的价值在于把“运行时计算卷积核”这件事变成可能。常见做法是构造基矩阵、求正规方程、解线性方程组，最后和单位脉冲向量相乘得出h。下面是一个典型的系数计算片段，与这个包的实现思路一致：

float* SGFilter_ComputeKernel(int window, int order) { int half = window / 2; Matrix A = Matrix_Create(window, order + 1); Matrix AtA = Matrix_Create(order + 1, order + 1); Matrix e0 = Matrix_Create(window, 1); Matrix coeff = Matrix_Create(order + 1, 1); float *kernel = (float*)malloc(window * sizeof(float)); /* 构造基矩阵：A[i][j] = (i - half)^j */ for (int i = 0; i < window; i++) { double x = i - half; A.data[i * A.cols + 0] = 1.0; for (int j = 1; j <= order; j++) { A.data[i * A.cols + j] = A.data[i * A.cols + j - 1] * x; } } /* AtA = A^T * A */ Matrix_MultiplyTranspose(&A, &AtA); /* 单位脉冲向量，只让窗口中心为1 */ e0.data[half * e0.cols + 0] = 1.0; /* 解正规方程 (A^T A) coeff = A^T e0 */ Matrix_Solve(&AtA, &A, &e0, &coeff); /* 卷积核 kernel[i] = A[i] * coeff */ Matrix_Multiply(&A, &coeff, kernel); return kernel; }

代码逻辑说明：先根据窗口中心建立多项式基矩阵，矩阵乘法函数求出AᵀA，随后用高斯消元解线性方程组。Matrix_Solve内部如果检测到主元为0，会返回错误码，说明窗口和阶数配置不合理。得到kernel后，在线滤波就是一组对称系数的乘加。如果窗口和阶数在编译期固定，我一般会把Init放在上电后执行一次，然后把系数打印出来固化成常量表，省去后续的矩阵运算。

3. 可移植到STM32的S-G滤波器：参数配置与内存估算

3.1 窗口、阶数与信号特征的匹配关系

参数选择直接决定滤波行为。窗口半宽m决定局部拟合范围，窗口越长，平滑越强，但对高频成分的保留越差。多项式阶数K决定拟合曲线追踪信号局部变化的能力，K越高，细节越完整，但也更容易把噪声拟合进去，同时数值条件变差。理论上必须满足K < window-1，工程上我一般限制K不超过窗口长度的一半。

下面这张表给出了我在项目中验证过的常见组合：

滑动平均本质上是0阶S-G，7点0阶就是普通7点平均。从滑动平均切到S-G后，最直观的区别是阶跃沿保住了，峰值不会被压低。如果信号本身有线性斜坡特征，用1阶拟合会比2阶更平顺，但响应速度稍慢。实际调参时，建议先用离线数据在PC上跑一遍不同组合，再定STM32上的最终值。

3.2 在STM32上的RAM与Flash占用评估

一个30点窗口的S-G滤波器，卷积核只需要30个float，120字节。真正占用RAM多的不是核，而是数据缓冲和矩阵运算临时数组。如果直接调用matrix.c的动态矩阵接口，要注意单片机的堆大小。在Cortex-M3/M4上我一般用预分配静态缓冲，避免malloc碎片：

#define SG_MAX_WINDOW 31 #define SG_MAX_ORDER 5 static float matrix_buffer[(SG_MAX_ORDER + 1) * (SG_MAX_ORDER + 1)]; static float ring_buffer[SG_MAX_WINDOW]; static SGFilter_Typedef sg; SGFilter_Init(&sg, 7, 3);

这段配置把窗口上限设到31，阶数上限设到5，matrix_buffer只需要36字节，ring_buffer约124字节，总共不到200字节RAM。STM32F103C8T6的20KB SRAM可以同时跑多组滤波器而毫无压力。需要注意Keil工程的浮点对齐，结构体内部指针编译器会自动处理，不需要手工字节对齐。

3.3 实时数据流处理与边界策略

S-G滤波器是非因果系统，输出y[n]依赖未来的x[n+m]。实时处理时，要么接受m个采样周期的输出延迟，要么对窗口边界做延拓。最常见的做法是维护一个长度为窗口的滑窗缓冲区：

static float data_buf[SG_MAX_WINDOW]; static float* pdata = data_buf; float SGFilter_NextSample(float newest) { /* 向前移动历史数据，把新值放到队尾 */ memmove(data_buf, &data_buf[1], (SG_MAX_WINDOW - 1) * sizeof(float)); data_buf[SG_MAX_WINDOW - 1] = newest; /* 只有窗口完全填充后才输出有效值 */ if (samples_collected < SG_MAX_WINDOW) { samples_collected++; return newest; /* 启动阶段直接旁路 */ } SGFilter_Process(&sg, data_buf, &output, SG_MAX_WINDOW); return output[SG_MAX_WINDOW / 2]; }

逻辑说明：每次新数据到达时，用memmove平移窗口，最新样本填到尾部，随后对整个窗口运行Process，取中心位置作为输出。由于Process是逐个点连续处理的，边界会按照内部策略处理，因此这里取center位置是安全的。如果要求输出与采样完全同步，启动阶段可以把第一个采样值重复填充整个窗口再开始滤波，避免前m个点产生暂态偏置。对大多数传感器信号，这个暂态可以忽略。

4. STM32 ADC采样数据平滑实战：从调用到调参

4.1 以STM32F407采集压力传感器为例

我在一套基于STM32F407的采集板上验证了这个库，信号源是桥式压力传感器经过仪表放大器后的电压，片内ADC配置为12位，采样率2kHz。原始数据里既有电源纹波，也有反激开关带来的窄脉冲干扰。采集链路用DMA双缓冲，每满256个样本进入中断，滤波放到主循环执行：

#include "SavitzkyGolayFilter.h" #define ADC_BUF_SIZE 512 static volatile uint16_t adc_buffer[ADC_BUF_SIZE]; static float input[ADC_BUF_SIZE], output[ADC_BUF_SIZE]; SGFilter_Typedef sg; int main(void) { MX_GPIO_Init(); MX_DMA_Init(); MX_ADC1_Init(); SGFilter_Init(&sg, 7, 3); HAL_ADC_Start_DMA(&hadc1, adc_buffer, ADC_BUF_SIZE); while (1) { if (dma_complete_flag) { dma_complete_flag = 0; for (int i = 0; i < ADC_BUF_SIZE; i++) { input[i] = adc_buffer[i] * 3.3f / 4096.0f; } SGFilter_Process(&sg, input, output, ADC_BUF_SIZE); } } }

逻辑说明：ADC原始值先换算成电压，再送滤波器。SGFilter_Process会对整块数据做连续处理，输出数组长度与输入一致。DMA回调里只置标志位，不直接做浮点转换，避免中断占用时间过长。如果你用定时器触发ADC，时序上要保证滤波频率与数据到达频率一致，不能让主循环里的Process覆盖新数据。

4.2 从采集波形判断参数是否合适

滤波好不好，不能只看标准差降了多少，还要看上升沿和过冲。我一般用原始信号、5点2阶、7点3阶、15点3阶四组数据对比，记录噪声标准差、阶跃延迟和过冲量。以下是一组实测结果：

从表里能看出，窗口越大噪声压得越低，但延迟和过冲都会变大。7点3阶在这个压力场景下是平衡点，噪声降到原来的三分之一，阶跃延迟只有3.5ms，过冲可接受。如果闭环控制带宽高，比如电机电流环或呼吸机压力环，建议不要超过11点；如果只是给上位机显示，15点3阶的表现更好。

4.3 在线调整系数的切换策略

这个库支持运行中改变窗口和阶数，Init会重新计算卷积核。切换时最需要注意的是暂态跳变：旧窗口数据还在，新卷积核直接作用会产生输出突变。我通常先调用SGFilter_Reset清空内部状态，再调用SGFilter_Init设置新参数，同时标记一个“静默期”，等新窗口填满后再启用输出。对于2kHz采样率，7点窗口静默期只有3.5ms，用户几乎感觉不到。

在线调整参数可以通过串口下发，比如用自定义协议帧包含window和order字段。接收到合法配置后，在空闲时刻执行初始化，这样不会干扰正在进行的ADC采集。需要注意的是，如果上位机在下发参数时把window设成偶数，Init会返回错误，必须在校验里拒绝。我在调Modbus类帧接收时也遇到过类似字段边界问题，解决办法是先解析到临时变量，确认合法后再写入运行结构体。

4.4 容易被忽略的坑

第一个坑是矩阵求逆失败。当order接近window-1时，Matrix_Solve返回值可能异常。Init阶段最好检查返回码，并打印行列式或条件数。第二个坑是STM32的FPU未启用。Cortex-M4F如果没有在启动代码里设置CPACR寄存器，浮点运算会全部走软件库，初始化一次Matrix_Solve可能要占几毫秒。第三个坑是边界效应。S-G在序列开头和结尾的处理策略直接影响波形，如果库没有延拓，开头几个点是无效的，需要丢弃。第四个坑是参数配置与采样率不匹配：窗口不变时，采样率越高，实际平滑时间越短，所以评估延迟要结合采样率一起看。

5. 用单位脉冲响应校验STM32上的S-G滤波系数精度

5.1 在主机端用Python生成理论核

矩阵求逆的精度难以靠肉眼观察来判断，最可靠的办法是让STM32处理一个单位脉冲序列，再把输出与离线计算的理论核对比。窗口7、阶数3的系数可以这样算出来：

import numpy as np window = 7 order = 3 half = window // 2 x = np.arange(window) - half A = np.vander(x, order + 1, increasing=True) e0 = np.zeros(window) e0[half] = 1.0 kernel = np.array([np.dot(A[i], np.linalg.lstsq(A, e0, rcond=None)[0]) for i in range(window)]) print(kernel)

这段代码求解的是窗口中心对单位脉冲的最小二乘响应，得到的kernel应当与STM32初始化后的系数一致。把STM32算出的kernel通过串口打印，逐项比较，误差小于1e-4说明矩阵运算正常。如果误差到1e-3级别，问题大多出在matrix.c用了float直接求逆，可以在初始化阶段改double，算完再转成float保存。

5.2 在Keil工程里嵌入自检用例

更实用的做法是把单位脉冲自检放到上电流程里，只在DEBUG模式下编译：

#ifdef SG_SELF_TEST static float test_in[7] = {0, 0, 0, 1.0f, 0, 0, 0}; static float test_out[7]; SGFilter_Process(&sg, test_in, test_out, 7); for (int i = 0; i < 7; i++) { /* 与预期kernel逐项比较，允许float误差1e-4 */ if (fabsf(test_out[i] - expected_kernel[i]) > 1e-4f) { /* 上报错误代码 */ } } #endif

代码逻辑：单位脉冲输入经过S-G滤波后，输出就是滤波器的脉冲响应，也就是卷积核本身。比较时要注意边界，因为数据序列长度等于窗口长度，输出前几个点受边界延拓影响，应该只拿窗口中心对齐后的有效k个点做比较。忽略这一点容易误报。

5.3 用DWT计数器测量耗时并做系数对称优化

在STM32上测滤波耗时，最简单的方法是用Cortex-M内核的DWT_CYCCNT。初始化DWT后，在Process前后读取计数差，就能得到精确到时钟周期的耗时。7点3阶S-G卷积每样本只做7次乘加，在168MHz的F407上大约几十个周期；如果实测几百周期，先查Keil优化是否在-O0，再把优化等级调到-O2。

更进一步，S-G卷积核是对称的，可以把首尾采样值先相加再乘同一系数。7点3阶的核是[0.2145, 0.0952, -0.1429, 0.3333, -0.1429, 0.0952, 0.2145]，利用对称性把乘法次数从7次降到4次，同时减少一半的内存读取。在滑窗实现里，用两个索引分别从首尾取数相加就行，代码改动不大。发布前建议把单位脉冲自检留在固件里，万一时钟或DMA配置异常导致系数偏差，自检能第一时间给出错误码，方便现场定位。

本文还有配套的精品资源，点击获取