DSP实验三和实验四

实验三离散时间系统及响应[实验目的]1、熟悉离散时间系统的时域和频域分析方法。

2、掌握稳定系统性的判别。

3、掌握利用MATLAB求解零极状态、差分方程和频率响应的方法。

[实验仪器]计算机、Matlab6.5（或更高版本）软件[实验参考书] ：自编实验指导书与本实验有关的MA TLAB函数：1. conv.m用来实现两个离散序列的线性卷积。

其调用格式是：y=conv(x,h)2．filter.m求离散系统的输出y(n) 。

若系统的h(n) 已知，可用conv.m文件可求出y(n) ；若系统的H(z) 已知，可用filter可求出y(n)，调用格式是:y=filter(b, a, x)；%其中x, y, a 和b都是向量。

3．impz.m在H(z) 已知情况下, 求系统的单位抽样响应h(n)。

调用格式是：h = impz(b, a, N)或[h,t]=impz(b,a,N)N是所需的的长度。

前者绘图时n从1开始，而后者从0开始。

4．freqz.m在H(z) 已知情况下, 求系统的频率响应。

基本的调用格式是：[H,w]=freqz(b,a,N,'whole',Fs)N是频率轴的分点数，建议N为2的整次幂；w是返回频率轴座标向量，绘图用；Fs是抽样频率，若Fs＝1，频率轴给出归一化频率；’whole’指定计算的频率范围是从0～FS,缺省时是从0～FS/2.5.zplane.m文件可用来显示离散系统的极－零图。

其调用格式是：zplane(z,p), 或zplane(b,a),前者是在已知系统零点的列向量z和极点的列向量p的情况下画出极－零图，后者是在仅已知H(Z)的A(z)、B(z) 的情况下画出极－零图。

[实验内容]1、序列的卷积运算。

下例是一个指数函数的卷积，结果如图所示。

clear all;fs=20;n =0:49;x=exp(-n/fs);subplot(211),stem(x)y=conv(x,x);subplot(212);stem(y);grid;05101520253035404550仿照上例，现在假设存在两个序列，()[,,,,,,]h n=54321,利用计x n=9121041215()[,,,,]算出两者的卷积，并画出图形。

２、时域离散系统和系统响应分析。

输入下列程序，观察实验结果。

N=16; n=0:N-1; %定义抽样的长度x=sin(2*pi*n/64)+sin(20*pi*n/64);%输入序列a=[1 -0.25]; b=[0.5 0.45 0.35]; %b 和a 分别是系统函数分子和分母多项式系数 subplot(221); zplane(b,a);%画出系统的零极点图 y=filter(b,a,x); %求系统的响应subplot(222 );stem(n,y) ; %画出系统响应的火柴杆图 [H,w]=freqz(b,a) ; %求出系统的频率响应H MagH=abs(H); %系统的幅频特性Subplot(223);plot(w/pi,MagH); %画出系统的幅频特性 PhaseH=angle(H); %系统的相频特性Subplot(224);plot(w/pi,PhaseH); %画出系统的相频特性 仿照此程序，完成下面问题的编程。

（1） 已知描述两个系统的系统函数分别为-1-2-3-4-1-2-3-40.001836+0.007344z +0.011016z +0.007374z +0.001836z ()1-3.0544+3.8291z -2.2925z +0.55075z H z =1-1-2-3-4-1-2-3-40.001836+0.007344z +0.011016z +0.007374z +0.001836z ()0.4-3.0544+3.8291z -2.2925z +0.55075z H z =1a. 确定系统的稳定性提示：零极点是分析系统频率响应的有力工具之一，在MALAB 中用zplane( )函数画出零点极点图，对于本例，利用零极点图分析系统是否稳定（分别画出两个系统的零极点图，说明其是否稳定）。

b. 假定输入为长度为100的矩形序列，试求出其对两系统的输出，绘出相应的图形。

进一步说明输入有界，输出是否有界。

提示：对于系统求解，在MATLAB 中可以调用filter( )函数完成。

