接上文：https://blog.csdn.net/weixin_53877178/article/details/122470759
目录
一、课题的任务
二、内容、步骤和要求
（1）语音信号的采集
（2）语音信号的分析
（3）语音信号的滤波
1）IIR滤波器的设计
1.标准理想滤波器的特性
2.数字滤波器的技术指标（以低通滤波器的性能指标为例）
3.常用的三种模拟滤波器
• 巴特沃斯滤波器（Butterworth）
巴特沃斯滤波器是指在通带内的幅度特性具有最大平坦特性，是一全极点型滤波器。
滤波器的平方幅频特性函数的前(2N－1)阶导数在模拟频率W=0处都为零。
• 切比雪夫滤波器（Chebyshev）
切比雪夫滤波器在过渡带比巴特沃斯滤波器的衰减快，但频率响应的幅频特性不如后者平坦。切比雪夫滤波器和理想滤波器的频率响应曲线之间的误差最小，但是在通频带内存在幅度波动。
在通带（或称“通频带”）上频率响应幅度等波纹波动的滤波器称为“I型切比雪夫滤波器”；
在阻带（或称“阻频带”）上频率响应幅度等波纹波动的滤波器称为“II型切比雪夫滤波器”。
• 椭圆滤波器（Ellipse）
在通带和阻带等波纹的一种滤波器。椭圆滤波器相比其他类型的滤波器，在阶数相同的条件下有着最小的通带和阻带波动。
4. 滤波器设计方法
• 脉冲响应不变法
优点： w = ΩT.
保证了模拟角频率与数字角频率之间的线性关系。
缺点：在 w＝±pi 附近产生混叠失真。
一般只适用于带限滤波器，高通、带阻滤波器不能直接进行变换。
• 双线性变换法
优点：不会出现由于高频部分超过折叠频率而混淆到低频部分去的现象。
缺点：使数字滤波器与模拟滤波器在响应与频率的对应关系上发生畸变。
以低通滤波器为例设计
相关指标
Wp：通带截止频率；
Ws：阻带起始频率；
Rp：通带最大衰减；
Rs：阻带最小衰减；
模拟滤波器部分
以巴特沃斯滤波器为例
[n, Wn] = buttord(Wp, Ws, Rp, Rs, 's');
获取模拟巴特沃斯滤波器的最小阶数 n 和截止频率 Wn.
以弧度/秒为单位指定频率 Wp 和 Ws. 带或阻带可以是无限的。
's'表示模拟滤波器。
Wp, Ws的取值（以下wp ws wp1 wp2 ws1 ws2 单位均为Hz）：
脉冲响应不变法：
低通 Wp = 2*pi*wp（通带截止频率）; Ws = 2*pi*ws（通带起始频率）;
高通 Wp = 2*pi*wp（通带起始频率）; Ws = 2*pi*ws（阻带截止频率）;
带通 Wp = 2*pi*[wp1, wp2]（低端通带起始频率 高端通带截止频率）;
Ws = 2*pi*[ws1, ws2]（低端阻带截止频率 高端阻带起始频率）;
带阻 Wp = 2*pi*[wp1, wp2]（低端通带截止频率 高端通带起始频率）;
Ws = 2*pi*[ws1, ws2]（低端阻带起始频率 高端阻带截止频率）;
双线性变换法：
在上述基础上先除以 Fs ，（采样频率）转化为数字滤波器的频率，设结果为 w，再利用2*Fs*tan(w/2) 转化为模拟滤波器指标。
[b, a] = butter(n, Wn, ftype, 's');根据 ftype 和 Wn ，设计低通、高通、带通或带阻巴特沃斯滤波器。
ftype 的取值：
低通 [b, a] = butter(n, Wn, 's');
高通 [b, a] = butter(n, Wn, ‘high’, 's');
带通 [b, a] = butter(n, Wn, 's');
带阻 [b, a] = butter(n, Wn, ‘stop’, 's');
其他模拟滤波器
切比雪夫Ⅰ型滤波器
[N, Wn] = cheb1ord(Wp, Ws, Ap, As, 's');
[b, a] = cheby1(N, Ap, Wc, 's');切比雪夫Ⅱ型滤波器
[N, Wn] = cheb2ord(Wp, Ws, Ap, As, 's');
[b, a] = cheby2(N, Ap, Wc, 's');椭圆滤波器
[N, Wn] = ellipord(Wp, Ws, Ap, As, 's');
[b, a] = ellip(N, Ap, As, Wc, 's');滤波器设计方法部分
脉冲响应不变法
[bz, az] = impinvar(b, a, fs);创建分别具有分子和分母系数bz和az的数字滤波器，其脉冲响应等于具有系数b和a的模拟滤波器的脉冲响应，按1/fs缩放，其中fs是采样率。b, a 为上述模拟滤波器设计所得。
双线性变换法
[bz, az] = bilinear(b, a, fs);将模拟域传递函数变换为数字域传递函数，Fs为取样频率。
滤波器效果展示：幅频响应特性曲线图
W = linspace(0, pi, 400*pi);
[H, W] = freqz(Bz, Az, W);
plot(app.filter_frequency, W*FS/(2*pi), 20*log10(abs(H)));y = linspace(x1, x2, n) 生成 n 个点。这些点的间距为 (x2-x1)/(n-1).
[h, w] = freqz(b, a, n)；返回数字滤波器的 n 点频率响应向量h和相应的角频率向量 w.
W*FS/(2*pi)：转换为Hz.
20*log10(abs(H))：转化为dB.
以低通滤波器为例设计的完整代码展示
%wp，ws自行设计，单位Hz
%一次只能出一个图，多个图同框可采用 subplot 函数
%巴特沃斯低通滤波器
%脉冲响应不变法
Wp = 2*pi*wp;
Ws = 2*pi*ws;
[N, Wc] = buttord(Wp, Ws, Ap, As, 's');
[B, A] = butter(N, Wc, 's');
[Bz, Az] = impinvar(B, A, FS);
W = linspace(0, pi,400*pi);
[H, W] = freqz(Bz, Az, W);
plot(W*FS/(2*pi), 20*log10(abs(H)));
%双线性变换法
Wp = 2*pi*wp/FS;
Ws = 2*pi*ws/FS;
[N, Wc] = buttord(2*FS*tan(Wp/2), 2*FS*tan(Ws/2), Ap, As, 's');
[B, A] = butter(N, Wc, 's');
[Bz, Az]=bilinear(B,A,FS);
W = linspace(0, pi,400*pi);
[H, W] = freqz(Bz, Az, W);
plot(W*FS/(2*pi), 20*log10(abs(H)));2）FIR滤波器的设计（窗函数设计法）
1.六种窗函数的基本参数比较
设计思想：
根据阻带最小衰减选择合适的窗函数进行设计。
①计算过渡带（过渡带见上述“数字滤波器的技术指标”）
低通 Wp = 2*pi*wp/Fs（通带截止频率）; Ws = 2*pi*ws/Fs（通带起始频率）;
W = Ws - Wp;
高通 Wp = 2*pi*wp/Fs（通带起始频率）; Ws = 2*pi*ws/Fs（阻带截止频率）;
W = Wp - Ws;
带通 Wp1 = 2*pi*wp1/Fs (低端通带起始频率);
Wp2 = 2*pi*wp2/Fs (高端通带截止频率);
Ws1 = 2*pi*ws1/Fs (低端阻带截止频率):
Ws2 = 2*pi*ws2/Fs (高端阻带起始频率);
W1 = Ws1 - Wp1;
W2 = Wp2 - Ws2;
W = min(W1, W2); （取过渡带的算术中心频率两者最小值）
带阻 Wp1 = 2*pi*wp1/Fs (低端通带截止频率);
Wp2 = 2*pi*wp2/Fs (高端通带起始频率);
Ws1 = 2*pi*ws1/Fs (低端阻带起始频率):
Ws2 = 2*pi*ws2/Fs (高端阻带截止频率);
W1 = Ws1 - Wp1;
W2 = Wp2 - Ws2;
W = min(W1, W2); （取过渡带的算术中心频率两者最小值）
②求阶数N
根据“六种窗函数的基本参数比较”求N.
N = ceil(x*pi/W);x 为相应窗对应的“过度带宽精确值”。
当 N 为偶数时，高通、带阻的设计会出现问题。于是利用下面代码将偶数变为奇数。
N = N + mod(N + 1, 2);mod 为求余函数。
③ 计算理想低通滤波器通带截止频率(关于pi归一化)
低通、高通：Wc = (Wp + Ws)/2;
带通、带阻：Wc1 = (Wp1 + Ws1)/2;
Wc2 = (wp2 + Ws2)/2;
Wc = [Wc1, Wc2];
b = fir1(N - 1, Wc, ftype, window);根据ftype,Wc的值和相应的窗函数设计滤波器。
ftype，window 的值：
%矩形窗
% 低通滤波器
hn = fir1(N - 1, Wc/pi, boxcar(N));
% 高通滤波器
hn = fir1(N - 1, Wc/pi, 'high', boxcar(N));
% 带通滤波器
hn = fir1(N - 1, Wc/pi, boxcar(N));
% 带阻滤波器
hn = fir1(N - 1, Wc/pi, 'stop', boxcar(N));
%三角窗
% 低通滤波器
hn = fir1(N - 1, Wc/pi, triang(N));
% 高通滤波器
hn = fir1(N - 1, Wc/pi, 'high', triang(N));
% 带通滤波器
hn = fir1(N - 1, Wc/pi, triang(N));
% 带阻滤波器
hn = fir1(N - 1, Wc/pi, 'stop', triang(N));
%汉宁窗
% 低通滤波器
hn = fir1(N - 1, Wc/pi, hanning(N));
% 高通滤波器
hn = fir1(N - 1, Wc/pi, 'high', hanning(N));
% 带通滤波器
hn = fir1(N - 1, Wc/pi, hanning(N));
% 带阻滤波器
hn = fir1(N - 1, Wc/pi, 'stop', hanning(N));
%哈明窗
% 低通滤波器
hn = fir1(N - 1, Wc/pi, hamming(N));
% 高通滤波器
hn = fir1(N - 1, Wc/pi, 'high', hamming(N));
% 带通滤波器
hn = fir1(N - 1, Wc/pi, hamming(N));
% 带阻滤波器
hn = fir1(N - 1, Wc/pi, 'stop', hamming(N));
%布莱克曼窗
% 低通滤波器
hn = fir1(N - 1, Wc/pi, blackman(N));
% 高通滤波器
hn = fir1(N - 1, Wc/pi, 'high', blackman(N));
% 带通滤波器
hn = fir1(N - 1, Wc/pi, blackman(N));
% 带阻滤波器
hn = fir1(N - 1, Wc/pi, 'stop', blackman(N));
%凯塞窗
% 低通滤波器
hn = fir1(N - 1, Wc/pi, kaiser(N,7.865));
% 高通滤波器
hn = fir1(N - 1, Wc/pi, 'high', kaiser(N,7.865));
% 带通滤波器
hn = fir1(N - 1, Wc/pi, kaiser(N,7.865));
% 带阻滤波器
hn = fir1(N - 1, Wc/pi, 'stop', kaiser(N,7.865));
滤波器效果展示：幅频响应特性曲线图
[H, W]= freqz(hn, 1, 256);
plot(W*FS/(2*pi), 20*log10(abs(H)));以矩形窗设计为例的完整代码
%Wp, Ws, Wp1, Wp2, Ws1, Ws2自行设计
%一次只能出一个图，多个图同框可采用 subplot 函数
%矩形窗
%低通滤波器
W = Ws-Wp; %计算过渡带
N = ceil(1.8*pi/W); %求阶数N
Wc = (Wp + Ws)/2; %计算理想低通滤波器通带截止频率(关于pi归一化)
hn = fir1(N - 1, Wc/pi, boxcar(N));
[H, W] = freqz(hn, 1, 500);
plot(W*FS/(2*pi), 20*log10(abs(H)));
%高通滤波器
W = Wp-Ws; %计算过渡带
N = ceil(1.8*pi/W); %求阶数N
N = N + mod(N + 1, 2);
Wc = (Wp + Ws)/2; %计算理想低通滤波器通带截止频率(关于pi归一化)
hn = fir1(N - 1, Wc/pi, 'high', boxcar(N));
[H, W] = freqz(hn, 1, 256);
plot(W*FS/(2*pi), 20*log10(abs(H)));
%带阻滤波器
W1 = Ws1 - Wp1; %计算过渡带
W2 = Wp2 - Ws2; %计算过渡带
W = min(W1,W2); %取过渡带的算术中心频率两者最小值
N = ceil(1.8*pi/W); %求阶数N
N = N + mod(N + 1, 2);
Wc1 = (Wp1 + Ws1)/2; %计算理想低通滤波器通带截止频率(关于pi归一化)
Wc2 = (Wp2 + Ws2)/2; %计算理想低通滤波器通带截止频率(关于pi归一化)
Wc = [Wc1, Wc2];
hn = fir1(N - 1, Wc/pi, 'stop', boxcar(N));
[H, W] = freqz(hn, 1, 256);
plot(W*FS/(2*pi), 20*log10(abs(H)));
%带通滤波器
W1 = Wp1 - Ws1; %计算过渡带
W2 = Ws2 - Wp2; %计算过渡带
W = min(W1, W2); %取过渡带的算术中心频率两者最小值
N = ceil(1.8*pi/W); %求阶数N
Wc1 = (Wp1 + Ws1)/2; %计算理想低通滤波器通带截止频率(关于pi归一化)
Wc2 = (Wp2 + Ws2)/2; %计算理想低通滤波器通带截止频率(关于pi归一化)
Wc = [Wc1, Wc2];
hn = fir1(N - 1, Wc/pi, boxcar(N));
[H, W] = freqz(hn, 1, 256);
plot(W*FS/(2*pi), 20*log10(abs(H)));3）滤波器对语音信号的处理
IIR滤波器：
y = filter(BZ, AZ, x);BZ，AZ 为 impinvar 或者 bilinear 函数返回值。
x 为 [x,Fs]=audioread(app.path); 中的x.
FIR滤波器：
y = fftfilt(hn, x);hn 为 fir1 函数返回值。
x 为 [x,Fs]=audioread(app.path); 中的x.
下文：基于Matlab App Designer的语音信号分析与处理（三）：App的设计_演员的必备修养的博客-CSDN博客
最终效果展示：
完整资源获取：
2876
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
