侯建
(北京電子工程總體研究所,北京 100854)
αβ濾波器是一種十分常見且應用廣泛的濾波器,在雷達數(shù)據(jù)處理中大量應用αβ濾波器,因此針對具體的應用需求進行αβ濾波器的參數(shù)優(yōu)化設計是經(jīng)常遇到的問題。本文選取了目前廣泛使用的2個優(yōu)化設計準則,建立了優(yōu)化設計模型,用Matlab編寫了優(yōu)化設計程序。
為了便于后面的論述,在這里先對αβ濾波器的有關(guān)基本概念進行一個簡要的梳理。
αβ濾波器的濾波算法如下[1]:
(1)
(2)
(3)

根據(jù)αβ濾波器的濾波方程,可以推出其傳遞函數(shù)為[1]
(4)
(5)
式中:Hf(z),Hp(z)分別為αβ濾波器的濾波估計和預測估計的傳遞函數(shù)。
可以求出αβ濾波器濾波估計的穩(wěn)態(tài)方差壓縮系數(shù)[1]:
(6)

可以推出αβ濾波器的等效噪聲帶寬為[2]
(7)
由傳遞函數(shù)可以看出αβ濾波器是一個二階離散系統(tǒng)。
對于一個穩(wěn)定的二階離散系統(tǒng),如果其階躍響應為衰減振蕩或等幅振蕩(具有一對共軛極點或一個二重極點),則傳遞函數(shù)可以寫成:
(8)

系統(tǒng)傳遞函數(shù)的極點為
z1=e-ξω0Tejωd T,
(9)
z2=e-ξω0Te-jωd T.
(10)
對比式(4),(5)和(8),可以得到:
1-α=e-2ξω0T,
(11)
2-α-β=2e-ξω0TcosωdT,
(12)
(13)
根據(jù)式(11)和(12)可以得到:
(14)
(15)
通常濾波周期T是由系統(tǒng)其他諸多因素確定的,因此這里僅討論在濾波周期已經(jīng)確定的情況下參數(shù)α和β的優(yōu)化設計問題。
從一般的意義上講,如果有n個參數(shù)需要確定,就必須有n個彼此獨立的方程,求解這n個方程構(gòu)成的方程組就可以得到n個參數(shù)。因此αβ濾波器的參數(shù)設計問題就轉(zhuǎn)化為如何確定2個相互獨立的方程。
在前面的論述中已經(jīng)給出了有關(guān)αβ濾波器的一些公式,如:濾波估計的穩(wěn)態(tài)方差壓縮系數(shù)公式(6)、濾波器的等效噪聲帶寬公式(7)、系統(tǒng)的阻尼系數(shù)和阻尼振蕩頻率公式(13)~(15)等。現(xiàn)在討論如何利用這些公式并結(jié)合一些準則來確定參數(shù)α和β。
αβ濾波器一般是用來對被測量進行跟蹤濾波的,因此跟蹤濾波精度是一個十分重要的考量因素,而穩(wěn)態(tài)方差則是衡量穩(wěn)態(tài)跟蹤精度的一個重要指標,穩(wěn)態(tài)方差壓縮系數(shù)表征了穩(wěn)態(tài)情況下濾波器輸出的方差相對輸入方差的比例系數(shù),該數(shù)值越小,則輸出方差就越小,因此通常情況下可以將式(6)作為一個方程;當然在某些情況下如果系統(tǒng)對等效噪聲帶寬有明確的要求也可以將式(7)作為一個方程。實際上式(6)和式(7)在很大程度上是等效的,因為如果帶寬越寬,則系統(tǒng)對起伏噪聲的濾波效果就越差,因此輸出方差就越大,穩(wěn)態(tài)方差壓縮系數(shù)就越大。
除了穩(wěn)態(tài)跟蹤精度之外,通常也關(guān)心系統(tǒng)的暫態(tài)特性,如一般希望系統(tǒng)的過渡過程盡可能短,超調(diào)盡可能小。前文已說明αβ濾波器是一個二階系統(tǒng),對于一個二階系統(tǒng),系統(tǒng)的阻尼系數(shù)ξ和自然諧振頻率ω0是與系統(tǒng)暫態(tài)特性密切相關(guān)的參數(shù);圖1給出了阻尼系數(shù)ξ取不同數(shù)值時二階系統(tǒng)的階躍響應。
式(16~18)給出了0<ξ<1時二階系統(tǒng)階躍響應暫態(tài)特性的計算公式[3]:
(16)
(17)
(18)
式中:Δ為二階系統(tǒng)階躍響應的超調(diào);Tr為二階系統(tǒng)階躍響應的上升時間;Ts為二階系統(tǒng)階躍響應的調(diào)整時間;ε為二階系統(tǒng)階躍響應的調(diào)整時間對應的歸一化誤差(階躍響應在Ts時刻的取值與穩(wěn)態(tài)值之差相對穩(wěn)態(tài)值歸一化)。
從圖1和式(16)~(18)可以看出:
(1) 當阻尼系數(shù)ξ=0時,系統(tǒng)的階躍響應是一個等幅振蕩,系統(tǒng)處于臨界穩(wěn)定狀態(tài),不能正常工作,因此要求ξ>0;
(2) 當0<ξ<1時,系統(tǒng)的階躍響應是一個衰減振蕩,阻尼系數(shù)ξ越大,系統(tǒng)的阻尼越大,衰減越快,超調(diào)也越小;當ξ一定時,要使上升時間縮短,則需要提高系統(tǒng)的自然諧振頻率ω0;
(3) 當ξ≥1時,系統(tǒng)的階躍響應是一個單調(diào)衰減(沒有振蕩),此時系統(tǒng)沒有超調(diào),阻尼系數(shù)ξ越大,階躍響應的上升越慢,調(diào)整時間也越長;雖然ξ≥1時,系統(tǒng)的階躍響應沒有振蕩(也沒有超調(diào)),但此時系統(tǒng)的階躍響應相對穩(wěn)態(tài)值的收斂速度要比0<ξ<1時慢。
通過上述分析可以看出如果從系統(tǒng)的暫態(tài)特性考慮,一般取0.7≤ξ≤1比較合適。因此確定參數(shù)α和β的第2個公式可以按下述原則確定:
(1) 取ξ=1,即系統(tǒng)設計工作在臨界阻尼狀態(tài);
(2) 在0.7≤ξ≤1范圍內(nèi)選取一個合適的阻尼系數(shù)。
下面分別討論上述2種情況。
當取ξ=1時,此時系統(tǒng)的特征方程有一對重根,這時可以推出參數(shù)α和β滿足[1]:
(19)
當在0.7≤ξ≤1范圍內(nèi)選取一個合適的阻尼系數(shù)ξ時,可以推出參數(shù)α和β滿足:
(20)
顯然,式(19)是式(20)在ξ=1時的特例。
此外,文獻[1]還給出了一種優(yōu)化設計情況,即在“給定速度估值暫態(tài)平方差條件下使速度估值穩(wěn)態(tài)方差壓縮系數(shù)最小”。根據(jù)該準則可以推出參數(shù)α和β滿足:
(21)
至此給出了確定參數(shù)α和β可以采用的公式。為便于描述,稱式(19)~(21)為優(yōu)化準則公式,稱式(6),(7)為設計指標公式。下面就可以利用所選用的公式,根據(jù)具體的設計要求(如要求的“穩(wěn)態(tài)方差壓縮系數(shù)/等效噪聲帶寬”、“阻尼系數(shù)——可以根據(jù)希望的系統(tǒng)暫態(tài)響應指標確定”)求解方程得到參數(shù)α和β。具體方法如下:

圖1 二階系統(tǒng)階躍響應Fig.1 Step Response of second-order system
首先從式(6)或式(7)中選擇一個作為設計指標公式并確定希望的指標設計值(如選擇式(6)希望的指標設計值就是方差壓縮系數(shù),如選擇式(7)希望的指標設計值就是等效噪聲帶寬);然后從式(19)~(21)中選擇一個作為優(yōu)化準則公式(當選擇式(20)時,需要同時確定一個合適的阻尼系數(shù)值ξ),這樣就構(gòu)成了一個由2個方程組成的二元方程組,求解這個方程組就得到了所要求的參數(shù)α和β。
需要說明的是:如果選擇了式(6)或式(7)再加上式(21)則可以得到一個關(guān)于α的二次方程,可以直接求解。但如果選擇了式(6)或式(7)再加上式(19)或式(20),則得到的是一個關(guān)于α的超越方程,只能采用數(shù)值方法求解。在進行數(shù)值解時需要首先確定一個初值,然后利用程序進行迭代求解。在這種情況下如果初值選取的不合適則有可能導致迭代運算不能收斂或收斂到錯誤的結(jié)果。因此,如何選定初值關(guān)系到迭代運算是否能夠正確地收斂。下面給出確定初值的方法。
首先需要確定α和β的取值范圍。

(4-2α-β)s2+2αs+β=0.
(22)
則根據(jù)系統(tǒng)穩(wěn)定性判定準則,必須滿足以下條件:

(23)


