# Joint vehicle state and parameters estimation
via Twin-in-the-Loop observers

F. Dettùa. S. Formentina and S. M. Savaresia CONTACT F. Dettù. email:
<EMAIL_ADDRESS>S. Formentin, email<EMAIL_ADDRESS>S.
M. Savaresi, email<EMAIL_ADDRESS>aDipartimento di Elettronica,
Informazione e Bioingegneria, Politecnico di Milano, Piazza Leonardo da Vinci
32, 20133, Milano, Italy.

###### Abstract

Vehicular control systems are required to be both extremely reliable and
robust to different environmental conditions, _e.g._ load or tire-road
friction. In this paper, we extend a new paradigm for state estimation, called
Twin-in-the-Loop filtering (TiL-F), to the estimation of the unknown
parameters describing the vehicle operating conditions. In such an approach, a
digital-twin of the vehicle (usually already available to the car
manufacturer) is employed on-board as a plant replica within a closed-loop
scheme, and the observer gains are tuned purely from experimental data. The
proposed approach is validated against experimental data, showing to
significantly outperform the state-of-the-art solutions.

###### keywords:

Parameter estimation; state estimation; driving simulator; sensors; signal
processing.

††articletype: RESEARCH ARTICLE

##  1 Introduction

The Twin-in-the-Loop (TiL) framework has been recently proposed for vehicle
dynamics control [2] and estimation of unknown states and variables [13]. In
this innovative framework, an high-fidelity off-the-shelf vehicle dynamics
simulator is employed in a control/estimation loop - thanks to the recent
developments in on-board hardware platforms 111For example, https://www.vi-
grade.com/en/products/autohawk/.. Being the simulator a commercial product, we
assume it to be a black-box, _i.e._ we can only extract its outputs given some
inputs, but we are not able to access its equations and update rules directly.
With regards to the TiL estimator, the authors in [13] showed how the
presented algorithm is able to estimate both unmeasured states (_e.g._ vehicle
sideslip) and other unknown variables, such as the tire-road forces. However,
they did not address the estimation of time-varying parameters; in fact,
modern vehicle dynamics controllers are strongly based on time-domain models
of the process [9, 19, 18]. Due to the variability of conditions during the
daily use of a vehicle - _e.g._ the most common case being the payload - such
controllers can underperform, or even destabilize the system, if such
variations are not properly addressed. Hence, following a well established
research trend on vehicle state and parameters estimation, we build this
research within the TiL framework, extending it to a new application.

Related works. Existing literature solutions for four-wheeled vehicles on-line
parameter estimation can be broadly split among direct sensing and software
sensing ones. In the first class of solutions, a sensor for the measurement of
the variable of interest is used (_e.g._ [26], where the mass is estimated via
strain gauges on the suspensions). The second class of solutions features
instead the estimation of unmeasurable variables via production available
sensors (_i.e._ without adding new piece of hardware). If we focus on software
sensing algorithms, we find many solutions in the literature - where we cite
here the most relevant contributions in the past years - [24, 15, 7, 25, 12,
27, 6, 14], which can be categorized according to four axes: estimated
quantities, employed sensors, used estimation algorithm, and employed model.
Other important contributions design software sensing with an highly specific
goal, like rollover risk evaluation [23]: we do not consider them in this
review.

On the estimated quantities side, the mass $M$ is indeed the most important
vehicle parameter, and is always considered in the cited works; vehicle yaw
inertia $J_{zz}$ is then considered in [24, 7, 25], while roll and pitch
inertia $J_{xx}$, $J_{yy}$ are considered in [15]. Other estimated parameters
include longitudinal [24, 15, 7], lateral [15] or vertical [23] center-of-mass
position variation. The tire-road friction coefficient is estimated in [14] \-
however, this parameter is different in nature from the other considered in
the literature, as its variations are not imputable to a variation of the
vehicle load conditions, but rather to the terrain.

With regards to the sensing layout, one common denominator is the use of
classical production vehicle sensors, like gyroscopes or Inertial Measurement
Units (IMU); yaw rate, which is a gyroscope output, is employed in all the
cited articles but [15]. Acceleration information is used in all references;
[15] use three acceleration sensors, to be placed above the unsprung masses.
Other considered sensors include wheel angular speed [27, 14], suspension
displacement [14], GPS position measurements [14] and vehicle longitudinal
speed [24, 25, 7].

Moving further to the considered estimation algorithm, most works are based on
Kalman Filtering (KF) [24, 25, 7, 12, 14], specifically declined - except for
[25, 12, 14] \- in its Dual formulation; Dual KF is based on two separated
observers, one devoted to classical state filtering, and the other one solely
devoted to parameter estimation, to be switched-off once a good estimate of
said parameters has been reached. [27] propose a similar formulation, where
the correction law is computed via a particle filter. On the other hand, [15]
propose a different approach, based on the use of free decay responses of the
vehicle - more similar to a classic identification algorithm rather than to a
state observer.

Finally, when coming to employed plant replica, the simplest model is the
single-track, which is however only used in [12], as it fails to capture most
of vehicle dynamics non-linearities. The most common model is the double-track
[24, 25, 7, 12, 27]; the latter can be enhanced with roll dynamics, as in [24,
7], or wheel dynamics [27], depending on the estimation target or the
available measurements. [15] consider a full vehicle model - _i.e._ featuring
also unsprung mass dynamics, roll and pitch dynamics. Finally, the most
complex model is the one proposed in [14]: the authors consider a multibody
model - featuring 14 states.

Contributions. From the literature review, an important point emerges. As also
pointed out in [13], classical vehicular state observers are based on
simplified models of the variables to be estimated: an ad-hoc model for each
variable is required, and each model has to be accurately calibrated. Given
that approximately $75\%$ of the time of a control system design project is
devoted to modelling, as indicated in [4], the resulting procedure turns out
to be cumbersome. Let this apart, consider that any estimator is based upon a
set of hyperparameters controlling its behaviour - _e.g._ the noise
covariances in Kalman Filters. Generally speaking, the more complex an
estimator is, the longer the hyperparameter selection is. Relying on simple
correction laws is thus beneficial in the implementation phase.

In TiL framework we employ a single model - a digital twin - for each variable
to be estimated: this solves the problem of designing ad-hoc models and of
having many models running at run-time on the same estimation algorithm. The
model calibration phase comes for free when using a digital twin, in the
assumption that a car manufacturer already has a library of simulators for
their production vehicles. On the other hand, the correction algorithm
complexity is now shifted to the model fidelity: the more faithful to the real
system the model is, the less intense - and thus simple - the correction need
to be. The authors in [13] actually show as a linear and time-invariant
correction law is sufficient to outperform the benchmark.

Matter-of-factly, the authors in [14] consider an high-fidelity model of the
plant as in TiL framework: however, being the considered model a white-box,
the authors are able to directly compute the update rules via a Kalman Filter-
like algorithm. This is not possible in TiL, as the underlying model is a
black-box, and the correction law, although linear, has to be tuned via model-
free approaches [13].

With respect to the research in [13], we extend the TiL framework to the
estimation of varying parameters, providing insights on how the correction law
has to be designed, and on the feasibility of the approach in terms of
robustness to measurement noise and uncertainty.

The paper is structured as follows. In Section 2 we provide to the reader some
preliminaries on the TiL estimator. Section 3 defines the problem of varying
mass and moment of inertia, and its estimation. Then, in Section 4 we conduct
an in-depth analysis of the parameters estimation problem by independently
considering set of variables, providing sensitivity analyses to uncertainty
and noise. Section 5 shows an application of the described and simulated
framework on experimental data, while Section 6 concludes the manuscript with
a few considerations.

##  2 Twin-in-the-Loop observer architecture

Figure 1: Twin-in-the-Loop estimation scheme. Figure 2: Twin-in-the-Loop
complete estimation architecture, featuring correction of simulator states and
augmented states.

