趙艷輝，張擁軍，李海峰，趙 霞

(中國空空導彈研究院，河南 洛陽 471009)

0 引 言

受氣動力和慣性力的作用，導彈在飛行過程中會產生彈性形變和模態(tài)振動。彈體的彈性振動信號通過彈上傳感器進入控制系統(tǒng)，導致控制面的附加輸出，使彈性振動加劇，嚴重時使控制回路失去穩(wěn)定性。因此,在設計飛行控制系統(tǒng)時,需對彈性振動進行抑制，使控制系統(tǒng)的頻率特性在彈體模態(tài)頻率處獲得適當的增益衰減，保證控制系統(tǒng)的伺服彈性穩(wěn)定性。文獻[1-2]指出，高超聲速飛行器在氣動/結構/推進系統(tǒng)間存在嚴重的耦合效應，彈性振動是控制系統(tǒng)必須解決的問題之一。文獻[3-5]對彈箭和飛機控制系統(tǒng)設計中彈性振動抑制問題的解決方案進行了全面的綜述,并指出合理設計陷波器是實現振動模態(tài)抑制的基本方法之一。文獻[6] 利用遞歸最小二乘算法在線估計彈體一階彎曲頻率，并以估計值作為結構濾波器的中心頻率，實現了對彈體彎曲模態(tài)振動的抑制。文獻[7]以氣動伺服彈性系統(tǒng)的彈性模態(tài)頻響峰值最小作為優(yōu)化目標，以剛體模態(tài)頻響特性作為設計約束，建立了一種基于多目標遺傳算法的結構陷幅濾波器優(yōu)化設計方法。文獻[8]提出一種基于修正反正切函數的陷波器自適應算法，克服了梯度自適應算法收斂速度慢、對噪聲敏感等缺點。文獻[9]提出一種并行主動噪聲濾波器自適應算法，可以獲得多個頻率的精確估計值，有效實現了主動噪聲控制。文獻[10] 通過粒子群多參數尋優(yōu)，以時域峭度最大原則對陷波器中心頻率及帶寬進行自適應選取，提高了陷波效率及準確性。文獻[11] 提出一種基于小波變換與奇異值分解的彈性自適應陷波濾波方法，在線估計彈性振型頻率，實現自適應陷波。上述文獻提供了多種結構濾波器的設計方法，部分方法在連續(xù)域設計結構濾波器，并分析其性能；在線辨識彈性模態(tài)頻率的方法對振蕩信號的信噪比有特定要求，限制了其在工程實際中的應用。

在工程實際中，結構濾波器的離散化過程是實現彈性模態(tài)振動抑制的關鍵環(huán)節(jié)，需要在嵌入飛控計算機的飛控軟件中實現與連續(xù)域結構濾波器性能一致的數字化濾波器。目前，結構濾波器的數字化基于控制系統(tǒng)采樣頻率遠高于濾波中心頻率這一假設展開。在飛控計算機采樣頻率一定的情況下，這種假設有時不能成立，甚至出現陷波中心頻率接近Nyquist頻率(0.5倍采樣頻率)的情況。采樣定理和數字控制系統(tǒng)設計理論表明：隨著濾波中心頻率接近Nyquist頻率，離散化后的數字濾波器會發(fā)生頻率特性畸變。從公開發(fā)表的文獻看，中心陷波頻率接近Nyquist頻率情況下的結構濾波器設計和離散化問題，目前仍是一項技術空白，如何解決這一問題是本文研究的重點。

1 常用結構濾波器離散化方法及其缺點

通常，數字化結構濾波器的設計步驟是：按照陷波中心頻率處的增益衰減要求，在連續(xù)域設計結構濾波器；通過合適的離散化方法，將連續(xù)域的濾波器傳遞函數離散化，并應用于數字化飛行控制系統(tǒng)。

復變量z和復變量s之間的轉換關系為

z=eTss

(1)

式中：Ts為控制系統(tǒng)的采樣頻率。

常用的離散化方法有雙線性變換法、預畸校正法和零極點匹配法[12]。

1.1 雙線性變換法(Tustin Method)

利用式(1)，在s=0的領域內進行1階Taylor展開，得到z的近似表達式：

(2)

利用式(2)反解連續(xù)域復變量s，得到雙線性變換公式：

(3)

雙線性變換方法簡單，應用廣泛，離散化后的濾波器一定是穩(wěn)定的。但由于連續(xù)域的頻率ωs和離散域的頻率ωd的非線性關系，造成ωd的扭曲現象：

在離散域看到的頻率ωd存在如下關系：

(4)

通過雙線性變換，將ωs映射到離散域后，得到的ωd小于ωs。隨著ωs接近于Nyquist頻率(0.5倍采樣頻率)，扭曲現象愈發(fā)嚴重。ωs和ωd之間的非線性關系如圖1所示。

