数字信号处理实验报告书

数字信号处理实验报告书班级：姓名：学号：联系方式：2019年12月实验一离散时间信号产生，图形显示及文件读存实验实验目的1、熟悉常用的MATLAB命令，掌握离散时间信号的时域表示方法；2、掌握用MATLAB描绘二维图像的方法；3、掌握用MATLAB存读文件的方法。实验原理及实验内容1、实验原理序列的表示方法有列举法、解析法和图形法，相应的用MATLAB也可以用以下几种表示方法，具体介绍如下：（1）列举法在MATLAB中，用一个列向量来表示一个有限长序列，由于一个列向量并不包含位置信息，因此需要用表示位置的n和表示量值的x两个向量来表示任意一个序列，如：例1.1：n=[-3,-2,-1,0,1,2,3,4];x=[251203-14];如果不对向量的位置进行定义，则MATLAB默认该序列的起始位置为n=0。由于内存有限，MATLAB不能表示一个无限序列。（2）解析法对于有解析表达式的确定信号，首先定义序列的范围即n的值，然后直接写出该序列的表达式，如：例1.2：实现实指数序列的MATLAB程序为：n=0:10;x=(0.5).^n;（3）图形法在MATLAB中用图形法表示一个序列，是在前两种表示方法的基础上将序列的各个量值描绘出来，即首先对序列进行定义，然后用相应的画图语句画图，如：例1.3：绘制在1中用列举法表示的序列的图形，则在向量定义之后加如下绘图语句：stem(n,x);MATLAB中数据格式主要有二进制格式和ASCII码格式，使用函数完成数据存取功能。2、实验内容（1）产生一个单位抽样序列和一个有延迟的单位抽样序列，并在同一界面绘制其序列图形。（2）产生一个向前时移3个时刻的单位阶跃序列：，其中，绘制其序列的图形。（3）产生一个复指数序列：，绘制其序列的图形，绘制其实部、虚部、幅度和相位波形。（4）将数据文件保存为.dat文件和ASCII码相应文件实验程序及实验结果1-1>>n=[-1:8];>>x=[zeros(1,1),1,zeros(1,8)];>>stem(n,x)>>axis([-5,10,0,2]);>>>>n=[-1:8];>>x=[zeros(1,6),1,zeros(1,3)];>>stem(n,x)>>axis([-5,10,0,2]);>>1-2>>n=[-5:5];>>x=[zeros(1,2),ones(1,9)];>>stem(n,x)>>1-3n0=-30;n2=30;n=n0:n2;x=5*exp(3j*n);figure(1)subplot(411)stem(n,real(x),'filled');axis([-30,30,min(real(x))-1,1.2*max(real(x))])title('复指数序列')ylabel('实部');grid;subplot(412)stem(n,imag(x),'filled');axis([-30,30,min(imag(x))-1,1.2*max(imag(x))])ylabel('虚部');xlabel('n');grid;subplot(413)stem(n,abs(x),'filled');axis([-30,30,min(abs(x))-1,1.2*max(abs(x))])ylabel('振幅');xlabel('n');grid;subplot(414)stem(n,angle(x),'filled');axis([-30,30,min(angle(x))-1,1.2*max(angle(x))])ylabel('相位');xlabel('n');grid;心得体会很久没用MATLAB了，不会新建脚本，在编辑器中写代码而在命令行窗口编辑代码，造成了很多不必要的麻烦。而后请教同学打开了编辑器。编写代码方便了许多。通过实验一，我学会了离散时间信号产生，图形显示及文件读存，在课本以及网上找类似的代码“临摹”，不仅有助于我学习使用MATLAB产生抽样序列和复指数序列的方法，更培养了我独立自学的能力。

