本文介绍Savitzky-Golay滤波器基本原理。

Savitzky-Golay滤波器（简称为S-G滤波器）被广泛地运用于数据平滑去噪，它是一种在时域内基于局域多项式最小二乘法拟合的滤波方法。这种滤波器最大的特点在于在滤除噪声的同时确保信号的形状，宽度不变。S-G滤波器滤波的效果和选取窗口宽度，多项式阶次有关。基于这种滤波器的特点，它在光谱分析（平滑滤波）中经常使用。

1.基本原理

S-G滤波器本质上属于FIR滤波器，我们知道设计FIR滤波器关键在于求得滤波器系数。对于S-G滤波器其系数即为多项式前的常系数，下面对其作简单推导。

设滤波器窗口宽度为L=2k+1，即以中心点s前后k个数据作为滤波器的窗口，多项式阶次为n。

则

进而，

H为范德蒙矩阵，经推导x的预测值

其中

1)B为与输入x无关的矩阵，也就是S-G滤波器系数

2)矩阵的行数由窗口宽度L决定，为L行

3)矩阵的列数由多项式阶次n决定，为n+1列

2.应用

1)滤波器系数求解

这里以窗口宽度为11，阶次为4为例。

在matlab中命令行输入：

order=4;
framelen=11;
v=-5:1:5;
A = fliplr(vander(v));
A=A(1:framelen,1:order+1);
B=A*inv(A'*A)*A';

在matlab中有专门设计S-G滤波器的命令，这里我们对比一下2个值，在命令行输入：

order=4;
framelen=11;
b = sgolay(order,framelen);

经过比较，这2个值是相等的。

注意：

a)矩阵B或b中间行向量（第6行）即为S-G滤波器的滤波器系数

b)矩阵B或b其他行在处理信号边缘（最左侧及最右侧）有用

2)信号边缘的采样值处理

a)采用sgolayfilt对信号进行处理

sgolayfilt是matlab内置的对输入信号进行S-G滤波的命令，这里以此为参考，对比信号边缘是否处理的差异。

order = 4;
framelen = 11;

lx = 36;
x = randn(lx,1);

sgf = sgolayfilt(x,order,framelen);

plot(x,':')
hold on
plot(sgf,'.-')
legend('signal','sgolay')

输出结果：

b)直接使用滤波器系数进行滤波

这里使输入信号和滤波器系数卷积进行滤波。

m = (framelen-1)/2;

B = sgolay(order,framelen);

steady = conv(x,B(m+1,:),'same');

plot(steady)
legend('signal','sgolay','steady')

输出结果：

可见，直接使用滤波器系数进行滤波，在信号边缘滤波效果并不好，和matlab内置的滤波器命令也有差异。

c)对信号边缘进行处理

靠近信号边缘的采样无法放在对称窗的中心，必须区别对待。

为了确定启动瞬变，对 B 的前 (framelen-1)/2 行与信号的前 framelen 个采样执行矩阵乘法。

ybeg = B(1:m,:)*x(1:framelen);

为了确定终止瞬变，对 B 的最后 (framelen-1)/2 行与信号的最后 framelen 个采样执行矩阵乘法。

yend = B(framelen-m+1:framelen,:)*x(lx-framelen+1:lx);

将瞬变部分与稳态部分连接起来以生成完整信号。

cmplt = steady;
cmplt(1:m) = ybeg;
cmplt(lx-m+1:lx) = yend;

plot(cmplt)
legend('signal','sgolay','steady','complete')
hold off

输出结果：

可见，经过信号边缘进行处理后，和sgolayfilt滤波后的曲线重合。

3.实现

采用Eigen库，可容易实现：

1)范德蒙矩阵生成及相关矩阵运算（主要求B）

2)卷积运算或FIR滤波

也可以使用matlab或其他工具生成滤波器系数，再使用FIR滤波器进行滤波器（注意信号边缘的处理）。这里就不做过多介绍。

本文介绍了Savitzky-Golay滤波器基本原理。