数字信号处理编程一
编程作业
实验:自己构造含高频噪声和低频有用信号的复杂信号,分析信号的频谱特征,据此给定一个低通滤波器的指标,分别借助于Butterworth、切比
雪夫1型、切比雪夫2型、椭圆模拟低通滤波器,用MATLAB设计此
低通数字滤波器并滤波。
实验要求:
1)给出完整的程序;
2)画出幅度谱和对应的增益图(dB);
3)分析不同方式所需的滤波器阶数;
4)画出并分析每种滤波器的滤波前后的结果;
实验过程:
1高频噪声信号和低频有用信号
程序:
x=0:100;
y=sin(0.1*pi*x);
noise=sin(0.5*pi*x);
y1=y+noise;
figure
plot(x,y,'r',x,y1,'b');
title('红色的为原始信号,蓝色为加噪声后的信号');
figure
stem(x/100,abs(fft(y1)));
xlabel('\omega/\pi');ylabel('幅度');
title('加了噪声信号的频谱');
运行结果图:
0102030405060708090100
-2
-1.5-1-0.500.511.5
2红色的为原始信号,蓝色为加噪声后的信号
ω/π
幅度
加了噪声信号的频谱
结果分析:
这次有用信号是y=sin(0.1*pi*x);频率为0.1π;
噪声信号是noise=sin(0.5*pi*x);频率是0.5π;
因此设置边界频率为wp=0.1π,通带的纹波不超过0.5dB,在阻带边
界的频率为ws=0.5π处的最小阻带衰减为15dB;
2、Butterworth滤波器的设计
程序:
x=0:100;
y=sin(0.1*pi*x);
noise=sin(0.5*pi*x);
y1=y+noise;
%自己设置的巴特沃斯滤波器的指标
wp=0.1;
ws=0.5;
rp=0.5;
rs=15;
%构造的巴特沃斯滤波器
[n,wn]=buttord(wp,ws,rp,rs);
[a,b]=butter(n,wn);
[h,omega] = freqz(a,b,256);
gain =20*log10(abs(h));
figure
plot (omega/pi,gain);grid;
xlabel('\omega/\pi'); ylabel('幅度, dB');
title('构造的巴特沃斯滤波器');
lvbo=filtfilt(a,b,y1);
figure
plot(x,y,'r',x,lvbo,'b');
title('红色的为原始信号,蓝色为用巴特沃斯滤波器过滤的信号'); disp('巴特沃斯滤波器所需要的阶数n='); n 程序运行结果:
00.10.20.30.4
0.50.60.70.80.91
-120
-100
-80
-60
-40
-20
ω/π
幅度, d B
构造的巴特沃斯滤波器
实验结果分析:巴特沃斯滤波器所需的阶数是n=3 3、切比雪夫1型滤波器 程序:
x=0:100;
y=sin(0.1*pi*x); noise=sin(0.5*pi*x); y1=y+noise;
%自己设置的切比雪夫1型波器的指标 wp=0.1; ws=0.5;
0102030405060708090100
-1.5
-1
-0.5
0.5
1
1.5
红色的为原始信号,蓝色为用巴特沃斯滤波器过滤的信
号
rp=0.5;
rs=15;
%构造的切比雪夫1型滤波器
[n,wn]=cheb1ord(wp,ws,rp,rs);
[a,b] = cheby1(n,rp,wn);
[h,omega] = freqz(a,b,256);
gain =20*log10(abs(h));
figure
plot (omega/pi,gain);grid;
xlabel('\omega/\pi'); ylabel('幅度, dB');
title('构造的切比雪夫1型滤波器');
lvbo=filtfilt(a,b,y1);
figure
plot(x,y,'r',x,lvbo,'b');
title('红色的为原始信号,蓝色为用切比雪夫1型波器过滤的信号');
disp('切比雪夫1型所需要的阶数n=');n
n
运行结果图:
00.10.20.30.4
0.50.60.70.80.91
-120
-100
-80
-60
-40
-20
ω/π
幅度, d B
构造的切比雪夫1型滤波器
10
20
30
40
50
60
70
80
90
100
-1-0.8-0.6-0.4-0.200.20.40.60.8
1红色的为原始信号,蓝色为用切比雪夫1型波器过滤的信号
运行结果分析:切比雪夫1型滤波器所需的阶数是2阶3、切比雪夫2型滤波器
程序:
x=0:100;
y=sin(0.1*pi*x);
noise=sin(0.5*pi*x);
y1=y+noise;
%自己设置的切比雪夫2型波器的指标
wp=0.1;
ws=0.5;
rp=0.5;
rs=15;
%构造的切比雪夫2型滤波器
[n,wn]=cheb2ord(wp,ws,rp,rs);
[a,b] = cheby2(n,rs,wn);
[h,omega] = freqz(a,b,256);
gain =20*log10(abs(h));
figure
plot (omega/pi,gain);grid;
xlabel('\omega/\pi'); ylabel('幅度, dB');
title('构造的切比雪夫2型滤波器');
lvbo=filtfilt(a,b,y1);
figure
plot(x,y,'r',x,lvbo,'b');
title('红色的为原始信号,蓝色为用切比雪夫2型波器过滤的信号');
disp('切比雪夫2型所需要的阶数n=');n
n
运行结果图:
00.10.20.30.4
0.50.60.70.80.91
-60
-50
-40
-30
-20
-10
ω/π
幅度, d B
构造的切比雪夫2型滤波器
0102030405060708090100
-1.5
-1
-0.5
0.5
1
1.5
红色的为原始信号,蓝色为用切比雪夫2型波器过滤的信号
运行结果分析:切比雪夫2型滤波器所需的阶数是2阶4、椭圆滤波器
程序:
x=0:200;
y=sin(0.1*pi*x);
noise=sin(0.5*pi*x);
y1=y+noise;
%自己设置的椭圆波器的指标
wp=0.1;
ws=0.5;
rp=0.5;
rs=15;
%构造的椭圆滤波器
[n,wn]=ellipord(wp,ws,rp,rs);
[a,b] = ellip(n,rp,rs,wn);
[h,omega] = freqz(a,b,256);
gain =20*log10(abs(h));
figure
plot (omega/pi,gain);grid;
xlabel('\omega/\pi'); ylabel('幅度, dB');
title('构造的椭圆滤波器');
lvbo=filtfilt(a,b,y1);
figure
plot(x,y,'r',x,lvbo,'b');
title('红色的为原始信号,蓝色为用切比雪夫2型波器过滤的信号'); disp('椭圆所需要的阶数n='); n
程序运行图:
00.10.20.30.4
0.50.60.70.80.91
-60
-50
-40
-30
-20
-10
ω/π
幅度, d B
构造的椭圆滤波器
运行结果分析:椭圆滤波器所需要的阶数是2阶 实验结果分析:
各个滤波器滤波所需的阶数
四种滤波器没有显示出太大的差别。 总结:
通过这次的实验编程,我掌握了如何编写巴特沃斯,切比雪夫1,,切比雪夫2型,椭圆滤波器。我进一步强化了自己学习的理论知识。
10
20
30
40
50
60
70
80
90
100
-1-0.8-0.6
-0.4-0.200.20.40.60.81红色的为原始信号,蓝色为用切比雪夫2型波器过滤的信号