基于UKF与修正Dugoff轮胎模型路面峰值附着系数估算方法

摘要

本发明公开一种基于UKF与修正Dugoff轮胎模型路面峰值附着系数估算方法，通过：实时采集车辆的各种传感器信号，利用车辆纵向动力学方程和模型的几何坐标关系对各轮滑移率和侧偏角估计；然后将估计的滑移率、垂向力、侧偏角等传给基于修正Dugoff模型的UKF系数计算模块，得到非线性系统的系数向量，将此向量与实时估计的纵向力发送到UKF路面峰值附着系数估计模块，求取峰值附着系数。本方法应用车辆状态观测系统实时采集信号，保证了计算的实时性，对于没有拟合过的路面情况估计准确度高。应用修正Dugoff轮胎模型和UKF理论，使得求解过程简单，运算量小、快捷，收敛时间短。本方法鲁棒性良好，能够较好的识别各轮的路面情况，适用于路面峰值附着系数的实时估计。

说明书

技术领域

本发明涉及一种分布式驱动电动汽车的路面峰值附着系数估算方法，特别是关于 一种基于UKF与修正Dugoff轮胎模型的路面峰值附着系数估算方法。

背景技术

路面峰值附着系数：是指轮胎与地面间作用的纵向力、侧向力的合力与垂向力之 比的最大值。目前国内外对于路面峰值附着系数实时估算方法已经进行了大量研究。 这些方法可以分为基于原因的方法和基于效果的方法两类。前者利用超声波传感器等 来检测路面状况来估算路面附着系数，该类方法需要外加昂贵的传感器，并且对于环境 的依赖程度较高。后者方法则是直接利用车辆与轮胎的动力学特性来估计路面附着系 数，例如用μ-s曲线斜率(附着系数与滑移率曲线)估算路面附着系数的方法、利用 侧向力或回正力矩与轮胎侧偏角的关系估计的方法。该类方法虽然不用加装额外的传 感器，但各方法使用的工况有局限性，不能实现纵、侧向力耦合存在的路面峰值附着 系数实时估算，且多数方法仅能实现转向轮或者驱动轮的路面峰值附着系数估计。

从观测器的角度分析，大部分学者采用基本卡尔曼滤波算法或者扩展卡尔曼算法 来解决小滑移率条件下(线性轮胎特性区域)的路面峰值附着系数实时估算，这两种 滤波器不适于应用于强非线性区下的路面峰值附着系数实时估算。而无迹卡尔曼滤波 (UKF)算法由于使用无迹(UT)变换来处理均值和协方差的非线性传递，对非线性函数的 概率密度分布进行近似，用一系列确定样本来逼近状态的后验概率密度，不需要求导 计算Jacobian矩阵。UKF没有线性化忽略高阶项，因此它特别适用于强非线性区下的 路面峰值附着系数实时估算。

从轮胎模型的角度分析，相关研究为保证求解的实时性，应用较为简单的轮胎模 型，如Dugoff轮胎模型、多项式轮胎模型，但这些模型在强非线性区的精度较低。而 在2014年美国汽车工程学会年会上发表的修正Dugoff轮胎模型能够在保持较简单的 形式的前提下，较大的提高了在强非线性区的精度。

分布式电驱动车辆区别于传统内燃机驱动的车辆，其四个车轮由四个轮毂电机直 接驱动，具有响应快速、传动链短等优点。电机的驱动转矩精确可知，这样就可以实 时准确估计出轮胎的纵向力。

发明内容

针对现有的路面峰值附着系数算法只适用工况单一的问题，本发明目的是提供一 种基于UKF与修正Dugoff轮胎模型路面峰值附着系数估算方法，用以实时准确的估计 各种工况下的路面情况。

为实现上述目的，本发明采取以下技术方案：一种基于UKF与修正Dugoff轮胎模 型路面峰值附着系数估算方法，适用于分布式驱动电动汽车每个车轮的路面峰值附着 系数估算，该方法包括以下步骤：

1)建立一种车辆路面峰值附着系数估算系统，包括用以获得参数的：一设置在车 辆驱动系统上的驱动力矩传感器、一设置在转向柱管处的方向盘转角传感器、一设置 在车轮处的轮速传感器、一设置在整车质心处的车速传感器和惯性传感器，以及包括 用以计算处理的：一车轮滑移率计算模块、一基于纵向动力学的轮胎力估计模块、一 侧偏角计算模块、一基于修正Dugoff模型的UKF系数计算模块和一UKF路面峰值附着 系数估计模块；

2)估算过程为：在车辆运行过程中，整车控制器取某两个相邻的采样时刻k-1和 k，分别从驱动力矩传感器接收车辆的驱动力矩信号T_i(k)，从轮速传感器接收实时的 轮速信号ω_i(k-1)和ω_i(k)，从车速传感器接收实时的纵向车速信号v_x(k)以及侧向车速 信号v_y(k)，从惯性传感器接收实时的纵向加速度信号a_x(k)和侧向加速度信号a_y(k)以 及横摆角速度信号发送到基于纵向动力学的轮胎力估计模块、侧偏角计算模块 和车轮滑移率计算模块；

3)基于纵向动力学的轮胎力估计模块根据实时采集到的驱动力矩信号T_i(k)、轮 速信号ω_i(k-1)和ω_i(k)，纵向加速度信号a_x(k)和侧向加速度信号a_y(k)，计算出各车轮 的纵向力F_xi(k)、垂向力F_zi(k)，将垂向力F_zi(k)发送到基于修正Dugoff模型的UKF系 数计算模块；将纵向力F_xi(k)发送到UKF路面峰值附着系数估计模块；

