数字信号处理MATLAB版上机实验操作.docx
《数字信号处理MATLAB版上机实验操作.docx》由会员分享，可在线阅读，更多相关《数字信号处理MATLAB版上机实验操作.docx（22页珍藏版）》请在冰豆网上搜索。
数字信号处理MATLAB版上机实验操作
实验一离散时间信号与系统
一、实验目的：
1、熟悉常见离散时间信号的产生方法；
2、熟悉离散时间系统的单位脉冲响应和单位阶跃响应的求解方法；
3、熟悉离散时间信号经过离散时间系统的响应的求解方法。
二、实验内容：
已知离散时间系统差分方程为y（n）-0.5y（n-1）+0.06y（n-2）=x（n）+x（n-1），求
1、该系统的单位脉冲响应并绘图；
2、该系统的单位阶跃响应并绘图；
3、已知x（n）=可自己指定用filter函数经过系统的响应并绘图；
4、用conv_m函数求系统响应并绘图。
三、实验平台：
MATLAB集成系统
四、设计流程：
此处写个人自己的设计流程
五、程序清单：
此处写程序内容
六、调试和测试结果：
此处写程序的执行结果和实验过程中的调试经过、出现的错误和对应的解决方法
七、教师评语与成绩评定
此处由老师填写
上机操作：
实验一离散时间信号与系统
实验内容：
1．脉冲响应
>>b=[1,1];a=[1,-0.5,0.06];n=[-10:
25];
>>impz（b,a,n）;
>>title（'ImpulseResponse'）;xlabel（'n'）;ylabel（'h（n）'）
2.单位阶跃响应
>>x=stepseq（0,-10,25）;s=filter（b,a,x）;
Warning:
Functioncallstepseqinvokesinexactmatchd:
\MATLAB7\work\STEPSEQ.M.
>>stem（n,s）
>>title（'StepResponse'）;xlabel（'n'）;ylabel（'s（n）'）
3.>>a=[1,-0.5,0.06];b=[1,1];
>>n=-20:
120;
>>x1=exp（-0.05*n）.*sin（0.1*pi*n+pi/3）;
>>s1=filter（b,a,x1）;
>>stem（n,s1）;;xlabel（'n'）;ylabel（'s1（n）'）;
4．>>a=[1,-0.5,0.06];b=[1,1];
>>n=-20:
120;
>>h=impz（b,a,n）;
>>x1=exp（-0.05*n）.*sin（0.1*pi*n+pi/3）;
>>[y,m]=conv_m（x1,n,h,n）;
Warning:
Functioncallconv_minvokesinexactmatchd:
\MATLAB7\work\CONV_M.M.
>>stem（m,y）;title（'系统响应'）;xlabel（'m'）;ylabel（'y（m）'）;
实验二离散信号与系统的连续频域分析
一、实验目的：
1、掌握离散时间信号的DTFT的MATLAB实现；
2、掌握离散时间系统的DTFT分析；
3、掌握系统函数和频率相应之间的关系。
二、实验内容：
1、自定义一个长度为8点的信号，信号幅度值也由自己任意指定，对该信号作DTFT，分别画出幅度谱和相位谱
2、已知离散时间系统差分方程为y（n）-0.5y（n-1）+0.06y（n-2）=x（n）+x（n-1），求出并画出其频率响应
3、求该系统系统函数，并画极零点图，并通过freqz函数求频率响应。
三、实验平台：
MATLAB集成系统
四、设计流程：
此处写个人自己的设计流程
五、程序清单：
此处写程序内容
六、调试和测试结果：
此处写程序的执行结果和实验过程中的调试经过、出现的错误和对应的解决方法
七、教师评语与成绩评定
此处由老师填写
上机操作
1
>>n=0:
7;
>>x=（0.9*exp（j*pi/3））.^n;
>>w=-pi:
pi/200:
pi;
>>X=dtft（x,n,w）;
>>magX=abs（X）;angX=angle（X）;
>>subplot（2,1,1）;plot（w/pi,magX）;xlabel（'w/pi'）;ylabel（'幅度|X|'）;
>>subplot（2,1,2）;plot（w/pi,angX）;xlabel（'w/pi'）;ylabel（'相位（rad/π）'）;
2
>>a=[1,-0.5,0.06];b=[1,1,0];
>>m=0:
length（b）-1;l=0:
length（a）-1;
>>w=0:
pi/500:
pi;
>>num=b*exp（-j*m'*w）;
>>den=a*exp（-j*l'*w）;
>>H=num./den;
>>magH=abs（H）;angH=angle（H）;
>>H1=freqz（b,a,w）;
>>magH1=abs（H1）;angH1=angle（H1）;
>>subplot（2,2,2）;plot（w/pi,angH/pi）;grid;xlabel（'w（frequencyinpiunits）'）;ylabel（'相位（rad/π）'）;
>>subplot（2,2,1）;plot（w/pi,magH）;grid;xlabel（'w（frequencyinpiunits）'）;ylabel（'幅度|H|'）;
>>subplot（2,2,3）;plot（w/pi,magH1）;grid;xlabel（'w（frequencyinpiunits）'）;ylabel（'幅度|H1|'）;
>>subplot（2,2,4）;plot（w/pi,angH1/pi）;grid;xlabel（'w（frequencyinpiunits）'）;ylabel（'相位（rad/π）'）;axis（[0,1,-0.8,0]）;
>>figure
（2）;zplane（b,a）;
实验3离散信号与系统的离散频域分析（DFT）
一、实验目的：
1、掌握离散时间系统的DFT的MATLAB实现；
2、熟悉DTFT和DFT之间的关系。
3、了解信号不同变形的DFT与原信号DFT之间的关系
二、实验内容：
自定义一个长为8点的信号；
1、对信号分别做8点、16点、32点DFT，分别与DTFT合并作图并比较DFT与DTFT之间的关系。
2、在信号后加零扩展至八点，每两个相邻样本之间插入一个零值，扩充为16点序列，作DFT，画出幅度谱和相位谱，并与原序列的DFT进行比较。
3、将信号以8为周期扩展，得到长为16的两个周期，作DFT，画出幅度谱和相位谱，并与原序列的DFT进行比较。
三、实验平台：
MATLAB集成系统
四、设计流程：
此处写个人自己的设计流程
五、程序清单：
此处写程序内容
六、调试和测试结果：
此处写程序的执行结果和实验过程中的调试经过、出现的错误和对应的解决方法
七、教师评语与成绩评定
此处由老师填写
上机操作
1.
>>n=0:
7;
>>x=（0.9*exp（j*pi/3））.^n;
>>x1=[xzeros（1,8）];x2=[x1zeros（1,16）];
>>Xk=dft（x,8）;Xk1=dft（x1,16）;Xk2=dft（x2,32）;
>>w=0:
pi/200:
2*pi;X=dtft（x,n,w）;
>>magX=abs（X）;angX=angle（X）;
>>magXk=abs（Xk）;angXk=angle（Xk）;magXk1=abs（Xk1）;angXk1=angle（Xk1）;
>>magXk2=abs（Xk2）;angXk2=angle（Xk2）;
>>subplot（4,2,1）;plot（w/pi,magX）;xlabel（'w/pi'）;ylabel（'幅度|X|'）;gridon;
>>subplot（4,2,2）;plot（w/pi,angX）;xlabel（'w/pi'）;ylabel（'相位（rad/π）'）;
>>subplot（4,2,3）;stem（0:
7,magXk）;xlabel（'k'）;ylabel（'幅度|X（k）|'）;axis（[0,8,0,6]）;
>>subplot（4,2,4）;stem（0:
7,angXk）;xlabel（'k'）;ylabel（'相位（rad/π）'）;axis（[0,8,-2,2]）;
>>subplot（4,2,5）;stem（0:
15,magXk1）;xlabel（'k'）;ylabel（'幅度|X1（k）|'）;axis（[0,16,0,6]）;
>>subplot（4,2,6）;stem（0:
15,angXk1）;xlabel（'k'）;ylabel（'相位（rad/π）'）;axis（[0,16,-2,2]）;
>>subplot（4,2,7）;stem（0:
31,magXk2）;xlabel（'k'）;ylabel（'幅度|X2（k）|'）;axis（[0,32,0,6]）;
>>subplot（4,2,8）;stem（0:
31,angXk2）;xlabel（'k'）;ylabel（'相位（rad/π）'）;axis（[0,32,-2,2]）;
2.
>>n=0:
7;
>>x=（0.9*exp（j*pi/3））.^n;
>>fori=1:
8%将原序列每两个相邻样本之间插入一个零值，扩展为长为16的序列
y1（2*i-1:
2*i）=[x（i）0];
end
>>Yk1=dft（y1,16）;
Warning:
FunctioncalldftinvokesinexactmatchD:
\MATLAB7\work\DFT.M.
>>magYk1=abs（Yk1）;angYk1=angle（Yk1）;
>>y2=[xx];%将信号以8为周期扩展，得到长为16的两个周期
>>Yk2=dft（y2,16）;
>>magYk2=abs（Yk2）;angYk2=angle（Yk2）;
>>Xk=dft（x,8）;
>>magXk=abs（Xk）;angXk=angle（Xk）;
>>subplot（3,2,1）;stem（0:
7,magXk）;xlabel（'k'）;ylabel（'幅度|X（k）|'）;
>>subplot（3,2,2）;stem（0:
7,angXk）;xlabel（'k'）;ylabel（'相位（rad/π）'）;
>>subplot（3,2,3）;stem（0:
15,magYk1）;xlabel（'k'）;ylabel（'幅度|Y1（k）|'）;
>>subplot（3,2,4）;stem（0:
15,angYk1）;xlabel（'k'）;ylabel（'相位（rad/π）'）;
>>subplot（3,2,5）;stem（0:
15,magYk2）;xlabel（'k'）;ylabel（'幅度|Y2（k）|'）;
>>subplot（3,2,6）;stem（0:
15,angYk2）;xlabel（'k'）;ylabel（'相位（rad/π）'）;
实验四IIR滤波器设计
一、实验目的：
1、掌握各种模拟原型滤波器的滤波特性；
2、掌握模数滤波器变换时的脉冲响应不变法和双线性变换法；
3、掌握低通滤波器变换成其他类型滤波器的方法；
4、能够根据指标选择合适的原型滤波器和合适的方法设计IIR滤波器。
二、实验内容：
1、自定设计指标（通带截止频率、通带波纹、阻带截止频率、阻带衰减），选择合适的模拟原型低通滤波器和合适的设计方法（脉冲响应不变法或双线性变换法），设计符合要求的IIR滤波器；并验证设计好的滤波器是否满足要求。
2、将此低通滤波器映射为高通、带通或带阻滤波器，并验证变换后的滤波器是否满足指标。
3、求输入x（n）=[cos（w1n）+cos（w2n）]u（n）经过系统后的输出y（n）。
其中w1位于通带内，w2位于阻带内，要求做一个两行两列的子图，第一个做x的时域波形，第二个做x的幅度谱，第三个做y的时域波形，第四个做y的幅度谱。
三、实验平台：
MATLAB集成系统
四、设计流程：
此处写个人自己的设计流程
五、程序清单：
此处写程序内容
六、调试和测试结果：
此处写程序的执行结果和实验过程中的调试经过、出现的错误和对应的解决方法
七、教师评语与成绩评定
此处由老师填写
上机操作
窗口设计法：
>>ws1=0.2*pi;wp1=0.35*pi;wp2=0.65*pi;ws2=0.8*pi;As=60;
>>tr_width=min（（wp1-ws1）,（ws2-wp2））;M=ceil（11*pi/tr_width）+1
M=
75
>>n=[0:
1:
M-1];wc1=（ws1+wp1）/2;wc2=（ws2+wp2）/2;
>>hd=ideal_lp（wc2,M）-ideal_lp（wc1,M）;
Warning:
Functioncallideal_lpinvokesinexactmatchD:
\MATLAB7\work\IDEAL_LP.M.
>>w_bla=（blackman（M））';h=hd.*w_bla;
Warning:
FunctioncallblackmaninvokesD:
\MATLAB7\toolbox\signal\signal\blackman.m,
however,functionD:
\MATLAB7\work\BLACKMAN.M
thatdiffersonlyincaseprecedesitonthepath.
>>[db,mag,pha,grd,w]=freqz_m（h,[1]）;delta_w=2*pi/1000;
Warning:
Functioncallfreqz_minvokesinexactmatchD:
\MATLAB7\work\FREQZ_M.M.
>>Rp=-min（db（wp1/delta_w+1:
1:
wp2/delta_w））
Rp=
0.0030
>>As=-round（max（db（ws2/delta_w+1:
1:
501）））
As=
75
>>subplot（2,2,1）;stem（n,hd）;xlabel（'n'）;ylabel（'hd（n）'）;axis（[0,M-1,-0.4,0.5]）;
>>subplot（2,2,2）;stem（n,w_bla）;xlabel（'n'）;ylabel（'w（n）'）;axis（[0,M-1,0,1.1]）;
>>subplot（2,2,3）;stem（n,h）;xlabel（'n'）;ylabel（'h（n）'）;axis（[0,M-1,-0.4,0.5）;
>>subplot（2,2,4）;plot（w/pi,db）;xlabel（'w'）;ylabel（'db'）;axis（[0,1,-139,10）;
接上面的程序：
>>w1=0.5*pi;w2=0.1*pi;
>>x=cos（w1*n）+cos（w2*n）;
>>[y,ny]=conv_m（x,n,h,n）;
>>X=dtft（x,n,w）;
>>Y=dtft（y,ny,w）;
>>magX=abs（X）;magY=abs（Y）;figure
（2）;
>>subplot（2,2,1）;stem（n,x）;xlabel（'n'）;ylabel（'x（n）'）;axis（[0,M-1,-2,2]）;
>>subplot（2,2,2）;plot（w/pi,magX）;xlabel（'w'）;ylabel（'幅度|X|'）;axis（[0,1,0,40]）;
>>subplot（2,2,3）;stem（ny,y）;xlabel（'n'）;ylabel（'y（n）'）;
>>subplot（2,2,4）;plot（w/pi,magY）;xlabel（'w'）;ylabel（'幅度|Y|'）;axis（[0,1,0,40]）;
频率采样法:
>>M=40;alpha=（M-1）/2;l=0:
M-1;wl=（2*pi/M）*l;
>>T1=0.109021;T2=0.59417456;
>>Hrs=[zeros（1,5）,T1,T2,ones（1,7）,T2,T1,zeros（1,9）,T1,T2,ones（1,7）,T2,T1,zeros（1,4）];
>>Hdr=[0,0,1,1,0,0];wdl=[0,0.2,0.35,0.65,0.8,1];
>>k1=0:
floor（（M-1）/2）;k2=floor（（M-1）/2）+1:
M-1;
>>angH=[-alpha*（2*pi）/M*k1,alpha*（2*pi）/M*（M-k2）];
>>H=Hrs.*exp（j*angH）;h=real（ifft（H,M））;
>>[db,mag,pha,grd,w]=freqz_m（h,1）;[Hr,ww,a,L]=Hr_Type2（h）;
Warning:
FunctioncallHr_Type2invokesinexactmatchD:
\MATLAB7\work\HR_TYPE2.M.
>>subplot（2,2,1）;plot（wl（1:
21）/pi,Hrs（1:
21）,'o',wdl,Hdr）
>>axis（[0,1,-0.1,1.1]）;xlabel（'w'）;ylabel（'Hr（k）'）
>>subplot（2,2,2）;stem（l,h）;xlabel（'n'）;ylabel（'h（n）'）;axis（[0,M,-0.3,0.3]）
>>subplot（2,2,3）;plot（ww/pi,Hr,wl（1:
21）/pi,Hrs（1:
21）,'o'）
>>axis（[0,1,-0.1,1.1]）;xlabel（'w'）;ylabel（'Hr（w）'）
>>subplot（2,2,4）;plot（w/pi,db）;axis（[0,1,-100,10]）;xlabel（'w'）;ylabel（'db'）
接上面的程序：
>>n=[0:
1:
M-1];
>>w1=0.5*pi;w2=0.1*pi;
>>x=cos（w1*n）+cos（w2*n）;
>>[y,ny]=conv_m（x,n,h,n）;
>>X=dtft（x,n,w）;
>>Y=dtft（y,ny,w）;
>>magX=abs（X）;magY=abs（Y）;figure
（2）;
>>subplot（2,2,1）;stem（n,x）;xlabel（'n'）;ylabel（'x（n）'）;axis（[0,M-1,-2,2]）;
>>subplot（2,2,2）;plot（w/pi,magX）;xlabel（'w'）;ylabel（'幅度|X|'）;axis（[0,1,0,21]）;
>>subplot（2,2,3）;stem（ny,y）;xlabel（'n'）;ylabel（'y（n）'）;
>>subplot（2,2,4）;plot（w/pi,magY）;xlabel（'w'）;ylabel（'幅度|Y|'）;axis（[0,1,0,21]）;
实验五FIR滤波器设计方法
一、实验目的：
1、掌握不同窗函数的性能指标；
2、根据指标选择窗函数设计FIR滤波器。
3、掌握频率采样法设计FIR滤波器
二、实验内容：
自定滤波器的类型（低通，高通或带通）、设计指标（通带截止频率、通带波纹、阻带截止频率、阻带衰减）
1、根据指标选择合适的窗函数，用窗口设计法设计符合指标的FIR滤波器；并验证其性能是否满足预定指标。
参考书上例题
2、根据指标选择合适的采样点数，用频率采样法设计符合指标的FIR滤波器；并验证其性能是否满足预定指标。
参考书上例题
求出通带内的Rp和阻带内的As，或者用两行两列的子图分别将结果列出来。
3、求输入x（n）=[cos（w1n）+cos（w2n）]u（n）经过系统后的输出y（n）。
其中w1位于通带内，w2位于阻带内，要求做一个两行两列的子图，第一个做x的时域波形，第二个做x的幅度谱，第三个做y的时域波形，第四个做y的幅度谱。
三、实验平台：
MATLAB集成系统
四、设计流程：
此处写个人自己的设计流程
五、程序清单：
此处写程序内容
六、调试和测试结果：
此处写程序的执行结果和实验过程中的调试经过、出现的错误和对应的解决方法
参考结果
设计指标：
wp=0.3*pi;ws=0.4*pi;Rp=1;As=40;
要点：
设计指标要与结果图一致
七、教师评语与成绩评定
此处由老师填写
上机操作
模拟原型Butterworth滤波器：
>>wp=0.4*pi;ws=0.6*pi;Rp=0.5;As=50;
>>[b,a]=afd_butt（wp,ws,Rp,As）;
***ButterworthFilterOrder=17
>>[db,mag,pha,w]=freqs_m（b,a,pi）;
>>[ha,x,t]=impulse（b,a,100）;
>>subplot（2,2,1）;plot（w/pi,mag）;xlabel（'w'）;ylabel（'|H|'）;gridon
>>subplot（2,2,2）;plot（w/pi,pha）;xlabel（'w'）;ylabel（'相位'）;gridon
>>subplot（2,2,3）;plot（w/pi,db）;xlabel（'w'）;ylabel（'db'）;gridon
>>subplot（2,2,4）;plot（t,ha）;xlabel（'t'）;ylabel（'ha（t）'）;axis（[0,100,-0.2,0.4]）;gridon
双线性变换法求得的数字Butterworth低通滤波器：
>>wp=0.4*pi;ws=0.6*pi;Rp=0.5;As=50;
>>T=2;OmegaP=（2/T）*tan（wp/2）;OmegaS=（2/T）*tan（ws/2）;
>>N=ceil（（log10（（10^（Rp/10）-1）/（10^（As/10）-1）））/（2*log10（OmegaP/OmegaS）））;
>>OmegaC=OmegaP/（（10^（Rp/10）-1）^（1/（2*N）））;wn=2/pi*atan（OmegaC*T/2）;
>>[b,a]=butter（N,wn）;
>>[db,mag,pha,grd,w]=freqz_m（b,a）;
>>n=0:
100;
>>h=impz（b,a,n）;
>>subplot（2,2,1）;plot（w/pi,mag）;xlabel（'w'）;ylabel（'|H|'）;axis（[0,1,0,1.15]）;gridon
>>subplot（2,2,2）;plot（w/pi,pha）;xlabel（'w'）;ylabel（'相位'）;axis（[0,1,-4,4]）;gridon
>>subplot（2,2,3）;plot（w/pi,db）;xlab
升级会员 