A Novel Kalman Filter Design and Analysis Method Considering Observability and Dominance Properties of Measurands Applied to Vehicle State Estimation

In Kalman filter design, the filter algorithm and prediction model design are the most discussed topics in research. Another fundamental but less investigated issue is the careful selection of measurands and their contribution to the estimation problem. This is often done purely on the basis of empirical values or by experiments. This paper presents a novel holistic method to design and assess Kalman filters in an automated way and to perform their analysis based on quantifiable parameters. The optimal filter parameters are computed with the help of a nonlinear optimization algorithm. To determine and analyze an optimal filter design, two novel quantitative nonlinear observability measures are presented along with a method to quantify the dominance contribution of a measurand to an estimate. As a result, different filter configurations can be specifically investigated and compared with respect to the selection of measurands and their influence on the estimation. An unscented Kalman filter algorithm is used to demonstrate the method’s capabilities to design and analyze the estimation problem parameters. For this purpose, an example of a vehicle state estimation with a focus on the tire-road friction coefficient is used, which represents a challenging problem for classical analysis and filter parameterization.


Introduction
The degree of automation and technical support for humans has increased rapidly in recent years. The basic requirement for any control system is the existence of measurable control variables. If they cannot be measured directly due to technical or economic reasons, state estimators are needed. Thus, state estimation forms the backbone of most modern control problems and is required for their implementation. One proven method for state estimation is Kalman filtering. The Kalman filter is an algorithm that provides optimal estimates for the states of a dynamical system sequentially in time. The disturbances in the underlying mathematical model of the system and measurement equations are assumed to be white noise [1]. This established and widely used method has been known and applied for more than 60 years now. Many different modifications of this method have been developed. Nevertheless, there is still no simple procedure for an optimal design or parametrization of a Kalman filter. This task is often executed by experienced control engineers based on empirical knowledge or by trial and error experiments through Monte-Carlo simulations, see, e.g., [2].
One of the most fundamental requirements for an estimation problem is its observability. This means that the estimated states, which are reconstructed, have to be contained in the measurements and uniquely extracted from them. In order to parametrize a filter, the covariance matrices need to be determined (when dealing with other filter types such as, e.g., an unscented Kalman filter, there might be even more parameters). This leads to the crucial question for the filter design as to which measurands should actually be used and how they influence the estimation problem. In many cases, the properties of the measurands for the estimation problem are only evaluated by simple plausibility checks. For instance, when estimating the state-of-charge (see e.g., [3]) of a battery, it is evident that the battery voltage should be used as a measurand. But how would other measurands affect the observability, and might they be even more dominant? Issues like these are addressed in this paper in a quantifiable way with the help of the novel design and analysis method presented below.

State of the Art
The authors of [4] present an automated method that allows the determination of proper Kalman filter parameterization. To evaluate the estimation accuracy, a performance index, which is analytically related to the filter parameterization, is introduced. By minimizing this index, optimal parameterization can be calculated. In [5], a relationship between the performance values and filter parameterization is investigated. However, neither an exact relationship between these values nor a calculation rule is presented. A similar relationship is presented in [6] using the example of a video tracking system. In [7], a method is shown where the filter is considered as a control system, thus allowing corresponding tuning criteria to be derived. In [8], the filter parameters are computed via optimization based on a genetic algorithm. The authors in [9] use Bayesian optimization for this purpose. Instead of minimizing a cost function, this approach tries to maximize the probability of improving the current best solution. In [10], a two-step method is presented in which particle swarm optimization (PSO) is used to tune both the filter parameterization and prediction model parameters.
Besides filter parameterization, the filter structure, i.e., the selection of state variables and measurands, is perhaps an even more important design issue. However, the methods mentioned so far cannot help to solve this problem. The so-called Programmable Kalman Filter Design Tool (PKFD) is shown in [11]. This tool provides both an optimal parameterization (system noise, measurement noise, and initial estimation error covariance matrices) and an optimal filter structure. Nevertheless, the tool rather serves as a rapid prototyping environment, allowing different Kalman filter setups to be compared, but not giving a precise analysis of the properties. Figure 1 shows a classification of the Kalman filter design methodologies mentioned above.

Contribution of This Paper
In the approaches mentioned so far, the design focus lies solely on the filter parameterization, but not on its structure. To fill this research gap, both design issues are addressed and evaluated in detail by the holistic method elaborated in this paper.
Using a novel state-specific nonlinear quantitative observability measure, the current "observability accuracy" can be determined at any time in an estimation problem. Additionally, a new "dominance analysis" method allows the percentage contribution of each measurand to an individual state to be calculated. Based on this knowledge about the observability and dominance properties of all possible measurands, different filter setups can be compared with each other. As a result, the benefit of a measurand can be evaluated in comparison with the effort required to provide it in reality. Furthermore, the filter covariance matrices (as well as other possible filter parameters) are determined by a nonlinear optimization algorithm considering the filter self-diagnosis and resulting in a minimum estimation error. Optimal filter structure using a quantitative nonlinear state-specific observability measure and a dominance analysis to evaluate the influence of the measurands' properties on the estimation problem.
This novel method is universally applicable. As an example, this paper considers its application to a vehicle state estimation problem using an unscented Kalman filter (UKF).

Fundamentals of State Estimation
A well-proven method for state estimation is Kalman filtering. Here, the system modeling, as well as the measurements, are described by their statistical characteristics, and an optimal estimation [1] is performed by an iterative procedure (prediction and measurement update). This section closely follows [3,12]; the interested reader is referred to these sources for more detail.

System Description
For many control system tasks, the plant model to be used in state estimation is naturally described as a nonlinear continuous-time state-space system: ̇= ( , ), = ( ), ( ) ∈ ℝ ×1 , ( ) ∈ ℝ n×1 , ( ) ∈ ℝ m×1 , ∈ ℝ. (1) where is the time, ( ) is the vector of inputs, ( ) is the vector of states, and ( ) is the vector of outputs. In a sampled data system (e.g., a microcontroller) the continuoustime model representation in Equation (1) cannot be used directly. Instead, a time-discrete representation is needed, and therefore, the time-discrete transformation of Equation (1) with additive Gaussian noise is used in the sequel: (2)
PKFD-Tool [11] (rather a rapid prototyping tool) In summary, the method enables a holistic filter design and provides quantitative criteria for an optimal filter configuration, namely:

•
Optimal filter parameterization using a nonlinear constrained optimization algorithm.

•
Optimal filter structure using a quantitative nonlinear state-specific observability measure and a dominance analysis to evaluate the influence of the measurands' properties on the estimation problem.
This novel method is universally applicable. As an example, this paper considers its application to a vehicle state estimation problem using an unscented Kalman filter (UKF).

Fundamentals of State Estimation
A well-proven method for state estimation is Kalman filtering. Here, the system modeling, as well as the measurements, are described by their statistical characteristics, and an optimal estimation [1] is performed by an iterative procedure (prediction and measurement update). This section closely follows [3,12]; the interested reader is referred to these sources for more detail.

System Description
For many control system tasks, the plant model to be used in state estimation is naturally described as a nonlinear continuous-time state-space system: . x = f(x, u), y = h(x), u(t) ∈ R s×1 , x(t) ∈ R n×1 , y(t) ∈ R m×1 , t ∈ R. (1) where t is the time, u(t) is the vector of inputs, x(t) is the vector of states, and y(t) is the vector of outputs. In a sampled data system (e.g., a microcontroller) the continuous-time model representation in Equation (1)  representation is needed, and therefore, the time-discrete transformation of Equation (1) with additive Gaussian noise is used in the sequel: Here, t k is the k-th sample time instant of a periodically sampled data system with u k = u(t k ), x k = x(t k ) and y k = y(t k ). The vectors w k and v k represent zero biased Gaussian white noise. The covariance matrices Q k and R k are defined as E w k w T Note that as a simplification, Q k is assumed to be a diagonal matrix, but in general, it may contain cross-correlation terms between the states (see, e.g., [13]). The operator E(·) calculates the expected value of a random variable [13]. The notation w k ∼ N(0, Q k ) indicates that w k is a Gaussian random variable with a mean vector of 0 and a covariance matrix of Q k = diag σ 2 1..n x , with the standard deviation σ.

Constrained State Estimation
Similar to an optimization problem, a state estimation problem can be advantageously simplified by introducing constraints on the possible solution space, e.g., by physical limits, thus reducing the complexity of the problem. There are several methods dealing with constraints for state estimators, see [14]. As they are beyond the scope of this paper, they will not be discussed in any further detail here. A summary of the methods with their advantages and disadvantages can be found in [12]. In the present paper, two linear constraints are applied (a lower and an upper bound for the maximum tire-road friction coefficient (TRFC)), see Section 4.1.3. For this purpose, the method described in [12] (p. 79) is used. By means of a root-finding problem, the feasible state variables are determined. For linear constraints, as is the case for this paper, the method provides an optimal solution to the problem.