4)车轮侧偏角计算模块根据接收到的实时方向盘转角信号ω_sw(k)、纵向车速信号 v_x(k)、侧向车速信号v_y(k)和横摆角速度信号计算各车轮的侧偏角α_i(k)，并将结 果发送给车轮滑移率计算模块和基于修正Dugoff模型的UKF系数计算模块；

5)车轮滑移率计算模块根据接收到的实时纵向车速信号v_x(k)、侧向车速信号v_y(k) 和横摆角速度信号各车轮的侧偏角α_i(k)和轮速信号ω_i(k)计算车轮的滑移率 s_i(k)，并将结果发送给基于修正Dugoff模型的UKF系数计算模块；

6)基于修正Dugoff模型的UKF系数计算模块根据侧偏角计算模块、滑移率计算 模块和基于纵向动力学的轮胎力估计模块发送过来的侧偏角α_i(k)、滑移率s_i(k)、车轮 的垂向力F_zi(k)，利用修正Dugoff模型算法计算得到建立UKF滤波器所需的非线性系 统的系数a_1i(k)、a_2i(k)、a_3i(k)、a_4i(k)、a_5i(k)，并将结果发送给UKF路面峰值附着系数估 计模块；

以上，i＝1、2、3、4，代表4个车轮。

7)UKF路面峰值附着系数估计模块根据接收到的实时纵向力F_xi(k)和系数 结合UKF观测方法得到各轮的路面峰值附着系数

所述步骤3)中，车轮纵向力F_xi(k)、垂向力F_zi(k)的计算方法为：

①基于纵向动力学的轮胎力估计模块根据前、后两个相邻采样时刻k-1和k实时 采集到的轮速信号ω_i(k-1)和ω_i(k)，首先计算得到在采样时刻k时的车轮角加速度

${\stackrel{·}{\omega }}_i\left(k\right)=\frac{{\omega }_i\left(k\right)-{\omega }_i(k-1)}{T}---\left(1\right)$

式中，T为采样步长；

②根据进一步计算得到车轮纵向力F_xi(k)：

$F_{\mathrm{xi}}\left(k\right)=\frac{T_i\left(k\right)-j{\stackrel{·}{\omega }}_i\left(k\right)}{R}---\left(2\right)$

式中，J为车轮转动惯量，R为车轮滚动半径；

③再根据惯性传感器在k时的车辆纵向加速度信号a_x(k)和侧向加速度信号a_y(k)， 进一步计算得到各车轮垂向力F_zi：

F_z1(k)＝mgl_r/(2l)-ma_x(k)h/(2l)-ma_y(k)hl_r/(lB)

F_z2(k)＝mgl_r/(2l)-ma_x(k)h/(2l)+ma_y(k)hl_r/(lB)

F_z3(k)＝mgl_f/(2l)+ma_x(k)h/(2l)-ma_y(k)hl_f/(lB)           (3)

F_z4(k)＝mgl_f/(2l)+ma_x(k)h/(2l)+ma_y(k)hl_f/(lB)

式中，l为轴距，B为轮距，l_r为质心到后轴距离，l_f为质心到前轴距离，h为质 心高度，m为车质量，g为重力加速度。

所述步骤4)中，各车轮的侧偏角α_i(k)的计算方法为：

①首先车轮侧偏角计算模块根据采样时刻k实时采集到的实时方向盘转角信号 ω_sw(k)，经转向传动比换算得到前面两转向轮的转角θ₁(k)和θ₂(k)；

②再根据实时纵向车速信号v_x(k)和侧向车速信号v_y(k)和横摆角速度信号计 算各车轮的侧偏角α_i(k)：

$\left(\begin{array}{c}{\alpha }_1\left(k\right)=-{\theta }_1\left(k\right)+\arctan \frac{v_y\left(k\right)+l_f\stackrel{·}{\psi }\left(k\right)}{v_x\left(k\right)-B\stackrel{·}{\psi }\left(k\right)}\\ {\alpha }_2\left(k\right)=-{\theta }_2\left(k\right)+\arctan \frac{v_y\left(k\right)+l_f\stackrel{·}{\psi }\left(k\right)}{v_x\left(k\right)+B\stackrel{·}{\psi }\left(k\right)}\\ {\alpha }_3\left(k\right)=\arctan \frac{v_y\left(k\right)-l_r\stackrel{·}{\psi }\left(k\right)}{v_x\left(k\right)-B\stackrel{·}{\psi }\left(k\right)}\\ {\alpha }_4\left(k\right)=\arctan \frac{v_y\left(k\right)-l_r\stackrel{·}{\psi }\left(k\right)}{v_x\left(k\right)+B\stackrel{·}{\psi }\left(k\right)}\end{array}\right)---\left(4\right)$

所述步骤5)中，计算车轮的滑移率s_i(k)，利用如下公式计算：