实验二离散信号的卷积和实验实验目的1、掌握两个离散信号卷积和的计算方法和编程技术；2、掌握用MATLAB描绘二维图像的方法。实验原理及实验内容1、卷积和的原理两个离散序列x(n)与y(n)的卷积和f(n)定义为：由于通常信号处理中所碰到的都是有始信号或有限时间信号，因此在实际计算卷积和时，求和是在有限范围内进行的。计算过程中上下限的选取和所得结果的分布区间取决于参与卷积的两个序列，下面将分别进行讨论：两个从n=0开始的序列和的卷积：(1)上式右边因子u(n)表示卷积和的结果也是一个从n=0开始的序列。从n=n1开始的序列和从n=n2开始的序列的卷积和，其中n1和n2为任意整数。(2)上式右边因子u(n-n1-n2)表示卷积和是一个从n=n1+n2开始的序列。从n=n1开始的长度为N1的加窗序列和从n=n2开始的长度为N2的加窗序列的卷积和，其中则(3)所得卷积和也是一个加窗序列，从n=n1+n2开始，长度为N1+N2-1。MATLAB提供了一个内部函数conv(x,h)用来计算两个有限长度序列的卷积，该函数得到的卷积结果默认从n=0开始，因此当参与卷积的两个序列的起始位置不是n=0时，则由该函数得到的计算结果将出现错误，此时需要重新定义结果的位置向量。由以上卷积运算的原理可知，两有限长序列卷积后仍为有限长序列，长度为两序列长度之和减1，结果的起始位置为两序列起始位置之和，截止位置为两序列截止位置之和。据此，可以得到卷积结果的位置向量。2、实验内容已知某系统的单位取样响应为，试用MATLAB求当激励信号为时，系统的零状态响应。系统的零状态响应是激励与系统的单位取样响应的卷积。实验程序及实验结果subplot(311)n=-2:8;y1=[zeros(1,2),ones(1,9)];y2=[zeros(1,6),ones(1,5)];x=y1-y2;stem(n,x)axis([-2,8,0,2]);subplot(312)n=-2:8;y1=[zeros(1,2),ones(1,9)];y3=[zeros(1,10),ones(1,1)];h=(0.8).^n.*(y1-y3);stem(n,h)axis([-2,8,0,2]);subplot(313)z=conv(x,h);plot(z);心得体会通过这次实验，我领会到了严谨治学，一丝不苟的态度，在科学研究中的重要性，有时科学研究会是极其枯燥的，可能是复杂的数学演算，或是一次次重复的实验，但不论是什么,我们都要抱以认真务实的态度对待科学，绝不可草草了事，我们现在是站在巨人的肩膀上学习，要抓重点，抓主线，把最本质的东西弄懂，弄通。比如傅立叶变换，才能在学习心知识的时候多思考，多总结发现内在规律才能有创造性的思考问题。

实验三吉布斯现象实验实验目的1、理解吉布斯现象的定义；2、掌握矩形波信号及傅里叶级数的表示方法；3、观察信号对比效果图，进一步认识吉布斯现象。实验原理及实验内容1、吉布斯现象介绍将具有不连续点的周期函数（如矩形脉冲)进行傅立叶级数展开后，选取有限项进行合成。当选取的项数越多，在所合成的波形中出现的峰起越靠近原信号的不连续点。当选取的项数很大时，该峰起值趋于一个常数，大约等于总跳变值的9%。2、实验内容（1）生成一个周期为T，幅值为F，占空比为50%的矩形波信号；【注】MATLAB中可能用到的函数：方波生成函数square()。一般形式为“y=A*square(B*t+C)”,其中A是方波信号的幅值，B是方波的角频率，C是方波的初相。（2）求矩形波信号的傅里叶展开式；（3）用plot函数绘制傅里叶级数展开式，观察求和项数变化对矩形波合成情况的影响。实验程序及实验结果t=-2:0.01:2;T=[152551101];fork=1:length(T)a=0;forn=1:2:T(k)b=exp(1i*n*pi*t)/n;a=a+b;endy=2*a/(1j*pi);x=0.5*square(pi*t);figure;p=plot(t,real(y),t,x);axis([-2,2,-0.8,0.8]);set(gca,'XTick',-2:1:2)set(gca,'YTick',-2:1:2)set(gca,'XTicklabel',{'-2','-1','0','1','2'})set(gca,'YTicklabel',{'-0.8','-0.4','0','0.4','0.8'})xlabel('自变量')ylabel('函数值')titlemsg=sprintf('吉布斯现象N=%d的合成波形',T(k));title(titlemsg)text(0,0.5,'\leftarrow方波函数','HorizontalAlignment','left')set(gcf,'Color','w')end心得体会通过此次课程设计，我对设计所用到的软件MATLAB有了更加深刻地了解，MATLAB不管在数值计算方面的功能很强大，而且其图形仿真功能更能满足各个领域的需要，因此我们以后更要经常运用MATLAB软件，使其成为自己不可或缺的工具。

