18、数字滤波器设计详解

最新推荐文章于 2025-11-10 10:54:02 发布

香菜滚出地球

最新推荐文章于 2025-11-10 10:54:02 发布

阅读量107

点赞数

CC 4.0 BY-SA版权

分类专栏： MATLAB助力数字信号处理文章标签：数字滤波器 IIR FIR

本文链接：https://blog.youkuaiyun.com/smartcontract5/article/details/151167324

MATLAB助力数字信号处理专栏收录该内容

19 篇文章 ¥499.90

订阅专栏¥69.90

会员秒杀 ¥9.9 重磅福利

超级会员免费看

数字滤波器设计详解

1. 数字滤波器设计基础

1.1 关键公式与参数计算

在数字滤波器设计中，有一系列关键公式用于计算相关参数。首先是公式（16.32）和（16.33）：
[
A = \frac{\omega_{u} - \omega_{l}}{\omega_{u} + \omega_{l}}
]
[
B = \frac{\omega_{u} \omega_{l}}{\omega_{u}^{2} - \omega_{l}^{2}}
]
接着，通过以下公式计算其他参数：
[
k_{1} = - 10 \log(A)
]
[
k_{2} = - 10 \log(1 + \epsilon^{2})
]
[
\gamma = \sqrt{\frac{A^{2} - 1}{\epsilon^{2}}}
]
[
n = \left\lceil\frac{\log\left(\frac{\gamma + \sqrt{\gamma^{2} - 1}}{\gamma - \sqrt{\gamma^{2} - 1}}\right)}{2 \log\left(\frac{\omega_{r1}}{\omega_{r2}}\right)}\right\rceil
]

1.2 函数 (B(s)) 与阶数 (n) 的关系

不同阶数 (n) 对应的函数 (B(s)) 如下表所示：
| (n) | (B(s)) |
| — | — |
| 1 | (s + 1) |
| 2 | (s^{2} + \sqrt{2}s + 1) |
| 3 | ((s + 1)(s^{2} + s + 1)) |
| 4 | ((s^{2} + 0.7653s + 1)(s^{2} + 1.8477s + 1)) |

1.3 (s) 平面传递函数

通过公式（16.38）可得到 (s) 平面的传递函数：
[
H(s) = \frac{B(s)}{s + \frac{\omega_{u} - \omega_{l}}{2}}
]

2. 双线性变换（IIR 数字滤波器）

2.1 双线性变换原理

设计 IIR 数字滤波器时，采用双线性变换将 (s) 平面的左半部分转换为 (z) 平面单位圆内部。模拟滤波器的传递函数 (H(s)) 通过以下变换：
[
s = 2f_{s} \frac{1 - z^{-1}}{1 + z^{-1}}
]
其中，(f_{s}) 为采样频率，(z^{-1}) 为延迟元素，延迟时间为 (T_{s})。若 (s = j\omega)，则 (z^{-1} = e^{-j\lambda})，(\lambda) 为数字频率（rad），(\omega) 为模拟频率（rad/sec），它们之间的关系为：
[
j\omega = 2f_{s} \frac{1 - e^{-j\lambda}}{1 + e^{-j\lambda}} = 2f_{s} \tan\left(\frac{\lambda}{2}\right)
]
或
[
\lambda = 2 \tan^{-1}\left(\frac{\omega}{2f_{s}}\right)
]

2.2 设计示例

以设计具有特定特性的数字低通滤波器（LPF）为例，步骤如下：
1. 频率转换 ：将数字频率 (\lambda) 转换为模拟频率 (\omega)，公式为 (\omega = 2f_{s} \tan\left(\frac{\lambda}{2}\right))。
2. 确定阶数 (n) 和函数 (B(s)) ：根据相关公式和表格确定阶数 (n) 和对应的 (B(s))。
3. 计算传递函数 (H(s)) ：根据 (B(s)) 和 (\omega_{c}) 计算 (H(s))。
4. 双线性变换得到 (H(z)) ：将 (s = 2f_{s} \frac{1 - z^{-1}}{1 + z^{-1}}) 代入 (H(s)) 得到 (H(z))。

以下是一个具体示例的 MATLAB 代码：

