多方法融合的汽車質心側偏角估計

高自群，謝桂芝，周兵，許艷，吳曉建，柴天

(1.湖南大學 汽車車身先進設計制造國家重點實驗室，湖南 長沙 410082；2.湖南大學 國家高效磨削工程技術研究中心，湖南 長沙 410082；3.南昌大學 先進制造學院，江西 南昌 330031)

隨著人們對汽車安全需求的不斷提升，越來越多的主動安全控制系統被應用于汽車[1-3]，質心側偏角是穩定性控制系統的關鍵參數，能夠及時準確地獲取該參數對于汽車安全駕駛的意義重大.由于成本原因，商用車上很少通過如光學傳感器、GNSS/INS 組合導航等檢測裝置直接獲取質心側偏角，而是采用軟測量技術及狀態觀測的方式獲取質心側偏角[4].質心側偏角觀測器的構建方式主要分為3 種：1)基于運動學、2)基于動力學和3)基于運動幾何.方式1）將傳感器量測量(橫擺角速度、側向加速度、縱向速度)與質心側偏角速度的微分關系式進行積分或者利用卡爾曼濾波算法得到質心側偏角[5].方式1)對車輛參數的不確定性、路面附著條件以及駕駛操縱輸入頻率的變化具有魯棒性，但對傳感器的測量精度要求很高，對環境噪聲、傳感器偏置以及標定誤差比較敏感[6].反饋補償可以消除該方式中由傳感器漂移、標定誤差引起的量測信號偏差[7].方式2)是基于車輛動力學模型和輪胎模型進行觀測器的構建，弊端在于估計結果的優劣依賴于觀測模型以及模型參數(車輛參數、路面附著系數、過程噪聲與量測噪聲協方差矩陣等)的準確性.為了簡化估計模型，建模時將如側向風引起的橫擺力矩與側向力的因素忽略，由此導致估計結果的偏差[8-9].Zhang 等[10]提出自適應擴展卡爾曼濾波器，并將其用于克服過程噪聲與環境噪聲統計特性的不確定性對于估計精度的不利影響.為了更好地逼近非線性輪胎模型，提高觀測器構建精度，Zhang 等[11]提出切換T-S 模糊觀測器，估計車輛質心側偏角.上述觀測方法在一定程度上改善了方式2)的估計效果，但該方式的估計精度嚴重依賴于影響估計算法實時性的觀測模型及參數的精確性.方式3）利用前、后輪側偏角與質心側偏角之間的運動幾何關系進行質心側偏角的估計[12]，特點是模型簡單、計算量小、無積分操作，對于量測噪聲、傳感器偏置和部分車輛參數的不確定性具有一定的魯棒性[13].方式3）所用運動幾何關系式中的輪胎側偏剛度取為定值，當輪胎側向力處于線性階段時，估計效果較好；當輪胎側向力處于非線性階段時，估計結果偏差較大[14].You 等[15]采用輪胎側偏剛度辨識方法來適應輪胎側偏剛度的非線性變化，該方法依賴路面附著系數的精確獲取.單一的估計方式都存在不足，因此學者們開始嘗試采用多方法融合估計方式.Nishio 等[16]將路面附著系數作為判斷條件，以實現運動學與動力學估計方式的邏輯切換.相比動力學方式，運動幾何方式結構更簡單，計算量更小，對環境噪聲、傳感器偏置具有很好的魯棒性[17].高博麟等[18]為運動學、動力學與運動幾何方式分別設計子濾波器，使用主濾波器對子濾波器結果進行信息融合，得到全局融合結果.該方法所用的動力學子濾波器對輪胎模型進行簡化，同時運動學子濾波器未考慮量測噪聲統計特性的不確定性.Li 等[5]以動態特性的強弱為融合指標，對運動幾何方式與運動學方式的估計結果進行擇優融合；為了克服質心側偏角的估計精度對于輪胎側偏剛度精確性的依賴，利用PI 控制器對輪胎側偏剛度進行實時反饋調節.該調整依賴于路面附著系數與輪胎模型的精確獲取，Li 等[5]還將運動學估計方式的量測噪聲統計特性(協方差矩陣R)預先設為定值，這與實際不符.

