CN103407451A - Method for estimating longitudinal adhesion coefficient of road - Google Patents

Abstract

The invention discloses a method for estimating the longitudinal adhesion coefficient of a road. The method comprises the following steps: preliminarily estimating a longitudinal road adhesion coefficient in real time on the road surface of a flat expressway by using a recursive least squares (RLS) method with a forgetting factor specific to a front-wheel steering four-wheel automobile on the basis of an overall longitudinal dynamics model and a simplified magic formula tire model; further filtering signal noise by using an extended Kalman filter (EKF) algorithm by taking the estimated longitudinal road adhesion coefficient and a tire model parameter as expansion states, and realizing self-adaptive adjustment of a tire model coefficient; finally, acquiring accurate and robust longitudinal road adhesion coefficient estimation in real time. The method can adapt to working conditions of high slip rate and low slip rate of flat expressways at the same time.

Description

The vertical adhesion value method of estimation of a kind of road

Technical field

The present invention relates to the vertical adhesion value method of estimation of a kind of road, its purpose is by the modeling of car load longitudinal dynamics and tire model, utilize recurrence least square and EKF method, under smooth express highway pavement environment, realization under high slip rate and low slip rate operating mode to vertical road-adhesion coefficient accurately, reliable estimation, the automobile longitudinal active safety system can be regulated control policy according to the estimated vertical road-adhesion coefficient gone out, improve vehicle safety, has precision high, cost is low, real-time and comformability wait well remarkable advantage, belonging to automobile active safety measures and control field.

Background technology

Along with socioeconomic development, the traffic safety problem becomes increasingly conspicuous, and has become a global difficult problem.The whole world all can cause a large amount of personal casualty and property damage because of traffic accident every year, and countries in the world are all in the generation of making great efforts to reduce traffic accident.In recent years, the automobile active safety technology has obtained development rapidly.The automobile active safety technology can prevent trouble before it happens, and initiatively avoids the generation of accident, has become one of topmost developing direction of modern automobile.Common active safety technology mainly comprises anti-skid brake system (ABS) at present, vehicle electric stability program (ESP), anti-slip regulation (TCS), automatically controlled driving skid control system (ASR), four-wheel steering stabilizing control system (4WS) etc.Can the quality of these system effects depends on to a great extent " road self adaptation ", if namely can estimate in real time road-adhesion coefficient, system just can be regulated control policy according to current road conditions, improves vehicle safety.Vertically road-adhesion coefficient is a key factor that affects automobile brake perfrmance, it is anti-collision warning/active collision avoidance (CW/CA), adaptive cruise (ACC), a requisite important parameter in the automobile longitudinal safety assisting systems such as anti-skid brake system (ABS), estimate that in real time vertical road-adhesion coefficient is the important technology prerequisite of vertical security algorithm of the realistic requirement of design, driving safety and the stability of its precision direct relation automobile, be vertical active safety control system can effectively be operated in depend on to a great extent road-adhesion coefficient can be by real time, measure accurately or estimate.

At present, in the automobile active safety field, road-adhesion coefficient mainly is divided into direct measurement and indirectly estimates two classes, direct measuring method is to utilize the sensor direct-detection road surfaces such as light, sound, image, radar, measure the larger factor of some road pavement adhesion value impacts, and according to previous experiences, predict the size of current road-adhesion coefficient, but these methods all need additionally to install additional sensor, and the sensor cost is all higher, be difficult to realize large-scale business application; Secondly need to carry out a large amount of test training, accuracy of identification depends on experience to a great extent, is difficult to accurately estimate the adhesion value of the road conditions that there is no test and trained.Indirect estimation methods is to carry out kinematics or Dynamic Modeling by the operational process to automobile, in conjunction with tire model, using relevant onboard sensor cheaply (as wheel speed sensors, gyroscope, accelerometer and GPS etc.) information as observation information, and then utilize suitable filtering algorithm for estimating to realize the estimation to road-adhesion coefficient.Existing road-adhesion coefficient indirect estimation methods comprises based on the vehicle lateral dynamics with based on two kinds of the researchs of longitudinal dynamics, wherein, the former is mainly used in the motor turning motion, and for the automobile longitudinal safety assisting system, the dynamic (dynamical) adhesion value of longitudinal direction of car that is based on of comparatively paying close attention to is estimated.On the road surface of different adhesion valuies, have different slip rates and normalization method tractive force relation, this is the basis that utilizes longitudinal direction of car dynam estimation adhesion value.Based on this relation, just can be by estimating the slippage slope, utilize sorting algorithm to estimate adhesion value, this is that the outer scholar of Present Domestic uses maximum methods, and the method thinks that the normalization method tractive force in little slip rate scope is directly proportional to slip rate, and its proportionate relationship is called the slippage slope, and road-adhesion coefficient is directly proportional to the slippage slope, therefore, by estimating the slippage slope, just can estimate the adhesion value on tire road surface.But this method can only be applied to low slip rate operating mode, and under high slip rate operating mode, the hypothesis that the normalization method tractive force is directly proportional to the slippage slope is no longer set up, and therefore can't utilize the method to carry out the road-adhesion coefficient estimation; Simultaneously on same road surface the slippage slope also there is some difference with the impact of other factors, brought certain difficulty for the identification of adhesion value.Another kind of method commonly used is to adopt multi-form non-linear formula in real time according to measure data fitting normalization method tractive force and slip rate curve, thereby according to the peak estimation adhesion value that calculates institute's matched curve, but this method need to just can make to estimate that result is more close to actual value under than giant tyre slip rate condition.Simultaneously, the method adopts nonlinear function to carry out match, and calculated amount is large, and real-time is poor, no matter and adopt any nonlinear function, the tire characteristics that all can't all coincide under all operating modes, error is larger.Visible, current on adhesion value is estimated existing subject matter be to be difficult under high slip rate and low slip rate operating mode, all to obtain estimation effect preferably simultaneously, limited its application on the automobile longitudinal active safety system.

Summary of the invention

For realizing under high slip rate and low slip rate operating mode, can carrying out accurately road-adhesion coefficient, reliably estimate, the present invention proposes the vertical adhesion value method of estimation of a kind of road.The present invention is directed to existing method is difficult under high slip rate and low slip rate operating mode, all obtain estimation effect preferably simultaneously, towards the front-wheel steering four-wheel automobile travelled on smooth express highway, adopt the simplification magic formula tire model of more realistic vehicle operating process, in conjunction with vertical Full Vehicle Dynamics model, utilize recurrent least square method to carry out the road-adhesion coefficient estimation, take full advantage of simultaneously the observed quantity of low cost vehicle-mounted sensor information as system, by the EKF method, realize the self adaptation adjustment of tire model coefficient, and remove vertical adhesion value estimated result institute Noise, export the vertical adhesion value estimated result of final road, has precision high, cost is low, real-time and comformability wait well remarkable advantage.

The vertical adhesion value method of estimation of a kind of road, it is characterized in that: on smooth express highway pavement, for the front-wheel steering four-wheel automobile, based on car load Longitudinal Dynamic Model and simplification magic formula tire model, utilization is with recurrence least square (the Recursive least squares of forgetting factor, RLS) method goes out the vertical adhesion value of road in real time according to a preliminary estimate, further using the estimated vertical adhesion value of road gone out and tire model parameter as extended mode, utilize EKF (Extended Kalman Filter, EKF) algorithm, the filtered signal noise, and realize the self adaptation adjustment of tire model coefficient, final Real-time Obtaining is accurate, the vertical adhesion value of the road of robust is estimated, high slip rate and low slip rate operating mode on can the smooth express highway of simultaneous adaptation,

Concrete steps comprise:

1) calculate the straight skidding rate

Use s _j(j=f, r) means the longitudinal direction of car slip rate, namely can be divided into again front wheel spindle straight skidding rate s _fWith hind axle straight skidding rate s _r, subscript j gets f or r, and f or r mean respectively front or rear wheel shaft, s _jMethod of calculating is:

In formula (1), R _eMean tire radius; v _TfAnd v _TrMean respectively on forward and backward wheel shaft the speed along the tire direction, all with vehicle absolute velocitye v, replace in the present invention, obtain by the GPS measurement; ω _fThe spin velocity on front wheel spindle is converted in the spin velocity equivalence that means two wheels on front wheel spindle; ω _rMean that on hind axle, the spin velocity on hind axle, ω are converted in two rotation of wheel cireular frequency equivalences _fAnd ω _rCan unify to be designated as ω _j(j=f, r) and