clc;
clear all; 
[a b] = butter (2, 0.125);
a = [0.0300 0.0599 0.0300];
b = [1.0000 -1.4542 0.5741];
H_z = (0.0300 + 0.0599*z^(-1) + 0.0300*z^(-2)) / (1 - 1.4542*z^(-1) + 0.5741*z^(-2));

3. 高阶 IIR 数字滤波器

3.1 高阶滤波器问题与解决方案

高阶递归滤波器对量化误差和溢出非常敏感，因此高阶 IIR 滤波器通常设计为级联的双二次部分。

3.2 设计示例

以设计四阶巴特沃斯型 IIR 低通数字滤波器为例，具体步骤如下：
1. 确定相对截止频率 ：如相对截止频率为 (\frac{\pi}{8})。
2. 使用 MATLAB 计算系数 ：

clc;
clear all;
[a b] = butter(4, 0.125);

转换为双二次部分 ：

[sos G] = tf2sos(a,b);

3.3 系统函数表示

最终系统函数可以表示为：
[
H(z) = G \frac{(1 + 2z^{-1} + z^{-2})}{(1 - 1.3651z^{-1} + 0.4776z^{-2})} \cdot \frac{(1 + 2z^{-1} + z^{-2})}{(1 - 1.6117z^{-1} + 0.7445z^{-2})}
]

4. IIR 数字高通、带通和带阻滤波器设计

4.1 设计原理

双线性变换同样适用于模拟高通、带通或带阻系统函数，以获得相应的数字滤波器。

4.2 设计示例

4.2.1 四阶高通 IIR 滤波器

clc;
clear all;
[a b] = butter(4,0.125,'high');
freqz(a,b);
[sos G] = tf2sos(a,b);
title('Response of HPF');

4.2.2 二阶（四阶）带通滤波器

clc;
clear all;
[a b] = butter(2,[0.5 0.75]);
freqz(a,b);
[sos G] = tf2sos(a,b);

4.2.3 四阶（八阶）带通滤波器

clc;
clear all;
[a b] = butter(4,[0.2 0.4]);
freqz(a,b);
axis([0 1 -40 0]);
[sos G] = tf2sos(a,b);

4.2.4 四阶（八阶）带阻滤波器

clc;
clear all;
[a b] = butter(4,[0.125 0.25], 'stop');
freqz(a,b);
axis([0 1 -40 0]);
[sos G] = tf2sos(a,b);

5. 使用 MATLAB 设计 IIR 滤波器

5.1 低通滤波器设计

clc;
clear all;
fs = 4000;
fc = 1500;
C = 2 - cos(2*pi*fc/fs);
A(2) = sqrt(C^2 - 1) - C;
A(1) = 1;
B = 1 + A(2);
[H, w] = freqz(B, A);
H = 20*log10(abs(H));
plot(w/(2*pi)*fs/1000, H, fc/1000, -3, 'o');
title('IIR - LPF Response');
xlabel('Frequency (kHz)');
ylabel('Magnitude (dB)');

5.2 高通滤波器设计

clc;
clear all;
fs = 20000;
fc = 1000;
C = 2 - cos(2*pi*fc/fs);
A(2) = sqrt(C^2 - 1);
A(1) = 1;
B = 1 - A(2);
[H, w] = freqz(B, A);
H = 20*log10(abs(H));
plot(w/(2*pi)*fs/1000, H);
title('IIR - HPF Response');
xlabel('Frequency (k Hz)');
ylabel('Magnitude (dB)');

5.3 带通滤波器设计

fs = 25000;
fc = 7000;
Bw = 4000;
A(3) = exp(-2*pi*(Bw/fs));
A(2) = -4*A(3)/(1 + A(3))*cos(2*pi*fc/fs);
A(1) = 1;
B = 1 - A(3) * sqrt(1 - A(2)*A(2)/(4*A(3)));
[H, w] = freqz(B, A);
H = 20*log10(abs(H));
plot(w/(2*pi)*fs/1000, H);
title('IIR BPF Response');
xlabel('Frequency (k Hz)');
ylabel('Magnitude (dB)');

6. FIR 滤波器设计

6.1 线性相位条件

对于某些 FIR 滤波器要实现线性相位，需满足以下条件：
[
h(n) = h(N - 1 - n), n = 0, 1, \cdots, N - 1
]
其中 (N) 为系数的奇数个数。

6.2 设计步骤

