序列逆Z变换的Matlab实现知识分享.ppt

《序列逆Z变换的Matlab实现知识分享.ppt》由会员分享，可在线阅读，更多相关《序列逆Z变换的Matlab实现知识分享.ppt（37页珍藏版）》请在淘文阁 - 分享文档赚钱的网站上搜索。

1、一、序列逆一、序列逆Z变换变换(binhun)的的Matlab实现实现 函数函数residuez:residuez:适合适合(shh)(shh)计算离散系统有理函数的留数计算离散系统有理函数的留数和极点，可以用于求解序列的逆和极点，可以用于求解序列的逆Z Z变换。变换。函数函数residuezresiduez基本调用方式基本调用方式:r,p,c=residuez(b,a);r,p,c=residuez(b,a);输入参数输入参数:b=b0:b=b0，b1 b1，bM bM为分子多项式的系数为分子多项式的系数,a=a0a=a0，a1 a1，aN aN为分母多项式的系数，这些多项式都按为分母多项式

2、的系数，这些多项式都按z z的降幂排列的降幂排列 输出输出(shch)(shch)参数参数:r:r是极点的留数，是极点的留数，p p是极点，是极点，c c是无穷项多是无穷项多项式的系数项，仅当项式的系数项，仅当MNMN时存在。时存在。第一页，共37页。例：计算例：计算(j sun)逆逆Z变换变换 例例 ：计算：计算(j sun)(j sun)的逆的逆Z Z变换。变换。解解:有理分式有理分式X(z)X(z)分子和分母多项分子和分母多项式都按式都按z z的降幂的降幂(jin m)(jin m)排列。排列。b=0,1;a=2,-3,1;%b=0,1;a=2,-3,1;%多项式的系数多项式的系数 r,

3、p,c=residuez(b,a);%r,p,c=residuez(b,a);%求留数、极点和系数项求留数、极点和系数项 disp(disp(留数留数:);disp(r);%:);disp(r);%显示输出参数显示输出参数 disp(disp(极点极点:);disp(p);:);disp(p);disp(disp(系数项系数项:);disp(c);:);disp(c);程序运行结果为程序运行结果为 留数留数:1 -1:1 -1 极点极点:1.0000 0.5000:1.0000 0.5000 系数项系数项:X X(z z)的部分分式形式的部分分式形式为为逆逆Z Z变换变换为为第二页，共37页。

4、二、二、DFT物理意义物理意义(yy)的的Matlab实现实现序列的N点DFT的物理意义:对X(ej)在0，2上进行N点的等间隔取样。函数fft用于快速计算离散傅里叶变换，调用方式为y=fft(x);y=fft(x,N);y=fft(x)利用FFT算法计算序列x的离散傅里叶变换。当x为矩阵时，y为矩阵x每一列的FFT。当x长度为2的整数次幂时，函数fft采用基-2的FFT算法，否则采用混合基算法。y=fft(x,N)采用N点FFT。当序列x长度小于N时，函数fft自动(zdng)对序列尾部补零，构成N点数据；当x长度大于N时，函数fft自动(zdng)截取序列前面N点数据进行FFT。第三页，共

5、37页。函数ifft用于快速计算向量或矩阵的离散(lsn)傅里叶逆变换，与函数fft的调用规则基本相同。调用方式为y=ifft(x);y=ifft(x,N);第四页，共37页。例：利用例：利用(lyng)FFT实现线性卷积实现线性卷积例：利用例：利用FFT实现线性卷积。已知序列实现线性卷积。已知序列x(n)=R4(n)，求：，求：（1）用）用conv函数函数(hnsh)求求x(n)与与x(n)的线性卷的线性卷积积y(n)，并绘出图形；，并绘出图形；（2）用）用FFT求求x(n)与与x(n)的的4点循环卷积点循环卷积y1(n)，并，并绘出图形；绘出图形；（3）用）用FFT求求x(n)与与x(n)

