風洞自由飛試驗中氣動參數辨識準度評價方法研究

張天姣，汪 清，何開鋒，錢煒祺

（1.空氣動力學國家重點實驗室，四川綿陽 621000；2.中國空氣動力研究與發展中心，四川綿陽 621000）

風洞自由飛試驗中氣動參數辨識準度評價方法研究

張天姣1，2，＊，汪 清1，2，何開鋒1，2，錢煒祺1，2

（1.空氣動力學國家重點實驗室，四川綿陽 621000；2.中國空氣動力研究與發展中心，四川綿陽 621000）

當高超聲速風洞自由飛試驗的測量數據被有色噪聲污染時，傳統的Cramér-Rao界作為參數估計準度的度量往往過于樂觀。文本采用一種修正協方差方法來處理傳統的最大似然估計的殘差，以便計算出有色殘差情況下精確的Cramér-Rao下界，對辨識參數結果進行不確定度評價。以10°半錐角尖錐模型為例，通過大量的Monte Carlo仿真試驗和風洞試驗驗證了修正協方差方法的有效性。結果表明，在風洞試驗測量存在有色噪聲情況下，修正協方差方法給出的標準差均值約為傳統的Cramér-Rao界方法給出的標準差的3～5倍，與參數估計的統計標準差一致，客觀反映了參數辨識結果的精準度。

自由飛；風洞試驗；有色噪聲；氣動參數辨識；準度評價方法

0 引 言

由于風洞自由飛試驗數據采集的特殊性，即測量數據只包含模型位置和姿態角信息，無法獲取線加速度、角速率及角加速度的測量量［1］，導致氣動力系數的辨識準度受位置信息等測量精度影響很大［2］，故待估氣動參數的散布較大，正確評價待估參數的辨識準度十分重要。目前，廣泛采用Cramér-Rao界作為參數估計準度的度量，但通常情況下，Cramér-Rao界過于樂觀［3-4］。文獻［4］的方法通過對殘差進行頻譜分析來對Cramér-Rao界進行修正，這需要分析者工程經驗十足，且工作量較大。文獻［5］首先提出一種方法來處理傳統的最大似然估計殘差，以便計算出有色殘差情況下精確的Cramér-Rao下界，并將其應用于F-18飛行器的縱向機動飛行中［6］。針對風洞自由飛試驗中測量數據被有色噪聲污染的實際情況，本文采用一種新的修正協方差方法替代傳統方法計算Cramér-Rao界，通過對多次Monte Carlo仿真試驗和風洞試驗實測數據進行氣動參數辨識結果的準度評價，驗證了修正協方差方法計算Cramér-Rao界的有效性，其客觀反映了參數辨識結果的精準度。

1 氣動參數辨識準度評價方法

1.1 氣動參數辨識原理

目前工程上應用最為廣泛的氣動力參數辨識方法是最大似然法（MLE：Maximum Likelihood Estimation），該方法將參數辨識問題轉化為一優化問題，通過優化選取氣動力模型參數值，使模型輸出與實測值間的偏差達極小［7］，其形式如下［8］：

式中：ν（i）為輸出誤差向量：

氣動參數辨識問題就是尋求待辨識參數θ的最大似然估計值，使準則函數J（θ）達最小值，即：

本文采用修正Newton-Raphson迭代算法求解此優化問題。其迭代公式為：

式中：

其中，定義靈敏度矩陣S，表示觀測量對待辨識參數的靈敏度，其表達式為：

1.2 Cramér-Rao不等式

Cramér-Rao不等式是準度估計理論中最重要的結果，它給出參數估計準度可達到的理論極限。在某種意義上，Cramér-Rao不等式給出了數據中所含信息量的度量。對于有偏估計，Cramér-Rao不等式為：

式中：!θtrb（θtr）表示估計的偏差對參數真值的一階梯度，M為信息矩陣。

對于無偏估計，b（θtr）＝0，式（8）化為：

式（9）表明，參數估計值與真值的方差總是大于信息矩陣之逆，只有當信息矩陣無窮大時，才能充分接近真值，即當試驗的數據足夠多時，式（9）的等式成立。經驗證明，只要數據含有系統的最大固有周期的幾倍長的數據，式（9）的等號就接近成立［9-11］。因此式（9）是最大似然參數估計準度的最好度量。在參數估計時，θtr是未知的，故以θ＾代替θtr。