$\left(\begin{array}{c}{s}_1\left(k\right)=\frac{|\sqrt{[{v_y\left(k\right)+l_f\stackrel{·}{\psi }\left(k\right)]}^2+{[v_x\left(k\right)-B\stackrel{·}{\psi }\left(k\right)]}^2}\cos {\alpha }_1\left(k\right)-R{\omega }_1\left(k\right)|}{\max |\sqrt{[{v_y\left(k\right)+l_f\stackrel{·}{\psi }\left(k\right)]}^2+{[v_x\left(k\right)-B\stackrel{·}{\psi }\left(k\right)]}^2}\cos {\alpha }_1\left(k\right),R{\omega }_1\left(k\right)|}\\ s_2\left(k\right)=\frac{|\sqrt{{\left[v_y\left(k\right)+l_f\stackrel{·}{\psi }\left(k\right)\right]}^2+{[v_x\left(k\right)+B\stackrel{·}{\psi }\left(k\right)]}^2}\cos {\alpha }_2\left(k\right)-R{\omega }_2\left(k\right)|}{\max |\sqrt{{[v_y\left(k\right)+l_f\stackrel{·}{\psi }\left(k\right)]}^2+{(v_x\left(k\right)+B\stackrel{·}{\psi }\left(k\right))}^2}\cos {\alpha }_2\left(k\right),R{\omega }_2\left(k\right)|}\\ s_3\left(k\right)=\frac{|v_x\left(k\right)-R{\omega }_3\left(k\right)|}{\max |v_x\left(k\right),R{\omega }_3\left(k\right)|}\\ s_4\left(k\right)=\frac{|v_x\left(k\right)-R{\omega }_4\left(k\right)|}{\max |v_x\left(k\right),R{\omega }_4\left(k\right)|}\end{array}\right)---\left(5\right)$

所述步骤6)中，基于修正Dugoff模型的UKF系数计算模块，利用修正Dugoff 模型算法计算UKF非线性系统的系数的方法如下：

$\left(\begin{array}{c}{a}_{1i}\left(k\right)=\frac{{3C}_s{s}_i^2\left(k\right)F_{\mathrm{zi}}^2\left(k\right)[-s_i\left(k\right)+1][1+s_i\left(k\right)]}{16\{{\left[C_s{s}_i\left(k\right)\right]}^2+{\left[C_{\alpha }\tan {\alpha }_i\left(k\right)\right]}^2\}}\\ a_{2i}\left(k\right)=\frac{-3C_s{s}_i\left(k\right)F_{\mathrm{zi}}\left(k\right)[-s_i^2\left(k\right)+s_i\left(k\right)][1+s_i\left(k\right)]}{4\sqrt{{\left[C_s{s}_i\left(k\right)\right]}^2+{\left[C_{\alpha }\tan {\alpha }_i\left(k\right)\right]}^2[1+s_i\left(k\right)]}}\\ -\frac{[1.15s_i^2\left(k\right)+1.63s_i\left(k\right)+1.27]{\left\{F_{\mathrm{zi}}\left(k\right)\right[1+s_i\left(k\right)\left]\right\}}^2}{4\{{\left[C_s{s}_i\left(k\right)\right]}^2+{\left[C_{\alpha }\tan {\alpha }_i\left(k\right)\right]}^2\}}\\ a_{3i}\left(k\right)=\frac{-C_s{s}_i\left(k\right)F_{\mathrm{zi}}\left(k\right)[1.15s_i^2\left(k\right)+1.63s_i\left(k\right)+1.27]}{\sqrt{{\left[C_s{s}_i\left(k\right)\right]}^2+{\left[C_{\alpha }\tan {\alpha }_i\left(k\right)\right]}^2}}\\ a_{4i}\left(k\right)=\frac{{3C}_s{s}_i^2\left(k\right)[-s_i\left(k\right)+1]}{4[1+s_i\left(k\right)]}\\ a_{5i}\left(k\right)=\frac{C_s{s}_i\left(k\right)[1.15s_i^2\left(k\right)+1.63s_i\left(k\right)+1.27]}{[1+s_i\left(k\right)]}\end{array}\right)---\left(6\right)$

式中，C_s、C_α分别是轮胎的纵向滑移刚度和侧偏刚度，是轮胎的参数。

步骤7)中，UKF路面峰值附着系数估计模块利用UKF观测方法得到各轮的路面峰 值附着系数具体实施过程如下：

①根据修正Dugoff轮胎模型，各轮胎的纵向力可表示为其所在路面峰值附着系数 μ_maxi和该车轮状态的函数，如下所示：

$F_{\mathrm{xi}}\left({\mu }_{\max i}\right)=\left(\begin{array}{cc}{a}_{1i}{\mu }_{\max i}^3+a_{2i}{\mu }_{\max i}^2+a_{3i}{\mu }_{\max i}& {\lambda }_i<1\\ a_{4i}{\mu }_{\max i}+a_{5i}& {\lambda }_i\ge 1\end{array}\right)---\left(7\right)$

${\lambda }_i\left(k\right)=\frac{{\mu }_{\max i}\left(k\right)·F_{\mathrm{zi}}\left(k\right)·[1+s_i\left(k\right)]}{2\sqrt{{\left[C_s{s}_i\left(k\right)\right]}^2+{\left[C_{\alpha }\tan {\alpha }_i\left(k\right)\right]}^2}}---\left(8\right)$

②以上式为基础建立无味卡尔曼滤波器(UKF)：

i)选取状态量和观测量如下：

x_i＝μ_maxi            (9)

y_i＝F_xi

ii)将二者关系表述为非线性系统，如下：

$\left(\begin{array}{c}{x}_i(k+1)=x_i\left(k\right)+w_i\left(k\right)\\ y_i\left(k\right)=\left(\begin{array}{cc}{a}_{1i}\left(k\right)x_i^3\left(k\right)+a_{2i}\left(k\right)x_i^2\left(k\right)+a_{3i}\left(k\right)x_i\left(k\right)+v_i\left(k\right)& {\lambda }_i\left(k\right)<1\\ a_{4i}{\left(k\right)x}_i\left(k\right)+a_{5i}\left(k\right)+v_i\left(k\right)& {\lambda }_i\left(k\right)\ge 1\end{array}\right)\\ {\lambda }_i\left(k\right)=\frac{x_i\left(k\right)·F_{\mathrm{zi}}\left(k\right)·[1+s_i\left(k\right)]}{2\sqrt{{\left[C_s{s}_i\left(k\right)\right]}^2+{\left[C_{\alpha }\tan {\alpha }_i\left(k\right)\right]}^2}}\end{array}\right)---\left(10\right)$

