Modelling and Wavelet-Based Identification of 3-DOF Vehicle Suspension System

In this paper, a three Degrees of Freedom (DOF) model of a quarter vehicle suspension system is proposed including the seat driver mass. The modal parameters of this system, which indicate the comfort and the safety of the suspension, are identified using Wavelet analysis. Two applications of wavelet analysis are presented: signal denoising based on the Discrete Wavelet Transform (DWT) and modal identification based on the Continuous Wavelet Transform (CWT). It is shown that the CWT analysis of the system response, initially denoised using DWT, allows the estimation of the natural pulsations and the damping ratios. The usefulness of the DWT in denoising and the accuracy of the CWT in modal identification are tested and confirmed by applying them to the proposed model. The complete modeling and identification of a 3-DOF vehicle suspension system is developed and the simulation results verify these statements and are satisfactory.


Introduction
The importance of the suspension system of a vehicle is that it is responsible of the comfort of the driver and the safety of the vehicle since it is the part of the vehicle which carries the body of the vehicle and which transmits all the forces between this body and the road.For that, the complete physical model of the suspension is often reduced under investigation of the vertical dynamics established on a quarter-vehicle model.The study of this model, and like any mechanical system, presents challenge for engineers and researchers.Noise reduction and parameters identification from system response are common problems in signal processing and vibrational applications [1][2][3].That's why numerous approaches have been developed recently and applied to denoise signals and to extract modal parameters of systems and structures in time and frequency domains [1].
The wavelet analysis is one of these new approaches that were revealed mainly because it represents an easy way of extracting time varying frequency components.Indeed, the Wavelet Transform technique is shown to be more effective than other denoising methods such the Fourier Transform (FT) or Compressed Sensing technique [1], and more useful than other identification methods such as Hilbert-Huang Transform [4][5][6].
In noise reduction, wavelets has been successfully used by Giaouris et al. [2] in order to propose a new pseudoadaptive denoising method based on the Wavelet Transform, this method adjust the level of signal decomposition and was applied effectively to an electric drive.
In systems identification, resonant pulsations and damping ratios are the most difficult quantity to determine and requires dynamic tests [3,4].Numerous approaches have been developed and applied to identify modal parameters of systems and structures in time-frequency domain such as the WT [5,6] and the Hilbert-Huang Transform (HHT) [7,8].
The objective of this study is focuses first to propose a model of a quarter vehicle suspense system and then to identify the parameters of the proposed model of the considered system using WT.
The different sections of this paper are organized as follows: in section II a complete modeling of a 3-DOF quarter suspension system is proposed with a dynamical study followed by a Simulink diagram of the considered vehicle model.Beginning with some mathematical fundaments of the discrete and the continuous WT, denoising and the identification methods are presented in section III.The simulation results are illustrated and discussed in section IV, and a conclusion is drawn in section V.

Proposed 3-DOF Model for the Vehicle Suspension
The suspension system has a big role in design and dynamic of the vehicle.Indeed its goal is to maintain the comfort of the passengers and the stability of the vehicle by supporting the weight of the vehicle, by reducing the disturbances of road and by avoiding the road excitations experienced by the tires from being transmitted to the passengers.
The suspension system is generally reduced to a quarter vehicle model and the majority of the suspension studies use a model with two DOF.In our study, we propose a three DOF model including the seat and drive mass ( s m ), chassis mass ( c ) and the tire mass ( t ).The proposed suspension system is modeled by interconnecting three mass-spring-damper systems as shown in Figure 1.

m m
Assuming that each sub-system is represented by a mass connected to parallel arrangement of a linear spring and a linear viscous damper, the three DOF are the vertical motions of: the seat, the chassis and the tire.
The model shown in Figure 1 includes for each mass a stiffness constant k, a damping constant c and a vertical displacement z.The vertical displacement r z denotes the change of road surface elevation.
Assuming that the tire does not break away from the road surface and that the displacements are measured form the static equilibrium positions, the linear equation of motion for the seat mass is: Similarly, the equation of motion of the chassis mass is given as: and that of the tire mass is : These three differential equations can be written in the matrix form as follows: M is the inertia matrix, D the damping matrix and S the stiffness matrix.
From these dynamical driving equations, we can derive the state space representation of the system by assuming the state space variables are as follows: As a result, the state space representation of the proposed model is: where   A is the state matrix,   B the input matrix,   C the output matrix,   D the transmission matrix, U the input of the system and Y its output.
1 0 0 0 0 0 0 1 0 0 0 0 where s represents the scale or the dilation parameter and  represents the time or the translation parameter.
Otherwise, the wavelet functions are able to adjust themselves as the transform is computed for each component of the analyzed signal: the scale parameter indicates the level of analysis.Small values of s provide a local (or high frequency) analysis while large values correspond to large scale (or low frequency) analysis.Changing the  parameter moves the time localization center of each wavelet.
From the equations of motion ( 1), ( 2) and ( 3) we can also deduce the following diagram block of the proposed model (Figure 2).
These scale and translation parameters can vary in a discrete or continuous way: it is what makes the difference between the Continues Wavelet Transform (CWT) and Discrete Wavelet Transform (DWT).

Wavelet Transform Analysis
Wavelet Transform (WT) is conducted in a manner similar to Short-Time Fourier Transform (STFT).However, Wavelet Transform is superior in the sense that its window function, the analyzing wavelet, is scaled (or dilated) in addition to being translated in time.This analyzing wavelet function is often called the mother wavelet because it gives rise to a family of wavelets through the dilations and translations.A generalized wavelet family In this study, we stated only some keys equations and concepts of each approach, more rigorous mathematical treatment of this subject can be found in [9,10].

DWT for Denoising
For WT in the discrete form, the widely used scale and translation parameter settings that create an orthonormal bases are: Analogous to the FT, the DWT calculates wavelet coefficients by taking the inner product of an input signal   f x with a function, that is in this case the wavelet family , j k  .DWT is designed to give good time resolution and poor frequency resolution at high frequencies and good frequency resolution and poor time resolution at low frequencies.A signal is decomposed using the DWT into two sets of coefficients called approximations coefficients represent low-frequency signal components defined by Vetterli et al. [11]: and details coefficients represent high-frequency signal components defined by: This decomposition, recursively repeated through filtering and down-sampling operations using low-pass and high-pass filters (Figure 3), breaks the signal into lower level coefficient sets.Without manipulating these coefficients, the original signal can be reconstructed exactly through the inverse DWT defined by: where   k g and  are respectively the high-pass and the low-pass filters (Figure 3).

 k h
Manipulating wavelet coefficients prior to signal reconstruction changes the original signal.The original signal can be modified, enhanced or denoised through various coefficient manipulation operations.In our work, the DWT will be used to denoise the response of the system on the road to identify his modal parameters using the CWT method.

CWT for Modal Identification
The Continuous Wavelet Transform (CWT) is defined as the convolution of the signal and the scaled-shifted versions where * , s   is the complex conjugate of the analyzing function , s   .
One of the most widely used mother wavelet for parameter identification is the well-known modified Morlet wavelet defined by [12][13][14] : where 0  is the central pulsation of the Morlet mother wavelet and N is parameter introduced to control the shape of the basic wavelet intend to offer a better compromise in terms of localization, in both time and frequency for a signal.
The optimal value of N is determined by minimization of the wavelet coefficients entropy [13] defined by: For a linear damped multi-DOF system with real modes, the free-decay response and impulse response can be used for modal identification.Hence, a simulated freedecay response with p modes is employed as the analytical signal to investigate the CWT based time-frequency decompositions.
Suppose the signal is given by: where 0k A is the amplitude of the kth mode, 0k  is the phase lag, k  is the damping ratio and k  is the natural pulsation.Since the wavelet transforms are linear representation of the signal, it follows that the CWT of this multi-components signal is: Thus, the CWT works as a time-frequency filter and for a fixed value k s of the scale parameter s, which maxi- , only the mode associated with k s gives a significant contribution in the CWT, while all others are negligible: These regions where the amplitude of the CWT is maximum are called the ridges.By extracting the values of the wavelet coefficients along the ridge yields the skeleton of the CWT.So, the CWT of each separated mode is defined by a skeleton expressed as follows: which implies that the CWT is able to decompose a multicomponent signal into separated modes and to represent them of a complex-valued signals each one defined by amplitude   k A  and a phase angle Note that the phase angle of this skeleton is: where a  is the slope of this linear curve, and by applying logarithmic operator to the amplitude of the skeleton, we obtain: Consequently, the natural pulsation k   and the damping ratio k   can be defined by combining the slopes of A specialized program has been developed on the MATLAB numeric computing environment for estimating these parameters from the system response.This program includes the following procedures: Step 1: Calculate the wavelet entropy in order to determinate the optimal value of N corresponding to the optimal Morlet wavelet.
Step 2: Transforming the time signal (response of the system) into time-scale domain using the CWT.
Step 3: Detecting ridges and extracting corresponding skeletons.
Step 4: Calculate modal parameters according to equation (21) by least-square interpolation

Simulation
In order to show the usefulness the WT in suspension study, the proposed 3-DOF suspension system was created using Simscape blockset library of Matlab software and the correspondent Simulink model shown on Figure 4 was implemented with the simple configuration given in Table 1.
In the simulation study, the proposed model is exited by an impulse signal and three cases of identification parameters are developed and illustrated.

Case 1: Identification from Original Response
The original response signal is shown in Figure 5.For this response he CWT of was firstly calculated for 1 N  as shown in Figure 6(a) and Figure 6(b) gives a crosssection of the CWT at different times.According to this time-scale representation, we cannot see clearly the three modes of the considered system.For this reason, we have to calculate the variation of the wavelet entropy.By increasing N from 1 to 50, the minimum value is obtained for N = 22 (Figure 6(c)), and the optimized CWT of the response corresponding to this value is shown in Figure 7(a) in which the three modes of the considered system can be easily observed.
According to the last study, each mode can be isolated by ridge extraction and used to estimate natural pulsation and damping ratio.The identification results are listed in Table 2 and shows that CWT method correctly identifies the natural pulsations and the damping ratios corresponding to the three      modes.We note that the damping ratio is slightly less accurate than the natural pulsation because its value is estimated from that previously estimated for the pulsation which causes a propagation of error.