[13] proposed an unified estimator for all the variables of interest in
vehicle dynamics control. In such an architecture (see Fig. 1), the two main
features are the following

  * •

The vehicle simulator, an high-fidelity plant replica, whose equations are a
black-box.

  * •

The closed-loop correction, based upon available measurements, perturbs the
model in order to properly estimate unknown variables.

The detailed architecture of the complete observer is given in Fig. 2.
Generically speaking, the simulator states at time $k$ are updated according
to the following equations

$\begin{split}\tilde{x}_{k+1}&=f(\hat{x}_{k},u_{k},\tilde{z}_{k},p),\\\
\tilde{y}_{k}&=g(\hat{x}_{k},u_{k},\tilde{z}_{k},p),\\\
\tilde{z}_{k}&=h(\hat{x}_{k},u_{k},\delta\hat{z}_{k-1},p)\end{split}$ (1)

where $f,\ g,\ h$ are unknown (possibly nonlinear) functions,
$y\in\mathbb{R}^{n_{y}}$ is the set of measurable outputs,
$x\in\mathbb{R}^{n_{x}}$ is the set of internal states of the simulator,
$z\in\mathbb{R}^{n}{{}_{z}}$ is a set of additional variables to be estimated
and $p$ is a set of constant parameters. We denote $\tilde{\nu}_{k}$ an
a-priori estimation for a certain variable $\nu_{k}$: then closed-loop
correction updates the a-priori estimates, and generates a-posteriori ones
$\hat{\nu}_{k}$ \- as in the widely known Kalman Filter formulation.

The authors in [13] extend the state vector to estimate tire-road contact
forces; more specifically, the additional variables $z$ are estimated by
correcting the nominal values provided by the simulator

$\hat{z}_{k}=\tilde{z}_{k}+\delta\hat{z}_{k}.$ (2)

The artificial variable $\delta\hat{z}_{k}$ is described by a fictitious
constant dynamics equation

$\delta\tilde{z}_{k+1}=\delta\hat{z}_{k}.$ (3)

The latter assumption is widely used in case the parameters to be estimated
are slowly varying. Overall, the set of states, including internal and
extended ones, now reads $x^{aug}=\begin{bmatrix}x^{T}&\delta
z^{T}\end{bmatrix}$. At this point, the innovation term at step $k$,
consisting in the mismatch among measured ($y_{k}$) and a-priori estimated
outputs ($\tilde{y}_{k}$) is to be used to correct the state vector $x^{aug}$.
The innovation is mapped onto the states via a linear law in [13]

$\begin{split}\Delta x_{k}^{aug}&=K\left(y_{k}-\tilde{y}_{k}\right)\\\
\hat{x}^{aug}_{k}&=\tilde{x}^{aug}_{k}+\Delta x^{aug}_{k}.\end{split}$ (4)

However, we will show in the following that the linear formulation is
inadequate for the estimation of certain parameters, _e.g._ the vehicle mass
or inertia. Let us thus define a mixed linear-nonlinear correction law

$\hat{x}_{k}^{aug}=\tilde{x}_{k}^{aug}+\begin{bmatrix}K&0\\\
0&\mathcal{K}\left(y_{k}\right)\end{bmatrix}\left(y_{k}-\tilde{y}_{k}\right),$
(5)

where $K\in\mathbb{R}^{n_{x}\times n_{y}}$ maps the innovation onto the
nominal states, and $\mathcal{K}\left(y_{k}\right)\in\mathbb{R}^{n_{z}\times
n_{y}}$ maps the innovation onto the extended states.

##  3 A parameter estimation case study: the uncertain load setting

Figure 3: Vehicle top view, considering differently added additional masses.

Among the physical parameters defining a vehicle, mass and inertia might
significantly change due to different load conditions, while being extremely
important in several on-board control systems. Hence, we consider in the
following the problem of estimating such parameters in a road vehicle via the
TiL estimator architecture. Consider the 5-seats vehicle in Fig. 3. Such
vehicle can be loaded with additional masses. On the mathematical perspective,
we assume each of the $m_{i},\ i=1,\ldots,\mathcal{M}$ to be a point mass,
characterized by a certain position
$\textrm{p}_{m_{i}}=\left(p_{m_{i,x}},p_{m_{i,y}},p_{m_{i,z}}\right)$ in the
reference frame defined by axes $X,\ Y,\ Z$ \- as displayed in Fig. 3. If we
assume that the vehicle chassis can be modeled as a rigid body with its mass
concentrated in a single point - a common assumption in vehicle dynamics - the
center-of-mass of such system is computed as [20]

$\begin{split}CM_{x}&=\sum_{i=0,\ldots,\mathcal{M}}m_{i}\cdot p_{m_{i,x}},\\\
CM_{y}&=\sum_{i=0,\ldots,\mathcal{M}}m_{i}\cdot p_{m_{i,y}},\\\
CM_{z}&=\sum_{i=0,\ldots,\mathcal{M}}m_{i}\cdot p_{m_{i,z}}.\end{split}$ (6)

whereas $m_{0}$ denotes the original lumped mass of the chassis. Indeed, the
center-of-mass including additional masses -
$CM=\left(d_{x},d_{y},d_{z}\right)$ \- is going to be different from the
original one - $CM_{0}=\left(d_{0,x},d_{0,y},d_{0,z}\right)$. The new lumped
mass is then simply defined as

$M_{tot}=m_{0}+\sum_{i=1,\ldots,\mathcal{M}}m_{i}.$ (7)

Consequently, the chassis moments and products of inertia change due to the
mass variation and the re-defined center-of-mass [20]. Considering the
Huygens-Steiner theorem [5], we obtain the new moments/products of inertia,
defined with respect to $CM$; this is done starting from the ones defined in
$CM_{0}$ and taking into account the newly added masses. _E.g._ , for the
moments of inertia - $J_{xx},\ J_{yy},\ J_{zz}$ \- one has

$\begin{split}J_{xx}=&J_{xx,0}-m_{0}\left[\delta\left(d_{y},d_{y,0}\right)+\delta\left(d_{z},d_{z,0}\right)\right]+\\\
&\sum_{i=1,\ldots,\mathcal{M}}m_{i}\left[\delta\left(p_{m_{i},y},d_{y}\right)+\delta\left(p_{m_{i},z},d_{z}\right)\right],\\\
J_{yy}=&J_{yy,0}-m_{0}\left[\delta\left(d_{x},d_{x,0}\right)+\delta\left(d_{z},d_{z,0}\right)\right]+\\\
&\sum_{i=1,\ldots,\mathcal{M}}m_{i}\left[\delta\left(p_{m_{i},x},d_{x}\right)+\delta\left(p_{m_{i},z},d_{z}\right)\right],\\\
J_{zz}=&J_{zz,0}-m_{0}\left[\delta\left(d_{x},d_{x,0}\right)+\delta\left(d_{y},d_{y,0}\right)\right]+\\\
&\sum_{i=1,\ldots,\mathcal{M}}m_{i}\left[\delta\left(p_{m_{i},x},d_{x}\right)+\delta\left(p_{m_{i},y},d_{y}\right)\right],\\\
\end{split}$ (8)

whereas $J_{xx,0},\ J_{yy,0},\ J_{zz,0}$ represent the chassis nominal moment
of inertia, and $\delta\left(x_{1},x_{2}\right)=\left(x_{1}-x_{2}\right)^{2}$.
Also product of inertia variations - _e.g._ $J_{xy}$\- can be quantified
analytically in a similar way; for the load conditions considered in the case
study - realistic ones for road vehicles - we don’t have significant
variations of these parameters.

Now, let us consider again scheme of Fig. 1: the high-fidelity digital twin is
modelled with a set of nominal parameters, _i.e._ lumped mass, inertia, and
center-of-mass location. The real vehicle is instead loaded with additional
masses, and thus ”perturbed” from the nominal parameters set.

