《数字信号处理》实验指导书学生版概述文档格式.docx

1、1数组的加减乘除和乘方运算输入A=1 2 3，B=4 5 6，求C=A+B，D=A-B，E=A.*B，F=A./B，G=A.B,并用stem画出A，B，C，D，E，F，G。再输入一些数组，进行类似运算。Stem函数格式可通过help stem命令查询2给出实现单位脉冲序列、单位阶跃序列的MATLAB函数，并用stem函数画出图形。注：（1） 单位脉冲序列可以通过借助MATLAB中的零矩阵函数zeros表示。全零矩阵zeros（1,N）产生一个由N个零组成的列向量，对于有限区间的可以通过以下MATLAB程序表示% 单位抽样序列实现程序k=-30:30;delta=zeros（1,30）,1,ze

2、ros（1,30）;stem（k,delta）（2） 单位阶跃序列可以通过借助MATLAB中的单位矩阵函数ones表示。单位矩阵ones（1,N）产生一个由N个1组成的列向量，对于有限区间的% 单位阶跃序列实现程序uk=zeros（1,30）,ones（1,31）;stem（k,uk）3已知两离散序列分别为和，试用MATLAB绘出它们的波形及的波形。对于离散序列来说，序列相加是将信号对应时间序号的值逐项相加，在这里不能象连续时间信号那样用符号运算来实现，而必须用向量表示的方法，即在MATLAB中离散序列的相加需表示成两个向量的相加，因而参加运算的两序列向量必须具有相同的维数。实现离散序列相加的

3、MATLAB实用子程序如下function f,k=lsxj（f1,f2,k1,k2）%实现f（k）=f1（k）+f2（k）,f1,f2,k1,k2是参加运算的二离散序列及其对应的时间序列向量,f和k为返回的和序列及其对应的时间序列向量k=min（min（k1）,min（k2）:max（max（k1）,max（k2）;%构造和序列长度s1=zeros（1,length（k）;s2=s1; %初始化新向量s1（find（k=min（k1）&（k=min（k2）&=max（k2）=1）=f2;%将f2中在和序列范围内但又无定义的点赋值为零f=s1+s2; %两长度相等序列求和stem（k,f,fi

4、lled）axis（min（min（k1）,min（k2）-1）,（max（max（k1）,max（k2）+1）,（min（f）-0.5）,（max（f）+0.5）%坐标轴显示范围该实验首先建立一个lsxj.m文件，即上段程序，然后存储在默认路径下。该实验参考程序：%求两离散序列之和实现程序f1=-2:2;k1=-2:f2=1 1 1;k2=-1:1;subplot 221;stem（k1,f1）,axis（-3 3 -2.5 2.5）;title（f1k）;subplot 222;stem（k2,f2）,axis（-3 3 -2.5 2.5）;f2ksubplot 223;f,k=lsxj（

5、f1,f2,k1,k2）;%调用自定义的lsxj函数stem（k,f）;fk=f1k+f2（k）4试用MATLAB绘出3中两离散序列乘法与离散序列加法相似，这里参加运算的两序列向量必须具有相同的维数。实现离散时间信号相乘的MATLAB实用子程序如下function f,k=lsxc（f1,f2,k1,k2）f=s1.*s2;同样，在该实验中现建立lsxj.m文件存放在默认路径下，然后再调用该文件实现两序列相乘。该实验参考程序如下%求两离散序列之积实现程序f,k=lsxc（f1,f2,k1,k2）;fk=f1k*f2（k） 实验报告要求1 简述实验目的。2 按实验步骤附上实验程序及结果。实验二

6、用MATLAB进行离散系统的Z域分析一、实验目的掌握应用MATLAB实现IIR、FIR数字滤波器结构的方法，即熟悉应用MATLAB实现IIR、FIR数字滤波器结构的相关函数。二、实验仪器PC机（装有MATLAB软件）三、MATLAB几个信号处理工具箱函数1impz：求H（z）的Z反变换h（n）h,T=impz（B,A,N）h为存放h（n）的列向量，时间变量N存放在列向量T中，当N为标量时，表示T=0：N-1，计算h（n），n=0，1，2，N-1；当N为向量时，T=N，仅计算N指定的整数点上的h（n）。2freqz：求数字滤波器H（z）的频率响应函数H=freqzB,A,w计算由向量w指定的数字

