哈尔滨工业大学试验方法与数字信号处理大作业精编版.docx
《哈尔滨工业大学试验方法与数字信号处理大作业精编版.docx》由会员分享，可在线阅读，更多相关《哈尔滨工业大学试验方法与数字信号处理大作业精编版.docx（16页珍藏版）》请在冰点文库上搜索。
哈尔滨工业大学试验方法与数字信号处理大作业精编版
HarbinInstituteofTechnology
大作业一
课程名称：
试验方法与数字信号处理
院系：
机械电子
班级：
15S0825
学号：
姓名：
哈尔滨工业大学
给出信号
1.绘出信号波形。
利用matla软件，绘制出的原信号波形如图1所示。
图1原波形信号
2.低通滤波，分别用FIR，IIR滤波器，保留10Hz，去除80Hz和200Hz，并画出波形，并与10Hz信号对比。
解：
原信号的最大Fmax=200Hz，取：
此时，满足采样定理。
（1）、用FIR滤波器（附录1）
选择低通滤波的截止频率为50Hz，滤波器项数为80，通过FIR滤波器公式，可得到滤波后的信号。
编写matlab程序，对比滤波后信号和10Hz信号，如图2所示。
图2FIR滤波后信号与10Hz信号对比
通过图2可以发现，滤波后的信号大致反应了10Hz信号的变化，相位一致，幅值衰减了一部分，说明滤波后，确实去除了80Hz，200Hz的信号。
为了进一步说明问题，绘制滤波后信号的频谱图，如图3所示。
从图3可以看出，随着N的增大，10Hz信号幅值衰减的程度变小，会趋于至原幅值的一半，其余信号幅值衰减的程度变大，滤波效果更加明显。
图3FIR滤波后频谱（N=8，30，80，800）10Hz
尝试用汉宁窗口对泄漏进行修正，修正前后的波形如图4所示。
图4采用汉宁窗口修正
（2）、用IIR滤波器（附录2）
选择低通滤波的截止频率为50Hz的二阶IIR滤波器，根据相关公式，可以得到IIR滤波器的滤波因子，进而可得到滤波后的信号。
编写matlab程序，对比滤波后信号和10Hz信号，如图5所示。
图5IIR滤波后信号与10Hz信号对比
通过图5可以发现，滤波后的信号大致反应了10Hz信号的变化，相位一致，幅值衰减了一部分，说明滤波后，确实去除了80Hz，200Hz的信号。
在滤波信号开始阶段，会出现一较大的波动，该波动会随滤波的进行而消失。
为了便于说明问题，绘制出滤波后信号的频谱，如图6所示。
从图6可以看出，滤波后的信号幅值基本与原幅值一样，且高频信号衰减幅度比较大，滤波效果比FIR滤波效果好。
图6IIR滤波后频谱
3、带通滤波，分别用FIR，IIR滤波器，保留80Hz，去除10Hz和200Hz，并画出波形，并与10Hz信号对比。
解：
原信号的最大Fmax=200Hz，取：
满足采样定理。
（1）、用FIR滤波器（附录3）
选择带通频率为40~120Hz，即F1=40Hz，F2=120Hz，滤波器项数为80，根据公式，可得相应的滤波因子，编写相应的程序，可得到滤波后的信号，如图7所示。
图7FIR滤波后信号与80Hz信号对比
通过图7可以发现，滤波后的信号大致反应了80Hz信号的变化，相位一致，幅值衰减了一部分，说明滤波后，确实去除了10Hz，200Hz的信号。
为了进一步说明问题，绘制滤波后信号的频谱图，如图8所示。
从图8可以看出，随着N的增大，80Hz幅值衰减的程度变小，会趋于至原幅值的一半，10Hz和200Hz信号幅值衰减程度变大，滤波效果更加明显。
图8FIR滤波后频谱（N=8，30，80，800）80Hz
（2）用FIR滤波（附录4）
选择带通频率为40~120Hz，即F1=40Hz，F2=120Hz，，根据公式，可得相应的滤波因子，编写相应的程序，可得到滤波后的信号，如图9所示。
图9IIR滤波后信号与80Hz信号对比
通过图9可以发现，滤波后的信号大致反应了80Hz信号的变化，相位一致，幅值衰减了一部分，说明滤波后，确实去除了10Hz，200Hz的信号。
在滤波信号开始阶段，会出现一较大的波动，该波动会随滤波的进行而消失。
为了便于说明问题，绘制出滤波后信号的频谱，如图10所示。
从图10可以看出，滤波后的信号幅值基本与原幅值一样，10Hz信号和200Hz信号的幅值衰减较大，滤波效果比FIR滤波效果好。
图10IIR滤波后频谱
（4）、原信号波形加5%的白噪声信号，进行滤波（附录5）
解：
利用matlab的awgn函数，对原信号添加50%的白噪声，命令如下：
y=awgn（x,SNR）在信号x中加入高斯白噪声。
信噪比SNR，本例中，SNR=2。
加入白噪声信号之后的信号波形如图11所示。
图11添加白噪声信号之后的信号波形
采用低通IIR滤波器，滤去80Hz,200Hz信号，保留10Hz信号，滤波后信号如图12所示。
图12加白噪声之后滤波信号与10Hz信号对比
为了便于分析，绘制滤波后的频谱，如图13所示。
图13加入白噪声滤波之后频谱
将该频谱与未加白噪声的滤波之后的信号的频谱（图6）对比可以发现，加入白噪声之后，滤波之后的信号同样被白噪声影响，并未滤去白噪声信号。
附录1
%采用FIR滤波器低通滤波器
%滤波效果和N，F有关
clc;clear;
Dt=0.0001;
t=0:
Dt:
0.5;
xt=@（t）sin（2*pi*10*t）+sin（2*pi*80*t）+sin（2*pi*200*t）;
F=50;%低通滤波的频率；
N=80;%滤波器项数；
fi_fir=sin（2*pi*F*（1:
N）*Dt）./（pi*（1:
N））;%滤波因子
f0_fir=2*F*Dt;
f_fir=[f0_firfi_fir];%得到的滤波因子序列
fork=1:
length（t）
k_t=Dt*（（k-N）:
k）;
x_k_t=xt（k_t）;
w=conv（f_fir,x_k_t）;
y（k）=w（length（f_fir））;
end
figure;
plot（t,y,'r'）;
holdon;
plot（t,sin（2*pi*10*t））;
title（'滤波后信号与10Hz信号对比'）;
xlabel（'时间t'）;
ylabel（'xt'）;
legend（'滤波后','y=sin（2*pi*10*t'）;
%%采用汉宁窗口对泄漏进行修正
holdon;
fi_hanning=0.5*fi_fir.*（1+cos（pi*（1:
N）/N））;
f_hanning=[f0_firfi_hanning];
fork=1:
length（t）
k_t=Dt*（（k-N）:
k）;
x_k_t=xt（k_t）;
w=conv（f_hanning,x_k_t）;
y_hanning（k）=w（length（f_hanning））;
end
figure;
holdon
plot（t,y,'b-',t,y_hanning,'g'）
title（'采用汉宁窗口修正对比'）;
xlabel（'时间t'）;
ylabel（'xt'）;
legend（'未修正','修正后'）;
%%频谱分析幅值频谱
subplot（4,1,4）;
N=length（t）;
Y=fft（y,N）/N*2;
ff=1/Dt/N*（0:
1:
N-1）;
plot（ff（1:
N/20）,abs（Y（1:
N/20）））;
title（'滤波后频谱N=800'）
xlabel（'频率（Hz）'）
ylabel（'H（f）'）;
附录2
%采用二阶IIR滤波器低通滤波器
clc;clear;
%绘制信号波形
Dt=1/1000;
t=0:
Dt:
0.5;
xt=sin（2*pi*10*t）+sin（2*pi*80*t）+sin（2*pi*200*t）;
F=50;%低通滤波的频率；
omega=tan（pi*F*Dt）;
f0=omega^2/（1+sqrt
（2）*omega+omega^2）;
f1=2*omega^2/（1+sqrt
（2）*omega+omega^2）;
f2=omega^2/（1+sqrt
（2）*omega+omega^2）;
g1=-2*（1-omega^2）/（1+sqrt
（2）*omega+omega^2）;
g2=（1-sqrt
（2）*omega+omega^2）/（1+sqrt
（2）*omega+omega^2）;
y
（1）=0;
y
（2）=xt
（2）;
fork=3:
length（t）
%x_k=xt（k）;x_k_1=xt（k-1）;x_k_2=xt（k-2）;
y（k）=f0*xt（k）+f1*xt（k-1）+f2*xt（k-2）-g1*y（k-1）-g2*y（k-2）;
end
plot（t,y）
holdon;
plot（t,sin（2*pi*10*t））;
title（'滤波后信号与10Hz信号对比'）;
xlabel（'时间t'）;
ylabel（'xt'）;
legend（'滤波后','y=sin（2*pi*10*t'）;
%%频谱分析幅值频谱
N=length（t）;
Y=fft（y,N）/N*2;
ff=1/Dt/N*（0:
1:
N-1）;
plot（ff（1:
N/2）,abs（Y（1:
N/2）））;
title（'滤波后频谱'）
xlabel（'频率（Hz）'）
ylabel（'H（f）'）;
附录3
%fir滤波器带通
clc;clear;
Dt=0.0001;
t=0:
Dt:
0.1;
xt=@（t）sin（2*pi*10*t）+sin（2*pi*80*t）+sin（2*pi*200*t）;
F1=40;F2=120;
N=800;%滤波器项数；
f0=2*Dt*（F2-F1）;
fi=2./（pi.*（1:
N））.*sin（pi*（F2-F1）.*（1:
N）*Dt）.*cos（pi*（F2+F1）.*（1:
N）*Dt）;
f=[f0fi];
fork=1:
length（t）
k_t=Dt*（（k-N）:
k）;
x_k_t=xt（k_t）;
w=conv（f,x_k_t）;
y（k）=w（length（f））;
end
plot（t,y）
holdon;
plot（t,sin（2*pi*80*t））;
title（'滤波后信号与80Hz信号对比'）;
xlabel（'时间t'）;
ylabel（'xt'）;
legend（'滤波后','y=sin（2*pi*80*t'）;
subplot（4,1,4）;
N=length（t）;
Y=fft（y,N）/N*2;
ff=1/Dt/N*（0:
1:
N-1）;
plot（ff（1:
N/20）,abs（Y（1:
N/20）））;
title（'滤波后频谱N=800'）
xlabel（'频率（Hz）'）
ylabel（'H（f）'）;
附录4
%iir滤波器带通
clc;clear;
%绘制信号波形
Dt=1/1000;
t=0:
Dt:
0.2;
xt=sin（2*pi*10*t）+sin（2*pi*80*t）+sin（2*pi*200*t）;
F1=40;F2=120;
omega=tan（pi*（F2-F1）*Dt）;
beta=cos（pi*（F2+F1）*Dt）/cos（pi*（F2-F1）*Dt）;
K=1+sqrt
（2）*omega+omega^2;
f0=omega^2/K;f1=0;f2=-2*f0;f3=0;f4=f0;
g1=-2*beta*（2+sqrt
（2）*omega）/K;
g2=2*（1+2*beta^2-omega^2）/K;
g3=-2*beta*（2-sqrt
（2）*omega）/K;
g4=（1-sqrt
（2）*omega+omega^2）/K;
y
（1）=xt
（1）;
y
（2）=xt
（2）;
y（3）=xt（3）;
y（4）=xt（4）;
fork=5:
length（t）
y（k）=f0*xt（k）+f1*xt（k-1）+f2*xt（k-2）+f3*xt（k-3）+f4*xt（k-4）-...
g1*y（k-1）-g2*y（k-2）-g3*y（k-3）-g4*y（k-4）;
end
plot（t,y）
holdon;
plot（t,sin（2*pi*80*t））;
title（'滤波后信号与80Hz信号对比'）;
xlabel（'时间t'）;
ylabel（'xt'）;
legend（'滤波后','y=sin（2*pi*80*t'）;
N=length（t）;
Y=fft（y,N）/N*2;
ff=1/Dt/N*（0:
1:
N-1）;
plot（ff（1:
N/2）,abs（Y（1:
N/2）））;
title（'滤波后频谱'）
xlabel（'频率（Hz）'）
ylabel（'H（f）'）;
附录5
%添加高斯白噪声，信号比为2
%采用二阶IIR滤波器
clc;clear;
%绘制信号波形
Dt=1/1000;
t=0:
Dt:
0.2;
xt=sin（2*pi*10*t）+sin（2*pi*80*t）+sin（2*pi*200*t）;
%添加白噪声
yt=awgn（xt,2）;
F=30;%低通滤波的频率；
omega=tan（pi*F*Dt）;
f0=omega^2/（1+sqrt
（2）*omega+omega^2）;
f1=2*omega^2/（1+sqrt
（2）*omega+omega^2）;
f2=omega^2/（1+sqrt
（2）*omega+omega^2）;
g1=-2*（1-omega^2）/（1+sqrt
（2）*omega+omega^2）;
g2=（1-sqrt
（2）*omega+omega^2）/（1+sqrt
（2）*omega+omega^2）;
%y
（1）,y
（2）如何赋值？
？
？
y
（1）=0;
%y
（2）=sin（2*pi*10*Dt）;
y
（2）=yt
（2）;
fork=3:
length（t）
%x_k=xt（k）;x_k_1=xt（k-1）;x_k_2=xt（k-2）;
y（k）=f0*yt（k）+f1*yt（k-1）+f2*yt（k-2）-g1*y（k-1）-g2*y（k-2）;
end
plot（t,y）
holdon;
plot（t,sin（2*pi*10*t））;
title（'滤波后信号与10Hz信号对比'）;
xlabel（'时间t'）;
ylabel（'xt'）;
legend（'滤波后','y=sin（2*pi*10*t'）;
N=length（t）;
Y=fft（y,N）/N*2;
ff=1/Dt/N*（0:
1:
N-1）;
plot（ff（1:
N/2）,abs（Y（1:
N/2）））;
title（'滤波后频谱'）
xlabel（'频率（Hz）'）
ylabel（'H（f）'）;
升级会员 