${\mathrm{\ω}}_f=\frac{1}{2}({\mathrm{\ω}}_{\mathrm{fR}}+{\mathrm{\ω}}_{\mathrm{fL}})$

(2)

${\mathrm{\ω}}_r=\frac{1}{2}({\mathrm{\ω}}_{\mathrm{rR}}+{\mathrm{\ω}}_{\mathrm{rL}})$

In formula (2), ω _FL, ω _FR, ω _RLAnd ω _RRThe spin velocity that means respectively the near front wheel, off front wheel, left rear wheel and off hind wheel, obtain by utilizing four wheel speed sensors to measure;

2) calculate the car load longitudinal force

With

Mean the normalization method tractive force, namely can be divided into again front wheel spindle normalization method tractive force With hind axle normalization method tractive force

Subscript j gets f or r, and f or r mean respectively front or rear wheel shaft:

In formula (3), F _XfAnd F _XrMean respectively longitudinal force on forward and backward wheel shaft, F _XfAnd F _XrCan unify to be designated as F _Xj(j=f, r), F _ZfAnd F _ZrVertical load before or after meaning respectively to be assigned on wheel shaft, can unify to be designated as F _Zj(j=f, r) and can be calculated as follows:

$F_{\mathrm{zf}}=\frac{\mathrm{mgb}}{(a+b)},F_{\mathrm{zr}}=\frac{\mathrm{mga}}{(a+b)}---\left(4\right)$

In formula (4), g means acceleration due to gravity, and m means vehicle mass, and a is the distance of vehicle front wheel shaft center to barycenter, and b is the distance of automobile back wheel wheel shaft center to barycenter;

The car load Longitudinal Dynamic Model is expressed as:

F _x=F _xf+F _xr=ma _x+D _av ²+C _rollmg (5)

In formula (5), a _xMean the real-time longitudinal acceleration of vehicle, by accelerometer measures, obtained D _aMean the air resistance constant, C _RollMean coefficient of rolling resistance, v means vehicle absolute velocitye, F _xMean the car load longitudinal force;

The present invention adopts and simplifies the magic formula tire model, and subscript j gets f or r, and f or r mean respectively front or rear wheel shaft:

In formula (6), μ means the Real-time Road adhesion value, and B, C mean model coefficient, supposes that the coefficient of road adhesion situation of each wheel is identical, and the normalization method tractive force of forward and backward wheel shaft is respectively:

By formula (5), (7), (8) Shi Kede

F _x=F _xf+F _xr=μ{F _zfsin[Carctan(Bs _f)]+F _zrsin[Carctan(Bs _r)]} (9)

3) based on vertical road-adhesion coefficient of method of least square according to a preliminary estimate

Formula (9) is expressed as to parameter criterion of identification form:

In formula (10), k means the discrete moment, y (k)=F _xExpression system output, can be calculated and be learnt by (5) formula; θ (k)=μ _RLSMean the solve for parameter vector, wherein μ _RLSExpression by the RLS method road-adhesion coefficient according to a preliminary estimate;

Mean the input regression vector, in the present invention, superscript ' expression is to matrix transpose; E (k) means identification error; Utilize with the RLS algorithm of forgetting factor and determine that in real time unique unknown quantity---the estimating step of road-adhesion coefficient is as follows:

(1) computing system output variable y (k), and calculate the iteration vector

(2) calculate identification error e (k);

(3) calculated gains vector N (k);

Wherein,

Wherein, parameter lambda is forgetting factor;

(4) calculate solve for parameter vector theta (k);

θ(k)=θ(k-1)+N(k)e(k)

So far, can go out according to a preliminary estimate in real time road-adhesion coefficient μ _RLS

4) vertical road-adhesion coefficient of extension-based Kalman filtering is estimated

By the estimated μ gone out of RLS method _RLSAs observed quantity, B, C and μ, as extended mode, utilize longitudinal dynamics equation (5) to set up the EKF model, and the equation of state after discretization and observational equation are suc as formula (11):

X(k)=f(X(k-1))+W(k-1)

(11)

Z(k)=h(X(k))+V(k)

In formula (11), state of the system vector X=[x ₁, x ₂, x ₃, x ₄] ', and x ₁=v, x ₂=μ, x ₃=B, x ₄=C; W means the system Gaussian white noise vector of zero-mean, W=[w ₁w ₂w ₃w ₄] ', be w wherein ₁, w ₂, w ₃And w ₄Mean respectively four systems Gaussian white noise component, the system noise covariance matrix Q (k-1) that W is corresponding is: $Q(k-1)=\left[\begin{array}{cccc}{{\mathrm{\&sigma;}}_{{\mathrm{<figure-callout\; id="1"\; label="w"\; filenames="130903154114.png"\; state="\{\{state\}\}">w</figure-callout>}}_1}}^2& 0& 0& 0\\ 0& {{\mathrm{\&sigma;}}_{w_2}}^2& 0& 0\\ 0& 0& {{\mathrm{\&sigma;}}_{{\mathrm{<figure-callout\; id="3"\; label="w"\; filenames="130903154105.png"\; state="\{\{state\}\}">w</figure-callout>}}_3}}^2& 0\\ 0& 0& 0& {{\mathrm{\&sigma;}}_{w_4}}^2\end{array}\right],$ Wherein And

Mean respectively system Gaussian white noise w ₁, w ₂, w ₃And w ₄Corresponding variance; Observation vector Z=[a _x-m, v _m, μ _RLS] ', a _x-mMean the longitudinal direction of car acceleration/accel that the acceleration pick-up measurement obtains; v _mThe car speed obtained is measured in expression by GPS; μ _RLSExpression estimates that by recurrence least square the road-adhesion coefficient of gained is worth according to a preliminary estimate; V means the mutual incoherent zero-mean observation white noise vector with W,

Expression by acceleration pick-up measure the longitudinal direction of car acceleration/accel obtained observation noise and

That average is 0, variance is

Gaussian white noise;

Expression by GPS measure the car speed obtained observation noise and

That average is 0, variance is

Gaussian white noise;

Mean by above-mentioned method of least square estimate acquisition road-adhesion coefficient observation noise and That average is 0, variance is

Gaussian white noise; The observation noise variance battle array R that V is corresponding can be expressed as $R=\left[\begin{array}{ccc}{{\mathrm{\&sigma;}}_{a_{x-m}}}^2& 0& 0\\ 0& {{\mathrm{\&sigma;}}_{v_m}}^2& 0\\ 0& 0& {{\mathrm{\&sigma;}}_{{\mathrm{\&mu;}}_{\mathrm{RLS}}}}^2\end{array}\right];$ F () and h () mean nonlinear system function vectorial sum observation function vector:

$f\left(X(k-1)\right)=\left[\begin{array}{c}{f}_1\left(X(k-1)\right)\\ f_2\left(X(k-1)\right)\\ f_3\left(X(k-1)\right)\\ f_4\left(X(k-1)\right)\end{array}\right],$ $h\left(X\left(k\right)\right)=\left[\begin{array}{c}{h}_1\left(X\left(k\right)\right)\\ h_2\left(X\left(k\right)\right)\\ h_3\left(X\left(k\right)\right)\end{array}\right],$ Wherein

$f_1\left(X(k-1)\right)=v(k-1)-\frac{T\{D_a{\left[v(k-1)\right]}^2+C_{\mathrm{roll}}\mathrm{mg}\}}{m}+$

$\frac{\mathrm{T\&mu;}(k-1)\left\{F_{\mathrm{zf}}\sin \right[C(k-1)\arctan \left(B(k-1)s_f\right)]+F_{\mathrm{zr}}\sin [C(k-1)\arctan \left(B(k-1)S_r\right)\left]\right\}}{m}$

$,f_2\left(X(k-1)\right)=\mathrm{\&mu;}(k-1),f_3\left(X(k-1)\right)=B(k-1),f_4\left(X(k-1)\right)=C(k-1)$

$h_1\left(X\left(k\right)\right)=\frac{\mathrm{\&mu;}\left\{F_{\mathrm{zf}}\sin \right[C\arctan \left({\mathrm{Bs}}_f\right)]+F_{\mathrm{zr}}\sin [C\arctan \left({\mathrm{Bs}}_r\right)\left]\right\}-D_a{\left[v(k-1)\right]}^2-C_{\mathrm{roll}}\mathrm{mg}}{m}$