iii)均值和方差的初始化：

$\left(\begin{array}{c}{\hat{x}}_i\left(0\right)=1\\ P_i\left(0\right)=0\end{array}\right)---\left(11\right)$

iv)无迹变换：

$\left(\begin{array}{c}{\hat{x}}_i^{\left(1\right)}(k-1)={\hat{x}}_i^{+}(k-1)+\sqrt{P_i^{+}(k-1)}\\ {\hat{x}}_i^2(k-1)={\hat{x}}_i^{+}(k-1)-\sqrt{P_i^{+}(k-1)}\end{array}\right)---\left(12\right)$

和是对k-1时刻的状态更新后的滤波值进行无迹变换得到的两个 点，P_i⁺(k-1)为k-1时刻的状态更新后的后验方差矩阵，为k-1时刻的状态 更新后的滤波值。

v)时间更新：

利用非线性方程进行非线性变换：

$\left(\begin{array}{c}{\hat{x}}_i^{\left(1\right)}\left(k\right)={\hat{x}}_i^{\left(1\right)}(k-1)\\ {\hat{x}}_i^{\left(2\right)}\left(k\right)={\hat{x}}_i^{\left(2\right)}(k-1)\end{array}\right)---\left(13\right)$

加权得到的预测的状态向量：

${\hat{x}}_i^{-}\left(k\right)={\hat{x}}_i^{\left(1\right)}\left(k\right)+{\hat{x}}_i^{\left(2\right)}\left(k\right)---\left(14\right)$

加权得到的预测的协方差矩阵：

$P_i^{-}\left(k\right)=\frac{1}{2}\{{[{\hat{x}}_i^{\left(1\right)}\left(k\right)-{\hat{x}}_i^{-}\left(k\right)]}^2+{[{\hat{x}}_i^{\left(2\right)}\left(k\right)-{\hat{x}}_i^{-}\left(k\right)]}^2\}+Q_i(k-1)---\left(15\right)$

Q_i(k-1)代表过程噪声矩阵。

利用非线性方程进行非线性变换：

先计算${\lambda }_i\left(k\right)=\frac{{\hat{x}}_i^{-}\left(k\right)·F_{\mathrm{zi}}\left(k\right)·[1+s_i\left(k\right)]}{2\sqrt{{\left[C_s{s}_i\left(k\right)\right]}^2+{\left[C_{\alpha }\tan {\alpha }_i\left(k\right)\right]}^2}}$

根据λ_i(k)计算

$\left(\begin{array}{c}{\hat{y}}_i^{\left(1\right)}\left(k\right)=\left(\begin{array}{cc}{a}_{1i}\left(k\right){\left[x_i^{\left(1\right)}\left(k\right)\right]}^3+a_{2i}\left(k\right)[x_i^{\left(1\right)}\left(k\right){]}^2+a_{3i}\left(k\right)x_i^{\left(1\right)}\left(k\right)+v_i\left(k\right)& {\lambda }_i\left(k\right)<1\\ a_{4i}\left(k\right)x_i^{\left(1\right)}{\left(k\right)}^3+a_{5i}\left(k\right)+v_i\left(k\right)& {\lambda }_i\left(k\right)\ge 1\end{array}\right)\\ {\hat{y}}_i^{\left(2\right)}\left(k\right)=\left(\begin{array}{cc}{a}_{1i}\left(k\right){\left[x_i^{\left(2\right)}\left(k\right)\right]}^3+a_{2i}\left(k\right){\left[x_i^{\left(2\right)}\left(k\right)\right]}^2+a_{3i}\left(k\right)x_i^{\left(2\right)}\left(k\right)+v_i\left(k\right)& {\lambda }_i\left(k\right)<1\\ a_{4i}\left(k\right)x_i^{\left(2\right)}{\left(k\right)}^3+a_{5i}\left(k\right)+v_i\left(k\right)& {\lambda }_i\left(k\right)\ge 1\end{array}\right)\end{array}\right)---\left(16\right)$

加权得到系统的预测值：

${\hat{y}}_i\left(k\right)=\frac{1}{2}[{\hat{y}}_i^{\left(1\right)}\left(k\right)+{\hat{y}}_i^{\left(2\right)}\left(k\right)]---\left(17\right)$

vi)量测更新：

后验估计的协方差矩阵：

$P_{y-i}\left(k\right)=\frac{1}{2}\{{[{\hat{y}}_i^{\left(1\right)}\left(k\right)-{\hat{y}}_i\left(k\right)]}^2+{[{\hat{y}}_i^{\left(2\right)}\left(k\right)-{\hat{y}}_i\left(k\right)]}^2\}+R_i\left(k\right)---\left(18\right)$

R_i(k)代表量测噪声矩阵。

先验估计的协方差矩阵：

$P_{\mathrm{xy}-i}\left(k\right)=\frac{1}{2}\left\{\right[{\hat{y}}_i^{\left(1\right)}\left(k\right)-{\hat{y}}_i\left(k\right)\left]\right[{\hat{x}}_i^{\left(1\right)}\left(k\right)-{\hat{x}}_i^{-}\left(k\right)]+[{\hat{y}}_i^{\left(2\right)}\left(k\right)-{\hat{y}}_i\left(k\right)\left]\right[{\hat{x}}_i^{\left(2\right)}\left(k\right)-{\hat{x}}_i^{-}\left(k\right)\left]\right\}---\left(19\right)$

滤波增益矩阵：