（２）一个特定的线性和时不变系统，描述它的差分方程如下：)3()1(2)()2(25.0)1(5.0)(-+-+=-+--n x n x n x n y n y n ya. 如果此系统的输入为)()]6.0sin(4)2.0cos(35[)(n u n n n x ππ++=。

在1000≤≤n 间求出)(n y 的响应。

b. 图示系统的频率特性提示：使用freqz()函数求解系统的频率特性，幅度响应A=abs(H),相位响应P=angle(H) 3、求下列两个序列的频谱，包括幅度谱和相位谱，分别画出相应的图形。

比较这两者的区别和联系()[,,,,,,,]()[,,,,,,,]h n h n ==124321123412344321提示：本题按离散傅立叶变换（DFT ）实现，DFT 是数字信号处理中最重要的变换之一。

其MATLAB 实现方法可以根据其定义式实现，也可以采用其矩阵形式式实现（实际中均以快速算法FFT 方法实现）。

下面给出其两种实现的MATLA 文件。

①、以定义形式实现%在此输入时域序列，可以采用直接输入的形式，也可采用在运行M 文件后在交互窗口输 %入的形式：xn=input(‘请输入序列x(n):’); N=length(xn);%获得输入序列的长度X=zeros(1,N);%定义信号谱序列X及长度for k=0:N-1for n=0:N-1X(k+1)=X(k+1)+xn(n+1)*exp(-j*2*pi*n*k/N);%按定义式计算序列的DFTend;end;②、以矩阵形式实现x n=input(‘请输入序列x(n):’);N=length(xn);n=0:N-1;k=n;nk=n’*k;WN= exp(-j*2*pi/N).^nk;Xk=xn*WN;[思考题]1、对于一个离散时间系统，如何利用零极点图分析系统是否稳定？2、试说明DTFT的理论分析与实验运算上的区别。

3、实验2中，信号的DTFT一个周期为7个点，相当于单位圆上多远距离取一个值？如果要修改一个周期内的抽样点数，怎么实现？实验四利用FFT和CZT对信号进行频谱分析[实验目的]1、加深对离散傅立叶变换（DFT）算法原理和性质的理解2、掌握离散傅立叶变换DFT的MATLAB实现，理解DFT的频谱泄漏和分辨率问题。

3、熟悉FFT算法原理和FFT子程序的应用4、掌握应用FFT对典型信号进行频谱分析的方法。

5、掌握应用FFT实现两个序列的线性卷积的方法。

[实验仪器]计算机、Matlab6.5（或更高版本）软件[实验参考书] ：自编实验指导书[原理简介]本实验是对信号的频谱或功率谱进行分析。

为了了解信号的特点和频谱分布，可以通过对信号进行谱分析，计算出信号的幅度谱、相位谱和功率谱来实现。

信号的谱分析可以用DFT（FFT）来实现。

在运用DFT进行频谱分析的过程中可能产生三种误差：(1)混叠序列的频谱时被采样信号的周期延拓，当采样速率不满足Nyquist定理时，就会发生频谱混叠，使得采样后的信号序列频谱不能真实的反映原信号的频谱。

避免混叠现象的唯一方法是保证采样速率足够高，使频谱混叠现象不致出现，即在确定采样频率之前，必须对频谱的性质有所了解，在一般情况下，为了保证高于折叠频率的分量不会出现，在采样前，先用低通模拟滤波器对信号进行滤波。

谱分析中的参数选择；A 若已知信号的最高频率c f ，为防止混叠，选定采样频率s f ：c s f f 2≥ (1) B 根据实际需要，选定频率分辨f ∆，一但选定后，即可确定FFT 所需的点数Nf f N s ∆=/ (2)我们希望f ∆越小越好，但f ∆越小，N 越大，计算量、存储量也随之增大。

一般取N 为2的整次幂，以便用FFT 计算，若已给定N ，可用补零方法便N 为2的整次幂。

C s f 和N 确定后，即可确定所需相应模拟信号)(t x 的长度s s NT f N T ==/ (3) 分辨率f ∆反比于T ，而不是N ，在给定的T 的情况下，靠减小s T 来增加N 是不能提高分辨率的，因为s NT T =为常数。

(2) 泄漏实际中我们往往用截短的序列来近似很长的甚至是无限长的序列，这样可以使用较短的DFT 来对信号进行频谱分析，这种截短等价于给原信号序列乘以一个矩形窗函数，也相当于在频域将信号的频谱和矩形窗函数的频谱卷积，所得的频谱是原序列频谱的扩展。

泄漏不能与混叠完全分开，因为泄漏导致频谱的扩展，从而造成混叠。

为了减少泄漏的影响，可以选择适当的窗函数使频谱的扩散减至最小。

(3) 栅栏效应DFT 是对单位圆上Z 变换的均匀采样，所以它不可能将频谱视为一个连续函数，就一定意义上看，用DFT 来观察频谱就好像通过一个栅栏来观看一个图景一样，只能在离散点上看到真实的频谱，这样就有可能发生一些频谱的峰点或谷点被“尖桩的栅栏”所拦住，不能别我们观察到。

减小栅栏效应的一个方法就是借助于在原序列的末端填补一些零值，从而变动DFT 的点数，这一方法实际上是人为地改变了对真实频谱采样的点数和位置，相当于搬动了每一根“尖桩栅栏”的位置，从而使得频谱的峰点或谷点暴露出来。

用FFT 计算线性卷积用FFT 可以实现两个序列的圆周卷积。

在一定的条件下，可以使圆周卷积等于线性卷积。

一般情况，设两个序列的长度分别为N1和N2，要使圆周卷积等于线性卷积的充要条件是FFT 的长度：N ≥N1＋N2－１对于长度不足N 的两个序列，分别将他们补零延长到N 。

当两个序列中有一个序列比较长的时候，我们可以采用分段卷积的方法。

有两种方法： ✧ 重叠相加法。

将长序列分成与短序列相仿的片段，分别用FFT 对它们作线性卷积，再将分段卷积各段重叠的部分相加构成总的卷积输出。

✧ 重叠保留法。

这种方法在长序列分段时，段与段之间保留有互相重叠的部分，在构成总的卷积输出时只需将各段线性卷积部分直接连接起来，省掉了输出段的直接相加。

[实验内容]实验中用到的信号序列a) Gaussian 序列2(),015()0n p q d e n x n --⎧⎪≤≤=⎨⎪⎩b) 衰减正弦序列s i n (2),0()0an b e fnn x n π⎧≤≤⎪=⎨⎪⎩c) 三角波序列1,03()8470c n n x n n n ⎧+≤≤⎪=-≤≤⎨⎪⎩d) 反三角波序列403()3470d n n x n n n -≤≤⎧⎪=-≤≤⎨⎪⎩e ）余弦序列和()cos(0.48)cos(0.52)e x n n n ππ=+实验例程：１、设x(n)是由两个正弦信号及白噪声叠加而成，试用 快速傅立叶变换（FFT ）对其进行频谱分析。