單一的估計方式向多方法融合估計方式的發展，一定程度上提高了質心側偏角的估計精度，但融合估計的方式仍然待深入研究：1)盡管有學者利用PI 控制器對輪胎側偏剛度進行了較準確的實時反饋調整，但該運動幾何估計方式依賴于輪胎模型及路面附著系數等參數已知.2)量測噪聲受環境影響較大，其統計特性具有不確定性.在融合估計過程中，量測噪聲的統計特性根據先驗值預先設為定值，但定值化處理不僅會降低估計精度，甚至會導致濾波發散.針對上述問題，本研究基于改進的質心側偏角運動學與運動幾何估計方法，設計質心側偏角融合估計算法.

1 車輛模型

如圖1 所示，傳統線性二自由度無法描述車輛在非線性狀態下的運動學關系，因此采用包含縱向、橫擺與橫向三自由度的車輛模型來表征車輛不同狀態與參數之間的數學關系.圖中， δf為前輪 轉 角， αf， αr分 別 為 前、后 輪 側 偏 角，Fyf、Fyr分別為前、后輪側向力，Fxf、Fxr分別為前、后輪縱向 力，vy、vx分 別 為 側 向、縱 向 速 度， ωr為 橫 擺 角速度， β 為質心側偏角，lf、lr為質心到前、后軸的距離，m為車輛質量，Iz為整車繞z軸轉動慣量，整車平衡方程由縱向、側向和橫擺運動組成，描述如下.

圖1 三自由度車輛模型Fig.1 Three degrees of freedom vehicle model

整車縱向運動：

整車側向運動：

整車橫擺運動：

為了后續觀測器設計，將式(1)～(3)變形得到

前后輪的輪胎側偏角描述為

考慮到縱向加速度對垂直載荷的影響，將垂直載荷描述為

式中：Fzf、Fzr分別為前、后軸垂直載荷，g為重力加速度，hg為質心高度.轉向時橫擺運動的產生伴隨著左、右車輪速度的不一致，前輪轉角也會影響車速與輪速之間的運動幾何關系.

式中： ωi j為車輪角速度，r車輪半徑，bf、br分別為前、后輪的輪距.

2 質心側偏角融合觀測器設計

2.1 質心側偏角融合觀測器整體結構

如圖2 所示，質心側偏角融合觀測器分為4 個模塊：1)自適應擴展卡爾曼濾波器(adaptive extended Kalman filter, AEKF)、2)運動幾何觀測器、3)融合權重計算、4)質心側偏角融合估計值的計算.模塊1）基于車輛縱、側向速度，橫擺角速度與縱、側向加速度之間的運動學方程；構建具有量測噪聲協方差R自適應調整功能的擴展卡爾曼觀測器(EKF)；實現對車輛縱、側向速度的估計.模塊2）由側向車速的代數計算式得到側向車速估計值vy_m，利用模糊控制器對輪胎側偏剛度進行實時校正，以提高側向車速的估計精度.模塊3）基于側向車速融合結果的變化率v˙y_fuse提取瞬態特征值，再利用提取的瞬態特性值 S tbvy進行融合權重分配計算.模塊4）依據融合權重，將模塊1）、2）的側向速度估計結果進行加權融合，得到側向速度融合估計值vy_fuse，再結合AEKF 估計的縱向車速vx_k，計算得到車輛質心側偏角的融合估計值 βfuse.融合觀測器利用車載傳感器獲取的信息有車輛縱向速度vx_s，橫擺角速度 ωr_s，車輪輪 速 ωi j，前輪 轉 角δf以 及車 輛 縱、側 向 加速 度ax_s、ay_s.

圖2 質心側偏角融合觀測器結構圖Fig.2 Structure diagram of sideslip angle fusion observer

2.2 自適應擴展卡爾曼濾波器的設計

化簡式(4)、(5)得到車輛運動學方程v˙x=ax+ωrvy，v˙y=ay-ωrvx，選取量測方程為式(11)～(14)，依據狀態估計流程，將運動學方程整理成狀態方程和測量方程標準形式：

取采樣時間為Ts=1 ms，離散化得到差分方程：

其 中Ak=I+A(t)Ts、Hk=H(t)，I為 與A(t) 同 維 的單位矩陣.

擴展卡爾曼濾波算法分為時間更新和量測更新2 個部分，算法流程如圖3 所示.EKF 將量測協方差矩陣R視為定值，但在車輛實際行駛過程中，外界的環境噪聲統計特性具有不確定性，將R視為定值會引起估計偏差.本研究在EKF 基礎上加入量測協方差矩陣自適應算法，Zhang 等[10]指出新息能夠很好地反應環境噪聲對估計誤差的影響，為此，可根據新息變化設計量測協方差矩陣自適應調整算法.定義新息 εk為量測量的實際值yk與預測值之間的誤差：

