《数字信号处理实验指导书.doc》由会员分享，可在线阅读，更多相关《数字信号处理实验指导书.doc（29页珍藏版）》请在taowenge.com淘文阁网|工程机械CAD图纸|机械工程制图|CAD装配图下载|SolidWorks_CaTia_CAD_UG_PROE_设计图分享下载上搜索。
1、数字信号处理实验指导书山东大学控制学院生物医学工程专业刘忠国2015-2-10数字信号处理实验目录实验一 离散时间信号与系统分析3实验二 离散时间信号与系统的Z变换分析 7实验三IIR滤波器的设计与信号滤波 13实验四 用窗函数法设计FIR数字滤波器15实验五 用FFT作谱分析17实验六综合实验19附 录： 各实验参考程序20实验一 离散时间信号与系统分析一、实验目的1掌握离散时间信号与系统的时域分析方法。2掌握序列傅氏变换的计算机实现方法，利用序列的傅氏变换对离散信号、系统及系统响应进行频域分析。 3熟悉理想采样的性质，了解信号采样前后的频谱变化，加深对采样定理的理解。 二、实验原理1离散时
2、间系统一个离散时间系统是将输入序列变换成输出序列的一种运算。若以来表示这种运算，则一个离散时间系统可由下图来表示： 图 离散时间系统即输出与输入之间关系用下式表示 （1）离散时间系统中最重要、最常用的是线性时不变系统。2离散时间系统的单位脉冲响应设系统输入，系统输出的初始状态为零，这时系统输出用表示，即，则称为系统的单位脉冲响应。可得到： （2）该式说明线性时不变系统的响应等于输入序列与单位脉冲序列的卷积。3连续时间信号的采样采样是从连续信号到离散时间信号的过渡桥梁，对采样过程的研究不仅可以了解采样前后信号时域和频域特性发生的变化以及信号内容不丢失的条件，而且有助于加深对拉氏变换、傅氏变换、Z
3、变换和序列傅氏变换之间关系的理解。对一个连续时间信号进行理想采样的过程可以表示为信号与一个周期冲激脉冲的乘积，即： （3）其中，是连续信号的理想采样，是周期冲激脉冲设模拟信号，冲激函数序列以及抽样信号的傅立叶变换分别为、和，即根据连续时间信号与系统中的频域卷积定理，式（3）表示的时域相乘，变换到频域为卷积运算，即 其中 由此可以推导出 由上式可知，信号理想采样后的频谱是原来信号频谱的周期延拓，其延拓周期等于采样频率。根据香农定理，如果原信号是带限信号，且采样频率高于原信号最高频率的2倍，则采样后的离散序列不会发生频谱混叠现象。4有限长序列的分析对于长度为N的有限长序列，我们只观察、分析在某些频
4、率点上的值。 一般只需要在之间均匀的取M个频率点，计算这些点上的序列傅立叶变换： 其中，。是一个复函数，它的模就是幅频特性曲线。三、主要实验仪器及材料微型计算机、Matlab软件。四、实验内容1知识准备认真复习离散信号与系统、单位脉冲响应、抽样定理等有关内容，阅读本实验原理与方法。a) 2编制信号产生子程序，用于产生实验中要用到的信号序列，并分析幅频响应（1）单位脉冲序列单位脉冲序列（2）系统单位脉冲响应序列（3）理想采样信号序列对信号进行理想采样，可以得到一个理想的采样信号序列，。其中为幅度因子，a是衰减因子，W是频率，为采样周期。这几个参数要在实验过程中输入，以产生不同的。首先产生理想采样
5、信号序列，使200, =50,=260。然后改变参数A1，=0.4, =700，产生理想采样信号序列。3离散信号、系统和系统响应的分析观察信号xb(n)和系统hb(n)的时域和频域特性；利用线性卷积求信号通过系统以后的响应。比较系统响应和信号的时域和幅频特性。绘出图形。4分析理想采样信号序列的特性产生理想采样信号序列，使：（1）首先选用采样频率为1000Hz，T=1/1000，观察所得理想采样信号的幅频特性，在折叠频率以内和给定的理想幅频特性无明显差异，并作记录。（2）改变采样频率为300Hz，T=1/300，观察所得理想采样信号的幅频特性曲线的变化，并作记录。（3）进一步减小采样频率为200
6、Hz，T=1/200，观察频谱混叠现象是否明显存在，说明原因，并记录此时的幅频特性曲线。5. 卷积定律的验证。采用参数444.128, =50,=50， T=1/1000，将和系统的傅氏变换相乘，直接求得，将得到的幅频特性曲线和先求后再求得的幅频特性曲线进行比较，观察二者有无差异。验证卷积定律。五、思考题1线性时不变系统的输出的长度与输入及系统的单位冲激响应的长度有什么关系？2. 对信号进行理想抽样时，抽样频率不同，相应理想采样序列傅立叶变换频谱的数字频率度量是否都相同？它们所对应的模拟频率是否相同? 为什么？六、实验报告要求 1简述实验原理及目的。2. 总结在上机实验内容中要求比较时域、幅频
7、曲线差异部分内容的结果。3.总结实验所得主要结论。4.简要回答思考题。实验二 离散时间信号与系统的Z变换分析一、 实验目的1、熟悉离散信号Z变换的原理及性质2、熟悉常见信号的Z变换3、了解正/反Z变换的MATLAB实现方法4、了解离散信号的Z变换与其对应的理想抽样信号的傅氏变换和拉氏变换之间的关系5、了解利用MATLAB实现离散系统的频率特性分析的方法二、实验原理1、正/反Z变换Z变换分析法是分析离散时间信号与系统的重要手段。如果以时间间隔对连续时间信号f(t)进行理想抽样，那么，所得的理想抽样信号为：理想抽样信号的双边拉普拉斯变换Fd (s)为：若令 ， ， 那么的双边拉普拉斯变换Fd (s
8、)为：则离散信号f(k)的Z变换定义为： 从上面关于Z变换的推导过程中可知，离散信号f(k)的Z变换F(z)与其对应的理想抽样信号的拉氏变换Fd (s)之间存在以下关系： 同理，可以推出离散信号f(k)的Z变换F(z)和它对应的理想抽样信号的傅里叶变换之间的关系为 如果已知信号的Z变换F(z)，要求出所对应的原离散序列f(k)，就需要进行反Z变换，反Z变换的定义为：其中，C为包围的所有极点的闭合积分路线。在MATLAB语言中有专门对信号进行正反Z变换的函数ztrans( ) 和itrans( )。其调用格式分别如下：l F=ztrans( f ) 对f(n)进行Z变换，其结果为F(z)l F=
9、ztrans(f,v) 对f(n)进行Z变换，其结果为F(v)l F=ztrans(f,u,v) 对f(u)进行Z变换，其结果为F(v)l f=iztrans ( F ) 对F(z)进行Z反变换，其结果为f(n)l f=iztrans(F,u) 对F(z)进行Z反变换，其结果为f(u)l f=iztrans(F,v,u ) 对F(v)进行Z反变换，其结果为f(u)注意： 在调用函数ztrans( )及iztrans( )之前，要用syms命令对所有需要用到的变量（如t,u,v,w）等进行说明，即要将这些变量说明成符号变量。例用MATLAB求出离散序列 的Z变换MATLAB程序如下：syms k
10、 zf=0.5k; %定义离散信号Fz=ztrans(f) %对离散信号进行Z变换运行结果如下：Fz = 2*z/(2*z-1)例已知一离散信号的Z变换式为 ，求出它所对应的离散信号f(k)MATLAB程序如下：syms k zFz=2* z/(2*z-1); %定义Z变换表达式fk=iztrans(Fz,k) %求反Z变换运行结果如下：fk = (1/2)k例：求序列的Z变换.clc;clear allsyms n hn=sym(kroneckerDelta(n, 1) + kroneckerDelta(n, 2)+ kroneckerDelta(n, 3)Hz=ztrans(hn)Hz=s
11、implify(Hz)2、离散系统的频率特性同连续系统的系统函数H(s)类似，离散系统的系统函数H(z)也反映了系统本身固有的特性。对于离散系统来说，如果把其系统函数H(z)中的复变量z换成（其中），那么所得的函数就是此离散系统的频率响应特性，即离散时间系统的频率响应为：其中， 称为离散系统的幅频特性，称为系统的相频特性。同连续系统一样，离散时间系统的幅频特性也是频率的偶函数，相频特性也是频率的齐函数。由于是频率w的周期函数，所以离散系统的频率响应特性也是频率w的周期函数，其周期为，或者角频率周期为。实际上，这就是抽样系统的抽样频率，而其中的Ts则是系统的抽样周期。频率响应呈现周期性是离散系统
12、特性区别于连续系统特性的重要特点。因此，只要分析在范围内的情况，便可分析出系统的整个频率特性。函数来表示离散系统的频率响应特性， 表示幅频特性，而相频特性仍用来表示。应该特别注意的是，虽然这里的变量w仍然称为频率变量，但是它已经不是原来意义上的角频率概念，而实际上是表示角度的概念。我们称之为数字频率。它与原来角频率的关系为：。也就是说，根据离散系统的系统函数H(z)，令其中的，并且代入0范围内不同的频率值（实际上是角度值），就可以逐个计算出不同频率时的响应，求出离散系统的频率响应特性。再利用离散系统频率特性的周期性特点（周期为2p），求出系统的整个频率特性。离散系统的幅频特性曲线和相频特性曲线
13、能够直观地反映出系统对不同频率的输入序列的处理情况。在函数随w的变换关系中，在w=0附近，反映了系统对输入信号低频部分的处理情况，而在w=p附近，则反映了系统对输入信号高频部分的处理情况。一般来说，分析离散系统频率响应特性就要绘制频率响应曲线，而这是相当麻烦的。虽然可以通过几何矢量法来定性画出频率响应特性曲线，但一般来说这也是很麻烦的。值得庆幸的是，MATLAB为我们提供了专门用于求解离散系统频率响应的函数freqz() ，其调用格式如下：l H,w=freqz(B,A,N) 其中，B和A分别是表示待分析的离散系统的系统函数的分子，分母多项式的向量，N为正整数，返回向量H则包含了离散系统频率响
14、应函数在范围内的N个频率等分点的值。向量w则包含范围内的N个频率等分点。在默认情况下N=512。l H,w=freqz(B,A,N,whole) 其中，B，A和N的意义同上，而返回向量H包含了频率响应函数在范围内N个频率等分点的值。由于调用freqz()函数只能求出离散系统频率响应的数值，不能直接绘制曲线图，因此，我们可以先用freqz()函数求出系统频率响应的值，然后再利用MATLAB的abs()和angle()函数以及plot()命令，即可绘制出系统在或范围内的幅频特性和相频特性曲线。例若离散系统的系统函数为，请用MATLAB计算频率范围内10个等分点的频率响应的样值。 MATLAB程序如
15、下：A=1 0; %分母多项式系数向量B=1 -0.5; %分子多项式系数向量H,w=freqz(B,A,10) %求出对应范围内10个频率点的频率响应样值运行结果如下：H = 0.5000 0.5245 + 0.1545i 0.5955 + 0.2939i 0.7061 + 0.4045i 0.8455 + 0.4755i 1.0000 + 0.5000i 1.1545 + 0.4755i 1.2939 + 0.4045i 1.4045 + 0.2939i 1.4755 + 0.1545iw = 0 0.3142 0.6283 0.9425 1.2566 1.5708 1.8850 2.19
16、91 2.5133 2.8274例用MATLAB计算前面离散系统在频率范围内200个频率等分点的频率响应值，并绘出相应的幅频特性和相频特性曲线。MATLAB程序如下：A=1 0;B=1 -0.5;H,w=freqz(B,A,200);H,w=freqz(B,A,200,whole); %求出对应范围内200个频率点的频率响应样值HF=abs(H); %求出幅频特性值HX=angle(H); %求出相频特性值subplot(2,1,1);plot(w,HF) %画出幅频特性曲线subplot(2,1,2);plot(w,HX) %画出相频特性曲线运行结果如下：运行结果分析：从该系统的幅频特性曲线
17、可以看出，该系统呈高通特性，是一阶高通滤波器。三、 实验内容1 求出下列离散序列的Z变换 2 已知下列单边离散序列的z变换表达式，求其对应的原离散序列。 3.已知离散系统的系统函数H (z)如下，请绘出系统的幅频和相频特性曲线，并说明系统的作用 4. 已知描述离散系统的差分方程为：请绘出系统的幅频和相频特性曲线，并说明系统的作用。四、 预习要求1、 熟悉正反z变换的意义及用MATLAB软件实现的方法2、 熟悉离散系统的频率响应特性及用MATLAB软件实现的方法3、 编写MATLAB程序五、 实验报告要求1、 简述实验目的及实验原理2、 计算相应z变换或反z变换的理论值，并与实验结果进行比较3、
18、 记录离散系统的频率响应特性曲线，分析系统作用4、 写出程序清单6、 收获与建议实验三IIR滤波器的设计与信号滤波、实验目的（）熟悉用双线性变换法设计IIR数字滤波器的原理与方法。（）掌握数字滤波器的计算机仿真方法。（）通过观察对实际心电图信号的滤波作用，获得数字滤波的感性知识。、实验原理利用双线性变换设计IIR滤波器（只介绍巴特沃斯数字低通滤波器的设计），首先要设计出满足指标要求的模拟滤波器的传递函数，然后由通过双线性变换可得所要设计的IIR滤波器的系统函数。如果给定的指标为数字滤波器的指标，则首先要转换成模拟滤波器的技术指标，这里主要是边界频率的转换，对指标不作变化。边界频率的转换关系为。
19、接着，按照模拟低通滤波器的技术指标根据相应设计公式求出滤波器的阶数和截止频率；根据阶数查巴特沃斯归一化低通滤波器参数表，得到归一化传输函数；最后，将代入去归一，得到实际的模拟滤波器传输函数。之后，通过双线性变换法转换公式，得到所要设计的IIR滤波器的系统函数。利用所设计的数字滤波器对实际的心电图采样信号进行数字滤波器。、实验步骤及内容（）复习有关巴特沃斯模拟滤波器的设计和用双线性变换法设计IIR数字滤波器的内容，用双线性变换法设计一个巴特沃斯IIR低通数字滤波器。设计指标参数为：在通带内频率低于时，最大衰减小于；在阻带内频率区间上，最小衰减大于。（）以为采样间隔，绘制出数字滤波器在频率区间上的
20、幅频响应特性曲线。（）用所设计的滤波器对实际心电图信号采样序列（实验数据在后面给出）进行仿真滤波处理，并分别绘制出滤波前后的心电图信号波形图，观察总结滤波作用与效果。（）编写程序完成各部分实验内容。4、实验用MATLAB函数介绍buttord(); butter(); bilinear(); freqz(); freqs(); filter(); figure(); plot(); stem(); abs();title(); xlabel(); ylabel(); text(); hold on; axis(); grid on; subplot();、思考题（）用双线性变换法设计数字滤波器
21、过程中，变换公式 中T的取值， 对设计结果有无影响? 为什么?（）如果用脉冲响应不变法设计该IIR数字低通滤波器，程序如何改动？ 、实验报告要求（）简述实验目的及实验原理。（）编程实现各实验内容，列出实验清单及说明。（）由绘制的特性曲线及设计过程简述双线性变换法的特点。（）对比滤波前后的心电图信号波形，说明数字滤波器的滤波过程与滤波作用。（）简要回答思考题。、心电图信号采样序列人体心电图信号在测量过程中往往受到工业高频干扰，所以必须经过低通滤波处理后，才能作为判断心脏功能的有用信息。下面给出的数据是一实际心电图信号采样序列样本，其中存在高频干扰。本实验中，以作为输入序列，滤除其中的干扰成分。实
22、验四 用窗函数法设计FIR数字滤波器、实验目的（）掌握用窗函数法设计FIR数字滤波器的原理和方法。（）熟悉线性相位FIR数字滤波器特性。（）了解各种窗函数对滤波特性的影响。、实验原理如果所希望的滤波器的理想频率响应函数为，则其对应的单位脉冲响应为窗函数设计法的基本原理是用有限长单位脉冲响应逼近。由于往往是无限长序列，且是非因果的，所以用窗函数将截断，并进行加权处理，得到：就作为实际设计的FIR数字滤波器的单位脉冲响应序列，其频率响应函数为式中，为所选窗函数的长度。这种对理想单位取样响应的加窗处理对滤波器的频率响应会产生以下三点影响：（1）使理想特性不连续的边沿加宽，形成一过渡带，过渡带的宽度取
23、决于窗函数频谱的主瓣宽度。（2）在过渡带两旁产生肩峰和余振，它们取决于窗函数频谱的旁瓣；旁瓣越多，余振也越多；旁瓣相对值越大，肩峰则越强。（3）增加截断长度，只能缩小窗函数频谱的主瓣宽度而不能改变旁瓣的相对值；旁瓣与主瓣的相对关系只决定于窗函数的形状。因此增加N，只能相对应减小过渡带宽。而不能改变肩峰值。肩峰值的大小直接决定通带内的平稳和阻带的衰减，对滤波器性能有很大关系。例如矩形窗的情况下，肩峰达8.95%，致使阻带最小衰减只有21分贝，这在工程上往往是不够的。怎样才能改善阻带的衰减特性呢？只能从改善窗函数的形状上找出路，所以希望的窗函数频谱中应该减少旁瓣，使能量集中在主瓣，这样可以减少肩峰
24、和余振，提高阻带衰减。而且要求主瓣宽度尽量窄，以获得较陡的过渡带，然而这两个要求总不能同时兼得，往往需要用增加主瓣宽度带换取较大的阻带衰减，于是提出了海明窗、汉宁窗、布莱克曼窗、凯塞窗、切比雪夫窗等窗函数。所以，用窗函数法设计的滤波器性能取决于窗函数的类型及窗口长度的取值。设计过程中，要根据对阻带最小衰减和过渡带宽度的要求选择合适的窗函数类型和窗口长度。设待求滤波器的过渡带用表示，它近似等于窗函数主瓣宽度。因过渡带近似与窗口长度成反比，决定于窗口形式。例如，矩形窗A=4，海明窗A=8等。按照过渡带及阻带衰减情况，选择窗函数形式。原则是在保证阻带衰减满足要求的情况下， 尽量选择主瓣窄的窗函数。这
25、样选定窗函数类型和窗口长度后，求出单位脉冲响应，再求出。是否满足要求，要进行验算。一般在的尾部加零使长度满足的整数次幂，以便用FFT计算。如果要观察细节，补零点数增多即可。如果不满足要求，则要重新选择窗函数类型和长度，再次验算，直至满足要求。如果要求线性相位特性，则还必须满足：根据上式中的正、负号和长度的奇偶性又将线性相位FIR滤波器分成四类。要根据所设计的滤波特性正确选择其中一类。例如，要设计线性相位低通特性，可选择一类， 而不能选一类。 、实验步骤及内容（）根据下列技术指标，设计一个线性相位的FIR数字低通滤波器。通带截止频率，通带允许波动；阻带截止频率，阻带衰减。 （）写出（）中理想低通
26、滤波器的频率响应和单位脉冲响应。（）写出（）中所设计的低通滤波器的单位脉冲响应；并调用fir1()函数得到所设计的低通滤波器的单位脉冲响应，调用fft()函数进行频响验证。打印输出各部分结果。（）编程验证窗长和窗形状对实际滤波器性能的影响。如要求用窗函数法设计一个线性相位FIR数字低通滤波器，用理想低通滤波器作为逼近滤波器，截止频率，用四种窗函数（矩形窗，汉宁窗（升余弦窗），哈明窗（改进的升余弦窗），布莱克曼窗）设计该滤波器，选择窗函数的长度两种情况。4、实验用MATLAB函数介绍fir1(); fft(); freqz(); boxcar(); hamming(); hanning(); b
27、lackman(); sin();figure(); plot(); stem(); abs();angle(); title(); xlabel(); ylabel(); text(); hold on; axis(); grid on; subplot();、思考题如果要求用窗函数法设计带通滤波器， 且给定上、 下边带截止频率为和，试求理想带通的单位脉冲响应。、实验报告要求（）简述实验目的及实验原理。（）编程实现各实验内容，列出实验清单及说明。（）总结窗函数法设计FIR数字滤波器的步骤。（）简要回答思考题。实验五 用FFT作谱分析、实验目的（）进一步加深DFT算法原理和基本性质的理解（因为
28、FFT只是DFT的一种快速算法，所以FFT的运算结果必然满足DFT的性质）（）熟悉FFT算法原理及子程序的应用。（）掌握用FFT对连续信号和时域离散信号进行频谱分析的基本方法。了解可能出现的分析误差和原因，以便在实际中正确应用FFT。、实验原理如果用FFT对模拟信号进行谱分析，首先要把模拟信号转换成数字信号，转换时要求知道模拟信号的最高截止频率，以便选择满足采样定理的采样频率。一般选择采样频率是模拟信号中最高频率的34倍。另外要选择对模拟信号的观测时间，如果采样频率和观测时间确定，则采样点数也确定了。这里观测时间和对模拟信号进行谱分析的分辨率有关，最小的观测时间和分辨率成倒数关系。要求选择的采
29、样点数和观测时间大于它的最小值。用FFT作谱分析时，要求做FFT的点数服从的整数幂，这一点在上面选择采样点数时可以考虑满足，即使满足不了，可以通过在序列尾部加完成。如果要进行谱分析的模拟信号是周期信号，最好选择观测时间是信号周期的整数倍。如果不知道信号的周期，要尽量选择观测时间长一些，以减少截断效应的影响。用FFT对模拟信号作谱分析是一种近似的谱分析。首先一般模拟信号（除周期信号外）的频谱是连续频谱，而用FFT作谱分析得到的是数字谱，因此应该取FFT的点数多一些，用它的包络作为模拟信号的近似谱。另外，如果模拟信号不是严格的带限信号，会因为频谱混叠现象引起谱分析的误差，这种情况下可以预先将模拟信
30、号进行预滤，或者尽量将采样频率取高一些。一般频率混叠发生在折叠频率附近，分析时要注意因频率混叠引起的误差。最后要注意一般模拟信号是无限长的，分析时要截断，截断的长度和分辨率有关，但也要尽量取长一些，取得太短因截断引起的误差会很大。举一个极端的例子，一个周期性正弦波，如果所取观察时间太短，例如取小于一个周期，它的波形和正弦波相差太大，肯定误差很大，但如果取得长一些，即使不是周期的整倍数，这种截断效应也会小一些。、实验步骤及内容（）复习DFT的定义、性质和用DFT作谱分析的有关内容。（）复习FFT算法原理与编程思想。（）编制信号产生程序，产生以下典型信号供谱分析用：式中频率自己选择；（）分别以变换
31、区间，对进行FFT，画出相应的幅频特性曲线。（）分别以变换区间，对，进行FFT，画出相应的幅频特性曲线。（）分别以变换区间，对进行FFT，画出相应的幅频特性曲线。（）分别对模拟信号选择采样频率和采样点数。对，周期，频率自己选择，采样频率，观测时间，采样点数用计算。对，选择采样频率，采样点数为，。（）分别将模拟信号转换成序列，用，表示，再分别对它们进行FFT，并画出相应的幅频特性曲线。4、实验用MATLAB函数介绍fft(); figure(); plot(); stem(); abs();title(); xlabel(); ylabel(); text(); hold on; axis();
32、 grid on; subplot(); sin(); cos(); 等。、思考题（）在N=8时，和的幅频特性会相同吗? 为什么? N=16呢?（）如果周期信号的周期预先不知道，如何用FFT进行谱分析?、实验报告要求（）简述实验目的及实验原理。（）编程实现各实验内容，列出实验清单及说明。（）将实验结果和理论分析结果进行比较，分析说明误差产生的原因以及用FFT作谱分析时有关参数的选择方法。并总结实验所得的主要结论。（）简要回答思考题。实验六综合实验、实验目的（）给定设计指标，能选用多种方法设计出符合指标要求的滤波器。（）掌握各种设计方法的原理，步骤以及特点。（）培养学生实际应用及综合设计的能力。
33、、实验内容（）设计一个数字带通滤波器，要求通带范围为0.25 rad到0.45 rad，通带最大衰减dB，0.15 rad以下和0.55 rad以上为阻带，阻带最小衰减dB。要求用IIR和FIR的各种方法设计数字滤波器，并绘制滤波器的幅频特性、相频特性、零极点分布，单位脉冲响应，阶跃响应，滤波器系数、群时延等特性。（）假设一信号包含有用信号和干扰信号，有用信号的频带范围是，干扰信号的频带范围是。现要求设计一数字滤波器，指标是在频率范围中幅度失真为；在时，衰减大于；分别用FIR和IIR两种滤波器进行滤除干扰，最后进行比较。、实验用MATLAB函数介绍（略）、实验报告要求（）简述实验目的。（）编程
34、实现各实验内容，列出实验清单及说明。（）总结实验结论并写出实验收获或感想。附：各实验参考程序实验一 离散时间系统分析参考程序四、2.（1）,（2）的参考程序：%单位脉冲序列的时域和幅频特性%在MatLab中，这一函数可以用zeros函数实现：n=1:50; x=zeros(1,50); %MatLab中数组的下标从1开始x(1)=1;close all;subplot(3,1,1);stem(x);title(单位冲击信号序列x(n);k=-24:25;X=x*(exp(-j*pi/25).(n-1)*k);magX=abs(X);subplot(3,1,2);stem(magX);title
35、(单位冲击信号序列的幅度谱);angX=angle(X);subplot(3,1,3);stem(angX);title(单位冲击信号序列的相位谱);%以下是的时域和幅频特性n=1:50;x=zeros(1,50);x(1)=1;x(2)=2.5;x(3)=2.5;x(4)=1;close all;subplot(3,1,1);stem(x);title(系统单位脉冲响应信号序列);k=-24:25;X=x*(exp(-j*pi/25).(n-1)*k);magX=abs(X);subplot(3,1,2);stem(magX);title(系统频率响应的幅度谱);angX=angle(X);
36、subplot(3,1,3);stem(angX);title(系统频率响应的相位谱)%以下矩形脉冲序列的时域和幅频特性(没有要求做，仅作参考)n=1:50; x=sign(sign(10-n)+1);close all;subplot(3,1,1);stem(x);title(矩形脉冲序列);k=-24:25;X=x*(exp(-j*pi/25).(n-1)*k);magX=abs(X);subplot(3,1,2);stem(magX);title(矩形脉冲序列傅立叶变换的幅度谱);angX=angle(X);subplot(3,1,3);stem(angX);title(矩形脉冲序列傅立
37、叶变换的相位谱);四、2.（3）的参考程序：n=1:50; %定义序列的长度是50A=444.128; %设置信号有关的参数a=50*sqrt(2.0)*pi;T=0.001; %采样率 w0=50*sqrt(2.0)*pi; %符号在MatLab中不能输入，用w代替x=A*exp(-a*(n-1)*T).*cos(w0*(n-1)*T); %pi是MatLab中定义的%close allfiguresubplot(3,1,1);stem(n-1,x); %绘制x（n）的图形title(理想采样信号序列x(n); %设置结果图形的标题k=-24:25;%W=(pi/25)*k;X=x*(exp
38、(-j*pi/25).(n-1)*k);magX=abs(X); %绘制x（n）的幅度谱subplot(3,1,2);stem(n-1,magX);title(理想采样信号序列的幅度谱);angX=angle(X); %绘制x（n）的相位谱subplot(3,1,3);stem(n-1,angX);title(理想采样信号序列的相位谱);四、3. 的参考程序（要用代替）：%卷积计算%在MatLab中提供了卷积函数conv，即y=conv（n，h），调用十分方便。%信号和系统单位脉冲响应的卷积n=1:50;hb=zeros(1,50);hb(1)=1;hb(2)=2.5;hb(3)=2.5;hb
39、(4)=1;close all;subplot(3,1,1);stem(hb);title(系统单位脉冲响应hbn);m=1:50;T=0.001;A=444.128;a=50*sqrt(2.0)*pi;w0=50*sqrt(2.0)*pi;x=A*exp(-a*(m-1)*T).*cos(w0*(m-1)*T);subplot(3,1,2);stem(x);title(输入信号xn);y=conv(x,hb);subplot(3,1,3);stem(y);title(输出信号yn);四、5. 的参考程序%参考程序clearhb=zeros(1,50);hb(1)=1;hb(2)=2.5;hb
40、(3)=2.5;hb(4)=1;m=1:50;T=0.001;A=444.128;a=50*sqrt(2.0)*pi;w0=50*sqrt(2.0)*pi;x=A*exp(-a*(m-1)*T).*cos(w0*(m-1)*T);n=1:50;k=-24:25;X=x*(exp(-j*pi/25).(n-1)*k);magX=abs(X);%绘制信号x(n)傅立叶变换的幅度谱subplot(3,2,1);stem(magX);title(输入信号傅立叶变换的幅度谱);angX=angle(X); %绘制x(n)的相位谱subplot(3,2,2);stem(angX);title(输入信号傅立
41、叶变换的相位谱);Hb=hb*(exp(-j*pi/25).(n-1)*k);magHb=abs(Hb); %绘制hb(n)的幅度谱subplot(3,2,3);stem(magHb);title(系统单位频率响应的幅度谱);angHb=angle(Hb); %绘制hb(n)的相位谱subplot(3,2,4);stem(angHb);title(系统频率响应的相位谱);n=1:99;k=1:99;y=conv(x,hb);Y=y*(exp(-j*pi/25).(n-1)*k);magY=abs(Y); %绘制y(n)的幅度谱subplot(3,2,5);stem(magY);title(输出
42、信号傅立叶变换的幅度谱);angY=angle(Y); %绘制y(n)的相位谱subplot(3,2,6);stem(angY);title(输出信号傅立叶变换的相位谱);%将以下验证的结果显示figureXHB=X.*Hb;subplot(4,1,1);stem(abs(XHB);title(x(n)的幅度谱与hb(n)的幅度谱相乘);axis(0,50,0,8000)subplot(4,1,2);stem(abs(Y);title(y(n) 傅立叶变换的幅度谱);axis(0,50,0,8000)angXHB=angle(XHB); %绘制y(n)的相位谱subplot(4,1,3);% stem(unwrap(angXHB);stem(angXHB);axis(0,50,0,6)title(输出x(n)信号的频谱与hb(n)的频谱相乘的相位谱);angY=angle(Y); %绘制y(n)的相位谱subplot(4,1,4);stem(angY);axis(0,50,0,6)% stem(unwrap(angY);title(输出y(n)信号傅立叶变换的相位谱);