$,h_3\left(X\left(k\right)\right)=v,h_3\left(X\left(k\right)\right)=\mathrm{\&mu;}$

T means the discrete cycle, and its representative value is 10 milliseconds, 20 milliseconds, 50 milliseconds or 100 milliseconds;

For the described equation of state of formula (11) and observational equation, use the EKF theory, Criterion filtering recursive process, this recursive process comprise that the time upgrades and measurement is upgraded:

Time upgrades:

State one-step prediction equation: $\hat{X}(k,k-1)=f\left(X(k-1)\right)$

One-step prediction error covariance matrix P (k, k-1):

P(k，k-1)=A(k，k-1)P(k-1)A′(k，k-1)+Q(k-1)

Measure and upgrade:

Filtering gain matrix K (k): K (k)=P (k, k-1) H ' is [H (k) P (k, k-1) H ' (k)+R (k)] (k) ^-1

State estimation: $\hat{X}\left(k\right)=\hat{X}(k,k-1)+K\left(k\right)[Z\left(k\right)-h[\hat{X}(k,k-1)]]$

Estimation error variance battle array P (k): P (k)=[I-K (k) H (k)] P (k, k-1) and I are 4 * 4 unit matrix;

Wherein, A is that state of the system functional vector f asks the Jacobi matrix of partial derivative to state vector X, and H is that observation function vector h asks the Jacobi matrix of partial derivative to state vector X, i.e. the capable j column element of the i A of matrix A and H _{[i, j]}(i=1,2,3,4 j=1,2,3,4) and H _{[i, j]}(i=1,2,3 j=1,2,3,4) can try to achieve by following formula respectively:

$A_{[i,j]}=\frac{{\&PartialD;f}_i}{{\&PartialD;x}_j}\left(\hat{X}(k,k-1)\right),(i=\mathrm{1,2,3,4},j=\mathrm{1,2,3,4})$

$H_{[i,j]}=\frac{{\&PartialD;h}_i}{{\&PartialD;x}_j}\left(\hat{X}(k,k-1)\right),(i=1,\mathrm{2,3},j=\mathrm{1,2,3,4})$

Particularly, the value of each element of a matrix is as follows:

A _[1，1]=1+T[(-2D _av)/m]

A _[1，2]=T{F _zfsin[Carctan(Bs _f)]+F _zrsin[Carctan(Bs _r)]}/m

$A_{\left[\mathrm{1,3}\right]}=\mathrm{T\&mu;}[\frac{F_{\mathrm{zf}}{\mathrm{Cs}}_f\cos \left(C\arctan \left({\mathrm{Bs}}_f\right)\right)}{1+{\left({\mathrm{Bs}}_f\right)}^2}+\frac{F_{\mathrm{zr}}{\mathrm{Cs}}_r\cos \left(C\arctan \left({\mathrm{Bs}}_r\right)\right)}{1+{\left({\mathrm{Bs}}_r\right)}^2}]/m$

A _[1，4]=Tμ[F _zf(arctan(Bs _f)cos(Carctan(Bs _f)))+F _zr(arctan(Bs _r)cos(Carctan(Bs _r)))]/m

A _[2，1]=A _[2，3]=A _[2，4]=A _[3，1]=A _[3，2]=A _[3，4]=A _[4，1]=A _[4，2]=A _[4，3]=0

A _[2，2]=A _[3，3]A _[4，4]=1

H _[1，1]=(-2D _av)/m H _[1，2]={F _zfsin[Carctan(Bs _f)]+F _zrsin[Carctan(Bs _r)]}/m

$H_{\left[\mathrm{1,3}\right]}=\mathrm{\&mu;}[\frac{F_{\mathrm{zf}}{\mathrm{Cs}}_f\cos \left(C\arctan \left({\mathrm{Bs}}_f\right)\right)}{1+{\left({\mathrm{Bs}}_f\right)}^2}+\frac{F_{\mathrm{zr}}{\mathrm{Cs}}_r\cos \left(C\arctan \left({\mathrm{Bs}}_r\right)\right)}{1+{\left({\mathrm{Bs}}_r\right)}^2}]/m$

H[ _1，4]=μ[F _zf(arctan(Bs _f)cos(Carctan(Bs _f)))+F _zr(arctan(Bs _r)cos(Carctan(Bs _r)))]/m

H _[2，2]＝H _[2，3]＝H _[2，4]=H _[3，1]=H _[3，3]=H _[3，4]=0

H _[2，1]=H _[3，2]＝1

The filtered μ value that the above-mentioned EKF filtering recursion of usining is exported is as final road-adhesion coefficient estimated result.

Advantage of the present invention and remarkable result:

1. the present invention proposes the vertical adhesion value method of estimation of a kind of road, can be used for smooth express highway pavement get on the car vertical active safety control to the measurement of the vertical adhesion value of road with estimate needs, have that precision is high, cost is low, real-time and an advantage such as comformability is good.

2. method of the present invention is can't the high slip rate of simultaneous adaptation and low slip rate operating mode for orthodox method, on the basis of car load Longitudinal Dynamic Model and simplified nonlinear magic formula tire model, propose, no matter at high slip rate, still under low slip rate operating mode, can obtain road-adhesion coefficient information comparatively accurately.

3. the vertical adhesion value method of estimation of road that proposes of the present invention has good comformability for the sudden change of road-adhesion coefficient, and response time is short, can meet the get on the car requirement of vertical active safety control of smooth express highway pavement.

The accompanying drawing explanation

Fig. 1 is method flow block diagram proposed by the invention;

Fig. 2 is the relation of simplifying magic formula tire model normalization method tractive force and slip rate;

Fig. 3 is slippage rate curve under the low single adhesion value operating mode of slip rate;

Fig. 4 is the adhesion value estimated result (in figure, the adhesion value estimated result is with black solid line representative, and Carsim output true value is with black dotted line representative) of the inventive method output under the operating mode corresponding with Fig. 3;

Fig. 5 is slippage rate curve under low slip rate adhesion value sudden load;

Fig. 6 is the adhesion value estimated result (in figure, the adhesion value estimated result is with black solid line representative, and Carsim output true value is with black dotted line representative) of the inventive method output under the operating mode corresponding with Fig. 5;

Fig. 7 is slippage rate curve under the single adhesion value operating mode of high slip rate;

Fig. 8 is the adhesion value estimated result (in figure, the adhesion value estimated result is with black solid line representative, and Carsim output true value is with black dotted line representative) of the inventive method output under the operating mode corresponding with Fig. 7;

Fig. 9 is the adhesion value estimated result (in figure, the adhesion value estimated result is with black solid line representative, and Carsim output true value is with black dotted line representative) of traditional slippage gradient method output under the operating mode corresponding with Fig. 7;

Figure 10 is slippage rate curve under high slip rate adhesion value sudden load;

Figure 11 is the adhesion value estimated result (in figure, the adhesion value estimated result is with black solid line representative, and Carsim output true value is with black dotted line representative) of the inventive method output under the operating mode corresponding with Figure 10;

Figure 12 is the adhesion value estimated result (in figure, the adhesion value estimated result is with black solid line representative, and Carsim output true value is with black dotted line representative) of traditional slippage gradient method output under the operating mode corresponding with Figure 10;

The specific embodiment

Embodiment 1

Along with socioeconomic development, the traffic safety problem becomes increasingly conspicuous, and has become a global difficult problem.The whole world all can cause a large amount of personal casualty and property damage because of traffic accident every year, and countries in the world are all in the generation of making great efforts to reduce traffic accident.In recent years, the automobile active safety technology has obtained development rapidly.The automobile active safety technology can prevent trouble before it happens, and initiatively avoids the generation of accident, has become one of topmost developing direction of modern automobile.Common active safety technology mainly comprises anti-lock braking system in automobiles (ABS) at present, vehicle electric stability program (ESP), anti-slip regulation (TCS), automatically controlled driving skid control system (ASR), four-wheel steering stabilizing control system (4WS) etc.Can the quality of these system effects depends on to a great extent " road self adaptation ", if namely can estimate in real time road-adhesion coefficient, system just can be regulated control policy according to current road conditions, improves vehicle safety.Vertically road-adhesion coefficient is a key factor that affects automobile brake perfrmance, it is anti-collision warning/active collision avoidance (CW/CA), adaptive cruise (ACC), a requisite important parameter in the automobile longitudinal safety assisting systems such as anti-skid brake system (ABS), estimate that in real time vertical road-adhesion coefficient is the important technology prerequisite of vertical security algorithm of the realistic requirement of design, driving safety and the stability of its precision direct relation automobile, be vertical active safety control system can effectively be operated in depend on to a great extent road-adhesion coefficient can be by real time, measure accurately or estimate.