DLR Kalman Filter Estimation Framework
For this research work, we use and extend the DLR Kalman Filter estimation framework [12], which uses prediction models based on continuous-time Modelica models and automatically generates model-based nonlinear state estimators. The approach is based on an extended FMI (Functional Mockup Interface) 2.0 co-simulation interface [15] that interacts with the state estimation algorithms implemented in the DLR Kalman Filter Library [16]. Starting from a multi-physical Modelica model (continuous-time, usually nonlinear), a nonlinear prediction model is automatically generated in the form of a sampled data system (cf. Equation (2)). The framework employs an intelligent separation of the model (encapsulated in a standardized FMI 2.0 for co-simulation [15]) and the estimation algorithm by utilizing modern computer technologies and recent developments in the Modelica language [17]. They enable automated discretization, integration, and derivative calculation of an object-oriented equation-based prediction model. The FMI defines a standardized interface to be used in computer simulations to develop complex cyber-physical systems. The following estimation algorithms are implemented reliably and efficiently in the DLR Kalman Filter Library: EKF, EKF SR (EKF square-root), EKF UD (EKF UD-decomposition), UKF (unscented Kalman filter), and UKF SR (UKF square-root), see [1,18]. Additionally, there are modified algorithms for parameter estimation as well as an extension to nonlinear moving horizon estimation (MHE) using a fast nonlinear gradient descent search, as is presented in [12]. Recently, the library's features were extended to meet the requirements of embedded targets. This was part of the ITEA EMPHYSIS project in which a new embedded FMI (eFMI) standard specification was designed [19].

Design Method for Kalman Filters
The structural filter design is based on the observability and dominance analysis of measurands. This section starts with the introduction of three different observability measures for nonlinear systems. In addition to the probably best-known one-rank analysis of the observability matrix, two new approaches for a quantitative statement about observability are presented. Next, these three methods are compared by means of an illustrative example. Moreover, the advantages of the two newly developed observability measures are discussed. Afterward, a new approach-the so-called dominance analysis-is described to quantify the contribution of a measurand to the estimation. Finally, the holistic filter design method is shown, which consists of the observability and dominance analysis embedded in an optimization framework.

Nonlinear Observability Measures
The observability problem together with its counterpart, the controllability problem, are important parts of the systems theory. Especially when designing observers, observability plays a key role. In this section, three methods for nonlinear observability analysis are presented. For reasons of clarity, the time dependence is not explicitly indicated by an argument for the respective variables.
The observability of a dynamic system is a property that is independent of the estimation method, but is solely determined by the • structure of the problem, i.e., which measurements are available and how they are linked to the states (measurement equation). • excitation of the system by the input u, which has to take a minimum value (persistence of excitation).

Observability via Rank Condition
Probably the best-known method for the observability analysis of nonlinear systems is the rank investigation of the observability matrix. This matrix is built using the n − 1 Lie derivatives y, . y, . . . y (n−1) [13], whereby y indicates the system output from Equation (1). Hence, the nonlinear observability matrix is: Sensors 2021, 21, x FOR PEER REVIEW 5 of 28 servability are presented. Next, these three methods are compared by means of an illustrative example. Moreover, the advantages of the two newly developed observability measures are discussed. Afterward, a new approach-the so-called dominance analysisis described to quantify the contribution of a measurand to the estimation. Finally, the holistic filter design method is shown, which consists of the observability and dominance analysis embedded in an optimization framework.

Nonlinear Observability Measures
The observability problem together with its counterpart, the controllability problem, are important parts of the systems theory. Especially when designing observers, observability plays a key role. In this section, three methods for nonlinear observability analysis are presented. For reasons of clarity, the time dependence is not explicitly indicated by an argument for the respective variables.
The observability of a dynamic system is a property that is independent of the estimation method, but is solely determined by the  structure of the problem, i.e., which measurements are available and how they are linked to the states (measurement equation).  excitation of the system by the input , which has to take a minimum value (persistence of excitation).

Observability via Rank Condition
Probably the best-known method for the observability analysis of nonlinear systems is the rank investigation of the observability matrix. This matrix is built using the − 1 Lie derivatives ,̇, … ( − ) [13], whereby indicates the system output from Equation (1). Hence, the nonlinear observability matrix is: The system is globally observable if a unique inverse function of Equation (3) can be found, which occurs, however, possible only in very rare cases in practical applications. Instead, Equation (3) making the observability investigation at certain points possible. A linear system of equations is obtained: If is invertible, i.e., has full rank , the states can be obtained via Δ = −1 ⋅ Δ . The nonlinear system in Equation (1) The criterion of a rank loss provides only a binary statement about the observability and is numerically highly sensitive. Quantitative information, regarding how good or bad the observability is, is not given. After the presented classical rank-based observability analysis, the next two sections introduce two new observability measures allowing a quantitative statement about the observability. The system is globally observable if a unique inverse function of Equation (3) can be found, which occurs, however, possible only in very rare cases in practical applications. Instead, Equation (3) can be linearized along the reference trajectories making the observability investigation at certain points possible. A linear system of equations is obtained: The criterion of a rank loss provides only a binary statement about the observability and is numerically highly sensitive. Quantitative information, regarding how good or bad the observability is, is not given. After the presented classical rank-based observability analysis, the next two sections introduce two new observability measures allowing a quantitative statement about the observability.

Quantitative Observability Measure Considering Numerical Condition Number
The crucial question for the observability quantification is how far O M is away from a rank loss. The geometric interpretation of the numerical condition number κ(·) of a matrix provides an answer to this question. Namely, the reciprocal condition number indicates the relative distance (w.r.t. the Euclidean norm) of a non-singular matrix to its nearest singular matrix [20] (p. 242): The "distance to the singularity", which can be quantified via κ(O M ), corresponds to the "distance to a rank loss" of O M and can, therefore, be interpreted as a quantitative observability measure.
The relationship between the numerical condition number and the observability properties can be clearly shown by considering the covariance-error-ellipsoid (a.k.a. confidence ellipsoid) of Equation (5). In Figure 2, this is shown as an example of a system with order n = 2. The crucial question for the observability quantification is how far is away from a rank loss. The geometric interpretation of the numerical condition number (•) of a matrix provides an answer to this question. Namely, the reciprocal condition number indicates the relative distance (w.r.t. the Euclidean norm) of a non-singular matrix to its nearest singular matrix [20] (p. 242): The "distance to the singularity", which can be quantified via ( ), corresponds to the "distance to a rank loss" of and can, therefore, be interpreted as a quantitative observability measure.
The relationship between the numerical condition number and the observability properties can be clearly shown by considering the covariance-error-ellipsoid (a.k.a. confidence ellipsoid) of Equation (5). In Figure 2, this is shown as an example of a system with order = 2.
providing a statement about the shape of the ellipsoid. For = 2, large values of the condition number mean a narrow ellipse. This means that existing small uncertainties of one state lead to large uncertainties of another state, which implies bad observability. A good condition number of , i.e., small values of ( ), is equivalent to good observability. One main disadvantage of the two methods presented above is the high computational effort required for the calculation of the Lie derivatives in the observability matrix.   One main disadvantage of the two methods presented above is the high computational effort required for the calculation of the Lie derivatives in the observability matrix. This effort grows with the numbers of states and outputs. Furthermore, the statements about observability are valid only for the whole system, but not for single states. This means the two presented measures show a loss of observability, even though only a single state is unobservable. Often, not all states are of equal interest, wherefore a state-specific quantitative observability measure is presented hereinafter.