圖2 αβ濾波器的參數(shù)取值范圍Fig.2 Effective parameter range of αβ filter
有了α和β的取值范圍和已經(jīng)選定的公式所確定的α和β的約束關(guān)系就可以按照下述方法確定α的計算初值,并進一步經(jīng)迭代運算得到最終設計值:
1) 在[0,1]范圍內(nèi)按照一定的間隔步長選取一系列的α值;
2) 對于每一個α,利用已經(jīng)選定的公式所確定的α和β的約束關(guān)系計算對應的β;
3) 根據(jù)前面兩步獲得的每一對α和β,利用公式(6),(7),(13)~(15)計算αβ濾波器的相關(guān)性能指標;
4) 取最接近所要求指標所對應的α值作為計算初值;
5) 利用所確定α計算初值啟動迭代運算,得到最終的α設計值;
6) 利用已經(jīng)選定的公式所確定的α和β的約束關(guān)系計算對應的β得到β的設計值。
圖3和圖4了給出了αβ濾波器的參數(shù)計算機優(yōu)化設計流程。關(guān)于參數(shù)α的迭代算法,可以根據(jù)情況選用各種合適的數(shù)值解法,本文選擇了最簡單的牛頓切線法。

圖3 αβ濾波器的參數(shù)計算機優(yōu)化設計流程Fig.3 Optimized design flat for alpha-beita filter parameter
根據(jù)本文所介紹的方法,我們按照圖3和圖4的流程用Matlab編寫了αβ濾波器的參數(shù)計算機優(yōu)化設計程序,這里給出2個設計實例。
實例1:假定要設計的αβ濾波器的采樣周期為T=36 ms,要求的穩(wěn)態(tài)方差壓縮系數(shù)為0.36,采用臨界阻尼準則。
設計結(jié)果為α=0.485 945,β=0.080 106,將參數(shù)α和β的設計結(jié)果代入式(6),(15)和(7)分別得到:
穩(wěn)態(tài)方差壓縮系數(shù):
阻尼系數(shù):
等效噪聲帶寬:
可見設計結(jié)果滿足設計要求(穩(wěn)態(tài)方差壓縮系數(shù)為0.36,阻尼系數(shù)為1——臨界阻尼)。
圖5~8給出了一些驗證曲線。
圖5為按照臨界阻尼準則,當參數(shù)α在[0,1]區(qū)間變化時,參數(shù)β的變化規(guī)律曲線,圖中的“*”點為本次設計結(jié)果(對應穩(wěn)態(tài)方差壓縮系數(shù)為0.36且臨界阻尼的參數(shù)α和β);
圖6為按照臨界阻尼準則,當參數(shù)α在[0,1]區(qū)間變化時,穩(wěn)態(tài)方差壓縮系數(shù)為的變化規(guī)律曲線,圖中的“*”點為本次設計結(jié)果(對應穩(wěn)態(tài)方差壓縮系數(shù)為0.36且臨界阻尼的參數(shù)α和β);

圖4 αβ濾波器的參數(shù)計算機優(yōu)化設計流程(續(xù))Fig.4 Optimized design flat for alpha-beita filter parameter(continue)
圖7為按照臨界阻尼準則,當參數(shù)α在[0,1]區(qū)間變化時,阻尼系數(shù)的變化規(guī)律曲線(曲線為常數(shù)1表明滿足臨界阻尼準則);

圖5 按照臨界阻尼準則,當α在[0,1]變化時, β的變化曲線Fig.5 Value distribution curve of parameter beita according to the rule of critical damping while the parameter alpha varies in range of [0,1]

圖6 按照臨界阻尼準則,當α在[0,1]變化時,穩(wěn)態(tài)方差壓縮系數(shù)的變化曲線Fig.6 Value distribution curve of variance compress coefficient according to the rule of critical damping while the parameter alpha varies in range of [0,1]
圖8為求解參數(shù)α時的超越方程f(x)=0對應的函數(shù)y=f(x)曲線。可見:由于y=f(x)曲線的一階導數(shù)符號不是恒定的,因此在用牛頓切線法求解方程時對迭代運算的初值必須正確選取,否則有可能得不到正確的收斂結(jié)果。

圖7 按照臨界阻尼準則,當α在[0,1]變化時, 阻尼系數(shù)的變化曲線Fig.7 Value distribution curve of damping coefficient according to the rule of critical damping while the parameter alpha varies in range of [0,1]

