ch4_4双线性变换法-北京交通大学陈后金教授信号处理课件

1、近代数字信号处理近代数字信号处理(Advanced Digital Signal Processing) 电子信息工程学院电子信息工程学院信号与图像处理研究室如何将模拟滤波器转变为数字滤波器如何将模拟滤波器转变为数字滤波器?W Wp,W Wsw wp,w wsH(s)H(z)频率频率变换变换设计模拟设计模拟滤波器滤波器AF到到DF的转换的转换利用利用AF-BW filter及脉冲响应不变法设计一及脉冲响应不变法设计一DF，满足，满足 W Wp=0.2p p, W Ws=0.6p p, Ap 2dB, As 15dB 。01-18-16-14

2、-12-10-8-6-4-20Normalized frequencyGain,dBAs = 14.2dB )(sH) (sH)(zHwwW/,/TT)2arctan(2TTwwT wW,)2tan(2WwT模拟频率与数字频率的关系为 )2tan(2WwT)2cos()2sin(2j)2tan(2jjWWWwTTwjsWjez11112zzTssTsTz/2/22j2j2j2jeeee2WWWWTWWjje1e12TWWwjje1e12jT双线性变换双线性变换 11112zzTssTsTz222222)/2()/2(wwTTz令令s= +jw w，则有，则有1) 1) 0,0,|z| 0,0,