$K_i\left(k\right)=P_{\mathrm{xy}-i}\left(k\right)P_{y-i}^{-1}\left(k\right)---\left(20\right)$

状态更新后的滤波值：

${\hat{x}}_i^{+}\left(k\right)={\hat{x}}_i^{-}\left(k\right)+K_i\left(k\right)·(y_i\left(k\right)-{\hat{y}}_i\left(k\right))---\left(21\right)$

状态更新后的后验方差矩阵：

$P_i^{+}\left(k\right)=P_i^{-}\left(k\right)-K_i^2\left(k\right)·P_{y-i}\left(k\right)---\left(22\right)$

到此为止，UKF滤波器建立完毕。滤波器输出的滤波值即为各轮处的附着系 数估计值

上述估算过程可以选取车辆4个车轮中的任意一个进行。

本发明由于采取以上技术方案，其具有以下优点：本发明通过建立一组估算系统， 首先采集驱动力矩信号、车辆的纵向速度信号、侧向速度信号、横摆加速度信号、方 向盘转角信号和轮速信号等，采取实时估计的方法，利用车辆纵向动力学方程和车辆 模型的几何坐标关系分别对各轮滑移率和侧偏角实时估计；然后将估计的滑移率值、 垂向力、侧偏角等传给车辆控制器中的基于修正Dugoff模型的UKF系数计算模块，计 算得到利用修正Dugoff轮胎模型推导出来的非线性系统的系数向量，将此向量与实时 与实时估计的纵向力发送到预先建好的UKF路面峰值附着系数估计模块，求取峰值附 着系数。本方法应用车辆状态观测系统实时采集信号，保证了计算的实时性，对于没有 拟合过的路面情况估计准确度高。本方法应用修正Dugoff轮胎模型和UKF理论，使得 求解过程简单，运算量小、快捷，收敛时间短；并且能在较宽滑移率范围内和纵向力 侧向力耦合的情况下计算，具有较高的估计准确性，适用范围广。本方法鲁棒性(在 不改变参数的情况下，能够对多种路况进行识别)良好，能够较好的识别各轮的路面 情况，适用于车辆在各种行驶过程中路面的峰值附着系数的实时估计。

附图说明

图1是本发明的系统关系示意图。

图2是本发明的车辆模型示意图。

具体实施方式

下面结合附图和实施例对本发明进行详细的描述。

图1所示是本发明的路面峰值附着系数实时估计方法所应用的系统图，外部硬件 包括一设置在车辆驱动系统上的驱动力矩传感器1、一设置在转向柱管处的方向盘转 角传感器2、一设置在车轮处的轮速传感器3，设置在整车质心处的一车速传感器和4 一惯性传感器5，设置于整车控制器中的处理模块包括一车轮滑移率计算模块6、一基 于纵向动力学的轮胎力估计模块7、一轮胎侧偏角计算模块8、一基于修正Dugoff模 型的UKF系数计算模块9和一UKF路面峰值附着系数估计模块10。上述驱动力矩传感 器1、轮速传感器3、惯性传感器5与基于纵向动力学的轮胎力估计模块7相连；方向 盘转角传感器2、车速传感器4、惯性传感器5与轮胎侧偏角计算模块8相连；车速传 感器4、惯性传感器5、轮速传感器3、轮胎侧偏角计算模块8与车轮滑移率计算模块 6相连；基于纵向动力学的轮胎力估计模块7、轮胎侧偏角计算模块8、车轮滑移率计 算模块6、车速传感器4，与基于修正Dugoff模型的UKF系数计算模块9相连；基于 纵向动力学的轮胎力估计模块7、基于修正Dugoff模型的UKF系数计算模块9与UKF 路面峰值附着系数估计模块10相连。

基于上述系统，本发明对车辆行驶过程中的路面峰值附着系数的实时估计方法包 括以下步骤：

1、在车辆运行过程中，整车控制器取某两个相邻的采样时刻k-1和k，分别从驱 动力矩传感器接收车辆的驱动力矩信号T_i(k)，从轮速传感器接收实时的轮速信号 ω_i(k-1)和ω_i(k)，从车速传感器接收实时的纵向车速信号v_x(k)以及侧向车速信号v_y(k)， 从惯性传感器接收实时的纵向加速度信号a_x(k)和侧向加速度信号a_y(k)以及横摆角速 度信号发送到基于纵向动力学的轮胎力估计模块、侧偏角计算模块和车轮滑移 率计算模块；

2、基于纵向动力学的轮胎力估计模块根据实时采集到的驱动力矩信号T_i(k)、轮 速信号ω_i(k-1)和ω_i(k)，纵向加速度信号a_x(k)和侧向加速度信号a_y(k)，计算出各车轮 的纵向力F_xi(k)、垂向力F_zi(k)，将垂向力F_zi(k)发送到基于修正Dugoff模型的UKF系 数计算模块；将纵向力F_xi(k)发送到UKF路面峰值附着系数估计模块；

3、轮胎侧偏角计算模块根据接收到的实时方向盘转角信号ω_sw(k)、纵向车速信号 v_x(k)、侧向车速信号v_y(k)和横摆角速度信号计算各车轮的侧偏角α_i(k)，并将结 果发送给车轮滑移率计算模块和基于修正Dugoff模型的UKF系数计算模块；

4、车轮滑移率计算模块根据接收到的实时纵向车速信号v_x(k)、侧向车速信号v_y(k) 和横摆角速度信号各车轮的侧偏角α_i(k)和轮速信号ω_i(k)计算车轮的滑移率 s_i(k)，并将结果发送给基于修正Dugoff模型的UKF系数计算模块；