State-Specific Quantitative Observability Measure Using Weighted Least Squares
The basic idea of the state-specific observability measure introduced in this section is assessment of the observability via the weighted least squares (WLS) solution for the states applied to the linearized measurement equation.
If the output in Equation (1) is linearized along the reference states x ref (denoted as the ground truth states; since the observability is independent of the used estimation algorithm and depends only on the underlying system structure and excitations (see Section 3.1), the linearization has to be performed around the reference states, but not around the estimated ones), the measurement sensitivity matrix H ref is obtained as and the linear measurement equation implies Remark: If an extended Kalman filter is used for the estimation, the measurement sensitivity matrix advantageously results from the filter algorithm as a by-product. However, instead of linearizing around the estimated states, the linearization has to be performed around the reference states.
Under the condition that at least as many measurements are available as states, i.e., in case of m ≥ n, and rank(H ref ) = n, Equation (10) is overdetermined and thus solvable. This means that with known measurements z, Equation (10) can be solved directly to find the states using a curve-fitting approach. A suitable method for this is a least-square approach. The measurements need to be weighted by the inverse measurement noise covariance matrix R, so that a weighted least square (WLS) problem formulation [23] can be used, leading to The basic idea of the state-specific observability measure introduced in this section is assessment of the observability via the weighted least squares (WLS) solution for the states applied to the linearized measurement equation.
If the output in Equation (1) is linearized along the reference states ref (denoted as the ground truth states; since the observability is independent of the used estimation algorithm and depends only on the underlying system structure and excitations (see Section 3.1), the linearization has to be performed around the reference states, but not around the estimated ones), the measurement sensitivity matrix ref is obtained as and the linear measurement equation implies Remark: If an extended Kalman filter is used for the estimation, the measurement sensitivity matrix advantageously results from the filter algorithm as a by-product. However, instead of linearizing around the estimated states, the linearization has to be performed around the reference states.
Under the condition that at least as many measurements are available as states, i.e., in case of ≥ , and rank( ref ) = , Equation (10) is overdetermined and thus solvable. This means that with known measurements , Equation (10) can be solved directly to find the states using a curve-fitting approach. A suitable method for this is a least-square approach. The measurements need to be weighted by the inverse measurement noise covariance matrix , so that a weighted least square (WLS) problem formulation [23] can be used, leading to with the matrix obs ∈ ℝ × . Remark: The calculation of obs via a singular value decomposition is computationally demanding (see [22]) and can be calculated by a QR decomposition in a more efficient way [12] (p. 39ff).
represent a quantitative observability measure. The covariance entries obs, 2 indicate the current observability of the ℎ state. Large variance values imply bad observability, while small values imply it is good. Compared to classical rank loss-based approaches, the presented method allows a state-specific quantitative statement about the observability. Furthermore, the covariance entries obs, 2 have the same physical units as the states, (11) with the matrix P obs ∈ R n×n .
Remark: The calculation of P obs via a singular value decomposition is computationally demanding (see [22]) and can be calculated by a QR decomposition in a more efficient way [12] (p. 39ff). P obs is the covariance matrix of the WLS estimator and indicates the impact of the uncertainty of the measurements z on the states x WLS . Exactly this quantification, namely, how well the states can be reconstructed from the measurements, corresponds to the definition of observability. Thus, the diagonal entries of P obs diag(P obs ) = σ 2 obs,x 1 · · · σ 2 obs,x n have the same physical units as the states, making their physical interpretability possible. Due to the simple WLS formulation, the method can be executed quickly and computationally efficiently.

Comparison of Observability Measures
This section compares the three observability measures using an illustrative example of the tire road friction coefficient (TRFC) estimation. For this purpose, the nonlinear two-track model from Section 4.1.2 is exploited (its exact knowledge is not yet required to understand the example). Since this section is only aimed at comparing the observability measures, the reference model corresponds to the filter prediction model. Sinusoidal steering with a constant vehicle speed is simulated as a maneuver. This implies that the system is excited solely by the lateral acceleration a C y . Therefore, observable and non-observable phases can be specifically generated for the TRFC µ max . Figure 3 shows a maneuver range without excitation between t = 7 s and t = 17 s, in which no acceleration, braking, or steering are present. It is obvious that in the area where the vehicle simply rolls straight ahead, the TRFC is not observable. steering with a constant vehicle speed is simulated as a maneuver. This implies that the system is excited solely by the lateral acceleration . Therefore, observable and non-observable phases can be specifically generated for the TRFC max . Figure 3 shows a maneuver range without excitation between = 7 and = 17 , in which no acceleration, braking, or steering are present. It is obvious that in the area where the vehicle simply rolls straight ahead, the TRFC is not observable. At the very top of Figure 3, the excitation by the lateral acceleration is given. The plots below show the TRFC max , its estimation, and the following three observability measures:


The rank of , which returns only a binary yes/no-observability assessment.  The reciprocal numerical condition number of (due to large condition numbers, the logarithm base 10 is taken for reasons of clarity). Small values correspond to poor observability, whereas large values indicate good observability.  The standard deviation of the TRFC obs, max using the WLS approach. Large values imply high uncertainty, i.e., poor observability, and small values vice versa. The specified uncertainty has the same physical unit as the TRFC, making its direct and clear interpretation possible. At the very top of Figure 3, the excitation by the lateral acceleration a C y is given. The plots below show the TRFC µ max , its estimation, and the following three observability measures:

•
The rank of O M , which returns only a binary yes/no-observability assessment.

•
The reciprocal numerical condition number of O M (due to large condition numbers, the logarithm base 10 is taken for reasons of clarity). Small values correspond to poor observability, whereas large values indicate good observability.

•
The standard deviation of the TRFC σ obs,µ max using the WLS approach. Large values imply high uncertainty, i.e., poor observability, and small values vice versa. The specified uncertainty has the same physical unit as the TRFC, making its direct and clear interpretation possible.
The jump of µ max to a low friction value can be correctly estimated shortly after its occurrence at t = 2 s if there is sufficient excitation. From t = 7 s on, there is no excitation anymore, so the observability is lost. Therefore, the estimated TRFC follows its first order lag (PT1) behavior given in the system model and tends towards the set value of µ max = 1.
The rank of the observability matrix decreases until t = 13 s, indicating a loss of observability with a delay of 6 s. The two quantitative measures indicate observability loss immediately after its occurrence. According to the WLS approach, the current uncertainty is σ obs,µ max ≈ 1.2, which corresponds to the total loss of observability. From t = 17 s on, there is again sufficient excitation, so that all three measures again show observability, and the TRFC can be estimated correctly. In the two quantitative measures, it can also be clearly seen that at each zero crossing of a C y , the observability becomes poor, as in this case, there is no excitation for a short time.
This example shows the advantages of the two introduced novel observability measures. While the well-known rank loss-based criterion indicates the loss of observability with insufficient accuracy, the measures κ(O M ) and σ obs,x immediately provide a quantified statement about how good or bad the current observability is. The condition number κ(O M ) only provides a statement about the whole system. However, in many cases, some states are observable, while others are not. Since the focus in this paper is directed towards the TRFC, i.e., the state µ max , the state-specific observability criterion σ obs,x i is used below.

Dominance Measure: Individual Contribution of Measurands to Estimated States
In the structural design of a Kalman filter, an essential question is which measurands contribute to the estimation of a state and how much. A simple and practical method to determine this contribution can be realized by considering the equation for the filter measurement update: clearly seen that at each zero crossing of , the observability becomes poor, as in this case, there is no excitation for a short time.
This example shows the advantages of the two introduced novel observability measures. While the well-known rank loss-based criterion indicates the loss of observability with insufficient accuracy, the measures ( ) and obs, immediately provide a quantified statement about how good or bad the current observability is. The condition number ( ) only provides a statement about the whole system. However, in many cases, some states are observable, while others are not. Since the focus in this paper is directed towards the TRFC, i.e., the state max , the state-specific observability criterion obs, is used below.

Dominance Measure: Individual Contribution of Measurands to Estimated States
In the structural design of a Kalman filter, an essential question is which measurands contribute to the estimation of a state and how much. A simple and practical method to determine this contribution can be realized by considering the equation for the filter measurement update: where ∈ ℝ ×1 is called the innovation. The second summand ∈ ℝ ×1 of Equation (13) indicates an incremental contribution of each measurand to the corresponding state: The row vectors ∈ ℝ 1× , ∈ {1, … , } of the Kalman gain matrix ∈ ℝ × contain the gains of the innovation for the ℎ state. Through the element-wise multiplication the vector , ∈ ℝ 1× is obtained, whose elements are the proportion of the measurements in the estimation of the ℎ state. For the investigation of a discrete time interval ∈ [ 0 , ], the absolute values of the measurements can be summed up. For better comparability, normalization is performed to quantify the relative shares in the overall estimation: Analysis of the measurands' influence makes sense only for time points in which the state is observable. Due to the indicator function (•) defined as where d k ∈ R m×1 is called the innovation. The second summand m k ∈ R n×1 of Equation (13) indicates an incremental contribution of each measurand to the corresponding state: Sensors 2021, 21, x FOR PEER REVIEW 9 of 28 clearly seen that at each zero crossing of , the observability becomes poor, as in this case, there is no excitation for a short time.
This example shows the advantages of the two introduced novel observability measures. While the well-known rank loss-based criterion indicates the loss of observability with insufficient accuracy, the measures ( ) and obs, immediately provide a quantified statement about how good or bad the current observability is. The condition number ( ) only provides a statement about the whole system. However, in many cases, some states are observable, while others are not. Since the focus in this paper is directed towards the TRFC, i.e., the state max , the state-specific observability criterion obs, is used below.

