数字信号处理实验指导书49418doc.docx

1、数字信号处理实验指导书49418doc数字信号处理实 验 指 导 书武汉理工大学教材中心2012年7月实验一 时域离散信号的产生一、实验目的1、了解常用时域离散信号及其特点；2、掌握MATLAB程序的编程方法；3、熟悉MATLAB函数的调用方法。二、实验原理在时间轴上的离散点取值的信号，称为离散时间信号。离散时间信号只在某些离散的瞬时给出函数值，而在其他时刻无定义。它是时间上不连续按一定先后次序排列的一组数的集合，称为时间序列，用x(n)表示，n取整数代表时间的离散时刻。在MATLAB中用向量来表示一个有限长度的序列。常用离散信号：1、单位抽样序列2、单位阶跃序列3、实指数序列4、复指数序列5

2、、正(余)弦序列6、随机序列 在利用计算机进行系统的研究时，经常需要产生随机信号，MATLAB提供一个工具函数rand来产生随机信号。7、周期序列 三、实验用函数1、stem功能：绘制二维图形。调用格式：stem(n,x)；n为横轴，x为纵轴的线性图形。2、length功能：计算某一变量的长度或采样点数。调用格式：N=length(t)；计算时间向量t的个数并赋给变量N。3、axis功能：限定图形坐标的范围。调用格式：axis(x1,x2,y1,y2)；横坐标从x1x2，纵坐标从y1y2。4、zeros功能：产生一个全0序列。调用格式：x=zeros(1,n)；产生n个0的序列。5、ones功

3、能：产生一个全1序列。调用格式：y=ones(1,n)；产生n个1的序列。四、参考实例例1.1 用Matlab产生单位抽样序列。%先建立函数impseq(n1,n2,n0)function x,n=impseq(n1,n2,n0) n=n1:n2; x=(n-n0)=0;%编写主程序调用该函数x,n=impseq(-2,8,2);stem(n,x)程序运行结果如图1-1所示：图1-1 单位抽样序列例1.2实数指数序列（运算符“.”）Matlab程序如下：n=0:10;x=0.9.n;stem(n,x)程序运行结果如图1-2所示图1-2 实数指数序列例1.3复数指数序列（）Matlab程序如下：

4、n=-10:10; alpha=-0.1+0.3*j; x=exp(alpha*n);real_x=real(x); image_x=imag(x);mag_x=abs(x); phase_x=angle(x);subplot(2,2,1); stem(n,real_x)subplot(2,2,2); stem(n,image_x)subplot(2,2,3); stem(n,mag_x)subplot(2,2,4); stem(n,phase_x)程序运行结果如图1-3所示图1-3 复数指数序列例1.4正、余弦序列（）Matlab程序如下：n=0:10;x=3*cos(0.1*pi*n+pi

5、/3);stem(n,x)程序运行结果如图1-4所示图1-4 正、余弦序列例1.5随机序列rand(1,N)产生其元素在0，1之间均匀分布长度为N的随机序列randn(1,N)产生均值为0，方差为1，长度为N的高斯随机序列例1.6周期序列如何生成周期序列 1、 将一个周期复制p次； 2、借助矩阵运算、matlab下标能力。先生成一个包含p列x(n)值的矩阵，然后用结构（:）来把p列串接成一个长周期序列。因为这个结构只能用于列向，最后还需要做矩阵转置获得所需序列。Matlab程序如下：x=1,2,3; %一个x(n)xn=x*ones(1,3) %生成p列x(n)xn=xn(:) %将p列串接成

6、长列序列并转置stem(xn)程序运行的结果如图1-5所示图1-5 周期序列五、实验任务1、调试部分例题程序，掌握Matlab基本操作方法。2、编写程序，完成下列函数波形： 1）利用zeros函数生成单位抽样序列； 2）利用zeros函数和ones函数生成单位阶跃序列； 六、实验报告1、简述实验目的、原理。2、写出上机调试通过的实验任务的程序并描述其图形曲线。实验二 离散序列的基本运算一、实验目的1、加强MATLAB运用。2、了解离散时间序列在时域中的基本运算。3、熟悉相关函数的使用方法，掌握离散序列运算程序的编写方法。二、实验原理离散序列的时域运算包括信号的相加、相乘，信号的时域变换包括信号

