数字信号处理_IIR及FIR设计
- 1、下载文档前请自行甄别文档内容的完整性,平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益,请联系客服反馈,我们会尽快为您处理(人工客服工作时间:9:00-18:30)。
IIR 数字滤波器设计作业
通信工程 201304015
7.12
设计一个数字切比雪夫I 型带通滤波器,给定指标为: (1) 波纹20dB p R ≤,当 200Hz 400Hz f ≤≤ (2) 衰减20dB s A ≥, 100Hz f ≤, 600Hz f ≥ (3) 抽样频率2kHz s f =
试用○
1冲激响应不变法,○2双线性变换法进行设计,最后写出()H z 的表达式,并画出系统的幅频响应特性(dB)。
解:
○
1冲激响应不变法设计程序如下 %冲激响应不变法,ex712.m clc;clear all
OmegaP1=2*pi*200;OmegaP2=2*pi*400;%带通截止频率 OmegaS1=2*pi*100;OmegaS2=2*pi*600;%1.5kHz 阻带 Rp=2;%波纹系数 As=20;%阻带衰减dB Fs=2*10^3;%抽样频率2khz
OmegaP=[OmegaP1,OmegaP2];OmegaS=[OmegaS1,OmegaS2]; wp=OmegaP/Fs;ws=OmegaS/Fs;%等效数字频率
[N,OmegaC]=cheb1ord(OmegaP,OmegaS,Rp,As,'s')%滤波器阶数截止频率 [b,a]=cheby1(N,Rp,OmegaC,'s');%AF 系统函数的分子 分母 [bz,az]=impinvar(b,a,Fs)%冲击不变法AF to DF w0=[wp,ws]%四个频点 Hx=freqz(bz,az,w0);%检验
[H,w]=freqz(bz,az);%计算0~pi 上的响应
dbHx=-20*log10(abs(Hx)/max(abs(H)))%归一化并求dB dbH=20*log10(abs(H)/max(abs(H)));%归一化的频率响应 plot(w/pi,dbH);%画图
xlabel('\Omega/\pi');ylabel('dB'); axis([0,1,-60,5]);grid
运行,得
N = 3;
OmegaC =1.0e+03 * 1.2566 2.5133 w0 = 0.6283 1.2566 0.3142 1.8850
bz = -0.0000 0.0272 -0.0581 0.0109 0.0437 -0.0237 0 az = 1.0000 -3.3030 6.0060 -6.7463 5.1356 -2.4093 0.6290 dbHx = 2.0022 2.0015 41.6739 30.7707
图 1冲激响应不变法设计IIR 带通滤波器
由程序返回得到的数值可以得知,这是一个3阶带通滤波器。
转换为模拟原型滤波器的频率见OmegaC 返回的值。
数字滤波器的边界频率见w0返回的值。
由bz,az 返回值可得
-1-2-3-4-5
-1-2-3-4-5-6
0.02720.05810.01090.01370.0237()1 3.303 6.006 6.7463 5.1356 2.40930.629z z z z z H z z z z z z z -++-=-+-+-+
由dbHx 返回值可知在四个边界频率处幅度响应大小: 对应到模拟频率,即有
满足题目要求的设计指标,通带波纹不大于2dB, 截止频率处衰减不小于20dB 。
○
2双线性变换法设计程序如下 %双线性变换法,ex7122m clc;clear all
OmegaP1=2*pi*200;OmegaP2=2*pi*400;%带通截止频率 OmegaS1=2*pi*100;OmegaS2=2*pi*600;%1.5kHz 阻带 Rp=2;%波纹系数 As=20;%阻带衰减dB Fs=2*10^3;%抽样频率2khz
OmegaP=[OmegaP1,OmegaP2];OmegaS=[OmegaS1,OmegaS2]; wp=OmegaP/Fs;ws=OmegaS/Fs;%等效数字频率
0.10.20.30.4
0.50.60.70.80.91
-60-40-20
Ω/π
d B
OmegaP_t=2*Fs*tan(wp/2);OmegaS_t=2*Fs*tan(ws/2); [N,OmegaC]=cheb1ord(OmegaP_t,OmegaS_t,Rp,As,'s') [b,a]=cheby1(N,Rp,OmegaC,'s');%AF 系统函数的分子 分母 [bz,az]=bilinear(b,a,Fs)%双线性变换法 AF to DF w0=[wp,ws];%四个频点
Hx=freqz(bz,az,w0);%计算四个频点上对应的幅度响应 [H,w]=freqz(bz,az);%计算0~pi 上的响应
dbHx=-20*log10(abs(Hx)/max(abs(H)))%归一化并求dB dbH=20*log10(abs(H)/max(abs(H)));%归一化的频率响应 plot(w/pi,dbH);%画图
xlabel('数字频率域频率\Omega/\pi');ylabel('dB'); axis([0,1,-300,5]); grid
运行,得
N = 2;
OmegaC =1.0e+03 *1.2997 2.9062 w0 = 0.6283 1.2566 0.3142 1.8850
bz = 0.0512 -0.0000 -0.1024 -0.0000 0.0512 az = 1.0000 -2.0733 2.4881 -1.5944 0.6125 dbHx = 1.9997 1.9997 24.1325 22.3051
图 2 双线性变换法设计IIR 带通滤波器
由程序返回得到的数值可以得知,这是一个2阶带通滤波器。
转换为模拟原型滤波器的频率见OmegaC 返回的值。
数字滤波器的边界频率见w0返回的值。
由bz,az 返回值可得-2-4
-1-2-3-4
0.05120.10240.0512()1 2.0733 2.4881 1.59440.6125z z H z z z z z -+=-+-+
由dbHx 返回值可知在四个边界频率处幅度响应大小:
对应到模拟频率,即有
00.10.20.3
0.40.50.60.70.80.91
-300
-200-100
数字频率域频率Ω/π
d B
满足题目要求的设计指标,通带波纹不大于2dB, 截止频率处衰减不小于20dB 。
7.14
设计一个数字切比雪夫I 型带阻滤波器,给定指标为: (1) 衰减30dB s A ≥,当 1kHz 2kHz f ≤≤ (2) 波纹3dB p R ≤,当 500Hz f ≤, 3kHz f ≥ (3) 抽样频率10kHz s f =
试用双线性变换法进行设计,最后写出()H z 的表达式,并画出系统的幅频响应特性
(dB)。
解:
双线性变换法设计程序如下
%冲激响应不变法,ex714.m %设计数字切比雪夫I 型带阻滤波器 %双线性变换法 clc;clear all
OmegaP1=2*pi*500;OmegaP2=2*pi*3000;%带通截止频率 OmegaS1=2*pi*1000;OmegaS2=2*pi*2000;%阻带截止频率 Rp=3;%通带波纹dB As=30;%阻带衰减dB
Fs=10*10^3;%抽样频率10khz
OmegaP=[OmegaP1,OmegaP2];OmegaS=[OmegaS1,OmegaS2]; wp=OmegaP/Fs;ws=OmegaS/Fs;%等效数字频率
OmegaP_t=2*Fs*tan(wp/2);OmegaS_t=2*Fs*tan(ws/2);
[N,OmegaC]=cheb1ord(OmegaP_t,OmegaS_t,Rp,As,'s')%AF 阶数和截至频率 [b,a]=cheby1(N,Rp,OmegaC,'stop','s');%AF 系统函数的分子 分母 [bz,az]=bilinear(b,a,Fs)%双线性变换法F to DF w0=[wp,ws];%四个频点
Hx=freqz(bz,az,w0);%计算两个频点上对应的幅度响应 [H,w]=freqz(bz,az);%计算0~pi 上的响应
dbHx=-20*log10(abs(Hx)/max(abs(H)))%归一化并求dB dbH=20*log10(abs(H)/max(abs(H)));%归一化的频率响应 plot(w/pi,dbH);%画图
xlabel('数字频率域频率\Omega/\pi');ylabel('dB');
axis([0,1,-150,5]); grid
运行,得
N = 3;
OmegaC = 1.0e+04 * 0.3430 2.7528 w0 = 0.3142 1.8850 0.6283 1.2566
bz = 0.0946 -0.3508 0.7174 -0.8802 0.7174 -0.3508 0.0946 az = 1.0000 -1.4600 0.3178 -0.3506 0.6885 0.2289 -0.3824 dbHx = 0.2548 3.0000 39.8925 39.8925
图 3双线性变换法设计IIR 带阻滤波器
由程序返回得到的数值可以得知,这是一个2阶带阻滤波器。
转换为模拟原型滤波器的频率见OmegaC 返回的值。
数字滤波器的边界频率见w0返回的值。
由bz,az 返回值可得
-1-2-3-4-5-6
-1-2-3-4-5-6
0.0946-0.3508+0.7174-0.88020.71740.35080.0946()1 1.46000.37180.35060.68850.22890.3824z z z z z z H z z z z z z z +-+=-+-+--
由dbHx 返回值可知在四个边界频率处幅度响应大小: 对应到模拟频率,即有
满足题目要求的设计指标,通带波纹不大于3dB, 截止频率处衰减大于30dB 。
7.15
00.10.20.3
0.4
0.50.60.70.80.91
-150
-100-50
数字频率域频率Ω/π
d B
设计一个数字切比雪夫Ⅱ型带阻滤波器,给定指标为: (4) 衰减30dB s A ≥,当 10kHz 20kHz f ≤≤ (5) 波纹3dB p R ≤,当 5kHz f ≤, 30kHz f ≥ (6) 抽样频率100kHz s f =
试用双线性变换法进行设计,最后写出()H z 的表达式,并画出系统的幅频响应特性
(dB)。
设计时请先想一想,这一题和上一题有什么相似处。
由此应该得出什么结
论。
解:
分析:这一题的各频率指标要求均是7.14的10倍,抽样频率也是10倍。
此时对应的数字频率是一样的。
如果用切比雪夫Ⅰ型滤波器,双线性变换法设计出来的与7.14中的肯定是一样的(见图4).这说明只要给定的各数字频率参数一样,用同种方法设计出来的数字滤波器是一样的。
程序如下:
%设计数字切比雪夫Ⅱ型带阻滤波器 %双线性变换法 clc;clear all
OmegaP1=2*pi*5000;OmegaP2=2*pi*30000;%带通截止频率 OmegaS1=2*pi*10000;OmegaS2=2*pi*20000;%阻带截止频率 Rp=3;%通带波纹dB As=30;%阻带衰减dB
Fs=100*10^3;%抽样频率10khz
OmegaP=[OmegaP1,OmegaP2];OmegaS=[OmegaS1,OmegaS2]; wp=OmegaP/Fs;ws=OmegaS/Fs;%等效数字频率
OmegaP_t=2*Fs*tan(wp/2);OmegaS_t=2*Fs*tan(ws/2);
[N,OmegaC]=cheb2ord(OmegaP_t,OmegaS_t,Rp,As,'s')%AF 阶数和截至频率 [b,a]=cheby2(N,As,OmegaS,'stop','s');%AF 系统函数的分子 分母 [bz,az]=bilinear(b,a,Fs)%双线性变换法F to DF w0=[wp,ws];%四个频点
Hx=freqz(bz,az,w0);%计算两个频点上对应的幅度响应 [H,w]=freqz(bz,az);%计算0~pi 上的响应
dbHx=-20*log10(abs(Hx)/max(abs(H)))%归一化并求dB dbH=20*log10(abs(H)/max(abs(H)));%归一化的频率响应 plot(w/pi,dbH);%画图
xlabel('数字频率域频率\Omega/\pi');ylabel('dB'); axis([0,1,-150,5]); grid
图 4双线性变换法设计IIR 带阻滤波器(fs=100kHz )
7.17
要求设计一个数字带通滤波器,其抽样频率25kHz s f =,通达截止频率为
15kHz p f =,27kHz p f =,通带衰减0.5dB p R =,阻带截止频率为1 3.5kHz st f =,
28.5kHz st f =,阻带衰减为45dB s A =。
(1) 利用MATLAB 工具箱中ellipord 及ellip 设计椭圆函数滤波器; (2) 利用MATLAB 工具箱中cheblord 及cheb1设计切比雪夫Ⅰ型滤波器; (3) 利用MATLAB 工具箱中cheb2ord 及cheb2设计切比雪夫Ⅱ型滤波器; (4) 利用buttord 及butter 设计巴特沃思型滤波器。
要求每种设计都给出系统函数并画出幅频特性(dB )、相频特性以及单位冲激响应。
解:
(1) 利用MATLAB 工具箱中ellipord 及ellip 设计椭圆函数滤波器:
核心语句为:
%设计椭圆带通滤波器 %双线性变换法 clc;clear all
OmegaP1=2*pi*5000;OmegaP2=2*pi*7000;%带通截止频率 OmegaS1=2*pi*3500;OmegaS2=2*pi*8500;%阻带截止频率 Rp=0.5;%通带波纹dB As=45;%阻带衰减dB
Fs=35*10^3;%抽样频率35khz
OmegaP=[OmegaP1,OmegaP2];OmegaS=[OmegaS1,OmegaS2]; wp=OmegaP/Fs/pi;ws=OmegaS/Fs/pi;%等效数字频率
00.10.20.3
0.40.50.60.70.80.91
-150
-100-50
数字频率域频率Ω/π
d B
w0=[wp,ws];%四个频点
[N,wc]=ellipord(wp,ws,Rp,As)%DF 阶数和截至频率 [bz,az]=ellip(N,Rp,As,wc)%DF 系统函数的分子 分母 Hx=freqz(bz,az,w0*pi);%计算四个频点上对应的幅度响应 [H,w]=freqz(bz,az);%计算0~1上的幅频响应
dbHx=-20*log10(abs(Hx)/max(abs(H)))%归一化并求dB
dbH=20*log10(abs(H)/max(abs(H)));%归一化的频率响应
pha=unwrap(angle(H));%计算0~1上的相频响应
figure(1)
plot(w/pi,dbH);%画幅频图
xlabel('数字频率域频率\Omega/\pi');ylabel('dB'); axis([0,1,-80,5]); grid figure(2)
plot(w/pi,pha);%画相频图
xlabel('数字频率域频率\Omega/\pi');ylabel('Phase\degree'); % axis([0,1,-150,5]); grid
disp('系统传递函数H(z)'); printsys(bz,az,'z'); figure(3)
h=dimpulse(bz,az); stem(h)
xlabel('n');ylabel('Impulse response');
图 5直接设计法设计IIR 椭圆带通滤波器幅度响应
0.10.20.3
0.40.50.60.70.80.91
数字频率域频率Ω/π
d B
图 6直接设计法设计IIR 椭圆带通滤波器相位响应
图 7直接设计法设计IIR 椭圆带通滤波器单位冲激响应
得到四个截止频率处的幅度响应如下:
满足题目设定的指标要求。
(2) 利用MATLAB 工具箱中cheb1ord 及cheby1设计切比雪夫Ⅰ型滤波器:
%设计chebshev I 带通滤波器 clc;clear all
OmegaP1=2*pi*5000;OmegaP2=2*pi*7000;%带通截止频率
0.1
0.2
0.3
0.40.50.60.70.8
0.9
1
-6-4-2024
6数字频率域频率Ω/π
P h a s e \d e g r e e
50
100
150
200
250
n
I m p u l s e r e s p o n s e
OmegaS1=2*pi*3500;OmegaS2=2*pi*8500;%阻带截止频率
Rp=0.5;%通带波纹dB
As=45;%阻带衰减dB
Fs=35*10^3;%抽样频率35khz
OmegaP=[OmegaP1,OmegaP2];OmegaS=[OmegaS1,OmegaS2];
wp=OmegaP/Fs/pi;ws=OmegaS/Fs/pi;%等效数字频率
w0=[wp,ws];%四个频点
[N,wc]=cheb1ord(wp,ws,Rp,As)%DF阶数和截至频率
[bz,az]=cheby1(N,Rp,wc);%DF系统函数的分子分母
Hx=freqz(bz,az,w0*pi);%计算四个频点上对应的幅度响应
[H,w]=freqz(bz,az);%计算0~1上的幅频响应
dbHx=-20*log10(abs(Hx)/max(abs(H)))%归一化并求dB
dbH=20*log10(abs(H)/max(abs(H)));%归一化的频率响应
pha=unwrap(angle(H));%计算0~1上的相频响应
figure(1)
plot(w/pi,dbH);%画幅频图
xlabel('数字频率域频率\Omega/\pi');ylabel('dB');
axis([0,1,-150,5]);
grid
figure(2)
plot(w/pi,pha);%画相频图
xlabel('数字频率域频率\Omega/\pi');ylabel('Phase\degree'); % axis([0,1,-150,5]);
grid
disp('系统传递函数H(z)');
printsys(bz,az,'z');
figure(3)
h=dimpulse(bz,az);
stem(h)
xlabel('n');ylabel('Impulse response');
grid
图 8直接设计法设计IIR 切比雪夫Ⅰ型带通滤波器幅度响应
图 9直接设计法设计IIR 切比雪夫Ⅰ型带通滤波器相位响应
图 10直接设计法设计IIR 切比雪夫Ⅰ型带通滤波单位冲激响应
得到四个截止频率处的幅度响应如下:
00.10.20.3
0.40.50.60.70.80.91
-150
-100-50
数字频率域频率Ω/π
d B
00.10.20.3
0.40.50.60.70.80.91
-15
-10
-5
5数字频率域频率Ω/π
P h a s e \d e g r e e
n
I m p u l s e r e s p o n s e
满足题目设定的指标要求。
(3)利用MATLAB工具箱中cheb2ord及cheby2设计切比雪夫Ⅱ型滤波器:
核心语句为:
%设计chebshev Ⅱ带通滤波器
clc;clear all
OmegaP1=2*pi*5000;OmegaP2=2*pi*7000;%带通截止频率
OmegaS1=2*pi*3500;OmegaS2=2*pi*8500;%阻带截止频率
Rp=0.5;%通带波纹dB
As=45;%阻带衰减dB
Fs=35*10^3;%抽样频率35khz
OmegaP=[OmegaP1,OmegaP2];OmegaS=[OmegaS1,OmegaS2];
wp=OmegaP/Fs/pi;ws=OmegaS/Fs/pi;%等效数字频率
[N,wc]=cheb2ord(wp,ws,Rp,As)%AF阶数和截至频率
[bz,az]=cheby2(N,As,wc);%AF系统函数的分子分母
w0=[wp,ws];%四个频点
Hx=freqz(bz,az,w0*pi);%计算两个频点上对应的幅度响应
[H,w]=freqz(bz,az);%计算0~1上的幅频响应
dbHx=-20*log10(abs(Hx)/max(abs(H)))%归一化并求dB
dbH=20*log10(abs(H)/max(abs(H)));%归一化的频率响应
pha=unwrap(angle(H));%计算0~1上的相频响应
figure(1)
plot(w/pi,dbH);%画幅频图
xlabel('数字频率域频率\Omega/\pi');ylabel('dB');
axis([0,1,-100,5]);
grid
figure(2)
plot(w/pi,pha);%画相频图
xlabel('数字频率域频率\Omega/\pi');ylabel('Phase\degree');
% axis([0,1,-150,5]);
grid
disp('系统传递函数H(z)');
printsys(bz,az,'z');
figure(3)
h=dimpulse(bz,az); stem(h)
xlabel('n');ylabel('Impulse response'); grid
图 11直接设计法设计IIR 切比雪夫Ⅱ型带通滤波器幅度响应
图 12 直接设计法设计IIR 切比雪夫Ⅱ型带通滤波器相位响应
00.10.20.3
0.40.50.60.70.80.91
-100
-80-60-40-20
0数字频率域频率Ω/π
d B
0.1
0.2
0.3
0.40.50.60.70.8
0.9
1
-6-4-2024
6数字频率域频率Ω/π
P h a s e \d e g r e e
-0.2
-0.1
0.1
0.2
n
I m p u l s e r e s p o n s e
图13直接设计法设计IIR切比雪夫Ⅱ型带通滤波器单位冲激响应
得到四个截止频率处的幅度响应如下:
满足题目设定的指标要求。
(4)利用buttord及butter设计巴特沃思型滤波器。
核心语句为:
%设计chebshev Ⅱ带通滤波器
clc;clear all
OmegaP1=2*pi*5000;OmegaP2=2*pi*7000;%带通截止频率
OmegaS1=2*pi*3500;OmegaS2=2*pi*8500;%阻带截止频率
Rp=0.5;%通带波纹dB
As=45;%阻带衰减dB
Fs=35*10^3;%抽样频率35khz
OmegaP=[OmegaP1,OmegaP2];OmegaS=[OmegaS1,OmegaS2];
wp=OmegaP/Fs/pi;ws=OmegaS/Fs/pi;%等效数字频率
w0=[wp,ws];%四个频点
[N,wc]=buttord(wp,ws,Rp,As)%AF阶数和截至频率
[bz,az]=butter(N,wc);%AF系统函数的分子分母
Hx=freqz(bz,az,w0*pi);%计算两个频点上对应的幅度响应
[H,w]=freqz(bz,az);%计算0~1上的幅频响应
dbHx=-20*log10(abs(Hx)/max(abs(H)))%归一化并求dB
dbH=20*log10(abs(H)/max(abs(H)));%归一化的频率响应
pha=unwrap(angle(H));%计算0~1上的相频响应
figure(1)
plot(w/pi,dbH);%画幅频图
xlabel('数字频率域频率\Omega/\pi');ylabel('dB');
axis([0,1,-200,5]);
grid
figure(2)
plot(w/pi,pha);%画相频图
xlabel('数字频率域频率\Omega/\pi');ylabel('Phase\degree');
% axis([0,1,-150,5]);
grid
disp('系统传递函数H(z)'); printsys(bz,az,'z'); figure(3)
h=dimpulse(bz,az); stem(h)
xlabel('n');ylabel('Impulse response'); grid
图 15直接设计法设计IIR 巴特沃思带通滤波器相位响应
0.1
0.2
0.3
0.40.50.60.70.8
0.9
1
-30-25-20-15-10-5
0数字频率域频率Ω/π
P h a s e \d e g r e e
0.10.20.3
0.40.50.60.70.80.91
数字频率域频率Ω/π
d B
图 14直接设计法设计IIR 巴特沃思带通滤波器幅度响应
图 16直接设计法设计IIR 巴特沃思带通滤波器单位冲激响应
得到四个截止频率处的幅度响应如下:
满足题目设定的指标要求。
20
40
60
80
100
120
140
160
180
n
I m p u l s e r e s p o n s e
FIR 数字滤波器设计作业
周昂 通信工程 201304015027
8.2
设计一个线性相位FIR 数字高通滤波器。
技术指标为p 0.7ωπ=,0.5st ωπ=,
55dB s A =。
画出()h n 及20log |()|(dB)j H e ω及arg |()|(dB)j H e ω曲线。
解:
分析:由于要求阻带最小衰减为55dB, 故考虑用blackman 窗进行设计。
程序如下:
%8.2程序,FIR 高通滤波器 clc,clear all; wp=0.7*pi;ws=0.5*pi; dw=abs(ws-wp);%过渡带宽
N0=ceil(11*pi/dw);%计算blackman 窗长 N=N0+mod(N0+1,2);%确保N 是奇数 n=[0:N-1];
wd=(blackman(N))';%求blakman 函数 wc=(wp+ws)/2/pi;%理想低通滤波器的截止频率 hn=fir1(N-1,wc,'high',blackman(N)); [H,w]=freqz(hn,1,1000);%求频率响应 dbH=20*log10(abs(H)); dw=2*pi/1000;
As=-max(dbH(1:wp/dw+1))%检查通带最大衰减 Rp=-min(dbH(ws/dw+1:501))%检查阻带最小衰减 pha=unwrap(angle(H));%计算0~1上的相频响应 figure(1);stem(n,hn);
title('单位冲激响应(db)');xlabel('n');ylabel('h(n)');grid figure(2)
subplot(2,2,[1 3]) plot(w/pi,dbH);
title('幅度响应(db)');xlabel('\omage/\pi'); ylabel('20log|H(e^j^/omega)|(dB)'); axis([0,1,-200,5]);grid; subplot(2,2,2) plot(w/pi,dbH);
title('幅度响应(db)');xlabel('\omage/\pi'); ylabel('20log|H(e^j^/omega)|(dB)'); axis([0.5,0.9,-10,10]);grid;
set(gca,'xtickmode','manual','xtick',[0.5:0.1:0.9]); set(gca,'ytickmode','manual','ytick',[-10,-3,0,3,10]); subplot(2,2,4)
plot(w/pi,dbH);
title('幅度响应(db)');xlabel('\omage/\pi'); ylabel('20log|H(e^j^/omega)|(dB)'); grid;
axis([0.3 0.6 -100 -50])
set(gca,'xtickmode','manual','xtick',[0.3:0.1:0.6]); set(gca,'ytickmode','manual','ytick',[-100,-80,-55,-50]); figure(3);
plot(w/pi,pha);%画幅频图
xlabel('数字频率域频率\Omega/\pi');ylabel('Phase/degree'); % axis([0,1,-150,5]); grid
运行,得
As =72.3744 Rp = 0.0023
图 17 加blackman 窗的FIR 高通滤波器单位冲激响应
10
20
3040
50
60
n
h (n )
图 18 加blakcman 窗的FIR 高通滤波器幅频响应
图 19加blackman 窗的FIR 高通滤波器相频响应
由运行结果可以看出,阻带最小衰减达到了72.37,满足设计指标。
由幅频图可以看出,高通滤波器在0.5π取最小衰减,已经达到了题目设计的55dB 的指标。
由相频响应可以看出,在通带满足线性相位。
信号得到无失真传输。
0.51
-200
-180-160-140-120
-100
-80-60-40-20
Ω/π
20l o g |H (e j /o m e g a )|(d B )
0.5
0.6
0.70.80.9
-10-3
03
10
Ω/π
20l o g |H (e j /o m e g a )|(d B )
0.3
0.40.50.6
-100
-80
-72-55
-50Ω/π
20l o g |H (e j /o m e g a )|(d B )
00.10.20.30.4
0.50.60.70.80.91
-50
-40-30-20-100
10Ω/π
P h a s e /d e g r e e
8.3
设计一个线性相位FIR 数字带通滤波器。
技术指标为10.4p ωπ=,20.5p ωπ=,
1
0.2st ωπ=,2
0.7st ωπ=,75dB s A =。
画出()h n 及20log |()|(dB)j H e ω及
arg |()|(dB)j H e ω曲线。
解:
分析:由于要求阻带最小衰减为75dB,故考虑用kaiser 窗进行设计。
用kaiser 窗设计时需要根据阻带最小衰减值As 计算出beta 值,这时kaiser 窗与其他窗在使用上的一点区别。
程序如下:
%8.3程序,FIR 带通滤波器 clc,clear all;
wp1=0.4*pi;wp2=0.5*pi; ws1=0.2*pi;ws2=0.7*pi; As0=75;
B=min(wp1-ws1,ws2-wp2);%过渡带宽
N0=ceil((As0-7.95)/(2.286*B))+1;%计算kaiser 窗长 N=N0+mod(N0+1,2);%确保N 是奇数 n=[0:N-1];
beta=0.1102*(As0-8.7);%beta wd=(kaiser(N,beta))';%求kaiser 函数 wc1=(wp1+ws1)/2/pi; wc2=(wp2+ws2)/2/pi; wc=[wc1,wc2];
hn=fir1(N-1,wc,'bandpass',kaiser(N,beta)); [H,w]=freqz(hn,1,1000);%求频率响应 dbH=20*log10(abs(H)); dw=pi/1000;
Rp=-min(dbH(wp1/dw+1:1:wp2/dw+1))%检查通带最大衰减
As=-max(max(dbH(1:1:ws1/dw+1)),max(dbH(ws2/dw+1:1:1000)))%检查
阻带最小衰减
pha=unwrap(angle(H));%计算0~1上的相频响应 figure(1);stem(n,hn);%画冲激响应图
title('单位冲激响应(db)');xlabel('n');ylabel('h(n)');grid figure(2)%画幅频图 subplot(2,2,[1 3]) plot(w/pi,dbH); xlabel('\Omega/\pi');
ylabel('20log|H(e^j^/omega)|(dB)');grid;
subplot(2,2,2) plot(w/pi,dbH); xlabel('\Omega/\pi');
ylabel('20log|H(e^j^/omega)|(dB)');grid; axis([0.3 0.6 -10 10])
set(gca,'xtickmode','manual','xtick',[0.3,0.4,0.5,0.6]); set(gca,'ytickmode','manual','ytick',[-10,-3,0,3,10]); subplot(2,2,4) plot(w/pi,dbH); xlabel('\Omega/\pi');
ylabel('20log|H(e^j^/omega)|(dB)');grid; axis([0.6 0.9 -100 -50])
set(gca,'xtickmode','manual','xtick',[0.6,0.7,0.8,0.9]); set(gca,'ytickmode','manual','ytick',[-100,-75,-50]); figure(3);
plot(w/pi,pha);%画相频图
xlabel('\Omega/\pi');ylabel('Phase/degree');grid
运行,得
As =72.8692 Rp = 0.0021
图 20 加blackman 窗的FIR 带通滤波器单位冲激响应
5101520
253035404550
n
h (n )
图 21 加blackman 窗的FIR 带通滤波器幅频响应
图 22加blackman 窗的FIR 带通滤波器相频响应
可以看出,设定阻带衰减为75dB 时,并不满足要求。
因此将衰减增至80dB,得到返回的阻带最大衰减为75.9251dB ,满足要求。
由相频响应可以看出,在通带满足线性相位。
信号得到无失真传输。
0.51
-160
-140-120-100-80-60-40-20
0Ω/π
20l o g |H (e j /o m e g a )|(d B )
0.3
0.40.50.6
-10-3
03
10
Ω/π
20l o g |H (e j /o m e g a )|(d B )
0.6
0.70.80.9
-100
-75
-50
Ω/π
20l o g |H (e j /o m e g a )|(d B )
0.1
0.2
0.3
0.4
0.50.6
0.7
0.8
0.9
1
-30-25-20-15-10-50
5
10Ω/π
P h a s e /d e g r e e
8.4
设计一个线性相位FIR 数字带阻滤波器。
技术指标为10.35p ωπ=,20.8p ωπ=,
1
0.5st ωπ=,2
0.65st ωπ=,80dB s A =。
画出()h n 及20log |()|(dB)j H e ω及
arg |()|(dB)j H e ω曲线。
解:
分析:由于要求阻带最小衰减为80dB,故考虑用kaiser 窗进行设计。
程序如下:
%8.4程序,FIR 带阻滤波器 %通带截止频率wp 0.35pi 0.8pi
%阻带截止频率ws 0.5pi 0.65pi 最小衰减As80dB clc,clear all;
wp1=0.35*pi;wp2=0.8*pi; ws1=0.5*pi;ws2=0.65*pi; As0=80;
B=min(ws1-wp1,wp2-ws2);%过渡带宽
N0=ceil((As0-7.95)/(2.286*B))+1;%计算kaiser 窗长 N=N0+mod(N0+1,2);%确保N 是奇数 n=[0:N-1];
beta=0.1102*(As0-8.7);%beta wd=(kaiser(N,beta))';%求kaiser 函数 wc1=(wp1+ws1)/2/pi; wc2=(wp2+ws2)/2/pi; wc=[wc1,wc2];
hn=fir1(N-1,wc,'stop',kaiser(N,beta)); [H,w]=freqz(hn,1,1000);%求频率响应 dbH=20*log10(abs(H)); dw=pi/1000;
As=-max(dbH(ws1/dw+1:1:ws2/dw+1))%检查阻带最小衰减
Rp=-min(min(dbH(1:1:wp1/dw+1)),min(dbH(wp2/dw+1:1:1000)))%检查
通带最大波纹
pha=unwrap(angle(H));%计算0~1上的相频响应 figure(1);stem(n,hn);
title('单位冲激响应(db)');xlabel('n');ylabel('h(n)');grid
figure(2)
subplot(2,2,[1 3]) plot(w/pi,dbH);%画幅频图 xlabel('\Omega/\pi');
ylabel('20log|H(e^j^/omega)|(dB)');grid;
set(gca,'ytickmode','manual','ytick',[-200,-150,-100,-75,-50,-2,0,2,50]);
subplot(2,2,2) plot(w/pi,dbH); xlabel('\Omega/\pi');
ylabel('20log|H(e^j^/omega)|(dB)');grid; axis([0.2 0.5 -10 10])
set(gca,'xtickmode','manual','xtick',[0.2,0.3,0.4,0.5]); set(gca,'ytickmode','manual','ytick',[-10,-3,0,3,10]); subplot(2,2,4) plot(w/pi,dbH); xlabel('\Omega/\pi');
ylabel('20log|H(e^j^/omega)|(dB)');grid; axis([0.4 0.7 -100 -50])
set(gca,'xtickmode','manual','xtick',[0.4,0.5,0.6,0.7]); set(gca,'ytickmode','manual','ytick',[-100,-80,-75,-50]); figure(3);
plot(w/pi,pha);%画相频图
xlabel('\Omega/\pi');ylabel('Phase/degree');grid subplot
运行,得
As =73.1121 Rp =0.0013
图 23 加kaiser 窗的FIR 带阻滤波器单位冲激响应
102030
40506070
n
h (n )
图 24 加kaiser 窗的FIR 带阻滤波器幅频响应
图 25加kaiser 窗的FIR 带阻滤波器相频响应
由程序计算返回的结果和幅频图可以看出,设定阻带衰减为75dB 时,在两个阻带截止频率点上衰减并不满足要求。
因此将衰减增至80dB,得到返回的阻带最大衰减为77.0511dB ,满足要求。
由相频响应可以看出,在通带满足线性相位。
信号得到无失真传输。
0.51
-150
-100
-75
-50
-2
2Ω/π
20l o g |H (e j /o m e g a )|(d B )
0.2
0.30.40.5
-10
-3
03
10
Ω/π
20l o g |H (e j /o m e g a )|(d B )
0.4
0.50.60.7
Ω/π
20l o g |H (e j /o m e g a )|(d B )
0.1
0.2
0.3
0.4
0.50.6
0.7
0.8
0.9
1
-80-70-60-50-40
-30-20
-100Ω/π
P h a s e /d e g r e e。