In [13], the authors design a closed-loop correction so as to force the
simulator states to be identical to the ones of the real vehicle, regardless
of the reasons leading to differences between the two instances of the
vehicle. On the other hand, in the present work we shift the philosophy,
trying to estimate also the exact differences among the two instances.

##  4 Simulation analysis

Nominal param. | Value
---|---
$m_{0}\ [kg]$ | $2125.8$
$J_{xx,0}\ [kg\cdot m^{2}]$ | $834.23$
$J_{yy,0}\ [kg\cdot m^{2}]$ | $3640.182$
$J_{zz,0}\ [kg\cdot m^{2}]$ | $3932.77$
$J_{xy,0}\ [kg\cdot m^{2}]$ | $0.14$
$J_{xz,0}\ [kg\cdot m^{2}]$ | $0.097$
$J_{yz,0}\ [kg\cdot m^{2}]$ | $3.86$
$d_{0,x}\ [cm]$ | $-125.0$
$d_{0,y}\ [cm]$ | $-0.003$
$d_{0,z}\ [cm]$ | $64.4$
Table 1: Nominal vehicle mass, inertia and CM parameters.

Additional load param. | Value
---|---
$m_{p,1}\ [kg]$ | $75\ kg$
$m_{p,2}\ [kg]$ | $80\ kg$
$m_{p,3}\ [kg]$ | $65\ kg$
$m_{p,4}\ [kg]$ | $75\ kg$
$m_{t,1}\ [kg]$ | $30\ kg$
$m_{t,2}\ [kg]$ | $30\ kg$
Table 2: Additional load parameters.

Perturbed param. | Value | Variation222Where the variation from value $v_{0}$ to value $v_{1}$ is defined as $\Delta=100(v_{0}/v_{1}-1)$.
---|---|---
$M_{tot}\ [kg]$ | $2480.8$ | $+16.70\ \%$
$J_{xx}\ [kg\cdot m^{2}]$ | $901.9$ | $+8.11\ \%$
$J_{yy}\ [kg\cdot m^{2}]$ | $4394.4$ | $+20.72\ \%$
$J_{zz}\ [kg\cdot m^{2}]$ | $4760.0$ | $+21.03\ \%$
$d_{x}\ [cm]$ | $-131.6$ | $+5\ \%$
$d_{y}\ [cm]$ | $1.6$ | $+100.2\ \%$
$d_{z}\ [cm]$ | $68.4$ | $+5.3\ \%$
Table 3: Perturbed vehicle parameters.

In this Section, we demonstrate the TiL methodology for the case study of
Section 3. For a better assessment of the potentialities and limitations of
the present method, we conduct a preliminary analysis in a controlled
simulation environment. The nominal chassis parameters of the considered
vehicle are given in Table 2. Said nominal parameters also include the
presence of a driver ($m_{d}$ in Fig. 3). The parameters are those of a
generic sport utility vehicle, which is modeled in VI-CarRealTime (CRT)
simulation environment [21]. The nominal vehicle model is perturbed by adding
masses according to the scheme of Fig. 3; the values of such loads are given
in Table 2. For this new configuration, the center-of-mass, mass and inertia
can be analytically computed according to what described in Section 3; the new
values for the most relevant parameters are in Table 3.

We apply the TIL architecture to the problem herein, defining the complete
extended state vector as

$x^{aug}=\begin{bmatrix}x^{T}&\delta M&\delta J_{xx}&\delta J_{yy}&\delta
J_{zz}\end{bmatrix}^{T},$ (9)

where $x$ is the set of VI-CarRealTime states, $x\in R$. In this simulation
analysis, only extended states - parameters - are going to be corrected, as we
want to focus on the feasibility of parameters correction. Each extended state
represents the deviation of said parameter from the nominal value provided by
the simulator as an output. Joint state and parameter correction will be
considered in an experimental case study in Section 5.

In order to reproduce a simulation environment the closer to a real scenario,
we introduce measurement noise. The noise is introduced on the signals coming
from an - simulated in CRT - inertial measurement unit (IMU) placed near the
vehicle CM. The following noise law is considered

$s_{n}=s+n_{s},\ n_{s}\sim WN(0,\sigma_{s}^{2}),$ (10)

where $s$ is a signal coming from an IMU, and $\sigma_{s}$ is the standard
deviation of a suitable white noise. In the simulation case study, $x$-axis
acceleration $a_{x}$ and $x,\ y,\ z$ axes angular rates $\omega_{x,y,z}$ are
employed as measurements in the observer, hence, only these signals are
perturbed with noise. When not specified, the $snr$ is tuned to be $\approx
10$: however, we also conduct sensitivity analyses to verify the impact of
measurement noise onto the estimator.

When quantifying the observer performance, we consider as a metric the root
mean square estimation error in the last $w=1$ seconds of test, as the
dynamics of the estimated parameters are rather slow, and considering larger
windows would take into account transient effects.

$rms_{v}=\dfrac{1}{wf_{s}}\sum_{k=N_{s}-wf_{s}}^{N_{s}}\left(\nu_{k}-\hat{\nu}_{k}\right)^{2},$
(11)

Where $N_{s}$ is the number of samples in an experiment, $f_{s}\ [Hz]$ is the
sampling frequency, and $\nu$ is a suitable parameter or variable to be
estimated, _e.g._ $\nu=\delta M$.

On the other hand, we also consider the percentage version of $rms_{w}$, in
order to quantify the amount of error within the parameters values. This
metric is defined as

$rms_{\nu,\%}=rms_{\nu}/\nu_{0},$ (12)

Where $\nu_{0}$ is the nominal value of a variable to be estimated; _e.g._ in
our case, the nominal value for $\delta M$ is the difference among the
perturbed value in Table 3 and the nominal value in Table 2.

In the following, we consider as acceptable a percentage $rms$ smaller than
the $10\ \%$ of a perturbated variable.

###  4.1 Center-of-mass position identification

The $CM$ position is perturbed from $CM_{0}$ in the modified layout, as
evident from Table 3. Most significant displacement occur along $x$ and $z$
axes, as $y$ axis variation are mostly due to extremely unbalanced - on the
left-right direction - loads, which are not common in road vehicles. Given
that the digital twin computes the torque and force balances referring to the
center-of-mass, a not exact knowledge of $CM$ might be detrimental for
estimating the other parameters.

Two possible ways of solving this problem exist:

  * •

Insert the variables $d_{x},\ d_{y},\ d_{z}$ onto the augmented state vector
$x^{aug}$, and consequently correct the mismatch from the nominal values;

  * •

Identify in a traditional way - see _e.g._ [17, 8] \- $d_{x},\ d_{y},\ d_{z}$
on-line via available data and simplified models.

In the present work, we decide to follow the second approach, as we want to
concentrate onto the estimation of mass and inertia: hence, we assume correct
knowledge of $CM$. However, we also perform sensitivity analyses to check the
effect of a wrong $CM$ onto the estimation performance: in general, we will
show as for the considered load conditions, which are reasonable for a
passenger car, the performance loss is negligible even in case of wrong $CM$
estimation.

###  4.2 Vehicle mass estimation

(a) Mass estimation in an urban driving context. (b) Urban driving commands.

Figure 4: Mass estimation in an urban driving like context, with noisy
measurements. Estimated mass for different initial conditions.

Figure 5: Mass estimation with varying noise levels.

Test conditions | $rms_{M}\ [\%]$
---|---
| • Noiseless data ($snr=\infty$);
---
• Correct $CM$
$0.64$
| • Noisy data ($snr=10$);
---
• Correct $CM$.
$1.61$
| • Noisy data ($snr=10$);
---
• Wrong $CM$.
$0.68$
Figure 6: Mass estimation performance in different conditions.

Vehicle chassis mass can be estimated by suitably exciting the longitudinal
dynamics. Consider the following simplified longitudinal dynamics balance for
a road vehicle

$M_{chassis}\cdot a_{x}=F_{trac}-F_{brake}-F_{roll}-F_{aero}.$ (13)

