Error model of geomagnetic-field measurement and extended Kalman-filter based compensation method

The real-time accurate measurement of the geomagnetic-field is the foundation to achieving high-precision geomagnetic navigation. The existing geomagnetic-field measurement models are essentially simplified models that cannot accurately describe the sources of measurement error. This paper, on the basis of systematically analyzing the source of geomagnetic-field measurement error, built a complete measurement model, into which the previously unconsidered geomagnetic daily variation field was introduced. This paper proposed an extended Kalman-filter based compensation method, which allows a large amount of measurement data to be used in estimating parameters to obtain the optimal solution in the sense of statistics. The experiment results showed that the compensated strength of the geomagnetic field remained close to the real value and the measurement error was basically controlled within 5nT. In addition, this compensation method has strong applicability due to its easy data collection and ability to remove the dependence on a high-precision measurement instrument.


Introduction
The geomagnetic-field is an inherent physical field of the earth. For navigation and orientation, the geomagnetic-field has several advantages. It is strongly self-contained and has antiinterference and all-weather features. Due to its small size, geomagnetic navigation equipment has low energy consumption and low cost. For these reasons, geomagnetic navigation is increasingly attracting the attention of domestic and foreign scholar and has become a highlight in the field of navigation around the world [1]. The accurate acquisition of geomagnetic information is the precondition for geomagnetic navigation and the basis for achieving highprecision geomagnetic navigation at the same time. Currently, a navigation carrier mostly measures magnetic information of the space via a strap-down magnetic sensor. However, as the spectral range of the geomagnetic field is very wide, the measurement of the geomagnetic field is subject to interference, causing the output of the magnetic sensor to not only contain the geomagnetic information needed for navigation and orientation but also contain various interferential information. To promote the accuracy of geomagnetic navigation and meet the requirement of high-precision geomagnetic navigation for highly accurate geomagnetic-field information, the output information of the magnetic sensor must be corrected [2][3]. a1111111111 a1111111111 a1111111111 a1111111111 a1111111111 By referring to the research findings of domestic and foreign scholars, this paper divided the sources of the measurement error of magnetic sensors into manufacturing error, installation error and interference from the external electromagnetic environment [4][5][6][7][8][9][10][11][12][13]. The study of an existing magnetic-interference model found that the simplified model cannot clearly describe the sources of measurement error of magnetic sensors. A 12-parameter error model including manufacturing error and installation error was built [4]. An error model including manufacturing error and interference from the external electromagnetic environment was built [5]. The error model built in [6][7][8] only mentioned the sensor's manufacturing error. The error model built in [9][10][11][12][13][14] only considered the interference of the external electromagnetic environment. It is necessary to note that these models do not address daily variations of the geomagnetic-field in the study of the interference of the external electromagnetic environment. As the major component of magnetic-field variation in quiet periods of geomagnetism, the accurate measurement of the geomagnetic daily variation field has a non-negligible effect on the accurate measurement of the geomagnetic-field [15]. To ensure the accuracy of measurement data, the model must be corrected based on the geomagnetic daily variation field.
Regarding compensation for magnetic-field interference, the existing methods include the multi-pose method, ellipse-fitting method, ellipsoid-fitting method, neural network method, deviation-compensation method, Tolles-Lawson equation-based magnetic-measurement compensation method and Kalman filtering. The "ellipse assumption method of two-dimensional magnetic field measuring track" was proposed [4]. The core idea of this method is to convert the problem of compensation of magnetic measurements to a problem of parameter estimation. However, two estimation processes contained in this algorithm were non-linear parameter estimations and such estimation processes are extremely complicated and involve large amounts of calculations. In the process of estimating the ellipse parameters in the first step, since the ellipticity of the conic section cannot be ensured, the ordinary least-square method is not applicable and this further intensifies the complexity of algorithm. In addition, this method can only solve for 9 of the 12 parameters of the constructed model; the estimation of the remainder of the parameters requires the help of accelerometer and reference information provided by GPS, which intensifies the complexity of the algorithm structure and increases the number of calculations.
A method based on ellipse fitting was proposed in [11][12] to use a measured value to fit an ellipse and find ellipse parameters, and obtain the error model's parameters in accordance with the relation between the ellipse parameters and error model. The algorithm of this method is rather complicated and it strictly requires the carrier to move only in the horizontal plane during data acquisition, therefore this method is only applicable to two-dimensional space and is greatly limited in practical application.
References [6,9,[16][17][18][19] put forward ellipsoid fitting algorithm, which calculates ellipse parameters with least-square method or iterative algorithm and indirectly performs error compensation based on the characteristic that the track of the geomagnetic field subjected to various interferences in a fixed space is ellipsoid. Although the ellipsoid fitting algorithm is able to obtain an ideal compensation effect, there is such a problem that iterative algorithm doesn't only require a large number of calculations but also relatively accurate initial conditions, and there is likely to be the problem of matrix singularity in the solving process of least-square method, making it unable to obtain the correct solution.
References [18][19] both proposed to use neural network algorithm for the compensation of carrier's magnetic field: reference [18] applied neural network algorithm in the parameter estimation of ellipsoidal model, and reference [19] compared neural network algorithm with ellipsoid fitting algorithm in terms of the result of the error compensation for magnetic measurement. The final results of both studies revealed that the neural network algorithm was able to achieve a good effect with respect to compensation accuracy, but its effect was not better than that of the ellipsoid fitting algorithm. Furthermore, in consideration of the long time it takes in model training, the application of this method is greatly limited.
The compensation effects of ellipse fitting algorithm, ellipsoid fitting algorithm and neural network algorithm are all sensitive to the measuring accuracy and noise of sensor, so additional measures are required to solve this problem; for example, reference [18] adopted wavelet de-noising method to eliminate the effect of noise, which undoubtedly increased the workload of compensation. As a recursive algorithm, Kalman-filter algorithm is able to achieve online operation and has a certain adaptability to the errors in initial states and a good ability to suppress noise; therefore it has been extensively used in the estimation and calibration of various models' parameters [20][21][22][23][24][25]. Although the extended Kalman-filter algorithm has been applied in geomagnetic navigation technology to a certain degree, the applications are mostly focused on the fusion of navigation information [26][27], and there lack the reports on error compensation for magnetic measurement.
During data acquisition, ellipse fitting method requires the carrier to move only in a twodimensional space [11]; multi-pose method requires rotating a specially made hexahedral apparatus to obtain 12 or 24 different poses [28]; despite the simplicity of data acquisition, ellipsoid fitting method is unable to reach a stable solution in the end because constraint matrix singularity will be caused when the ellipsoid where the data acquired are approaches a sphere [16][17]. By contrast, the method proposed in the paper requires the sensor only to rotate in a three-dimensional space during data acquisition, and there are no additional constraints, so it has better adaptability.
On the basis of study done by others, this paper improved the existing magnetic-sensor measurement-error model, introduced the geomagnetic daily variation field and built a new model that includes the magnetic sensor's manufacturing error, installation error and the interference of the external electromagnetic environment. The new model is able to indicate the source of the measurement error. With respect to the estimation of the model parameter, an extended Kalman-filter with the model parameter as the state variable and the magnetic sensor's output quantity as the observation quantity was designed according to the model's characteristics. The new filter can fully make use of all measured data to resolve the problem of the poor utilization rate of the measured data. In addition, with the increase in the amount of measured data, the estimated value of the model parameter will well converge to the real value to obtain the optimal solution in the sense of statistics. In addition, due to the Kalman filter's good capability to constrain noise, this method can effectively prevent strange solutions in the process of seeking the solution and eventually realize high-precision compensation of measured data.

Magnetic sensor's measurement-error model
This paper will proceed from the magnetic sensor's manufacturing error, installation error and interference of the external electromagnetic environment to study each of these three factors in detail and build their measurement-error models. Finally, a complete model containing the above three factors will be built.

Manufacturing-error model
Manufacturing error, caused by the level of production technology, material properties and many other factors that make the sensor unable to reach the ideal working conditions, mainly includes sensitivity error, zero error and non-orthogonal error. Sensitivity error is the result of the different sensitivities of the magnetic sensor's three measuring axes and can be represented using a diagonal matrix. Zero error occurs because the zero points of the geomagnetic sensor, analog circuit and analog digital converter are not zero and can be represented using a fixed vector form. Non-orthogonal error occurs because the actual directions of the magnetic sensor's three measuring axes are not completely orthogonal and its principle is shown in Fig 1 oxyz represents the ideal coordinates of the geomagnetic sensor's measuring axes. This set of coordinates is an orthogonal one, while the actual coordinates of the geomagnetic sensor's measuring axes, o-x 0 y 0 z 0 , are non-orthogonal. The origin o of the coordinates is the center of the geomagnetic sensor. The ox-axis coincides with the ox 0 -axis; the o-xy plane coincides with the o-x 0 y 0 plane; and the included angle of the oy-axis and oy 0 -axis is α. The included angle between the projection of the o-z 0 -axis on the o-yz plane and the o-z-axis is β, and the included angle between the o-z 0 -axis and its projection is γ [5][6].
Assuming that the magnetic sensor rotates in a certain fixed area and under the condition of only considering the effect of its manufacturing error, the relation between the magnetic output and the real magnetic field of the area is [7][8] where H m represents the output of the magnetic sensor's three axes and H i represents the magnetic-field vector in the sensor's coordinates. respectively, represent the sensitivity matrix and non-orthogonal matrix. s x and s y and s z respectively, represent the sensitivity coefficients of the geomagnetic sensor's three axes and are con- and it can be regarded as a fixed vector. Error model of geomagnetic-field and compensation method

Magnetic sensor's installation-error model
The difference between the actual installation position and the ideal installation position is the main cause of the magnetic sensor's installation error. The transformation process between the two coordinates is shown in Fig 2: o-x 0 y 0 z 0 represents the carrier's coordinates and o-xyz represents the magnetic sensor's coordinates. The three transformation processes completely reflect the transformation from the magnetic sensor's coordinates to the carrier's coordinates [4][5].
Under the condition of only considering the effect of the installation error, the transformation of the geomagnetic field from the magnetic sensor's coordinates to the carrier's coordinates is where H m represents the magnetic-field vector in the carrier's coordinates, H i represents the magnetic-field vector in the sensor's coordinates and E p represents the installation-error matrix. For the same carrier, the strap-down magnetic sensor will not be adjusted in position once it is fixed, therefore, the installation-error matrix E p will not change under normal conditions.

Interference-error model of the external electromagnetic environment
The magnetic sensor's measurement of the geomagnetic-field is likely to be subject to interference from the surrounding electromagnetic environment, which is determined by the characteristic wide spectral range of the geomagnetic-field. The interference from the external electromagnetic environment includes hard magnetic interference, soft magnetic interference, random magnetic-field interference and geomagnetic daily variation interference. Hard magnetic interference refers to the interfering magnetic field that is formed as the hard magnetic material inside the carrier is magnetized by the external magnetic field. As hard magnetic material has high coercive force and a remanence value that can remain unchanged over a long period, the interfering field formed by hard magnetic material remains unchanged in the carrier's fixed coordinates and it can be regarded as a fixed vector. Soft magnetic interference refers to the induced magnetic field that is formed as soft magnetic material inside the carrier is magnetized by the external magnetic field. Different from hard magnetic interference, soft magnetic interference will vary with the variation of the external magnetic field, and the carrier's soft magnetic interfering field can be regarded as the summation of a number of magnetized magnetic dipoles and magnetic-dipole moments [12]. The induced magnetic field is described in detail as follows. Assume that M x ,M y ,M z are the three components of the induced magnetic moment produced by H ix ,H iy ,H iz , respectively, that magnetize soft iron in the sensor's coordinates. The induced magnetic moment in the three axial directions in the carrier's coordinates is where K x ,K y ,K z are the magnetized coefficients, which are non-dimensional constants. At random attitude, the direction of the carrier's induced magnetic moment aligns with the direction of the magnetic field and its three-axis direction coincides with the sensor's threeaxis direction. As shown in Fig 3(A), assume that the point p is the sensor's position, UGV represents the rectangular coordinates of p's plane, og is the vertical line between the x 0 oy 0 and ugv planes, the U-axis is parallel to the x 0 -axis, and the V-axis is parallel to the y 0 -axis. As the sensor is installed in the carrier in a strap-down configuration, the relative positions of the coordinates UGV and the sensor's coordinates ox 0 y 0 z 0 remain unchanged. Meanwhile, spherical coordinates are created with the equivalent magnetic-moment point o as the center and the magnetized magnetic moment M z as the reference to coincide with the z i -axis and z 0 -axis. Thus, the M z direction is aligned with the z i direction. Then, according to the magnetic-field theory of the magnetic dipole [14], we have where r ¼ ffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi and μ is the dielectric permeability, which is an inherent dielectric property.
With the magnetized magnetic moment M x as the reference, keeping coordinates UGV and coordinates ox 0 y 0 z 0 motionless and rotating coordinates ox i y i z i along y i , we obtain coordinates ox 0 i y 0 i z 0 i , making the z 0 i -axis coincide with the x 0 -axis and the M x direction align with the z 0 i direction, as shown in Fig 3(B). Then, where a x ¼ arccos À 1 u r À Á and b x ¼ arctan v d À Á , With the magnetized magnetic moment M y as the reference, keeping coordinates UGV and coordinates ox 0 y 0 z 0 motionless and rotating coordinates ox i y i z i along x i on the basis of the previous rotation of coordinates ox i y i z i , we obtain ox 00 i y 00 i z 00 i , making the z 00 i -axis coincide with the y 0 -axis and the M y direction align with the z 00 i direction, as shown in Fig 3(C). Then, where a y ¼ arccos À 1 v r À Á and b y ¼ arctan u d À Á , Then, at point p, the induced magnetic field F m formed by the carrier's soft magnetic material can be represented as Substituting Eqs (4), (5) and (6) into Eq (7) gives where the coefficient matrix on the right contains 6 unknowns: u, v, d, K x , K y and K 2. These 6 unknowns are all constants for a carrier installed with a strap-down magnetic sensor. Eq (9) can be simplified to the following.
The generation of random magnetic field contributes to the imbalanced variation of current inside the carrier and interference in the operation of onboard radio-transmission equipment. This random interfering field is not a major component of the carrier's magnetic field and can be reduced or eliminated through a rational wiring plan. Under the condition of only considering the interference from the electromagnetic environment inside the carrier, the measured value of the geomagnetic-field vector in the carrier's coordinates can be represented as follows: where matrix A soft represents the interference generated by the soft magnetic material's induction field on the measurement, and it is specifically expressed as follows: where a ij is the interference coefficient of the induced magnetic field in which the i direction is acted on by the j direction; 1+a ij represents the increased multiple of strength of the magnetic field in the magnetic sensor's i-axis direction under the influence of the induced magnetic field; H hard represents the hard magnetic interference and is specifically expressed as H hard = [H hardx H hardy H hardz ] T ; and H w represents the random magnetic field inside the carrier. Eq (11) shows the interference model of the external electromagnetic environment that is generally adopted. This model only considers the electromagnetic-interference field inside the carrier and neglects the geomagnetic daily variation field, which has significant influence on geomagnetic-field measurement. It was noted that the geomagnetic daily variation field is one of the important sources of error for geomagnetic-field measurement and has significant influence on geomagnetic navigation [15]. The geomagnetic daily variation field component can be divided into two parts: the regular variation field and the irregular variation field. The regular variation field includes daily variations and annual variations, and its variation is continuous and independent of time; the irregular variation field includes magnetic storm and substorm, and its variation is sporadic and discontinuous. The regular variation field's interference imposed on the magnetic field can be regarded as a fixed vector. The irregular variation field is regarded as noise in this paper. Therefore, by comprehensively considering the influence of the magnetic sensor's external electromagnetic environment on magnetic-field measurement, Eq (11) can be written as follows.
where H r represents the regular variation field and H ir the irregular variation field in the geomagnetic daily variation field.

Measurement-error model
In summary, comprehensively considering the magnetic sensor's manufacturing error, installation error and the error of the magnetic sensor's external magnetic interference can produce the magnetic sensor's geomagnetic-field measurement-error model: Normally, C is an invertible matrix and the noise ε can be eliminated with the use of a certain device and technology [4], thus the magnetic sensor's measurement-error compensation model can be described as where the inverse of C can be described as C À 1 ¼

Filter design
It has been previously noted that the sensitivity-error matrix S, installation-error matrix E p , hard magnetic interference H hard and regular variation field H r in the geomagnetic daily variation field can be regarded as fixed values for handling. However, only when the aircraft's attitude and geomagnetic field remain unchanged will the interference-error matrix A soft be fixed. Therefore, the equation does not satisfy the linearity condition.
In designing the filter, choose all elements in the matrix C and vector b's three components as the system's state variables, totaling 12 dimensions: x ¼ ½c xx c xy c xz c yx c yy c yz c zx c zy c zz b As the parameters in the state variables are constant values, the filter system's state equation can be expressed as where w(t) is process noise. Assume that this noise is zero-mean Gaussian noise. Normally there are two options for the selection of the observable, namely the strength and vector of the geomagnetic-field. The algorithm is much simpler when the strength of the geomagnetic-field is used as the observable; therefore, this paper uses the total strength of the geomagnetic-field, including the noise measured using a magnetic sensor as the observable: Further reorganizing Eq (14) gives the following equation.
Eq (18) can be transformed into the following equation.
The above equation can be further specified as follows.
Eqs (16) and (21) constitute the filter mathematical model estimated using the magnetic sensor's error parameters. As the observation equation does not satisfy the linearity condition, it is not possible to directly perform parameter estimation using the Kalman-filter algorithm [21]. This paper used an extended Kalman-filter algorithm to perform parameter estimation. Based on Eqs (16) and (21), the filter model can be expressed as follows.
( _ xðtÞ ¼ wðtÞ where f( ) represents the non-linear function of x(t). The system model is linearized to as follows.
where the state-transition matrix F k,k−1 is a 12-dimensional unit matrix and the observation matrix H k is the Jacobi matrix of partial derivatives of with respect to the state variable x(t).  In summary, the parameter-estimation algorithm of the magnetic sensor's measurementerror model based on an extended Kalman-filter structure is as follows: The model's error parameter can be estimated by the above filtering process, and the compensation of magnetic measurement data can thus be made according to Eq (14).

Experimental verification 4.1. Simulation and analysis
This paper adopts a numerical simulation to verify the validity of the compensation algorithm. First, assuming the strength of the clean geomagnetic-field of the location to be 45306nT, 1,000 groups of data of the geomagnetic-field under ideal conditions were generated; second, the author further generated output data, including interfering items in the light of the conditions set in Table 1, according to an already-built error model as observables for filter estimation. On the basis of the previous two steps, the author performed an estimation of the error model's parameters with the designed filter. After obtaining the estimated values of the parameters, the author compensated the existing parameters and examined the result of the compensation. In the meantime, the author used the RLS (recursive least-square) method and ellipsoid-fitting method to compensate the actual output data and compared the effects of the three compensation methods. The parameters' initial states in the three compensation methods were the same.

Condition settings of simulation.
Settings of the relative values, including the error-coefficient matrix C, are as shown in Table 1.
To test the filter ability, take into account some errors in initial states during the setting of initial state vector, and set the initial state vector as: In general, the measuring noise covariance matrix can be obtained by means of some offline samples, and it represents the measuring noise level of sensor; the process noise covariance matrix is used to describe the uncertainty of model. According to experience, the performance and convergence of estimation algorithm are subject to the impact of parameter values in the covariance matrix; thus it's necessary to conduct some simulations to determine the parameter values. In view of the fact that the values of all principal diagonal elements of error-coefficient matrix C are all around 1 and those of other elements close to 0, and the values of three elements of vector b are relatively large, the initial covariance matrix and the process noise covariance matrix are respectively set as: P 0 ¼ diagð½1; 1 À 4 ; 1 À 4 ; 1 À 4 ; 1; 1 À 4 ; 1 À 4 ; 1 À 4 ; 1; 10 4 ; 10 4 ; 10 4 Þ; Q ¼ diagð½10 À 4 ; 10 À 4 ; 10 À 4 ; 10 À 4 ; 10 À 4 ; 10 À 4 ; 10 À 4 ; 10 À 4 ; 10 À 4 ; 10 4 ; 10 4 ; 10 4 Þ: Since the standard deviation of measuring noise is set as 100, the initial value of observation noise covariance matrix is set as: Comparison of strengths of geomagnetic field before and after compensation. After obtaining the estimation results of filtering, the result of compensating magnetic-measurement data is as shown in Fig 8: the strength of the geomagnetic-field fluctuates greatly before compensation and the fluctuation is obviously reduced after compensation using the extended Kalmanfilter method and the RLS method. However, the former's compensation effect is much better, as the fluctuation of the strength of the geomagnetic-field basically disappears and remains close to the real value after being compensated by the former. Although the fluctuation of the strength of the geomagnetic-field is obviously reduced after being compensated by the latter, its fluctuation margin is larger than that of the former and deviates from the real value. It is known from calculation that the compensation error of the extended Kalman-filter can mostly remain within 5nT, with a mean value of 0.49nT and standard deviation of 16.23nT; however, Error model of geomagnetic-field and compensation method the mean value of the compensation error of RLS is -297.99nT and the standard deviation is 611.90nT.

Results of comparison with ellipsoid fitting algorithm.
The least-square method based on the ellipsoid constraint is used for comparison with the method proposed in the paper. The ellipsoid fitting method is unable to directly work out the parameters of model, so the paper only lists the comparison diagram of geomagnetic field's differential signals after compensations with the two methods as shown in Fig 9. It is known from calculation that the compensation error of the extended Kalman-filter can mostly remain within 5nT, with a mean value of 0.49nT and standard deviation of 16.23nT; however, the mean value of the compensation error of ellipsoid fitting is -0.37nT and the standard deviation is 27.63nT. The method proposed in the paper achieves better compensation effect than the ellipsoid fitting algorithm.

Experiment and analysis
To examine the effect of the method proposed in this paper, an experiment was conducted. First, a clean place with respect to the electromagnetic environment was selected and the strength of the clean geomagnetic field of this place was measured as 48476nT. Second, a three-axis sensor and ferromagnetic materials on a three-axis non-magnetic rotary table were fixed, and the non-magnetic rotary table was randomly rotated in space. The sensor's output data was recorded and 50 groups of data were collected. The strength of the data is as shown in Fig 10, and compensation error is as shown in Fig 11. It is known from calculation that the mean value of the compensation error of the extended Kalman filter is -1.96nT and the standard deviation is 19.45nT; however, the mean value of the compensation error of RLS is -189.39nT and the standard deviation is 245.55nT. The mean value of the compensation error of ellipsoid fitting is -0.37nT and the standard deviation is 27.63nT. It's thus clear that the method proposed in the paper is obviously better than the least-square method with respect to compensation effect and meanwhile avoids the problem of matrix singularity likely to arise in the solving process of the latter. The method proposed in the paper achieves better compensation effect than the ellipsoid fitting algorithm and also avoids the problems, such as heavy computation and the need for good initial conditions, in the solving process of the latter. Error model of geomagnetic-field and compensation method

Conclusion
The measurement of the geomagnetic field using a magnetic sensor is affected by multiple factors. To address these factors, this paper built a brand-new geomagnetic-measurement model that can accurately describe the source of measurement error. It not only considered the magnetic sensor's manufacturing error, installation error and the interference from the external  Error model of geomagnetic-field and compensation method electromagnetic environment but also introduced the previously neglected geomagnetic daily variation field. The author designed an extended Kalman-filter according to the characteristics of the newly built model and performed experiment. The experimental results showed that parameter estimation with this method has good convergence, and that the strength of the geomagnetic-field after being compensated always remains close to the real value. Thus, this method is able to meet the demands of geomagnetic navigation for high accuracy of geomagnetic-field measurement. This compensation method has strong applicability due to its easy data collection and its ability to get rid of the dependence on high-precision measurement instruments.