5、基于修正Dugoff模型的UKF系数计算模块根据侧偏角计算模块、滑移率计算 模块和基于纵向动力学的轮胎力估计模块发送过来的侧偏角α_i(k)、滑移率s_i(k)、车轮 的垂向力F_zi(k)，利用修正Dugoff模型算法计算得到建立UKF滤波器所需的非线性系 统的系数a_1i(k)、a_2i(k)、a_3i(k)、a_4i(k)、a_5i(k)，并将结果发送给UKF路面峰值附着系数估 计模块；

以上，i＝1、2、3、4，代表4个车轮。

6、UKF路面峰值附着系数估计模块根据接收到的实时纵向力F_xi(k)和系数 a_1i(k)、(，结合UKF观测方法得到各轮的路面峰值附着系数 ${\hat{\mu }}_{\max i}.$

上述步骤2中，车轮纵向力F_xi(k)、垂向力F_zi(k)的计算方法为：

①基于纵向动力学的轮胎力估计模块根据前、后两个相邻采样时刻k-1和k实时 采集到的轮速信号ω_i(k-1)和ω_i(k)，首先计算得到在采样时刻k时的车轮角加速度 ${\stackrel{·}{\omega }}_i\left(k\right):$

${\stackrel{·}{\omega }}_i\left(k\right)=\frac{{\omega }_i\left(k\right)-{\omega }_i(k-1)}{T}---\left(\text{1}\right)$

式中，T为采样步长；

②根据进一步计算得到车轮纵向力F_xi(k)：

$F_{\mathrm{xi}}\left(k\right)=\frac{T_i\left(k\right)-J{\stackrel{·}{\omega }}_i\left(k\right)}{R}---\left(2\right)$

式中，J为车轮转动惯量，R为车轮滚动半径；

③再根据惯性传感器在k时的车辆纵向加速度信号a_x(k)和侧向加速度信号a_y(k)， 进一步计算得到各车轮垂向力F_zi：

F_z1(k)＝mgl_r/(2l)-ma_x(k)h/(2l)-ma_y(k)hl_r/(lB)

F_z2(k)＝mgl_r/(2l)-ma_x(k)h/(2l)+ma_y(k)hl_r/(lB)

F_z3(k)＝mgl_f/(2l)+ma_x(k)h/(2l)-ma_y(k)hl_f/(lB)  (3)

F_z4(k)＝mgl_f/(2l)+ma_x(k)h/(2l)+ma_y(k)hl_f/(lB)

式中，l为轴距，B为轮距，l_r为质心到后轴距离，l_f为质心到前轴距离，h为质 心高度，m为车质量，g为重力加速度，相关参数指示见图2。

上述步骤3中，各车轮的侧偏角α_i(k)的计算方法为：

①首先车轮侧偏角计算模块根据采样时刻k实时采集到的实时方向盘转角信号 ω_sw(k)，经转向传动比换算得到前面两转向轮的转角θ₁(k)和θ₂(k)，相关参数指示见图 2；

②再根据实时纵向车速信号v_x(k)和侧向车速信号v_y(k)和横摆角速度信号计 算各车轮的侧偏角α_i(k)：

$\left(\begin{array}{c}{\alpha }_1\left(k\right)=-{\theta }_1\left(k\right)+\arctan \frac{v_y\left(k\right)+l_f\stackrel{·}{\psi }\left(k\right)}{v_x\left(k\right)-B\stackrel{·}{\psi }\left(k\right)}\\ {\alpha }_2\left(k\right)=-{\theta }_2\left(k\right)+\arctan \frac{v_t\left(k\right)+l_f\stackrel{·}{\psi }\left(k\right)}{v_x\left(k\right)+B\stackrel{·}{\psi }\left(k\right)}\\ {\alpha }_3\left(k\right)=\arctan \frac{v_y\left(k\right)-l_r\stackrel{·}{\psi }\left(k\right)}{v_x\left(k\right)-B\stackrel{·}{\psi }\left(k\right)}\\ {\alpha }_4\left(k\right)=\arctan \frac{v_y\left(k\right)-l_r\stackrel{·}{\psi }\left(k\right)}{v_x\left(k\right)+B\stackrel{·}{\psi }\left(k\right)}\end{array}\right)---\left(4\right)$

上述步骤4中，计算车轮的滑移率s_i(k)，利用如下公式计算：

$\left(\begin{array}{c}{s}_1\left(k\right)=\frac{|\sqrt{{(v_y\left(k\right)+l_f\stackrel{·}{\psi }\left(k\right))}^2+{(v_x\left(k\right)-B\stackrel{·}{\psi \left(k\right)})}^2}\cos {\alpha }_1\left(k\right)-R{\omega }_1\left(k\right)|}{\max |\sqrt{{(v_y\left(k\right)+l_f\stackrel{·}{\psi }\left(k\right))}^2+{(v_x\left(k\right)-B\stackrel{·}{\psi \left(k\right)})}^2}\cos {\alpha }_1\left(k\right),R{\omega }_1\left(k\right)|}\\ s_2\left(k\right)=\frac{|\sqrt{{(v_y\left(k\right)+l_f\stackrel{·}{\psi \left(k\right)})}^2+{(v_x\left(k\right)+B\stackrel{·}{\psi }\left(k\right))}^2}\cos {\alpha }_2\left(k\right)-R{\omega }_2\left(k\right)|}{\max |\sqrt{{(v_y\left(k\right)+l_f\stackrel{·}{\psi }\left(k\right))}^2+{(v_x\left(k\right)+B\stackrel{·}{\psi \left(k\right)})}^2}\cos {\alpha }_2\left(k\right),R{\omega }_2\left(k\right)|}\\ s_3\left(k\right)=\frac{|v_x\left(k\right)-R{\omega }_3\left(k\right)|}{\max |v_x\left(k\right),R{\omega }_3\left(k\right)|}\\ s_4\left(k\right)=\frac{|v_x\left(k\right)-R{\omega }_4\left(k\right)|}{\max |v_x\left(k\right),R{\omega }_4\left(k\right)|}\end{array}\right)---\left(5\right)$

上述步骤5中，基于修正Dugoff模型的UKF系数计算模块，利用修正Dugoff模 型算法计算UKF非线性系统的系数的方法如下：

$\left(\begin{array}{c}{a}_{1i}\left(k\right)=\frac{{3C}_s{s}_i^2\left(k\right)F_{\mathrm{zi}}^2\left(k\right)[-s_i\left(k\right)+1][1+s_i\left(k\right)]}{16\{{\left[C_s{s}_i\left(k\right)\right]}^2+{\left[C_{\alpha }\tan {\alpha }_i\left(k\right)\right]}^2\}}\\ a_{2i}\left(k\right)=\frac{-3C_s{s}_i\left(k\right)F_{\mathrm{zi}}\left(k\right)[-s_i^2\left(k\right)+s_i\left(k\right)][1+s_i\left(k\right)]}{4\sqrt{{\left[C_s{s}_i\left(k\right)\right]}^2+{\left[C_{\alpha }\tan {\alpha }_i\left(k\right)\right]}^2[1+s_i\left(k\right)]}}\\ -\frac{[1.15s_i^2\left(k\right)+1.63s_i\left(k\right)+1.27]{\left\{F_{\mathrm{zi}}\left(k\right)\right[1+s_i\left(k\right)\left]\right\}}^2}{4\{{\left[C_s{s}_i\left(k\right)\right]}^2+{\left[C_{\alpha }\tan {\alpha }_i\left(k\right)\right]}^2\}}\\ a_{3i}\left(k\right)=\frac{-C_s{s}_i\left(k\right)F_{\mathrm{zi}}\left(k\right)[1.15s_i^2\left(k\right)+1.63s_i\left(k\right)+1.27]}{\sqrt{{\left[C_s{s}_i\left(k\right)\right]}^2+{\left[C_{\alpha }\tan {\alpha }_i\left(k\right)\right]}^2}}\\ a_{4i}\left(k\right)=\frac{{3C}_s{s}_i^2\left(k\right)[-s_i\left(k\right)+1]}{4[1+s_i\left(k\right)]}\\ a_{5i}\left(k\right)=\frac{C_s{s}_i\left(k\right)[1.15s_i^2\left(k\right)+1.63s_i\left(k\right)+1.27}{[1+s_i\left(k\right)]}\end{array}\right)---\left(6\right)$

式中，C_s、C_α分别是轮胎的纵向滑移刚度和侧偏刚度，是轮胎的参数。

上述步骤6中，UKF路面峰值附着系数估计模块利用UKF观测方法得到各轮的路 面峰值附着系数具体实施过程如下：

①根据修正Dugoff轮胎模型，各轮胎的纵向力可表示为其所在路面峰值附着系数 μ_maxi和该车轮状态的函数，如下所示：

$F_{\mathrm{xi}}\left({\mu }_{\max i}\right)=\left(\begin{array}{cc}{a}_{1i}{\mu }_{\max i}^3+a_{2i}{\mu }_{\max i}^2+a_{3i}{\mu }_{\max i}& {\lambda }_i<1\\ a_{4i}{\mu }_{\max i}+a_{5i}& {\lambda }_i\ge 1\end{array}\right)---\left(7\right)$

${\lambda }_i\left(k\right)=\frac{{\mu }_{\max i}\left(k\right)·F_{\mathrm{zi}}\left(k\right)·[1+s_i\left(k\right)]}{2\sqrt{{\left[C_s{s}_i\left(k\right)\right]}^2+{\left[C_{\alpha }\tan {\alpha }_i\left(k\right)\right]}^2}}---\left(8\right)$

②以上式为基础建立无味卡尔曼滤波器(UKF)：

i)选取状态量和观测量如下：

$\left(\begin{array}{c}{x}_i={\mu }_{\max i}\\ y_i=F_{\mathrm{xi}}\end{array}\right)---\left(9\right)$

ii)将二者关系表述为非线性系统，如下：