Dominance Measure: Individual Contribution of Measurands to Estimated States
In the structural design of a Kalman filter, an essential question is which measurands contribute to the estimation of a state and how much. A simple and practical method to determine this contribution can be realized by considering the equation for the filter measurement update: where ∈ ℝ ×1 is called the innovation. The second summand ∈ ℝ ×1 of Equation (13) indicates an incremental contribution of each measurand to the corresponding state: The row vectors ∈ ℝ 1× , ∈ {1, … , } of the Kalman gain matrix ∈ ℝ × contain the gains of the innovation for the ℎ state. Through the element-wise multiplication the vector , ∈ ℝ 1× is obtained, whose elements are the proportion of the measurements in the estimation of the ℎ state. For the investigation of a discrete time interval ∈ [ 0 , ], the absolute values of the measurements can be summed up. For better comparability, normalization is performed to quantify the relative shares in the overall estimation: Analysis of the measurands' influence makes sense only for time points in which the state is observable. Due to the indicator function (•) defined as The row vectors k x l ∈ R 1×m , l ∈ {1, . . . , n} of the Kalman gain matrix K k ∈ R n×m contain the m gains of the innovation d k for the l th state. Through the element-wise multiplication the vector γ k,x l ∈ R 1×m is obtained, whose elements are the proportion of the measurements in the estimation of the l th state. For the investigation of a discrete time interval k ∈ [k 0 , k end ], the absolute values of the measurements can be summed up. For better comparability, normalization is performed to quantify the relative shares in the overall estimation: The row vectors ∈ ℝ 1× , ∈ {1, … , } of the Kalman gain matrix ∈ ℝ × contain the gains of the innovation for the ℎ state. Through the element-wise multiplication the vector , ∈ ℝ 1× is obtained, whose elements are the proportion of the measurements in the estimation of the ℎ state. For the investigation of a discrete time interval ∈ [ 0 , ], the absolute values of the measurements can be summed up. For better comparability, normalization is performed to quantify the relative shares in the overall estimation: Analysis of the measurands' influence makes sense only for time points in which the state is observable. Due to the indicator function (•) defined as only those time points are considered, in which the observability measure for the -th state ,obs, does not exceed a fixed upper bound obs, . This results in a scalar value for each measurand. The larger the value, the more dominant the sensor is in the overall estimation. Hereinafter, the method is, therefore, referred to as dominance analysis.
Analysis of the measurands' influence makes sense only for time points in which the state is observable. Due to the indicator function The row vectors ∈ ℝ 1× , ∈ {1, … , } of the Kalman gain matrix ∈ ℝ × contain the gains of the innovation for the ℎ state. Through the element-wise multiplication the vector , ∈ ℝ 1× is obtained, whose elements are the proportion of the measurements in the estimation of the ℎ state. For the investigation of a discrete time interval ∈ [ 0 , ], the absolute values of the measurements can be summed up. For better comparability, normalization is performed to quantify the relative shares in the overall estimation: Analysis of the measurands' influence makes sense only for time points in which the state is observable. Due to the indicator function (•) defined as only those time points are considered, in which the observability measure for the -th state ,obs, does not exceed a fixed upper bound obs, . This results in a scalar value for each measurand. The larger the value, the more dominant the sensor is in the overall estimation. Hereinafter, the method is, therefore, referred to as dominance analysis.
defined as

∶=
The row vectors ∈ ℝ 1× , ∈ {1, … , } of the Kalman gain matrix ∈ ℝ × contain the gains of the innovation for the ℎ state. Through the element-wise multiplication the vector , ∈ ℝ 1× is obtained, whose elements are the proportion of the measurements in the estimation of the ℎ state. For the investigation of a discrete time interval ∈ [ 0 , ], the absolute values of the measurements can be summed up. For better comparability, normalization is performed to quantify the relative shares in the overall estimation: Analysis of the measurands' influence makes sense only for time points in which the state is observable. Due to the indicator function (•) defined as only those time points are considered, in which the observability measure for the -th state ,obs, does not exceed a fixed upper bound obs, . This results in a scalar value for each measurand. The larger the value, the more dominant the sensor is in the overall estimation. Hereinafter, the method is, therefore, referred to as dominance analysis. (17) only those time points k are considered, in which the observability measure for the l-th state σ k,obs,x l does not exceed a fixed upper bound σ obs,x l . This results in a scalar value for each measurand. The larger the value, the more dominant the sensor is in the overall estimation. Hereinafter, the method is, therefore, referred to as dominance analysis.

Novel Holistic Method for Kalman Filter Design
Based on the concepts described in the previous sections, a holistic method is presented hereinafter to optimally design a Kalman filter with respect to structure and parameterization.
Optimal filter structure: To tackle the fundamental issue of appropriate measurands' selection, their respective influence on the estimation problem has to be quantified. On the one hand, the state-individual quantitative observability measure presented in Section 3.1.3 is used to analyze the properties of a measurand with respect to observability. On the other hand, the dominance analysis according to Section 3.2 allows the individual contribution of each measurand to the estimation to be evaluated. This enables the identification of sensor variables with a high or low information content and to accordingly adjust their selection.
Optimal filter parameterization: The DLR's Multi-Objective Parameter Synthesis (MOPS) optimization tool [24] is used to determine the optimal filter parameters. MOPS provides a variety of optimization methods. For the present design method, a pattern search approach is exploited, being able to deal with constrained nonlinear problems.
The holistic Kalman filter design method is an iterative process combining the concepts presented above as shown in Figure 4.
Starting from an arbitrary sensor configuration j, j ∈ N, the system and measurement noise covariance matrices Q j and R j , as well as the UKF parameters (α, β, κ) j,UKF , are determined with the help of MOPS. Considering the requirement for an optimal tracking behavior with the estimation error x to be minimized, the optimization problem, i.e., the cost function, can be formulated as

Novel Holistic Method for Kalman Filter Design
Based on the concepts described in the previous sections, a holistic method is presented hereinafter to optimally design a Kalman filter with respect to structure and parameterization.
Optimal filter structure: To tackle the fundamental issue of appropriate measurands' selection, their respective influence on the estimation problem has to be quantified. On the one hand, the state-individual quantitative observability measure presented in Section 3.1.3 is used to analyze the properties of a measurand with respect to observability. On the other hand, the dominance analysis according to Section 3.2 allows the individual contribution of each measurand to the estimation to be evaluated. This enables the identification of sensor variables with a high or low information content and to accordingly adjust their selection.
Optimal filter parameterization: The DLR's Multi-Objective Parameter Synthesis (MOPS) optimization tool [24] is used to determine the optimal filter parameters. MOPS provides a variety of optimization methods. For the present design method, a pattern search approach is exploited, being able to deal with constrained nonlinear problems.
The holistic Kalman filter design method is an iterative process combining the concepts presented above as shown in Figure 4. Starting from an arbitrary sensor configuration , ∈ ℕ, the system and measurement noise covariance matrices and , as well as the UKF parameters ( , , ) ,UKF , are determined with the help of MOPS. Considering the requirement for an optimal tracking behavior with the estimation error to be minimized, the optimization problem, i.e., the cost function, can be formulated as s.t. inequality constraints

Optimal Parameterization
Optimal Structure (18) s.t. inequality constraints where c l , l {1, . . . , n}, is a specific weight factor for each state estimation error. The inequality constraint considers the filter self-diagnostics, according to which the filter only works correctly if the estimation errors lie within the confidence interval, i.e., if the standard deviation of the estimation errors is ± diag P + k . According to the Gaussian distribution, ≈68% of all estimation errors lie within the confidence interval corresponding to ± diag P + k . Thus, about 30% of the values are not considered, although the filter would work correctly. For that reason, the constraints appear to be rather restrictive.
3.1.3 is used to analyze the properties of a measurand with respect to observability. On the other hand, the dominance analysis according to Section 3.2 allows the individual contribution of each measurand to the estimation to be evaluated. This enables the identification of sensor variables with a high or low information content and to accordingly adjust their selection.
Optimal filter parameterization: The DLR's Multi-Objective Parameter Synthesis (MOPS) optimization tool [24] is used to determine the optimal filter parameters. MOPS provides a variety of optimization methods. For the present design method, a pattern search approach is exploited, being able to deal with constrained nonlinear problems.
The holistic Kalman filter design method is an iterative process combining the concepts presented above as shown in Figure 4. Starting from an arbitrary sensor configuration , ∈ ℕ, the system and measurement noise covariance matrices and , as well as the UKF parameters ( , , ) ,UKF , are determined with the help of MOPS. Considering the requirement for an optimal tracking behavior with the estimation error to be minimized, the optimization problem, i.e., the cost function, can be formulated as s.t. inequality constraints  The filter optimally parameterized with MOPS can then be evaluated for a maneuver by means of the observability analysis with the criterion σ obs,µ max presented in Section 3.1.3, as well as the dominance analysis presented in Section 3.2. This loop can be run through for a filter configuration j. By comparing different sensor settings, a quantitative statement of how much a sensor contributes to the observability (observability analysis) and to the estimation (dominance analysis) can be made. Thus, the effort to provide a measurand, in reality, can be compared to its benefit for the estimation.