It is clear that, provided that roughly the same forces act on the vehicle, an
increased mass yields a lower longitudinal acceleration. Hence, if the digital
twin is characterized by a smaller mass than the real vehicle, a suitable
correction can be based upon the difference between longitudinal
accelerations. Indeed, if the terms on the right hand side of Eq. (13) are
significantly different, we are not able to determine whether an acceleration
variation is due to a mass difference or rather to an error in modelling the
brake system or the tire-road interaction. For this reason, in this simulation
analysis the sole modelling difference among the two simulated instances of
the vehicle lie in the mass perturbation. In Section 5 we will test the mass
estimation on real data, where other modelling differences might exist.

The correction applied on $\delta M$ at time step $k$ reads

$\Delta\delta M_{k}=K_{a_{x}-\delta
M}\cdot\textrm{sgn}(a_{x,k})\cdot(a_{x,k}-\tilde{a}_{x,k}),$ (14)

Where $K_{a_{x},\delta M}$ is a suitable gain to be tuned and
$\textrm{sgn}\left(\cdot\right)$ extracts the sign of the argument. Note as
$\mathcal{K}\left(y_{k}\right)=K_{a_{x}-\delta M}\cdot\textrm{sgn}(a_{x,k})$
in this context.

The importance of the sign operator is easily understood with a simple
example. Consider a coasting down - deceleration due to inertia $a_{x}<0$ \-
event: if the vehicle has an higher mass than the simulator, it will coast
down faster, thus $\left(a_{x}-\tilde{a}_{x}\right)<0$. Correcting $\delta M$
proportionally to $\left(a_{x}-\tilde{a}_{x}\right)$, with the same sign,
would lead to a decrement of the extended state: however, $\delta M$ should
instead increase to compensate for the greater mass on the vehicle. Thus, the
introduction of the sign operator accounts for this issue, an introduces a
necessary non-linearity in the correction law.

The considered experiment for the mass estimation is an urban driving one -
see Fig. 4b: in said test, only the correction of Eq. (14) is applied.

The results - in case of noisy data and correct center-of-mass - can be
appreciated in Fig. 4a: the mass is correctly estimated even by starting from
different initial conditions, and the estimate eventually converges to the
correct value. Table 6 reports the estimation performance in different
conditions: as one can notice, the mass estimate is highly robust both to
measurement noise and to imperfect knowledge of the center-of-mass. In the
worst case, $rms_{M}=0.68\ \%$: since $\delta M=355\ kg$, the error is smaller
than $3\ kg$.

Figure 6 compares the estimator performance for an $snr$ growing from $5$ to
$20$: we can conclude that for reasonable noise levels, the estimator
performance does not drop below the $10\%$, which was originally set to be our
limit.

Finally, note that $\delta M$ converges to the real value even if we are not
correcting $\delta J_{xx},\ \delta J_{yy},\delta J_{zz}$: this highlights the
good decoupling among the estimated variables, without the need for designing
a separated model for each of them.

###  4.3 Yaw inertia and roll inertia estimation

(a) Correct center-of-mass estimation case. (b) Wrong center-of-mass
estimation case.

Figure 7: Roll and yaw inertia estimation by means of a swept steer
experiment, with noisy measurements. The upper plot depicts the estimated
parameters in case the center-of-mass is exactly known, while the lower one
depicts them in case the CM is not exactly known.

Figure 8: Roll inertia estimation with varying noise levels.

Figure 9: Yaw inertia estimation with varying noise levels.

Roll and yaw inertia can be estimated by exciting lateral dynamics and left-
to-right load transfer, which often come together, _e.g._ while negotiating a
curve. A proper experiment to estimate both variables at the same time is a
swept steer maneuver, at constant speed, from $0$ to $4$ $Hz$. The non-linear
correction law to be applied onto $\delta J_{xx}$ and $\delta J_{zz}$ follows
from what described in Section 4.2 for the mass

$\begin{split}\Delta\delta J_{xx,k_{k}}&=K_{\dot{\omega}_{x}-\delta
J_{xx}}(\dot{\omega}_{x,k}-\tilde{\dot{\omega}}_{x,k})\cdot\textrm{sgn}(\dot{\omega}_{x,k}),\\\
\Delta\delta J_{zz,k_{k}}&=K_{\dot{\omega}_{z}-\delta
J_{xx}}(\dot{\omega}_{z,k}-\tilde{\dot{\omega}}_{z,k})\cdot\textrm{sgn}(\dot{\omega}_{z,k}).\end{split}$
(15)

Where the angular accelerations $\dot{\omega}$ are employed - such signals are
easily obtained from the angular rates via differentiation. Note that while
estimating $J_{xx}$ and $J_{zz}$, $\delta M$ could have been already estimated
as described before. Same goes for the center-of-mass. Figure 7a depicts the
estimation results in case of noisy measurements, wrong mass estimate, and
correct center-of-mass; only corrections in Eq. (15) are applied. The
estimates converge to the true value, even starting from different initial
conditions. On the other hand, Fig. 7b depicts the estimation of the same
variable, in case the nominal center-of-mass is kept within the digital twin:
as expected, we note a steady-state error in the estimated inertia - mostly
noticeable for $J_{zz}$.

Table 4 quantifies the estimation performance in different cases. One can
notice as the error on $J_{xx}$ is generally higher and grows faster as the
uncertainty increases. This is easily explained as the nominal value and the
variation of $J_{xx}$ are much smaller than for the other inertia. Finally,
Fig. 9 and Fig. 9 show the estimator performance when considering different
noise levels on the employed measurements. As expected, $J_{xx}$ estimation
performance significantly drops, and is never within the $10\%$ limit. On the
other hand, $J_{zz}$ estimation shows good performance up to $snr\approx 7.5$,
which is a reasonable bound for yaw acceleration noise.

Table 4: Roll and yaw inertia estimation performance in different conditions. Test conditions | $rms_{J_{xx}}\ [\%]$ | $rms_{J_{zz}}\ [\%]$
---|---|---
| Noiseless data ($snr=\infty$);
---
Correct mass;
Correct CM.
$0.38$ | $0.96$
| Noiseless data ($snr=\infty$);
---
Wrong mass;
Correct CM.
$4.22$ | $0.25$
| Noisy data ($snr=10$);
---
Correct mass;
Correct CM.
$23.91$ | $6.18$
| Noisy data ($snr=10$);
---
Wrong mass;
Correct CM.
$26.50$ | $5.38$
| Noisy data ($snr=10$);
---
Wrong mass;
Wrong CM.
$45.43$ | $26.50$

###  4.4 Pitch inertia estimation

Figure 10: Pitch inertia estimation in case the road is known without
uncertainty, with noisy measurements.

Figure 11: Pitch inertia estimation with varying noise levels. Oscillation in
$rms_{\%}$ is intrinsically due to the complex nature of the Twin-in-the-Loop
observer.

(a) Sensitivity to $J_{yy}$ estimation error, based on the amount of noise in
the road profile. (b) Noisy versus real road profile, in case
$\varepsilon_{z}=0.1$.

Figure 12: Analysis on the effect of road profile noise onto the estimation of
$J_{yy}$. (a) shows the sensitivity to $J_{yy}$ estimation error for different
noise levels, while (b) compares noisy and real measurements for a specific
noise level.

As done above, pitch inertia $J_{yy}$ can be estimated by suitably exciting
the corresponding dynamics. The most straightforward way of exciting pitch
dynamics is the excitation of the vehicle body with a non-zero road profile:
hence, we assume that the vehicle is being driven at constant speed on a paved
road. The road profile at front wheels is modeled as [16]

$z_{f}=\int_{0}^{t}\eta_{z},\ \eta_{z}\sim WN(0,\sigma_{z}^{2}),$ (16)

Where $\sigma_{z}=0.01\ m$, and $WN$ is the realization of a white noise.
Consequently, the profile on the rear wheels is equal to
$z_{r}(t)=z_{f}(t-v/wb)$, with $v$ the vehicle longitudinal speed, and $wb$
the wheelbase.