圖3 擴展卡爾曼濾波算法流程圖Fig.3 Flowchart of extended Kalman filter algorithm

新息的理論協方差矩陣

式中：E為期望值，為k時刻誤差協方差矩陣先驗值.新息的實際協方差矩陣

式中：M為滑動窗口的長度，取M=50.依據系統的實際新息協方差與理論協方差的相對大小，實時更新R，調整因子 αk的表達式[19]為

式中：tr() 為矩陣的跡，n為量測向量的維度.當Sk>S?k時，應 將R調 小，此 時 αk<1 ；當Sk≤S?k時，理論上應將R調大，考慮到將R值增大可能導致濾波發散，保持R不變，R=1.利用調整因子調整卡爾曼增益計算式：其中 αk可對R進行 自適應修正.

2.3 運動幾何觀測器設計

假設輪胎側向力處于線性區域，則前、后輪胎側向力分別為Fyf=Cfαf、Fyr=Crαr，其中Cf、Cr分別為前、后輪的側偏剛度.定義Cf=kaCr，其中ka=Fzf/Fzr[14]，聯合式(7)、(8)可以得到

由may_s=Fyf+Fyr、vx_s=vx、ωr_s=ωr，聯立式(23)、(24)得到側向速度代數計算式為

式中：ay_s、ωr_s、vx_s分別為傳感器測得的側向加速度、橫擺角速度和縱向速度.式(25)只在推導過程中涉及輪胎力，實際應用時無需計算輪胎力.式(25)為代數計算式，無積分累計過程，即使量測信號存在偏差，估計結果也不會因累計誤差而發散；式中除去可通過車載傳感器與CAN 總線獲取的狀態量與輸入量、整車質量m與質心到前后軸距離外，精確地獲取輪胎側偏剛度值是得到準確的側向速度估計值vy_m的關鍵，在實際應用中，輪胎側偏剛度受到垂直載荷、輪胎氣壓和路面附著條件等因素影響，難以獲取精確的輪胎側偏剛度值.如圖4 所示為帶有輪胎側偏剛度自適應反饋調節的運動幾何觀測器.圖中，側向速度估計值vy_m通過式(25)計算得到，聯合車輛側向速度估計值、傳感器測得的車輛縱向速度的變化率v˙x_s、橫 擺 角速 度 ωr_s， 獲取車 輛縱向 加速度 的計算 值ax，模糊控制器利用傳感器測得的車輛縱向加速度ax_s與ax的差值e及其變化率er，對輪胎側偏剛度進行實時反饋調節.該調節方式未涉及輪胎模型，因此對于路面附著系數變化具有魯棒性.模糊控制器設計過程如下.1）確定模糊控制的輸入與輸出：當輪胎側偏剛度越準確，縱向加速度量測值ax_s與理論計算值ax之差越接近零；反之，其絕對值將越大.將e、er作為模糊控制的輸入，以后輪胎側偏剛度Cr作為模糊控制的輸出.2）定義模糊控制的輸入輸出模糊集：NVB、NB、NVM、NM、NVS、NS、Z、PS、PM、PB 分別表示負很大、負大、負很中、負中、負很小、負小、零、正小、正中、正大.3）定義輸入輸出的隸屬度函數Df，如圖5 所示.4）通過大量的仿真試驗與理論分析，積累足夠的專家經驗后，制定如表1 所示的模糊規則.利用Matlab 的模糊邏輯工具箱生成Simulink 中模糊控制器模塊的可執行文件，其中模糊控制器的解模糊化方法為重心法.

表1 后輪胎側偏剛度的模糊規則Tab.1 Fuzzy rule of rear tyre lateral stiffness

圖4 運動幾何觀測器流程圖Fig.4 Flowchart of kinematic-geometric observer

圖5 隸屬度函數Fig.5 Membership function

2.4 2 種估計方式的對比分析

運動學方式在方向盤處于機動狀態(車輛狀態的變化率較大)時有較好的估計效果[20]，為了分析所提運動幾何估計方式在不同車輛狀態變化率下的具體表現，參考Li 等[5]探究不同估計方式適用范圍的試驗方案，設置階躍-回正試驗進行仿真分析.