對于飛行器動力學系統模型，當過程噪聲分布矩陣為零，觀測噪聲分布矩陣已知時，信息矩陣M的表達式為：

定義信息矩陣M的逆矩陣為散布矩陣D：

那么，矩陣D的對角線元素djj給出第j個參數估計值方差的下界，Cramér-Rao不等式及Cramér-Rao下界表達式為：

Cramér-Rao界是預測準度的方法，這已被氣動參數辨識的大量數字仿真所證實。但是對于飛行試驗實測數據的參數辨識實踐卻表明，重復多次飛行試驗所得參數估計結果的標準差比Cramér-Rao界高5～10倍［12-14］。原因在于，采用式（13）估算Cramér-Rao界時，假定了系統噪聲是平均分布于直至Nyquist頻率的白噪聲，而對氣動參數辨識有重要影響的噪聲卻是局限于低頻區域的有色噪聲。噪聲特性的差異導致的真實標準差與Cramér-Rao界的比值為Nyquist頻率與系統截止頻率之比，約為5～10。因此，采用Cramér-Rao界作為參數估計準度的度量往往偏小，過于樂觀，在工程上參數估計的不確定度通常取Cramér-Rao界的5～10倍。

1.3 修正協方差方法

在有色測量噪聲下，參數估計的協方差為［15］：

當最大似然估計算法收斂時，離散噪聲協方差陣和觀測靈敏度獨立于參數估計，于是上式可改寫為：

對于有色殘差，有：

將式（16）代入式（15），得：

式（18）給出了有色測量噪聲情況下參數估計的協方差陣，其中，和R均采用式（17）進行估計。

參數估計的標準差為：

2 仿真辨識算例

2.1 白噪聲對辨識結果的影響

本節以某10°半錐角尖錐模型為例，具體說明氣動參數辨識算法在高超聲速風洞自由飛試驗中的應用。模型及試驗工況的基本參數為：m＝0.2217kg，l＝0.168 96m，底部直徑D＝0.06m，總壓p0＝514 520Pa，總溫T0＝410K，來流馬赫數Ma＝5，采樣頻率f＝2000Hz。給定模型仿真初始條件：V＝7m／s、θ＝5°、ψ＝5°、γ＝0°和1組氣動參數，在6個觀測量上添加白噪聲，其中，位置坐標、俯仰角、偏航角和滾轉角的噪聲標準差分別為1mm、0.1°、0.1°和0.02°。進行300次Monte Carlo仿真試驗，采用最大似然法對其進行氣動參數辨識，待辨識參數為。

表1列出了Monte Carlo試驗氣動參數估計的統計結果，文中涉及到的角度單位以弧度計。表1中，第2列“真值”為仿真時所用參數值，第3～4列為300次試驗的參數估計值的均值和標準差s，第5～ 6列為傳統的Cramér-Rao界方法給出的標準差的均值及其與s的比值s／，第7～8列為修正協方差方法給出的標準差的均值及其與s的比值s／。

表1 白噪聲對氣動參數辨識結果的影響Table 1 Influence of white noise on aerodynamic parameter estimation results

由表1可見，參數估計的統計標準差s與Cramér-Rao界方法、修正協方差方法給出的標準差均值一致。結果表明，在測量噪聲為白噪聲情況下，2種方法均較好地反映了辨識結果的分布，均可以作為參數估計準度的有效度量。

2.2 有色噪聲對辨識結果的影響

由于對氣動參數辨識有重要影響的噪聲主要是低頻區域的有色噪聲，故在6個觀測量上添加不同水平的有色噪聲。下面簡要介紹以氣動參數辨識為目的添加的有色測量噪聲的生成過程。由于高超聲速風洞自由飛試驗空間較小，試驗時間較短等固有特點，故本文采用的軸對稱尖錐模型的縱向短周期運動為主要運動模態，有色噪聲的主要功率譜應分布于［0，ωn］上，大于Nyquist頻率的區間上功率譜為0，其余部分采用高斯白噪聲描述。ωn為縱向短周期模態近似頻率［17］，表達式如下：

式中：