The real vehicle, whose load is unknown and has to be estimated, is fed with
the generated road profile. The digital twin is then fed with the same road
profile. In pratice, we are considering that the road profile at the current
time is known; this information can be easily extracted _e.g._ via Kalman
Filtering [11] or sensing devices [19]. To achieve a realistic simulation, we
assume the presence of noise in such information

$z_{f,k}^{n}=z_{f,k}\cdot u_{z,k},\
u_{z,k}\sim\mathcal{U}\left(1-\varepsilon_{z},1+\varepsilon_{z}\right).$ (17)

Equation (17) models a multiplicative uniformly distributed noise onto
$z_{f}$.
The correction law to be applied is derived as usual

$\Delta\delta J_{yy,k_{k}}=K_{\dot{\omega}_{y}-\delta
J_{yy}}(\dot{\omega}_{y,k}-\tilde{\dot{\omega}}_{y,k})\cdot\textrm{sgn}(\dot{\omega}_{y,k}),$
(18)

Where the $y$-axis angular acceleration $\dot{\omega}_{y}$ is employed as a
measurement. In the following, only correction in Eq. (18) is applied.

As also discussed for roll and yaw inertia, a preliminary estimation of the
vehicle mass can yield benefits also to inertia estimation; the same goes for
the center-of-mass position. The inertia estimation - considering noisy
measurements - are showed in Fig. 11. The estimate converges to the correct
value even by starting from different initial conditions. Also, note that the
observer is robust with respect to imperfect knowledge of other inertia
parameters $J_{xx}$ and $J_{zz}$, which are not being corrected.

Then, Table 5 quantifies the estimation performance for different simulation
layouts: even in the worst scenario, _i.e._ in case of noisy measurements and
road profile information, and wrong knowledge of mass and center-of-mass
location, the percentage $rms$ settles at $\approx 10\%$.

Finally, Fig. 12a shows the variation of $rms_{\%}$ for increasing road
profile noise ($\varepsilon_{z}$): a significant performance drop arises for
$\varepsilon_{z}>0.1$.

For a better understanding of the introduced noise level, Fig. 12b depicts the
real versus noisy road profiles. The $rms$ among the two signals is $\approx
0.2\ cm$, which is comparable to what observed in [11].

Table 5: Pitch inertia estimation performance in different conditions of measurement noise, road profile noise, and with or without the mass estimation. Test conditions | $rms_{J_{yy}}\ [\%]$
---|---
| Noiseless data ($snr=\infty$);
---
Noiseless road profile ($\varepsilon_{z}=0$);
Correct mass;
Correct CM.
$0.11$
| Noiseless data ($snr=\infty$);
---
Noisy road profile ($\varepsilon_{z}=0.1$);
Correct mass;
Correct CM.
$3.19$
| Noisy data ($snr\approx 10$);
---
Noiseless road profile ($\varepsilon_{z}=0$);
Correct mass;
Correct CM
$0.33$
| Noisy data ($snr\approx 10$);
---
Noisy road profile ($\varepsilon_{z}=0.1$);
Correct mass;
Correct CM
$3.24$
| Noisy data ($snr\approx 10$);
---
Noisy road profile ($\varepsilon_{z}=0.1$);
Wrong mass;
Correct CM.
$2.38$
| Noisy data ($snr\approx 10$);
---
Noisy road profile ($\varepsilon_{z}=0.1$);
Correct mass;
Wrong CM.
$5.51$
| Noisy data ($snr\approx 10$);
---
Noisy road profile ($\varepsilon_{z}=0.1$);
Wrong mass;
Wrong CM.
$10.20$

###  4.5 Implementation details

Figure 13: Schematic representation of the operations flow of a possible
implementation of the TiL parameters estimator.

In the Sections above, we showed how the TIL approach can properly estimate
mass and inertia under different noise and uncertainty conditions; this
allowed us to focus on each problem independently. In an actual implementation
of the full parameter estimator, one could consider simple logic rules to
prioritize estimation of one or more variables. A possibility is that of Fig.
13: the dashed arrows represent outputs/input of a given module of the
estimator, while solid arrows describe the flow of operations.

The first stage consists in estimating the vehicle mass. We verified (see
Table 6) as this part of the algorithm is extremely robust to knowledge of the
center-of-mass and other inertial parameters; hence, this module is the first
one to be executed. Estimation during straight driving can be enforced to
avoid too highly dynamic conditions: to do so, a simple threshold on vehicle
yaw-rate is sufficient. The transition to the next module can be based upon a
time threshold: setting a reasonably high threshold ensures that the estimator
has converged to the true value when transitioning.

Now, one could estimate the pitch inertia $J_{yy}$, which has been showed to
be robust to center-of-mass knowledge (Table 5); as discussed, this step
requires knowledge of the road profile. Also in this case, enforcing straight
driving can enhance robustness by avoiding parameter corrections in highly
dynamic conditions.

On the other hand, in order to estimate other inertia parameters $J_{xx},\
J_{zz}$ we need to know the center-of-mass variation beforehand, so as to have
reasonable estimation errors (within the $10\ \%$ bound). A suitable
estimation procedure can thus be applied to solve for this problem; this
enables estimation of roll and yaw inertia, for which a data selection to
isolate cornering events is necessary - a threshold on yaw-rate or lateral
acceleration solves the last point.

##  5 Experimental validation

In the following, we test the parameter estimator onto real data. Data has
been collected in a proving ground with an high-performance car, and provided
by a partner car manufacturer whose name cannot be disclosed for
confidentiality reasons.

The vehicle is provided with a set of sensors, yielding the following
measurements

  * •

Center-of-mass 3D accelerations and angular rates (IMU);

  * •

Four wheel angular rates (encoders).

For validation purposes, a double antenna GPS is used to measure longitudinal
and lateral speed ($v_{x}$, $v_{y}$). Consequently, the vehicle sideslip
$\beta$ can be computed as $\beta=atan(v_{y}/v_{x})$. We consider joint
estimation of parameters and states: namely, we want to estimate vehicle
sideslip and mass. Hence, the augmented state vector is

$x^{aug}=\begin{bmatrix}x^{T}&\delta M\end{bmatrix}^{T},$ (19)

where $x\in\mathbb{R}^{28\times 1}$ is the set of simulator states.
Measurements vector $y\in\mathbb{R}^{10\times 1}$ reads

$y=\begin{bmatrix}a_{x}&a_{y}&a_{z}&\omega_{x}&\omega_{y}&\omega_{z}&\omega_{fl}&\omega_{fr}&\omega_{rl}&\omega_{rr}\end{bmatrix}^{T}$
(20)

The most relevant vehicle nominal parameters are collected in Table 6.

Parameter
---
$M_{tot}$ | $J_{xx}$ | $J_{yy}$ | $J_{zz}$ |  $d_{x}$ [$cm$] |  $d_{y}$ [$cm$] |  $d_{z}$ [$cm$]
$1391.29$ | $328.281$ | $1698.21$ | $1864.95$ | $157.14$ | $0$ | $46.22$
Table 6: Test vehicle nominal parameters, as modeled on the digital twin.

###  5.1 Tuning the TiL observer

1 Select an experiment, consisting of measured states, input, measurements
$\left(x_{k},u_{k},y_{k}\right)$

2 Select $\mathcal{N}$ number of total iterations, $n_{seed}$ number of
initial points ($n_{seed}<\mathcal{N}$)

3 Evaluate the objective function in $n_{seed}$ random initial points

4 $i\leftarrow n_{seed}$

5 while  _$i <\mathcal{N}$_ do

6 Update surrogate function $\hat{f}_{i}(\theta)$

7 Compute acquisition function $a_{i}(\theta)$

8 Next point to evaluate is $\theta_{i+1}=\arg\min_{{\theta}}a_{i}(\theta)$

9 $i\leftarrow i+1$

10 Evaluate $y_{i}=f(\theta_{i})$

11

12 end while