At present, in the automobile active safety field, road-adhesion coefficient mainly is divided into direct measurement and indirectly estimates two classes, direct measuring method is to utilize the sensor direct-detection road surfaces such as light, sound, image, radar, measure the larger factor of some road pavement adhesion value impacts, and according to previous experiences, predict the size of current road-adhesion coefficient, but these methods all need additionally to install additional sensor, and the sensor cost is all higher, be difficult to realize large-scale business application; Secondly need to carry out a large amount of test training, accuracy of identification depends on experience to a great extent, is difficult to accurately estimate the adhesion value of the road conditions that there is no test and trained.Indirect estimation methods is to carry out kinematics or Dynamic Modeling by the operational process to automobile, in conjunction with tire model, using relevant onboard sensor cheaply (as wheel speed sensors, gyroscope, accelerometer and GPS etc.) information as observation information, and then utilize suitable filtering algorithm for estimating to realize the estimation to road-adhesion coefficient.Existing indirect method comprises based on the vehicle lateral dynamics with based on two kinds of the researchs of longitudinal dynamics, wherein, the former is mainly used in the motor turning motion, and for the automobile longitudinal safety assisting system, the dynamic (dynamical) adhesion value of longitudinal direction of car that is based on of comparatively paying close attention to is estimated.On the road surface of different adhesion valuies, have different slip rates and normalization method tractive force relation, this is the basis that utilizes longitudinal direction of car dynam estimation adhesion value.Based on this relation, just can utilize sorting algorithm to estimate adhesion value by estimating the slippage slope, this is that the outer scholar of Present Domestic uses maximum methods.The method thinks that the normalization method tractive force in little slip rate scope is directly proportional to slip rate, and its proportionate relationship is called the slippage slope, and road-adhesion coefficient is directly proportional to the slippage slope, therefore, by estimating the slippage slope, just can estimate the adhesion value on tire road surface.But this method can only be applied to low slip rate operating mode, and under high slip rate operating mode, the hypothesis that the normalization method tractive force is directly proportional to the slippage slope is no longer set up, and therefore can't carry out the road-adhesion coefficient estimation.Simultaneously on same road surface the slippage slope also there is some difference with the impact of other factors, brought certain difficulty for the identification of adhesion value.Another kind of method commonly used is to adopt power function, exponential function, the non-linear formula that logarithmic function etc. are multi-form and combination thereof come in real time according to measure data fitting normalization method tractive force and slip rate curve, thereby according to the peak estimation road-adhesion coefficient that calculates institute's matched curve, but this method need to just can make to estimate that result is more close to actual value under than giant tyre slip rate condition, simultaneously, the method is based on the non-linear tire model supposed to the match of normalization method tractive force and slip rate curve, if select tire model too complicated, real-time is difficult to guarantee, if the tire model of choosing is too simple, can't guarantee accuracy again.Visible, current on adhesion value is estimated existing subject matter be to be difficult under high slip rate and low slip rate operating mode, all to obtain estimation effect preferably simultaneously, limited its application on the automobile longitudinal active safety system.

For realizing can carrying out accurately road-adhesion coefficient under the high slip rate of smooth express highway pavement and low slip rate operating mode, reliable estimation, the present invention proposes the vertical adhesion value method of estimation of a kind of road.This method is based on based on recurrence least square (Recursive least squares, RLS) and EKF (Extended Kalman Filter, EKF).The present invention is directed to existing method and can't be applicable to simultaneously high slip rate and low slip rate operating mode, adopt the simplification magic formula tire model of more realistic vehicle operating process, in conjunction with vertical Full Vehicle Dynamics model, utilize recurrent least square method to carry out the road-adhesion coefficient estimation, take full advantage of simultaneously the observed quantity of low cost vehicle-mounted sensor information as system, by the EKF method, realize the self adaptation adjustment of tire model coefficient, and remove vertical road-adhesion coefficient estimated result institute Noise, export the vertical adhesion value estimated result of final road, has precision high, cost is low, real-time and comformability wait well remarkable advantage, the estimated road-adhesion coefficient gone out can be used for the automobile longitudinal active safety control travelled on smooth express highway pavement.Concrete thought of the present invention is as follows:

Recurrence least square is the iterative algorithm to unknown vector, and the minimum variance of model error of take is target, for each sampling period, uses existing sampled data by iterating, to calculate unknown vector.Kalman filter is to take the optimal State Estimation filter of Minimum Mean Square Error as criterion, it does not need to store observed reading in the past, only, according to the estimated valve of current observed value and previous moment, utilize computing machine to carry out recursion calculating, just can realize the estimation to live signal.Recurrence least square and Kalman filtering all have the characteristics that memory data output is little, algorithm is easy.

For adapt under the high slip rate of smooth express highway pavement and low slip rate operating mode automobile longitudinal active safety control to the measurement of vertical road-adhesion coefficient with estimate requirement, at first automobile and tire are carried out to suitable Dynamic Modeling.For application of the present invention, the present invention (should have the widest situation at present for the four wheeler of the front-wheel steering on common road traffic environment that travels, the car of exemplary such as front-wheel steering), ignore side direction, yaw, pitching, inclination and upper and lower bounce motion, ignore left and right wheels difference, because scope of the present invention is flat road surface, therefore can ignore road grade, the coefficient of road adhesion situation of supposing each wheel is identical, and the car load Longitudinal Dynamic Model is expressed as:

ma _x=F _x-D _av ²-C _rollmg (1)

In formula (1), m means vehicular gross combined weight, a _xMean the real-time longitudinal acceleration of vehicle, by acceleration pick-up, obtained F _xMean the car load longitudinal force, D _aMean the air resistance constant, C _RollMean coefficient of rolling resistance, v means vehicle absolute velocitye, is obtained by GPS, and g means acceleration due to gravity.

Formula (1) deformable is:

F _x=F _xf+F _xr=ma _x+D _av ²+C _rollmg (2)

In formula (2), F _XfAnd F _XrMean respectively longitudinal force on forward and backward wheel shaft, F _XfAnd F _XrCan unify to be designated as F _Xj(j=f, r).