根據式（21）計算得出ωn約為25Hz。有色噪聲的窄帶部分先由零均值高斯白噪聲通過一個如圖1所示的通帶邊界頻率為25Hz、阻帶截止頻率為50Hz的6階Chebyshev II型低通濾波器生成，再與一個獨立的寬帶零均值高斯白噪聲相疊加。通過對噪聲幅值進行調節，達到所需信噪比。在進行Monte Carlo仿真時，對所有仿真輸出量均采用此種方式獨立地添加噪聲。圖2給出了第300次仿真中疊加到θ上的有色噪聲時間序列和頻譜圖。采用此種方法生成的有色噪聲較好地模擬了真實飛行試驗數據殘差序列的特點，這也是本文采用此有色噪聲模型的原因所在。

圖1 6階Chebyshev II型低通濾波器Fig.1 Sixth-order Chebyshev type II low-pass filter

圖2 有色噪聲時間歷程和頻譜圖（θ，RUN＝300）Fig.2 Time history curve and power spectrum of example simulated colored noise（θ，RUN＝300）

在6個觀測量上添加上文所述形式的有色噪聲，噪聲水平與2.1節白噪聲水平一致，進行300次Monte Carlo仿真試驗，并對其進行參數辨識。圖3示出了和的估計結果分布，可見對于最大似然估計，當試驗次數足夠多時，參數估計值的概率分布趨近于以真值為中心的高斯正態分布。

表2列出了Monte Carlo試驗氣動參數估計的統計結果。由表2可見，參數估計的統計標準差s約為傳統的Cramér-Rao界方法給出的標準差均值的 5倍，表明傳統的Cramér-Rao界方法過于樂觀。而參數估計的統計標準差s與修正協方差方法給出的標準差均值大致相當，表明修正協方差方法較為合理，可以作為參數估計準度的度量。

圖3 Monte Carlo仿真試驗的參數估計值分布圖Fig.3 Aerodynamic parameter estimation results from Monte Carlo simulation

表2 有色噪聲對氣動參數辨識結果的影響Table 2 Influence of colored noise on aerodynamic parameter estimation results

3 實測數據辨識算例

3.1 實測數據算例

上文中仿真算例是將各觀測量的仿真時間歷程作為實測值，即使對仿真值疊加白噪聲或有色噪聲，也并不能完全反映真實試驗飛行歷程，所以本節針對10°半錐角尖錐模型的風洞自由飛試驗實測數據開展氣動參數辨識研究。

試驗是在中國空氣動力研究與發展中心超高速空氣動力研究所Φ1m高超聲速風洞中進行的。共進行5次試驗，試驗模型構型及工況與2.1節近似，各具體參數如表3所示。觀測量與待辨識氣動參數的選取與2.1節一致，主要氣動參數的估計結果及其Cramér-Rao界如表4和圖4所示。圖4中符號“○”表示參數估計值，分別將采用傳統方法和修正協方差方法計算的待估參數的標準差σ和σc作為誤差帶，以紅色（左側）和藍色（右側）實線表示。

表3 各尖錐模型及試驗工況的基本參數Table 3 Basic parameters of pointed cone models and operating conditions

3.2 結果分析

由圖4可見，各參數采用修正協方差方法計算的標準差σc遠大于傳統方法計算結果σ，說明采用傳統方法計算的Cramér-Rao界過于樂觀。圖4中，、和的散布較大，相應地，表4中對應氣動參數的標準差亦很大，這是由力系數對位置坐標的測量準度十分敏感造成的，與文獻［18］中的計算結果一致。實際上，由于力矩導數主要影響短周期運動，力導數影響長周期運動［19］，而風洞自由飛試驗的空間較小，試驗時間較短，長周期運動體現不充分，位置坐標隨迎角、側滑角的變化不明顯，因此，力導數的辨識存在一定難度。只有在模型上安裝過載測量設備的情況下，風洞自由飛試驗才能對力導數進行有效辨識。

圖4 5次風洞試驗主要氣動參數辨識結果Fig.4 Aerodynamic parameter estimation results from 5wind tunnel free-flight tests

表4 5次風洞試驗氣動參數辨識統計結果Table 4 Statistics of aerodynamic parameter estimation results from 5wind tunnel free-flight tests

從表4可見，采用傳統方法計算的Cramér-Rao界過于樂觀，遠小于由5次飛行試驗得到的待估氣動參數的標準差，而采用修正協方差方法計算的Cramér-Rao界是傳統方法計算結果的3～5倍，且與5次飛行試驗得到參數的標準差接近。因此，在測量數據被有色噪聲污染情況下，采用修正協方差方法計算得到的Cramér-Rao界是待估氣動參數準度的一個有效度量。