圖1 ωs和ωd之間的非線性關系Fig.1 Nonlinear relation between ωs and ωd

這種頻率扭曲現象隨著頻率逐漸接近Nyquist頻率越發(fā)嚴重，甚至會使結構濾波器失效。

1.2 預畸校正法(Prewaping Method)

為改善雙線性變換離散化方法的缺陷，需要對結構濾波器的陷波中心頻率進行預畸校正。根據離散域的中心頻率的期望值，對連續(xù)域濾波器的中心頻率進行預畸校正，使ωd與ωs一致，預畸校正過程如下：

簡化得

(5)

通過預畸校正可以使離散化的濾波器和連續(xù)域濾波器的中心頻率和衰減深度一致。但隨著濾波中心頻率接近Nyquist頻率，結構濾波器在中心頻率附近的陷波寬度逐漸減小。這是由于預畸校正方法僅能對某一個特定頻率進行預畸校正所造成的。

1.3 零極點匹配法(Matched Pole-Zero Method)

利用式(1)可以將連續(xù)域濾波器的零極點映射為離散域濾波器的零極點。按照如下步驟進行連續(xù)域濾波器的零極點匹配離散化過程：

(1) 將連續(xù)域傳遞函數Gc(s)的有限極點在z域中用z=eTss表示。

(2) 將連續(xù)域傳遞函數Gc(s)的有限零點在z域中用z=eTss表示。

(6)

(7)

式中：m≤n；zi=-e-μiTs，i=1,…,m；pk=-e-λkTs，k=1,…,n。

(4) 調整Gd(z)的增益Kz。

對于結構濾波器設計，應滿足約束條件：

Gd(z)|z=1=Gc(s)|s=0

(8)

零極點匹配方法避免了預畸校正方法帶來的陷波寬度變窄的問題，但是隨著濾波中心頻率接近Nyquist頻率，濾波器的相位滯后和幅值衰減程度遠超過連續(xù)域傳遞函數的設計結果，需要通過合理的設計方法加以克服。

2 數字化結構濾波器設計

2.1 數字化結構濾波器解析表達式

結構濾波器連續(xù)域傳遞函數：

(9)

雙線性變換后得到的數字濾波器：

(10)

預畸校正后，再經雙線性變換得到數字濾波器：

(11)

零極點匹配離散化后得到數字濾波器：

(12)

利用歐拉公式得到化簡后的離散化傳遞函數：

(13)

考慮到控制系統(tǒng)的有限帶寬特性，利用低頻匹配關系式(8)確定結構濾波器歸一化增益Kz：

(14)

式(10)～(11)和式(13)為不同離散化方法對應的數字化結構濾波器解析表達式，可應用于在線實時更新濾波器參數的情形，也可以作為眾多結構濾波器設計方法離散化過程的有益補充。

2.2 近Nyquist頻率情況下的濾波性能改進方法

當結構濾波器中心頻率接近Nyquist頻率時，現有離散化方法存在缺點，本文給出2種通過修正連續(xù)域結構濾波器分母多項式阻尼比ξ2至一個期望值ξ2c來改進零極點匹配離散化設計結果的優(yōu)化方法。

2.2.1 方法1：衰減增益約束條件下的修正方法

假設修正后的濾波器離散化結果如下：

(15)

若Fnf(z)在指定的關鍵頻率ωc處的增益和式(9)所示的連續(xù)域濾波器的增益相同，即滿足約束關系式：

|Fnf(z)|z=ejωcTs|=|Fnf(s)|s=jωc|

(16)

式中：ωc=2π(fm-Δf)，fm為濾波器中心頻率，Δf為彈性模態(tài)頻率的最大攝動量。

式(16)為超越方程，難以獲得解析解，可以通過數值解法獲得ξ2c的數值解。式(16)右端的模為常值，左端的模為隨ξ2c變化的函數，二者交點對應的橫坐標即為ξ2c的數值解。

2.2.2 方法2：相位滯后約束條件下的修正方法

對于式(15)，若Fnf(z)在指定的角頻率ωl處的相位和式(9)所示的連續(xù)域濾波器的相位相同，即滿足約束關系式：

∠Fnf(z)|z=ejωlTs=∠Fnf(s)|s=jωl

(17)

式中：ωl=2πfl，fl為控制系統(tǒng)設計者所關心的與控制系統(tǒng)回路帶寬有關的低頻段頻率值。

式(17)為超越方程，難以獲得解析解，可以通過數值解法獲得ξ2c的數值解。式(17)右端的相位為常值，左端的相位為隨ξ2c變化的函數，二者交點對應的橫坐標即為ξ2c的數值解。