The magic formula tire model is the highest Empirical tire model of fitting precision of generally acknowledging, that high slip rate or low slip rate operating mode are all applicable equally, but it is the nonlinear function of the complexity that combined by trigonometric function, and in model, need the unknown factor of determining also more, calculated amount is larger, is unwell to real-time use.Therefore, this paper adopt its simplified model [but list of references: Bian Mingyuan, for the vertical simplification tire model [J] of road-adhesion coefficient assessment. Chongqing University of Technology's journal (natural science), 2012,26 (1): 1-5.].With

Mean the normalization method tractive force, namely can be divided into again front wheel spindle normalization method tractive force

With hind axle normalization method tractive force

Subscript j gets f or r, and f or r mean respectively front or rear wheel shaft:

(3)

In formula (3), F _XfAnd F _XrMean respectively longitudinal force on forward and backward wheel shaft, F _XfAnd F _XrCan unify to be designated as F _Xj(j=f, r), F _ZfAnd F _ZrVertical load before or after meaning respectively to be assigned on wheel shaft, can unify to be designated as F _Zj(j=f, r) and can be calculated as follows:

$\begin{array}{cc}{F}_{\mathrm{zf}}=\frac{\mathrm{mgb}}{(a+b)},& F_{\mathrm{zr}}=\frac{\mathrm{mga}}{(a+b)}\end{array}$

(4)

In formula (4), a is the distance of vehicle front wheel shaft center to barycenter, and b is the distance of automobile back wheel wheel shaft center to barycenter.

Use s _j(j=f, r) means the longitudinal direction of car slip rate, namely can be divided into again front wheel spindle straight skidding rate s _fWith hind axle straight skidding rate s _r, subscript j gets f or r, and f or r mean respectively front or rear wheel shaft, s _jMethod of calculating is:

(5)

In formula (5), R _eMean the wheel tyre radius; ν _TfAnd ν _TrMean respectively on forward and backward wheel shaft the speed along the tire direction, v _TfAnd v _TrCan unify to be designated as v _Tj(j=f, r), because scene that the present invention was suitable for is express highway, the road curvature of express highway is less, Vehicle Driving Cycle is on express highway the time, its yaw velocity is also less, and the present invention's method used is based on longitudinal dynamics to carry out vertical road-adhesion coefficient estimation, so v in the present invention _TjCan be approximately equal to vehicle absolute velocitye v, obtain by the GPS measurement; ω _fThe spin velocity on front wheel spindle is converted in the spin velocity equivalence that means two wheels on front wheel spindle; ω _rMean that on hind axle, the spin velocity on hind axle, ω are converted in two rotation of wheel cireular frequency equivalences _fAnd ω _rCan unify to be designated as ω _j(j=f, r) and

${\mathrm{\&omega;}}_f=\frac{1}{2}({\mathrm{\&omega;}}_{\mathrm{fR}}+{\mathrm{\&omega;}}_{\mathrm{fL}})$

(6)

${\mathrm{\&omega;}}_r=\frac{1}{2}({\mathrm{\&omega;}}_{\mathrm{rR}}+{\mathrm{\&omega;}}_{\mathrm{rL}})$

In formula (6), ω _FL, ω _FR, ω _RLAnd ω _RRThe spin velocity that means respectively the near front wheel, off front wheel, left rear wheel and off hind wheel, obtain by utilizing four wheel speed sensors to measure;

Simplifying the magic formula tire model as shown in Figure 2, is μ=0.2,0.4 in Fig. 2,0.6, the relation of 0.8,1.0 o'clock normalization method tractive force and slip rate, when low slip rate, become linear approximate relationship as seen, this is also the basis that traditional slippage gradient method utilizes slippage slop estimation road-adhesion coefficient.According to simplifying the magic formula tire model, have

In formula (7), μ means the Real-time Road adhesion value, and B, C mean model coefficient, and the normalization method tractive force of forward and backward wheel shaft is respectively:

By formula (2), (8), (9) Shi Kede

F _x=F _xf+F _xr=μ{F _zfsin[Carctan(Bs _f)]+F _zrsin[Carctan(Bs _r)]} (10)

Step is as follows according to a preliminary estimate based on vertical road-adhesion coefficient of method of least square:

Formula (10) is expressed as to parameter criterion of identification form:

In formula (11), k means the discrete moment, y (k)=F _xExpression system output, can be calculated and be learnt by (2) formula, θ (k)=μ _RLSMean the solve for parameter vector, wherein μ _RLSExpression by the RLS method road-adhesion coefficient according to a preliminary estimate;

Mean the input regression vector, in the present invention, superscript ' expression is to matrix transpose; E (k) means identification error, utilizes with the RLS algorithm of forgetting factor and determines that in real time unique unknown quantity---the estimating step of road-adhesion coefficient is as follows:

(1) computing system output variable y (k), and calculate the iteration vector

(2) calculate identification error e (k):

(3) calculated gains vector N (k):

Wherein,

Wherein, parameter lambda is forgetting factor, can effectively reduce the impact of no longer relevant to model legacy data, and prevents that covariance from dispersing, and span is in [0.9,1] usually, and the present invention gets 0.995.

(4) calculate solve for parameter vector theta (k)

θ(k)＝θ(k-1)+N(k)e(k)

So far, can go out according to a preliminary estimate in real time road-adhesion coefficient μ _RLS.

In the vehicle operating process, tire model B parameter and C value be not constant always, affected by the factor such as the speed of a motor vehicle and have certain change, and simultaneously, the vertical speed of a motor vehicle of measured car load also contains larger noise, all can be on the certain impact of estimation generation of adhesion value.Based on such consideration, by the estimated μ gone out of RLS _RLSAs observed quantity, B, C and μ, as extended mode, utilize longitudinal dynamics equation (1) to set up the EKF model, realize the real-time update of B, C value and the filtering of v, μ, and the equation of state after discretization and observational equation are suc as formula (12):

X(k)＝f(X(k-1))+W(k-1)

(12)

Z(k)＝h(X(k))+V(k)

In formula (12), state of the system vector X=[x ₁, x ₂, x ₃, x ₄] ', and x ₁=v, x ₂=μ, x ₃=B, x ₄=C; In the present invention, superscript ' expression is to matrix transpose; W means the system Gaussian white noise vector of zero-mean, W=[w ₁w ₂w ₃w ₄] ', be w wherein ₁, w ₂, w ₃And w ₄Mean respectively four systems Gaussian white noise component, the system noise covariance matrix Q (k-1) that W is corresponding is: $Q(k-1)=\left[\begin{array}{cccc}{{\mathrm{\&sigma;}}_{w_1}}^2& 0& 0& 0\\ 0& {{\mathrm{\&sigma;}}_{w_2}}^2& 0& 0\\ 0& 0& {{\mathrm{\&sigma;}}_{w_3}}^2& 0\\ 0& 0& 0& {{\mathrm{\&sigma;}}_{w_4}}^2\end{array}\right],$ Wherein And

Mean respectively system Gaussian white noise w ₁, w ₂, w ₃And w ₄Corresponding variance; Observation vector Z=[a _x-m, v _m, μ _RLS] ', a _x-mMean the longitudinal direction of car acceleration/accel that the acceleration pick-up measurement obtains; v _mThe car speed obtained is measured in expression by GPS; μ _RLSExpression estimates that by recurrence least square the road-adhesion coefficient of gained is worth according to a preliminary estimate; V means to observe white noise inwards with the mutual incoherent zero-mean of W,

Expression by acceleration pick-up measure the longitudinal direction of car acceleration/accel obtained observation noise and

That average is 0, variance is

Gaussian white noise; Expression by GPS measure the car speed obtained observation noise and

That average is 0, variance is

Gaussian white noise;

Mean by above-mentioned method of least square estimate acquisition road-adhesion coefficient observation noise and

That average is 0, variance is

Gaussian white noise; The observation noise variance battle array R that V is corresponding can be expressed as $R=\left[\begin{array}{ccc}{{\mathrm{\&sigma;}}_{a_{x-m}}}^2& 0& 0\\ 0& {{\mathrm{\&sigma;}}_{v_m}}^2& 0\\ 0& 0& {{\mathrm{\&sigma;}}_{{\mathrm{\&mu;}}_{\mathrm{RLS}}}}^2\end{array}\right];$ F () and h () mean nonlinear system function vectorial sum observation function vector:

$f\left(X(k-1)\right)=\left[\begin{array}{c}{f}_1\left(x(k-1)\right)\\ f_2\left(x(k-1)\right)\\ f_3\left(x(k-1)\right)\\ f_4\left(x(k-1)\right)\end{array}\right],$ $h\left(X\left(k\right)\right)=\left[\begin{array}{c}{h}_1\left(X\left(k\right)\right)\\ h_2\left(X\left(k\right)\right)\\ h_3\left(X\left(k\right)\right)\end{array}\right],$ Wherein

$f_1\left(X(k-1)\right)=v(k-1)-\frac{T\{D_a{\left[v(k-1)\right]}^2+C_{\mathrm{roll}}\mathrm{mg}\}}{m}+$

$\frac{\mathrm{T\&mu;}(k-1)\left\{F_{\mathrm{zf}}\sin \right[C(k-1)\arctan \left(B(k-1)s_f\right)]+F_{\mathrm{zr}}\sin [C(k-1)\arctan \left(B(k-1)s_r\right)\left]\right\}}{m}$

$,f_2\left(X(k-1)\right)=\mathrm{\&mu;}(k-1),f_3\left(X(k-1)\right)=B(k-1),f_4\left(X(k-1)\right)=C(k-1)$

$h_1\left(X\left(k\right)\right)=\frac{\mathrm{\&mu;}\left\{F_{\mathrm{zf}}\sin \right[C\arctan \left({\mathrm{Bs}}_f\right)]+F_{\mathrm{zr}}\sin [C\arctan \left({\mathrm{Bs}}_r\right)\left]\right\}-D_a{\left[v(k-1)\right]}^2-C_{\mathrm{roll}}\mathrm{mg}}{m}$

$,h_2\left(X\left(k\right)\right)=v,h_3\left(X\left(k\right)\right)=\mathrm{\&mu;}$

T means the discrete cycle, and its representative value is 10 milliseconds, 20 milliseconds, 50 milliseconds or 100 milliseconds.

For the described equation of state of formula (12) and observational equation, use the EKF theory, Criterion filtering recursive process, this recursive process comprise that the time upgrades and measurement is upgraded:

Time upgrades:

State one-step prediction equation: $\hat{X}(k,k-1)=f\left(X(k-1)\right)$

One-step prediction error covariance matrix P (k, k-1):

P(k，k-1)＝A(k，k-1)P(k-1)A′(k，k-1)+Q(k-1)

Measure and upgrade:

Filtering gain matrix K (k): K (k)=P (k, k-1) H ' is [H (k) P (k, k-1) H ' (k)+R (k)] (k) ^-1