7、的移位、反折、倒相及尺度变换等。在MATLAB中，序列的相加和相乘运算是两个向量之间的运算，因此参加运算的两个序列必须具有相同的长度，否则不能直接进行运算，需要进行相应的处理后再进行运算。三、实验用函数1、find功能：寻找非零元素的索引号。调用格式：find(n=min(n1)&(n0时，原序列向右移m位，当m=min(n1)&(n=min(n2)&(n=max(n2)=1)=x2;y=y1+y2;1）序列相乘 序列相加是指两个序列中相同序号的序列值逐项对应相乘，形成新的序列。 参加运算的两个序列的维数不同时处理方法与序列相加相同。五、实验任务1、理解序列运算的性质，了解函数语句的意义。2、

8、利用例题函数完成下列序列运算 1）已知x1(n)=u(n+1) (-3n5); x2(n)=u(n-3) (-4n7) 求：x(n)=x1(n)+x2(n) 2）已知x1(n)=3e-0.25n (-2n8) x2(n)=u(n+1) (-3n6) 求：x(n)=x1(n)*x2(n)六、实验报告1、简述实验目的和原理。2、列写上机调试通过的程序，并描绘其波形曲线。实验三 离散卷积的原理及应用一、实验目的1、通过实验进一步理解卷积定理，了解卷积过程；2、掌握应用线性卷积求解离散时间系统响应的基本方法。二、实验原理对于线性移不变离散系统，任意的输入信号x(n)可以用及其位移的线性组合来表示，即

9、当输入为时，系统的输出y(n)=h(n)，由系统的线性移不变性质可以得到系统对x(n)的响应y(n)为 称为离散系统的线性卷积，简写为 y(n)=x(n)*h(n)也就是说，如果已知系统的冲激响应，将输入信号与系统的冲激响应进行卷积运算，即可求得系统的响应。三、实验用函数1、卷积函数conv功能：进行两个序列的卷积运算。调用格式：y=conv(x,h)；用于求解两有限长序列的卷积。2、sum功能：求各元素之和。调用格式：y=sum(x)；求序列x中各元素之和。3、hold功能：控制当前图形窗口是否刷新的双向切换开关。调用格式：hold on：使当前图形窗口中的图形保持且不被刷新，准备接受绘制新

10、的图形。hold off：使当前图形窗口中的图形不具备不被刷新的性质。4、impz功能：求解数字系统的冲激响应。调用格式：h,t=impz(b,a)；求解数字系统的冲激响应h，取样点数为缺省值。h,t=impz(b,a,n)；求解数字系统的冲激响应h，取样点数由n确定。impz(b,a)；在当前窗口用stem(t,h)函数绘制图形。5、dstep功能：求解数字系统的阶跃响应。调用格式：h,t=dstep(b,a)；求解数字系统的冲激响应h，取样点数为缺省值。h,t=dstep(b,a,n)；求解数字系统的冲激响应h，取样点数由n确定。dstep(b,a)；在当前窗口用stairs(t,h)函数

11、绘制图形。四、参考实例在利用Matlab提供的卷积函数进行卷积运算时，主要是确定卷积结果的时间区间。conv函数默认两信号的时间序列从n=0开始，卷积结果对应的时间序列也从n=0开始。如果信号不是从0开始，则编程时必须用两个数组确定一个信号，其中，一个数组是信号波形的幅度样值，另一个数组是其对应的时间向量。%建立一个适用于信号从任意时间开始的通用函数function y,ny=sconv(x,h,nx,nh,p) y=conv(x,h); n1=nx(1)+nh(1); %计算y的非零样值的起点位置 n2=nx(length(x)+nh(length(h); %计算y的非零样值的宽度 ny=n

12、1:p:n2; %确定y的非零样值的时间向量五、实验任务已知一个IIR数字低通滤波器的系统函数公式为输入一个矩形信号序列 x=square(n/5) (-2n10)，求该系统的响应。六、实验报告1、简述实验目的、原理。2、写出上机调试通过的实验任务的程序并描述其图形曲线。实验四 离散傅立叶级数一、实验目的1、加深对离散周期序列傅里叶级数基本概念的理解。2、掌握MATLAB求解周期序列傅里叶级数变换和逆变换的方法。3、观察离散周期序列的重复周期数对频谱特性的影响。二、实验原理离散时间序列x(n)满足x(n)=x(n+rN)，称为离散周期序列，其中N为周期，x(n)为主值序列。周期序列可用离散傅里

13、叶级数表示成 n=0,1,N-1其中，是周期序列离散傅里叶级数第K次谐波分量的系数，也称为周期序列的频谱，可表示为 k=0,1,N-1上面两式是周期序列的一对傅里叶级数变换对。令，以上两式可简写为：三、实验用函数1、mod功能：模除求余。调用格式：mod(x,m)：x整除m取正余数。2、floor功能：向-舍入为整数。调用格式：floor(x)：将x向-舍入为整数。四、实例1、周期序列的傅里叶变换和逆变换 依据变换公式编写通用函数 1）离散傅里叶级数正变换通用函数function xk=dfs(xn,N) n=0:1:N-1; %n的行向量 k=n; %k的行向量 WN=exp(-j*2*pi

14、/N); %WN因子 nk=n*k; %产生一个含nk值的N乘N维矩阵 WNnk=WN.nk; %DFS矩阵 xk=xn* WNnk; %DFS系数行向量2）离散傅里叶级数逆变换通用函数function xn=idfs(xk,N) n=0:1:N-1; %n的行向量 k=n; %k的行向量 WN=exp(-j*2*pi/N); %WN因子 nk=n*k; %产生一个含nk值的N乘N维矩阵 WNnk=WN.(-nk); %DFS矩阵 xn=xk* WNnk/N; %DFS系数行向量例：已知一个周期性矩形序列的脉冲宽度占整个周期的1/4，一个周期的采样点为16点。用傅里叶级数求信号的幅度和相位频谱

15、；求傅里叶级数逆变换的图形，与原信号图形进行比较。MATLAB程序N=16;xn=ones(1,N/4),zeros(1,3*N/4);n=0:N-1;xk=dfs(xn,N);xn1=idfs(xk,N);subplot(2,2,1); stem(n,xn); title(x(n);subplot(2,2,2); stem(n,abs(xn1); title(idfs(|X(k)|);subplot(2,2,3); stem(n,abs(xk); title(|X(k)|);subplot(2,2,4); stem(n,angle(xk); title(arg|X(k)|);程序运行结果如图

16、4-1 图4-12、周期重复次数对序列频谱的影响理论上讲，周期序列不满足绝对可积条件，要对周期序列进行分析，可以先取K个周期进行处理，然后让K无限增大，研究其极限情况。这样可以观察信号序列由非周期到周期变换时，频谱由连续谱逐渐向离散谱过渡的过程。例：已知一个矩形序列的脉冲宽度占整个周期的1/2，一个周期的采样点数为10，用傅立叶级数变换求信号的重复周期数分别为1、4、7、10时的幅度频谱。MATLAB程序：xn=ones(1,5),zeros(1,5);Nx=length(xn);Nw=1000;dw=2*pi/Nw;k=floor(-Nw/2+0.5):(Nw/2+0.5);for r=0:

17、3; K=3*r+1; nx=0:(K*Nx-1); x=xn(mod(nx,Nx)+1); Xk=x*(exp(-j*dw*nx*k)/K; subplot(4,2,2*r+1); stem(nx,x) axis(0,K*Nx-1,0,1.1); ylabel(x(n); subplot(4,2,2*r+2); plot(k*dw,abs(Xk) axis(-4,4,0,1.1*max(abs(Xk); ylabel(X(k);end程序运行结果如图4-2图4-2从上图可以看出，信号序列的周期数越多，则频谱越是向几个频点集中，当信号周期数趋于无穷大时 ，频谱转化为离散谱。五、实验任务1、输入

18、并运行例题程序，熟悉基本指令的使用。2、已知一个信号序列的主值为x(n)=0,1,2,3,2,1,0，显示两个周期的信号序列波形。求解：用傅里叶级数求信号的幅度和相位频谱；求傅里叶级数逆变换的图形，与原信号图形进行比较。六、实验报告1、简述实验目的、原理。2、写出上机调试通过的实验任务的程序并描述其图形曲线。实验五 离散傅里叶变换一、实验目的1、加深对离散傅里叶变换基本概念的理解。2、了解有限长序列傅里叶变换与周期序列傅里叶级数的联系。3、熟悉相关函数的使用方法。二、实验原理有限长序列的傅里叶变换和逆变换对于非周期序列，在实际中常常使用有限长序列。有限长序列x(n)表示为x(n)是非周期序列，

19、但可以理解为某一周期序列的主值序列。由离散傅立叶级数DFS和IDFS引出有限长序列的离散傅立叶正、逆变换关系式。DFT与DFS的关系比较两者的变换对，可以看出两者的区别仅仅是将周期序列换成了有限长序列。有限长序列x(n)可以看作是周期序列的一个周期；反之周期序列可以看作是有限长序列x(n)以N为周期的周期延拓。由于公式非常相似，在程序编写上也基本一致。三、实例1、已知有限长序列x(n)为：x(n)=0,1,2,3,4,5,6,7,8,9，求x(n)的DFT和IDFT。要求1）画出序列傅里叶变换对应的|X(k)|和argX(k)图形。2）画出原信号与傅里叶逆变换IDFTX(k)图形进行比较。MA

20、TLAB程序：xn=0,1,2,3,4,5,6,7,8,9; N=length(xn);n=0:N-1; k=0:N-1;xk=xn*exp(-j*2*pi/N).(n*k);x=(xk*exp(j*2*pi/N).(n*k)/N;subplot(2,2,1); stem(n,xn); title(x(n);subplot(2,2,2); stem(n,abs(x); title(IDFT|X(k)|);subplot(2,2,3); stem(k,abs(xk); title(|X(k)|);subplot(2,2,4); stem(k,angle(Xk); title(arg|X(k)|)

21、;程序运行结果如图5-1：图5-1由上图可看出，与周期序列不同的是，有限长序列本身是仅有N点的离散序列，相当于周期序列的主值部分。因此，其频谱也对应序列的主值部分，是含N点的离散序列。2、有限长序列DFT与周期序列DFS的联系已知周期序列的主值x(n)=0,1,2,3,4,5，求x(n)周期重复次数为4次时的DFS。要求1）画出原主值序列和信号周期序列；2）画出序列傅里叶变换对的图形。MATLAB程序：xn=0,1,2,3,4,5; N=length(xn);n=0:4*N-1; k=0:4*N-1;xn1=xn(mod(n,N)+1);xk=xn1*exp(-j*2*pi/N).(n*k);

22、subplot(2,2,1); stem(xn); title(原主值信号x(n);subplot(2,2,2); stem(n,xn1); title(周期序列信号);subplot(2,2,3); stem(k,abs(xk); title(|X(k)|);subplot(2,2,4); stem(k,angle(xk); title(arg|X(k)|);程序运行结果如图5-2：图5-2 与上一个例题比较，有限长序列x(n)可以看成是周期序列的一个周期，反之，周期序列可以看成是有限长序列以N为周期的周期延拓。频域上的情况也是相同的。从这个意义上说，周期序列只是有限个序列值有意义。四、实验

23、任务1、输入并运行例题程序，熟悉基本指令的使用。2、验证离散傅里叶变换的线性性质。有两个有限长序列分别为x1(n)和x2(n)，长度分别为N1和N2，且y(n)=ax1(n)+bx2(n)，(a，b均为常数)，则该y(n)的N点DFT为 Y(k)=DFTy(n)=aX1(k)+bX2(k) (0=k=N-1) 其中：N=max(N1,N2)，X1(k)和X2(k)分别为x1(n)和x2(n)的N点DFT。已知序列： x1(n)=0,1,2,4 x2(n)=1,0,1,0,1五、实验报告1、简述实验目的、原理。2、写出上机调试通过的实验任务的程序并描述其图形曲线。实验六 快速傅里叶变换一、实验目

24、的1、加深对快速傅里叶变换基本理论的理解。2、了解用MATLAB语言进行快速傅里叶变换的方法。3、掌握常用函数的调用方法。二、实验原理DFT是在时域和频域均为离散序列的变换方法，它适用于有限长序列。但如果按照变换公式进行运算的话，当序列长度很大时，将占用很大的内存空间，运算时间也会很长，无法实时处理问题。快速傅里叶变换是用于提高DFT运算的高速运算方法的统称，FFT是其中的一种，FFT不是一种新的变换形式，它仅仅只是一种快速算法。FFT主要有时域抽取算法和频域抽取算法，基本思想是将一个长度为N的序列分解成多个短序列再进行运算，如基2算法、基4算法等等，从而可以大大缩短运算时间。三、实验用函数1

25、、fft功能：一维基2快速傅里叶变换。调用格式：y=fft(x)：利用FFT算法计算矢量x的离散傅里叶变换。当x为2的幂次方时，采用高速基2FFT算法，否则为稍慢的混合算法。y=fft(x,N)：采用n点FFT。当x的长度小于N时，FFT函数在x的尾部补零，以构成N点数据，当x的长度大于N时，FFT函数在x的尾部截断，以构成N点数据。2、ifft功能：一维基2快速傅里叶逆变换。调用格式：与FFT调用方法相同，只需改换函数名。3、fftshift功能：对FFT的输出进行重新排列，将零频分量移到频谱的中心。调用格式：y=fftshift(x)：对FFT的输出进行重新排列，将零频分量移到频谱的中心。

26、四、实例1、用MATLAB工具箱函数FFT进行频谱分析时需要注意：1）函数fft的返回值y的数据结构的对称性若已知序列x=4,3,2,6,7,8,9,0，求X(k)=DFTx(n)利用函数fft计算，其MATLAB程序如下：N=8;n=0:N-1;xn=4,3,2,6,7,8,9,0;xk=fft(xn)程序运行结果如下：xk = 39.0000 -10.7782 - 6.2929i 0 + 5.0000i 4.7782 + 7.7071i 5.0000 4.7782 - 7.7071i 0 - 5.0000i -10.7782 + 6.2929i由程序运行结果可见，xk的第一行元素对应频率值

27、为0，第五行元素对应频率为莱奎斯特（Nyquist）频率，即标准频率值为1。因此第一行至第五行对应的标准频率为01。而第五行至第八行对应的是负频率，其K(x)值是以Nyquist频率为轴对称。一般情况，对于N点的x(n)序列的FFT是N点的复数序列，其点n=N/2+1对应Nyquist频率，作谱分析时仅取序列X(k)的前一半即可，其后一半序列和前一半是对称的。2）频率计算若N点序列x(n)(n=0,1,N-1)是在采样频率fs(Hz)下获得的。它的FFT也是N点序列，即X(k)（k=0,1,N-1），则第K点对应实际频率值为: f=k*fs/N3）作FFT分析时，幅值大小与FFT选择的点数有关，但不影响分析结果。2、已知信号由15Hz幅值0.5的正弦信号和40Hz幅值2的正弦信号组成，数据采样频率为100Hz,试绘制N=128点DFT的幅频图。MATLAB程序如下：fs=100;N=128;n=0:N-1;t=n/fs;x=0.5*si