《数字信号处理》实验指导书(完整)

精品文档数字信号处理实验指导书通信教研室安阳工学院二零零九年三月40欢迎下载40欢迎下载。第1章 系统响应及系统稳定性1.1 实验目的l 学会运用MATLAB求解离散时间系统的零状态响应；l 学会运用MATLAB求解离散时间系统的单位取样响应；l 学会运用MATLAB求解离散时间系统的卷积和。1.2 实验原理及实例分析1.2.1 离散时间系统的响应离散时间LTI系统可用线性常系数差分方程来描述，即 （1-1）其中，（，1，N）和（，1，M）为实常数。MATLAB中函数filter可对式（13-1）的差分方程在指定时间范围内的输入序列所产生的响应进行求解。函数filter的语句格式为y=filter(b,a,x)其中，x为输入的离散序列；y为输出的离散序列；y的长度与x的长度一样；b与a分别为差分方程右端与左端的系数向量。【实例1-1】 已知某LTI系统的差分方程为试用MATLAB命令绘出当激励信号为时，该系统的零状态响应。解：MATLAB源程序为a=3 -4 2;b=1 2;n=0:30;x=(1/2).n;y=filter(b,a,x);stem(n,y,fill),grid onxlabel(n),title(系统响应y(n)程序运行结果如图1-1所示。图1-1 实例1-1系统的零状态响应1.2.2 离散时间系统的单位取样响应系统的单位取样响应定义为系统在激励下系统的零状态响应，用表示。MATLAB求解单位取样响应可利用函数filter，并将激励设为前面所定义的impDT函数。例如，求解实例1-1中系统的单位取样响应时，MATLAB源程序为a=3 -4 2;b=1 2;n=0:30;x=impDT(n);h=filter(b,a,x);stem(n,h,fill),grid onxlabel(n),title(系统单位取样响应h(n)程序运行结果如图1-2所示。图1-2 实例1-1的系统单位取样响应MATLAB另一种求单位取样响应的方法是利用控制系统工具箱提供的函数impz来实现。impz函数的常用语句格式为impz(b,a,N)其中，参数N通常为正整数，代表计算单位取样响应的样值个数。【实例1-2】 已知某LTI系统的差分方程为利用MATLAB的impz函数绘出该系统的单位取样响应。解：MATLAB源程序为a=3 -4 2;b=1 2;n=0:30;impz(b,a,30),grid ontitle(系统单位取样响应h(n)程序运行结果如图1-3所示，比较图1-2和图1-3，不难发现结果相同。图1-3 系统单位取样响应1.2.3 离散时间信号的卷积和运算由于系统的零状态响应是激励与系统的单位取样响应的卷积，因此卷积运算在离散时间信号处理领域被广泛应用。离散时间信号的卷积定义为 （1-2）可见，离散时间信号的卷积运算是求和运算，因而常称为“卷积和”。MATLAB求离散时间信号卷积和的命令为conv，其语句格式为y=conv(x,h)其中，x与h表示离散时间信号值的向量；y为卷积结果。用MATLAB进行卷积和运算时，无法实现无限的累加，只能计算时限信号的卷积。例如，利用MALAB的conv命令求两个长为4的矩形序列的卷积和，即，其结果应是长为7（4+4-1=7）的三角序列。用向量1 1 1 1表示矩形序列，MATLAB源程序为x1=1 1 1 1;x2=1 1 1 1;g=conv(x1,x2)g=1 2 3 4 3 2 1如果要绘出图形来，则利用stem命令，即n=1:7;stem(n,g,fill),grid on,xlabel(n)程序运行结果如图1-4所示。图1-4 卷积结果图对于给定函数的卷积和，我们应计算卷积结果的起始点及其长度。两个时限序列的卷积和长度一般等于两个序列长度的和减1。【实例1-3】 已知某系统的单位取样响应为，试用MATLAB求当激励信号为时，系统的零状态响应。解：MATLAB中可通过卷积求解零状态响应，即。由题意可知，描述向量的长度至少为8，描述向量的长度至少为4，因此为了图形完整美观，我们将向量和向量加上一些附加的零值。MATLAB源程序为nx=-1:5; %x(n)向量显示范围(添加了附加的零值)nh=-2:10; %h(n)向量显示范围(添加了附加的零值)x=uDT(nx)-uDT(nx-4);h=0.8.nh.*(uDT(nh)-uDT(nh-8);y=conv(x,h);ny1=nx(1)+nh(1); %卷积结果起始点%卷积结果长度为两序列长度之和减1,即0到(length(nx)+length(nh)-2)%因此卷积结果的时间范围是将上述长度加上起始点的偏移值ny=ny1+(0:(length(nx)+length(nh)-2);subplot(311)stem(nx,x,fill),grid onxlabel(n),title(x(n)axis(-4 16 0 3)subplot(312)stem(nh,h,fill),grid onxlabel(n),title(h(n)axis(-4 16 0 3)subplot(313)stem(ny,y,fill),grid onxlabel(n),title(y(n)=x(n)*h(n)图1-5 利用卷积和法求解系统的零状态响应axis(-4 16 0 3)程序运行结果如图1-5所示。1.3 编程练习1. 试用MATLAB命令求解以下离散时间系统的单位取样响应，并判断系统的稳定性。（1）（2）2. 已知某系统的单位取样响应为，试用MATLAB求当激励信号为时，系统的零状态响应。附：1. 单位取样序列单位取样序列，也称为单位冲激序列，定义为 （12-1）要注意，单位冲激序列不是单位冲激函数的简单离散抽样，它在n=0处是取确定的值1。在MATLAB中，冲激序列可以通过编写以下的impDT.m文件来实现，即function y=impDT(n)y=(n=0); %当参数为0时冲激为1,否则为0调用该函数时n必须为整数或整数向量。【实例2-1】 利用MATLAB的impDT函数绘出单位冲激序列的波形图。解：MATLAB源程序为n=-3:3;x=impDT(n);stem(n,x,fill),xlabel(n),grid ontitle(单位冲激序列)axis(-3 3 -0.1 1.1)图2-1 单位冲激序列程序运行结果如图12-1所示。2. 单位阶跃序列单位阶跃序列定义为 （12-2）在MATLAB中，冲激序列可以通过编写uDT.m文件来实现，即function y=uDT(n)y=n=0; %当参数为非负时输出1调用该函数时n也同样必须为整数或整数向量。【实例2-2】 利用MATLAB的uDT函数绘出单位阶跃序列的波形图。解：MATLAB源程序为n=-3:5;x=uDT(n);stem(n,x,fill),xlabel(n),grid ontitle(单位阶跃序列)axis(-3 5 -0.1 1.1)图2-2 单位阶跃序列程序运行结果如图12-2所示。第2章 z变换及离散时间LTI系统的z域分析2.1 实验目的l 学会运用MATLAB求离散时间信号的z变换和z反变换；l 学会运用MATLAB分析离散时间系统的系统函数的零极点；l 学会运用MATLAB分析系统函数的零极点分布与其时域特性的关系；l 学会运用MATLAB进行离散时间系统的频率特性分析。2.2 实验原理及实例分析2.2.1 z正反变换序列的z变换定义为 （2-1）其中，符号表示取z变换，z是复变量。相应地，单边z变换定义为 （2-2）MATLAB符号数学工具箱提供了计算离散时间信号单边z变换的函数ztrans和z反变换函数iztrans，其语句格式分别为Z=ztrans(x)x=iztrans(z)上式中的x和Z分别为时域表达式和z域表达式的符号表示，可通过sym函数来定义。【实例2-1】 试用ztrans函数求下列函数的z变换。（1）； （2）。解：（1）z变换MATLAB源程序为x=sym(an*cos(pi*n);Z=ztrans(x);simplify(Z)ans=z/(z+a)（2）z变换MATLAB源程序为x=sym(2(n-1)-(-2)(n-1);Z=ztrans(x);simplify(Z)ans=z2/(z-2)/(z+2)【实例2-2】 试用iztrans函数求下列函数的z反变换。（1） （2）解：（1）z反变换MATLAB源程序为Z=sym(8*z-19)/(z2-5*z+6);x=iztrans(Z);simplify(x)ans=-19/6*charfcn0(n)+5*3(n-1)+3*2(n-1)其中，charfcn0(n)是函数在MATLAB符号工具箱中的表示，反变换后的函数形式为。（2）z反变换MATLAB源程序为Z=sym(z*(2*z2-11*z+12)/(z-1)/(z-2)3);x=iztrans(Z);simplify(x)ans=-3+3*2n-1/4*2n*n-1/4*2n*n2其函数形式为。如果信号的z域表示式是有理函数，进行z反变换的另一个方法是对进行部分分式展开，然后求各简单分式的z反变换。设的有理分式表示为 （2-3）MATLAB信号处理工具箱提供了一个对进行部分分式展开的函数residuez，其语句格式为R,P,K=residuez(B,A)其中，B，A分别表示X(z)的分子与分母多项式的系数向量；R为部分分式的系数向量；P为极点向量；K为多项式的系数。若X(z)为有理真分式，则K为零。【实例2-3】 试用MATLAB命令对函数进行部分分式展开，并求出其z反变换。解：MATLAB源程序为B=18;A=18,3,-4,-1;R,P,K=residuez(B,A)R= 0.3600 0.2400 0.4000P= 0.5000 -0.3333 -0.3333K= 从运行结果可知，表示系统有一个二重极点。所以，X(z)的部分分式展开为因此，其z反变换为2.2.2 系统函数的零极点分析离散时间系统的系统函数定义为系统零状态响应的z变换与激励的z变换之比，即 （2-4）如果系统函数的有理函数表示式为 （2-5）那么，在MATLAB中系统函数的零极点就可通过函数roots得到，也可借助函数tf2zp得到，tf2zp的语句格式为Z,P,K=tf2zp(B,A)其中，B与A分别表示的分子与分母多项式的系数向量。它的作用是将的有理分式表示式转换为零极点增益形式，即 （2-6）【实例2-4】 已知一离散因果LTI系统的系统函数为试用MATLAB命令求该系统的零极点。解：用tf2zp函数求系统的零极点，MATLAB源程序为B=1,0.32;A=1,1,0.16;R,P,K=tf2zp(B,A)R= -0.3200P= -0.8000 -0.2000K= 1因此，零点为，极点为与。若要获得系统函数的零极点分布图，可直接应用zplane函数，其语句格式为zplane(B,A)其中，B与A分别表示的分子和分母多项式的系数向量。它的作用是在Z平面上画出单位圆、零点与极点。【实例2-5】 已知一离散因果LTI系统的系统函数为试用MATLAB命令绘出该系统的零极点分布图。解：用zplane函数求系统的零极点，MATLAB源程序为B=1,0,-0.36;A=1,-1.52,0.68;zplane(B,A),grid onlegend(零点,极点)title(零极点分布图)程序运行结果如图14-1所示。可见，该因果系统的极点全部在单位圆内，故系统是稳定的。图2-1 零极点分布图2.2.3 系统函数的零极点分布与其时域特性的关系与拉氏变换在连续系统中的作用类似，在离散系统中，z变换建立了时域函数与z域函数之间的对应关系。因此，z变换的函数从形式可以反映的部分内在性质。我们仍旧通过讨论的一阶极点情况，来说明系统函数的零极点分布与系统时域特性的关系。【实例2-6】 试用MATLAB命令画出现下列系统函数的零极点分布图、以及对应的时域单位取样响应的波形，并分析系统函数的极点对时域波形的影响。（1） （2） （3）（4） （5） （6）（7）解：MATLAB源程序为b1=1,0;a1=1,-0.8;subplot(121)zplane(b1,a1)title(极点在单位圆内的正实数)subplot(122)impz(b1,a1,30);grid on;figureb2=1,0;a2=1,0.8;subplot(121)zplane(b2,a2)title(极点在单位圆内的负实数)subplot(122)impz(b2,a2,30);grid on;figureb3=1,0;a3=1,-1.2,0.72;subplot(121)zplane(b3,a3)title(极点在单位圆内的共轭复数)subplot(122)impz(b3,a3,30);grid on;figureb4=1,0;a4=1,-1;subplot(121)zplane(b4,a4)title(极点在单位圆上为实数1)subplot(122)impz(b4,a4);grid on;figureb5=1,0;a5=1,-1.6,1;subplot(121)zplane(b5,a5)title(极点在单位圆上的共轭复数)subplot(122)impz(b5,a5,30);grid on;figureb6=1,0;a6=1,-1.2;subplot(121)zplane(b6,a6)title(极点在单位圆外的正实数)subplot(122)impz(b6,a6,30);grid on;figureb7=1,0;a7=1,-2,1.36;subplot(121)zplane(b7,a7)title(极点在单位圆外的共轭复数)subplot(122)impz(b7,a7,30);grid on;程序运行结果分别如图2-2的（a）、（b）、（c）、（d）、（e）、（f）、（g）所示。(a)(b)(c)(d)(e)(f)图2-2 系统函数的零极点分布与其时域特性的关系(g)从图14-2可知，当极点位于单位圆内时，为衰减序列；当极点位于单位圆上时，为等幅序列；当极点位于单位圆外时，为增幅序列。若有一阶实数极点，则为指数序列；若有一阶共轭极点，则为指数振荡序列；若的极点位于虚轴左边，则序列按一正一负的规律交替变化。2.2.4 离散时间LTI系统的频率特性分析对于因果稳定的离散时间系统，如果激励序列为正弦序列，则系统的稳态响应为。其中，通常是复数。离散时间系统的频率响应定义为 （2-7）其中，称为离散时间系统的幅频特性；称为离散时间系统的相频特性；是以（，若零，）为周期的周期函数。因此，只要分析在范围内的情况，便可分析出系统的整个频率特性。MATLAB提供了求离散时间系统频响特性的函数freqz，调用freqz的格式主要有两种。一种形式为H,w=freqz(B,A,N)其中，B与A分别表示的分子和分母多项式的系数向量；N为正整数，默认值为512；返回值w包含范围内的N个频率等分点；返回值H则是离散时间系统频率响应在范围内N个频率处的值。另一种形式为H,w=freqz(B,A,N,whole)与第一种方式不同之处在于角频率的范围由扩展到。【实例2-6】 用MATLAB命令绘制系统的频率响应曲线。 解：利用函数freqz计算出，然后利用函数abs和angle分别求出幅频特性与相频特性，最后利用plot命令绘出曲线。MATLAB源程序为b=1 -0.96 0.9028;a=1 -1.56 0.8109;H,w=freqz(b,a,400,whole);Hm=abs(H); Hp=angle(H);subplot(211)plot(w,Hm),grid onxlabel(omega(rad/s),ylabel(Magnitude)title(离散系统幅频特性曲线)subplot(212)plot(w,Hp),grid onxlabel(omega(rad/s),ylabel(Phase)title(离散系统相频特性曲线)程序运行结果如图2-3所示。图4-3 离散系统频响特性曲线2.3 编程练习1. 试用MATLAB的residuez函数，求出的部分分式展开和。结果：MATLAB源程序为B=2 16 44 56 32;A=3,3，-15,18，-12;R,P,K=residuez(B,A)R = -0.0177 9.4914 -3.0702 + 2.3398i -3.0702 - 2.3398iP = -3.2361 1.2361 0.5000 + 0.8660i 0.5000 - 0.8660iK = -2.66672. 试用MATLAB画出下列因果系统的系统函数零极点分布图，并判断系统的稳定性。（1） 结果：程序：B=2,-1.6,-0.9;A=1,-2.5,1.96,-0.48;R,P,K=tf2zp(B,A)zplane(B,A) （2）结果：程序：B=1,-1;A=1,-0.9,-0.65,0.873;R,P,K=tf2zp(B,A)zplane(B,A)3. 试用MATLAB绘制系统的频率响应曲线。结果：程序： b=1; a=1,-0.75,0.125; H,w=freqz(b,a,400,whole); Hm=abs(H); Hp=angle(H); subplot(211) plot(w,Hm),grid on xlabel(omega(rad/s),ylabel(Magnitude) title(离散系统幅频特性曲线) subplot(212) plot(w,Hp),grid on xlabel(omega(rad/s),ylabel(Phase) title(离散系统相频特性曲线)实验三 时域采样与频域采样1. 实验目的 时域采样理论与频域采样理论是数字信号处理中的重要理论。要求掌握模拟信号采样前后频谱的变化，以及如何选择采样频率才能使采样后的信号不丢失信息；要求掌握频率域采样会引起时域周期化的概念，以及频率域采样定理及其对频域采样点数选择的指导作用。2. 实验原理与方法 时域采样定理的要点是：a) 对模拟信号以间隔T进行时域等间隔理想采样，形成的采样信号的频谱是原模拟信号频谱以采样角频率（）为周期进行周期延拓。公式为： b) 采样频率必须大于等于模拟信号最高频率的两倍以上，才能使采样信号的 频谱不产生频谱混叠。 利用计算机计算上式并不方便，下面我们导出另外一个公式，以便用计算机上进行实验。 理想采样信号和模拟信号之间的关系为： 对上式进行傅立叶变换，得到： 在上式的积分号内只有当时，才有非零值，因此： 上式中，在数值上，再将代入，得到： 上式的右边就是序列的傅立叶变换，即 上式说明理想采样信号的傅立叶变换可用相应的采样序列的傅立叶变换得到，只要将自变量用代替即可。 频域采样定理的要点是：a) 对信号x(n)的频谱函数X(ej)在0，2上等间隔采样N点，得到则N点IDFT得到的序列就是原序列x(n)以N为周期进行周期延拓后的主值区序列，公式为： b) 由上式可知，频域采样点数N必须大于等于时域离散信号的长度M(即NM)，才能使时域不产生混叠，则N点IDFT得到的序列就是原序列x(n),即=x(n)。如果NM，比原序列尾部多N-M个零点；如果NM，z则=IDFT发生了时域混叠失真，而且的长度N也比x(n)的长度M短，因此。与x(n)不相同。 在数字信号处理的应用中，只要涉及时域或者频域采样，都必须服从这两个采样理论的要点。 对比上面叙述的时域采样原理和频域采样原理，得到一个有用的结论，这两个采样理论具有对偶性：“时域采样频谱周期延拓，频域采样时域信号周期延拓”。因此放在一起进行实验。3. 实验内容及步骤（1）时域采样理论的验证。给定模拟信号， 式中A=444.128，=50，=50rad/s，它的幅频特性曲线如图10.2.1 图3.1 的幅频特性曲线 现用DFT(FFT)求该模拟信号的幅频特性，以验证时域采样理论。 安照的幅频特性曲线，选取三种采样频率，即=1kHz，300Hz，200Hz。观测时间选。 为使用DFT，首先用下面公式产生时域离散信号，对三种采样频率，采样序列按顺序用，表示。 因为采样频率不同，得到的，的长度不同， 长度（点数）用公式计算。选FFT的变换点数为M=64，序列长度不够64的尾部加零。X(k)=FFTx(n) ， k=0,1,2,3,-,M-1式中k代表的频率为 。要求： 编写实验程序，计算、和的幅度特性，并绘图显示。观察分析频谱混叠失真。（2）频域采样理论的验证。给定信号如下： 编写程序分别对频谱函数在区间上等间隔采样32和16点，得到： 再分别对进行32点和16点IFFT，得到： 分别画出、的幅度谱，并绘图显示x(n)、的波形，进行对比和分析，验证总结频域采样理论。提示：频域采样用以下方法容易变程序实现。 直接调用MATLAB函数fft计算就得到在的32点频率域采样 抽取的偶数点即可得到在的16点频率域采样，即。 当然也可以按照频域采样理论，先将信号x(n)以16为周期进行周期延拓，取其主值区（16点），再对其进行16点DFT(FFT),得到的就是在的16点频率域采样。 4思考题： 如果序列x(n)的长度为M，希望得到其频谱在上的N点等间隔采样，当NM时， 如何用一次最少点数的DFT得到该频谱采样？5. 实验报告及要求a) 运行程序打印要求显示的图形。 b) 分析比较实验结果，简述由实验得到的主要结论c) 简要回答思考题d) 附上程序清单和有关曲线。4、实验程序清单1） 时域采样理论的验证程序清单% 时域采样理论验证程序exp2a.mTp=64/1000;%观察时间Tp=64微秒%产生M长采样序列x(n)% Fs=1000;T=1/Fs;Fs=1000;T=1/Fs;M=Tp*Fs;n=0:M-1;A=444.128;alph=pi*50*20.5;omega=pi*50*20.5;xnt=A*exp(-alph*n*T).*sin(omega*n*T);Xk=T*fft(xnt,M); %M点FFTxnt)subplot(3,2,1);stem(n,xnt,.);grid on;%调用绘图函数stem绘制序列图title(a) Fs=1000Hz);k=0:M-1;fk=k/Tp;subplot(3,2,2);plot(fk,abs(Xk);title(a) T*FTxa(nT),Fs=1000Hz);xlabel(f(Hz);ylabel(幅度);axis(0,Fs,0,1.2*max(abs(Xk)%=% Fs=300Hz和 Fs=200Hz的程序与上面Fs=1000Hz完全相同。(1) Fs=1000Hz的图片：如a所示(2) Fs=200Hz的图片：如b所示(3) Fs=300Hz的图片：如c所示(4) 综合图片如e所示2 频域采样理论的验证程序清单%频域采样理论验证程序exp2b.mM=27;N=32;n=0:M;%产生M长三角波序列x(n)xa=0:floor(M/2); xb= ceil(M/2)-1:-1:0; xn=xa,xb;Xk=fft(xn,1024);%1024点FFTx(n), 用于近似序列x(n)的TFX32k=fft(xn,32);%32点FFTx(n)x32n=ifft(X32k);%32点IFFTX32(k)得到x32(n)X16k=X32k(1:2:N);%隔点抽取X32k得到X16(K)x16n=ifft(X16k,N/2);%16点IFFTX16(k)得到x16(n)subplot(3,2,2);stem(n,xn,.);box ontitle(b) 三角波序列x(n);xlabel(n);ylabel(x(n);axis(0,32,0,20)k=0:1023;wk=k/1024;%subplot(3,2,1);plot(wk,abs(Xk);title(a)FTx(n);xlabel(omega/pi);ylabel(|X(ejomega)|);axis(0,1,0,200)k=0:N/2-1;subplot(3,2,3);stem(k,abs(X16k),.);box ontitle(c) 16点频域采样);xlabel(k);ylabel(|X_1_6(k)|);axis(0,8,0,200)n1=0:N/2-1;subplot(3,2,4);stem(n1,x16n,.);box ontitle(d) 16点IDFTX_1_6(k);xlabel(n);ylabel(x_1_6(n);axis(0,32,0,20)k=0:N-1;subplot(3,2,5);stem(k,abs(X32k),.);box ontitle(e) 32点频域采样);xlabel(k);ylabel(|X_3_2(k)|);axis(0,16,0,200)n1=0:N-1;subplot(3,2,6);stem(n1,x32n,.);box ontitle(f) 32点IDFTX_3_2(k);xlabel(n);ylabel(x_3_2(n);axis(0,32,0,20)结果如图d所示实验四 用FFT对信号作频谱分析1实验目的 学习用FFT对连续信号和时域离散信号进行谱分析的方法，了解可能出现的分析 误差及其原因，以便正确应用FFT。2. 实验原理用FFT对信号作频谱分析是学习数字信号处理的重要内容。经常需要进行谱分析的信号是模拟信号和时域离散信号。对信号进行谱分析的重要问题是频谱分辨率D和分析误差。频谱分辨率直接和FFT的变换区间N有关，因为FFT能够实现的频率分辨率是，因此要求。可以根据此式选择FFT的变换区间N。误差主要来自于用FFT作频谱分析时，得到的是离散谱，而信号（周期信号除外）是连续谱，只有当N较大时离散谱的包络才能逼近于连续谱，因此N要适当选择大一些。周期信号的频谱是离散谱，只有用整数倍周期的长度作FFT，得到的离散谱才能代表周期信号的频谱。如果不知道信号周期，可以尽量选择信号的观察时间长一些。对模拟信号进行谱分析时，首先要按照采样定理将其变成时域离散信号。如果是模拟周期信号，也应该选取整数倍周期的长度，经过采样后形成周期序列，按照周期序列的谱分析进行。3实验步骤及内容（1）对以下序列进行谱分析。 选择FFT的变换区间N为8和16 两种情况进行频谱分析。分别打印其幅频特性曲线。 并进行对比、分析和讨论。（2）对以下周期序列进行谱分析。 选择FFT的变换区间N为8和16 两种情况分别对以上序列进行频谱分析。分别打印其幅频特性曲线。并进行对比、分析和讨论。（3）对模拟周期信号进行谱分析 选择 采样频率，变换区间N=16,32,64 三种情况进行谱分析。分别打印其幅频特性，并进行分析和讨论。 4思考题（1）对于周期序列，如果周期不知道，如何用FFT进行谱分析？（2）如何选择FFT的变换区间？（包括非周期信号和周期信号）（3）当N=8时，和的幅频特性会相同吗？为什么？N=16 呢？5实验报告要求（1）完成各个实验任务和要求。附上程序清单和有关曲线。（2）简要回答思考题。附录： 实验程序清单function y=mstem(x,N)k=0:N-1;wk=2*k/N;stem(wk,abs(x),.);% 用FFT对信号作频谱分析clear all;close all%实验内容(1)=x1n=ones(1,4); %产生序列向量x1(n)=R4(n)M=8;xa=1:(M/2); xb=(M/2):-1:1; x2n=xa,xb; %产生长度为8的三角波序列x2(n)x3n=xb,xa;X1k8=fft(x1n,8); %计算x1n的8点DFTX1k16=fft(x1n,16); %计算x1n的16点DFTX2k8=fft(x2n,8); %计算x1n的8点DFTX2k16=fft(x2n,16); %计算x1n的16点DFTX3k8=fft(x3n,8); %计算x1n的8点DFTX3k16=fft(x3n,16); %计算x1n的16点DFT%以下绘制幅频特性曲线subplot(2,2,1);mstem(X1k8); %绘制8点DFT的幅频特性图title(1a) 8点DFTx_1(n);xlabel(/);ylabel(幅度);axis(0,2,0,1.2*max(abs(X1k8)subplot(2,2,3);mstem(X1k16); %绘制16点DFT的幅频特性图title(1b)16点DFTx_1(n);xlabel(/);ylabel(幅度);axis(0,2,0,1.2*max(abs(X1k16)figure(2)subplot(2,2,1);mstem(X2k8); %绘制8点DFT的幅频特性图title(2a) 8点DFTx_2(n);xlabel(/);ylabel(幅度);axis(0,2,0,1.2*max(abs(X2k8)subplot(2,2,2);mstem(X2k16); %绘制16点DFT的幅频特性图title(2b)16点DFTx_2(n);xlabel(/);ylabel(幅度);axis(0,2,0,1.2*max(abs(X2k16)subplot(2,2,3);mstem(X3k8); %绘制8点DFT的幅频特性图title(3a) 8点DFTx_3(n);xlabel(/);ylabel(幅度);axis(0,2,0,1.2*max(abs(X3k8)subplot(2,2,4);mstem(X3k16); %绘制16点DFT的幅频特性图title(3b)16点DFTx_3(n);xlabel(/);ylabel(幅度);axis(0,2,0,1.2*max(abs(X3k16)%实验内容(2) 周期序列谱分析=N=8;n=0:N-1; %FFT的变换区间N=8x4n=cos(pi*n/4);x5n=cos(pi*n/4)+cos(pi*n/8);X4k8=fft(x4n); %计算x4n的8点DFTX5k8=fft(x5n); %计算x5n的8点DFTN=16;n=0:N-1; %FFT的变换区间N=16x4n=cos(pi*n/4);x5n=cos(pi*n/4)+cos(pi*n/8);X4k16=fft(x4n); %计算x4n的16点DFTX5k16=fft(x5n); %计算x5n的16点DFTfigure(3)subplot(2,2,1);mstem(X4k8); %绘制8点DFT的幅频特性图title(4a) 8点DFTx_4(n);xlabel(/);ylabel(幅度);axis(0,2,0,1.2*max(abs(X4k8)subplot(2,2,3);mstem(X4k16); %绘制16点DFT的幅频特性图title(4b)16点DFTx_4(n);xlabel(/);ylabel(幅度);axis(0,2,0,1.2*max(abs(X4k16)