Vehicle State Estimator
The presented universally applicable design and analysis method is exemplarily applied to a vehicle state estimation problem using an unscented Kalman filter (UKF). The vehicle considered in this paper is the ROboMObil (ROMO)-an innovative robotic electric vehicle concept developed at the DLR Robotics and Mechatronics Center. The design of the ROMO is based on the so-called "wheel robot" concept [25] with all wheel by-wire steering capabilities, where the drivetrain, brakes, steering system, spring, and dampers are integrated into each of the four wheels. The dissemination of electric vehicles in the automotive market results in a variety of estimation problems related, e.g., to the battery state [3] as well as to the vehicle dynamics, see, e.g., [26].
The most important variables for describing driving stability are probably the vehicle side-slip angle as a measure of the current vehicle stability, and the maximum coefficient of friction between the tire and the road (TRFC) as a measure of the stability limit. Although the vehicle side-slip angle can be directly measured by a high-precision IMU (inertial measurement unit) or an optical road sensor [12], it is only estimated in production vehicles due to the high costs of such sensors. Methods for estimating the TRFC have been intensely researched in science and technology for decades [27,28]. The reconstruction of the TRFC from the measurements is quite complex since a certain excitation threshold has to be exceeded in order to distinguish between the different coefficients of friction. Knowledge about the TRFC is of paramount importance for manned driving since the road conditions are often underestimated even by experienced drivers. At the same time, the TRFC also plays a key role for autonomous driving functions because computers, unlike human drivers, have no intuition that, for example, the speed should be reduced on a wet road covered with leaves. In addition to the vehicle side-slip angle and the TRFC, the vehicle velocity over ground and the vehicle yaw rate are also estimated. However, the focus will lie on the TRFC due to its growing importance and the complexity of its estimation.

Vehicle Models
In this section, two vehicle models are presented. The first one is a detailed and highly accurate full vehicle model, which serves as the reference. Since the model is simulated only once per maneuver to generate the reference data, the high effort required for this computation is acceptable. The second model is the filter prediction model, which is permanently simulated during the design and analysis process. Therefore, this model has to be computationally more efficient than the first one.

High Fidelity Reference Model
The full vehicle model is a detailed, high-precision validated multiphysics Modelica model, which, as a reference, represents real vehicle behavior in a very accurate manner [29]. In addition to the vehicle model, there is an extensive environmental model in which parameters such as the road gradient, the ambient temperature, or the coefficient of friction between tire and road can be defined. Moreover, there is a driver model that allows individual trajectories to be driven in addition to predefined standard maneuvers. The chassis is a kinematic multi-body model. Pacejka's Magic Formula 5.2 [30] is used as the tire model. The spring and damper characteristics are described by nonlinear characteristic curves. The overall model has over 100 states and a nonlinear system of equations with a dimension >500.

Filter Prediction Model
A nonlinear two-track model is used for the prediction model in the filter, see Figure 5. The quantities marked by (·) C are expressed in the car coordinate system with an origin in the center of gravity (CoG), while those indicated by (·) W i are in the i th wheel robot coordinate system. The vehicle state is described by four states: side-slip angle β C , vehicle velocity over ground v C , maximum tire-road friction coefficient µ max , and yaw rate . ψ C .
The equations of motion are thus obtained as: . . ..
For the tire model representing the main nonlinearity in the vehicle model, a slightly simplified version of Pacejka's Magic Formula 5.2 [30] is used. The detailed tire model can be found in Appendix B. The maximum tire-road friction coefficient TRFC, which is tire-specific in reality, is described as a whole-vehicle equal parameter in this paper. The TRFC is a parameter in the tire model, which is often described in the literature as . µ max = 0 (representing a random walk process under the presence of white noise). Another possibility is to model . µ max as a first-order filter (PT1). When the excitation is insufficient, the friction value would converge to a predetermined value [31], e.g., µ max = 1 (see Equation (22) with T µ = 1s), which corresponds to a high friction value. This state description is also called artificial stabilization [32]. An additional benefit of the PT1 formulation is its anti-windup effect in a constrained estimation. As presented in the next section, the TRFC is constrained in the measurement update to a maximum value of 1. In a random walk state description, for example, this would cause the integrators to accumulate the TRFC even during saturation (wind-up). The system would need some iteration steps to "unload" again. A detailed presentation of the model equations complementing Equations (20)-(23) can be found in Appendix A.

Constraint of the Tire Road Friction Coefficient
For the estimation problem presented here, it is useful to integrate prior physical knowledge about the state TRFC µ max into the estimation via the constraints. The currently utilized frictional force between the tire and the road, i.e., the instantaneous friction value µ act , is calculated via the relationship of the Kamm circle (see Figure 6): Obviously, the measured longitudinal and lateral vehicle accelerations can be exploited to approximate the currently utilized friction value. With the assumption that the vehicle is in a stable state of motion, the friction value at time step k for the dynamic lower bound is determined to be: µ max,low,k ≥ µ act,k .
Sensors 2021, 21, x FOR PEER REVIEW 13 of 28 (representing a random walk process under the presence of white noise). Another possibility is to model ̇m ax as a first-order filter (PT1). When the excitation is insufficient, the friction value would converge to a predetermined value [31], e.g., max = 1 (see Equation (22) with = 1 ), which corresponds to a high friction value. This state description is also called artificial stabilization [32]. An additional benefit of the PT1 formulation is its anti-windup effect in a constrained estimation. As presented in the next section, the TRFC is constrained in the measurement update to a maximum value of 1. In a random walk state description, for example, this would cause the integrators to accumulate the TRFC even during saturation (wind-up). The system would need some iteration steps to "unload" again. A detailed presentation of the model equations complementing Equations (20)-(23) can be found in Appendix A.

Constraint of the Tire Road Friction Coefficient
For the estimation problem presented here, it is useful to integrate prior physical knowledge about the state TRFC max into the estimation via the constraints. The currently utilized frictional force between the tire and the road, i.e., the instantaneous friction value act , is calculated via the relationship of the Kamm circle (see Figure 6): Obviously, the measured longitudinal and lateral vehicle accelerations can be exploited to approximate the currently utilized friction value. With the assumption that the vehicle is in a stable state of motion, the friction value at time step for the dynamic lower bound is determined to be: max,low, ≥ act, .
(25) Figure 6. TRFC constraint based on the Kamm circle and its limits.
For the maximum possible TRFC, a value of 1 is assumed (depending on the tire type, however, values for max up to 1.2 are also possible on dry surfaces [28]). For the upper bound, this provides:  For the maximum possible TRFC, a value of 1 is assumed (depending on the tire type, however, values for µ max up to 1.2 are also possible on dry surfaces [28]). For the upper bound, this provides: µ max,up,k ≤ 1.
The presented constraints primarily concern the Kalman filter measurement update step. The values of the prediction step are also indirectly limited, see Equation (22). Without additional limitation of the predictor, this fact would lead to a wind-up effect of the integral terms when the system is saturated. This has to be avoided in the interest of an optimal estimation (in this paper, implemented by a PT1 modeling approach, Section 4.1.2).

Sensors and Measurands
Any information about the vehicle state is provided by the measurands. Consequently, the selection of the measured variables is a central issue in the filter design. The ROMO is equipped with sensors to determine the following variables: with the index i ∈ { f l, f r, rl, rr} for the tire position. The ROMO is equipped with several other sensors that are not relevant for this research work, see [12]. All signals are provided by the ROMO's vehicle dynamics control (VDC) [33]. Most of the listed variables are also used by the ESC (electronic stability control), and are thus, in principle, also available in a production vehicle. The ROMO's sensors are analyzed in [12] with respect to the noise properties, bias, delay, etc. As a result, the sensor measurements can be represented as signals with realistic properties to synthetically generate the measurement data by simulating the reference vehicle model.
Virtual Measurands: In addition to the directly available measurements, the so-called virtual measurands can be used. They are calculated from an advantageous conversion of directly measurable quantities. One of the benefits of using virtual measurands is the wheelrelated consideration, providing a more precise physical description of the TRFC. Another advantage is the inclusion of additional measurands in the estimation problem (driving and braking torques as well as wheel speeds). The concept of using virtual measurands for vehicle state estimation is based on [32]. For a detailed derivation, reference is, therefore, made to this. To distinguish the virtual measurands from the real ones, they are marked by a tilde.
For the calculation of the virtual measurands, the vehicle is considered as a single-track model (STM). For that reason, the axle-wise mean value of the two steering angles (left and right) can be used as a simplification δ W p STM = δ W i + δ W i /2, with p ∈ {front, rear}, as well as the mean value of the wheel speeds on the left-and right-hand sides of the front

Virtual vehicle velocity:
The virtual speed v C is obtained from the wheel speed ω W f STM and the effective tire radius r eff at the front axle: Virtual longitudinal axle force: A virtual longitudinal axle force can be determined from the driving and braking torques at the wheels and the tire radius: Virtual lateral axle force: The lateral acceleration at the center of gravity, the yaw rate, the virtual longitudinal axle force, and the effective steering angle are used to calculate a virtual lateral axle force: Virtual side-slip angle: The calculation of the virtual side-slip angle β C is based on the assumption that the slip angle at the front axle α f can be approximated by the associated cornering stiffness c α, f as α f ≈ F W f y / c α, f . The identification of the cornering stiffness c α, f is performed via a parameter optimization for representative driving dynamics maneuvers with DLR MOPS. With the help of the approximated tire-slip angle α f , the virtual sideslip angle is obtained using virtual longitudinal and lateral vehicle velocities v C x and v C y , respectively, [32]: (30)