7、频率点上数字滤波器H（z）的频率响应，结果存在H向量中。H,w=freqzB,A,M计算出M个频率点上的频率响应，存放在H向量中，M个频率存放在w向量中，freqz函数自动将这M个频率点均匀设置在频率范围0，之间。若缺省w和M时，函数自动选取512个频率点计算。不带输出向量的freqz函数将自动绘制幅频和相频曲线。3zplane:绘制H（z）的零极点图zplane（z,p）绘制出列向量z中的零点（以符号“”表示）和列向量p中的极点（以符号“”表示）以及参考单位圆，并在多阶零点和极点的右上角标出其阶数。如果z和p为矩阵，则会以不同颜色绘出z和p各列中的零点和极点。zplane（B,A）绘制出系统

8、函数H（z）的零极点图。四、实验内容1利用MATLAB绘制离散系统的零极图对于离散系统其系统函数可由式表示，则系统函数的零点和极点可以用MATLAB的多项式求根函数roots（）来实现，调用该函数的命令格式为：p=roots（A）其中A为待求根的多项式的系数构成的行向量，返回向量p则包含该多项式所有的根位置列向量。值得注意的是：求系统函数零极点时，离散系统的系统函数可能有两种形式，一种是分子和分母多项式按Z的降幂次序排列，另一种是分子和分母多项式按的升幂次序排列。若是以Z的降幂次序排列，则系数向量一定要由多项式的最高幂次开始，一直到常数项，缺项要用0补齐; 若是以按的升幂次序排列，则分子和分母

9、多项式系数向量的维数一定要相同，不足的要用0补齐，不则的零点或极点就可能被漏掉。下面是求系统函数零极点，并绘制其零极点图的MATLAB实用函数ljdt（）,该函数在绘制系统零极点图的同时，还绘出了Z平面的单位圆。function ljdt（A,B）% The function to draw the pole-zero diagram for discrete systemp=roots（A）; %求系统极点q=roots（B）; %求系统零点p=p; %将极点列向量转置为行向量q=q %将零点列向量转置为行向量x=max（abs（p q 1）; %确定纵坐标范围x=x+0.1;y=x; %确

10、定横坐标范围clfhold onaxis（-x x -y y） %确定坐标轴显示范围w=0:pi/300:2*pi;t=exp（i*w）;plot（t） %画单位园axis（squareplot（-x x,0 0） %画横坐标轴plot（0 0,-y y） %画纵坐标轴text（0.1,x,jImztext（y,1/10,Rezplot（real（p）,imag（p）,x） %画极点plot（real（q）,imag（q）,o） %画零点pole-zero diagram for discrete system） %标注标题hold off完成下列实验内容：（1）已知某离散系统的系统函数为试用

11、MATLAB求出该系统的零极点，并画出零极点分布图，判断系统是否稳定。（2） 已知某离散系统的系统函数为试用MATLAB求出该系统的零极点，并画出零极点分布图，求系统的单位冲激响应和幅频响应，并判断系统的是否稳定。2. 利用MATLAB分析离散系统的零极图分布与系统单位响应时域特性的关系已知离散系统的零极分布分别如图2-1所示，其中虚线表示单位圆，试用MATLAB分析系统单位响应的时域特性。图2-1 离散系统的零极点分布图因系统的零极点分布图已知，则系统的系统函数就可知，故可以利用MATLAB函数impz（）函数可求出系统单位响应。对于图2-1所示的系统，其系统函数分别为、实验报告要求2 按实