3、|z|1S域虚轴映射到域虚轴映射到z域单位圆上域单位圆上S域右半平面映射到域右半平面映射到z域单位圆外域单位圆外) 2/tan(2WTwWp)(WjeH)(wjHpWsWWpwsww)2(tan2Tw幅度响应不是常数时会产生幅度失真幅度响应不是常数时会产生幅度失真无混叠无混叠1. 1. 将数字滤波器的频率指标将数字滤波器的频率指标 W Wk 转换为转换为 模拟滤波器的频率指标模拟滤波器的频率指标 w wk 2. . 由模拟滤波器的指标设计模拟滤波器的由模拟滤波器的指标设计模拟滤波器的H(s)。3. 利用双线性变换法，将利用双线性变换法，将H(s)转换转换H(z)。)2tan(2kkTWw111

4、12)()(zzTssHzHW Wp,W Wsw wp,w wsH(s)H(z)设计模拟设计模拟滤波器滤波器双线性变换双线性变换：AF分子、分母多项式的系数向量分子、分母多项式的系数向量Fs：抽样频率：抽样频率：DF分子、分母多项式的系数向量分子、分母多项式的系数向量利用利用MATLAB)2tan(2WwT11112)()(zzTssHzH 利用利用BW型模拟低通滤波器和双线性变换法设计满足指型模拟低通滤波器和双线性变换法设计满足指标标W Wp=p p/3，Ap=3dB，N=1的数字低通滤波器，并与的数字低通滤波器，并与脉冲响应不变法设计的脉冲响应不变法设计的DF比较。比较。设双线性变换中的参

5、数为设双线性变换中的参数为T(1) 将将DF的频率指标转换为的频率指标转换为AF的频率指标的频率指标)2tan(2ppWwT(2) 设计设计3dB截频为截频为w wp的一阶的一阶BW型模拟低通滤波器，即型模拟低通滤波器，即N=1, w wc = w wp1/1)(cwssH1/1pws故故1)2/tan(21pWsT 利用利用BW型模拟低通滤波器和双线性变换法设计满足指型模拟低通滤波器和双线性变换法设计满足指标标W Wp=p p/3，Ap=3dB，N=1的数字低通滤波器，并与的数字低通滤波器，并与脉冲响应不变法设计的脉冲响应不变法设计的DF比较。比较。设双线性变换中的参数为设双线性变换中的参数

6、为T(3) 用双线性变换法将模拟滤波器转换为数字滤波器用双线性变换法将模拟滤波器转换为数字滤波器参数参数T的取值和最终的设计结果无关。的取值和最终的设计结果无关。 为简单起见，一般取为简单起见，一般取T=2 1)2/tan(21)(pWsTsH1pp1p) 1)2(tan()2tan(1)1)(2tan()(z/z/zHWWW11112zzTs112697. 01366. 0366. 0zz 利用利用BW型模拟低通滤波器和双线性变换法设计满足指型模拟低通滤波器和双线性变换法设计满足指标标W Wp=p p/3，Ap=3dB，N=1的数字低通滤波器，并与的数字低通滤波器，并与脉冲响应不变法设计的脉

7、冲响应不变法设计的DF比较。比较。双线性变换法设计的双线性变换法设计的DF的系统函数为的系统函数为112697. 01)366. 0366. 0)(zzzH双脉冲响应不变法设计的脉冲响应不变法设计的DF的系统函数为的系统函数为1ppe1e1)(zzHWW脉令令z=ejW W ，可分别获得两者的幅度响应。，可分别获得两者的幅度响应。1/3e1/3z 利用利用BW型模拟低通滤波器和双线性变换法设计满足型模拟低通滤波器和双线性变换法设计满足指标指标W Wp=p p/3，Ap=3dB，N=1的数字低通滤波器，的数字低通滤波器，并与脉冲响应不变法设计的并与脉冲响应不变法设计的DF比较。比较。脉冲响应不变

8、法脉冲响应不变法双线性变换法双线性变换法 脉冲响应不变法存在频脉冲响应不变法存在频谱混叠，所设计的谱混叠，所设计的DF不满不满足给定指标。而双线性变足给定指标。而双线性变换法不存在频谱混叠，所换法不存在频谱混叠，所设计的设计的DF满足给定指标。满足给定指标。 利用利用BW型模拟低通滤波器和双线性变换法设计满足指型模拟低通滤波器和双线性变换法设计满足指标标W Wp=p p/3，Ap=3dB，N=1的数字低通滤波器，并与的数字低通滤波器，并与脉冲响应不变法设计的脉冲响应不变法设计的DF比较。比较。H双双(z)和和H脉脉(z)幅度响应比较的幅度响应比较的MATLAB实现实现Wp=pi/3;b=1-e

9、xp(-Wp);b1=tan(Wp/2)*1 1;a=1 -exp(-Wp);a1=1+tan(Wp/2) tan(Wp/2)-1;w=linspace(0,pi,512);h=freqz(b,a,w);h1=freqz(b1,a1,w);plot(w/pi,(abs(h),w/pi,(abs(h1) );xlabel(Normalized frequency);ylabel(Amplitude);set(gca,ytick,0 0.7 1);set(gca,xtick,0 Wp/pi 1);grid;利用利用AF-BW filter及及双线性变换法双线性变换法设计一设计一DF，满足，满足 W

10、 Wp=0.2p p, W Ws=0.6p p, Ap 2dB, As 15dB (1) 将数字低通指标转换成模拟低通指标，取将数字低通指标转换成模拟低通指标，取T=2 Ap2db, As15db (2) 设计模拟低通滤波器设计模拟低通滤波器 （BW型）型）)/lg(2)110110(lgps1 . 01 . 0pswwAAN=2)2/(11 . 0sc) 110(sNAww=0.585 1 12)(1)(c2cLwwsssH9324. 0)2tan(2ppWwT4376. 1)2tan(2ssWwT3342. 05827. 03342. 02ss利用利用AF-BW filter及及双线性变换

11、法双线性变换法设计一设计一DF，满足，满足 W Wp=0.2p p, W Ws=0.6p p, Ap 2dB, As 15dB (3) 用用双线性变换法将模拟低通滤波器转换成数字低通滤波器双线性变换法将模拟低通滤波器转换成数字低通滤波器 11112)()(zzTssHzH21213237. 02606. 018157. 05315. 08157. 0zzzz利用利用AF-BW filter及及双线性变换法双线性变换法设计一设计一DF，满足，满足 W Wp=0.2p p, W Ws=0.6p p, Ap 2dB, As 15dB %Design DF BW low-pass filter usi

12、ng impulse invariance%DF BW LP specficationWp=0.2*pi; Ws=0.6*pi; Ap=2; As=15;T=2;Fs=1/T; %Sampling frequency(Hz)%Analog Butterworth specficationwp=2*tan(Wp/2)/T;ws=2*tan(Ws/2)/T;%determine the order of AF filter and the 3-dB cutoff frequency N,wc=buttord(wp,ws,Ap,As,s)%determine the AF-BW filternuma

13、,dena=butter(N,wc,s)利用利用AF-BW filter及及双线性变换法双线性变换法设计一设计一DF，满足，满足 W Wp=0.2p p, W Ws=0.6p p, Ap 2dB, As 15dB %determine the DF filternumd,dend=bilinear(numa,dena,Fs)%plot the frequency responsew=linspace(0,pi,1024);h=freqz(numd,dend,w);plot(w/pi,20*log10(abs(h);axis(0 1 -50 0);grid;xlabel(Normalized f

14、requency);ylabel(Gain,dB);%computer Ap As of the designed filterw=Wp Ws;h=freqz(numd,dend,w);fprintf(Ap= %.4fn,-20*log10( abs(h(1);fprintf(As= %.4fn,-20*log10( abs(h(2);利用利用AF-BW filter及及双线性变换法双线性变换法设计一设计一DF，满足，满足 W Wp=0.2p p, W Ws=0.6p p, Ap 2dB, As 15dB 01-50-45-40-35-30

15、-25-20-15-10-50Normalized frequencyGain,dbAp= 0.3945As= 15.0000利用利用AF-BW filter及及双线性变换法双线性变换法设计一设计一DF，满足，满足 W Wp=0.2p p, W Ws=0.6p p, Ap 2dB, As 15dB 01-50-45-40-35-30-25-20-15-10-50Normalized frequencyGain,dB将将双线性变换法与脉冲响应不变法所双线性变换法与脉冲响应不变法所设计设计DF的结果比较。的结果比较。双线性变换双线性变换Ap=

16、0.3945As= 15.0000脉冲响应不变法脉冲响应不变法双线性变换法双线性变换法脉冲响应不变脉冲响应不变Ap= 1.1187As= 12.3628W Wp,W Wsw wp,w wsH(z)w =W/T脉冲响应不变法脉冲响应不变法双线性变换法双线性变换法H(s)模拟频模拟频率变换率变换设计原型设计原型低通滤波器低通滤波器sp,ww)(LsH复频率变换复频率变换)2tan(2Tw 脉冲响应不脉冲响应不变法不能设计变法不能设计和和数字滤数字滤波器波器1e111zpsTpll11112zzTsW Wp,W WsH(z)w =W/T脉冲响应不变法脉冲响应不变法双线性变换法双线性变换法数字频数字频

17、率变换率变换设计原型设计原型低通滤波器低通滤波器sp,ww)(LsH)2tan(2WwT1e111zpsTpll11112zzTssp,WW)(LzHz域域变换变换试设计满足下列指标的试设计满足下列指标的BW型数字带阻滤波器型数字带阻滤波器 W Wp1=2.8113rad/s, W Wp2=2.9880rad/s, Ap 1dB , W Ws1=2.9203rad/s, W Ws2=2.9603rad/s, As 10dB 。 脉冲响应不变法不适合设计数字带阻滤波器，脉冲响应不变法不适合设计数字带阻滤波器，因此采用双线性变换法设计。因此采用双线性变换法设计。 (1) 将数字带阻滤波器指标转换成

18、模拟带阻滤波器指标将数字带阻滤波器指标转换成模拟带阻滤波器指标 )2tan(2WwT取取T=2，利用，利用得模拟带阻指标为得模拟带阻指标为w wp1=6rad, w wp2=13rad, w ws1=9rad, w ws2=1rad1,Ap 1dB, As 10dB 试设计满足下列指标的试设计满足下列指标的BW型数字带阻滤波器型数字带阻滤波器 W Wp1=2.8113rad/s, W Wp2=2.9880rad/s, Ap 1dB , W Ws1=2.9203rad/s, W Ws2=2.9603rad/s, As 10dB 。(2) 将模拟带阻滤波器指标转换成模拟低通滤波器指标将模拟带阻滤波

19、器指标转换成模拟低通滤波器指标 21s2swwB9499. 92s1s0www3714. 0,max202p2p2202p1p1pwwwwwwwBBAp 1dB, As 10dB 1sww wp1=6rad, w wp2=13rad, w ws1=9rad, w ws2=1rad1,Ap 1dB, As 10dB 模拟带阻指标模拟带阻指标试设计满足下列指标的试设计满足下列指标的BW型数字带阻滤波器型数字带阻滤波器 W Wp1=2.8113rad/s, W Wp2=2.9880rad/s, Ap 1dB , W Ws1=2.9203rad/s, W Ws2=2.9603rad/s, As 10d

20、B 。(3) 设计原型设计原型BW型模拟低通滤波器型模拟低通滤波器 Ap 1dB, As 10dB 1,3714. 0spww原型模拟低通指标原型模拟低通指标2)/(lg2)110110lg(ps1 . 01 . 0pswwAAN5774. 0) 110()2/(11 . 0scsNAww原型模拟低通滤波器的系统函数为原型模拟低通滤波器的系统函数为 12)(1)(c2cLwwsssH3333. 08165. 03333. 02ss试设计满足下列指标的试设计满足下列指标的BW型数字带阻滤波器型数字带阻滤波器 W Wp1=2.8113rad/s, W Wp2=2.9880rad/s, Ap 1dB

21、 , W Ws1=2.9203rad/s, W Ws2=2.9603rad/s, As 10dB (4) 由复频率变换将原型模拟低通转换为模拟带阻滤波器由复频率变换将原型模拟低通转换为模拟带阻滤波器 3333. 08165. 03333. 0)(2LsssH202)()(LBSwsBsssHsH9801485210899. 4980119823424ssssss(5) 由双线性变换模拟带阻滤波器转换成数字带阻滤波器由双线性变换模拟带阻滤波器转换成数字带阻滤波器 11112)()(zzTssHzH432143219067. 06412. 35601. 58241. 319522. 07327.

22、35624. 57327. 39522. 0zzzzzzzz 其中其中Wp, Ws为归一化角频率。为归一化角频率。 例例W Wp=0.1p p, 则则Wp=0.1 若为若为或或滤波器，则滤波器，则Wp=Wp1, Wp2; Ws=Ws1, Ws2 其中其中Wc=W1, W2 其中其中Wc=W1, W2 * 带通带通、带阻带阻中中W=W1, W2 * 带通带通、带阻带阻中中W=W1, W2 * 带通带通、带阻带阻中中W=W1, W2模拟低通原型模拟低通原型 (prototype) 频率变换频率变换综合设计综合设计IIR 滤波器最小阶数的计算滤波器最小阶数的计算离散化离散化利用利用MATLAB实现数

23、字带阻滤波器实现数字带阻滤波器 W Wp1=2.8113rad/s, W Wp2=2.9880rad/s, Ap 1dB , W Ws1=2.9203rad/s, W Ws2=2.9603rad/s, As 10dB。Wp=2.813,2.9880;Ws=2.9203,2.9603;Ap=1;As=10;N,Wc=buttord(Wp/pi,Ws/pi,Ap,As)numd,dend=butter(N,Wc,stop)运行结果运行结果N =2numd =0.9522, 3. 7327, 5.5624, 3.7327, 0.9522dend =1.0000, 3.8242 , 5.5601, 3

24、.6412, 0.9067432143219067. 06412. 35601. 58241. 319522. 07327. 35624. 57327. 39522. 0)(zzzzzzzzzH利用IIR数字滤波器处理脉搏信号 中医或西医领域的多种疾病诊断分析（心血管疾病等） 疲劳状态预测（司机疲劳驾驶等） 情绪变化分析（学生、孕妇情绪监测） 心理学分析（测谎、犯罪等心理研究）采集脉搏信号，经数字滤波、特征提取等处理与分析，可实现多种应用：清醒状态脉搏信号及频谱疲劳状态脉搏信号及频谱时域：波形密集频域：波动较大时域：波形稀疏频域：波动平缓疲劳检测利用IIR数字滤波器处理脉搏信号人体脉搏脉搏信号

25、采集系统输出信号往往含有噪声基线漂移、人体呼吸等低频干扰50Hz工频干扰其他高频干扰信号谱分析脉搏信号处理设计滤波器去噪特征检测相关应用 需要对采集的脉搏信号进行去噪处理!如何处理?利用IIR数字滤波器处理脉搏信号信号谱分析实现脉搏信号频谱%信号采样率参数fs=200; %A/D设备T=1/fs; %脉搏时域信号x=load(x.txt);M=length(x);t=0:M-1*T;figure(1);plot(t,x);xlabel(time/(s);ylabel(x(t);title(脉搏时域信号);grid on; %FFT频谱分析N=M*2-1;X=fft(x,N);figure(2)

26、;f=-M+1:M-1*fs/N;plot(f,abs(fftshift(X);xlabel(frequence/(Hz);ylabel(amplitude);title(脉搏频谱);grid on;脉搏时域波形利用IIR数字滤波器处理脉搏信号信号频谱脉搏信号频谱分析：结论：脉搏信号频率范围： 1-20Hz噪声信号频率范围： 0-1Hz，20Hz以上方案一： 带通滤波器方案二： 高通和低通滤波器级联基线漂移等低频干扰，频率小于1Hz工频干扰，频率为50Hz其他高频干扰20Hz以上如何去噪如何去噪?时域信号正常人的脉搏信号在120Hz的频率范围内利用IIR数字滤波器处理脉搏信号滤波器指标确定：滤

27、波器指标确定：滤除0-1Hz和20Hz以上频率范围噪声过渡带小，滤波器阶数高，实现成本高。1.带通滤波器指标2.高通与低通滤波器级联指标实际工程中常采用高通滤波器与低通滤波器级联方案。10.9Hzpw10.1Hzsw220.1Hzpw220.9Hzsw0.9Hzpw0.1Hzsw50Hzsw30Hzpw高通滤波器低通滤波器利用IIR数字滤波器处理脉搏信号数字滤波器实现数字滤波器实现%Butterworth高通滤波器滤波；Wp1= 0.9/(fs/2); Ws1= 0.1/(fs/2);n,Wn=buttord(Wp1,Ws1,0.5,40);b,a = butter(n,Wn,high);fi

28、gure(3);freqz(b,a,512,200); title(Butterworth Highpass Filter)y=filtfilt(b,a,x);figure(4);plot(t,y); xlabel(time/(s);ylabel(x(t);title(Butterworth高通滤波后脉搏信号);grid on;Y=fft(y,N);figure(5);plot(f,abs(fftshift(Y);xlabel(frequence/(Hz);ylabel(amplitude);title(FFT(幅频曲线);grid on;% Butterworth低通滤波器滤波；Wp = 3

29、0/(fs/2);Ws = 50/(fs/2);Rp = 3; Rs = 60;n,Wc = buttord(Wp,Ws,Rp,Rs);b,a = butter(n,Wc,low);figure(6);freqz(b,a,512,200); title(Butterworth Lowpass Filter)q=filter(b,a,y);figure(7);plot(t,q);xlabel(time/(s);ylabel(x(t);title(Butterworth低通滤波脉搏信号);grid on;Y=fft(q,N);figure(8);plot(f,abs(fftshift(Y);xla

30、bel(frequence/(Hz);ylabel(amplitude);title(FFT(幅频曲线);grid on;利用IIR数字滤波器处理脉搏信号ButterwothChebyshev-IEllipticN=2N=3N=3阶数最小相位近似线性低频噪声滤除效果好利用IIR数字滤波器处理脉搏信号N=3ButterwothChebyshev-IElliptic高频噪声滤除效果好N=2N=3利用IIR数字滤波器处理脉搏信号时域信号对比信号频谱对比原始信号滤波后信号原始信号频谱滤波后信号频谱利用IIR数字滤波器处理脉搏信号信号的时域波形加噪信号的时域波形问题研讨滤波器的幅度响应滤波后信号的频谱滤波后信号的波形问题研讨加噪信号的频谱信号的频谱问题研讨原信号时域波形含噪声信号时域波形问题