3 設計實例與仿真驗證

通過設計實例演示現有離散化方法的缺點，并證明近Nyquist頻率情況下濾波性能改進方法的有效性。設控制系統(tǒng)的離散化周期Ts為0.004 s，分別在陷波中心頻率遠小于及接近Nyquist頻率的情況下進行仿真驗證。

陷波器中心頻率40 Hz及110 Hz情況下3種離散化方法的頻率特性比較結果如圖2所示。

圖2 陷波器中心頻率40 Hz,110 Hz情況下不同離散化方法的頻率特性比較結果Fig.2 Frequency characteristic comparison of different discretization methods when center frequency is 40 Hz and 110 Hz

從圖2 (a) 可知,濾波器中心頻率遠小于Nyquist頻率的情況下，雙線性變換形成的離散化濾波器的中心頻率為37 Hz,偏離了設計值40 Hz；通過預畸校正可以使離散化后的濾波器中心頻率與連續(xù)域設計結果一致；零極點匹配離散化結果的幅頻特性更接近連續(xù)域的設計結果。

從圖2 (b) 可知,濾波器中心頻率接近Nyquist頻率的情況下，通過預畸校正仍可以使離散化后的濾波器中心頻率與連續(xù)域設計結果一致，但濾波器的陷波寬度變窄；彈體彈性模態(tài)頻率的攝動超過一定邊界的情況下，經預畸校正獲得的離散化濾波器將失去彈性濾波作用；采用零極點匹配方法獲得的數字化結構濾波器最為有效，但需要克服低頻相位滯后偏大的缺點。

通過方法1得到修正阻尼比ξ2c的數值解如圖3所示，數字化結構濾波器的性能如圖4所示。

圖3 采用方法1獲得的ξ2c的數值解Fig.3 Numerical solution of ξ2c obtained by method 1

圖4 采用方法1獲得的數字化結構濾波器設計結果Fig.4 Design results of digital notch filter obtained by method 1

從圖3～4可以看出，方法1通過將分母多項式的阻尼比ξ2c的數值修正為0.4，可以使得近Nyquist頻率情況下的離散域濾波器的陷波特性更接近連續(xù)域設計結果，同時相位滯后也相應減小。

通過方法2得到修正阻尼比ξ2c的數值解如圖5所示，數字化結構濾波器的性能如圖6所示。

圖6 采用方法2獲得的數字化結構濾波器設計結果Fig.6 Design results of digital notch filter obtained by method 2

從圖5～6可以看出,方法2通過將分母多項式的阻尼比ξ2c的數值修正為0.35，可以使得近Nyquist頻率情況下的離散域濾波器低頻相位特性幾乎與連續(xù)域設計結果重合，陷波特性也更接近連續(xù)域設計結果。

圖5 采用方法2獲得的ξ2c的數值解Fig.5 Numerical solution of ξ2c obtained by method 2

假設濾波器輸入信號為幅值1°、頻率110 Hz的正弦信號，對上述濾波器的陷波效果進行檢驗。彈體彈性模態(tài)頻率為名義值時的仿真結果如圖7所示，彈體彈性模態(tài)頻率攝動情況下的仿真結果如圖8～9所示。

圖7 彈性模態(tài)頻率取名義值110Hz時的仿真結果Fig.7 Simulationresultswhenelasticmodalfrequencyis110Hz圖8 彈性模態(tài)頻率取邊界值105Hz時的仿真結果Fig.8 Simulationresultswhenelasticmodalfrequencyis105Hz圖9 彈性模態(tài)頻率取邊界值115Hz時的仿真結果Fig.9 Simulationresultswhenelasticmodalfrequencyis115Hz

從圖8～9可以看出,在彈體彈性模態(tài)頻率攝動的情況下，采用預畸校正離散化方法形成的數字化結構濾波器已經不能滿足彈性濾波需求；而本文給出的2種修正零極點匹配方法的濾波性能依舊能夠和連續(xù)域設計結果保持一致，證明了本文提出的近Nyquist頻率情況下修正零極點匹配離散方法的有效性。

4 結 論

本文回顧了常用的結構濾波器離散化方法，給出了3種數字化結構濾波器的解析表達式，并分析了陷波中心頻率接近Nyquist頻率的情況下存在的缺點。通過分析離散域濾波器的頻率特性，得出了陷波中心頻率接近Nyquist頻率的情況下適合采用零極點匹配方法離散化方法的結論，并給出了2種修正分母多項式阻尼比的數值方法，使修正零極點匹配方法獲得的數字化結構濾波器性能和連續(xù)域設計結果一致；最后，通過輸入信號頻率存在攝動的濾波仿真試驗驗證了本文提出的近Nyquist頻率情況下修正零極點匹配離散化方法的有效性。