State estimation: $\hat{X}\left(k\right)=\hat{X}(k,k-1)+K\left(k\right)[Z\left(k\right)-h[\hat{X}(k,k-1)\left]\right]$

Estimation error variance battle array P (k): P (k)=[I-K (k) H (k)] P (k, k-1) and I are 4 * 4 unit matrix;

Wherein, A is that state of the system functional vector f asks the Jacobi matrix of partial derivative to state vector X, and H is that observation function vector h asks the Jacobi matrix of partial derivative to state vector X, i.e. the capable j column element of the i A of matrix A and H _{[i, j]}(i=1,2,3,4 j=1,2,3,4) and H _{[i, j]}(i=1,2,3 j=1,2,3,4) can try to achieve by following formula respectively:

$A_{[i,j]}=\frac{{\&PartialD;f}_i}{{\&PartialD;x}_j}\left(\hat{X}(k,k-1)\right),(i=\mathrm{1,2,3,4},j=\mathrm{1,2,3,4})$

$H_{[i,j]}=\frac{{\&PartialD;h}_i}{{\&PartialD;x}_j}\left(\hat{X}(k,k-1)\right),(i=\mathrm{1,2,3},j=\mathrm{1,2,3,4})$

Particularly, the value of each element of a matrix is as follows:

A _[1，1]＝1+T[(-2D _av)/m]

A _[1，2]＝T{F _zfsin[Carctan(Bs _f)]+F _zrsin[Carctan(Bs _r)]}/m

$A_{\left[\mathrm{1,3}\right]}=\mathrm{T\&mu;}[\frac{F_{\mathrm{zf}}{\mathrm{Cs}}_f\cos \left(C\arctan \left({\mathrm{Bs}}_f\right)\right)}{1+{\left({\mathrm{Bs}}_f\right)}^2}+\frac{F_{\mathrm{zr}}{\mathrm{Cs}}_r\cos \left(C\arctan \left({\mathrm{Bs}}_r\right)\right)}{1+{\left({\mathrm{Bs}}_r\right)}^2}]/m$

A _[1，4]＝Tμ[F _zf(arctan(Bs _f)cos(Carctan(Bs _f)))+F _zr(arctan(Bs _r)cos(Carctan(Bs _r)))]/m

A _[2，1]＝A _[2，3]＝A _[2，4]＝A _[3，1]＝A _[3，2]＝A _[3，4]＝A _[4，1]＝A _[4，2]＝A _[4，3]＝0

A _[2，2]＝A _[3，3]＝A _[4，4]＝1

H _[1，1]＝(-2D _av)/m H _[1，2]＝{F _zfsin[Carctan(Bs _f)]+F _zrsin[Carctan(Bs _r)]}/m

$H_{\left[\mathrm{1,3}\right]}=\mathrm{\&mu;}[\frac{F_{\mathrm{zf}}{\mathrm{Cs}}_f\cos \left(C\arctan \left({\mathrm{Bs}}_f\right)\right)}{1+{\left({\mathrm{Bs}}_f\right)}^2}+\frac{F_{\mathrm{zr}}{\mathrm{Cs}}_r\cos \left(C\arctan \left({\mathrm{Bs}}_r\right)\right)}{1+{\left({\mathrm{Bs}}_r\right)}^2}]/m$

H _[1，4]＝μ[F _zf(arctan(Bs _f)cos(Carctan(Bs _f)))+F _zr(arctan(Bs _r)cos(Carctan(Bs _r)))]/m

H _[2，2]＝H _[2，3]＝H _[2，4]＝H _[3，1]＝H _[3，3]＝H _[3，4]＝0

H _[2，1]＝H _[3，2]＝1

In above-mentioned filtering recursion computation process, can determine that automobile is at each tire model B parameter, C, speed v and road-adhesion coefficient μ constantly, quantity of state B, C that EKF is exported feed back to tire model, form closed loop feedback, the real-time update recurrence least square is estimated, can further improve the adhesion value estimated accuracy, respond in time the change of travel conditions, the road-adhesion coefficient μ value that EKF was exported of usining is as final estimated result.

Embodiment 2

The actual effect of the method for estimation of the travel condition of vehicle based on improved EKF proposed for check the present invention, carried out the simulating, verifying experiment on professional vehicle dynamics simulation software CarSim.

CarSim is the special simulation software for vehicle dynamics by the exploitation of U.S. MSC (Mechanical Simulation Corporation) company, by numerous in the world automakers, components supplying business, adopted at present, be widely used in the business development of modern automobile control system, become the standard software of auto trade, enjoyed a very good reputation.Vehicle dynamic model in Carsim be by respectively to car body, the suspension of automobile, turn to, the height modeling true to nature of each subsystem such as braking and each tire realizes, has very high degree of freedom, the actual information of travel condition of vehicle accurately that approaches very much can be provided, therefore, the travel condition of vehicle information of Carsim output can be used as the reference output of vehicle.

The estimation effect of algorithm under low slip rate and high slip rate operating mode for check the present invention proposition, and the algorithm of the present invention's proposition is to the response of road-adhesion coefficient sudden change, in emulation experiment, be provided with respectively high slip rate and low slip rate operating mode, and the sudden change of adhesion value, and with traditional slippage gradient method, contrast under high slip rate operating mode.Based on the method for slippage slope, think that the normalization method tractive force in little slip rate scope is directly proportional to slip rate, its proportionate relationship is called the slippage slope, and road-adhesion coefficient is directly proportional to the slippage slope, therefore, by estimating the slippage slope, just can estimate the adhesion value on tire road surface.Under low slip rate operating mode, simulation time is set to 45s, and under high slip rate operating mode, simulation time is set to 50s.Sampling frequency is 10Hz, and namely cycle T is 0.1s.Vehicle used is the four-wheeled of a front-wheel steering, and principal parameter is as follows: m=1220 (kilogram), R _e=310.8 (millimeters), a=1.040 (rice), b=1.560 (rice).The measurement noise of setting the linear velocity (cireular frequency recorded by wheel speed sensors is multiplied by tire radius and obtains) of four wheels is that average is 0, standard deviation is the Gaussian white noise of 0.04 (meter per second).But it is 14 and the 1.3[list of references that tire model B parameter, C get respectively initial value: Fredrik Gustafsson, Automotive Safety Systems.Replacing costly sensors with software algorithms, 2009, IEEE Signal Processing Magazine, 2009,26 (4): 32-47.], gravity acceleration g gets 9.78, the forgetting factor of recurrence least square is set to 0.995, and the standard deviation of the system zero average Gaussian white noise of Kalman filtering is respectively

(meter per second), And

The standard deviation of the zero-mean Gaussian white noise of three observed quantities of Kalman filtering is respectively

(meter per second ²),

(meter per second) reaches

Related results such as Fig. 3～shown in Figure 12.

Fig. 3-Fig. 6 is simulation result under low slip rate operating mode, wherein Fig. 3 is the slip rate that adhesion value was set to 0.2 o'clock, Fig. 4 utilizes algorithm of the present invention to carry out result (the black solid line representative of adhesion value estimated result in figure of adhesion value estimation under the corresponding operating mode of Fig. 3, Carsim output true value represents with black dotted line), Fig. 5 is that adhesion value is set to from the slip rate of 0.8 transition to 0.2, Fig. 6 utilizes algorithm of the present invention to carry out result (the black solid line representative of adhesion value estimated result in figure of adhesion value estimation under the corresponding operating mode of Fig. 5, Carsim output true value represents with black dotted line), by Fig. 3, visible its slip rate of Fig. 5 is all less, the method that visible this paper proposes is under low slip rate operating mode, can identify rapidly road-adhesion coefficient, accuracy is higher, error is less than 0.1, as seen from Figure 6, when road-adhesion coefficient transition, algorithm can converge to rapidly near true value in 2s, the sudden change of response road-adhesion coefficient fast.