Estimator Setups
The first step of estimator design is the selection of a specific filter type. In intensive investigations, a UKF with a sampling time of 20 ms has proven to be a good compromise between computational effort and accuracy. A linear Kalman filter is not an option because of strong nonlinearities. An MHE is also not taken into account due to the enormous computation time required. An extended Kalman filter (EKF) basically needs less computational effort than a UKF, but it delivers good results only for small sampling times. At the same time, a UKF provides satisfactory estimates, i.e., with a lower computational cost, even for higher sampling times.
The filter parameters, namely, the system noise and measurement noise covariance matrices Q and R, respectively, as well as the factors α UKF , β UKF , γ UKF for the approximation of the probability densities, are determined with the help of the design method presented in Section 3.3. The entries of the state vector x ∈ R 4×1 , which have to be estimated, are the side-slip angle, the vehicle velocity over ground, the maximum tire-road friction coefficient, and the yaw rate (see Equation (20)- (23)): The inputs u ∈ R 11×1 consist of four-wheel steering angles, four-wheel speeds, the longitudinal and lateral acceleration, and the yaw rate: The state µ max is constrained by the method presented in Section 2.2. According to the available sensor signals presented in the previous section, 12 variables can be measured. The measured variables y are divided into four setups for further analysis and compared with each other. The influence of different measurands on the observability of the states is assessed by the nonlinear quantitative observability measure according to Equation (12). The contribution of different measurands to the estimation is quantified in the dominance analysis according to Section 3.2 Table 1 shows the four measurement setups. The green color indicates that the respective sensor signal is used as a measurement. Setup 1 represents the maximum configuration with all 12 available measurands. In comparison, the real side-slip angle and the real vehicle velocity are no longer available for setup 2. In setup 3, the tire self-aligning torque (SAT) is additionally removed. In Setup 4, the virtual longitudinal and lateral axle forces are not considered anymore (no driving and braking torque sensors at wheels available). With the remaining five measurands, this setup represents the minimum configuration.

Test Track and Maneuver
The selection of suitable excitations is of paramount importance for the analysis of an estimation filter. The excitations should contain the broadest possible spectrum of components that can arise during filter operation. In the present case, which deals with vehicle state estimation, a suitable test track or driving maneuver has to be chosen. A maneuver on a lying and crossing eight, a so-called figure-eight maneuver, is used. It covers many characteristic properties of the vehicle dynamics, and therefore, enables an optimal filter analysis: The track is divided into segments with different TRFCs in the range of 0.5 ≤ µ max ≤ 1 (see Figure 7), which roughly correspond to driving on a road with a thin layer of ice and a dry road with a high coefficient of friction, respectively. Further specifications of the track and the test maneuver can be found in Table 2.

Results
This section presents the results of the four filter setups mentioned above, which are parameterized and analyzed with the help of the design method according to Section 3.3. In the beginning, the results of the estimated states by a UKF are shown, i.e., the tracking performance of the filters (Section 5.1). On the one hand, the fit value is used to assess the estimation quality, i.e., to evaluate the error between estimates and true trajectories [34]. This gives a percentage fit value (100% perfect fit) and is thus a well-interpretable metric. On the other hand, the root-mean-square error (RMSE) is used as a physically interpretable error measure, see Appendix C.
Using the state-specific quantitative observability measure, observability analysis is performed in Section 5.2 for every state in each of the four setups. To compare the respective setups with each other, the following evaluation measure is defined: Here, the curves of the observability measure of the i th state σ obs,x i are summed up over all time steps k for every setup j and divided by the corresponding value of reference configuration 4. This provides a percentage value showing the observability improvement of setup j compared to setup 4 (minimum configuration). Next, in Section 5.3, dominance analysis is performed for each setup using the state µ max as an example to investigate the influence of the respective measurands on the estimation. The results of the dominance analysis are presented in Section 5.4 and used to rank the measurands according to their importance for the estimation. Finally, the most significant insights of the design analysis are summarized in Section 5.5.

Vehicle State Estimation
The estimation results are presented in Table 3 using the fit value and the RMSE for the four different filter setups for each of the estimated states. Very good or very bad estimation results are represented by the colors green or red, respectively. Figure 8 shows the estimated trajectories corresponding to the performance values listed in Table 3 as well as the reference values.
Setup 1 (maximum configuration) can estimate all four states with a very good accuracy.
Since the side-slip angle and the vehicle velocity are measured directly, the fit values are over 90%. The TRFC can be reconstructed from the measurements very well resulting in a fit value of almost 80% or, equivalently, an RMSE of 0.04. The loss of the direct side-slip angle and vehicle velocity measurements in the second setup has only a minor effect on the estimation performance. The coefficient of friction loses almost six percentage points in the fit value compared to setup 1. However, it can be well-reconstructed thanks to the sensors, which are still available. In the third setup, in which the measurement of the tire selfaligning torque (SAT) is removed, larger changes in the estimation accuracy are observed. Nevertheless, the TRFC can be still estimated satisfactorily with a fit value of almost 60% or, equivalently, an RMSE of 0.08. For the fourth setup (minimum configuration), where the driving and braking torque information at the wheels is no longer available, the TRFC cannot be estimated satisfactorily, i.e., with a fit value of less than 15% or an RMSE of 0.15. It is noteworthy that the side-slip angle, the estimation of which is based only on the virtual measurand β C , can still be well-reconstructed with a fit value of more than 80%. The vehicle velocity can be estimated by all setups identically well (fit value >97%). The yaw rate is available as a direct measurand and does not have to be reconstructed. Therefore, this state is not further examined and its time course is not shown in Figure 8. After having discussed the estimation quality in general, some specific estimation segments are analyzed in more detail below.
The individual state variables are obviously strongly coupled with each other, so that the estimation deviations over time are seen not only in a single state. In 1 , the side-slip angle is heavily underestimated for setup 4, caused by an overestimation of the coefficient of friction in the same time range 4 . The opposite case can be seen in 2 , where the side-slip angle for setups 3 and 4 is overestimated because the corresponding coefficients of friction in 6 are considerably underestimated. Since the side-slip angle and the vehicle velocity are measured directly, the fit values are over 90%. The TRFC can be reconstructed from the measurements very well resulting The observed strong coupling between the deviations in the estimated side-slip anglê β C and the coefficient of frictionμ max is also indicated, without going into further detail, by high correlation coefficients between the estimation errors of these variables. For that reason, it can be stated that a good estimate of the coefficient of friction can only be obtained by a good estimate of the side-slip angle.
In 3 and 5 , the coefficient of friction shows for setup 4 a step-like behavior at the beginning of isolated longitudinal excitation phases (see velocity plot), converging then to a value of 1. An isolated braking maneuver 7 demonstrates similar behavior. For setup 4, the excitation is obviously insufficient in these regions to reconstruct the coefficient of friction from the measurements. The state is no longer observable in this case, which is also confirmed by the quantitative observability measure in the following section, and tends towards the preset friction value of 1 due to the modeled first-order dynamics of µ max (artificial stabilization, see Section 4.1.2).

Observability Analysis
The comparison of the observability properties of different setups (see Equation (33)) can be found in Table 4. The percentage improvement of the observability related to setup 4, with its minimum sensor equipment, is shown here for all states. No analysis is performed for the estimated yaw rate since this state is directly measured for each setup and is thus fully observable per se. In general, it can be seen from Table 4 that the more measurands there are available (setup 4 = minimum number of measurands, setup 1 = maximum number), the better the observability. The trajectories of the quantitative observability measure σ obs,x i over time are shown in Figure 9.

Side-slip angleβ C
Although the side-slip angle is directly measured in setup 1, its observability is hardly any better compared to setup 2 measuring this angle only virtually. Setup 2 leads to a large increase in observability compared to setup 3. This means that the measurand tire self-aligning torque (SAT) M W f z of setup 2 provides highly relevant information for the side-slip angle. The fact that the observability of setup 3 is by over 40% better than that of setup 4 implies that the driving and braking torques at the wheels also represent a good information source for the side-slip angle estimation.
Generally, it can be noted that setups 1 and 2 have similar and the best observability accuracies of the vehicle side-slip angle. Setup 3 provides significantly worse observability, while setup 4 delivers the worst one.

Vehicle velocityv C
The observability measure σ obs,v (see center of Figure 9) of setup 1 is with a few exceptions the smallest, which can be anticipated, since the vehicle velocity is directly measured. Peaks with poor observability are present in each setup. They appear when there is a sign change in the longitudinal vehicle acceleration (i.e., switching from acceleration to deceleration and vice versa). In these moments, there is no excitation for a short period because of the acceleration reversal. This phenomenon makes the observability worsen and is captured by the quantitative criterion. The fact that, in these moments, setup 2 is better than setup 3 implies that the additional measurand SAT contains information also about the vehicle velocity and additionally supports the estimation ofv C in such moments. Overall, the observability is very good for all setups, even for the minimally equipped setup 4, meaning that not many sensors are required when only the vehicle velocity is estimated.
Sensors 2021, 21, x FOR PEER REVIEW 20 of 28 Table 4 that the more measurands there are available (setup 4 = minimum number of measurands, setup 1 = maximum number), the better the observability. The trajectories of the quantitative observability measure , over time are shown in Figure 9.