12、验步骤编写相应的程序，有实验报告中附上实验程序及结果。实验三 傅立叶变换1 熟悉傅立叶变换的各种性质2 熟悉基本信号的频域转换3 熟悉应用FFT对典型信号进行频谱分析的方法4 了解应用FFT进行信号频谱分析过程中可能出现的问题，以便在实际中正确应用FFT实验原理时域信号处理我们已经比较熟悉，信号的频谱函数对于我们却是一个全新的概念。一个信号的时域转换可以通过傅立叶变换（DFT）来完成。有限长序列的DFT是其Z变换在单位圆上的等距采样，或者说是序列傅立叶变换的等距采样，因此可以用于序列的频谱分析。FFT并不是与DFT不同的另一种变换，而是为了减少DFT运算次数的一种快速算法。它是对变换式进行一次

13、次分解，使其成为若干小数点的组合，从而减少运算量。常用的FFT是以2为基数的，其长度N=2L。从频率采样定理知道，N点序列x（n）的N个离散时间傅立叶变换等间隔样本能唯一地重构x（n）。这些单位圆上的N个样本叫做离散傅立叶系数。设，为一周期（具有无限持续时间）序列，则它的主周期为具有有限持续时间的离散傅立叶变换，N点序列的离散傅立叶变换由下式给出：N点序列的离散傅立叶反变换为：在MATLAB中用fft函数实现快速傅立叶变换，如fft（x）实现2点fft变换，fft（x,N） 实现N点fft变换，N必须为2n。实验内容1、MATLAB实现基本序列的离散傅里叶变换计算已知复正弦序列，余弦序列，分别

14、对序列求当时的DFT，并绘出幅频特性曲线，并分析两种N值下DFT是否有差别，及产生原因。程序清单如下：N=16;N1=8;n=0:N-1;k=0:N1-1;x1n=exp（j*pi*n/8）; %产生x1（n）X1k=fft（x1n,N）; %计算N点DFTx1（n）X2k=fft（x1n,N1）; %计算N1点DFTx1（n）x2n=cos（pi*n/8）; %产生x2（n）X3k=fft（x2n,N）; %计算N点DFTx2（n）X4k=fft（x2n,N1）; %计算N1点DFTx2（n）Subplot（2,2,1）;stem（n,abs（X1k）,.axis（0,20,0,20）;yl

15、abel（|X1（k）|16点的DFTx1（n）subplot（2,2,2）;stem（n,abs（X3k）,|X2（k）|16点的DFTx2（n）subplot（2,2,3）;stem（k,abs（X2k）,8点的DFTx1（n）subplot（2,2,4）;stem（k,abs（X4k）,8点的DFTx2（n）N点离散傅里叶变换的一种物理解释就是是以N为周期的周期延拓序列的离散傅立叶级数系数的主值区间序列，即当N=16时， 正好分别是的一个周期，所以的周期延拓序列就是这两个单一频率的正弦序列，其离散傅里叶级数的系数分别如图3.14中的16点的DFTx1（n）和16点的DFTx2（n）所示。

16、而当N=8时， 的半个周期，所以的周期延拓序列就不再是单一频率的正弦序列，而是含有丰富的谐波成分，其离散傅立叶级数的系数与N=16时的差别很大，因此对信号进行谱分析的时候，一定要截取整个周期，否则得到错误的频谱。2、MATLAB验证N点DFT的物理意义已知，绘制相应的幅频和相频曲线，并计算图示N=8和N=16时的DFT。N2=16; % 两种FFT的变换长度k1=0: k2=0:N2-1;w=2*pi*（0:2047）/2048;Xw=（1-exp（-j*4*w）./（1-exp（-j*w）; %对x（n）的频谱函数采样2048个点可以近似的看作是连续的频谱xn=（n=0）&（n4）; %产生