13Update surrogate function $\hat{f}_{N}(\theta)$

14 return _1) best evaluated point
$\overline{\theta}=\arg\min_{i}f(\theta_{i})$_

2) best predicted feasible point
$\hat{\overline{\theta}}=\arg\min_{\theta}\hat{f}_{\mathcal{N}}(\theta)$

Algorithm 1 BO for TiL calibration: pseudo-code

Applying the correction in Eq. (5) with set of states and measurements defined
above would yield to state-output mapping matrix $K\in\mathbb{R}^{28\times
10}$. With regards to the non-linear correction $\mathcal{K}$, to be applied
for the estimation, the expression in Eq. (14) is used. All considered
literature solutions rely on Kalman Filters to compute the correction gain to
be applied: unfortunately, given that the digital twin is a black-box, we
cannot access its equations, and thus we cannot employ standard techniques -
_e.g._ the system linearization, as in Extended Kalman Filter. A data-driven
approach is instead here proposed, following what showed in [13]. An
optimization procedure is performed based on experimental data: a set of
states, driver inputs and measurements is necessary $\langle
x_{k},u_{k},y_{k}\rangle$. Such procedure calibrates the parameters off-line
by minimizing a certain cost function - depending on the estimator objective.

Given that the cost metric to be optimized cannot be written in closed-form -
being the simulator equation unknown - we need to rely on a black-box
optimizer. These methods suffer of huge scalability problems, and the number
of optimization variables shall be reduced as much as possible; this can be
enforced by removing variables from the measurements and states vector,
reducing the matrix order. Specifically, we consider the following subset of
$y$

$y_{ss}=\begin{bmatrix}a_{x}&a_{y}&\omega_{y}&\omega_{z}&\omega_{fl}&\omega_{fr}&\omega_{rl}&\omega_{rr}\end{bmatrix}^{T}.$
(21)

$a_{z}$ and $\omega_{x}$ have been removed from the set, as we are interested
in estimating planar dynamics; on the other hand, one can verify that
including the pitch rate $\omega_{y}$ among the measurements can enhance the
estimation performance, as the variable is highly related to the road profile
- and thus to suspension forces.

The same considerations can be applied to the state vector, thus obtaining a
reduced version

$x_{ss}=\begin{bmatrix}v_{x}&v_{y}&\omega_{y}&\omega_{z}&\omega_{fl}&\omega_{fr}&\omega_{rl}&\omega_{rr}&\delta
M\end{bmatrix}^{T}.$ (22)

Where the planar velocities allow to algebraically estimate $\beta$; wheel,
pitch and yaw angular rates $\omega_{fl,fr,rl,rr},\ \omega_{y},\ \omega_{z}$
are corrected with the corresponding measurements in order to guarantee
simulator stability. At this point, mapping matrix has dimensions $9\times 8$,
meaning that $72$ parameters shall be calibrated. $72$ parameters are still
computationally intractable, hence, we promote physics-inspired sparsity in
the matrix [1]. In details, we consider the following corrections

  * •

$K_{\omega-\omega}$. Correction onto the wheel angular rates via the
corresponding measurements;

  * •

$K_{a_{x}-v_{x}}$ and $K_{a_{y}\rightarrow v_{y}}$. Correction onto the
vehicle longitudinal/lateral speed via the corresponding accelerations;

  * •

$K_{\omega_{z}-\omega_{z}}$. Correction onto the yaw rate via the
corresponding measurement;

  * •

$K_{\omega_{y}-\omega_{y}}$. Correction onto the pitch rate via the
corresponding measurement.

Hence, the set of parameters to be tuned reads

$\theta=\begin{bmatrix}K_{\omega-\omega}&K_{a_{x}-v_{x}}&K_{a_{y}-v_{y}}&K_{\omega_{y}-\omega_{y}}&K_{\omega_{z}-\omega_{z}}&K_{a_{x}-\delta
M}\end{bmatrix}$ (23)

Where the term $K_{a_{x}-\delta M}$ is the correction on the vehicle mass; we
apply it through the switching law

$\begin{split}\Delta\delta M_{k}&=K_{a_{x}-\delta
M}\cdot\epsilon\left(a_{x,k},\omega_{z,k}\right)\cdot(a_{x,k}-\tilde{a}_{x,k}),\\\
\epsilon&=\begin{array}[]{c}\begin{cases}0&,\
\left|\omega_{z,k}\right|>\bar{\omega}_{z},\\\ 1&,\
\textrm{sgn}\left(a_{x,k}\right)\geq
0\wedge\left|\omega_{z,k}\right|<\bar{\omega}_{z},\\\ -1&,\
\textrm{sgn}\left(a_{x,k}\right)<0\wedge\left|\omega_{z,k}\right|<\bar{\omega}_{z}.\end{cases}\end{array}\end{split}$
(24)

Equation (24) is very similar to Eq. (14): however, we introduce a scheduling
law based on the yaw-rate $\omega_{z}$, with $\bar{\omega}_{z}=3\ deg/s$ \- as
also discussed in Section 4.5 \- in order to perform the estimation during
straight driving. Now, the optimal observer calibration $\theta^{*}$ can be
found by solving the following optimization problem

$\begin{split}\min_{{\theta}}\ &J\left(\theta\right)\\\ \textrm{subject to}&\
\theta\subseteq\Theta.\end{split}$ (25)

The cost function has to be selected based on the estimation targets. In this
case, we have

$J(\theta)=\sqrt{\dfrac{1}{N_{s}}\sum_{i=1}^{N_{s}}k_{\beta}\left(\beta_{k}-\hat{\beta}_{k}\right)^{2}}+\sqrt{\dfrac{1}{10f_{s}}\sum_{i=N_{s}-10f_{s}}^{N_{s}}\left(\delta
M_{k}-\hat{\delta M}_{k}\right)^{2}}.$ (26)

Whereas $k_{\beta}=100$ accounts for numeric differences among the two cost
function terms, and $N_{s}$ is the number of samples in the experiment. Note
that the second term of the cost function weights the mass in the last $10$
seconds of experiment. This is because we don’t need and don’t expect the mass
estimation to be as fast as the sideslip one, hence we don’t give great
relevance to the first part of experiment.

We select Bayesian Optimization (BO) to solve the problem in Eq. (25): BO is a
model-free optimizer employed when the cost function to be minimized cannot be
written down in form of equation, but we are able to evaluate for some values
of $\theta$. A Gaussian Process proxy of the cost function is then found by
regression on said points, and used to estimate the global optimum of the
problem in Eq. (25). The BO procedure is schematized in Algorithm 1, and more
information can be found in one of the many papers about this topic, _e.g._
[3].

###  5.2 Benchmark estimator

Figure 14: Double-track vehicle model scheme.

In order to fairly compare the TiL estimator with a benchmark, we build an
observer with a double-track planar model as a system replica (see Fig. 14);
the same correction architecture as for the TiL, based on constant gains, is
retained. This model has been frequently used in the state and parameter
estimation literature ([7, 27]). Furthermore, we use double-track model for
consistency to our previous research on TiL topic [13].
The planar vehicle dynamics are written - in discrete time - as

$\begin{split}v_{x,k+1}&=v_{x,k}+T_{s}\left(\dfrac{F_{x}^{T}}{M_{tot}+\delta
M_{k}}+v_{y,k}\omega_{z,k}\right),\\\
v_{y,k+1}&=v_{y,k}+T_{s}\left(\dfrac{F_{y}^{T}}{M_{tot}+\delta
M_{k}}-v_{x,k}\omega_{z,k}\right),\\\
\omega_{z,k+1}&=\omega_{z,k}+T_{s}\left(\dfrac{M_{z}^{T}}{J_{zz}}\right),\\\
\delta M_{k+1}&=\delta M_{k}.\end{split}$ (27)

Where $T_{s}$ is the sampling time. The equations are enhanced with the
extended state for the mass estimation, $\delta M$. The total forces
$F_{x}^{T},\ F_{y}^{T}$ and moment $M_{z}^{T}$ are respectively written as a
function of the tire-road forces