实验四离散信号FFT的实现实验实验目的1、了解离散信号FFT算法的思想2、了解DFT，DIF算法的思想3、学习并使用MATLAB掌握FFT算法。实验原理及实验内容1.简介在信号处理中，常常遇到的不是一个连续的函数，而是离散的数列。以等间隔时间取样的时间序列{，，，.......，}为例，这里的数据个数为N。对于有限长离散数字信号{x[n]}，0nN-1，其离散谱{x[k]}可以由离散傅氏变换（DFT）求得。DFT的定义为：上式也可表示为：2.快速傅里叶变换快速傅立叶变换是傅立叶变换的一种快速算法，简称FFT，采用FFT算法能够大大减小计算离散傅立叶变换所需要的乘法次数，特别是被变换的抽样点数N越多，FFT算法计算量的节省就越显著。由4.1中式X[k]易知一次复数惩罚相当于4此=次实数乘法和2次实数加法，而一次复数加法相当于两次实数加法。因此，每计算一个X[k]需要4N次实数乘法和2N+2（N-1）=2（2N-1）次实数加法。DFT计算共需要计算N个X[k]，因此，完成整个DFT算法共需要4N次实数乘法和4N2-2N次实数加法。可见直接计算N点DFT所需要的时间量和N是成正比的的，当N非常大时，计算量会显著提升。快速傅立叶算法的思路，就是利用复指数函数的性质减小计算量。由于是一个周期函数，具有以下性质。，它可利用周期性与对称性改进运算，使DFT运算中部分项合并，提高计算效率。同时也可以利用的周期性和对称性，把长度为N点的大点数的DFT运算依次分解为若干的小点数DFT。因为DFT的计算量正比与，N越小，计算量也越小。快速傅立叶算法的思路正是基于以上两种想法发展起来的，它有多种形式，但基本上可大致分为两类，第一种是按n域(时域)划分的DIT-FFT算法(decimation-in-timeFFT，库利·图基算法)。第二种则是按k域(频域)划分的DIF-FFT算法(decimation-in-frequencyFFT，桑德·图基算法)。3.快速傅立叶变换的MATLAB实现在MATLAB中，函数fft用于快速计算DFT，调用格式为x是取样的样本，可以是一个向量或者矩阵，y是x的快速傅里叶变换。在实际操作中会对x进行补0操作，使x的长度等于2的整数次幂，这样能够显著提高程序的计算速度。实验程序及实验结果fs=200;figure(2)subplot(411)N_data=64;N=64;n=0:N_data-1;t=n/fs;x=2*sin(50*2*pi*t)+10*sin(20*2*pi*t);y=fft(x,N);mag=abs(y);f=(0:N-1)*fs/N;plot(f(1:N/2),mag(1:N/2)*2/N);xlabel('频率/Hz');ylabel('振幅');title('Ndata=64;FFT采样点个数=64');gridon;subplot(412)N_data=64;N=128;n=0:N_data-1;t=n/fs;x=2*sin(50*2*pi*t)+10*sin(20*2*pi*t);y=fft(x,N);mag=abs(y);f=(0:N-1)*fs/N;plot(f(1:N/2),mag(1:N/2)*2/N);xlabel('频率/Hz');ylabel('振幅');title('Ndata=64;FFT采样点个数=128');gridon;subplot(413)N_data=128;N=128;n=0:N_data-1;t=n/fs;x=2*sin(50*2*pi*t)+10*sin(20*2*pi*t);y=fft(x,N);mag=abs(y);f=(0:N-1)*fs/N;plot(f(1:N/2),mag(1:N/2)*2/N);xlabel('频率/Hz');ylabel('振幅');title('Ndata=128;FFT采样点个数=128');gridon;subplot(414)N_data=128;N=512;n=0:N_data-1;t=n/fs;x=2*sin(50*2*pi*t)+10*sin(20*2*pi*t);y=fft(x,N);mag=abs(y);f=(0:N-1)*fs/N;plot(f(1:N/2),mag(1:N/2)*2/N);xlabel('频率/Hz');ylabel('振幅');title('Ndata=128;FFT采样点个数=512');gridon;

心得体会通过本实验加深了对FFT的理解，更熟练的应用Matlab软件，提高了自己编写相关程序的能力，在实验的过程中出现了些许错误，有操作上也有程序上的，通过本实验加强了自己动手纠错的能力。