$\left(\begin{array}{c}{x}_i(k+1)=x_i\left(k\right)+w_i\left(k\right)\\ y_i\left(k\right)=\left(\begin{array}{cc}{a}_{1i}\left(k\right)x_i^3\left(k\right)+a_{2i}\left(k\right)x_i^2\left(k\right)+a_{3i}\left(k\right)x_i\left(k\right)+v_i\left(k\right)& {\lambda }_i\left(k\right)<1\\ a_{4i}\left(k\right)x_i\left(k\right)+a_{5i}\left(k\right)+v_i\left(k\right)& {\lambda }_i\left(k\right)\ge 1\end{array}\right)\\ {\lambda }_i\left(k\right)=\frac{x_i\left(k\right)·F_{\mathrm{zi}}\left(k\right)·[1+s_i\left(k\right)]}{2\sqrt{{\left[C_s{s}_i\left(k\right)\right]}^2+{\left[C_{\alpha }\tan {\alpha }_i\left(k\right)\right]}^2}}\end{array}\right)---\left(10\right)$

iii)均值和方差的初始化：

$\left(\begin{array}{c}{\hat{x}}_i\left(0\right)=1\\ P_i\left(0\right)=0\end{array}\right)---\left(11\right)$

iv)无迹变换：

$\left(\begin{array}{c}{\hat{x}}_i^{\left(1\right)}(k-1)={\hat{x}}_i^{+}(k-1)+\sqrt{P_i^{+}(k-1)}\\ {\hat{x}}_i^2(k-1)={\hat{x}}_i^{+}(k-1)-\sqrt{P_i^{+}(k-1)}\end{array}\right)---\left(12\right)$