6、的的8点循环卷积点循环卷积y2(n)，并，并将结果与（将结果与（1）比较，说明线性卷积与循环卷积之间）比较，说明线性卷积与循环卷积之间的关系。的关系。解解 程序程序(chngx)如下：如下：第五页，共37页。N1=4;N2=8;n1=0:1:N1-1;n2=0:1:N2-1;x=1,1,1,1;%构造构造(guzo)序列序列x(n)x1=1,1,1,1,0,0,0,0;%在序列在序列x(n)后补后补4个零个零figure(1)subplot(2,2,1)stem(n1,x),grid on;title(序列序列x(n)y1=conv(x,x);%y1为为x(n)与与x(n)的线性卷积的线性卷积

7、subplot(2,2,2)stem(0:1:length(y1)-1,y1),grid on;title(x(n)与与x(n)线性卷积线性卷积)第六页，共37页。X2=fft(x);%计算计算(j sun)x(n)与与x(n)的的4点循环卷积点循环卷积Y2=X2.*X2;y2=ifft(Y2);subplot(2,2,3)stem(n1,y2),grid on;title(x(n)与与x(n)的的4点循环卷积点循环卷积)X3=fft(x1);%计算计算(j sun)x(n)与与x(n)的的8点循环卷积点循环卷积Y3=X3.*X3;y3=ifft(Y3)subplot(2,2,4)stem(n

8、2,y3),grid on;title(x(n)与与x(n)的的8点循环卷积点循环卷积)第七页，共37页。程序运行结果程序运行结果(ji gu)图图第八页，共37页。三、三、频域取样频域取样(qyng)定理的定理的Matlab实现实现 例：例：设设x(n)=(0.7)nu(n)，在单位圆上以，在单位圆上以M=5和和M=20，对其，对其Z变换取样，研究时域信号受变换取样，研究时域信号受M变化的影响。变化的影响。（1）对）对x(n)进行进行(jnxng)Z变换；变换；（2）对）对X(z)进行进行(jnxng)等角取样，取样点数为等角取样，取样点数为M，求，求X(k)；（3）对）对X(k)进行进行(

9、jnxng)IDFT变化，得到变化，得到M点序列，请比点序列，请比较几个序列，并作分析。较几个序列，并作分析。解解 x(n)=(0.7)nu(n)的的Z变换变换(binhun)为为第九页，共37页。程序清单程序清单n=0:19;x=0.7.n;na=0:4;za=exp(j*2*pi*na/5);%在在z平面的单位平面的单位(dnwi)圆上对其进行圆上对其进行5点的等角距取样点的等角距取样Xa=za./(za-0.7);xa=abs(ifft(Xa);nb=0:19;zb=exp(j*2*nb*pi/20);%在在z平面的单位平面的单位(dnwi)圆上对其进行圆上对其进行20点的等角距取样点的

10、等角距取样Xb=zb./(zb-0.7);xb=abs(ifft(Xb);第十页，共37页。figure(1)subplot(2,2,1);%画出原始时域信号画出原始时域信号(xnho)stem(n,x)title(时域信号时域信号(xnho)x(n)subplot(2,2,2);xa=xa,xa,xa,xa;stem(n,xa)title(5点取样恢复的序列点取样恢复的序列)subplot(2,2,3);stem(n,xb)title(20点取样恢复的序列点取样恢复的序列)第十一页，共37页。程序运行结果程序运行结果(ji gu)第十二页，共37页。四、四、用用FFT进行进行(jnxng)谱

11、分析的谱分析的Matlab实现实现设模拟信号设模拟信号 ，以，以 t=0.01n(n=0:N-1)进行取样，试用进行取样，试用(shyng)fft函数对其做频谱分析。函数对其做频谱分析。N分别为：分别为：(1)N=45；(2)N=50；(3)N=55；(2)N=60。程序清单如下程序清单如下(rxi)%计算计算N=45的的FFT并绘出其幅频曲线并绘出其幅频曲线N=45;n=0:N-1;t=0.01*n;q=n*2*pi/N;x=2*sin(4*pi*t)+5*cos(8*pi*t);y=fft(x,N);figure(1)subplot(2,2,1)plot(q,abs(y)title(FFT

12、 N=45)第十三页，共37页。程序清单程序清单%计算计算(j sun)N=50的的FFT并绘出其幅频曲线并绘出其幅频曲线N=50;n=0:N-1;t=0.01*n;q=n*2*pi/N;x=2*sin(4*pi*t)+5*cos(8*pi*t);y=fft(x,N);figure(1)subplot(2,2,2)plot(q,abs(y)title(FFT N=50)第十四页，共37页。%计算计算(j sun)N=55的的FFT并绘出其幅频曲线并绘出其幅频曲线N=55;n=0:N-1;t=0.01*n;q=n*2*pi/N;x=2*sin(4*pi*t)+5*cos(8*pi*t);y=ff

13、t(x,N);figure(1)subplot(2,2,3)plot(q,abs(y)title(FFT N=55)第十五页，共37页。%计算计算(j sun)N=60的的FFT并绘出其幅频曲线并绘出其幅频曲线N=60;n=0:N-1;t=0.01*n;q=n*2*pi/N;x=2*sin(4*pi*t)+5*cos(8*pi*t);y=fft(x,N);figure(1)subplot(2,2,4)plot(q,abs(y)title(FFT N=60)第十六页，共37页。程序运行结果程序运行结果(ji gu)从图中可以看出从图中可以看出(kn ch)(kn ch)，这几种情况下均有较好，这

14、几种情况下均有较好的精度。的精度。第十七页，共37页。程序运行结果程序运行结果(ji gu)分析分析分析分析(fnx)：由：由t=0.01n进行取样可得，采样频率进行取样可得，采样频率fs=100Hz。而连续信号的最高模拟角频率为。而连续信号的最高模拟角频率为8 ，由由2 f可得，最高频率为可得，最高频率为8 /2=4Hz。因此，满足。因此，满足采样定理的要求。采样定理的要求。采样序列为采样序列为即即为周期为周期(zhuq)序列，周期序列，周期(zhuq)N=50。将程序中将程序中plot改为改为stem函数，则可以更清楚地看出频谱。函数，则可以更清楚地看出频谱。第十八页，共37页。修改修改(

15、xigi)程序运行结果程序运行结果第十九页，共37页。6.5 IIR数字滤波器的数字滤波器的Matlab仿真仿真(fn zhn)实现实现IIR数字(shz)滤波器设计 模拟滤波器到数字(shz)滤波器的转换 第二十页，共37页。五、五、IIR数字数字(shz)滤波器设计滤波器设计设数字滤波器系统函数(hnsh)为模拟滤波器的系统函数(hnsh)为函数(hnsh)butter和cheby1可以确定Butterworth和Chebyshev I型滤波器的系统函数(hnsh)。第二十一页，共37页。函数函数(hnsh)butter的调用格式的调用格式 函数butter的调用格式为b,a=butter

16、(n,Wc,)%设计数字(shz)Butterworth滤波器b,a=butter(n,Wc,ftype)%设计模拟Butterworth滤波器其中，n为滤波器阶数，Wc为截止频率。第二十二页，共37页。函数函数(hnsh)cheby1的调用格式的调用格式 函数cheby1的调用格式为b,a=cheby1(n,Rp,Wc)%设计数字Chebyshev滤波器b,a=cheby1(n,Rp,Wc,ftype)%设计模拟Chebyshev滤波器其中(qzhng)，n为滤波器阶数，Rp为通带内的纹波系数，Wc为截止频率。第二十三页，共37页。例：设计例：设计(shj)butterworth(shj)b

17、utterworth低通滤低通滤波器波器 设设计计一一模模拟拟(mn)butterworth(mn)butterworth低低通通滤滤波波器器，通通带带截截止止频频率率300Hz300Hz，通通带带最最大大衰衰减减2dB2dB，阻阻带带截截止止频频率率800Hz800Hz，阻阻带带最最小小衰衰减减30dB30dB。解解 滤波器的阶数和截止频率可由式和确定滤波器的阶数和截止频率可由式和确定(qudng)(qudng)，程序段为，程序段为 Wp=2*pi*300;Ws=2*pi*800;Rp=2;Rs=30;Wp=2*pi*300;Ws=2*pi*800;Rp=2;Rs=30;N=ceil(log

18、10(10(0.1*Rs)-1)/(10(0.1*Rp)-N=ceil(log10(10(0.1*Rs)-1)/(10(0.1*Rp)-1)/(2*log10(Ws/Wp);1)/(2*log10(Ws/Wp);Wc=Wp/(10(Rp/10)-1)(1/(2*N);Wc=Wp/(10(Rp/10)-1)(1/(2*N);b,a=butter(N,Wc,s);b,a=butter(N,Wc,s);freqs(b,a)freqs(b,a)第二十四页，共37页。程序运行结果程序运行结果(ji gu)运运行行(ynxng)(ynxng)程程序序，得得到到N=N=4 4，Wc=Wc=2.0157e+0

19、032.0157e+003。幅频特性和相频特性如图幅频特性和相频特性如图7.167.16所示。所示。第二十五页，共37页。模拟模拟(mn)滤波器到数字滤波器的转换滤波器到数字滤波器的转换设模拟滤波器系统函数为数字滤波器的系统函数为从模拟滤波器到数字滤波器的转换有两种方法(fngf)，即脉冲响应不变法和双线性变换法。(b,a分别为模拟滤波器的系统函数的系数，bz,az分别为数字滤波器的系统函数的系数)第二十六页，共37页。脉冲响应不变法脉冲响应不变法(bin f)脉冲响应不变法:用代换Ha(s)中的(s-sk)即可得到H(z)，从而将模拟滤波器转换(zhunhun)为数字滤波器格式。可用函数im

20、pinvar实现，调用格式为bz,az=impinvar(b,a,fs)其中，fs为取样频率。第二十七页，共37页。双线性变换法双线性变换法双线性变换法:用代换Ha(s)中的s即可得到H(z)，从而将模拟滤波器转换为数字滤波器格式。可用函数bilinear实现，调用格式为zd,pd,kd=bilinear(z,p,k,fs)其中，z，p，k和zd，pd，kd分别为s域和z域系统函数的零点(ln din)、极点和增益。第二十八页，共37页。例：模拟例：模拟(mn)(mn)滤波器转换数字滤波滤波器转换数字滤波器器利用(lyng)impinvar将一模拟低通滤波器变换成数字滤波器(取样频率为10Hz

21、)，程序段为b,a=butter(4,.3,s);bz,az=impinvar(b,a,10);程序运行结果(ji gu)为bz=1.0e-006*-0.0000 0.1324 0.5192 0.1273 0az=1.0000 -3.9216 5.7679 -3.7709 0.9246第二十九页，共37页。六、六、FIR数字数字(shz)滤波器的滤波器的Matlab仿仿真实现真实现窗函数(hnsh)法设计FIR滤波器FIR滤波器的优化设计第三十页，共37页。窗函数窗函数(hnsh)法设计法设计FIR滤波器滤波器窗函数(hnsh)法是通过对理想滤波器的单位取样响应加窗来逼近理想滤波器的。函数(h

22、nsh)fir1用于设计标准的低通、带通、高通和带阻滤波器。第三十一页，共37页。函数函数fir1的调用的调用(dioyng)格式格式 函数fir1的调用(dioyng)格式为 b=fir1(n,Wc,ftype,Windows)其中，n为滤波器阶数，Wc为截止频率ftype决定滤波器类型，ftype=high，设计高通FIR滤波器，ftype=stop，设计带阻FIR滤波器。Windows指定窗函数类型，默认为Hamming窗；可选Hanning、Hamming、Blackman、triangle、bartlett和boxcar窗，每种窗都可以由Matlab的相应函数生成。第三十二页，共37

23、页。例：设计例：设计(shj)butterworth(shj)butterworth低通滤波低通滤波器器 设计(shj)一个15阶FIR低通滤波器，截止频率为0.2。b=fir1(15,0.2);freqz(b,1,512);函数freqz(b,a,N)用于计算由a和b构成的数字滤波器的频率响应，并用图形(txng)方式分别表示其幅度响应和相位响应。第三十三页，共37页。程序运行结果程序运行结果(ji gu)程序运行结果(ji gu)如图所示。第三十四页，共37页。FIR滤波器的优化设计滤波器的优化设计(shj)前面介绍了FIR滤波器的优化设计方法，通过迭代的方法求解FIR滤波器，过程十分复杂

24、。在Matlab中，可以调用函数remez实现滤波器的设计。b=remez(n,f,m)函数remez采用Parks-McClellan算法设计线性相位(xingwi)FIR滤波器，n为滤波器阶数，其幅频特性由f和m指定。第三十五页，共37页。例：设计例：设计(shj)(shj)带通滤波器带通滤波器采采用用Parks-McClellanParks-McClellan算算法法设设计计一一个个1717阶阶带带通通滤滤波波器器，并并画画出出期期望望的的幅幅频频特特性性曲曲线线(qxin)(qxin)和和实实际际的的幅幅频频特特性曲线性曲线(qxin)(qxin)，程序段为f=0 0.3 0.4 0.6 0.7 1;m=0 0 1 1 0 0;b=remez(17,f,m);h,w=freqz(b,1,512);plot(f,m,w/pi,abs(h);第三十六页，共37页。程序运行结果程序运行结果(ji gu)程序运行结果(ji gu)如图所示。第三十七页，共37页。