By being set, vehicle take comparatively strong brake operating to produce high slip rate operating mode.And estimate that with traditional method based on the slippage slope road-adhesion coefficient is to compare, based on the method for slippage slope, think that the normalization method tractive force in little slip rate scope is directly proportional to slip rate, its proportionate relationship is called the slippage slope, and road-adhesion coefficient is directly proportional to the slippage slope, therefore, by estimating the slippage slope, just can estimate the adhesion value on tire road surface.Fig. 7-Figure 12 is the simulation result under high slip rate operating mode.Wherein Fig. 7 is the slip rate that adhesion value was set to 0.2 o'clock, Fig. 8 utilizes algorithm of the present invention to carry out result (the black solid line representative of adhesion value estimated result in figure of adhesion value estimation under the corresponding operating mode of Fig. 7, Carsim output true value represents with black dotted line), Fig. 9 utilizes traditional slippage gradient method to carry out the result of adhesion value estimation (in figure, the adhesion value estimated result is with deceiving the solid line representative, and Carsim output true value represents with black dotted line) under the corresponding operating mode of Fig. 7, Figure 10 is that adhesion value is set to from the slip rate of 0.8 transition to 0.2, Figure 11 utilizes algorithm of the present invention to carry out result (the black solid line representative of adhesion value estimated result in figure of adhesion value estimation under the corresponding operating mode of Figure 10, Carsim output true value represents with black dotted line), Figure 12 utilizes traditional slippage gradient method to carry out (the black solid line representative of adhesion value estimated result in figure of the result of adhesion value estimation under the corresponding operating mode of Figure 10, Carsim output true value represents with black dotted line), by Fig. 7, the visible straight skidding rate of Figure 10 is all larger, because the straight skidding rate is no longer approximate linear with longitudinal force of tire, tradition slippage gradient method has produced great error, and the accuracy that institute of the present invention employing method is being estimated, the promptness aspect of road pavement sudden change response has good effect.

Claims (1)

1. vertical adhesion value method of estimation of road, it is characterized in that: on smooth express highway pavement, for the front-wheel steering four-wheel automobile, based on car load Longitudinal Dynamic Model and simplification magic formula tire model, utilization is with recurrence least square (the Recursive least squares of forgetting factor, RLS) method goes out the vertical adhesion value of road in real time according to a preliminary estimate, further using the estimated vertical adhesion value of road gone out and tire model parameter as extended mode, utilize EKF (Extended Kalman Filter, EKF) algorithm, the filtered signal noise, and realize the self adaptation adjustment of tire model coefficient, final Real-time Obtaining is accurate, the vertical adhesion value of the road of robust is estimated, high slip rate and low slip rate operating mode on can the smooth express highway of simultaneous adaptation,

Concrete steps comprise:

1) calculate the straight skidding rate

Use S _j(j=f, r) means the longitudinal direction of car slip rate, namely can be divided into again front wheel spindle straight skidding rate s _fWith hind axle straight skidding rate s _r, subscript j gets f or r, and f or r mean respectively front or rear wheel shaft, s _jMethod of calculating is:

In formula (1), R _eMean tire radius; v _TfAnd v _TrMean respectively on forward and backward wheel shaft the speed along the tire direction, all with vehicle absolute velocitye v, replace in the present invention, obtain by the GPS measurement; ω _fThe spin velocity on front wheel spindle is converted in the spin velocity equivalence that means two wheels on front wheel spindle; ω _rMean that on hind axle, the spin velocity on hind axle, ω are converted in two rotation of wheel cireular frequency equivalences _fAnd ω _rCan unify to be designated as ω _j(j=f, r) and

${\mathrm{\&omega;}}_f=\frac{1}{2}({\mathrm{\&omega;}}_{\mathrm{fR}}+{\mathrm{\&omega;}}_{\mathrm{fL}})$

(2)

${\mathrm{\&omega;}}_r=\frac{1}{2}({\mathrm{\&omega;}}_{\mathrm{rR}}+{\mathrm{\&omega;}}_{\mathrm{rL}})$

In formula (2), ω _FL, ω _FR, ω _RLAnd ω _RRThe spin velocity that means respectively the near front wheel, off front wheel, left rear wheel and off hind wheel, obtain by utilizing four wheel speed sensors to measure;

2) calculate the car load longitudinal force

With

Mean the normalization method tractive force, namely can be divided into again front wheel spindle normalization method tractive force

With hind axle normalization method tractive force

Subscript j gets f or r, and f or r mean respectively front or rear wheel shaft:

In formula (3), F _XfAnd F _XrMean respectively longitudinal force on forward and backward wheel shaft, F _XfAnd F _XrCan unify to be designated as F _Xy(j=f, r), F _ZfAnd F _ZrVertical load before or after meaning respectively to be assigned on wheel shaft, can unify to be designated as F _Zf(j=f, r) and can be calculated as follows:

$F_{\mathrm{zf}}=\frac{\mathrm{mgb}}{(a+b)},F_{\mathrm{zr}}=\frac{\mathrm{mga}}{(a+b)}---\left(4\right)$

In formula (4), g means acceleration due to gravity, and m means vehicle mass, and a is the distance of vehicle front wheel shaft center to barycenter, and b is the distance of automobile back wheel wheel shaft center to barycenter;

The car load Longitudinal Dynamic Model is expressed as:

F _x＝F _xf+F _xr＝ma _x+D _av ²+C _rollmg (5)

In formula (5), a _xMean the real-time longitudinal acceleration of vehicle, by accelerometer measures, obtained D _aMean the air resistance constant, C _RollMean coefficient of rolling resistance, v means vehicle absolute velocitye, F _xMean the car load longitudinal force;

The present invention adopts and simplifies the magic formula tire model, and subscript j gets f or r, and f or r mean respectively front or rear wheel shaft:

In formula (6), μ means the Real-time Road adhesion value, and B, C mean model coefficient, supposes that the coefficient of road adhesion situation of each wheel is identical, and the normalization method tractive force of forward and backward wheel shaft is respectively:

By formula (5), (7), (8) Shi Kede

F _x＝F _xf+F _xr=μ{F _zfsin[Carctan(Bs _f)]+F _zrsin[Carctan(Bs _r)]} (9)

3) based on vertical road-adhesion coefficient of method of least square according to a preliminary estimate

Formula (9) is expressed as to parameter criterion of identification form:

In formula (10), k means the discrete moment, y (k)=F _xExpression system output, can be calculated and be learnt by (5) formula; θ (k)=μ _RLSMean the solve for parameter vector, wherein μ _RLSExpression by the RLS method road-adhesion coefficient according to a preliminary estimate;

Mean the input regression vector, in the present invention, superscript ' expression is to matrix transpose; E (k) means identification error; Utilize with the RLS algorithm of forgetting factor and determine that in real time unique unknown quantity---the estimating step of road-adhesion coefficient is as follows:

(1) computing system output variable y (k), and calculate the iteration vector

(2) calculate identification error e (k);

(3) calculated gains vector N (k);

Wherein,

Wherein, parameter lambda is forgetting factor;

(4) calculate solve for parameter vector theta (k);

θ(k)=θ(k-1)+N(k)e(k)

So far, can go out according to a preliminary estimate in real time road-adhesion coefficient μ _RLS

4) vertical road-adhesion coefficient of extension-based Kalman filtering is estimated

By the estimated μ gone out of RLS method _RLSAs observed quantity, B, C and μ, as extended mode, utilize longitudinal dynamics equation (5) to set up the EKF model, and the equation of state after discretization and observational equation are suc as formula (11):

X(k)＝f(X(k-1))+W(k-1)

(11)

Z(k)＝h(X(k))+V(k)