实验五离散信号IFFT的实现实验实验目的1．学习并掌握IFFT算法；2．加深对DIF算法和IFFT算法的理解。实验原理及实验内容1.实验原理（1）DFT的定义：对于有限长离散数字信号{x[n]}，0≤n≤N-1，其离散谱{x[k]}可以由离散傅氏变换（DFT）求得。DFT的定义为：，k=0,1…,N-1通常令，称为旋转因子。直接计算DFT的问题及FFT的基本思想由DFT的定义可以看出，在x[n]为复数序列的情况下，完全直接运算N点DFT需要(N-1)2次复数乘法和N(N-1)次加法。因此，对于一些相当大的N值(1024)来说，直接计算它的DFT所作的计算量是很大的。FFT的基本思想在于，将原有的N点序列分成两个较短的序列，这些序列的DFT可以很简单的组合起来得到原序列的DFT。例如，若N为偶数，将原有的N点序列分成两个(N/2)点序列，那么计算N点DFT将只需要约[(N/2)2·2]=N2/2次复数乘法。即比直接计算少作一半乘法。因子(N/2)2表示直接计算（N/2）点DFT所需要的乘法次数，而乘数2代表必须完成两个DFT。上述处理方法可以反复使用，即(N/2)点的DFT计算也可以化成两个(N/4)点的DFT（假定N/2为偶数），从而又少作一半的乘法。这样一级一级的划分下去一直到最后就划分成两点的FFT运算的情况。（3）基2按频率抽取的（DIF）的FFT算法思想对于情况下的一种普遍使用FFT的结构是频率抽取法。对于频率抽取法，输入序列不是按奇偶数分解，而是按前后对半分开，这样便将N点DFT写成前后两部分，即，进一步分解偶数组和奇数组，即令于是有 这正是两个N/2点的DFT运算，即将一个N点的DFT分解为两个N/2点的DFT。与时间抽取法一样。由于，N/2仍是一个偶数，这样，一个点的DFT通过M次分解后，最后剩下的全部是2点的DFT，2点DFT实际上只有加减运算。（4）IFFT算法思想FFT变换是将信号从时域转换到频域，这样在时域复杂的信号转换到频域看起来就方便容易了很多。但有时候也需要将频域信号转换到时域，所以这时运用到IFFT变换。（详见书本）（5）在FFT，IFFT计算中使用到的MATLAB命令函数fft(x)可以计算R点序列的R点DFT值；而fft(x,N)则计算R点序列的N点DFT，若R>N，则直接截取R点DFT的前N点，若R<N，则x先进行补零扩展为N点序列再求N点DFT。函数ifft(X)可以计算R点的谱序列的R点IDFT值；而ifft(X,N)同fft(x,N)的情况。2、实验内容已知信号，其中，，采样频率200Hz，fft采点个数为128，绘制出原始信号，原始信号快速傅里叶变换的频谱图，最后用ifft处理频域信号（由fft得到的频域信号）返回时间域。实验程序及实验结果t=-1:0.001:1;x=0.7*sin(2*pi*20*t)+3*sin(80*3*pi*t);subplot(3,1,1)plot(t,x)xlabel('时间t');ylabel('函数值x');title('原始函数');fs=200;N_data=64;N=128;n=0:N_data-1;t=n/fs;y=fft(x,N);mag=abs(y);f=(0:N-1)*fs/N;subplot(3,1,2)plot(f(1:N/2),mag(1:N/2)*2/N);xlabel('频率/hz');ylabel('振幅');title('FFT变换');subplot(3,1,3)z=ifft(y,N);plot(z);xlabel('时间t');ylabel('函数值');title('IFFT逆变换');gridon;心得体会通过本实验加深了对IFFT的理解，更熟练的应用Matlab软件，提高了自己编写相关程序的能力，在实验的过程中出现了些许错误，有操作上也有程序上的，通过本实验加强了自己动手纠错的能力。