如圖6 所示為在vx=90 km/h，路面附著系數μ=0.85條件下，階躍-回正仿真試驗結果對比圖.稱圖6(b)中側向速度變化率較大的填充區域為瞬態區域，該區域車輛側向速度變化較劇烈；根據有無轉向輸入，將剩余側向速度變化率v˙y較小的區域分別稱為非零、零穩態區域.依據圖6 分析2 種估計方式的特點：在瞬態區域t=4～7 s 和t=15～17 s，運動學方式相比于運動幾何方式具有更好的實時跟隨估計效果；在非零穩態區域t=7～15 s，運動學方式由于積分累計誤差的影響，估計偏差呈現逐漸增大的趨勢，運動幾何方式的估計結果雖存在一定偏差，但偏差值相對較小，且并未隨時間增大；在方向盤回正后的零穩態區域t=17～20 s，運動學方式的估計結果依然存在偏差，運動幾何方式的估計偏差接近于零.經上述分析可知，運動學方式在瞬態區域估計效果相對較好，運動幾何方式在零、非零穩態區域內有可接受的估計精度.

圖6 運動學與運動幾何方式的估計結果對比Fig.6 Comparison of estimation results between kinematic method and kinematic-geometry method

2.5 融合估計算法設計與驗證

2.5.1 融合估計算法設計 本研究提出以狀態量的瞬態特性強弱（用于表征狀態量變化率的大小）為融合指標的融合算法，該算法的流程圖如圖7 所示.Cheli 等[20]將側向加速度、橫擺角速度和方向盤轉角在0～0.1 s 的均方根加權之和作為衡量狀態量穩態特性強弱的指標，用于融合算法的設計.借鑒該思想，將側向速度變化率v˙y_fuse作為瞬態特性指標的提取對象.圖中緩存器輸出當前時刻前推0.1 s 的側向速度變化率，接著利用側向速度變化率進行均方根計算，計算結果通過一階低通濾波器去除噪聲，得到瞬態特性指標值Stbvy.將S tbvy輸入融合權重計算模塊，該模塊利用Stbvy與前輪轉角 δf進行融合權重的分配計算.融合權重計算式為

圖7 融合算法流程Fig.7 Fusion algorithm process

1) 當 S tbvy小 于 閾 值 S tbvy_t1， δf小 于 閾 值 δf_t1時，車輛狀態處于零穩態區域，依據2.4 節的結論，該區域運動幾何方式的估計結果偏差接近于零，因此，給與運動幾何方式的融合權重w=1.00.2)當Stbvy小 于 閾 值 S tbvy_t2時，如 果 δf小 于 閾 值 δf_t1，融合權重依據式(27)計算，該區域為瞬態與零穩態區 域 的 過 渡 區 間.如 果 δf大 于 閾 值 δf_t1，融 合 權 重依據式(28)計算，該區域為瞬態與非零穩態區域的過渡區間.設置過渡區間可避免融合權重的突變導致融合估計結果的突變.3)當 S tbvy小于設定的 閾 值 S tbvy_t1， δf大 于 設 定 的 閾 值 δf_t1，表 明 車 輛狀態處于非零穩態區域，以運動幾何方式的估計結果為主，融合權重為ks(依據運動幾何方式穩態偏差大小而定)，本研究取0.75.4)當 S tbvy大于設定的閾值 S tbvy_t2，此時車輛狀態處于瞬態區域，以運動學方式的估計結果為主，融合權重w=0.06.

2.5.2 融合估計算法驗證 階躍-回正試驗進行融合算法的驗證，取 S tbvy_t1=0.01, Stbvy_t2=0.1, δf_t1=0.3°，ks=0.75.如圖8 所示為融合權重分配結果圖.在零穩態區域(t=0～4 s、t=17～20 s)與非零穩態區域(t=7～15 s)，瞬態特性值 S tbvy較小，運動幾何法的權重為1；瞬態特性值 S tbvy高 的區域(t=4～7 s、t=15～17 s)運動幾何的融合權重為0.06，不同區域之間的融合權重也能平滑過渡；輸出的融合權重能依據瞬態特性值大小進行動態變化，符合設計要求.

圖8 融合權重分配圖Fig.8 Fusion weight distribution chart

3 仿真驗證