$\begin{split}F_{x}^{T}&=\left(F_{x,fl}+F_{x,fr}\right)\cos(\delta)-(F_{y,fl}+F_{y,fr})\sin(\delta)+\\\
&F_{x,rl}+F_{x,rr},\\\
F_{y}^{T}&=\left(F_{x,fl}+F_{x,fr}\right)\sin(\delta)+(F_{y,fl}+F_{y,fr})\cos(\delta)+\\\
&F_{y,rl}+F_{y,rr},\\\
M_{z}^{T}&=l_{f}\left(F_{y,fl}+F_{y,fr}\right)\cos(\delta)+\dfrac{t}{2}(F_{y,fl}-F_{y,fr})\sin(\delta)-\\\
&\dfrac{t}{2}(F_{x,fl}-F_{x,fr})\cos(\delta)+l_{f}\left(F_{x,fl}+F_{x,fr}\right)\sin(\delta)-\\\
&l_{r}\left(F_{y,rl}+F_{y,rr}\right)-\dfrac{t}{2}\left(F_{x,rl}-F_{x,rr}\right).\end{split}$
(28)

Where $l_{f}$, $l_{r}$, $t$ are parameters influenced by the center-of-mass
position, which could thus change over time. However, given that the
considered case study is a 2-seats car, the expected load variations are not
of significant influence to the CM position.

The wheel forces $F_{x,ij},\ F_{y,ij}$ are modeled through a simplified
Pacejka model [10]

$\begin{split}F_{x}^{ij}&=F_{z,ij}D_{x}\sin\left(C_{x}\arctan\left(B_{x}\lambda_{ij}-E_{x}\left(B_{x}\lambda_{ij}-\arctan\left(B_{x}\lambda_{ij}\right)\right)\right)\right),\\\
F_{y}^{ij}&=F_{z,ij}D_{y}\sin\left(C_{y}\arctan\left(B_{y}\alpha_{ij}-E_{y}\left(B_{y}\alpha_{ij}-\arctan\left(B_{y}\alpha_{ij}\right)\right)\right)\right).\end{split}$
(29)

Where $F_{z}^{ij}$ is the normal force at each wheel; the same can be easily
estimated by using center-of-mass measured accelerations [22]. The wheel
longitudinal and lateral slips $\lambda_{ij},\ \alpha_{ij}$ are computed as

$\begin{split}\lambda_{ij}&=\dfrac{R_{ij}\omega_{ij}-v_{x,ij}}{\max\left(R_{ij}\omega_{ij},v_{x,ij}\right)},\\\
\alpha_{ij}&=\arctan\left(v_{y}^{ij}/v_{x}^{ij}\right).\end{split}$ (30)

Where $v_{x,ij},\ v_{y,ij}$ are obtained from simple kinematic considerations.
In principle, one could include the wheel dynamics as in [13], in order to
filter the wheel angular rate states. However, since the estimation target is
different, and we assume the angular rate measurements to be reliable, we
directly consider measured $\omega_{ij}$ in Eq. (30). This also avoids us to
model the braking and traction torques effect onto the wheels, a potential
further source of error. Overall, Eq. (27), Eq. (29), Eq. (30) can be combined
in order to obtain the double-track dynamics

$\begin{split}x_{bench,k+1}&=x_{bench,k}+T_{s}\cdot
f_{bench}\left(x_{bench,k},u_{bench,k}\right),\\\
y_{bench,k}&=\begin{bmatrix}a_{x,k}\\\ a_{y,k}\\\
\omega_{z,k}\end{bmatrix}=\begin{bmatrix}\dfrac{F_{x}^{T}}{M_{tot}+\delta
M_{k}}+v_{y,k}\omega_{z,k}\\\ \dfrac{F_{y}^{T}}{M_{tot}+\delta
M_{k}}-v_{x,k}\omega_{z,k}\\\ \omega_{z,k}\end{bmatrix}.\end{split}$ (31)

Whereas $x_{bench}=\begin{bmatrix}v_{x}&v_{y}&\dot{\psi}&\delta
M\end{bmatrix}^{T}\in\mathbb{R}^{4}$, and
$u_{bench,k}=\begin{bmatrix}\delta&\omega_{fl}&\omega_{fr}&\omega_{rl}&\omega_{rr}\end{bmatrix}\in\mathbb{R}^{5}$.
$f_{bench}$ is a suitable non-linear function.
We apply on the benchmark estimator the same correction law proposed in
Section 2. The correction gains are stored in vector $\theta_{bench}$

$\theta_{bench}=\begin{bmatrix}K_{a_{x}-v_{x}}&K_{\omega_{z}-\omega_{z}}&K_{a_{x}-\delta
M}&K_{a_{y}-v_{y}}\end{bmatrix}.$ (32)

Bayesian Optimization is employed to tune $\theta_{bench}$, as for the TiL
observer.

######  Remark.

Let us remark that the double-track model - especially the Pacejka parameters
- need to be identified based on experimental data: this increases the number
of tunable parameters and the calibration complexity of the benchmark. On the
other hand, a digital twin of the vehicle is usually already well calibrated,
as the same is used for simulation purposes by the car manufacturer.
Furthermore, the double-track model is not valid anymore if we want to
estimate roll or pitch inertia: we would need to design another ad-hoc model,
or to enhance the double-track with non-planar dynamics.

###  5.3 Final results

Figure 15: Speed, yaw-rate and lateral acceleration profiles in a series of
circuit laps (optimization experiment).

Figure 16: Speed, yaw-rate and lateral acceleration profiles in a series of
double-lane-change maneuvers (validation experiment).

Figure 17: Joint mass and sideslip estimation in a series of circuit laps (optimization experiment). The TIL estimator is compared with a benchmark estimator. Figure 18: Joint mass and sideslip estimation in a series of circuit laps (optimization experiment). The TIL estimator is compared with a benchmark estimator - highlighted portion. Figure 19: Joint mass and sideslip estimation in a series of circuit laps (validation experiment). The TIL estimator is compared with a benchmark estimator. Parameter | Lower bound | Upper bound | Optimized - TiL | Optimized - bench
---|---|---|---|---
$k_{\omega-\omega}$ | $0$ | $1.5$ | $0.259$ | /
$k_{a_{x}-v_{x}}$ | $0.01$ | $0.01$ | $3.287\cdot 10^{-3}$ | $9.282\cdot 10^{-2}$
$k_{a_{y}-v_{y}}$ | $0.01$ | $0.01$ | $7.430\cdot 10^{-2}$ | $5.271\cdot 10^{-2}$
$k_{\omega_{y}-\omega_{y}}$ | $0.2$ | $1.5$ | $0.913$ | /
$k_{\omega_{z}-\omega_{z}}$ | $0.2$ | $1.5$ | $1.046$ | $0.314$
$k_{a_{x}-\delta M}$ | $0$ | $1000$ | $248.17$ | $2.06$
Table 7: Upper, lower bounds and optimal parameters found for the TiL
estimator via BO.

The following results are obtained by testing the TiL estimator and the
benchmark defined in Section 5.2. The experiments have been performed on the
vehicle having a nominal mass, hence, in order to test the algorithms, we
modify the initialization of the estimator internal models, so to have a wrong
mass value $\delta_{M,0}=-350\ kg$.

Figure 16 depicts vehicle speed, steering wheel angle and lateral acceleration
in the experiment used for solving optimization problem in Eq. (25) - the test
consists in a series of circuit laps. As the reader can note, the experiment
is extremely challenging, in that accelerations up to $1.5\ g$ and speed up to
$230\ km/h$ are reached.

For this test, we apply the BO procedure solving the problem of Eq. (25): the
optimal parameters in Table 7 are obtained. The variable bounds can be found
via trial-and-error and prior physical knowledge, as discussed in [13] The
estimation results are displayed in Fig. 17. An highlight on the last part of
the test is showed in Fig. 18. As one can notice, TiL is able to outperform
the double-track based observer both on mass and sideslip estimation. The mass
estimate convergence is indeed slower than in the simulation test of Section
4.2: this is however expected, as the real world experiment is significantly
more challenging.

