在数字信号处理中,滤波器是用于去除信号中不需要的频率成分的重要工具。然而,当使用滤波器进行多次滤波时,由于相位延迟的累积,会出现截止频率降低的问题。本文将探讨这一问题,并介绍如何正确地校正滤波器的截止频率。
根据David A. Winter在其著作《Biomechanics and Motor Control of Human Movement》第四版第69页的描述,以及D. Gordon E. Robertson等人在其著作《Research Methods in Biomechanics》第二版第288页的描述,对于零相位延迟的Butterworth滤波器,需要对截止频率进行校正。这是因为每次滤波过程都会引入相位延迟,为了补偿这种延迟,滤波器需要反向运行,每次滤波过程都会使得-3dB截止频率降低。Winter展示了如何通过应用校正因子来补偿二阶Butterworth滤波器的截止频率衰减。
Winter使用的调整后的角截止频率的变量名可能会引起误解,因此本文采用了与Robertson等人相似的变量名(仅将Ω替换为U)。
C = (((2^(1/filter_passes))-1)^(1/4)); % David A. Winter Butterworth校正因子
Wc = 2*pi*cutoff_frequency; % 角截止频率
Uc = tan(Wc/(2*sampling_frequency)); % 调整后的角截止频率
Un = Uc / C; % David A. Winter校正
请注意,这两位作者的早期版本的书籍和论文可能没有使用这种更新的校正方法。
如果在互联网上搜索关于截止频率调整的信息,会发现许多人都错误地实现了这种方法。他们没有将校正因子应用于调整后的角截止频率(Uc),而是直接应用于截止频率。例如,对于双次滤波,有些情况下校正被错误地执行如下:
f_corrected = f_cutoff / 0.802
虽然这种方法给出的结果并不完全错误,但在较低频率下仍然有偏差。为了研究如果直接将校正因子应用于截止频率,可以从校正后的参数Un(见结果部分如何实现)推导出截止频率。结果表明,这是一个频率依赖的因子。将其称为频率依赖因子Cf:
Cf = f_cutoff / f_corrected
如果有人知道如何直接得到这个因子,请告知。
根据Winter的方法,获取二阶低通Butterworth系数的正确方式是:
C = (((2^(1/filter_passes))-1)^(1/4)); % David A. Winter butterworth校正因子
Wc = 2*pi*cutoff_frequency; % 角截止频率
Uc = tan(Wc/(2*sampling_frequency)); % 调整后的角截止频率
Un = Uc / C; % David A. Winter校正
k1 = sqrt(2) * Un;
k2 = Un^2;
a0 = k2 / (1 + k1 + k2);
a1 = 2 * a0;
a2 = a0;
k3 = a1 / k2;
b1 = k3 - a1;
b2 = 1 - a1 - k3;
b = [ a0 a1 a2 ];
a = [ 1 -b1 -b2 ];
data_filtered = myfilter(data,b,a,'low',filter_passes);
如果想使用一个只接受截止频率作为输入的现有函数(如MATLAB的butter函数),那么可以将校正后的参数Un转换回频率。
filter_passes = 2; % filtfilt使用2次通过
C = (((2^(1/filter_passes))-1)^(1/4)); % David A. Winter butterworth校正因子
Wc = 2*pi*cutoff_frequency; % 角截止频率
Uc = tan(Wc/(2*sampling_frequency)); % 调整后的角截止频率
Un = Uc / C; % David A. Winter校正
f_corrected = atan(Un)*sampling_frequency/pi; % 校正后的截止频率
[b, a] = butter(2, f_corrected / (sampling_frequency / 2), 'low');
data_filtered = filtfilt(b, a, data);
现在可以测试参考单次滤波、未校正的多次滤波、错误的频率校正的多次滤波以及Winter的校正的多次滤波之间的差异。为了测试多次通过的影响,使用了一个以滤波次数为输入的滤波器(见附件文件)。下图显示的是滤波噪声的频谱。
在两次通过时,未校正的截止频率变为134Hz而不是150Hz,错误校正的截止频率变得过高,而Winter提出的校正截止频率非常准确,如下图所示:
在更高的通过次数下,差异变得更大(在4次通过时,错误校正的截止频率在fs/3时已经达到fs/2)。
同样的原则也适用于高通滤波器。