频率单位转换 ：将 (\omega) 的单位从 rad/sec 转换为 rad。
[
\omega_{1_{new}} = \frac{\pi}{2} f_{L}T_{s}
]
[
\omega_{2_{new}} = \frac{\pi}{2} f_{H}T_{s}
]
[
\omega_{1_{new}} = \omega_{c}
]
计算样本数 (N) ：(N) 必须为奇数。
[
N = \frac{2k}{\pi(\omega_{2_{new}} - \omega_{1_{new}})}
]
计算参数 (a) ：
[
a = \frac{N - 1}{2}
]
选择窗口类型和滤波器类型 。

6.3 使用窗口设计 FIR 滤波器

6.3.1 原理

通过截断无限长脉冲响应序列来获得有限长脉冲响应。若 (H_{d}(e^{j\omega})) 是理想期望频率响应，则：
[
H_{d}(e^{j\omega}) = \sum_{n = -\infty}^{\infty} h_{d}(n) e^{-jn\omega}
]
其中 (h_{d}(n)) 是相应的脉冲响应序列：
[
h_{d}(n) = \frac{1}{2\pi} \int_{-\pi}^{\pi} H_{d}(e^{j\omega}) e^{jn\omega} d\omega
]
通过窗口函数 (w(n)) 截断 (h_{d}(n)) 得到 (h(n))：
[
h(n) = h_{d}(n) w(n)
]
例如，矩形窗口：
[
w(n) = \begin{cases}
1, & 0 \leq n \leq N - 1 \
0, & \text{otherwise}
\end{cases}
]

6.3.2 设计示例

以设计一个阶数 (N = 15) 的低通滤波器为例，采样频率 (f_{s} = 4kHz)，通带截止频率 (f_{0} = 1kHz)。

fS = 4000;
f0 = 1000;
omega = pi * f0 / fS;
H_ideal = @(w) (abs(w) <= omega);
n = -7:7;
hd = (sin(omega * n))./(pi * n);
hd(8) = omega / pi;
h = hd .* rectwin(15)';

7. FIR 滤波器的输入输出方程

FIR 滤波器长度为 (L)（阶数为 (L - 1)）的输入输出方程为：
[
y(n) = \sum_{i = 0}^{L - 1} b_{i} x(n - i)
]
其中 (b = [b_{0}, b_{1}, \cdots, b_{L - 1}]^{T}) 为滤波器系数向量，(x(n) = [x(n), x(n - 1), \cdots, x(n - L + 1)]^{T}) 为信号向量。

7.1 设计示例

以设计一个 5 抽头 FIR 带通滤波器为例，步骤如下：
1. 确定归一化截止频率 ：

fs = 8000;
fL = 2000;
fH = 2400;
WL = 2*pi*fL/fs;
WH = 2*pi*fH/fs;

计算非因果 FIR 系数 ：

M = 2;
n = -M:M;
h = (sin(WH * n) - sin(WL * n))./(pi * n);
h(M + 1) = (WH - WL) / pi;

延迟得到滤波器系数 ：

b = h;

计算传递函数和差分方程 ：

H_z = b(1) + b(2)*z^(-1) + b(3)*z^(-2) + b(4)*z^(-3) + b(5)*z^(-4);
y_n = b(1)*x_n + b(2)*x_n_1 + b(3)*x_n_2 + b(4)*x_n_3 + b(5)*x_n_4;

计算频率响应 ：

W = 0:0.01:2*pi;
H = b(1) + b(2)*exp(-1j*W) + b(3)*exp(-2j*W) + b(4)*exp(-3j*W) + b(5)*exp(-4j*W);
plot(W/(2*pi)*fs, 20*log10(abs(H)));

8. 数字滤波器设计总结

8.1 IIR 与 FIR 滤波器对比

滤波器类型	优点	缺点	适用场景
IIR	阶数低，计算量小	相位特性差	窄带、锐截止且对相位要求不高的场合
FIR	相位特性好	阶数高，计算量大	对波形形状要求高的场合
### 8.2 设计流程总结
数字滤波器设计流程如下：
1. 规格制定：在频域中确定所需的频率响应。
2. 设计标准选择：选择衡量滤波器质量的标准。
3. 实现：选择滤波器类型、结构和长度，并选择算法计算所需的滤波器系数。
4. 量化：对给定字长的滤波器系数进行量化，分析有限精度效应，并优化量化滤波器以减少量化误差。
5. 实现：使用仿真验证设计的滤波器，并在软件（或硬件）中实现。