Figure 16 shows speed, yaw-rate and lateral acceleration for another
experiment, used for validating the estimators - this test consists in a
series of double lane change and braking maneuvers, and is thus different in
nature from the optimization one. The estimation results are given in Fig. 19.
Also in this case, both filters are able to properly estimate the sideslip
angle - with a comparable performance. On the other hand, the double-track
model is not able to correctly estimate the mass, and seems to be diverging
over time.

##  6 Conclusions

In this manuscript, we show that the mass and inertia of a vehicle can be
estimated - jointly with classical vehicle dynamics states - via the recently
introduced Twin-in-the-Loop filtering approach. The method is validated
extensively in simulation, considering independent estimation of each
parameter of interest, and performing sensitivity analyses to noise and
uncertainty. Then, the TiL estimator is tested against experimental data
collected on a high-performance car, showing high performance in joint
estimation of mass and sideslip. A comparison with another observer based on
state-of-the-art double-track vehicle modelling provides further insights into
the potential of the TiL architecture.

Future work will be dedicated to the study of the real-time estimation of the
tire-road friction coefficient.

## Conflicts of Interest

The authors declare no conflict of interest.

## References

  * Delcaro et al. [2023] Delcaro, G., Dettù, F., Formentin, S., Savaresi, S., 2023\.  Dealing with the curse of dimensionality in twin-in-the-loop observer design, in: 2023 IFAC World Congress, IFAC.
  * Dettù et al. [2022] Dettù, F., Formentin, S., Savaresi, S.M., 2022.  The twin-in-the-loop approach for vehicle dynamics control.  arXiv:2209.02263 .
  * Gelbart et al. [2014] Gelbart, M.A., Snoek, J., Adams, R.P., 2014.  Bayesian optimization with unknown constraints.  arXiv preprint arXiv:1403.5607 .
  * Gevers [2005] Gevers, M., 2005.  Identification for control: From the early achievements to the revival of experiment design.  European journal of control 11, 335–352.
  * Gobbi et al. [2011] Gobbi, M., Mastinu, G., Previati, G., 2011.  A method for measuring the inertia properties of rigid bodies.  Mechanical Systems and Signal Processing 25, 305–318.
  * Gong et al. [2020] Gong, X., Suh, J., Lin, C., 2020.  A novel method for identifying inertial parameters of electric vehicles based on the dual h infinity filter.  Vehicle System Dynamics 58, 28–48.
  * Hong et al. [2014] Hong, S., Lee, C., Borrelli, F., Hedrick, J.K., 2014\.  A novel approach for vehicle inertial parameter identification using a dual kalman filter.  IEEE Transactions on Intelligent Transportation Systems 16, 151–161.
  * Huang and Wang [2014] Huang, X., Wang, J., 2014.  Real-time estimation of center of gravity position for lightweight vehicles using combined akf–ekf method.  IEEE Transactions on Vehicular Technology 63, 4221–4231.  doi:10.1109/TVT.2014.2312195.
  * Li et al. [2016] Li, L., Lu, Y., Wang, R., Chen, J., 2016.  A three-dimensional dynamics control framework of vehicle lateral stability and rollover prevention via active braking with mpc.  IEEE Transactions on Industrial Electronics 64, 3389–3401.
  * Pacejka [2005] Pacejka, H., 2005.  Tire and vehicle dynamics.  Elsevier.
  * Qin et al. [2017] Qin, Y., Langari, R., Wang, Z., Xiang, C., Dong, M., 2017.  Road profile estimation for semi-active suspension using an adaptive kalman filter and an adaptive super-twisting observer, in: 2017 American Control Conference (ACC), IEEE. pp. 973–978.
  * Reina et al. [2017] Reina, G., Paiano, M., Blanco-Claraco, J.L., 2017.  Vehicle parameter estimation using a model-based estimator.  Mechanical Systems and Signal Processing 87, 227–241.
  * Riva et al. [2022] Riva, G., Formentin, S., Corno, M., Savaresi, S.M., 2022\.  Twin-in-the-loop state estimation for vehicle dynamics control: theory and experiments.  arXiv:2204.06259 .
  * Rodríguez et al. [2021] Rodríguez, A.J., Sanjurjo, E., Pastorino, R., Naya, M.Á., 2021\.  State, parameter and input observers based on multibody models and kalman filters for vehicle dynamics.  Mechanical Systems and Signal Processing 155, 107544.
  * Rozyn and Zhang [2010] Rozyn, M., Zhang, N., 2010.  A method for estimation of vehicle inertial parameters.  Vehicle system dynamics 48, 547–565.
  * Savaresi et al. [2010] Savaresi, S.M., Poussot-Vassal, C., Spelta, C., Sename, O., Dugard, L., 2010.  Semi-active suspension control design for vehicles.  Elsevier.
  * Solmaz et al. [2008] Solmaz, S., Akar, M., Shorten, R., Kalkkuhl, J., 2008\.  Real-time multiple-model estimation of centre of gravity position in automotive vehicles.  Vehicle System Dynamics 46, 763–788.
  * Tavernini et al. [2020] Tavernini, D., Vacca, F., Metzler, M., Savitski, D., Ivanov, V., Gruber, P., Hartavi, A.E., Dhaens, M., Sorniotti, A., 2020.  An explicit nonlinear model predictive abs controller for electro-hydraulic braking systems.  IEEE Transactions on Industrial Electronics 67, 3990–4001.  doi:10.1109/TIE.2019.2916387.
  * Theunissen et al. [2020] Theunissen, J., Sorniotti, A., Gruber, P., Fallah, S., Ricco, M., Kvasnica, M., Dhaens, M., 2020.  Regionless explicit model predictive control of active suspension systems with preview.  IEEE Transactions on Industrial Electronics 67, 4877–4888.
  * Thornton and Marion [2004] Thornton, S., Marion, J., 2004\.  Classical Dynamics of Particles and Systems.  Brooks/Cole.
  * VI-Grade [2022] VI-Grade, 2022.  Vi-carrealtime.  https://www.vi-grade.com/en/products/vi-carrealtime/.  [Online; accessed 23-March-2022].
  * Viehweger et al. [2021] Viehweger, M., Vaseur, C., van Aalst, S., Acosta, M., Regolin, E., Alatorre, A., Desmet, W., Naets, F., Ivanov, V., Ferrara, A., et al., 2021\.  Vehicle state and tyre force estimation: demonstrations and guidelines.  Vehicle system dynamics 59, 675–702.
  * Wang et al. [2020] Wang, C., Wang, Z., Zhang, L., Cao, D., Dorrell, D.G., 2020\.  A vehicle rollover evaluation system based on enabling state and parameter estimation.  IEEE Transactions on Industrial Informatics 17, 4003–4013.
  * Wenzel et al. [2006] Wenzel, T.A., Burnham, K., Blundell, M., Williams, R., 2006\.  Dual extended kalman filter for vehicle state and parameter estimation.  Vehicle system dynamics 44, 153–171.
  * Wielitzka et al. [2015] Wielitzka, M., Dagen, M., Ortmaier, T., 2015.  Joint unscented kalman filter for state and parameter estimation in vehicle dynamics, in: 2015 IEEE Conference on Control Applications (CCA), IEEE. pp. 1945–1950.
  * Yang et al. [2008] Yang, S., Liu, T., Cheng, Y., 2008.  Automatic measurement of payload for heavy vehicles using strain gages.  Measurement 41, 491–502.
  * Zhu et al. [2019] Zhu, J., Wang, Z., Zhang, L., Zhang, W., 2019.  State and parameter estimation based on a modified particle filter for an in-wheel-motor-drive electric vehicle.  Mechanism and Machine Theory 133, 606–624.