龔 濤
·(中國艦船研究設(shè)計中心,武漢 430064)
基于Matlab的IIR數(shù)字低通濾波器的設(shè)計方法
龔 濤
·(中國艦船研究設(shè)計中心,武漢 430064)
詳細的闡述了IIR數(shù)字濾波器設(shè)計的一般步驟,首先根據(jù)給定的目標數(shù)字濾波器的性能指標,將其轉(zhuǎn)換為模擬濾波器的指標參數(shù),再根據(jù)巴特沃思(Butterworth)低通模擬濾波器的設(shè)計原理設(shè)計出相應技術(shù)指標下的模擬低通濾波器,最后采用雙線性變換將設(shè)計出的模擬低通濾波器轉(zhuǎn)換為數(shù)字濾波器,即得到了期望的IIR數(shù)字濾波器。針對實例,在Matlab環(huán)境下分別設(shè)計出了4階、11階、25階的IIR數(shù)字濾波器,通過仿真驗證了濾波器的速度與濾波效果之間不是完全的統(tǒng)一。最后,分析了雙線性變換頻率畸變現(xiàn)象,造成相同頻率下離散后的幅值變小。通過頻率預曲折方法對雙線性變換進行相應的補償。以25階濾波器為例,仿真結(jié)果顯示,在指定頻率處實現(xiàn)了幅值在轉(zhuǎn)換前后大小不變的效果。
數(shù)字濾波器,雙線性變換,頻率畸變,預曲折
實現(xiàn)數(shù)字化是控制系統(tǒng)的重要發(fā)展方向,而數(shù)字信號處理已在通信、語音、圖像、自動控制、雷達、軍事、航空航天等領(lǐng)域廣泛應用。數(shù)字信號處理方法通常涉及變換、濾波、頻譜分析、編碼解碼等處理。數(shù)字濾波是重要環(huán)節(jié),它能滿足濾波器對幅度和相位特性的嚴格要求,能克服模擬濾波器所無法解決的電壓和溫度漂移以及噪聲等問題。
數(shù)字濾波器從實現(xiàn)方法上有無限脈沖響應(IIR)濾波器和有限脈沖響應(FIR)濾波器之分,它們的設(shè)計一般都包括如下3個步驟:
(1)給出所需要的濾波器技術(shù)指標;
(2)設(shè)計離散系統(tǒng)傳遞函數(shù)H(z)使其逼近所需要的技術(shù)指標;
(3)實現(xiàn)傳遞函數(shù)H(z)設(shè)計系統(tǒng)。
FIR濾波器可以對給定的頻率特性直接進行設(shè)計,而IIR濾波器目前最通用的方法是利用已經(jīng)很成熟的模擬濾波器的設(shè)計方法來進行設(shè)計。而模擬濾波器的設(shè)計方法又有Butterworth、Chebyshev(I型、II型)、橢圓濾波器等不同的設(shè)計方法。模擬濾波器的設(shè)計已經(jīng)有了一套相當成熟的方法,它不但有完整的設(shè)計公式,而且還有較為完整的圖表供查詢,因此,充分利用這些已有的資源將會給數(shù)字濾波器的設(shè)計帶來很大的方便。
本文針對一個實例采用雙線性z變換設(shè)計IIR低通濾波器,具體設(shè)計步驟如下:
(1)根據(jù)雙線性變換s平面到z平面的映射關(guān)系,將給出的數(shù)字濾波器的性能指標轉(zhuǎn)換為模擬濾波器的技術(shù)指標;
(2)依據(jù)轉(zhuǎn)換后的技術(shù)指標設(shè)計巴特沃斯(Butterworth)模擬低通濾波器G(s);
(3)再按一定規(guī)則將G(s)轉(zhuǎn)換為H(z),采用雙線性Z變換。
設(shè)計流程圖如下所示:
圖1 雙線性Z變換IIR數(shù)字低通濾波器設(shè)計流程圖
數(shù)字濾波器設(shè)計首先要找到合適的映射關(guān)系,將其離散時間傳遞函數(shù)H(z)從z平面映射到連續(xù)時間傳遞函數(shù)的s復平面上,這種關(guān)系保證:
(1)s平面的整個jΩ軸只映射到z平面的單位圓上;
(2)若G(s)是穩(wěn)定的,由G(s)映射得到的H(z)也應該是穩(wěn)定的;
(3)這種映射是可逆的,既能由G(s)得到H(z),也能由H(z)得到G(s);
(4)如果G(j0)=1,那么H(ej0)也應等于1。滿足以上4個條件的映射關(guān)系為:
此關(guān)系稱為雙線性z變換,由此關(guān)系不難求出:
用s平面和z平面中頻率表示,由(1)中映射關(guān)系可以得到以下頻率映射關(guān)系:
式(4)給出了s平面和z平面之間的頻率映射關(guān)系。由此可知,s平面的虛軸對應z平面的單位圓。當Ω由0變至+∞時,z平面中ω相應的由0變至π;當Ω由0變至-∞時,z平面ω則相應的由0變至-π。這種映射正是利用了正切函數(shù)的非線性特點,把整個jΩ軸壓縮到了單位圓的一周上。
這樣,當給定了數(shù)字濾波器的技術(shù)指標:通帶截止頻率ωp、阻帶截止頻率ωs、峰值通帶波紋αp、最小阻帶衰減αs后,根據(jù)(4)式,有:
Butterworth低通濾波器的設(shè)計可按以下3個步驟來進行。
(1)將實際頻率Ω歸一化得到歸一化的幅平方特性。
由此可以看出,在|G(jλ)|2中只有兩個參數(shù)C和N,N是濾波器的階次。
(2)依據(jù)由目標數(shù)字濾波器技術(shù)指標轉(zhuǎn)換后得到的模擬濾波器技術(shù)指標,求出參數(shù)C和N。
由式(9)可得:
這樣C和N便可求出。
3)確定G(s)
若令αp=3dB,則C=1,則Butterworth濾波器形式如下:
因為p=jλ,有
由上式解可以看出,2N個極點將均勻分布在s平面的半徑為1的圓上,相距(π/Ν)rad。為了保證所設(shè)計的濾波器是穩(wěn)定的,應將前N個極點,即位于s平面左半平面的N個極點賦予G(p),即
這樣,得到了歸一化的轉(zhuǎn)移函數(shù)G(p)。因為p=jλ=jΩ/Ωp=s/Ωp,所以,在求得G(p)后,用s/Ωp代替變量p,即得到實際需要的G(s)。
上述Butterworth低通模擬濾波器G(s),即是由期望的數(shù)字濾波器技術(shù)指標轉(zhuǎn)換后而求出的相應的模擬濾波器。利用雙線性z變換,即可獲得具備期望技術(shù)指標的IIR數(shù)字低通濾波器H(z),即
由于雙線性變換采用了正切函數(shù)的非線性特點,將s平面jΩ軸一一映射到z平面的單位圓周上,高度壓縮必然導致頻率映射的非線性,即出現(xiàn)畸變現(xiàn)象。
以一些具體的數(shù)據(jù)可以看出,s平面及z平面中出現(xiàn)的頻率畸變?nèi)缦卤硭荆?/p>
表1 具體數(shù)據(jù)舉例
由上表數(shù)據(jù)可知,
圖2 雙線性變換出現(xiàn)的頻率畸變
為了減小畸變的影響,就應修正雙線性變換,以確保在某一頻率處變換前后的幅值相同,即。具體步驟如下:
1)計算預修正頻率ω?
2)將D(s)修正為
3)采用頻率預曲折的方法對雙線性變換進行相應的補償,修正后的雙線性變換表達式如下:
此節(jié),針對具體實例進 IIR數(shù)字濾波器的設(shè)計。對一個含有2、3、5次諧波和高斯噪聲的交流電流進行濾波,濾除3次及5次諧波,得到只含有基波和2次諧波的電流,采樣頻率是1kHz,具體參數(shù)如下表。
表2 含諧波、噪聲電流具體參數(shù)
詳細程序設(shè)計:
運行程序,仿真波形如下所示:
圖3 輸入信號的時域波形與頻譜
圖4 fp=105Hz,fs=145Hz,αp=3dB,αs=12dB 技術(shù)指標下的4階IIR數(shù)字濾波器濾波效果
圖5 fp=105Hz,fs=145Hz,αp=3dB,αs=12dB 技術(shù)指標下數(shù)字濾波器的幅頻特性
在以上給定技術(shù)指標下,得到的4階IIR數(shù)字濾波器已經(jīng)將5次諧波完全濾除,3次諧波大大衰減,但是仍含有部分3次諧波。為了增強濾波效果,獲得更加陡的幅頻特性,縮短了fp與fs之間的過渡區(qū)域,并增大了阻帶頻率開始處的下降增益αs,由濾波器的設(shè)計過程可知,所獲得濾波器階數(shù)將會相應增大。技術(shù)指標修改如下:fp=110Hz,fs=135Hz,αp=3dB,αs=20dB獲得11階的IIR數(shù)字濾波器,濾波效果如下所示。
圖6 11階數(shù)字濾波器濾波效果
可以看出,信號中3次、5次諧波已經(jīng)完全被消除,從其頻譜圖中可以看出,在通帶0-100Hz以內(nèi),卻出現(xiàn)了幅值比較低的連續(xù)頻率信號。此時濾波器的幅頻特性如下所示。由幅頻特性可以看出,所設(shè)計的數(shù)字濾波器出現(xiàn)比較大的負的相位裕度,此濾波器如果加入實際系統(tǒng)將會造成系統(tǒng)的不穩(wěn)定。
圖7 11階濾波器幅頻特性
繼續(xù)增大濾波器階數(shù),給定數(shù)字濾波器技術(shù)指標如下:fp=110Hz,fs=125Hz,αp=3dB,αs=30dB,得到25階數(shù)字濾波器,其濾波效果及頻譜特性如下圖所示??芍?,出現(xiàn)的雜頻相對更加明顯,幅頻特性也出現(xiàn)更大的負的相位裕度。
圖8 25階數(shù)字濾波器濾波效果
圖9 25階濾波器幅頻特性
由以上仿真可知,濾波器的階數(shù)越高,濾波的幅頻特性越陡,即濾波速度越快。但是卻帶來了更大的負的相位裕度,同時也引入了在通帶0-100Hz內(nèi)的其他雜頻。可知,在濾波速度與效果上不是完全統(tǒng)一的。
由以上仿真結(jié)果可知,經(jīng)過濾波后,雖已濾除了3次、5次諧波,但是由于雙線性變換的頻率畸變導致2次諧波幅值有些變小。此時,要保證基波及2次諧波處幅值在離散前后不發(fā)生變化,需要利用頻率預曲折方法對雙線性變換進行相應的補償。指令如下:
[num1 den1]=C2DM(num,den,Ts,'prewarp',100);
即在頻率100Hz處保證離散前后幅值大小不發(fā)生變化。
圖10 修正后的25階濾波器濾波效果
可知,經(jīng)過頻率預曲折方法修正后的雙線性變換,在指定頻率處的幅值大小與原信號2次諧波幅值大小保持不變。達到補償?shù)男Ч?/p>
本文詳細的闡述了IIR數(shù)字濾波器設(shè)計的一般步驟,即首先根據(jù)給定的目標數(shù)字濾波器的技術(shù)指標,將其轉(zhuǎn)換為模擬濾波器的技術(shù)指標,再根據(jù)巴特沃思(Butterworth)低通模擬濾波器的設(shè)計原理設(shè)計出相應技術(shù)指標下的模擬低通濾波器,最后采用雙線性變換將設(shè)計出的模擬低通濾波器轉(zhuǎn)換為數(shù)字濾波器,即得到了期望技術(shù)指標的IIR數(shù)字濾波器。
針對實例,設(shè)計出了4階IIR數(shù)字濾波器,達到較為理想的濾波效果,3次、5次諧波已基本被濾除。通過改變技術(shù)指標,分別設(shè)計了11階、25階數(shù)字濾波器,通過仿真,驗證了濾波器的速度與濾波效果之間不是完全的統(tǒng)一。階數(shù)越高,濾波特性越陡峭,速度越快,但同時造成較大的負的相位裕度,也引入通帶的其他雜頻。
最后,分析了雙線性變換頻率畸變現(xiàn)象,在z平面上相映射的頻率要小于其在s平面上所對應的頻率,并且相同頻率下z平面與s平面中所對應的幅值也不相同,造成相同頻率下離散后的幅值變小。通過頻率預曲折方法對雙線性變換進行相應的補償。以25階濾波器為例,仿真結(jié)果顯示,在指定頻率處實現(xiàn)了幅值在轉(zhuǎn)換前后大小不變的效果。
[1]丁美玉.數(shù)字信號處理.西安:西安電子科技大學出版社,2005
[2]S.K.Mitra.Digital Signal Processing–A Computer-Based Approach,Third Ed.,McGraw-Hill,2005
[3]鄭南寧.數(shù)字信號處理.北京:清華大學出版社,2007
[4]張德豐.MATLAB數(shù)字信號處理與應用.北京:清華大學出版社,2010
[5]陳懷琛.數(shù)字信號處理教程.北京:電子工業(yè)出版社,2009
[6]葛哲學.精通MATLAB.北京:國防工業(yè)出版社,2008
Designing and Implementation of Lowpass IIR Digital Filter
GONG Tao
(China Ship Development and Design Center,Wuhan 430064,China)
This paper described the general steps of IIR digital filter design.Firstly,the given specifications of the digital filter were converted to those of the equivalent analog filter.Then,a low pass analog filter with the specifications was implemented according to the Butterworth low-pass analog filter designing method.Finally,the analog filter was converted to the correspondent digital filter by using bilinear transform so that the desired IIR low-pass digital filter was derived.For verification,in the MATLAB environment the 4-order,11-order and 25-order IIR low-pass digital filters were designed and verified by simulation.The results have shown that the speed and the filtering performance cannot be balanced.Finally,the frequency distortion phenomenon of the bilinear transform was analyzed and it is demonstrated that the amplitude in the same frequency will be small after discretization.Therefore,the Prewarping method can be used to compensate the bilinear transform appropriately.Simulations of a 25-order digital filter verify the validity of the compensation effect of the frequency prewarping.
digital filter,bilinear transform,frequency distortion,prewarping
U665.14
A
龔濤(1975-), ,碩士,高級工程師,研究方向:電氣工程。