通过以上内容，我们详细介绍了数字滤波器的设计方法，包括 IIR 和 FIR 滤波器的设计原理、步骤以及使用 MATLAB 进行设计的示例。希望这些内容能帮助你更好地理解和应用数字滤波器设计。

9. 数字滤波器设计案例分析

9.1 案例一：消除特定频率干扰的滤波器设计

需求是设计一个数字滤波器，采样率为 1000 Hz，用于消除输入信号中 105 Hz 的正弦分量，且 3 dB 带宽不大于 10 Hz。

操作步骤

初始化参数 ：

clc;
clear all;
FS = 1000;
FL = 100;
FU = 110;

设计滤波器 ：

[a b] = butter(2, [FL,FU]/(FS/2), 'stop');

绘制频率响应 ：

freqz(a,b);
title('Response of SPF');

转换为双二次部分 ：

[SOS G] = tf2sos(a,b);

9.2 案例二：二阶低通数字滤波器设计

已知模拟滤波器传递函数 (H(s) = \frac{\omega_{c}^{2}}{s^{2} + \sqrt{2}\omega_{c}s + \omega_{c}^{2}})，采样频率为 5 Hz，期望低通截止频率为 0.318 Hz（2 rad/sec），期望数字截止频率为 0.4 rad。

操作步骤

计算数字截止频率对应的模拟截止频率 ：

omega_c = 2 * tan(0.4 / 2);

双线性变换 ：

s = tf('s');
H_s = omega_c^2 / (s^2 + sqrt(2)*omega_c*s + omega_c^2);
fs = 5;
z = tf('z', 1/fs);
s_sub = 2*fs*(1 - z^(-1))/(1 + z^(-1));
H_z = subs(H_s, s, s_sub);

10. 数字滤波器设计的优化与注意事项

10.1 优化策略

减少计算量 ：对于 FIR 滤波器，可采用快速卷积算法减少计算量；对于 IIR 滤波器，合理选择滤波器阶数和结构。
提高精度 ：在量化过程中，选择合适的字长，减少量化误差。可使用 MATLAB 等工具进行量化滤波器设计。
改善相位特性 ：对于对相位要求高的应用，优先选择 FIR 滤波器；对于 IIR 滤波器，可采用全通滤波器进行相位校正。

10.2 注意事项

采样频率选择 ：采样频率应满足奈奎斯特采样定理，避免信号混叠。
滤波器阶数选择 ：阶数过高会增加计算量和延迟，阶数过低可能无法满足滤波器性能要求。
量化误差分析 ：在固定点 DSP 处理器上实现滤波器时，要充分考虑量化误差对滤波器性能的影响。

11. 数字滤波器设计的未来发展趋势

11.1 智能化设计

随着人工智能和机器学习技术的发展，未来数字滤波器设计可能会采用智能化方法。例如，使用神经网络自动优化滤波器系数，以满足不同的应用需求。

11.2 多学科融合

数字滤波器设计将与其他学科如通信、控制、生物医学等深度融合，产生更多创新应用。例如，在生物医学信号处理中，设计专门的滤波器用于提取微弱的生物电信号。

11.3 硬件集成化

为了提高滤波器的处理速度和降低功耗，未来数字滤波器可能会更多地集成到硬件芯片中，实现片上系统（SoC）。

12. 总结与展望

12.1 总结

本文全面介绍了数字滤波器的设计方法，涵盖了 IIR 和 FIR 滤波器的设计原理、步骤、MATLAB 实现示例以及案例分析。同时，讨论了数字滤波器设计的优化策略、注意事项和未来发展趋势。

12.2 展望

数字滤波器在现代信号处理中具有重要地位，随着技术的不断发展，其设计方法和应用领域也将不断拓展。希望读者通过本文的学习，能够掌握数字滤波器的设计技巧，并在实际应用中灵活运用。

以下是数字滤波器设计流程的 mermaid 流程图：

graph LR
    A[规格制定] --> B[设计标准选择]
    B --> C[实现]
    C --> D[量化]
    D --> E[验证与实现]

通过以上内容，我们对数字滤波器设计有了更深入的了解，从基础理论到实际案例，再到未来发展趋势，希望能为读者在数字滤波器设计领域的学习和实践提供有价值的参考。