圖8 求解參數(shù)α時的超越方程f(x)=0的 函數(shù)y=f(x)曲線Fig.8 Value curve of function y=f(x) where f(x) is the transcendental equation f(x)=0 to solve to determine the parameter alpha
圖9為所設計的αβ濾波器對給定輸入的濾波輸出。圖中曲線“in”,“p”,“f”分別表示輸入、預測輸出、濾波輸出。輸入信號是一個附加了噪聲的勻加速變化的信號。從圖中可以看出所設計的αβ濾波器的預測和濾波輸出能夠跟隨輸入變化,并且從濾波器的輸入噪聲和輸出起伏誤差的比較可以看出預測和濾波輸出的起伏要比輸入有一定程度的抑制,經(jīng)過對輸入和輸出的統(tǒng)計,輸入中加入的噪聲方差為1.105 547,預報輸出的起伏方差為0.604 702,濾波輸出的起伏方差為0.448 041;濾波輸出的方差壓縮倍數(shù)為0.405 266(設計值為0.36)。
實例2:假定要設計的αβ濾波器的采樣周期為T=36 ms,要求的穩(wěn)態(tài)方差壓縮系數(shù)為0.36,采用給定阻尼準則,要求阻尼系數(shù)為0.78。

圖9 所設計濾波器的輸入和輸出Fig.9 Input and output of the filter designed
設計結(jié)果為α=0.457 137,β=0.113 303,將參數(shù)α和β的設計結(jié)果代入式(6),(15)和(7)分別得到:
穩(wěn)態(tài)方差壓縮系數(shù):
阻尼系數(shù):
等效噪聲帶寬:
可見設計結(jié)果滿足設計要求(穩(wěn)態(tài)方差壓縮系數(shù)為0.36,阻尼系數(shù)為0.78)。
圖10~13給出了一些驗證曲線。
圖10為按照給定阻尼系數(shù)為0.78作為設計準則,當參數(shù)α在[0,1]區(qū)間變化時,參數(shù)β的變化規(guī)律曲線,圖中的“*”點為本次設計結(jié)果(對應穩(wěn)態(tài)方差壓縮系數(shù)為0.36且阻尼系數(shù)為0.78的參數(shù)α和β);
圖11為按照給定阻尼系數(shù)為0.78作為設計準則,當參數(shù)α在[0,1]區(qū)間變化時,穩(wěn)態(tài)方差壓縮系數(shù)為的變化規(guī)律曲線,圖中的“*”點為本次設計結(jié)果(對應穩(wěn)態(tài)方差壓縮系數(shù)為0.36且阻尼系數(shù)為0.78的參數(shù)α和β);
圖12為按照給定阻尼系數(shù)為0.78作為設計準則,當參數(shù)α在[0,1]區(qū)間變化時,阻尼系數(shù)的變化規(guī)律曲線(曲線為常數(shù)0.78表明滿足給定阻尼準則);
圖13為求解參數(shù)α時的超越方程f(x)=0對應的函數(shù)y=f(x)曲線??梢姡河捎趛=f(x)曲線的一階導數(shù)符號不是恒定的,因此在用牛頓切線法求解方程時對迭代運算的初值必須正確選取,否則有可能得不到正確的收斂結(jié)果。

圖10 按照給定阻尼0.78準則,當α在[0,1]變化時, β的變化曲線Fig.10 Value distribution curve of parameter beita according to the rule of given damping with value of 0.78 while the parameter alpha varies in range of [0,1]

圖11 按照給定阻尼0.78準則,α在[0,1]區(qū)間變化時,穩(wěn)態(tài)方差壓縮系數(shù)的變化曲線Fig.11 Value distribution curve of variance compress coefficient according to the rule of given damping with value of 0.78 while the parameter alpha varies in range of [0,1]

圖12 按照給定阻尼0.78準則,當α在[0,1] 區(qū)間變化時,阻尼系數(shù)的變化曲線Fig.12 Value distribution curve of damping coefficient according to the rule of given damping with value of 0.78 while the parameter alpha varies in range of [0,1]

圖13 求解參數(shù)α時的超越方程f(x)=0的函數(shù)y=f(x)曲線Fig.13 Value curve of function y=f(x) where f(x) is the transcendental equation f(x)=0 to solve to determine the parameter alpha
本文給出了3個優(yōu)化設計準則:
準則1:臨界阻尼(ξ=1);
準則2:給定阻尼(0<ξ<1);
準則3:給定速度估值暫態(tài)平方差條件下使速度估值穩(wěn)態(tài)方差壓縮系數(shù)最小。
準則1和準則2都已經(jīng)過實例驗證,得到了預期效果,表明了設計模型的正確性。由于采用準則3時得到一個關(guān)于α的二次方程,可以直接求解,比較簡單,所以本文沒有給出設計算例。
利用本文所建立的模型可以實現(xiàn)對αβ濾波器的參數(shù)優(yōu)化設計。