第四章线性时不变离散时间系统的频域分析一、传输函数和频率响应例4.1传输函数分析Q4.1clear;M=input('EnterthefilterlengthM:');w=0:2*pi/1023:2*pi;num=(1/M)*ones(1,M);den=[1];h=freqz(num,den,w);subplot(2,1,1)plot(w/pi,abs(h));gridtitle('MagnitudeSpectrum|H(e^{j\omega})|')xlabel('\omega/\pi');ylabel('Amplitude');subplot(2,1,2)plot(w/pi,angle(h));gridtitle('PhaseSpectrumarg[H(e^{j\omega})]')xlabel('\omega/\pi');ylabel('Phaseinradians');M=2M=10M=15幅度谱为偶对称,相位谱为奇对称,这是一个低通滤波器。M越大,通带越窄且过渡带越陡峭。Q4.2使用修改后的程序P3.1,计算并画出当w=[0,pi]时传输函数的因果线性时不变离散时间系统的频率响应。它表示哪种类型的滤波器?w=0:pi/511:pi;num=[0.150-0.15];den=[1-0.50.7];如下图1这是一个带通滤波器。图1图2Q4.3对下面的传输函数重做习题Q4.2:,式(4.36)和式(4.37)给出的两个滤波器之间的区别是什么?你将选择哪一个滤波器来滤波,为什么?w=0:pi/511:pi;num=[0.150-0.15];den=[0.7-0.51];如上图2也是一个带通滤波器,这两个滤波器的幅度谱是一样的,相位谱不太一样,我会选择第一个带通滤波器,因为它的相位谱更加平滑,相位失真小。Q4.4使用MATLAB计算并画出当w=[0,pi]时因果线性时不变离散时间系统的群延迟。系统的传输函数为。clf;w=0:pi/511:pi;num=[1-1.21];den=[1-1.31.04-0.222];h=grpdelay(num,den,w);plot(w/pi,h);xlabel('w/pi');ylabel('群延迟');Q4.5使用Q3.50中编写的程序,分别计算并画出式(4.36)和式(4.37)确定的两个滤波器的冲激响应中的前一百个样本。讨论你的结果。clf;num=[0.150-0.15];den=[0.7-0.51];L=input('输入样本数L:');[gt]=impz(num,den,L);stem(t,g);title(['前',num2str(L),'脉冲响应的样本']);xlabel('时间序号n');ylabel('h[n]');(4.36)式(4.37)式由图可知:这些情节由impz给生成的因果的脉冲响应实现的H(z)。我们观察到Q4.3因果滤波器与H(z)在(4.36)稳定,这意味着H[n]是绝对可和,我们看到交替和指数衰减的脉冲响应。在另一方面,因果编档人员与H(z)在(4.37)极点以外的单位圆,是不稳定的。不足为奇的是,相应的h[n]上图显示与n指数增长。Q4.6传输函数的极零点图同样能分析线性时不变离散时间系统的性质。使用命令zplane可以很容易地得到系统的极零点图。使用zplane分别生成式(4.36)和式(4.37)确定的两个滤波器的极零点图。讨论你的结果。clf;num=[0.150-0.15];den=[1-0.50.7];[zpk]=tf2zpk(num,den);disp('Zeros:');disp(z);disp('Poles:');disp(p);input('Hitreturntocontinue...');[sosk]=zp2sos(z,p,k)input('Hitreturntocontinue...');zplane(z,p);式(4.36)式(4.37)由图可知:过滤器在(4.36)在单位圆和两极因此它的因果实现稳定;较低的图显示过滤器(4.37)极点在单位圆外,其因果关系的实现是不稳定的。二、传输函数的类型例4.2滤波器Q4.7clf;fc=0.25;n=[-6.5:1:6.5];y=2*fc*sinc(2*fc*n);k=n+6.5;stem(k,y);title('N=14');axis([013-0.20.6]);xlabel('Timeindexn');ylabel('Amplitude');grid;图1图2如图1低通有限冲激滤波器的长度为14,决定滤波器长度的语句为n=[-6.5:1:6.5],而控制截止频率的参数是fc=0.25。Q4.8fc=0.45;n=[-9.5:1:9.5];y=2*fc*sinc(2*fc*n);k=n+9.5;stem(k,y);title('N=20');axis([019-0.20.7]);xlabel('Timeindexn');ylabel('Amplitude');grid;修改参数fc和n,得到如上图2,可知低通有限冲激滤波器的长度变为20.Q4.9clf;fc=0.65;n=[-7.0:1:7.0];y=2*fc*sinc(2*fc*n);k=n+7.0;stem(k,y);title('N=14');axis([014-0.41.4]);xlabel('Timeindexn');ylabel('Amplitude');grid;Q4.10clear;N=input('EnterthefiltertimeshiftN:');No2=N/2;fc=0.25;n=[-No2:1:No2];y=2*fc*sinc(2*fc*n);w=0:pi/511:pi;h=freqz(y,[1],w);plot(w/pi,abs(h));grid;title(strcat('|H(e^{j\omega})|,N=',num2str(N)));xlabel('\omega/\pi');ylabel('Amplitude');上图依次分别为N=5,10,30,100的四幅图,从这四幅图可以看出随着阶数N的增大,低通滤波器的过渡带越来越窄,阻带衰减越来越快,滤波器越来越接近理想低通滤波器。Q4.11clf;M=2;num=ones(1,M)/M;[g,w]=gain(num,1);plot(w/pi,g);gridaxis([01-500.5])xlabel('\omega/\pi');ylabel('GainindB');title(['M=',num2str(M)])可以验证3dB截止频率在π/2处。Q4.12clear;K=input('EnterthenumberofsectionsK:');Hz=[1];fori=1:K;Hz=conv(Hz,[11]);end;Hz=(0.5)^K*Hz;[g,w]=gain(Hz,1);ThreedB=-3*ones(1,length(g));t1=2*acos((0.5)^(1/(2*K)))*ones(1,512)/pi;t2=-50:50.5/511:0.5;plot(w/pi,g,w/pi,ThreedB,t1,t2);grid;axis([01-500.5])xlabel('\omega/\pi');ylabel('GainindB');title(['K=',num2str(K),';Theoretical\omega_{c}=',num2str(t1(1))]);Q4.13clear;M=input('EnterthefilterlengthM:');n=0:M-1;num=(-1).^n.*ones(1,M)/M;[g,w]=gain(num,1);plot(w/pi,g);grid;axis([01-500.5]);xlabel('\omega/\pi');ylabel('GainindB');title(['M=',num2str(M)]);其3dB截止频率约为0.82piQ4.14设计一个在0.45pi处具有3dB截止频率wc的一阶无限冲激响应低通滤波器和一阶无限冲激响应高通滤波器。用MATLAB计算并画出它们的增益响应,验证设计的滤波器是否满足指标。用MATLAB证明两个滤波器是全通互补和功率互补的。Q4.15级联10个式(4.15)所示一阶无限冲激响应低通滤波器,设计一个在0.3pi处具有3dB截止频率wc的无限冲激响应低通滤波器。把它与一个具有相同截止频率的一阶无限冲激响应低通滤波器的增益响应作比较。Q4.16设计一个中心频率wo在0.61pi处、3dB带宽为0.51pi的二阶无限冲激响应带通滤波器。由于式(4.20)是α的二次方程,为了产生相同的3dB带宽,参数α将有两个数值,得到的传输函数HBP(z)也会有两个不同的表达式。使用函数zplane可产生两个传输函数的极零点图,从中可以选择一个稳定的传输函数。用MATLAB计算并画出你所设计的滤波器的增益响应,并验证它确实满足给定的条件。用设计的稳定无限冲激响应带通滤波器的传输函数的参数α和β,生成一个二阶无限冲激响应带阻滤波器的传输函数HBS(z)。用MATLAB证明HBP(z)和HBS(z)都是全通互补和功率互补的。Q4.17用MATLAB计算并画出一个梳状滤波器的幅度响应,该梳状滤波器是在L取不同值的情况下,由式(4.40)给出的原型有限冲激响应低通滤波器得到的。证明新滤波器的幅度响应在处有L个极小值,在处有L个极大值,k=0,1,2,3......,L-1.Q4.18用MATLAB计算并画出一个梳状滤波器的幅度响应,该梳状滤波器是在L取不同值的情况下,由式(4.42)在M=2时给出的原型有限冲激响应低通滤波器得到的。确定这种梳状滤波器冲激响应的极大值和极小值的位置。从这些情节我们观察,梳状滤波器极距为1kπ/L,山峰为(2k+1)π/L.Q4.19clf;b=[1-8.530.5-63];num1=[b81fliplr(b)];num2=[b8181fliplr(b)];num3=[b0-fliplr(b)];num4=[b81-81-fliplr(b)];n1=0:length(num1)-1;n2=0:length(num2)-1;subplot(2,2,1);stem(n1,num1);xlabel('Timeindexn');ylabel('Amplitude');grid;title('Type1FIRFilter');subplot(2,2,2);stem(n2,num2);xlabel('Timeindexn');ylabel('Amplitude');grid;title('Type2FIRFilter');subplot(2,2,3);stem(n1,num3);xlabel('Timeindexn');ylabel('Amplitude');grid;title('Type3FIRFilter');subplot(2,2,4);stem(n2,num4);xlabel('Timeindexn');ylabel('Amplitude');grid;title('Type4FIRFilter');pausesubplot(2,2,1);zplane(num1,1);title('Type1FIRFilter');subplot(2,2,2);zplane(num2,1);title('Type2FIRFilter');subplot(2,2,3);zplane(num3,1);title('Type3FIRFilter');subplot(2,2,4);zplane(num4,1);title('Type4FIRFilter');disp('ZerosofType1FIRFilterare');disp(roots(num1));disp('ZerosofType2FIRFilterare');disp(roots(num2));disp