Side-slip angle ̂
Although the side-slip angle is directly measured in setup 1, its observability is hardly any better compared to setup 2 measuring this angle only virtually. Setup 2 leads to a large increase in observability compared to setup 3. This means that the measurand tire self-aligning torque (SAT) of setup 2 provides highly relevant information for the side-slip angle. The fact that the observability of setup 3 is by over 40% better than that of

Tire road friction coefficientμ max
The observabilities of setups 1 and 2 are nearly identical. In both cases, the observability is improved by over 96% compared to the reference setup. This fact implies that the measurands' side-slip angle and vehicle velocity, which can only be found in setup 1, do not yield in this case any benefits in terms of observability.
The observability of setup 3 is worse than that of setup 1 or 2 because of missing measurands required for maneuver segments with lateral excitation (namely, β C , v C , and M W f z ). In contrast, segments with longitudinal excitation-such as, e.g., 9 -demonstrate a similarly good observability accuracy as setups 1 and 2. This indicates that the sensor signals driving and braking torques M In about the vehicle velocity and additionally supports the estimation of ̂ in such moments. Overall, the observability is very good for all setups, even for the minimally equipped setup 4, meaning that not many sensors are required when only the vehicle velocity is estimated.
Tire road friction coefficient ̂ The observabilities of setups 1 and 2 are nearly identical. In both cases, the observability is improved by over 96% compared to the reference setup. This fact implies that the measurands' side-slip angle and vehicle velocity, which can only be found in setup 1, do not yield in this case any benefits in terms of observability.
The observability of setup 3 is worse than that of setup 1 or 2 because of missing measurands required for maneuver segments with lateral excitation (namely, , , and ). In contrast, segments with longitudinal excitation-such as, e.g., ⑨-demonstrate a similarly good observability accuracy as setups 1 and 2. This indicates that the sensor signals driving and braking torques and , respectively, contain enough information about the TRFC in this case.
In ⑩, setups 3 and 4 deliver worse observability than setups 1 and 2. This phenomenon is also reflected in the estimation accuracy shown in ⑥ (see Figure 8), where significant underestimation is present in both setups. In maneuver segments with isolated longitudinal excitation, namely ⑧, ⑨ (accelerating), and ⑪ (braking), the total loss of observability for setup 4 is captured by the observability measure when obs, max ≥ 1. In the corresponding time plots of the estimated TRFC in Figure 8, namely ③, ⑤, and ⑦, max converges to its default stabilized value of 1 due to lack of observability.

Dominance Analysis
The measurands were analyzed with respect to their observability property in the previous subsection, and their contributions to the estimation are presented now in what follows (dominance properties). Due to the focus of the paper, the dominance analysis is limited to the state TRFC max . However, the methodology would also allow analysis to be performed of other states as well. Table 5 shows the contributions of the measurands to the estimation of the state max for different setups.
Although a total of 12 measurands are available in setup 1, only three measurandsnamely, ̃a nd ̃( i.e., driving and braking torque and , respectively) and the SAT -are used by the filter to reconstruct max and account for almost 90% of the dominance. The measurements of the side-slip angle and the vehicle velocity appear, at first glance, to be very valuable, with a negligible contribution of less than 5% in total. The dominance of the three measurands mentioned above does not change for setup 2. The removal of and seems to be compensated by the virtual axle lateral force ̃.
, setups 3 and 4 deliver worse observability than setups 1 and 2. This phenomenon is also reflected in the estimation accuracy shown in 6 (see Figure 8), where significant underestimation is present in both setups. In maneuver segments with isolated longitudinal excitation, namely about the vehicle velocity and additionally supports the estimation of ̂ in such moments. Overall, the observability is very good for all setups, even for the minimally equipped setup 4, meaning that not many sensors are required when only the vehicle velocity is estimated.

Tire road friction coefficient ̂
The observabilities of setups 1 and 2 are nearly identical. In both cases, the observability is improved by over 96% compared to the reference setup. This fact implies that the measurands' side-slip angle and vehicle velocity, which can only be found in setup 1, do not yield in this case any benefits in terms of observability.
The observability of setup 3 is worse than that of setup 1 or 2 because of missing measurands required for maneuver segments with lateral excitation (namely, , , and ). In contrast, segments with longitudinal excitation-such as, e.g., ⑨-demonstrate a similarly good observability accuracy as setups 1 and 2. This indicates that the sensor signals driving and braking torques and , respectively, contain enough information about the TRFC in this case.
In ⑩, setups 3 and 4 deliver worse observability than setups 1 and 2. This phenomenon is also reflected in the estimation accuracy shown in ⑥ (see Figure 8), where significant underestimation is present in both setups. In maneuver segments with isolated longitudinal excitation, namely ⑧, ⑨ (accelerating), and ⑪ (braking), the total loss of observability for setup 4 is captured by the observability measure when obs, max ≥ 1. In the corresponding time plots of the estimated TRFC in Figure 8, namely ③, ⑤, and ⑦, max converges to its default stabilized value of 1 due to lack of observability.

Dominance Analysis
The measurands were analyzed with respect to their observability property in the previous subsection, and their contributions to the estimation are presented now in what follows (dominance properties). Due to the focus of the paper, the dominance analysis is limited to the state TRFC max . However, the methodology would also allow analysis to be performed of other states as well. Table 5 shows the contributions of the measurands to the estimation of the state max for different setups.
Although a total of 12 measurands are available in setup 1, only three measurandsnamely, ̃a nd ̃( i.e., driving and braking torque and , respectively) and the SAT -are used by the filter to reconstruct max and account for almost 90% of the dominance. The measurements of the side-slip angle and the vehicle velocity appear, at first glance, to be very valuable, with a negligible contribution of less than 5% in total. The dominance of the three measurands mentioned above does not change for setup 2. The removal of and seems to be compensated by the virtual axle lateral force ̃.
, about the vehicle velocity and additionally supports the estimation of ̂ in such moments. Overall, the observability is very good for all setups, even for the minimally equipped setup 4, meaning that not many sensors are required when only the vehicle velocity is estimated.

Tire road friction coefficient ̂
The observabilities of setups 1 and 2 are nearly identical. In both cases, the observability is improved by over 96% compared to the reference setup. This fact implies that the measurands' side-slip angle and vehicle velocity, which can only be found in setup 1, do not yield in this case any benefits in terms of observability.
The observability of setup 3 is worse than that of setup 1 or 2 because of missing measurands required for maneuver segments with lateral excitation (namely, , , and ). In contrast, segments with longitudinal excitation-such as, e.g., ⑨-demonstrate a similarly good observability accuracy as setups 1 and 2. This indicates that the sensor signals driving and braking torques and , respectively, contain enough information about the TRFC in this case.
In ⑩, setups 3 and 4 deliver worse observability than setups 1 and 2. This phenomenon is also reflected in the estimation accuracy shown in ⑥ (see Figure 8), where significant underestimation is present in both setups. In maneuver segments with isolated longitudinal excitation, namely ⑧, ⑨ (accelerating), and ⑪ (braking), the total loss of observability for setup 4 is captured by the observability measure when obs, max ≥ 1. In the corresponding time plots of the estimated TRFC in Figure 8, namely ③, ⑤, and ⑦, max converges to its default stabilized value of 1 due to lack of observability.

Dominance Analysis
The measurands were analyzed with respect to their observability property in the previous subsection, and their contributions to the estimation are presented now in what follows (dominance properties). Due to the focus of the paper, the dominance analysis is limited to the state TRFC max . However, the methodology would also allow analysis to be performed of other states as well. Table 5 shows the contributions of the measurands to the estimation of the state max for different setups.
Although a total of 12 measurands are available in setup 1, only three measurandsnamely, ̃a nd ̃( i.e., driving and braking torque and , respectively) and the SAT -are used by the filter to reconstruct max and account for almost 90% of the dominance. The measurements of the side-slip angle and the vehicle velocity appear, at first glance, to be very valuable, with a negligible contribution of less than 5% in total. The dominance of the three measurands mentioned above does not change for setup 2. The removal of and seems to be compensated by the virtual axle lateral force ̃.
(accelerating), and about the vehicle velocity and additionally supports the estimation of ̂ in such moments. Overall, the observability is very good for all setups, even for the minimally equipped setup 4, meaning that not many sensors are required when only the vehicle velocity is estimated.

