数字信号处理课程设计--使用Matlab工具进行数字信号处理.docx

课程设计报告( 2009 - 2010年度第2 学期)名称：数字信号处理题目：使用Matlab工具进行数字信号处理院系班级：学号：学生姓名指导教师：设计周数：一周成绩：日期：2010年7 月8日课程课程设计报告数字信号处理（自）课程设计任务书一、目的与要求是使学生通过上机使用Matlab工具进行数字信号处理技术的仿真练习，加深对数字信号处理（自）课程所学基本理论和概念的理解，培养学生应用Matlab等工具进行数字信号处理的基本技能和实践能力，为工程应用打下良好基础。二、主要内容1了解Matlab基本使用方法，掌握Matlab数字信号处理的基本编程技术。掌握数字信号的基本概念。1）使用Matlab（生成几种典型数字信号（正弦信号、周期信号、高斯随机信号等）。2）编程计算离散信号的特征值（均值、方差等）。3）进行信号加减运算。2Matlab编程实现典型离散信号（正弦信号、周期信号、随机信号）的离散傅立叶变换，显示时域信号和频谱图形（幅值谱和相位谱）；以周期信号为例，观察讨论基本概念（混叠、泄漏、整周期截取、频率分辨率等）。3设计任意数字滤波器，并对某类型信号进行滤波，并对结果进行显示和分析。三、进度计划序号设计(实验)内容完成时间备注1布置课程设计任务、示范讲解Matlab程序使用方法0.5天2上机进行课程设计3天3整理课程设计报告1天4演示课程设计内容并答辩0.5天四、设计成果要求1 提交完成设计内容的程序2 提交设计报告五、考核方式课程设计报告、设计内容演示和答辩相结合。考核内容：考勤、纪律、课程设计报告、实际编程能力和基本概念掌握程度等。学生姓名：指导教师 2010 年 7 月2日一、 课程设计的目的与要求通过上机使用Matlab工具进行数字信号处理技术的仿真练习，加深对数字信号处理（自）课程所学基本理论和概念的理解，培养学生应用Matlab等工具进行数字信号处理的基本技能和实践能力，为工程应用打下良好基础。二、 设计正文1. 了解Matlab基本使用方法，掌握Matlab数字信号处理的基本编程技术。掌握数字信号的基本概念。1.1使用Matlab（生成几种典型数字信号（正弦信号、周期信号、高斯随机信号等）。%生成正弦信号omega=pi/8 ;%设置数字角频率ns1=0;nf1=32;n1=ns1:nf1;xn1=sin(omega*n1);subplot(3,1,1);stem(n1,xn1,.);axis(0,35,-1.2,1.2);xlabel(n);ylabel(xn);title(a)正弦信号);grid%生成周期信号x=1 1 1 1 0 0 0;%周期为7xn2=x*ones(1,3);%产生3个周期长度的序列xn2=xn2(:);xn2=xn2;n2=0:length(xn2)-1;subplot(3,1,2);stem(n2,xn2,.);axis(0,20,-1.5,1.5);xlabel(n);ylabel(xn);title(b)周期信号);grid%生成高斯随机信号n3=40;%设置序列长度为40xn3=randn(1,n3);%产生均值为0方差为1的高斯随机序列subplot(3,1,3);stem(xn3);xlabel(n);ylabel(xn);title(c)高斯随机信号);grid1.2编程计算离散信号的特征值（均值、方差等）。%生成正弦信号omega=pi/8;%设置数字角频率clcns1=0;nf1=32;n1=ns1:nf1;xn1=sin(omega*n1);%生成周期信号x=1 1 1 1 0 0 0;%周期为7xn2=x*ones(1,3);%产生3个周期长度的序列xn2=xn2(:);xn2=xn2;n2=0:length(xn2)-1;%生成高斯随机信号n3=40;%设置序列长度为40xn3=randn(1,n3);%产生均值为0方差为1的高斯随机序列fprintf(正弦信号的均值为%.4f 方差分别为%.4fn,mean(yn1),var(yn1,1);fprintf(周期信号的均值为%.4f 方差分别为%.4fn,mean(x),var(x,1);fprintf(高斯随机信号的均值为%.4f 方差分别为%.4fn,mean(xn),var(xn,1);程序运行结果如下正弦信号的均值为-0.0000 方差分别为0.4848周期信号的均值为0.5714 方差分别为0.2449高斯随机信号的均值为-0.0024 方差分别为0.95061.3进行信号加减运算。%生成正弦信号clcomega=pi/8;%设置数字角频率ns1=0;nf1=32;n1=ns1:nf1;xn1=sin(omega*n1);%生成周期信号x=1 1 1 1 0 0 0;%周期为7xn2=x*ones(1,3);%产生3个周期长度的序列xn2=xn2(:);xn2=xn2;ns2=0;nf2=length(xn2)-1;n2=0:nf2;ny=0:max(nf1,nf2);%y的位置向量y1=zeros(1,length(ny);y2=y1;%y1,y2序列初始化y1(find(ny=nf1)=xn1;y2(find(ny=nf2)=xn2;ya=y1+y2;%两序列相加ys=y1-y2;%两序列相减%画图subplot(4,2,1);stem(n1,xn1,.);xlabel(n);ylabel(xn);title(正弦信号);gridsubplot(4,2,2);stem(ny,y1,.);xlabel(n);ylabel(xn);title(修正后的正弦信号);gridsubplot(4,2,3);stem(n2,xn2,.);xlabel(n);ylabel(xn);title(周期信号);gridsubplot(4,2,4);stem(ny,y2,.);xlabel(n);ylabel(xn);title(修正后的周期信号);gridsubplot(4,2,6);stem(ny,ya,.);xlabel(n);ylabel(xn);title(序列相加);gridsubplot(4,2,8);stem(ny,ys,.);xlabel(n);ylabel(xn);title(序列相减);grid2Matlab编程实现典型离散信号（正弦信号、周期信号、随机信号）的离散傅立叶变换，显示时域信号和频谱图形（幅值谱和相位谱）；以周期信号为例，观察讨论基本概念（混叠、泄漏、整周期截取、频率分辨率等）。2.1 Matlab编程实现典型离散信号的离散傅立叶变换 %生成正弦信号clcomega=pi/8;%设置数字角频率ns1=0;nf1=32;n1=ns1:nf1;xn1=sin(omega*n1);xk164=fft(xn1,128);%计算xn1的128点dft%以下为正弦信号的绘图k=0:127;wk=2*k/128;%产生128点dft对应的采样点频率（关于pi归一化值）subplot(3,1,1);stem(n1,xn1,.);xlabel(n);ylabel(xn);title(正弦信号);gridsubplot(3,1,2);stem(wk,abs(xk164),.);%绘制128点dft的幅频特性图title(128点dft的幅频特性图);xlabel(omega/pi);ylabel(幅度);subplot(3,1,3);stem(wk,angle(xk164),.);%绘制128点dft的相频特性图title(128点dft的相频特性图);xlabel(omega/pi);ylabel(相位);%生成周期信号clcx=1 1 1 1 0 0 0;%周期为7xn=x*ones(1,10);%产生10个周期长度的序列xn=xn(:);xn=xn;n=0:length(xn)-1;xk=fft(xn,128);%计算xn的128点dft%以下为周期信号的绘图k=0:127;wk=2*k/128;%产生128点dft对应的采样点频率（关于pi归一化值）subplot(3,1,1);stem(n,xn,.);xlabel(n);ylabel(xn);title(周期信号);gridsubplot(3,1,2);stem(wk,abs(xk),.);%绘制128点dft的幅频特性图title(128点dft的幅频特性图);xlabel(omega/pi);ylabel(幅度);subplot(3,1,3);stem(wk,angle(xk),.);%绘制128点dft的相频特性图title(128点dft的相频特性图);xlabel(omega/pi);ylabel(相位);%生成高斯随机信号clcn3=40;%设置序列长度为40xn3=randn(1,n3);%产生均值为0方差为1的高斯随机序列xk364=fft(xn3,128);%计算xn3的128点dft%以下为高斯随机信号的绘图k=0:127;wk=2*k/128;%产生128点dft对应的采样点频率（关于pi归一化值）subplot(3,1,1);stem(xn3,.);xlabel(n);ylabel(xn);title(高斯随机信号);gridsubplot(3,1,2);stem(wk,abs(xk364),.);%绘制128点dft的幅频特性图title(128点dft的幅频特性图);xlabel(omega/pi);ylabel(幅度);subplot(3,1,3);stem(wk,angle(xk364),.);%绘制128点dft的相频特性图title(128点dft的相频特性图);xlabel(omega/pi);ylabel(相位);2.2时域混叠的分析以周期信号为例，一个周期信号长度为7，截取10个周期。对其进行128点DFT，如图所示（a）、（b），然后抽取其中的64点进行IDFT 还原，如图（d）。因为采样点数N小于序列长度M，故无法由频域采样恢复原信号，在两周期处发生时域混叠。对128点DFT进行IDFT 还原，可以无失真的恢复原信号，如图（c）所示。上述混叠实质上是频域采样定理所描述的情形。%生成周期信号clcx=1 1 1 1 0 0 0;%周期为7xn=x*ones(1,10);%产生10个周期长度的序列xn=xn(:);xn=xn;n=0:length(xn)-1;xk=fft(xn,128);%计算xn的128点dftxn1=ifft(xk,128);xk1=xk(1:2:128);%抽取64点dft信号xn2=ifft(xk1,128);%以下为周期信号的绘图k=0:127;wk=2*k/128;%产生128点dft对应的采样点频率（关于pi归一化值）subplot(4,1,1);stem(n,xn,.);%绘制原始信号xlabel(n);ylabel(xn);title(a)周期信号);gridsubplot(4,1,2);stem(wk,abs(xk),.);%绘制128点dft的幅频特性图title(b)128点dft的幅频特性图);xlabel(omega/pi);ylabel(幅度);gridsubplot(4,1,3);stem(0:length(xn1)-1,xn1,.);title(c)抽取128点dft还原的信号);xlabel(n);ylabel(xn);gridsubplot(4,1,4);stem(0:length(xn2)-1,xn2,.);title(d)抽取64点dft还原的信号);xlabel(n);ylabel(xn);grid2.3频谱泄露的分析周期信号是无限长的，要用DFT分析其频谱，必须对周期信号进行截断。用矩形窗截断后信号的频谱向两边展宽，使频谱变模糊，谱分辨率降低。这种现象称为频谱泄露。一般矩形窗的长度越长，展宽就越窄。%生成周期信号clcx=1 1 1 1 0 0 0;%周期为7for i=1:2;switch i;case 1xn=x*ones(1,3);%产生3个周期长度的序列case 2xn=x*ones(1,20);%产生20个周期长度的序列 endfigure(i)xn=xn(:);xn=xn; n=0:length(xn)-1;xk=fft(xn,128);%计算xn的128点dft %以下为周期信号的绘图 k=0:127;wk=2*k/128;%产生128点dft对应的采样点频率（关于pi归一化值）subplot(2,1,1);stem(n,xn,.);xlabel(n);ylabel(xn);title(周期信号);grid subplot(2,1,2);stem(wk,abs(xk),.);%绘制128点dft的幅频特性图 title(128点dft的幅频特性图);xlabel(omega/pi);ylabel(幅度);end前两图是矩形窗为3个周期长度截取的信号及其频谱，后两图是矩形窗为20个周期长度截取的信号及其频谱。观察两者频谱，很容易理解频谱泄露的含义。2.4整周期截取的分析%生成周期信号clcx=1 1 1 1 0 0 0;%周期为7xn=x*ones(1,10);%产生10个周期长度的序列xn=xn(:);xn=xn;n=0:length(xn)-1;xk=fft(xn,128);%计算xn的128点dftxn1=xn(1:45);%截取非整周期序列xk1=fft(xn1,128);%以下为周期信号的绘图k=0:127;wk=2*k/128;%产生128点dft对应的采样点频率（关于pi归一化值）subplot(4,1,1);stem(n,xn,.);%绘制原始信号xlabel(n);ylabel(xn);title(a)整周期截取的周期信号);gridsubplot(4,1,2);stem(wk,abs(xk),.);%绘制128点dft的幅频特性图title(b)整周期截取128点dft的幅频特性图);xlabel(omega/pi);ylabel(幅度);gridsubplot(4,1,3);stem(0:length(xn1)-1,xn1,.);%绘制非整周期截取信号xlabel(n);ylabel(xn);title(c)非整周期截取周期信号);gridsubplot(4,1,4);stem(wk,abs(xk1),.);%绘制128点dft的幅频特性图title(d)非整周期截取128点dft的幅频特性图);xlabel(omega/pi);ylabel(幅度);grid很明显，非正整周期截取时它的最大数字频率不在0.20.4之间，这与理论值相差是比较大的。在截取长度不是很长时必须要整周期截取，这样才能保证截断误差不是很大。2.5频率分辨率的分析频谱分辨率是指能够分辨的两个频率分量最小的间隔。3设计任意数字滤波器，并对某类型信号进行滤波，并对结果进行显示和分析。%设计巴特沃兹滤波器clcwp=input(请输入通带边界频率wp=);ws=input(请输入阻带边界频率ws=);rp=input(请输入通带最大衰减rp=);rs=input(请输入阻带最小衰减rs=);N,wc=buttord(wp,ws,rp,rs,s)%计算相应的模拟低通滤波器阶数N和3dB截止频率B,A=butter(N,wc,s);%计算相应的模拟滤波器的系统函数tf(B,A)Bz,Az=impinvar(B,A,8000);%用脉冲响应不变法将模拟滤波器转换成数字滤波器H,w=freqz(Bz,Az);%计算数字滤波器的频率响应db=20*log10(abs(H);figure(1);hold onplot(w/pi,db);gridxlabel(omega/pi);ylabel(幅度（db）);title(sprintf(%d阶数字低通滤波器的损耗函数,N);%以下滤波过程x=0:0.1:40;y1=10*sin(pi*x/8);figure(2);subplot(2,2,1);plot(x,y1);gridtitle(正弦信号)y2=randn(1,length(x);subplot(2,2,2);plot(x,y2);gridtitle(高斯随机信号)y3=y1+y2;subplot(2,2,3);plot(x,y3);gridtitle(滤波前信号);y=filter(Bz,Az,y3);