Case 2: Identification from Noised Response
In this case, two white Gaussian noises with Signal-to-Noise Ratio (SNR) equal to 20 dB and 2 dB were added to the original response in order to simulate CWT performances under noisy conditions.Signals with the two levels of noise are shown on Figure 11 and Figure 12.
The same developed identification procedure is applied (Figure 13).By minimization of the wavelet entropy, the optimal value of N is the same as in the first case.
The CWT has been applied successfully at the two noise levels.Figure 14 gives the identification of the first mode in the case of SNR equal to 2 dB.The results are not as good as in the case without noise, but they are acceptable as shown in Table 2.

Case 3: Identification from Denoised Response
In this case, the DWT-based denoising process was applied to the previous two noised signals using Symlet wavelet and the troubles components affecting the signal were deleted as shown in Figure 15.Also in this case the identification was made following the same procedure.The results are similar to the previous and that is why only the identification of the first mode is represented in Figure 16.
Table 2 presents a comparison between actual and estimated parameters.In order to highlight the results obtained in each case, it useful to introduce a simple measure of the error estimation for a given modal parameter estimated estimated via the CWT and compared to the exact value actual .Therefore, the Mean Absolute Error (MAE) which is an average of the absolute errors actual estimated i deduced in a case with the three modal parameters of the system.This MAE is given by ( 22) and is introduced at the last column         The simulation results obtained give a very accurate estimation of the actual parameters.This verifies the effectiveness of the denoising and identification method based on DWT and CWT for the proposed model of suspension system.

Conclusions
We have proposed in this paper a 3-DOF model of a quarter-vehicle suspension.After a dynamic study and the modeling of the system, we have investigated the perform-ances of the Wavelet analysis in the identification of pulsation and damping parameters of the system.
Due to its time-scale representation of signals, and the location of a mode by a ridge and the corresponding skeleton, the CWT was able to identify the three modes of suspension system in all cases.
In the case of poor resolution, the first two modes were not distinguishable, for this a compromise between frequency and temporal resolutions has been established by calculating the entropy of wavelet coefficients, and was able to optimize the CWT to better locate these particular modes.
The important advantage of this approach is that it is not too sensitive to noise.Indeed, for all different noise levels, the modal identification was made with low error rates.The simulation results obtained underlined the accuracy and the efficiency of the developed method in the cases where the response signals are denoised by the DWT, and even in the presence of noise.
Finally, we can conclude that wavelets are a powerful tool for the modal analysis in vibration applications especially for suspension systems.

Figure 1 .
Figure 1.Proposed 3-DOF model of a quarter-vehicle suspension system.


wavelet family (7) then becomes: , j k   s   described in the normalized form is:

Figure 2 .
Figure 2. Diagram Block of the proposed 3-DOF Suspension System.

Figure 3 .
Figure 3. Signal decomposition using DWT and reconstruction using inverse DWT.

Figure 7 ( 1 ,
b) gives the three related ridges representation at different times, from this plot we can deduce the scale parameters corresponding to each mode: the first mode is localized in the neighborhood of wavelet coefficients to each obtained ridge gives the corresponding skeleton.

Figure 8
shows the plots of the real part (Figure 8(a)) and the imaginary part (Figure 8(b)) of the skeleton related to the first mode.And as previously established by Equation (21), an interpolation of the phase (Figure 8(c)) and the modulus (Figure 8(d)) of this first skeleton allows the estimated natural pulsation and the damping ratio of the first mode using a linear least-square fit procedure.Similarly, Fig- ure 9 and Figure 10 present the identification of the other modal parameters.104 s 

Figure 5 .
Figure 5. Original response of the suspension system.

Figure 6 .
Figure 6.CWT of the original response signal.a) CWT for N = 1; b) Cross section at different instants for N = 1; c) Minimization of the wavelet entropy.

Figure 7 .
Figure 7. CWT of the original response signal.a) CWT for N = 22; b) Cross section at different instants.

Figure 8 .
Figure 8. Identification of the skeleton corresponding to the 1 st mode.a) Real part; b) Imaginary part; c) Phase; d) Envelop.

Figure 9 .
Figure 9. Identification of the skeleton corresponding to the 2 nd mode.a) Real part; b) Imaginary part; c) Phase; d) Envelop.

Figure 10 .
Figure 10.Identification of the skeleton corresponding to the 3 rd mode.a) Real part; b) Imaginary part; c) Phase; d) Envelop.

Figure 14 .
Figure 14.Identification of the skeleton corresponding to the first mode (SNR = 2 dB).a) Real part; b) Imaginary part; c) Phase; d) Envelop.

Figure 16 .
Figure 16.Identification of the skeleton corresponding to the first mode of the denoised response (SNR = 2 dB).a) Real part; b) Imaginary part; c) Phase; d) Envelop.