可注意到，表4中部分氣動參數采用修正協方差方法計算的Cramér-Rao界雖然優于傳統方法計算結果，但仍與參數估計的統計標準差略有差別，這是由于本文中高超聲速風洞自由飛試驗僅進行了5次，樣本量較少，且每次試驗所用模型及試驗工況等各項條件不可能做到完全一致造成的，故可通過增加風洞試驗次數和提高樣本量來修正。

4 結 論

本文將基于最大似然理論的氣動參數辨識方法推廣應用于高超聲速風洞自由飛試驗中。多次Monte Carlo仿真試驗和風洞實測數據的辨識結果表明：

（1）在風洞試驗測量數據被白噪聲污染情況下，傳統方法和修正協方差方法計算的Cramér-Rao界接近，與多次重復試驗的參數估計的統計標準差一致，2種方法均較好地反映了辨識結果的散布程度，均可以作為參數估計準度的有效度量。

（2）當試驗測量數據被有色噪聲污染情況下，傳統方法計算的Cramér-Rao界遠小于參數估計的統計標準差，將其作為參數估計準度的度量過于樂觀，而采用修正協方差方法計算的Cramér-Rao界與參數估計的統計標準差十分接近，將其作為辨識參數結果的準度評價，更能準確地反映試驗中待估氣動參數的散布程度。

（3）采用Chebyshev II型低通濾波器生成的有色噪聲信號能量主要集中于低頻段，而在高頻段能量較低，此種有色噪聲信號較好地仿真了真實的高超聲速風洞試驗測量噪聲，并將其應用于Monte Carlo仿真試驗中。

（4）采用修正協方差方法的Cramér-Rao界的全部計算在時域上進行，避免了對有色殘差進行頻域分析，因此該方法適用范圍廣，對使用者的工程經驗不作要求。

［1］張守言.模型自由飛試驗［M］.北京：國防工業出版社，2002.

［2］Prislin R H.High amplitude dynamic stability characteristics of a blunt 10-degree cones［R］.AIAA-66-465，1966.

［3］Klein V.Determination of stability and control parameters of a light airplane from flight data two estimation methods［R］.NASA TP-1306，1979.

［4］Maine R E，Iliff K W.The theory and practice of estimating the accuracy of dynamic flight-determined coefficients［R］.NASA RP-1077，1981.

［5］Morelli E A，Klein V.Determining the accuracy of maximum likelihood parameter estimates with colored residuals［R］.NASA CR-194893，1994.

［6］Morelli E A，Klein V.Accuracy of aerodynamic model parameters estimated from flight test data［J］.Journal of Guidance，Control and Dynamics，1997，20（1）：74-80.

［7］蔡金獅.飛行器系統辨識學［M］.北京：國防工業出版社，2003.

［8］汪清，錢煒祺，何開鋒.導彈氣動參數辨識與優化輸入設計［J］.宇航學報，2008，29（3）：789-793.Wang Q，Qian W Q，He K F.Aerodynamic parameter identification and optimal input design for missile［J］.Journal of Astrinautics，2008，29（3）：789-793.

［9］Murphy P C.A methodology for airplane parameter estimation and confidence interval determination in nonlinear estimation problems［R］.NASA RP-1077，1981.

［10］Murphy P C，Klein V.Maximum likelihood algorithm using an efficient scheme for computing sensitivities and parameter confidence intervals［R］.AIAA-84-2084，1984.

［11］Murphy P C.Efficient computation of confidence intervals of parameters［R］.AIAA-87-2624，1987.

［12］Iliff K W，Maine R E.Further observations on maximum likelihood estimates of stability and control characteristics obtained from flight data［R］.AIAA-77-1133，1977.

［13］Maine R E，Iliff K W.Use of Cramér-Rao bounds on flight data with colored residuals［J］.Journal of Guidance，Control and Dynamics，1981，4（2）：207-213.

［14］Balakrishnan A V，Maine R E.Improvements in aircraft extraction programs［R］.NASA CR-145090，1975.

［15］Morelli E A，Klein V.Determing the accuracy of aerodynamic model parameters estimated from flight test data［R］.NASA CR-203351，1995.

