数字信号处理第8章.ppt
第8章上 机 实 验,8.1实验一:系统响应及系统稳定性8.2 实验二:时域采样与频域采样8.3实验三:用FFT对信号作频谱分析8.4实验四:IIR数字滤波器设计及软件实现8.5实验五:FIR数字滤波器设计与软件实现8.6 实验六:数字信号处理在双音多频拨号系统中的应用,8.1 实验一:系统响应及系统稳定性8.1.1 实验指导1.实验目的(1)掌握求系统响应的方法。(2)掌握时域离散系统的时域特性。(3)分析、观察及检验系统的稳定性。,2.实验原理与方法在时域中,描写系统特性的方法是差分方程和单位脉冲响应,在频域可以用系统函数描述系统特性。已知输入信号可以由差分方程、单位脉冲响应或系统函数求出系统对于该输入信号的响应。本实验仅在时域求解。在计算机上适合用递推法求差分方程的解,最简单的方法是采用MATLAB语言的工具箱函数filter函数。也可以用MATLAB语言的工具箱函数conv函数计算输入信号和系统的单位脉冲响应的线性卷积,求出系统的响应。,系统的时域特性指的是系统的线性时不变性质、因果性和稳定性。重点分析实验系统的稳定性,包括观察系统的暂态响应和稳定响应。系统的稳定性是指对任意有界的输入信号,系统都能得到有界的系统响应。或者系统的单位脉冲响应满足绝对可和的条件。系统的稳定性由其差分方程的系数决定。,实际中检查系统是否稳定,不可能检查系统对所有有界的输入信号,输出是否都是有界输出,或者检查系统的单位脉冲响应满足绝对可和的条件。可行的方法是在系统的输入端加入单位阶跃序列,如果系统的输出趋近一个常数(包括零),就可以断定系统是稳定的12。系统的稳态输出是指当n时,系统的输出。如果系统稳定,则信号加入系统后,系统输出的开始一段称为暂态效应,随着n的加大,幅度趋于稳定,达到稳态输出。注意在以下实验中均假设系统的初始状态为零,3 实验内容及步骤(1)编制程序,包括产生输入信号、单位脉冲响应序列的子程序,用filter函数或conv函数求解系统输出响应的主程序。程序中要有绘制信号波形的功能。(2)给定一个低通滤波器的差分方程为y(n)=0.05x(n)+0.05x(n1)+0.9y(n1)输入信号x1(n)=R8(n)x2(n)=u(n),分别求出x1(n)=R8(n)和x2(n)=u(n)的系统响应y1(n)和y2(n),并画出其波形。求出系统的单位脉冲响应,画出其波形。(3)给定系统的单位脉冲响应为h1(n)=R10(n)h2(n)=(n)+2.5(n1)+2.5(n2)+(n3)用线性卷积法求x1(n)=R8(n)分别对系统h1(n)和h2(n)的输出响应y21(n)和y22(n),并画出波形。,(4)给定一谐振器的差分方程为y(n)=1.8237y(n1)0.9801y(n2)+b0 x(n)b0 x(n2)令b0=1/100.49,谐振器的谐振频率为0.4 rad。用实验方法检查系统是否稳定。输入信号为u(n)时,画出系统输出波形y31(n)。给定输入信号为x(n)=sin(0.014n)+sin(0.4n)求出系统的输出响应y32(n),并画出其波形。,4 思考题(1)如果输入信号为无限长序列,系统的单位脉冲响应是有限长序列,可否用线性卷积法求系统的响应?如何求(2)如果信号经过低通滤波器,信号的高频分量被滤掉,时域信号会有何变化?用前面第一个实验的结果进行分析说明。5 实验报告要求(1)简述在时域求系统响应的方法。(2)简述通过实验判断系统稳定性的方法。分析上面第三个实验的稳定输出的波形。(3)对各实验所得结果进行简单分析和解释。(4)简要回答思考题。(5)打印程序清单和要求的各信号波形。,8.1.2 实验参考程序实验1程序:exp1.m%实验1:系统响应及系统稳定性close all;clear all%=%内容1:调用filter解差分方程,由系统对u(n)的响应判断稳定性A=1,-0.9;B=0.05,0.05;%系统差分方程系数向量B和Ax1n=1 1 1 1 1 1 1 1 zeros(1,50);%产生信号x1n=R8n,x2n=ones(1,128);%产生信号x2n=unhn=impz(B,A,58);%求系统单位脉冲响应h(n)subplot(2,2,1);y=h(n);tstem(hn,y);%调用函数tstem绘图title(a)系统单位脉冲响应h(n)y1n=filter(B,A,x1n);%求系统对x1n的响应y1nsubplot(2,2,2);y=y1(n);tstem(y1n,y);title(b)系统对R8(n)的响应y1(n)y2n=filter(B,A,x2n);%求系统对x2n的响应y2nsubplot(2,2,4);y=y2(n);tstem(y2n,y);title(c)系统对u(n)的响应y2(n)%=,%内容2:调用conv函数计算卷积x1n=1 1 1 1 1 1 1 1;%产生信号x1n=R8nh1n=ones(1,10)zeros(1,10);h2n=1 2.5 2.5 1 zeros(1,10);y21n=conv(h1n,x1n);y22n=conv(h2n,x1n);figure(2)subplot(2,2,1);y=h1(n);tstem(h1n,y);%调用函数tstem绘图title(d)系统单位脉冲响应h1(n)subplot(2,2,2);y=y21(n);tstem(y21n,y);,title(e)h1(n)与R8(n)的卷积y21(n)subplot(2,2,3);y=h2(n);tstem(h2n,y);%调用函数tstem绘图title(f)系统单位脉冲响应h2(n)subplot(2,2,4);y=y22(n);tstem(y22n,y);title(g)h2(n)与R8(n)的卷积y22(n)%=%内容3:谐振器分析un=ones(1,256);%产生信号unn=0:255;xsin=sin(0.014*n)+sin(0.4*n);%产生正弦信号,A=1,1.8237,0.9801;B=1/100.49,0,1/100.49;%系统差分方程系数向量B和Ay31n=filter(B,A,un);%谐振器对un的响应y31ny32n=filter(B,A,xsin);%谐振器对正弦信号的响应y32nfigure(3)subplot(2,1,1);y=y31(n);tstem(y31n,y)title(h)谐振器对u(n)的响应y31(n)subplot(2,1,2);y=y32(n);tstem(y32n,y);title(i)谐振器对正弦信号的响应y32(n),8.1.3 实验结果与波形实验结果与波形如图所示。,图,8.1.4 分析与讨论(1)综合起来,在时域求系统响应的方法有两种,第一种是通过解差分方程求得系统输出,注意要合理地选择初始条件;第二种是已知系统的单位脉冲响应,通过求输入信号和系统单位脉冲响应的线性卷积求得系统输出。用计算机求解时最好使用MATLAB语言进行。(2)实际中要检验系统的稳定性,其方法是在输入端加入单位阶跃序列,观察输出波形,如果波形稳定在一个常数值上,系统稳定,否则不稳定。上面第三个实验是稳定的。(3)谐振器具有对某个频率进行谐振的性质,本实验中的谐振器的谐振频率是0.4 rad,因此稳定波形为sin(0.4n)。,(4)如果输入信号为无限长序列,系统的单位脉冲响应是有限长序列,可用分段线性卷积法求系统的响应,具体方法请参考DFT一章的内容。如果信号经过低通滤波器,则信号的高频分量被滤掉,时域信号的变化减缓,在有阶跃处附近产生过渡带。因此,当输入矩形序列时,输出序列的开始和终了都产生了明显的过渡带,见第一个实验结果的波形。,8.2 实验二:时域采样与频域采样8.2.1 实验指导1.实验目的时域采样理论与频域采样理论是数字信号处理中的重要理论。要求掌握模拟信号采样前后频谱的变化,以及如何选择采样频率才能使采样后的信号不丢失信息;要求掌握频域采样会引起时域周期化的概念,以及频率域采样定理及其对频域采样点数选择的指导作用。,2.实验原理与方法1)时域采样定理的要点 时域采样定理的要点是:(1)对模拟信号xa(t)以T进行时域等间隔理想采样,形成的采样信号的频谱会以采样角频率s(s=2/T)为周期进行周期延拓。公式为,(2)采样频率s必须大于等于模拟信号最高频率的两倍以上,才能使采样信号的频谱不产生频谱混叠。利用计算机计算并不方便,下面我们导出另外一个公式,以便在计算机上进行实验。理想采样信号和模拟信号xa(t)之间的关系为,对上式进行傅里叶变换,得到,在上式的积分号内只有当t=nT时,才有非零值,因此,上式中,在数值上xa(nT)x(n),再将=T代入,得到,上式的右边就是序列的傅里叶变换X(ej),即,上式说明理想采样信号的傅里叶变换可用相应的采样序列的傅里叶变换得到,只要将自变量用T代替即可。,2)频域采样定理的要点频域采样定理的要点是:(1)对信号x(n)的频谱函数X(ej)在0,2上等间隔采样N点,得到,则N点IDFTXN(k)得到的序列就是原序列x(n)以N为周期进行周期延拓后的主值区序列,公式为,(2)由上式可知,频域采样点数N必须大于等于时域离散信号的长度M(即NM),才能使时域不产生混叠,这时N点IDFTXN(k)得到的序列xN(n)就是原序列x(n),即xN(n)=x(n)。如果NM,则xN(n)比原序列尾部多NM个零点;如果NM,则xN(n)=IDFTXN(k)发生了时域混叠失真,而且xN(n)的长度N也比x(n)的长度M短,因此,xN(n)与x(n)不相同。,在数字信号处理的应用中,只要涉及时域采样或者频域采样,都必须服从这两个采样理论的要点。对比上面叙述的时域采样原理和频域采样原理,得到一个有用的结论,即两个采样理论具有对偶性:“时域采样频谱周期延拓,频域采样时域信号周期延拓”。因此把这两部分内容放在一起进行实验。,3.实验内容及步骤1)时域采样理论的验证给定模拟信号xa(t)=Aet sin(0t)u(t)式中,A=444.128,=50,0=50 rad/s,它的幅频特性曲线如图所示。现用DFT(FFT)求该模拟信号的幅频特性,以验证时域采样理论。按照xa(t)的幅频特性曲线,选取三种采样频率,即Fs=1 kHz,300 Hz,200 Hz。观测时间选Tp=50 ms。,图8.2.1 xa(t)的幅频特性曲线,为使用DFT,首先用下面公式产生时域离散信号,对三种采样频率,采样序列按顺序用x1(n)、x2(n)、x3(n)表示。x(n)=xa(nT)=AenTsin(0nT)u(nT)因为采样频率不同,得到的x1(n)、x2(n)、x3(n)的长度不同,长度(点数)用公式N=TpFs计算。选FFT的变换点数为M=64,序列长度不够64的尾部加零。X(k)=FFTx(n)k=0,1,2,3,M1式中,k代表的频率为,要求:编写实验程序,计算x1(n)、x2(n)和x3(n)的幅度特性,并绘图显示。观察分析频谱混叠失真。2)频域采样理论的验证给定信号如下:,0n13,14n26,其它n,编写程序分别对频谱函数X(ej)=FTx(n)在区间0,2上等间隔采样32点和16点,得到X32(k)和X16(k):,再分别对X32(k)和X16(k)进行32点和16点IFFT,得到x32(n)和x16(n):,x32(n)=IFFTX32(k)32 n=0,1,2,31x16(n)=IFFTX16(k)16 n=0,1,2,15分别画出X(ej)、X32(k)和X16(k)的幅度谱,并绘图显示x(n)、x32(n)和x16(n)的波形,进行对比和分析,验证和总结频域采样理论。提示:频域采样用以下方法容易编程序实现。(1)直接调用MATLAB函数fft计算X32(k)=FFTx(n)32,就得到X(ej)在0,2的32点频率域采样。,(2)抽取X32(k)的偶数点即可得到X(ej)在0,2的16点频率域采样X16(k),即X16(k)=X32(2k)k=0,1,2,15。(3)也可以按照频域采样理论,先将信号x(n)以16为周期进行周期延拓,取其主值区(16点),再对其进行16点DFT(FFT),得到的就是X(ej)在0,2的16点频率域采样X16(k)。,4 思考题 如果序列x(n)的长度为M,希望得到其频谱X(ej)在0,2上的N点等间隔采样,当NM时,如何用一次最少点数的DFT得到该频谱采样?5.实验报告及要求(1)运行程序,打印要求显示的图形。(2)分析比较实验结果,简述由实验得到的主要结论。(3)简要回答思考题。(4)附上程序清单和有关曲线。,8.2.2 实验程序清单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:M1;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),yn=xa(nT);subplot(3,2,1);tstem(xnt,yn);%调用自编绘图函数tstem绘制序列图box on;title(a)Fs=1000Hz);k=0:M1;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的程序完全相同。,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=2*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/21;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/21;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:N1;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:N1;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),8.2.3 实验程序运行结果(1)时域采样理论的验证程序exp2a.m的运行结果如图所示。由图可见,当采样频率为1000 Hz时,频谱混叠很小;当采样频率为300 Hz时,频谱混叠很严重;当采样频率为200 Hz时,频谱混叠更很严重。,图,(2)频域采样理论的验证程序exp2b.m的运行结果如图所示。该图验证了频域采样理论和频域采样定理。对信号x(n)的频谱函数X(ej)在0,2上等间隔采样N=16时,N点IDFTXN(k)得到的序列正是原序列x(n)以16为周期进行周期延拓后的主值区序列:,由于NM,所以发生了时域混叠失真,因此,xN(n)与x(n)不相同。,图,8.2.4 简答思考题对实验中的思考题简答如下:先对原序列x(n)以N为周期进行周期延拓后取主值区序列,,再计算N点DFT,则得到N点频域采样:,8.3 实验三:用FFT对信号作频谱分析8.3.1 实验指导1 实验目的学习用FFT对连续信号和时域离散信号进行频谱分析(也称谱分析)的方法,了解可能出现的分析误差及其原因,以便正确应用FFT。2.实验原理用FFT对信号作频谱分析是学习数字信号处理的重要内容。经常需要进行谱分析的信号是模拟信号和时域离散信号。对信号进行谱分析的重要问题是频谱分辨率D和分析误差。,频谱分辨率直接和FFT的变换区间N有关,因为FFT能够实现的频率分辨率是2/N,因此要求2/ND。可以根据此式选择FFT的变换区间N。误差主要来自于用FFT作频谱分析时,得到的是离散谱,而信号(周期信号除外)是连续谱,只有当N较大时,离散谱的包络才能逼近于连续谱,因此N要适当选择大一些。周期信号的频谱是离散谱,只有用整数倍周期的长度作FFT,得到的离散谱才能代表周期信号的频谱。如果不知道信号周期,可以尽量选择信号的观察时间长一些。,对模拟信号进行谱分析时,首先要按照采样定理将其变成时域离散信号。如果是模拟周期信号,也应该选取整数倍周期的长度,经过采样后形成周期序列,按照周期序列的谱分析进行。3 实验内容及步骤(1)对以下序列进行谱分析:x1(n)=R4(n)n+1 0n38n 4n70 其它n 4n 0n3 n3 4n7 0 其它n,x2(n)=,x3(n)=,选择FFT的变换区间N为8和16的两种情况进行频谱分析。分别打印其幅频特性曲线,并进行对比、分析和讨论。(2)对以下周期序列进行谱分析:,选择FFT的变换区间N为8和16的两种情况分别对以上序列进行频谱分析。分别打印其幅频特性曲线。并进行对比、分析和讨论。(3)对模拟周期信号进行谱分析:x6(t)=cos8t+cos16t+cos20t选择采样频率Fs=64 Hz,变换区间N=16,32,64的三种情况进行谱分析。分别打印其幅频特性,并进行分析和讨论。,4 思考题(1)对于周期序列,如果周期不知道,如何用FFT进行谱分析?(2)如何选择FFT的变换区间(包括非周期信号和周期信号)?(3)当N=8时,x2(n)和x3(n)的幅频特性会相同吗?为什么?N=16时呢?5 实验报告要求(1)完成各个实验任务和要求,附上程序清单和有关曲线。(2)简要回答思考题。,8.3.2 实验程序清单实验三程序exp3.m%用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点DFT,X2k8=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:N1;%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:N1;%FFT的变换区间N=16x4n=cos(pi*n/4);x5n=cos(pi*n/4)+cos(pi*n/8);X4k16=fft(x4n);%计算x4n的16点DFT,X5k16=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),subplot(2,2,2);mstem(X5k8);%绘制8点DFT的幅频特性图title(5a)8点DFTx_5(n);xlabel(/);ylabel(幅度);axis(0,2,0,1.2*max(abs(X5k8)subplot(2,2,4);mstem(X5k16);%绘制16点DFT的幅频特性图title(5b)16点DFTx_5(n);xlabel(/);ylabel(幅度);axis(0,2,0,1.2*max(abs(X5k16)%实验内容(3)模拟周期信号谱分析=figure(4),Fs=64;T=1/Fs;N=16;n=0:N1;%FFT的变换区间N=16x6nT=cos(8*pi*n*T)+cos(16*pi*n*T)+cos(20*pi*n*T);%对x6(t)16点采样X6k16=fft(x6nT);%计算x6nT的16点DFTX6k16=fftshift(X6k16);%将零频率移到频谱中心 Tp=N*T;F=1/Tp;%频率分辨率Fk=N/2:N/21;fk=k*F;%产生16点DFT对应的采样点频率(以零频率为中心)subplot(3,1,1);stem(fk,abs(X6k16),.);box on%绘制8点DFT的幅频特性图,title(6a)16点|DFTx_6(nT)|);xlabel(f(Hz);ylabel(幅度);axis(N*F/21,N*F/21,0,1.2*max(abs(X6k16)N=32;n=0:N-1;%FFT的变换区间N=16x6nT=cos(8*pi*n*T)+cos(16*pi*n*T)+cos(20*pi*n*T);%对x6(t)32点采样X6k32=fft(x6nT);%计算x6nT的32点DFTX6k32=fftshift(X6k32);%将零频率移到频谱中心 Tp=N*T;F=1/Tp;%频率分辨率Fk=N/2:N/21;fk=k*F;%产生16点DFT对应的采样点频率(以零频率为中心),subplot(3,1,2);stem(fk,abs(X6k32),.);box on%绘制8点DFT的幅频特性图title(6b)32点|DFTx_6(nT)|);xlabel(f(Hz);ylabel(幅度);axis(N*F/21,N*F/21,0,1.2*max(abs(X6k32)N=64;n=0:N1;%FFT的变换区间N=16x6nT=cos(8*pi*n*T)+cos(16*pi*n*T)+cos(20*pi*n*T);%对x6(t)64点采样X6k64=fft(x6nT);%计算x6nT的64点DFT,X6k64=fftshift(X6k64);%将零频率移到频谱中心 Tp=N*T;F=1/Tp;%频率分辨率Fk=N/2:N/21;fk=k*F;%产生16点DFT对应的采样点频率(以零频率为中心)subplot(3,1,3);stem(fk,abs(X6k64),.);box on%绘制8点DFT的幅频特性图title(6a)64点|DFTx_6(nT)|);xlabel(f(Hz);ylabel(幅度);axis(N*F/21,N*F/21,0,1.2*max(abs(X6k64),8.3.3 实验程序运行结果实验三程序exp3.m运行结果如图所示。,图,8.3.4 分析与讨论请读者注意,用DFT(或FFT)分析频谱和绘制频谱图时,最好将X(k)的自变量k换算成对应的频率并作为横坐标,便于观察频谱。,为了便于读取频率值,最好关于归一化,即以/作为横坐标。,1.实验内容(1)图(1a)和(1b)说明x1(n)=R4(n)的8点DFT和16点DFT分别是x1(n)的频谱函数的8点和16点采样;因为x3(n)=x2(n+3)8R8(n),所以,x3(n)与x2(n)的8点DFT的模相等,如图(2a)和(3a)所示。但是,当N=16时,x3(n)与x2(n)不满足循环移位关系,所以图(2b)和(3b)的模不同。,2.实验内容(2)对周期序列谱分析。x4(n)=cosn的周期为8,所以N=8和N=16均是其周期的整数倍,得到正确的单一频率正弦波的频谱,仅在0.25处有1根单一谱线。如图(4a)和(4b)所示。x5(n)=cos(n/4)+cos(n/8)的周期为16,所以N=8不是其周期的整数倍,得到的频谱不正确,如图(5a)所示。N=16是其一个周期,得到正确的频谱,仅在0.25和0.125处有2根单一谱线,如图(5b)所示。,3.实验内容(3)对模拟周期信号谱分析。x6(t)=cos8t+cos16t+cos20t x6(t)有3个频率成分,f1=4 Hz,f2=8 Hz,f3=10 Hz。所以x6(t)的周期为0.5 s。采样频率Fs=64 Hz=16f1=8f2=6.4f3。变换区间N=16时,观察时间Tp=16T=0.25 s,不是x6(t)的整数倍周期,所以所得频谱不正确,如图(6a)所示。变换区间N=32,64 时,观察时间Tp=0.5 s,1 s,是x6(t)的整数周期,所以所得频谱正确,如图(6b)和(6c)所示。图中3根谱线正好分别位于4、8、10 Hz处。,变换区间N=64 时,频谱幅度是变换区间N=32时的2倍,这种结果正好验证了用DFT对中期序列谱分析的理论。注意:(1)用DFT(或FFT)对模拟信号分析频谱时,最好将X(k)的自变量k换算成对应的模拟频率fk并作为横坐标绘图,以便于观察频谱。这样,不管变换区间N取信号周期的几倍,画出的频谱图中有效离散谐波谱线所在的频率值不变,如图(6b)和(6c)所示。,(2)本程序直接画出采样序列N点DFT的模值,实际上分析频谱时最好画出归一化幅度谱,这样就避免了幅度值随变换区间N变化的缺点。本实验程序这样绘图只是为了验证用DFT对序列谱分析的理论。思考题(1)和(2)的答案请读者在教材节中去找,思考题(3)的答案在程序运行结果分析讨论中已经详细回答。,8.4 实验四:IIR数字滤波器设计及软件实现8.4.1 实验指导1 实验目的(1)熟悉用双线性变换法设计IIR数字滤波器的原理与方法。(2)学会调用MATLAB信号处理工具箱中滤波器设计函数(或滤波器设计分析工具fdatool)设计各种IIR数字滤波器,学会根据滤波需求确定滤波器指标参数。(3)掌握IIR数字滤波器的MATLAB实现方法。(4)通过观察滤波器输入输出信号的时域波形及其频谱,建立数字滤波的概念。,2 实验原理设计IIR数字滤波器一般采用间接法(脉冲响应不变法和双线性变换法),应用最广泛的是双线性变换法。其基本设计过程是:先将给定的数字滤波器的指标转换成过渡模拟滤波器的指标;设计过渡模拟滤波器;将过渡模拟滤波器的系统函数转换成数字滤波器的系统函数。MATLAB信号处理工具箱中的各种IIR数字滤波器设计函数都是采用双线性变换法。,教材第6章介绍的滤波器设计函数butter、cheby1、cheby2 和ellip可以分别被调用来直接设计巴特沃斯、切比雪夫1、切比雪夫2以及椭圆模拟和数字滤波器。本实验要求读者调用如上函数直接设计IIR数字滤波器。本实验的数字滤波器的MATLAB实现是指调用MATLAB信号处理工具箱函数filter对给定的输入信号x(n)进行滤波,得到滤波后的输出信号y(n)。,3.实验内容及步骤(1)调用信号产生函数mstg产生由三路抑制载波调幅信号相加构成的复合信号s(t),该函数还会自动绘图显示s(t)的时域波形和幅频特性曲线,如图所示。由图可见,三路信号时域混叠无法在时域分离。但频域是分离的,所以可以通过滤波的方法在频域分离,这就是本实验的目的。,图8.4.1 三路调幅信号s(t)的时域波形和幅频特性曲线,(2)要求将s(t)中三路调幅信号分离,通过观察s(t)的幅频特性曲线,分别确定可以分离s(t)中三路抑制载波单频调幅信号的三个滤波器(低通滤波器、带通滤波器、高通滤波器)的通带截止频率和阻带截止频率。要求滤波器的通带最大衰减为0.1 dB,阻带最小衰减为60 dB。提示:抑制载波单频调幅信号的数学表示式为,其中,cos(2fct)称为载波;fc为载波频率;cos(2f0t)称为单频调制信号;f0为调制正弦波信号频率,且满足fcf0。由上式可见,所谓抑制载波单频调幅信号,就是两个正弦信号相乘,它有两个频率成分:和频fc+f0和差频fcf0,这两个频率成分关于载波频率fc对称。所以,1路抑制载波单频调幅信号的频谱图是关于载波频率fc对称的2根谱线。容易看出,图中三路调幅信号的载波频率分别为250 Hz、500 Hz、1000 Hz。,(3)编程序调用MATLAB滤波器设计函数ellipord和ellip分别设计这三个椭圆滤波器,并绘图显示其幅频响应特性曲线。(4)调用滤波器实现函数filter,用三个滤波器分别对信号产生函数mstg产生的信号s(t)进行滤波,分离出s(t)中的三路不同载波频率的调幅信号y1(n)、y2(n)和y3(n),并绘图显示y1(n)、y2(n)和y3(n)的时域波形,观察分离效果。,4 信号产生函数mstg清单function st=mstg%产生信号序列向量st,并显示st的时域波形和频谱%st=mstg 返回三路调幅信号相加形成的混合信号,长度N=1600N=1600%N为信号st的长度Fs=10000;T=1/Fs;Tp=N*T;%采样频率Fs=10 kHz,Tp为采样时间t=0:T:(N1)*T;k=0:N1;f=k/Tp;fc1=Fs/10;%第1路调幅信号的载波频率fc1=1000 Hzfm1=fc1/10;%第1路调幅信号的调制信号频率fm1=100 Hz,fc2=Fs/20;%第2路调幅信号的载波频率fc2=500 Hzfm2=fc2/10;%第2路调幅信号的调制信号频率fm2=50 Hzfc3=Fs/40;%第3路调幅信号的载波频率fc3=250 Hzfm3=fc3/10;%第3路调幅信号的调制信号频率fm3=25 Hzxt1=cos(2*pi*fm1*t).*cos(2*pi*fc1*t);%产生第1路调幅信号xt2=cos(2*pi*fm2*t).*cos(2*pi*fc2*t);%产生第2路调幅信号xt3=cos(2*pi*fm3*t).*cos(2*pi*fc3*t);%产生第3路调幅信号,st=xt1+xt2+xt3;%三路调幅信号相加fxt=fft(st,N);%计算信号st的频谱%以下为绘图部分,绘制st的时域波形和幅频特性曲线subplot(3,1,1)plot(t,st);grid;xlabel(t/s);ylabel(s(t);axis(0,Tp/8,min(st),max(st);title(a)s(t)的波形)subplot(3,1,2),stem(f,abs(fxt)/max(abs(fxt),.);grid;title(b)s(t)的频谱)axis(0,Fs/5,0,1.2);xlabel(f/Hz);ylabel(幅度)5 实验程序框图实验程序框图如图所示,供读者参考。,图8.4.2,6 思考题(1)请阅读信号产生函数mstg,确定三路调幅信号的载波频率和调制信号频率。(2)信号产生函数mstg中采样点数N=800,对st进行N点FFT可以得到6根理想谱线。如果取N=1000,可否得到6根理想谱线?为什么?N=2000呢?请改变函数mstg中采样点数N的值,观察频谱图,验证您的判断是否正确。(3)修改信号产生函数mstg,给每路调幅信号加入载波成分,产生调幅(AM)信号,重复本实验,观察AM信号与抑制载波调幅信号的时域波形及其频谱的差别。提示:AM信号表示式:s(t)=1+cos(2f0t)cos(2fct),7 实验报告要求(1)简述实验目的及原理。(2)画出实验主程序框图,打印程序清单。(3)绘制三个分离滤波器的损耗函数曲线。(4)绘制经过滤波分离出的三路调幅信号的时域波形。(5)简要回答思考题。,8.4.2 滤波器参数及实验程序清单1.滤波器参数选取观察图可知,三路调幅信号的载波频率分别为250 Hz、500 Hz、1000 Hz。带宽(也可以由信号产生函数mstg清单看出)分别为50 Hz、100 Hz、200 Hz。所以,分离混合信号st中三路抑制载波单频调幅信号的三个滤波器(低通滤波器、带通滤波器、高通滤波器)的指标参数选取如下:对载波频率为250 Hz的调幅信号,可以用低通滤波器分离,其指标为通带截止频率fp=280 Hz,通带最大衰减ap=0.1 dB;阻带截止频率fs=450 Hz,阻带最小衰减as=60 dB。,对载波频率为500 Hz的调幅信号,可以用带通滤波器分离,其指标为通带截止频率fpl=440 Hz,fpu=560 Hz,通带最大衰减ap=0.1 dB;阻带截止频率fsl=275 Hz,fsu=900 Hz,阻带最小衰减as=60 dB。对载波频为1000 Hz的调幅信号,可以用高通滤波器分离,其指标为通带截止频率fp=890 Hz,通带