摘 要:在使用定點DSP芯片實現(xiàn)IIR數(shù)字濾波器的工程實踐中,對于給定的濾波器設計指標,當取樣頻率提高時,在Matlab輔助設計中能夠?qū)崿F(xiàn)的IIR濾波器,在TMS320C54x上卻特性改變,甚至無法使用。對此進行了較深入的仿真分析和研究,找出了在定點DSP上實現(xiàn)IIR數(shù)字濾波器時,濾波器系數(shù)取值的約束條件,摸索了解決問題的基本途徑。
關(guān)鍵詞:IIR數(shù)字濾波器;濾波器系數(shù);Matlab;DSP
中圖分類號:TN713文獻標識碼:B文章編號:1004373X(2008)1911904
Simulation Analysis of IIR Digital Filter in Fixed-point DSP
LIU Shufan,ZHANG Xiaohong,REN Shujie
(College of Science,PLA Science and Engineering University,Nanjing,211101,China)
Abstract:
In practical of realizing IIR digital filter by using fixed-point DSP chip,in regard to the design of filter indicator,when the sampling rate increase,it can realize IIR filter in the Matlab aided design,characteristics on TMS320C54x is actuallychanged,even unable to use to realize.IIR digital filter on fixed-point DSP,the constraint of filter coefficient value and methods to solve the problem is given.
Keywords:IIR digital filter;filter coefficient;Matlab;DSP
1 引 言
IIR數(shù)字濾波器是數(shù)字信號處理中最常用的濾波器類型之一。由于它可用較低的階數(shù)獲得較高的選擇性,所用的存儲單元少,計算量小,同時其設計方法沿用了模擬濾波器設計的成熟技術(shù),因而獲得較為廣泛的應用。同時,在許多實際應用中,模擬信號是實際存在的信號,使用先進的數(shù)字濾波器進行處理時,往往要考慮模擬和數(shù)字信號及系統(tǒng)之間的關(guān)系。因此如何正確進行IIR數(shù)字濾波器的設計是一項很有意義的研究工作。
由奈奎斯特(Nyquist)取樣定理我們知道,只有當模擬濾波器的頻率響應是有限帶寬的,且頻帶限于折疊頻率ΩS/2以內(nèi),等價的數(shù)字濾波器的頻率響應才能避免產(chǎn)生混疊。即有限帶寬與取樣頻率應符合|Ω|≤ΩS/2的關(guān)系。
現(xiàn)在的問題是,根據(jù)給定模擬濾波器的技術(shù)指標進行IIR數(shù)字濾波器的設計時,采樣頻率取值能夠無限制的增大嗎? 是不是只要滿足奈奎斯特取樣定理,設計的IIR數(shù)字濾波器在定點的DSP上就能實現(xiàn)?
在使用定點的DSP芯片TMS320C54x實現(xiàn)IIR數(shù)字濾波器的工程實踐中,我們注意到當取樣頻率提高時,原來在Matlab輔助設計中能夠?qū)崿F(xiàn)的IIR濾波器,在TMS320C54x上卻變得特性改變,甚至無法使用。利用Matlab軟件和DSP調(diào)試軟件CCS,對此進行了較深入的研究。分別觀察兩種軟件仿真的結(jié)果,找出了IIR數(shù)字濾波器在定點DSP實現(xiàn)中出現(xiàn)問題的基本原因,摸索了一些IIR數(shù)字濾波器取樣頻率選擇的規(guī)律。
2 Matlab輔助設計IIR數(shù)字濾波器
在用DSP芯片實現(xiàn)IIR數(shù)字濾波器之前,首先應根據(jù)設計技術(shù)指標用Matlab進行輔助設計。選取取樣頻率時,注意滿足取樣定理,即假定濾波器的有限帶寬最高頻率為fm,則選擇取樣頻率Fs≥2fm(Ωs≥2Ωm)。
在下面的例中選取了兩個不同的取樣頻率,以便進行設計結(jié)果的比較。
以設計一個巴特沃斯IIR數(shù)字低通濾波器為例,要求其通帶截止頻率fp=1 500 Hz,Rp≥3 dB;阻帶邊界頻率fst=2 500 Hz,As≥20 dB。這里假定該濾波器的有限帶寬最高頻率為fm=fst=2 500 Hz。
(1) 選擇濾波器取樣頻率Fs=7 500 Hz,通過Matlab程序可以觀察到該數(shù)字濾波器的幅度頻率特性和零極點分布圖如圖1,并得到設計數(shù)據(jù):
N=3Rp0=1.94As0=20.1
b=
+1.178 903 170 636 770 3e-001
+3.536 709 511 910 310 9e-001
+3.536 709 511 910 310 9e-001
+1.178 903 170 636 770 3e-001
a=+1.000 000 000 000 000 0e+000
-3.943 481 425 156 598 6e-001
+3.745 323 506 447 410 0e-001
-3.706 167 161 966 492 1e-002
圖1 Fs=7 500 Hz時的數(shù)字低通濾波器特性
(2) 將取樣頻率提高到Fs=12 000 Hz,再觀察數(shù)字濾波器的幅度頻率特性和零極點分布圖如圖2所示,并得到設計數(shù)據(jù):
N=4Rp0=2.5As0=20.4
b=
+1.156 831 136 894 841 8e-002
+4.627 324 547 579 367 3e-002
+6.940 986 821 369 050 9e-002
+4.627 324 547 579 367 3e-002
+1.156 831 136 894 841 8e-002
a=+1.000 000 000 000 000 0e+000
-1.891 564 214 634 106 7e+000
+1.638 343 396 929 202 3e+000
-6.720 580 318 738 741 4e-001
+1.103 718 314 819 532 9e-001
由上面兩組圖可見,其取樣頻率均滿足取樣定理,即Fs≥2fm≥5 000 Hz;其極點均分布在單位圓內(nèi),為穩(wěn)定系統(tǒng);實際的通帶衰減和阻帶衰減均在設計指標范圍內(nèi),幅頻響應滿足設計要求。
圖2 Fs=12 000 Hz時的數(shù)字低通濾波器特性
所不同的是,取樣頻率Fs=7 500 Hz時,數(shù)字濾波器的階數(shù)N=3,濾波器的b,a系數(shù)中除a0外均小于1,所有b,a系數(shù)間的絕對值差別較?。蝗宇l率Fs=12 000 Hz時,數(shù)字濾波器的階數(shù)N=4,濾波器b,a系數(shù)中除a0外出現(xiàn)了其他大于1的數(shù)據(jù),所有b,a系數(shù)間的絕對值差別較大。
3 CCS的數(shù)字濾波器仿真測試
在Matlab設計完成后,由Matlab圖形分析,兩組數(shù)據(jù)似乎都是穩(wěn)定的滿足設計要求的設計指標。在DSP硬件系統(tǒng)中能否實現(xiàn),則需進一步進行檢驗。
使用DSP實現(xiàn)IIR數(shù)字濾波器,需先在DSP調(diào)試軟件CCS中用Simulator進行軟件仿真,如果Simulator調(diào)試通過,則再將數(shù)字濾波器程序及數(shù)據(jù)輸入TMS320C54x硬件,進行硬件系統(tǒng)調(diào)試。
進行Simulator軟件仿真時,還需進行以下幾步工作:
(1) 為了測試的需要,先編制程序產(chǎn)生一個測試信號。本例的測試信號是2個正弦信號疊加的周期信號,圖3(a)、(b)分別為輸入時間信號的波形和其對應的頻譜圖。其中一個頻率選自通帶區(qū)域f1=1 000 Hz,另一個頻率選自阻帶區(qū)域f2=2 600 Hz。將該測試信號轉(zhuǎn)換成定點Q0格式后,調(diào)入IIR數(shù)字濾波程序。
(2) 分別將Matlab生成的數(shù)字濾波器的系數(shù)轉(zhuǎn)換成定點Q0格式,調(diào)入數(shù)字濾波的運算程序。
(3) 建立并運行IIR數(shù)字濾波程序,利用CCS的圖形窗觀察濾波結(jié)果。
圖3為Fs=7 500 Hz時的輸入、輸出信號的時域和頻域特性。其中圖3(c)、(d)為濾波后輸出時間信號的波形和其對應的頻譜圖。
從輸出信號的時域波形圖看,基本上是一個正弦信號。從輸出信號的頻譜圖中可見,處于通帶內(nèi)的1 000 Hz頻率能夠順利通過,阻帶內(nèi)2 600 Hz頻率基本被濾除。
圖3 Fs=7 500 Hz時的輸入輸出時域和頻域特性
根據(jù)輸出信號的頻譜圖圖3(d),可以測出通帶與阻帶的相對幅度分別為1 012 270和46 999.6,可以計算出阻帶與通帶之間的相對衰減為:
As0=-20lgU阻帶U通帶=-20lg46 999.61 012 270=26.664 dB
按通帶衰減Rp≥3 dB,阻帶衰減As≥20 dB的設計指標,可知通阻帶相對衰減應該在≥17 dB的范圍內(nèi)。可見,阻帶內(nèi)2 600 Hz頻率處,滿足通阻帶相對衰減的指標要求。
圖4為Fs=12 000 Hz時的輸入輸出時域和頻域特性。其中圖4(a)、(b)為輸入時間信號的波形和其對應的頻譜圖;圖4(c)、(d)為濾波后輸出時間信號的波形和其對應的頻譜圖。
圖4 Fs=12 000 Hz時的輸入輸出時域和頻域特性
從輸出信號的時域波形圖看,仍然是一個合成信號。從輸出信號的頻譜圖中可見,處于通帶內(nèi)的1 000 Hz頻率能夠順利通過;阻帶內(nèi)2 600 Hz頻率存在著不小的幅度,沒有被完全濾除。
根據(jù)輸出信號的頻譜圖圖4(d),可以測出通帶與阻帶的相對幅度分別為143 941和65 236.5,可以計算出阻帶與通帶之間的相對衰減為:
As0=-20lgU阻帶U通帶=-20lg65 236.5143 941=6.87 dB
可見,阻帶內(nèi)2 600 Hz頻率處,不滿足通阻帶相對衰減≥17 dB的指標要求。與在Matlab中設計的結(jié)果As0=20.4 dB有了很大的差別,設計結(jié)果由合格變?yōu)椴缓细瘛?/p>
4 誤差原因的分析與證明
這兩組數(shù)據(jù)設計的結(jié)果,表現(xiàn)是對于相同的設計指標,取樣頻率提高時,不能實現(xiàn)的可能性增大,實質(zhì)問題主要出在哪里呢?
一個數(shù)字濾波器的實現(xiàn)需要進行圖5所示的一系列處理。
因此,進行Simulator軟件仿真時,與實際在DSP上實現(xiàn)基本相同,在下面4個方面將引入誤差:
(1) 對IIR數(shù)字濾波器的系數(shù),采用定點Q0格式進行數(shù)據(jù)的處理,將引入系數(shù)量化誤差;
(2) 對輸入的測試信號,采用定點Q0格式進行數(shù)據(jù)的處理,將引入輸入信號的量化誤差;
(3) 由于使用的A/D,D/A轉(zhuǎn)換器件為有限字長,如D/A轉(zhuǎn)換器的字長僅10 b或12 b(本次仿真輸入信號仍保持16 b)。在軟件程序中事先要把數(shù)據(jù)化為相同的字長位數(shù),再送到D/A轉(zhuǎn)換器,也將引入量化誤差。
(4) 對數(shù)字濾波的運算過程編制相應的程序,其中濾波器選用直接型、級聯(lián)型還是并聯(lián)型,將產(chǎn)生不同的運算量化誤差。
圖5 在DSP中實現(xiàn)一個數(shù)字濾波器
以上幾個方面的運算都會引入一定的量化誤差。但這些問題對取樣頻率高與低,同樣存在。為什么有的設計結(jié)果能實現(xiàn),有的就無法實現(xiàn)?產(chǎn)生誤差的主要原因在哪里?有沒有可能事先初步預見一個數(shù)字濾波器比較難以實現(xiàn)?
(1) 濾波器系數(shù)中出現(xiàn)大于1的系數(shù)時,在定點DSP上較難實現(xiàn)
根據(jù)Matlab運算的數(shù)據(jù)結(jié)果可知,其數(shù)據(jù)處理采用的是單精度或者雙精度的浮點數(shù),而TMS320C54x是一個16 b定點運算芯片。在實現(xiàn)數(shù)字濾波器的過程中,對原來的一組輸入信號或者濾波器系數(shù)的浮點數(shù),必須先轉(zhuǎn)換為定點數(shù)。一般是先進行歸一化,然后再轉(zhuǎn)換成一定格式(如Q0格式)的定點數(shù)。
如對上述選擇取樣頻率Fs=7 500 Hz時的一組濾波器系數(shù), 由于本組數(shù)據(jù)負數(shù)的絕對值沒有大于1,因此可以直接在整組濾波器系數(shù)中選絕對值最大的數(shù)進行處理。將數(shù)據(jù)進行歸一化處理再轉(zhuǎn)換為Q0格式,得到:
bq=3 863, 11 589, 11 589,3 863
aq=32 767,-12 922, 12 273, -1 214
將化為Q0格式后的數(shù)據(jù)與原數(shù)據(jù)相比,求(b-bq./bqmax)./b,可以計算出兩組數(shù)據(jù)之間的相對誤差為:
1.0e-003 *
-0.022 8-0.022 8-0.022 8 -0.022 8
0 -0.030 5-0.056 9 0.329 2
可以看出,當Fs=7 500 Hz時,濾波器系數(shù)數(shù)據(jù)之間的相對誤差均小于10-3,基本不影響在DSP硬件環(huán)境中的實現(xiàn)。
而對上述選擇取樣頻率FS=12 000 Hz時的一組濾波器系數(shù)轉(zhuǎn)換為Q0格式,由于本組數(shù)據(jù)負數(shù)的絕對值有大于1的數(shù)據(jù),因此將正負數(shù)分別進行處理。得到:
bq=379,1 516, 2 274, 1 516,379
aq=32 767,-32 768,32 767,-22 022,3 617
計算化為Q0格式后的數(shù)據(jù)與原數(shù)據(jù)之間的相對誤差為:
0.000 20.000 20.000 20.000 20.000 2
0 0.471 30.389 6-0.000 0-0.000 1
即,當Fs=12 000 Hz時,兩組數(shù)據(jù)之間的相對誤差不均勻,最大相對誤差竟達到47.13%,引起在DSP硬件環(huán)境中不能達到設計要求。這是由于Fs=12 000 Hz時原數(shù)據(jù)組中有兩個絕對值大于1的數(shù)據(jù)。由于a0必須為1,絕對值大于1的數(shù)在處理時將被歸一,因此在這些絕對值大于1數(shù)上產(chǎn)生了較大的誤差。
(2) 改變數(shù)字濾波器的結(jié)構(gòu),是減小誤差的根本途徑
上述低通濾波器采用了IIR數(shù)字濾波器直接Ⅱ型進行仿真,因此效果較差。如果采用級聯(lián)型,情況應有所改善。如圖6所示,采用了級聯(lián)型實現(xiàn)上述IIR低通濾波器后,F(xiàn)s=12 000 Hz時的輸入輸出時域和頻域特性得到了明顯的提高。
圖6 采用級聯(lián)型實現(xiàn)IIR低通濾波器,
Fs=12 000 Hz時的輸入輸出時域和頻域特性
此時設計數(shù)據(jù)為:
SOS matrix:
1211-1.072 066 107 391 336
0.564 141 651 976 095
1211-0.819 498 107 242 771
0.195 645 599 106 783
Scale Factors:
0.123 018 886 146 190
0.094 036 872 966 003
1.000 000 000 000 000
測出的通帶與阻帶的相對幅度分別為161 098和18 911.2,可以計算出阻帶與通帶之間的相對衰減:
As0=-20lgU阻帶U通帶=-20lg18 911.2161 098=18.6 dB
可見,阻帶內(nèi)2 600 Hz頻率處,已經(jīng)滿足通阻帶相對衰減≥17 dB的指標要求,比原來使用直接型結(jié)構(gòu)的6.87 dB有了較大的改善。觀察上述設計數(shù)據(jù)可以看到仍然存在著系數(shù)大于1的情況。不過這個大于1的數(shù)據(jù)為-1.072 066 107 391 336,比較接近1。化為Q0格式后的數(shù)據(jù)與原數(shù)據(jù)之間的相對誤差為6.72%,比采用直接型結(jié)構(gòu)時的47.13%誤差有了明顯的改善。如果采用并聯(lián)型結(jié)構(gòu),該組數(shù)字濾波器的設計指標應該能夠得到更大地改善,更加容易在DSP上實現(xiàn),并獲得更加富余的相對衰減數(shù)值。
5 結(jié) 語
通過尋找多組數(shù)據(jù)反復仿真分析,基本可以得出這樣的結(jié)論:
(1) 在16位的定點DSP中實現(xiàn)IIR數(shù)字濾波器,除a0外,一般不應出現(xiàn)其他系數(shù)絕對值大于1(特別是大于1.1)的情況,否則該濾波器將可能不能實現(xiàn)。絕對值超過1的程度越大,引起的誤差也越大,不能實現(xiàn)的可能性越大。
(2) 在滿足采樣定理的條件下,設計指標相同,采樣頻率選得越大,越有可能出現(xiàn)系數(shù)絕對值大于1的情況,IIR數(shù)字濾波器也越難以實現(xiàn)。
(3) 當直接型不能實現(xiàn)IIR數(shù)字濾波器時,可以試用并聯(lián)型或級聯(lián)型來實現(xiàn)。用并聯(lián)型或級聯(lián)型實現(xiàn)IIR數(shù)字濾波器時,第(1)條同樣須滿足。
參考文獻
[1]陳懷琛.數(shù)字信號處理教程——Matlab釋義與實現(xiàn)[M].北京:電子工業(yè)出版社,2005.
[2]張雄偉,陳亮,徐光輝.DSP集成開發(fā)與應用實例[M].北京:電子工業(yè)出版社,2002.
[3]唐向宏.數(shù)字信號處理——原理、實現(xiàn)與仿真[M].北京:高等教育出版社,2006.
[4]蘇濤,藺麗華.DSP實用技術(shù)[M].西安:電子科技大學出版社,2002.
[5]樓順天,李博菡.基于Matlab的系統(tǒng)分析與設計——信號處理[M].西安:電子科技大學出版社,2000.
注:本文中所涉及到的圖表、注解、公式等內(nèi)容請以PDF格式閱讀原文