实验六IIR滤波器的设计及应用实验实验目的1、熟悉用双线性变换法设计IIR数字滤波器的原理与方法；2、学会调用MATLAB信号处理工具箱中滤波器设计函数（或滤波器设计分析工具fdatool）设计各种IIR数字滤波器，学会根据滤波需求确定滤波器指标参数。3、掌握IIR数字滤波器的MATLAB实现方法。4、通过观察滤波器输入输出信号的时域波形及其频谱，建立数字滤波的概念。实验原理及实验内容1、原理分析设计IIR数字滤波器一般采用间接法（脉冲响应不变法和双线性变换法），应用最广泛的是双线性变换法。基本设计过程是：①先将给定的数字滤波器的指标转换成过渡模拟滤波器的指标；②设计过渡模拟滤波器；③将过渡模拟滤波器系统函数转换成数字滤波器的系统函数。MATLAB信号处理工具箱中的各种IIR数字滤波器设计函数都是采用双线性变换法。滤波器设计函数butter、cheby1、cheby2和ellip可以分别被调用来直接设计巴特沃斯、切比雪夫1、切比雪夫2和椭圆模拟和数字滤波器。本实验要求读者调用如上函数直接设计IIR数字滤波器。本实验的数字滤波器的MATLAB实现是指调用MATLAB信号处理工具箱函数filter对给定的输入信号x(n)进行滤波，得到滤波后的输出信号y(n）。2、实验内容（1）调用信号产生函数mstg产生由三路抑制载波调幅信号相加构成的复合信号st，该函数还会自动绘图显示st的时域波形和幅频特性曲线，如图1所示。由图可见，三路信号时域混叠无法在时域分离。但频域是分离的，所以可以通过滤波的方法在频域分离，这就是本实验的目的。图1三路调幅信号st的时域波形和幅频特性曲线（2）要求将st中三路调幅信号分离，通过观察st的幅频特性曲线，分别确定可以分离st中三路抑制载波单频调幅信号的三个滤波器（低通滤波器、带通滤波器、高通滤波器）的通带截止频率和阻带截止频率。要求滤波器的通带最大衰减为0.1dB,阻带最小衰减为60dB。提示：抑制载波单频调幅信号的数学表示式为其中，称为载波，fc为载波频率，称为单频调制信号，f0为调制正弦波信号频率，且满足。由上式可见，所谓抑制载波单频调幅信号，就是2个正弦信号相乘，它有2个频率成分：和频和差频，这2个频率成分关于载波频率fc对称。所以，1路抑制载波单频调幅信号的频谱图是关于载波频率fc对称的2根谱线，其中没有载频成分，故取名为抑制载波单频调幅信号。容易看出，图1中三路调幅信号的载波频率分别为250Hz、500Hz、1000Hz。如果调制信号m(t)具有带限连续频谱，无直流成分，则就是一般的抑制载波调幅信号。其频谱图是关于载波频率fc对称的2个边带（上下边带），在专业课通信原理中称为双边带抑制载波(DSB-SC)调幅信号,简称双边带(DSB)信号。如果调制信号m(t)有直流成分，则就是一般的双边带调幅信号。其频谱图是关于载波频率fc对称的2个边带（上下边带），并包含载频成分。（3）编程序调用MATLAB滤波器设计函数ellipord和ellip分别设计这三个椭圆滤波器，并绘图显示其幅频响应特性曲线。（4）调用滤波器实现函数filter，用三个滤波器分别对信号产生函数mstg产生的信号st进行滤波，分离出st中的三路不同载波频率的调幅信号y1(n)、y2(n)和y3(n)，并绘图显示y1(n)、y2(n)和y3(n)的时域波形，观察分离效果。（3）信号产生函数mstgfunctionst=mstg%产生信号序列向量st,并显示st的时域波形和频谱%st=mstg返回三路调幅信号相加形成的混合信号，长度N=1600N=1600%N为信号st的长度。Fs=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=1000Hz,fm1=fc1/10; %第1路调幅信号的调制信号频率fm1=100Hzfc2=Fs/20; %第2路调幅信号的载波频率fc2=500Hzfm2=fc2/10; %第2路调幅信号的调制信号频率fm2=50Hzfc3=Fs/40; %第3路调幅信号的载波频率fc3=250Hz,fm3=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的频谱%====以下为绘图部分，绘制st的时域波形和幅频特性曲线====================subplot(3,1,1)plot(t,st);grid;xlabel('t/s');ylabel('s(t)');axis([0,Tp/8,min(st),max(st)]);title('(a)s(t)的波形')subplot(3,1,2)stem(f,abs(fxt)/max(abs(fxt)),'.');grid;title('(b)s(t)的频谱')axis([0,Fs/5,0,1.2]);xlabel('f/Hz');ylabel('幅度')实验程序及实验结果clearall;clearall;st=mstg;Fs=10000;fp=280;fs=450;wp=2*fp/Fs;ws=2*fs/Fs;rp=0.1;rs=60;[N,wp0]=ellipord(wp,ws,rp,rs);[B,A]=ellip(N,rp,rs,wp0);y1t=filter(B,A,st);figure(1);subplot(2,1,1);[H1,w]=freqz(B,A,1000);m=abs