采用Carsim-Simulink 的仿真與駕駛員在環(DIL)試驗平臺，以Carsim 輸出的值為參考值，驗證本研究提出的融合觀測器的有效性，整車參數如表2 所示.模擬真實的傳感器量測信號，須為Carsim 的理想輸出添加合適的零均值高斯白噪聲，各量測輸入所加噪聲的方差[21]取值如下：Varωi j= 0.002 rad/s、Varax= 0.001 m/s2、Varay= 0.002 m/s2、Varωr= 0.004 rad/s、 V arvx= 0.001m/s.設計角階躍、雙移線、變附著系數對接路面-雙移線等仿真工況以及駕駛員在環試驗，對所提融合觀測器的估計準確性進行驗證.為了體現矩陣R自適應對于運動學方式估計精度的提高效果，仿真結果對有無矩陣R自適應的運動學方式進行比較.為了定量分析估計誤差，選用均方根誤差作為評價指標：

表2 Carsim 車輛參數Tab.2 Parameters of Carsim vehicle

式中： βi、^βi為質心側偏角的真實值和估計值.

1) 雙移線、角階躍工況.初始車速vx= 100 km/h，μ=0.85，仿真時間t=10 s，方向盤轉角如圖9 所示.依據仿真結果圖10 和表3 可知，有自適應相比于無自適應的運動學方式，雙移線工況、角階躍工況的估計結果的均方根誤差值分別減少了37%、51%；融合方式相比于自適應運動學方式與運動幾何方式，估計結果的均方根誤差值在雙移線工況中分別減少了24%、58%，在角階躍工況中分別減少了56%、38%.因此加入自適應可提高運動學方式的估計精度，與此同時，融合運動學方式（自適應）與運動幾何方式的估計結果可進一步提高估計精度.

表3 2 種工況估計結果的均方根誤差Tab.3 Root mean squared error of estimation results under two conditions

圖9 2 種工況的方向盤轉角Fig.9 Steering angles of two conditions

圖10 2 種工況的質心側偏角估計結果Fig.10 Results of sideslip angle estimation for two conditions

2) 變路面附著系數對接路面-雙移線工況.為了驗證所提觀測器對路面附著系數變化的魯棒性，選取文獻[5]中提出的融合估計方式作為類比對象，記為融合方式B，該方式的估計精度依賴于路面附著系數的準確性.所提觀測器記為融合方式A.利用變路面附著系數對接路面-雙移線工況進行驗證，在Carsim 中設置長直道，縱向0～300 m路面附著系數 μ=0.85 ，300～600 m 路面附著系數μ=0.5， 橫向方向路面附著系數一致.車速vx= 90 km/h，方向盤轉角如圖11(a)，仿真時間t=25 s，仿真結果如圖11(b)所示，融合方式A、B 估計結果的均方根誤差分別為0.033、0.126.當t=0～15 s，融合方式A 的估計精度略高于融合方式B；t> 15 s，路面附著系數由0.85 突變到0.50，融合方式A 依然能夠得到較好的估計結果；由于算法中的路面附著系數不能及時更新為0.5，融合方法B 的估計精度明顯降低，融合方法B 總體估計結果的均方根誤差值相比于融合方法A 增加了73%.

3) DIL 試驗.為了進一步驗證所提觀測器的有效性與實時性，開展DIL，如圖12 所示為DIL仿 真 平 臺.初 始 車 速vx= 90 km/h ， μ =0.85 ，分 別 進行雙移線、角階躍試驗，試驗結果如圖13 所示，均方根誤差如表4 所示.相比于自適應運動學估計方式與運動幾何方式，本研究所提融合估計方式的均方差誤差值均有所減少，雙移線試驗中分別減少了34%、56%，角階躍試驗中分別減少了61%、47%.在試驗所用美國國家儀器(national instruments, NI)機箱的數據輸出頻率為1 000 Hz，即時間間隔為1 ms，觀測器每一步仿真平均用時(simulation time)分別為350、321 μs ，均小于設置的仿真步長，實時性滿足要求.

圖12 駕駛員在環仿真平臺Fig.12 Driver in loop simulation platform

圖13 2 種工況的駕駛員在環試驗質心側偏角估計結果Fig.13 Results of sideslip angle estimation of two conditions for driver in loop simulation test

4 結 語

設計可實現量測噪聲協方差矩陣自適應以及輪胎側偏剛度模糊自適應的質心側偏角融合估計算法.在Carsim-Simulink 的仿真與駕駛員在環仿真平臺展開的多工況試驗結果顯示，本研究設計的質心側偏角融合觀測器可實時獲取較準確的質心側偏角估計結果，即使在路面附著系數突變條件下，也能實現較準確估計.所提融合估計算法的驗證僅在Carsim-Simulink 以及駕駛員在環仿真平臺進行，后續考慮進行實車試驗，對所提算法做進一步驗證.