［16］Bendat J S，Piersol A G.Random data analysis and measurement procedures［M］.2nd ed.New York：Wiley，1986.

［17］方振平，陳萬春，張曙光.航空飛行器飛行動力學［M］.北京：北京航空航天大學出版社，2005.

［18］張天姣，錢煒祺，何開鋒，等.基于最大似然法的風洞自由飛試驗氣動力參數辨識技術研究［J］.實驗流體力學，2015，29（5）：8-14.Zhang T J，Qian W Q，He K F，et al.Research on aerodynamic parameter identification technology in wind tunnel free-flight test based on Maximum Likelihood Estimation［J］.Journal of Experiments in Fluid Mechanics，2015，29（5）：8-14.

［19］顧均曉.飛行穩定性和自動控制［M］.北京：國防工業出版社，2008.

［20］賈區耀，楊益農，陳農.風洞自由飛實驗測量M＝5，10°半錐角尖錐、鈍錐氣動動穩特性［J］.空氣動力學學報，2007，25（2）：24-28.Jia Q Y，Yang Y N，Chen N.Measurement of aerodynamic dynamic stability characters ofM＝5pointed and blunted cones with semi-cone angles equal 10°in wind tunnel free flight experiment［J］.Acta Aerodynamica Sinica，2007，25（2）：24-28.

［21］賈區耀，楊益農.錐動失穩與無控飛行器彈著點精度［J］.氣動物理—理論與應用，2007，2（2）：183-186.Jia Q Y，Yang Y N.Conical lose of stability and unguided aircraft impact accuracy［J］.Physics of Gases：Theory and Applications，2007，2（2）：183-186.

［22］陳農，楊益農，賈區耀.小迎角下鈍錐氣動動態穩定性［C］／／全國第九屆分離流，旋渦和流動控制會議論文集，2002：74-83.

［23］馬家驩，潘文欣，翟曼玲，等.10°尖錐標模高超聲速動導數的實驗測量［J］.空氣動力學學報，1997，15（4）：452-457.Ma J H，Pan W X，Zhai M L，et al.10°Cone model free flight experiment in hypersonic impulse type wind tunnel for dynamic stability measurement［J］.Acta Aerodynamica Sinica，1997，15（4）：452-457.

Research on accuracy assessment method of aerodynamic parameters identified from wind tunnel free-flight test data

Zhang Tianjiao1，2，＊，Wang Qing1，2，He Kaifeng1，2，Qian Weiqi1，2
（1.State Key Laboratory of Aerodynamics，Mianyang Sichuan 621000，China；2.China Aerodynamics Research and Development Center，Mianyang Sichuan 621000，China）

The conventional Cramér-Rao lower bounds method is too optimistic to be a good quantitative assessment of the accuracy of aerodynamic parameters identified from the wind tunnel free-flight test data，considering the colored noise in the measurement data.This paper introduces a technique，that modified covariance matrix method，to process the residuals from a conventional maximum likelihood estimation to compute the accurate Cramér-Rao lower bounds for colored residuals.The modified accuracy assessment method is validated by Monte Carlo simulation and wind tunnel tests of pointed cone models with the semi-cone angle being 10°.The identified results indicate that the Cramér-Rao lower bounds calculated by the modified covariance matrix method are 3～5times the quantity of the conventional.The modified results can be used as an accurate and impersonal assessment of the aerodynamic parameters estimated，which are consistent with the sample standard errors for the estimated parameters for colored residuals.

free-flight；wind tunnel test；colored noise；aerodynamic parameter identification；accuracy assessment method

V212

A

（編輯：楊 娟）

2016-10-08；

2016-11-30

＊通信作者E-mail：ruisazheng＠sina.com

ZhangTJ，WangQ，HeKF，etal.Researchonaccuracyassessmentmethodofaerodynamicparametersidentifiedfromwindtunnel free-flighttestdata.JournalofExperimentsinFluidMechanics，2017，31（1）：39-46.張天姣，汪 清，何開鋒，等.風洞自由飛試驗中氣動參數辨識準度評價方法研究.實驗流體力學，2017，31（1）：39-46.

1672-9897（2017）01-0039-08

10.11729／syltlx20160148

張天姣（1985-），女，天津人，助理研究員。研究方向：飛行器氣動參數辨識。通信地址：四川省綿陽市二環路南段6號13信箱03分箱（621000）。E-mail：ruisazheng＠sina.com