Tire road friction coefficient ̂
The observabilities of setups 1 and 2 are nearly identical. In both cases, the observability is improved by over 96% compared to the reference setup. This fact implies that the measurands' side-slip angle and vehicle velocity, which can only be found in setup 1, do not yield in this case any benefits in terms of observability.
The observability of setup 3 is worse than that of setup 1 or 2 because of missing measurands required for maneuver segments with lateral excitation (namely, , , and ). In contrast, segments with longitudinal excitation-such as, e.g., ⑨-demonstrate a similarly good observability accuracy as setups 1 and 2. This indicates that the sensor signals driving and braking torques and , respectively, contain enough information about the TRFC in this case.
In ⑩, setups 3 and 4 deliver worse observability than setups 1 and 2. This phenomenon is also reflected in the estimation accuracy shown in ⑥ (see Figure 8), where significant underestimation is present in both setups. In maneuver segments with isolated longitudinal excitation, namely ⑧, ⑨ (accelerating), and ⑪ (braking), the total loss of observability for setup 4 is captured by the observability measure when obs, max ≥ 1. In the corresponding time plots of the estimated TRFC in Figure 8, namely ③, ⑤, and ⑦, max converges to its default stabilized value of 1 due to lack of observability.

Dominance Analysis
The measurands were analyzed with respect to their observability property in the previous subsection, and their contributions to the estimation are presented now in what follows (dominance properties). Due to the focus of the paper, the dominance analysis is limited to the state TRFC max . However, the methodology would also allow analysis to be performed of other states as well. Table 5 shows the contributions of the measurands to the estimation of the state max for different setups.
Although a total of 12 measurands are available in setup 1, only three measurandsnamely, ̃a nd ̃( i.e., driving and braking torque and , respectively) and the SAT -are used by the filter to reconstruct max and account for almost 90% of the dominance. The measurements of the side-slip angle and the vehicle velocity appear, at first glance, to be very valuable, with a negligible contribution of less than 5% in total. The dominance of the three measurands mentioned above does not change for setup 2. The removal of and seems to be compensated by the virtual axle lateral force ̃.
(braking), the total loss of observability for setup 4 is captured by the observability measure when σ obs,µ max ≥ 1. In the corresponding time plots of the estimated TRFC in Figure 8, namely 3 , 5 , and 7 , µ max converges to its default stabilized value of 1 due to lack of observability.

Dominance Analysis
The measurands were analyzed with respect to their observability property in the previous subsection, and their contributions to the estimation are presented now in what follows (dominance properties). Due to the focus of the paper, the dominance analysis is limited to the state TRFC µ max . However, the methodology would also allow analysis to be performed of other states as well. Table 5 shows the contributions of the measurands to the estimation of the state µ max for different setups. z are used by the filter to reconstruct µ max and account for almost 90% of the dominance. The measurements of the side-slip angle β C and the vehicle velocity v C appear, at first glance, to be very valuable, with a negligible contribution of less than 5% in total. The dominance of the three measurands mentioned above does not change for setup 2. The removal of β C and v C seems to be compensated by the virtual axle lateral force F y W f .
In setup 3, where one of the dominant measurands, the SAT M W f z (30% of the total dominance), is no longer available, the contribution of the remaining sensors changes. One part of the estimation contribution of M W f z seems to be compensated by the virtual axle lateral forces, and another part by the measurement of the vehicle lateral acceleration a C y . At the same time, part of the contribution of the virtual axle forces F x W j is taken over by the longitudinal acceleration a C x . However, the measurands F x W j still remain dominant. The virtual side-slip angle makes a small contribution to the TRFC estimate for the first time in setup 3.
In setup 4, all previously dominant measurands are no longer available. The virtual side-slip angle β C and the vehicle lateral acceleration a C y take over the dominance of the sensors with a total contribution of more than 80%. The virtual vehicle velocity v C as well as the longitudinal acceleration a C x account for less than 10% of the overall dominance. As shown in the previous section, the problem is not observable for large parts of longitudinal maneuvers in setup 4. For that reason, it is completely plausible that no information can be extracted from the corresponding measurands.
For all setups, it seems that the yaw rate . ψ C is involved in the estimation of the TRFC with a relatively small contribution. However, it should be noted that . ψ C is used to calculate the virtual lateral axle forces F y W j as well as the virtual side-slip angle β C resulting in its additional indirect contribution through these variables.

Measurands Ranking
Based on the dominance analysis of the four different setups performed in the previous section (see Table 5), the influence and thus the "value" of each measurand for estimating the TRFC can be assessed. In Table 6, a ranking of the measurands is given. 1 2 The measurands with the biggest contribution to the TRFC estimation are the virtual longitudinal axle forces calculated from the wheel driving and braking torques. They are followed by the tire self-aligning torque (SAT). The sensor signals from the IMU, i.e., accelerations and velocities, are ranked next. The virtual side-slip angle and virtual vehicle velocity share last place in the ranking. Due to their purely approximate calculation ( v C with uncertainty in the tire radius, see Equation (27), β C with uncertainty in the cornering stiffness approximation, see Equation (30)), they are used more as additional sources. The information from directly available sensor signals, i.e., real ones, should be preferred. When looking at the table sequence, it is noticeable that information about the coefficient of friction is extracted from (listed in descending order) • forces or moments • accelerations (including virtual lateral axle forces since they are calculated from these measurands) • velocity or rotational speeds Since the tire forces and moments represent the most precise physical descriptor of the TRFC in terms of causality, they contain the most information about this state. For other descriptors, (inaccurate) conversion is necessary. It could be a purely analytical relationship (e.g., a = F/m) or an integration (e.g., v = a).

Summary of the Design Analysis
The side-slip angle could be estimated quite well with a total fit value of over 80% with all setups, even with those that only have virtual measurands and IMU information. The vehicle velocity could be reconstructed very well by all setups with a fit value of 97%. Therefore, it can be stated that no additional complex sensor technology is required to estimate this state. The coefficient of friction between the tire and the road is the state that is most difficult to reconstruct. The removal of sensors, starting from setup 1, shows large, direct effects on the estimation accuracy and observability properties. In this case, measurements of the side-slip angle and vehicle velocity do not provide any significant advantage. Therefore, the considerable effort required to provide these measurands can be saved in the case of production vehicles.
A setup that only uses the basic information from the IMU (setup 4) provides unsatisfactory results for TRFC estimation (approx. 15% fit value). For maneuvers with higher excitation, better estimation results are possible since the problem is then more observable. The estimation errors of the side-slip angle and TRFC are correlated, so that a good estimate of the TRFC is only possible if there is a good estimate of the side-slip angle.

Conclusions
In this paper, a novel and universal design and analysis method for nonlinear Kalman filters was presented. The method allows a systematic investigation of the measurands' influence on the estimation problem in terms of

•
Observability properties: Two novel quantitative nonlinear observability measures were presented.
Evaluation of the overall system via the numerical condition number of the observability matrix. A state-specific and physically interpretable observability measure via a weighted least-squares approach.
• Dominance properties: A new method quantifying the contribution and information content of a measurand for the state reconstruction.
To determine an optimal filter parameterization, the method uses an optimization algorithm. As an example, the method was applied to a vehicle state estimation problem focusing on the coefficient of friction between the tire and the road. For this purpose, an unscented Kalman filter with constraints was used. A nonlinear two-track model served as the prediction model, while a high-accuracy Modelica multi-body model was used as the reference model. The TRFC, whose estimation was the most challenging one among all states, could be reconstructed with a fit value lying between 15% and 80%. Using the design method, four filter setups with different available sensors (from a maximum to a minimum number of them) were analyzed and compared. Here, it was shown that for the structural design of an estimation filter, it is worth doing preliminary investigations about the influence of the measurands. For this estimation problem, it could be identified, for example, that the measurements of the side-slip angle and vehicle velocity over ground, which are costly to provide in a production vehicle, have a negligible influence on the investigated estimation properties for certain setups. Both sensors possess a total contribution to the estimation (dominance) of less than 5%. Thus, a sensor for measuring these signals is not necessary.
We plan to perform the analysis with additional sensors in the future, especially with new types of sensors such as camera, radar, and lidar. The vehicle state estimation problem should be implemented with the most promising setups in real driving tests, both prototypically and on an ECU (embedded system). Acknowledgments: The authors' thanks go to Andreas Pfeiffer for his valuable support.

Conflicts of Interest:
The authors declare no conflict of interest. Quantity expressed in the i-th wheel robot coordinate system (·) C Quantity expressed in the car coordinate system with origin in CoG (·) Virtual measurand (·) Estimated state

Appendix A. Two-Track Model
For the equations of the two-track model given below, reference is also made to Figure 5. The air resistance forces are given by: The longitudinal and lateral forces at the center of gravity are determined by: The accelerations at the vehicle center of gravity are calculated as follows: a C x = The quasi-stationary wheel loads are geometrically calculated from the lateral and longitudinal acceleration at the vehicle center of gravity: y . (A5)

Appendix B. Tire Model
For the tire model, a slightly simplified version of Pacejka's Magic Formula 5.2 [30] is used. As a simplification, among others, the wheel camber is neglected (γ W i = 0), all lambda scaling values are chosen to be λ = 1 (except for the coefficient of friction dependent factors λ µx = λ µy = µ max ), and the transient tire behavior is neglected. The tire position is described by the index i ∈ { f l, f r, rl, rr}.
The wheel load-dependent tire radius is: The longitudinal tire slip: The wheel load increment: The lateral tire slip: