Three-Axes Mems Calibration Using Kalman Filter and Delaunay Triangulation Algorithm

. MEMS-IMUs are widely used in research, industry, and commerce. A proper calibration technique must reduce their innate errors. In this study, a turntable-based IMU calibration approach was presented. Parameters such as the bias, lever arm, and scale factor, in addition to misalignment, are included in the general nonlinear model of the IMU output. Accelerometer error parameters were estimated using the transformed unscented Kalman flter (TUKF) with triangulation algorithm is suggested for calibrating inertial measurement unit (MPU6050) three-axes accelerometer. In contrast to the present methods, the suggested method uses the gravitational signal as a constant reference and necessitates no external equipment. Te technique requires that the sensor be positioned in a rough orientation and that basic rotations be adopted. Tis technology also ofers a quicker and easier calibration. Comparing the experimental fndings with other works, Allan deviation shows signifcant improvements for the bias instability, where a bias instability of (0.116 μ g) is achieved at temperatures between ( − 15 ° C) and (80 ° C).


Introduction
MEMS-IMUs have found widespread use in a variety of felds, including research, industry, and business.A suitable calibration method needs to decrease the mistakes that are inherent in them.Tese sensors found usage in a wide range of applications, including navigation and positioning systems, educational software, medical electronics, and more [1][2][3].Researchers have also recently examined the usage of MEMS devices in user interfaces [4][5][6].Despite the benefts mentioned above, MEMS sensors' measurements are less precise than those from tactical-grade sensors because of biases, noise, scalability errors, nonlinearity, thermal changes, and misalignment errors.It is clear that these sensors can only produce signifcant precision if a calibration technique is used to correctly ofset their inaccuracies.
When calibrating an instrument, it is necessary to compare its outputs to well-known reference data in order to identify the coefcients that, for a variety of output values, cause the output to coincide with the reference data [7,8].
Traditionally, mechanical platforms and rate tables are used to spin the sensor to accurately control orientation and rotation rates during the calibration of the inertial sensors [9].However, due to their high cost and need for a controlled environment, these are neither practical, economically, nor appropriate for calibrating MEMS sensors.However, the aim of this paper is to calibrate and minimize the MEMS sensor output error for the three axes (X, Y, and Z) to be with less bias error.

Literature Review
Alternate calibration techniques have been investigated over time by researchers to get around the challenges of employing pricey mechanical platforms.An accelerometer and gyroscope calibration method that has been fled was suggested by Ferraris et al. [10].Using a gravity-based 6position approach, they calibrated the accelerometer.For the sensor to be perfectly aligned with gravity, the sensor enclosure was meticulously built.Static intervals were used to calibrate the gyro's bias, and full rotations around each axis were done to determine the scale and misalignments.
Te rotations could only be applied in 2 separate ways because their approach needed an external reference block.A calibration technique using a gyroscope and an accelerometer was put out by Skog and Händel [11].Te accelerometer is calibrated using a technique that does not need the sensor to be positioned in a specifc way and instead uses gravity as a trustworthy external reference.With the goal of calibrating the inertial measurement unit (IMU) in the feld, Jurman et al. [12] combined the gyro calibration method from [9] with the gravity-based accelerometer calibration strategy.Olivares et al. [13] ofered a modifed approach that would do away with the need of using a turntable for such gyroscope calibration.First, in order to calibrate the accelerometer, they employed a gravity-based 6-position technique, and then, in order to calibrate the gyroscope, they used a previously calibrated accelerometer in conjunction with a bicycle wheel.
Ma et al. suggest a brand-new, used the applicable modulated feedback method to reduce the overall bias instability of the closed-loop system, MEMS capacitive accelerometers [14].Using modulated feedback allows for the separation of the 1/f noise that is brought into the system by the electrical interface.Comparative experiments have been set up, and a signifcant improvement in the bias instability is indicated by the Allan variance.Te bias instability has been reduced to under 13 g [14].A pair of double-ended tuning forks (DETFs), a pair of pair beams, two microlevers, proof mass, and even a quartz frame are all components of the centrally symmetric quartz accelerometer assembly that was designed by Chao Han et al. [15].In order to maximize the utilization of sensor space and enable sensor downsizing, microleverages and DETFs are placed around the chip perpendicular to one another.Because of this, the velocity random variable is calculated to be 0.84 g/Hz, and the bias instability is calculated to be 3.05 g.
Zhao et al. [16] proposed an automatic gain-control (AGC) for self-compensation method that combines a polarization source and reference.DC and AC polarization scenarios are covered, and it eliminates any 1/f noise and drift that the polarization source itself introduces.According to the experimental fndings, the suggested AC polarizing and self-compensation technique can increase the bias stability of the investigated SOA from 2.58 g to 0.51 g.
Gravity-based techniques can be used in the feld with ease because they do away with the need for additional equipment in the case of an accelerometer.In contrast, to hand, calibrating a gyroscope with a turntable or other external device restricts the sensor's ability to be reliably calibrated inside the feld and increases overall costs.An accelerometer, gyroscope, and the magnetometer are all integrated onto a single chip in the MPU6050 sensor that we used as an example.Inertial sensor businesses are increasingly attempting to merge numerous sensors onto the single chip in order to further decrease both the size and cost of their products as MEMS technology continues to progress.
In this research, we suggest a calibration algorithm for an inexpensive IMU (MPU6050) that may be used in the feld without any additional hardware.Te process also requires less computation and can be applied physically.Te document is structured in the manner described below.Te output model for the accelerometer is developed in Section 2. Te calibration method is discussed in sections 3-5.Te Delaunay triangulation technique is demonstrated in Section 6, the experiment setup is illustrated in Section 7, and the conclusion is presented in Section 8.

Output Model for the Accelerometer
Te modelling of the accelerometer output is carried out in the same manner as the modelling of the gyroscope output.Te scale factor, bias, and misalignment errors are also present in the gyro error components, which are comparable to those found in the static error components.In addition to this, the accelerometer has a problem with the error caused by the lever arm.
Te following is one way to explain the information provided by the accelerometers when M a ≜ G a − E a ′ , G a ≜ W T a are defned, respectively: and where M a and G a , respectively, indicate the actual and nominal confgurations of the sensors, Λ a scale factor matrix of errors, W a is the nominal direction matrix, V a and U a the orthogonal orientation vectors, b a is the vector of bias, and η a is the measured noise.Te ideal acceleration vector may be represented as follows [17,18] if we assume that the accelerometer was positioned at the distance from the center of the turntable. where T indicates the acceleration that would be taking place at the table's center of rotation in relation to the body frame T is a measure of the angular velocity of the turntable's acceleration, and d is defned as d � [d x d y d z ] T is the vector that is used to represent the lever arm.Now, this (3) is recast in a more succinct form as seen in the following equation: where 2 Applied Computational Intelligence and Soft Computing (5) If we substitute ( 4) for ( 1), then the outputs of the accelerometer will be in which (7)

TUKF Calibration Method
Te sensor output is nonlinear.As a result, an estimation of the error parameters has to be carried out using a nonlinear flter [19].TUKF, a new nonlinear flter, was recently created and published in [20].In contrast to the transformed unscented Kalman flter (TUKF), the TUKF does not sufer from linearization error since the TUKF approach is founded on the sigma point [20].Tis allows the TUKF to perform more accurately.In addition, in contrast to other flters that are based on sigma points, such as the unscented Kalman flter, this flter does not need any tweaking (UKF).As a result, the TUKF was used in this investigation in the capacity of an estimator (observer).For TUKF to be able to execute the estimate of unknown parameters, it requires a dynamic model in addition to observations.All readings and process model of the calibrating flters will be provided in the sections that follow.

Model of Measurement and Process for the Accelerometer
For the purpose of representing measurement, which is contingent just on the output of the system of the accelerometer (1), the following equation may be used: (8) Te Markovian model of the frst order is a reasonable assumption for η a [21], given that where X ai stands for the state variable of a frst order Gauss-Markov process, C xai stands for the coefcient, and v ai may be thought of as white Gaussian noise with a mean of zero and a covariance of E v a i v a i T   � R a i .If we assume that the calibration was fnished in a reasonably short amount of time, then we will be able to describe the error parameters as constant processes.Tis is because we will have assumed that the time it took to complete the calibration was relatively short.Because of this, it would seem that the error parameters are immune to the efects of variations in temperature.Te dynamism of the process in the state space may be described as follows: where the calculation for A a , B a an is as follows: Moreover, there is an angular acceleration α indicated in the form of (8).In general, there are two diferent methods that may be used to calculate this signal.Te frst problem is associated with the twice diferentiated output of the encoder.Te existence of noise, on the other hand, causes this method to result in an increase in the magnitude of the inaccuracy.Another potential strategy involves the angular acceleration of the system as an additional unknown variable in the state space.In order to accomplish this goal, a dynamic model for angular acceleration should be explored.For the purpose of this investigation, a Gauss-Markov dynamic of the frst order was used to model angular acceleration.An angular acceleration is basically just a straightforward function of an angular velocity like that.
where a, b, and c are constant coefcients that are estimated in the same way as any other coefcients.Te process noise's covariance matrix is denoted by the equation Te order of 27 is assigned to the process dynamics of accelerometers (9), on the other hand.Te complete order of process dynamics is increased to 30, which are employed in the TUKF for the purpose of computing the overall error characteristics of accelerometers.Tis is done via including its angular acceleration (11), which brings the total up to 30.

Formulation of TUKF
In [20], a comprehensive description of TUKF is provided.Just its mechanism is described in this section, while the Applied Computational Intelligence and Soft Computing algorithm that governs it is given out in Table 1.Take into consideration a process model of the nth order, the measurement models for which would be as follows: where u represents the input and x k represents the k th sequence of state variables, ω and σ are the measurement Gaussian white noises process, each one with covariance matrices values Q k and R k .Te TUKF approach for the described process and also the measurement dynamics consists of three components: initialization, prediction, then the measurement update.Sigma values, an estimate of the state variables ( x − o ), and an error covariance matrix (  P − o ) are all set to zero at the beginning.In order to make inferences about the future, the process model is used to compute estimates for state variables ( x − k ), as well as the square root of the covariance matrix.Te square root of the covariance matrix and the state variables ( x + k ) are then updated based on the actual values of the measurement and its model.
According to the information shown in Table 1, E [.] designates the expected value operator, qr (.) stands for the QR decomposition, chol (.) is regarded here as the Cholesky factor decomposition process, and I n is an identity matrix of the nth order.Te calibration procedure was carried out on an MPU6050 sensor in order to determine whether or not the approach described in this article is indeed practicable.
IMU and Arduino both contributed to the data collection process.Te inertial measurement unit (IMU), which consists of three orthogonal gyroscopes and accelerometers, is housed in the turntable in such a manner that allows us to take into account the following matrices while confguring it: To fgure out how well the proposed method works, the results of calibration using least squares (LS) [22] were compared to those of calibration using the proposed method.

Delaunay Triangulation Algorithm
Tis section examines how the transformed unscented Kalman flter's output is impacted by the Delaunay triangulation process.Figure 1 shows the ideal inputs to the IMU in three axes.Figure 2 shows the corrupted outputs at the zero acceleration condition of IMU, but there is a bias in the (Z) axis because of the 1 g gravity applied on the Z-axis in addition to the noise source such as Brownian noise.With the help of this algorithm, we can improve the accuracy of the Kalman flter's output values.A triangulation divides a polygon neatly into triangles, making it possible to; for example, quickly calculate the area or even a guarding of a polygon.Utilizing trigonometry (T) for interpolation is another popular application scenario: Consider an extending a function (f) "fairly" and continuously with (P), defned upon those vertices of the polygon (P • ).
Afterward, locate the triangle (t) with a (T) which contains a point p ϵ (P • ).As (p) is a convex expression.With the vertices v 1 , v 2 , v 3 of t combined ( 3 i λ i v i ), we can simply interpolate the function values using the same coefcients (f(p) �  3 i λ i (fv i )).Te procedure is to regard each of the three Kalman flter, output values as a triangle composed of three vertices, with one serving as the triangle's head and the other two as its base, as shown in Figure 3, if we do this, the center of the triangle, which will serve as such an alternative output value since it is closest to the real value as it can be seen in Figure 4. Figure 4 illustrates how, as a result of a Kalman flter output enhancement, the Z-axis now has the lowest amount of noise.

Experimental Results
Te TUKF and maximum likelihood (ML) errors of the calibrated accelerometer signal are shown for just one direction in Figure 5, which was acquired during the validation manoeuvre.
Te Allan variance approach may be used to ascertain the properties that are underlying random processes which are accountable for the generation of the data noise.By carrying out certain operations throughout the whole of the data, this method may be used to describe a wide variety of error terms that may be present in the readings from the inertial sensors.
For describing the methodology for evaluating the velocity random walk (VRW), bias stability, and acceleration random walk (ARW), suppose that there are N successive data points and that each of them has a sampling time of t o .Each of the n successive data points (with n < N/2) that come together to form a group is considered to be a cluster, as is seen in Figure 6.
A time, denoted by the letter T, that is equivalent to nt 0 is connected to each cluster.In the event when the instantaneous output rate of the inertial sensor is denoted by Ω(t) equal to, the cluster average is defned as [23]: where Ω k (T) stands for the cluster average of the output rate for a cluster that begins at the k th data point and has a total of n data points.Te succeeding cluster average may be defned as the following equation: where t k+1 � t k+1 + T. Tis Allan variance with length T may thus be defned as [23] σ It is abundantly clear that a limited number of clusters of a predetermined length T may be generated from any given fnite number of data points N. As a result, equation provides an estimates of the quantity σ 2 (T), the accuracy of which is 4 Applied Computational Intelligence and Soft Computing Table 1: Te outline of transformed unscented Kalman flter.
(1) Initializing Variables pertaining to the state and the square root (sr) of the error covariance matrix (2) Prediction (time update) TUKF cubature points after transformation

Sample propagation points
For i� 1, 2, . . ., 2n Update of Sr-innovation covariance matrix Update of cross-covariance matrix Updating the Kalman gain Updating the states Update Sr-error covariance matrix Applied Computational Intelligence and Soft Computing dependent on the number of separate clusters of a constant length that can be constructed.A sensor's capacity to maintain a consistent output is quantifed by the Allan variance.As a result, it must be connected to the statistical features of the inherent random processes, which infuence how well the sensor works.Te power spectral density (PSD) of the noise components in the original data are connected to the Allan variance that is acquired by conducting the procedures that have been detailed.Te Allan variance is related to the two-sided PSD, which is given by where S Ω (f) is the power spectrum density of the random process.
Te white noise spectrum of the accelerometer output is characteristic of the angle (velocity) random walk that is being performed.Te PSD is exemplifed by the following equation [23]: N � coefcient of the angle random walk.Te equation that was just presented yields a line that has a slope of (−1/2) when it is plotted on a log-log graph of σ 2 (τ) vs τ.Te units of N are (m/s/√h).Applied Computational Intelligence and Soft Computing Putting these values into the original equation for PSD and integrating them gives the following results: Te red noise (also known as Brownian noise) spectrum of the accelerometer output is a characteristic feature of the rate random walk.Te following elements constitute the PSD: K is the coefcient of the rate random walk.Putting these values into the original equation for PSD and integrating them gives the following results: Te equation that was just presented yields a line that has a slope of (1/2) when it is plotted on a log-log graph of σ 2 (τ) vs τ.Te units of K are (m/s/s 2 ).
Te pink noise (also known as ficker noise) spectrum of the accelerometer output is a characteristic feature of the bias instability.Te following elements constitute the PSD: where B is the coefcient of the bias instability, f o is the ficker noise corner frequency.Putting these values into the original equation for PSD and integrating them gives the following results: where x � πf o τ, Ci � is the cosine-integral function.

Applied Computational Intelligence and Soft Computing
When plotted on a log-log graph of σ 2 (τ) vs τ, the equation that was just discussed results in a line that has a slope of 0 degrees.B has the units of (m/s 2 ).
Figure 7 shows the setup for the data reading from the MPU6050 put in the oven chamber.50000 pieces of data then recorded after that.Tree measurements from the (IMU), taken while the sensor was kept still in a fxed location, made up each set.Calculated Allan deviation for the accelerometer's three axes is displayed in Figure 8 which give us information about MPU6050 sensor to be more stable after applying the proposed algorithm.As illustrated in Table 2, the bias instability has improved compared with the previous works [24][25][26] for three axes, which shows a very competitive outcome for the performances of bias stability and acceleration random walk in the temperature range (−15 °C to 80 °C) is obtained in this work.

Comparison with Previous Works
In this section the comparison with other works has been presented.In [24], an enhanced robust flter with a double state model by using the chi-square distribution of the square of the Mahalanobis distance as our theoretical foundation is presented.By using the upgraded robust Kalman flter, achieve position gains of up to 33 and 30 percent, respectively, in the north and east components has been achieved.Te robust flter has given a 57% reduction in the amount of azimuth error, and then the upgraded robust flter has the potential to provide greater performance.Te azimuth error has been decreased in loosely coupled and tightly coupled systems, however, the improvement in pitch and roll accuracy was smaller than the increase in azimuth accuracy owing to limited observability with the PC for Data Oven Chamber experimental trajectories being performed on a stationary track.Ten, the performance of this modifed robust flter was satisfactory, and it was suitable for use in the MEMS application.
Figure 9 shows the pitch and roll for the accelerometer while Figure 10 is for the azimuth in reference [24] after utilizing Kalman flter and the robust flter.
In [25], a transformed unscented Kalman flter TUKFbased calibration approach is used and a nonlinear model was suggested in this work as a means of estimating the sensor model parameters of a three-axes gyroscope and an accelerometer.Te model was developed using TUKF.Te TUKF dynamics were constructed using an expanded state space model of the sensor that was described here.Te real signals were collected from a tri-axis turntable as well, so that the fndings of the experiment could be compared to the real ones.Te estimated parameters obtained from the TUKF and the calibrated signals obtained from the least square (LS) technique are used in the process of adjusting the sensor output signals.
Figure 11 is the X axis error for the accelerometer for the reference after utilizing TUKF [25].
In [26], an adaptive and robust theory was used in the development of the Modifed Sage Husa Adaptive Robust Kalman Filter (MSHARKF), which is presented in this study.A new adaptive scale factor is added to the Adaptive Robust Kalman Filter (ARKF) algorithm in the MSHARKF algorithm, and the ARKF algorithm's state is modifed after each iteration.MSHARKF is implemented in the measurements of the MEMS IMU in order to reduce the amount of random noise and drift inaccuracy.It can be seen from the Allan variance analysis that both angle random walk (ARW) and bias instability (BI) are two dominating random noises in the raw MEMS inertial measurement unit (IMU) data.
Figure 12 shows the Allan deviation for the reference [26].

. Conclusion
Tis study investigated the calibration of devices such as the microelectromechanical IMU MPU6050 accelerometer, also aiming to develop a TUKF using triangular calibration algorithm that can be used in the feld, outside of a lab, without the need for additional equipment.Te method described removes restrictions on the orientation of the sensor during calibration by using just the gravitational feld of the earth as a reference.Te practical readings taken from the MPU6050 accelerometer at each (1) second from all the axis, which have been achieved using just a standard IMU, a custombuilt IMU, as well as an aviation graded table, supported the practicality of the proposed method which gave us acceleration random walk of 0.71 μg and a bias of 0.116 μg.As can be seen from the results, the bias error minimization for the three axes (X, Y, and Z) has been reached, in addition, the Markov or the sinusoidal noise and exponentially correlated noise can be calculated in the future work.

Figure 2 :
Figure 2: IMU corrupted and biased outputs of the three axes (X, Y, and Z) at zero acceleration condition.

Figure 1 :
Figure 1: IMU ideal inputs of the three axes (X, Y, and Z).

τ = nt 0 τ = nt 0 Figure 6 :Figure 5 :
Figure 6: Diagram of the data structure that is used in the Allan variance algorithm.

Figure 7 :
Figure 7: Photographs of the experiment setup.

Figure 8 :
Figure 8: Allan deviation after the (TUKF) and the Delaunay triangulation of the three axes (X, Y, and Z).