信息处理方向课程设计.docx

江南大学 物联网工程学院 信息处理方向课程设计实验报告班级: 姓名: 学号: 一、课程设计目的及要求数字信号处理是一门理论性和实践性都很强的学科，通过课程设计可以加深理解掌握基本理论，培养学生分析问题和解决问题的综合能力，为将来走向工作岗位奠定坚实的基础，因此做好课程设计是学好本课程的重要教学辅助环节。本指导书结合教材数字信号处理教程的内容，基于MATLAB程序语言提出课程设计的题目及要求，在做课程设计之前要求学生要尽快熟悉MATLAB语言，充分预习相关理论知识，独立编写程序，以便顺利完成课程设计。二、课程设计任务课程设计的过程是综合运用所学知识的过程。课程设计主要任务是围绕数字信号的频谱分析、特征提取和数字滤波器的设计来安排的。根据设计题目的具体要求，运用MATLAB语言完成题目所规定的任务及功能。设计任务包括：查阅专业资料、工具书或参考文献，了解设计课题的原理及算法、编写程序并在计算机上调试，最后写出完整、规范的课程设计报告书。课程设计地点在信息学院机房，一人一机，在教师统一安排下独立完成规定的设计任务。三、课程设计题目根据大纲要求提供以下四个课程设计题目供学生选择，根据实际情况也可做其它相关课题。1. DFT在信号频谱分析中的应用2. 数字滤波器的单位脉冲响应输出及其稳定性3. 离散时间系统频域分析4.数字滤波器的设计实现四、实验内容及步骤. 设计一 DFT在信号频谱分析中的应用1. 用MATLAB语言编写计算序列x(n)的N点DFT的m函数文件dft.m。并与MATLAB中的内部函数文件fft.m作比较。dft.m程序：function Xk=dft(xn,N) if length(xn)N xn=xn,zeros(1,N-length(xn); end n=0:N-1; for k=0:N-1 Xk(1,k+1)=sum(xn.*exp(-1)*j*n*k*(2*pi/N); End 比较后发现与fft.m结果一样2. 对离散确定信号 作如下谱分析：(1) 截取使成为有限长序列N()，(长度N自己选)写程序计算出的N点DFT ,并画出相应的幅频图。 程序：n=0:11; xn=cos(0.48*pi*n)+cos(0.52*pi*n); Xk=fft(xn,12); subplot(2,1,1);stem(n,xn);grid;title(xn波形图); subplot(2,1,2);stem(n,abs(Xk);grid;title(幅频图); (2) 将 (1)中补零加长至M点（长度M自己选），编写程序计算的M点DFT ,并画出相应的图。 程序：n=0:11; xn=cos(0.48*pi*n)+cos(0.52*pi*n); n1=0:29; xn1=xn,zeros(1,18); Xk1=fft(xn1,30); subplot(2,1,1);stem(n,xn);grid; title(xn);subplot(2,1,2);stem(n1,abs(Xk1);grid;title(|Xk1|k);(3) 利用补零DFT计算 (1)中N点有限长序列频谱并画出相应的幅频图。 程序：n=0:11; xn=cos(0.48*pi*n)+cos(0.52*pi*n); n1=0:999; xn1=xn,zeros(1,988); Xk1=fft(xn1,1000); plot(n1,abs(Xk1);3. 在计算机上用DFT对模拟信号进行谱分析时, 只能以有限大的采样频率对模拟信号的近似频谱。其误差主要来源于截断效应(频谱泄漏和谱间干扰)和频谱混叠失真。前者使谱分辨率降低，产生谱间干扰；后者使折叠频率附近的频谱产生较大失真。 实践证明，加大截取长度可提高频率分辨率；选择合适的窗函数可降低谱间干扰；而频谱混叠失真要通过提高采样频率和预滤波来改善。按题目要求编写程序，验证截断效应及加窗的改善作用，参数选取如下： (1) 采样频率fs=400Hz,T=1/fs； (2) 对作4096点DFT作为的近似频谱； (3) 取三种截取时间长度分别： Tp=0.04s, Tp=4*0.04s,Tp=8*0.04s;两种窗函数分别为矩形窗和Hamming窗。实验结果记录：程序：1fs=800; T=1/fs; Tp=0.04; N=Tp*fs; N1=N,4*N,8*N; for m=1:3; n=1:N1(m); xn=cos(200*pi*n*T)+sin(100*pi*n*T)+cos(50*pi*n*T); Xk=fft(xn,4096); fk=0:4095/4096/T; subplot(3,2,2*m-1);%plot(xn);plot(fk,abs(Xk)/max(abs(Xk); if m=1 title(矩形窗截取); endend%加hamming窗改变谱间干扰for m=1:3; n=1:N1(m); wn2=hamming(N1(m); xn=wn2*(cos(200*pi*n*T)+sin(100*pi*n*T)+cos(50*pi*n*T); Xk=fft(xn,4096); fk=0:4095/4096/T; subplot(3,2,2*m);%plot(xn); plot(fk,abs(Xk)/max(abs(Xk); if m=1 title(hamming窗截取); endend 图形如下：设计二 数字滤波器的单位脉冲响应输出及其稳定性1. 自行编制程序计算数字滤波器的单位脉冲响应，给定差分方程：y(n)-y(n-1)+0.9y(n-2)=x(n)编制文件himp.m，实现数字滤波器的单位脉冲响应程序：A=1,-1,0.9;B=1; %系统差分方程系数向量B和A impz(B,A,70); %求系统单位脉冲响应h(n) title(系统单位脉冲响应h(n)的图像);2.给定计算此数字滤波器的单位阶跃响应，并判断系统稳定性。用stem(n,y)画出相应的图形。实验结果记录：程序：A=1,-1,0.9;B=1;hn=impz(B,A); xn=ones(1,20); yn=conv(hn,xn)%求系统阶跃响应stem(yn);title(系统单位阶跃响应的图形); z=roots(A);magz=abs(z) 结果：magz =0.9487 和0.9487在单位圆内，所以稳定。 图形：设计三 离散时间系统频域分析1. 用MATLAB语言编写计算N阶差分方程所描述系统频响函数的m函数文件fr.m。fr.m程序：function H=fr(b,a,w);%计算N阶差分方程所描述系统频响函数m=0:length(b)-1;l=0:length(a)-1; num=b*exp(-j*m*w);den=a*exp(-j*l*w); H=num./den2. 根据频响特性与系统零极点的关系，自己构造一个N阶差分方程，使该差分方程为数字低通滤波器。利用MATLAB程序画出相应的幅频图。 程序：b=1,1;a=1,-1/2,1/6;u=0:999;w=(pi/1000)*u; H=fr(b,a,w);plot(w,abs(H),b);title(低通滤波器的幅频图);ylabel(|H(ejw)|); xlabel(w);3. 改变2.中差分方程的系数，使该差分方程分别为数字高通及全通滤波器。利用MATLAB程序画出相应的幅频图。数字高通 程序：b=1/7,-2/7,1/7;a=1,3/4,1/4;u=0:999;w=(pi/1000)*u; H=fr(b,a,w);plot(w,abs(H),b);title(高通滤波器的幅频图);ylabel(|H(ejw)|); xlabel(w);数字全通 程序：b=1,1,1;a=1,1,1;u=0:999;w=(pi/1000)*u; H=fr(b,a,w);plot(w,abs(H),b);title(全通滤波器的幅频图);ylabel(|H(ejw)|); xlabel(w); 设计四 数字滤波器的设计及实现抑制载波单频调幅信号的数学表达式为 其中，称为载波，为载波频率，称为单频调制信号，为调制正弦波信号频率，且满足。由式可见，所谓抑制载波单频调制信号，就是两个正弦信号相乘，它有2个频率成分：和频，差频，这两个频率成分关于载波频率对称。所以，1路抑制载波单频调幅信号的频谱图是关于载波频率对称的两根谱线。显然，当调制频率和（或）载波频率不同时，可以得到包含不同频率成分的单频调幅信号，将几路不同频率成分的单频调幅信号相加后形成混合信号，产生复合信号的函数，三路调幅信号的载波频率分别为250Hz、500Hz和1000Hz。1. 调用信号产生函数mstg产生由三路抑制载波调幅信号相加构成的复合信号st，观察st的时域波形和幅频特性曲线；程序：function st=mstg%产生信号序列st，并显示st的时域波形和频谱%st=mstg返回三路调幅信号相加形成的混合信号，长度N=800N=800; %信号长度N为800Fs=10000;T=1/Fs;Tp=N*T; %采样频率Fs=10kHz，Tp为采样时间t=0:T:(N-1)*T;k=0:N-1;f=k/Tp; fc1=Fs/10; %第1路调幅信号载波频率fc1=1000Hzfm1=fc1/10; %第1路调幅信号的调制信号频率fm1=100Hzfc2=Fs/20; %第2路调幅信号载波频率fc2=500Hzfm2=fc2/10; %第2路调幅信号的调制信号频率fm2=50Hzfc3=Fs/40; %第3路调幅信号载波频率fc3=250Hzfm3=fc3/10; %第3路调幅信号的调制信号频率fm3=25Hzxt1=cos(2*pi*fm1*t).*cos(2*pi*fc1*t); %产生第1路调幅信号xt2=cos(2*pi*fm2*t).*cos(2*pi*fc2*t); %产生第2路调幅信号xt3=cos(2*pi*fm3*t).*cos(2*pi*fc3*t); %产生第3路调幅信号st=xt1+xt2+xt3; %三路信号相加，得到复合信号fxt=fft(st,N); %计算信号st的频谱%以下为绘图命令subplot(2,1,1);plot(t,st);grid;xlabel(t/s);ylabel(s(t);axis(0,Tp,min(st),max(st);title(s(t)的波形)subplot(2,1,2);stem(f,abs(fxt)/max(abs(fxt),.);grid;title(s(t)的频谱)axis(0,Fs/8,0,1.2);xlabel(f/Hz);ylabel(幅度);波形及曲线：2. 通过观察st的幅频特性曲线，分别确定可以分离st中三路抑制载波单频调幅信号的三个滤波器（低通滤波器、带通滤波器、高通滤波器）的通带截止频率和阻带截止频率；假定要求滤波器的通带最大衰减为0.1dB，阻带最小衰减为60dB，编程调用MATLAB滤波器设计函数分别设计这三个数字滤波器，并绘图显示其幅频特性曲线。3. 用所设计的三个滤波器分别对复合信号st进行滤波，分离出st中的三路不同载波频率的调幅信号，并绘图显示滤波后信号的时域波形和频谱，观察分离效果。 实验结果记录:N=800; %信号长度N为800Fs=10000;T=1/Fs;Tp=N*T; %采样频率Fs=10kHz，Tp为采样时间t=0:T:(N-1)*T;k=0:N-1;f=k/Tp; fc1=Fs/10; %第1路调幅信号载波频率fc1=1000Hzfm1=fc1/10; %第1路调幅信号的调制信号频率fm1=100Hzfc2=Fs/20; %第2路调幅信号载波频率fc2=500Hzfm2=fc2/10; %第2路调幅信号的调制信号频率fm2=50Hzfc3=Fs/40; %第3路调幅信号载波频率fc3=250Hzfm3=fc3/10; %第3路调幅信号的调制信号频率fm3=25Hzxt1=cos(2*pi*fm1*t).*cos(2*pi*fc1*t); %产生第1路调幅信号xt2=cos(2*pi*fm2*t).*cos(2*pi*fc2*t); %产生第2路调幅信号xt3=cos(2*pi*fm3*t).*cos(2*pi*fc3*t); %产生第3路调幅信号st=xt1+xt2+xt3; %三路信号相加，得到复合信号低通：fp1=100;fs1=200; wp1=2*fp1*pi/Fs;ws1=2*fs1*pi/Fs;rp=0.1;rs=60; %滤波器边界频率（关于Pi归一化）N1,wc1=buttord(wp1,ws1,rp,rs); B1,A1=butter(N1,wc1,low); h1,w1=freqz(B1,A1,N); x1=filter(B1,A1,st); %滤波器软件实现，即对信号st进行滤波X1=fft(x1,N); %x1的N点离散傅里叶变换%以下为绘图部分subplot(3,1,1); %低通滤波器的幅频特性曲线plot(w1/pi,abs(h1) title(低通滤波器的幅频特性); xlabel(X轴(w/pi); ylabel(幅度); subplot(3,1,2); %绘制滤波后信号的时域波形plot(t,x1); title(滤波后信号的时域波形); xlabel(t/s); ylabel(xt1); axis(0 0.1 min(x1),1.2*max(x1); subplot(3,1,3); %绘制滤波后信号的频谱stem(f,abs(X1)/max(abs(X1),.); title(滤波后信号的频谱); xlabel(f/Hz); ylabel(幅度); axis(0 1500 0 2);带通：fp2=440;fpu=560;fs2=400;fsu=700; wp2=2*fp2/Fs,2*fpu/Fs;ws2=2*fs2/Fs,2*fsu/Fs; rp=0.1;rs=60; N2,wpo2=ellipord(wp2,ws2,rp,rs); B2,A2=ellip(N2,rp,rs,wpo2); h2,w2=freqz(B2,A2,N); x2=filter(B2,A2,st); %滤波器软件实现，即对信号st进行滤波X2=fft(x2,N);%x2的N点离散傅里叶变换%以下为绘图部分subplot(3,1,1);%带通滤波器的幅频特性曲线plot(w2/pi,abs(h2) title(带通滤波器的幅频特性); xlabel(X轴（w/pi）); ylabel(幅度); subplot(3,1,2);%绘制滤波后信号的时域波形plot(t,x2); title(滤波后信号的时域波形); xlabel(t/s); ylabel(xt2); axis(0 0.1 min(x2),1.2*max(x2); subplot(3,1,3);%绘制滤波后信号的频谱stem(f,abs(X2)/max(abs(X2),.); title(滤波后信号的频谱); xlabel(f/Hz); ylabel(幅度); axis(0 1500 0 2);高通：fp3=850;fs3=500;%由图取高通滤波器的通、阻带截止频率为850和500 wp3=2*fp3/Fs;ws3=2*fs3/Fs;rp=0.1;rs=60;%高通滤波器的指标参数N3,wc3=buttord(wp3,ws3,rp,rs); %调用buttord和butter直接设计数字滤波器B3,A3=butter(N3,wc3,high); h3,w3=freqz(B3,A3,N);% N个点的频率响应值x3=filter(B3,A3,st);%滤波器软件实现，即对信号st进行滤波X3=fft(x3,N);%x3的N点离散傅里叶变换%以下为绘图部分subplot(3,1,1);%绘制高通滤波器的幅频特性曲线plot(w3/pi,abs(h3) title(高通滤波器的幅频特性); xlabel(X轴（w/pi）); ylabel(幅度);subplot(3,1,2);%绘制滤波后信号的时域波形plot(t,x3); title(滤波后信号的时域波形); xlabel(t/s); ylabel(xt3); axis(0 0.01 min(x3),1.2*max(x3); subplot(3,1,3);%绘制滤波后信号的频谱stem(f,abs(X3)/max(abs(X3),.); title(滤波后信号的频谱); xlabel(f/Hz); ylabel(幅度); axis(0 1500 0 2);五、实验思考题分析（1）说明补零DFT的作用，补零DFT能否提高信号的频谱分辨率，说明降低频谱泄露、频谱分辨率的措施各是什么？答：作用：DF