In formula (11), state of the system vector X=[x ₁, x ₂, x ₃, x ₄] ', and x ₁=v, x ₂=μ, x ₃=B, x ₄=C; W means the system Gaussian white noise vector of zero-mean, W=[w ₁w ₂w ₃w ₄] ', be w wherein ₁, w ₂, w ₃And w ₄Mean respectively four systems Gaussian white noise component, the system noise covariance matrix Q (k-1) that W is corresponding is: $Q(k-1)=\left[\begin{array}{cccc}{{\mathrm{\&sigma;}}_{w_1}}^2& 0& 0& 0\\ 0& {{\mathrm{\&sigma;}}_{w_2}}^2& 0& 0\\ 0& 0& {{\mathrm{\&sigma;}}_{w_3}}^2& 0\\ 0& 0& 0& {{\mathrm{\&sigma;}}_{w_4}}^2\end{array}\right],$ Wherein

And

Mean respectively system Gaussian white noise w ₁, w ₂, w ₃And w ₄Corresponding variance; Observation vector Z=[a _x-m, v _m, μ _RLS] ', a _x-mMean the longitudinal direction of car acceleration/accel that the acceleration pick-up measurement obtains; V _MThe car speed obtained is measured in expression by GPS; μ _RLSExpression estimates that by recurrence least square the road-adhesion coefficient of gained is worth according to a preliminary estimate; V means the mutual incoherent zero-mean observation white noise vector with W,

Expression by acceleration pick-up measure the longitudinal direction of car acceleration/accel obtained observation noise and

That average is 0, variance is Gaussian white noise;

Expression by GPS measure the car speed obtained observation noise and

That average is 0, variance is

Gaussian white noise;

Mean by above-mentioned method of least square estimate acquisition road-adhesion coefficient observation noise and

That average is 0, variance is

Gaussian white noise; The observation noise variance battle array R that V is corresponding can be expressed as $R=\left[\begin{array}{ccc}{{\mathrm{\&sigma;}}_{a_{x-m}}}^2& 0& 0\\ 0& {{\mathrm{\&sigma;}}_{v_m}}^2& 0\\ 0& 0& {{\mathrm{\&sigma;}}_{{\mathrm{\&mu;}}_{\mathrm{RLS}}}}^2\end{array}\right];$ F () and h () mean nonlinear system function vectorial sum observation function vector:

$f\left(X(k-1)\right)=\left[\begin{array}{c}{f}_1\left(X(k-1)\right)\\ f_2\left(X(k-1)\right)\\ f_3\left(X(k-1)\right)\\ f_4\left(X(k-1)\right)\end{array}\right],$ $h\left(X\left(k\right)\right)=\left[\begin{array}{c}{h}_1\left(X\left(k\right)\right)\\ h_2\left(X\left(k\right)\right)\\ h_3\left(X\left(k\right)\right)\end{array}\right],$ Wherein

$f_1\left(X(k-1)\right)=v(k-1)-\frac{T\{D_a{\left[v(k-1)\right]}^2+C_{\mathrm{roll}}\mathrm{mg}\}}{m}+$

$\frac{\mathrm{T\&mu;}(k-1)\left\{F_{\mathrm{zf}}\sin \right[C(k-1)\arctan \left(B(k-1)s_f\right)]+F_{\mathrm{zr}}\sin [C(k-1)\arctan \left(B(k-1)s_r\right)\left]\right\}}{m}$

$,f_2\left(X(k-1)\right)=\mathrm{\&mu;}(k-1),f_3\left(X(k-1)\right)=B(k-1),f_4\left(X(k-1)\right)=C(k-1)$

$h_1\left(X\left(k\right)\right)=\frac{\mathrm{\&mu;}\left\{F_{\mathrm{zf}}\sin \right[C\arctan \left({\mathrm{Bs}}_f\right)]+F_{\mathrm{zr}}\sin [C\arctan \left({\mathrm{Bs}}_r\right)\left]\right\}-D_a{\left[v(k-1)\right]}^2-C_{\mathrm{roll}}\mathrm{mg}}{m}$

$,h_2\left(X\left(k\right)\right)=v,h_3\left(X\left(k\right)\right)=\mathrm{\&mu;}$

T means the discrete cycle, and its representative value is 10 milliseconds, 20 milliseconds, 50 milliseconds or 100 milliseconds;

For the described equation of state of formula (11) and observational equation, use the EKF theory, Criterion filtering recursive process, this recursive process comprise that the time upgrades and measurement is upgraded:

Time upgrades:

State one-step prediction equation: $\hat{X}(k,k-1)=f\left(X(k-1)\right)$

One-step prediction error covariance matrix P (k, k-1):

P(k，k-1)=A(k，k-1)P(k-1)A’(k，k-1)+Q(k-1)

Measure and upgrade:

Filtering gain matrix K (k): K (k)=P (k, k-1) H ' is [H (k) P (k, k-1) H ' (k)+R (k)] (k) ^-1

State estimation: $\hat{X}\left(k\right)=\hat{X}(k,k-1)+K\left(k\right)[Z\left(k\right)-h[\hat{X}(k,k-1)\left]\right]$

Estimation error variance battle array P (k): P (k)=[I-K (k) H (k)] P (k, k-1) and I are 4 * 4 unit matrix;

Wherein, A is that state of the system functional vector f asks the Jacobi matrix of partial derivative to state vector X, and H is that observation function vector h asks the Jacobi matrix of partial derivative to state vector X, i.e. the capable j column element of the i A of matrix A and H _{[i, j]}(i=1,2,3,4 j=1,2,3,4) and H _{[i, j]}(i=1,2,3 j=1,2,3,4) can try to achieve by following formula respectively:

$\begin{array}{c}{A}_{[i,j]}=\frac{\&PartialD;f_i}{{\&PartialD;x}_j}\left(\hat{X}(k,k-1)\right)\end{array},(i=\mathrm{1,2,3,4},j=\mathrm{1,2,3,4})$

$\begin{array}{c}{H}_{[i,j]}=\frac{\&PartialD;h_i}{{\&PartialD;x}_j}\left(\hat{X}(k,k-1)\right),(i=\mathrm{1,2,3},j=\mathrm{1,2,3,4})\end{array}$

Particularly, the value of each element of a matrix is as follows:

A _[1，1]＝1+T[(-2D _av)/m]

A _[1，2]=T{F _zfsin[Carctan(Bs _f)]+F _zrsin[Carctan(Bs _r)]}/m

$A_{\left[\mathrm{1,3}\right]}=\mathrm{T\&mu;}[\frac{F_{\mathrm{zf}}C{s}_f\cos \left(C\arctan \left({\mathrm{Bs}}_f\right)\right)}{1+{\left({\mathrm{Bs}}_f\right)}^2}+\frac{F_{\mathrm{zr}}{\mathrm{Cs}}_r\cos \left(\mathrm{Catc}\tan \left({\mathrm{Bs}}_r\right)\right)}{1+{\left(B{s}_r\right)}^2}]/m$

A _[1,4]＝Tμ[F _zf(arctan(Bs _f)cos(Carctan(Bs _f)))+F _zr(arctan(Bs _r)cos(Carctan(Bs _r)))]/m

A _[2，1]=A _[2，3]=A _[2，4]=A _[3，1]=A _[3，2]=A _[3，4]=A _[4，1]=A _[4，2]=A _[4，3]=0

A _[2，2]=A _[3，3]=A _[4，4]=1

H _[1，1]＝(-2D _av)/m H _[1，2]＝{F _zfsin[Carctan(Bs _f)]+F _zrsin[Carctan(Bs _r)]}/m

$H_{\left[\mathrm{1,3}\right]}=\mathrm{\&mu;}[\frac{F_{\mathrm{zf}}C{s}_f\cos \left(C\arctan \left({\mathrm{Bs}}_f\right)\right)}{1+{\left({\mathrm{Bs}}_f\right)}^2}+\frac{F_{\mathrm{zr}}{\mathrm{Cs}}_r\cos \left(C\arctan \left({\mathrm{Bs}}_r\right)\right)}{1+{\left(B{s}_r\right)}^2}]/m$

H _[1，4]＝μ[F _zf(arctan(Bs _f)cos(Carctan(Bs _f)))+F _zr(arctan(Bs _r)cos(Carctan(Bs _r)))]/m

H _[2，2]=H _[2，3]=H _[2，4]=H _[3，1]=H _[3，3]＝H _[3，4]=0

H _[2，1]=H _[3，2]＝1

The filtered μ value that the above-mentioned EKF filtering recursion of usining is exported is as final road-adhesion coefficient estimated result.

Images