17、x（n）X1k=fft（xn,N1）; %计算N1=8点的X1（k）X2k=fft（xn,N2）; %计算N2=16点的X2（k）subplot（3,2,1）;plot（w/pi,abs（Xw）;xlabel（w/subplot（3,2,2）;plot（w/pi,angle（Xw）;axis（0,2,-pi,pi）;line（0,2,0,0）;subplot（3,2,3）;stem（k1,abs（X1k）,k（=2k/N1）plot（N1/2*w/pi,abs（Xw） %图形上叠加连续频谱的幅度曲线subplot（3,2,4）;stem（k1,angle（X1k）;axis（0,N1,-pi,

18、pi）;line（0,N1,0,0）;） ;ArgX1（k）plot（N1/2*w/pi,angle（Xw） %图形上叠加连续频谱的相位曲线subplot（3,2,5）;stem（k2,abs（X2k）,k（=2k/N2）plot（N2/2*w/pi,abs（Xw）subplot（3,2,6）;stem（k2,angle（X2k）,axis（0,N2,-pi,pi）;line（0,N2,0,0）;ArgX2（k）plot（N2/2*w/pi,angle（Xw）前面的理论已经知道，序列x（n）的N点DFT的物理意义是对在0，2上进行N点的等间隔采样。由结果可以直观的看到X（k）与之间的采样关系。

19、3、MATLAB实现快速卷积用快速卷积法计算下面两个序列的卷积，给出实验结果并进行分析。快速卷积计算框图如图所示：图 快速卷积框图FFT的变换长度L必须满足LN+M-1，输出y（n）才等于x（n）和h（n）的线性卷积。计算两个序列的卷积时，也可以直接调用函数conv来计算，因为MATLAB中的计时比较粗糙，所以只有当N和M较大的时候，才能比较两种方法的执行时间快慢。1. 简述实验目的。2. 按实验步骤附上实验程序及结果，并给出相应的分析。实验四 IIR及FIR滤波器的MATLAB实现 一、实验目的三、实验内容编写MATLAB 程序，实现以下题目：（一）IIR滤波器结构实现1在MATLAB中，用

20、filter函数实现IIR的直接形式用直接型实现系统函数为的IIR数字滤波器，求单位脉冲响应和单位阶跃信号的输出。2在MATLAB中给定级联型系统函数，由扩展函数casfiltr实现IIR的级联形式系统函数，用级联型结构实现。b0=4;B=1,1,0;1,-1.4142136,1;A=1,-0.5,0;1,0.9,0.81;N=60;delta=impseq（0,0,N）;h=casfiltr（b0,B,A,delta）;x=ones（1,5）,zeros（1,N-5）;y=casfiltr（b0,B,A,x）;subplot（1,2,1）;stem（h）;直接型h（n）subplot（1,2

21、,2）;stem（y）;直接型y（n）扩展函数casfiltr程序：function y=casfiltr（b0,B,A,x）;% IIR各FIR滤波器的级联型的实现% % y=casfiltr（b0,B,A,x）;% y=输出序列% b0=级联型的增益系数% B=包含各bk和K乘3维实系数矩阵% A=包含各ak的K乘3维实系数矩阵% x=输入序列%K,L=size（B）;N=length（x）;w=zeros（K+1,N）;w（1,:）=x;for i=1:1:K w（i+1,:）=filter（B（i,:）,A（i,:）,w（i,:）;endy=b0 * w（K+1,:单位抽样序列的生成函

22、数impseq.mfunction x,n=impseq（n0,ns,nf）n=ns:nf;x=（n-n0）=0; %序列的起点为ns，终点为nf，在n=n0点处生成一个单位脉冲。3在MATLAB中用给定并联系统函数，由扩展函数parfiltr实现IIR的并联形式。用并联型实现系统函数为的滤波器。程度清单如下：C=0;B=-14.75,-12.90;24.50,26.82;A=1,-7/8,3/32;1,-1,0.5;h=parfiltr（C,B,A,delta）;y=parfiltr（C,B,A,x）;并联型h（n）并联型y（n）扩展函数parfiltr:function y=parfiltr（C,B,A,x）;% IIR滤波器的并联型实现% % y=parfiltr（C,B,A,x）;%C=