v)时间更新：

利用非线性方程进行非线性变换：

$\left(\begin{array}{c}{\hat{x}}_i^{\left(1\right)}\left(k\right)={\hat{x}}_i^{\left(1\right)}(k-1)\\ {\hat{x}}_i^{\left(2\right)}\left(k\right)={\hat{x}}_i^{\left(2\right)}(k-1)\end{array}\right)---\left(13\right)$

加权得到的预测的状态向量：

${\hat{x}}_i^{-}\left(k\right)={\hat{x}}_i^{\left(1\right)}\left(k\right)+{\hat{x}}_i^{\left(2\right)}\left(k\right)---\left(14\right)$

加权得到的预测的协方差矩阵：

$P_i^{-}\left(k\right)=\frac{1}{2}\{{[{\hat{x}}_i^{\left(1\right)}\left(k\right)-{\hat{x}}_i^{-}\left(k\right)]}^2+{[{\hat{x}}_i^{\left(2\right)}\left(k\right)-{\hat{x}}_i^{-}\left(k\right)]}^2\}+Q_i(k-1)---\left(15\right)$

利用非线性方程进行非线性变换：

先计算${\lambda }_i\left(k\right)=\frac{{\hat{x}}_i^{-}\left(k\right)·F_{\mathrm{zi}}\left(k\right)·[1+s_i\left(k\right)]}{2\sqrt{{\left[C_s{s}_i\left(k\right)\right]}^2+{\left[C_{\alpha }\tan {\alpha }_i\left(k\right)\right]}^2}}$

根据λ_i(k)计算

$\left(\begin{array}{c}{\hat{y}}_i^{\left(1\right)}\left(k\right)=\left(\begin{array}{cc}{a}_{1i}\left(k\right){\left[x_i^{\left(1\right)}\left(k\right)\right]}^3+a_{2i}\left(k\right){\left[x_i^{\left(1\right)}\left(k\right)\right]}^2+a_{3i}\left(k\right)x_i^{\left(1\right)}\left(k\right)+v_i\left(k\right)& {\lambda }_i\left(k\right)<1\\ a_{4i}\left(k\right)x_i^{\left(1\right)}{\left(k\right)}^3+a_{5i}\left(k\right)+v_i\left(k\right)& {\lambda }_i\left(k\right)\ge 1\end{array}\right)\\ {\hat{y}}_i^{\left(2\right)}\left(k\right)=\left(\begin{array}{cc}{a}_{1i}\left(k\right){\left[x_i^{\left(2\right)}\left(k\right)\right]}^3+a_{2i}\left(k\right){\left[x_i^{\left(2\right)}\left(k\right)\right]}^2+a_{3i}\left(k\right)x_i^{\left(2\right)}\left(k\right)+v_i\left(k\right)& {\lambda }_i\left(k\right)<1\\ a_{4i}\left(k\right)x_i^{\left(2\right)}{\left(k\right)}^3+a_{5i}\left(k\right)+v_i\left(k\right)& {\lambda }_i\left(k\right)\ge 1\end{array}\right)\end{array}\right)---\left(16\right)$

加权得到系统的预测值：

${\hat{y}}_i\left(k\right)=\frac{1}{2}[{\hat{y}}_i^{\left(1\right)}\left(k\right)+{\hat{y}}_i^{\left(2\right)}\left(k\right)]---\left(17\right)$

vi)量测更新：

后验估计的协方差矩阵：

$P_{y-i}\left(k\right)=\frac{1}{2}\{{[{\hat{y}}_i^{\left(1\right)}\left(k\right)-{\hat{y}}_i\left(k\right)]}^2+{[{\hat{y}}_i^{\left(2\right)}\left(k\right)-{\hat{y}}_i\left(k\right)]}^2\}+R_i\left(k\right)---\left(18\right)$

先验估计的协方差矩阵：

$P_{\mathrm{xy}-i}\left(k\right)=\frac{1}{2}\left\{\right[{\hat{y}}_i^{\left(1\right)}\left(k\right)-{\hat{y}}_i\left(k\right)\left]\right[{\hat{x}}_i^{\left(1\right)}\left(k\right)-{\hat{x}}_i^{-}\left(k\right)]+[{\hat{y}}_i^{\left(2\right)}\left(k\right)-{\hat{y}}_i\left(k\right)\left]\right[{\hat{x}}_i^{\left(2\right)}\left(k\right)-{\hat{x}}_i^{-}\left(k\right)\left]\right\}---\left(19\right)$

滤波增益矩阵：

$K_i\left(k\right)=P_{\mathrm{xy}-i}\left(k\right)P_{y-i}^{-1}\left(k\right)---\left(20\right)$

状态更新后的滤波值：

${\hat{x}}_i^{+}\left(k\right)={\hat{x}}_i^{-}\left(k\right)+K_i\left(k\right)·(y_i\left(k\right)-{\hat{y}}_i\left(k\right))---\left(21\right)$

状态更新后的后验方差矩阵：

$P_i^{+}\left(k\right)=P_i^{-}\left(k\right)-K_i^2\left(k\right)·P_{y-i}\left(k\right)---\left(22\right)$

到此为止，UKF滤波器建立完毕。滤波器输出的滤波值即为各轮处的附着系 数估计值