Next Article in Journal
Theoretical and Experimental: The Synthetic and Anion-Binding Properties of Tripodal Salicylaldehyde Derivatives
Next Article in Special Issue
Optimal Magnetic Sensor Vests for Cardiac Source Imaging
Previous Article in Journal
A Continuous Liquid-Level Sensor for Fuel Tanks Based on Surface Plasmon Resonance
Previous Article in Special Issue
A Magnetoresistive Tactile Sensor for Harsh Environment Applications
 
 
Font Type:
Arial Georgia Verdana
Font Size:
Aa Aa Aa
Line Spacing:
Column Width:
Background:
Article

A Novel Attitude Estimation Algorithm Based on the Non-Orthogonal Magnetic Sensors

Department of Automation, Nanjing University of Science and Technology, No. 200, Xiaolingwei Street, Xuanwu District, Nanjing 210094, China
*
Author to whom correspondence should be addressed.
Sensors 2016, 16(5), 730; https://doi.org/10.3390/s16050730
Submission received: 8 March 2016 / Revised: 13 May 2016 / Accepted: 15 May 2016 / Published: 19 May 2016
(This article belongs to the Special Issue Magnetoresistive Sensors under Extreme Conditions)

Abstract

:
Because the existing extremum ratio method for projectile attitude measurement is vulnerable to random disturbance, a novel integral ratio method is proposed to calculate the projectile attitude. First, the non-orthogonal measurement theory of the magnetic sensors is analyzed. It is found that the projectile rotating velocity is constant in one spinning circle and the attitude error is actually the pitch error. Next, by investigating the model of the extremum ratio method, an integral ratio mathematical model is established to improve the anti-disturbance performance. Finally, by combining the preprocessed magnetic sensor data based on the least-square method and the rotating extremum features in one cycle, the analytical expression of the proposed integral ratio algorithm is derived with respect to the pitch angle. The simulation results show that the proposed integral ratio method gives more accurate attitude calculations than does the extremum ratio method, and that the attitude error variance can decrease by more than 90%. Compared to the extremum ratio method (which collects only a single data point in one rotation cycle), the proposed integral ratio method can utilize all of the data collected in the high spin environment, which is a clearly superior calculation approach, and can be applied to the actual projectile environment disturbance.

1. Introduction

The geomagnetic field is the earth’s natural resource, which can provide a natural coordinate system for aeronautics, astronautics and marine applications. The vehicle’s attitude information can be achieved by using magnetic detection technology to measure the components of the geomagnetic field. A projectile body rotates about its longitudinal axis at very high speed. Due to the limitation of the measuring range of the triaxial gyroscope, it is difficult to use the gyroscope in an actual high-speed, high-spin projectile environment. In contrast, a magnetic sensor has advantages such as a rapid response speed, small size, high resistance to overload, and being cumulatively error free [1,2,3]. A magnetic sensor is applicable to the high-speed, high-spin general projectile attitude measurement [4,5]. Projectile attitude information includes heading angle, pitch angle and rolling angle, and its characteristic is an approximately constant heading angle. Therefore, the projectile attitude information can be estimated using roll angle and pitch angle.
There are several methods for measuring the projectile attitude by unitizing magnetic sensors, including the two-axis non-orthogonal magnetic sensor, three-axis orthogonal magnetic sensor [6,7], and four-axis magnetic sensor. The non-orthogonal magnetic sensors [8] utilize two mounted non-orthogonal magnetic sensors to collect geomagnetic field data during one projectile spin cycle. Because the heading angle and pitch angle are nearly invariable when a projectile rotates one cycle, the roll angle of the projectile can be calculated directly. There are two non-orthogonal magnetic sensor methods for calculating the pitch angle of spinning projectiles: the zero crossing method [9,10] and the extremum ratio method [11]. Both methods utilize the relationship between the eigenvalue ratios of the sensor output and the pitch angle to calculate the projectile attitude. Both of these methods can obtain the estimated value of projectile attitude at a special point. The difference between these methods is as follows: the zero crossing method utilizes phase information from two sensors to calculate the projectile attitude [12], whereas the extremum ratio method utilizes extremum value information from two sensors to calculate the projectile attitude. However, because those two methods use the special point of the sampled magnetic data to calculate the attitude value as the projectile rotates one cycle, the two methods are vulnerable to random disturbance. As the error of sampling value becomes smaller, the error of the estimated projectile attitude becomes smaller, and vice versa [13]. The twin-channel tracking differentiator [14] method utilizes a differential filter to process data from the triaxial magnetic sensors, thereby reducing the random disturbance encountered by the projectile attitude.
To reduce the random disturbance that exists in the extremum ratio method, an integral ratio method is proposed in this paper. This method can reduce the random disturbance by establishing the integral model and can thus achieve greater attitude calculation accuracy than the extremum ratio method.
The remainder of this paper is organized as follows. In Section 2, the projectile attitude measurement principle of magnetic sensor is introduced. In Section 3, the mathematical model of integral ratio method is established. The analytical expression of integral ratio method is derived in Section 4. In Section 5, the simulation results demonstrate that the accuracy of the proposed algorithm is significantly higher than that of the extremum ratio method. Conclusions are drawn in Section 6.

2. Attitude Measurement Principle of the Magnetic Sensors

The geographical coordinate is O x n y n z n , the xn axis points to the magnetic north, the yn axis points to the sky, and the Zn axis points to the magnetic east. As shown in Figure 1, vector Hn represents the magnitude and direction of the magnetic field in the geographical coordinate, Ψ is the heading angle, θ is the pitch angle, and γ is the roll angle.
The intersection angle between vector Hn and its projection in horizontal plane Pxnzn is geomagnetic inclination I. The carrier coordinate Oxcyczc can be achieved by rotating geographical coordinate Oxnynzn I degrees around the zn axis. The xc axis in Oxcyczc coincides with the magnetic field vector Hc, and Hc represents the geomagnetic magnitude and direction of Oxcyczc.
The instantaneous centroid of projectile is chosen as the origin of the projectile coordinate Oxbybzb. The Oxb axis coincides with the lengthwise axis of projectile, pointing to the head, which is positive. Oyb is perpendicular to Oxb, and upward is positive. Ozb is perpendicular to Oxbyb, and the direction is determined by using the right hand rule. Vector Hb represents the magnitude and direction of the magnetic field in the projectile coordinate Oxbybzb.
The coordinates Hbx, Hby, Hbz of vector Hb in the projectile coordinate Oxbybzb can be represented as the transformation from carrier coordinates to projectile coordinates O x c y c z c O x n y n z n O x b y b z b .
H b = L ( γ , ψ , θ ) L ( 0 , 0 , I ) T H c = L ( γ , ψ , θ I ) H c = L ( γ , ψ , θ m ) H c
where θ m is the pitch angle, which includes magnetic dip, θ m = θ I ; L ( γ , ψ , θ ) is the relationship that can be described by the direction cosine matrix.
L ( γ , ψ , θ ) = [ cos θ cos ψ cos ψ sin θ sin ψ sin γ sin ψ cos θ cos γ sin θ sin γ sin θ sin ψ + cos γ cos θ sin γ cos ψ cos γ sin ψ cos θ + sin γ sin θ cos γ sin θ sin ψ sin γ cos θ cos γ cos ψ ]
Vector Hc in carrier coordinates Oxcyczc is
H c = [ h 0 0 ] T
where h is the magnitude of the vector Hc.
The coordinates Hbx, Hby, Hbz of vector Hb in the projectile coordinates Oxbybzb can be represented as [11]
[ H b x H b y H b z ] = [ h cos θ m cos ψ h ( sin γ sin ψ cos θ m cos γ sin θ m ) h ( cos γ sin ψ cos θ m + sin γ sin θ m ) ]
Two non-orthogonal single-axis magnetic sensors, S1 and S2, are assembled separately on A and B along with axis Oxb of the projectile. As shown in Figure 2, the axis Oxb coincides with the longitudinal axis of the projectile. Two sensitive axes are in the plane Oxbyb. The sensitive axis of S1 is assembled along with axis Ozb, and the sensitive axis of S2 is assembled in the angle of λ with axis Oxb.
According to Equation (4), the measured data of magnetic sensor S1 form Hbz, and the measured data of magnetic sensor S2 are the combination of Hbx and Hbz As a result, the measured data Hs1 and Hs2, of two magnetic sensors, can be expressed by θ m , γ and h:
H S 1 = h ( cos γ sin ψ cos θ m + sin γ sin θ m )
H S 2 = h ( cos θ m cos ψ cos λ + cos γ sin ψ cos θ m sin λ + sin γ sin θ m sin λ )
The projectile rotates at nearly constant speed during one spinning circle, and the heading angle Ψ and pitch angle θ are nearly invariable. The roll angles during one spinning cycle can be calculated by obtaining one roll angle at a particular time. Assuming that the particular time is the zero point time of Hs1 and Hs2, the unknown magnetic field intensity scalar h in Equations (5) and (6) can be eliminated when the measured data of Hs1 or Hs2 are zero. The effects of environmental disturbance can thus be reduced.
If the measured data Hs1 = 0, then Equation (5) can be rewritten as
cos γ sin ψ cos θ m + sin γ sin θ m = 0
Furthermore, we can obtain
γ = arctan 2 ( ( 1 ) j + 1 sin ψ cos θ m , ( 1 ) j sin θ m )
where j = 1,2. The two components of Equation (8) cannot be zero simultaneously.
If the measured data Hs2 = 0, then Equation (6) can be rewritten as
cos θ m cos ψ cos λ + cos γ sin ψ cos θ m sin λ + sin γ sin θ m sin λ = 0
If
{ a = cos θ m cos ψ cos λ b = sin ψ cos θ m sin λ c = sin θ m sin λ
then
a + b cos γ + c sin γ = 0
When sin θ m 0 , because sin λ 0 , then c 0 . Thus, Equation (11) can be described as
sin γ = b cos γ a c
Substituting Equation (12) into sin 2 γ + cos 2 γ = 1 gives
( b 2 + c 2 ) cos 2 γ + 2 a b cos γ + a 2 c 2 = 0
When | sin λ | | cos θ m cos ψ | , ( 2 a b ) 2 4 ( b 2 + c 2 ) ( a 2 c 2 ) 0 , Equation (13) has real solution cos γ
cos γ = a 1 b 1 ± c 1 2 ( b 1 2 + c 1 2 a 1 2 ) b 1 2 + c 1 2
where a 1 = cos θ m cos ψ cos λ , b 1 = sin ψ cos θ m sin λ , and c 1 = sin θ m sin λ .
When | sin λ | > | cos θ m cos ψ | , Hs2 has two zero points.
When | sin λ | = | cos θ m cos ψ | , Hs2 has one zero point.
When | sin λ | < | cos θ m cos ψ | , Hs2 does not have a zero point.
From the above, we can see that the value of roll angle γ can be calculated from Equations (8) and (14) when the outputs of two magnetic sensors, S1 and S2, are equal to zero. There are four zero-crossing points that can be used individually to calculate the corresponding roll angle γ. Consequently, the roll angles at any time during one cycle can be calculated. Assuming that the heading angle Ψ does not change [15], when the included angle λ between magnetic sensor S2 and axis Oxb is given, the value of γ is related to the pitch angle θ m . As a result, to further ensure the value of the roll angle γ, it is important to utilize the magnetic data obtained as the projectile rotates one cycle to calculate θ m . Thus, assuming that the rotation speed remains constant, the calculation error of rolling angle is determined by the calculation error of the pitch angle as the projectile rotates one cycle.

3. Mathematical Model

3.1. Mathematical Model of the Extremum Ratio Method

Assume that the heading angle and the pitch angle are invariable in one cycle. The extremum ratio method utilizes the characteristics that there exist one maximum value and one minimum value in every cycle. The pitch angle θ m after the projectile rotates one cycle can be calculated. During one projectile rotation cycle, the roll angle changes between 0 and 2 π . When the output curves Hs1 and Hs2 are at their extremum values, Equations (5) and (6) should meet the following conditions [6]:
d H S 1 d γ = h ( sin γ sin ψ cos θ m + cos γ sin θ m ) = 0
d H S 2 d γ = h ( sin γ sin ψ cos θ m + cos γ sin θ m ) = 0
The mathematic model of the integral ratio method can be described as
g ( θ m ) = H S 2 m H S 1 m = sin λ ± cos θ m cos ψ cos λ sin 2 θ m + sin 2 ψ cos 2 θ m
where Hs1m and Hs2m are the maximum or minimum values achieved by Hs1 and Hs2, respectively. g ( θ m ) represents the ratio of Hs1m and Hs2m. If the heading angle Ψ and the included angle λ between S2 and axis O x b are known, then the estimated pitch angle θ m can be calculated by Equation (17).

3.2. Mathematical Model of the Integral Ratio Method

The mathematical model g ( θ m ) of the extremum ratio method utilizes the ratio between two particular points on magnetic sensor curves to calculate the pitch angle θ m . Its computational accuracy depends on the numerical precision of the sampling values Hs1 and Hs2 at the extreme point. Under practical conditions, due to various types of disturbance, some sampling values may have high error rates. In particular, if impulsive interference exists at any non-extreme point, a rather large sampling value will be generated. As a result, that value will be treated as the largest sampling value to calculate, thereby generating large errors.
A novel integral model is derived in this paper. During the rotation process, all of the samples of Hs1 and Hs2 are used to perform the integral calculation; subsequently, the ratio calculation is performed. The integral model is described as:
f ( θ m ) = 0 2 π H S 2 2 ( γ ) d γ 0 2 π H S 1 2 ( γ ) d γ = f 2 ( θ m ) f 1 ( θ m )
The estimated value f ^ ( θ m ) of f ( θ m ) can be expressed as
f ^ ( θ m ) = 2 π N n = 1 N H ^ S 2 2 ( n ) 2 π N n = 1 N H ^ S 1 2 ( n ) = f ^ 2 ( θ m ) f ^ 1 ( θ m )
where N is the total sampling number during one projectile rotation cycle and n is the sampling point. H ^ S 1 ( n ) = H S 1 ( n ) + v 1 ( n ) , H ^ S 2 ( n ) = H S 2 ( n ) + v 2 ( n ) , H S 1 ( n ) and H S 2 ( n ) are true value of magnetic sensors S1 and S2, respectively. H ^ S 1 ( n ) and H ^ S 2 ( n ) are sampled value of S1 and S2, respectively. v 1 ( n ) and v 2 ( n ) are both Gaussian white noise with a mean of zero and a variance of σ 2 . In Equation (19), one value f ^ ( θ m ) can be calculated in every cycle using the quadratic sum and mean calculation with the sampling values from two magnetic sensors.
The mathematical expectation of H ^ S 1 2 ( n ) can be described as
E [ H ^ S 1 2 ( n ) ] = H S 1 2 ( n ) + E [ v 1 ( n ) 2 ]
where v 1 ~ N ( 0 , σ 2 ) ; according to the properties of normal distribution, we can get
( v 1 ( n ) σ ) 2 ~ χ 2 ( 1 )
Therefore, E [ ( v 1 ( n ) σ ) 2 ] = 1 , and
E [ ( v 1 ( n ) ) 2 ] = σ 2
Substituting Equation (22) into Equation (20) gives
E [ H ^ S 1 2 ( n ) ] = H S 1 2 ( n ) + σ 2
D [ H ^ S 1 2 ( n ) ] = E [ H ^ S 1 4 ( n ) ] E [ H ^ S 1 2 ( n ) ] 2 = 4 H S 1 2 ( n ) σ 2
The mathematical expectation of f ^ 1 ( θ m ) can be given by
E [ f ^ 1 ( θ m ) ] = 2 π N n = 1 N ( H S 1 2 ( n ) + σ 2 ) = 2 π N n = 1 N H S 1 2 ( n ) + 2 π σ 2
Similarly, the mathematical expectation of f ^ 2 ( θ m ) can be given by
E [ f ^ 2 ( θ m ) ] = 2 π N n = 1 N ( H S 2 2 ( n ) + σ 2 ) = 2 π N n = 1 N H S 2 2 ( n ) + 2 π σ 2
Therefore, f ^ 1 ( θ m ) and f ^ 2 ( θ m ) are not unbiased estimators. Compared the mathematical expectation with unbiased estimation of f ^ 1 ( θ m ) and f ^ 2 ( θ m ) , the magnitude of single sampling error variance is amplified 2 π times. So the effect on parameter estimation is small and acceptable. f ^ 1 ( θ m ) and f ^ 2 ( θ m ) can be regarded as powers of two magnetic sensors’ sampled signal, and the power of the magnetic sensor signal is actually sampled as a characteristic value in the integral model to calculate the projectile attitude.
Thus, Equations (5) and (6) can be rewritten as
H S 1 = h ( d 2 cos γ + e 2 sin γ )
H S 2 = h ( a 2 + b 2 cos γ + c 2 sin γ )
where a 2 = cos θ m cos ψ cos λ , b 2 = sin ψ cos θ m sin λ , c 2 = sin θ m sin λ , d 2 = sin ψ cos θ m , and e 2 = sin θ m .
The integral expression of H S 1 can be described as
0 2 π H S 1 2 ( γ ) d γ = 0 2 π ( h ( d 2 cos γ + e 2 sin γ ) ) 2 d γ = h 2 0 2 π ( d 2 2 cos 2 γ + e 2 2 sin 2 γ + 2 d 2 e 2 sin γ cos γ ) d γ
Applying integral calculation on each items in Equation (29) separately gives
0 2 π cos 2 γ d γ = π
0 2 π sin 2 γ d γ = π
0 2 π sin γ cos γ d γ = 0
Substituting Equations (30)–(32) into Equation (29) gives
0 2 π H S 1 2 ( γ ) d γ = h 2 2 ( π d 2 2 + π e 2 2 )
The integral expression of H S 2 can be described as
0 2 π H S 2 2 ( γ ) d γ = 0 2 π ( h ( a 2 + b 2 cos γ + c 2 sin γ ) ) 2 d γ = h 2 0 2 π ( a 2 2 + b 2 2 cos 2 γ + c 2 2 sin 2 γ + 2 b 2 c 2 sin γ cos γ + 2 a 2 b 2 cos γ + 2 a 2 c 2 sin γ ) d γ
Applying the integral calculation on each items in Equation (34) separately gives
0 2 π cos γ d γ = 0
0 2 π sin γ d γ = 0
Substituting Equations (30)–(32), (35) and (36) into Equation (34) gives
0 2 π H S 2 2 ( γ ) d γ = h 2 ( 2 π a 2 2 + π b 2 2 + π c 2 2 )
Using Equations (33) and (37), the integral model f ( θ m ) can be described as
f ( θ m ) = 0 2 π H S 2 2 ( γ ) d γ 0 2 π H S 1 2 ( γ ) d γ = 2 cos 2 θ m cos 2 ψ cos 2 λ + sin 2 ψ cos 2 θ m sin 2 λ + sin 2 θ m sin 2 λ sin 2 ψ cos 2 θ m + sin 2 θ m
Substituting sin 2 θ m = 1 cos 2 θ m into Equation (38) gives
cos 2 θ m ( 2 cos 2 ψ cos 2 λ + sin 2 ψ sin 2 λ sin 2 λ f ( θ m ) sin 2 ψ + f ( θ m ) ) = f ( θ m ) sin 2 λ
We define u = 2 cos 2 ψ cos 2 λ + sin 2 ψ sin 2 λ sin 2 λ f ( θ m ) sin 2 ψ + f ( θ m ) and v = f ( θ m ) sin 2 λ . Then, Equation (39) can be rewritten as
u v   sin 2 θ m = u ( u v ) cos 2 θ m
The pitch angle θ m is therefore given by
θ m = arctan 2 ( ± | u ( u v ) | , ± | u v | )
In Equation (41), the two components cannot be zero simultaneously.
The period of pitch angle θ m is 2 π . The period of function f ( θ m ) with respect to pitch angle θ m is π . As the projectile rotates one cycle, one f ( θ m ) value and the corresponding four θ m values can be obtained. Assuming that h = 1 , Figure 3 gives the comparison between f ( θ m ) and the two integral outputs of f ( θ m ) and f 2 ( θ m ) . It can be seen that when the pitch angle of projectile takes any value, feature values ( f 1 ( θ m ) and f 2 ( θ m ) ) can be got from integral models and the estimated value of θ m can also be obtained throughout missile flight. However, one f ( θ m ) may correspond to two or four possible values of θ m , and the value of the pitch angle θ m can be chosen based on the quadrant in which the pitch angle θ m is located during the actual flight.

4. Novel Integral Ratio Method

A novel integral ratio method is proposed in this paper, and the flow chart of this method is shown in Figure 4. Because this method utilizes all of the data collected during the projectile rotation process, it is convenient to preprocess those data.

4.1. Data Preprocessing

Noise error exists in the measured data from the magnetic sensor. To reduce the noise disturbance, the filtering algorithm is usually adopted to preprocess the data smoothly. A maximum least squares filtering algorithm is used.
(1)
When the samples number is less than the designated value p, the increasing memory filter is successively used to perform the filtering. p is chosen as 10 in this paper.
r ^ ( k ) = 2 ( 1 k ) ( k + 1 ) ( k + 2 ) i = 0 k r ( k i ) + 6 ( k + 1 ) ( k + 2 ) i = 0 k ( k i ) r ( k i )
(2)
When the sample number is greater than p, the fixed memory filter can be used. p is the filter order of moving average part.
r ^ ( k ) = 2 ( 1 p ) ( p + 1 ) ( p + 2 ) i = k p k r ( k i ) + 6 ( p + 1 ) ( p + 2 ) i = k p k ( k i ) r ( k i )
where r ( k ) is the sampling value at time k and r ^ ( k ) is the least squares filtering value at time k .

4.2. The Principle of the Integral Ratio Method

Extreme ratio method only using a single extremum of sample data as characteristic value for processing, which is susceptible to random interference. But the improved integration ratio method utilizes the sum of all data squares of projectile spin as a characteristic value for processing, which can reduce the random noise on the estimated values.
The implementation procedure of the proposed integral ratio method is given as follows:
Step 1: Preprocessing the data from the magnetic sensor using Equations (42) and (43).
The sampling values from the magnetic sensor are filtered by utilizing the maximum least squares filtering algorithm. After filtering the noise error, the value can be input into the integral model to calculate the model function f ^ ( θ m ) .
Step 2: Calculating the integral value of H S 1 and H S 2 .
1)
Determine the integral expression 2 π N n = 1 N H ^ S 1 2 ( n ) .
2)
Determine the integral expression 2 π N n = 1 N H ^ S 2 2 ( n ) .
Step 3: Calculating the integral ratio . f ^ ( θ m ) .
One f ^ ( θ m ) value can be calculated according to Equation (19) by using the quadratic sum and the integral of the samples of two magnetic sensors.
Step 4: Calculating the pitch angle θ m .
During one projectile rotation cycle, assuming that both the included angle λ and the heading angle ψ are invariable, the pitch angle θ m can be determined by Equation (41).
Step 5: Calculating the roll angle γ using Equations (8) and (14).

5. Simulation and Results

5.1. Comparison of the Algorithms

In this section, the noise restraining abilities of the extremum ratio method and the proposed integral ratio method are compared to each other. Assuming that the projectile’s heading angle is ψ = π / 6 , the included angle between magnetic sensor S2 and axis O x b is λ = π / 4 , and the data range of pitch angle is θ m ( 0 , π / 2 ) . The noise of the magnetic sensor is zero mean white noise. Figure 5 and Figure 6 show a comparison of the pitch angle error Δ θ m for different noise intensities ( σ 2 = 0.001 and σ 2 = 0.01 , respectively) after applying the extremum ratio and the proposed integral ratio methods.
Figure 5 and Figure 6 show that both algorithms have rather large errors around θ m = 0 . The reason for the error is that the numerical difference of the magnetic component is too small near the zero point and the error of extremum is too large. The calculation error of the pitch angle also depends on the noise intensity. The error during the variance 0.01 is one order of magnitude greater than the error observed with a variance of 0.001. In the range of the pitch angle, the calculation error of the integral ratio method is much smaller than that of the extremum ratio method.
Table 1 shows the mean and variance of the pitch angle errors, which were calculated by applying the extremum ratio method and the integral ratio method, as well as the error comparison between these methods under different noise intensities. When the noise variance is 0.001, the mean and variance of the errors by applying integral ratio method can be reduced by 84.6% and 89.2%, respectively, compared with that of the errors obtained by applying the extremum ratio method. When the variance of noise is 0.01, the mean and variance of the errors by applying integral ratio method can be reduced by 90.2% and 96.1%, respectively, compared with that of the errors obtained by applying the extremum ratio method.
Table 1 indicates that the mean and variance of the errors obtained by applying the integral ratio method are significantly smaller than those of the errors obtained by applying the extremum ratio method. In addition, as the random errors from the magnetic sensor become larger, the integral ratio method becomes clearly superior in terms of calculation accuracy. Thus, the integral ratio method is clearly superior when the magnetic sensor operates under an environment with a high level of disturbance.
The extremum ratio method utilizes the extremum to calculate the pitch angle as the projectile rotates one cycle. If the extremum encounters any disturbance, the calculation will produce errors. However, the integral ratio method uses the result of the square of the magnetic sensor output to perform the integral calculation. Because the integral calculation performs the function of a filter when the data encounters any white noise, the integral ratio method has strong noise suppression ability.
From the above-mentioned information, we can see that the proposed integral ratio method is clearly superior in both the mean value and the variance of the calculation error compared to the extremum ratio method. When the variance of random noise is 0.01, the error of the integral ratio method is only 10% of the error of the extremum ratio method.

5.2. Ballistic Simulation

The dynamic model of the projectile with magnetic sensor assembled is described as
{ d x d t = v x d y d t = v y d z d t = v z d v x d t = c H τ ( y ) G ( v τ ) v x d v y d t = c H τ ( y ) G ( v τ ) v y g d v z d t = c H τ ( y ) G ( v τ ) v z
where x, y and z are the position components and v x , v y and v z are the velocity components.
The ballistic coefficient c = ( i D 2 / m ) × 10 3 , D is the diameter of projectile, i is the elastic coefficient, and m is the mass of projectile.
H τ ( y ) = H ρ ( y ) τ / τ o n , τ o n is the standard virtual temperature. τ = τ o n 0.0065 h f is the virtual temperature, H ρ ( y ) = ρ / ρ o n is the function of air density, ρ o n = 1.225   kg/m3, ρ = ρ o n exp ( g h f / ( 287.15 τ ) ) , and h f is the flight height.
G ( v τ ) is the resistance function, and v τ = ( v x 2 + v y 2 + v z 2 ) ( τ o n / τ ) .
g = g 0 ( 1 2 h f / ( R 0 + h f ) ) is the gravitational acceleration, g 0 = 9.780   m/s2, and R 0 = 6378160   m.
The equations of the projectile attitude is denoted as
{ ψ = arctan v z v x θ = arctan v y v x γ = ω g J 0.4 L D 3 exp ( 0.4 L D 3 J t )
The initial conditions of the flight trajectory are set as follows: D = 152 mm, L = 1300 mm, m = 52.8 kg, ω g = 40 π   rad/s, and J = 7.5 kg·m2. The muzzle velocity is 550m/s, the initial heading angle is ψ = π / 6 , the pitch angle is θ m = π / 4 , and the rolling angle is γ = 0 . The included angle between the pitch angle of magnetic sensor S2 and axis O x b is λ = π / 4 . The noise of the magnetic sensor is white noise, with a mean of 0 and a variance of 0.01 or 0.001. Both the magnetic sensor sampling period and the updating period of the projectile attitude are 1 ms.
Figure 7 shows the curves of the heading angle, pitch angle and rolling angle. The simulation curve of the rolling angle runs only 1 s, from 50 s to 51 s.
Figure 8 shows the curve of the projectile attitude obtained by utilizing the proposed integral ratio method. The simulation curve of the rolling angle runs only 1 s. The result shows that the proposed integral ratio method can produce an estimated value of the projectile attitude during the whole trajectory. Figure 9 shows the error curve of the projectile attitude by utilizing the integral ratio method.
Figure 10 indicates that the calculation error of the pitch angle obtained by utilizing the integral ratio method is much smaller than that obtained by utilizing the extremum ratio method.
Table 2 lists the variance of the calculation error of the pitch angle and the rolling angle obtained by utilizing two different methods. The data show that the variance of calculation error of pitch angle obtained by utilizing the integral ratio method is one order of magnitude smaller than that obtained by utilizing the extremum ratio method. The variance of the calculation error of the rolling angle obtained by utilizing the integral ratio method is only half that obtained by utilizing the extremum ratio method.

6. Conclusions

To solve the computational problem of a projectile’s pitch angle, the extremum ratio method can be used to select the extreme point to perform the calculation. According to the problem relating to the fact that the extremum ratio method will result in a large error under a noise disturbance, a novel integral ratio method was proposed, and the corresponding mathematical model was derived. In this paper, the computational expression of pitch angle was derived by utilizing the integral ratio method. The simulation results provided comparison diagrams between the extreme ratio method and the integral ratio method under different noise intensities. Compared to the extremum ratio method, the integral ratio method has several advantages. When the variance of noise becomes larger, the integral ratio method is clearly superior in terms of the calculation error. In addition, the extremum ratio method only selects one data point in every spin cycle, whereas the integral ratio method selects multi-group data in every spin cycle; as a result, multiple data processing algorithms can be performed to get more accurate attitude estimation. Furthermore, the integral model allows for digital processing. An existing digital device can be utilized to perform high-speed, high-precision AD (analog digital conversion) sampling.

Acknowledgments

This work is partially supported by National Nature Science Foundation of China (61473153), Project Funded by the Priority Academic Program Development of Jiangsu Higher Education Institutions (PAPD), Six Talent Peaks Project in Jiangsu Province (2015-XXRJ-006) and Natural Science Foundation of Jiangsu Province (BK20131352).

Author Contributions

All authors contributed extensively to the study presented in this manuscript. Jianliang Zhu and Panlong Wu designed the software, performed the experiment and analyzed the measurements. Yuming Bo supervised the work and contributed with valuable discussions and scientific advice. All authors contributed in writing this manuscript.

Conflicts of Interest

The authors declare no conflict of interest.

References

  1. Park, S.G.; Jeong, H.C.; Kim, J.W.; Hwang, D.; Lee, S.J. Magnetic compass fault detection method for GPS/INS/magnetic compass integrated navigation systems. Int. J. Control Autom. Syst. 2011, 9, 276–284. [Google Scholar] [CrossRef]
  2. Rogers, J.; Costello, M.; Harkins, T.; Hamaoui, M. Effective Use of Magnetometer Feedback for Smart Projectile Applications. Navigation 2011, 58, 203–219. [Google Scholar] [CrossRef]
  3. Zhao, W.; Bu, X.; Yu, G.; Xiang, C. Feedback-type giant magneto-impedance sensor based on longitudinal excitation. J. Magn. Magn. Mater. 2012, 324, 3073–3077. [Google Scholar] [CrossRef]
  4. Včelák, J.; Kubík, J. Influence of sensor imperfections to electronic compass attitude accuracy. Sensors Actuators A Phys. 2009, 15, 233–240. [Google Scholar] [CrossRef]
  5. Long, D.F.; Lin, J.; Zhang, X.M.; Li, J. Orientation estimation algorithm applied to high-spin projectiles. Meas. Sci. Technol. 2014, 25, 233–243. [Google Scholar] [CrossRef]
  6. Xiang, C.; Bu, X.; Zhu, Y. The Design of New Attitude Measuring Method of Non-Spinning Missile Based on Magnetic Sensors. Appl. Mech. Mater. 2012, 226, 1825–1828. [Google Scholar] [CrossRef]
  7. Xiang, C.; Bu, X.; Yang, B. Three Different Attitude Measurements of Spinning Projectile Based on Magnetic Sensors. Measurement 2014, 47, 331–340. [Google Scholar] [CrossRef]
  8. Hepner, D.; Harkins, T. Determining inertial orientation of a spinning body with body-fixed sensors. In Proceedings of the SPIE in Acquisition, Tracking, and Pointing XIV, Orlando, FL, USA, 26–27 April 2000; pp. 68–78.
  9. Harkins, T.E.; Hepner, D.J. MAGSONDE: A Device for Making Angular Measurements on Spinning Projectiles Using Magnetic Sensors. Available online: www.dtic.mil/cgi-bin/GetTRDoc?Location=U2&doc=GetTRDoc.pdf&AD=ADA386941 (accessed on 17 May 2016).
  10. Harkins, T.E.; Davis, B.S.; Hepner, D.J. Novel on-board sensor systems for making angular measurements on spinning projectiles. Proc. SPIE 2001, 4365, 176–187. [Google Scholar]
  11. Li, D.; Bu, X. Attitude measurement on high-spinning projectile using magnetic sensors and accelerometers. Trans. Nanjing Univ. Aeronaut. Astronaut. 2008, 25, 106–112. [Google Scholar]
  12. Liu, X. An Improved Interpolation Method for Wind Power Curves. IEEE Trans. Sustain. Energy 2012, 3, 528–534. [Google Scholar]
  13. Li, X.; Li, Z. A new calibration method for tri-axial field sensors in strap-down navigation systems. Meas. Sci. Technol. 2012, 23, 2852–2855. [Google Scholar] [CrossRef]
  14. Huang, Y.; Wang, H. Application of tracking-differentiator in angular measurements on spinning projectiles using magnetic sensors. In Proceedings of the IHMSC 2015, Hangzhou, China, 26–27 August 2015; Volume 25, pp. 433–436.
  15. Changey, S.; Fleck, V.; Beauvois, D. Projectile attitude and position determination using magnetic meter sensor only. Proc. SPIE 2005, 5803, 49–58. [Google Scholar] [CrossRef]
Figure 1. Coordinate transformation schematic diagram.
Figure 1. Coordinate transformation schematic diagram.
Sensors 16 00730 g001
Figure 2. Schematic diagram of mounted magnetic sensors.
Figure 2. Schematic diagram of mounted magnetic sensors.
Sensors 16 00730 g002
Figure 3. Comparison between the period of function ( f ( θ m ) ) and the two integral outputs of feature values f 1 ( θ m ) and f 2 ( θ m ) .
Figure 3. Comparison between the period of function ( f ( θ m ) ) and the two integral outputs of feature values f 1 ( θ m ) and f 2 ( θ m ) .
Sensors 16 00730 g003
Figure 4. The flow chart of integral ratio method.
Figure 4. The flow chart of integral ratio method.
Sensors 16 00730 g004
Figure 5. Comparison of the pitch angle error with σ 2 = 0.001 .
Figure 5. Comparison of the pitch angle error with σ 2 = 0.001 .
Sensors 16 00730 g005
Figure 6. Comparison of the pitch angle error with σ 2 = 0.001 .
Figure 6. Comparison of the pitch angle error with σ 2 = 0.001 .
Sensors 16 00730 g006
Figure 7. The curves of the attitude angle during projectile flight.
Figure 7. The curves of the attitude angle during projectile flight.
Sensors 16 00730 g007
Figure 8. The attitude angle calculated by the integral ratio method.
Figure 8. The attitude angle calculated by the integral ratio method.
Sensors 16 00730 g008
Figure 9. Attitude angle calculating error of the integral ratio method.
Figure 9. Attitude angle calculating error of the integral ratio method.
Sensors 16 00730 g009
Figure 10. The pitch angle calculated by the extremum ratio method and the integral ratio method.
Figure 10. The pitch angle calculated by the extremum ratio method and the integral ratio method.
Sensors 16 00730 g010
Table 1. The comparison of the pitch angle error.
Table 1. The comparison of the pitch angle error.
MethodError MeanError Variance
(unit: 0.001 rad)(unit: 0.001)
σ2 = 0.001σ2 = 0.01σ2 = 0.001σ2 = 0.01
Extremum ratio method0.4538911.20.0177220.87621
Integral ratio method0.0700041.10.00191820.034168
Error reduction84.6%90.2%89.2%96.1%
Table 2. The error variance comparison between the extremum ratio method and the proposed integral ratio method.
Table 2. The error variance comparison between the extremum ratio method and the proposed integral ratio method.
MethodPitch AngleRoll Angle
(unit: 0.001)(unit: 0.001)
Extremum ratio method1.45260.20321
Integral ratio method0.186520.093807
Error reduction87.2%53.8%

Share and Cite

MDPI and ACS Style

Zhu, J.; Wu, P.; Bo, Y. A Novel Attitude Estimation Algorithm Based on the Non-Orthogonal Magnetic Sensors. Sensors 2016, 16, 730. https://doi.org/10.3390/s16050730

AMA Style

Zhu J, Wu P, Bo Y. A Novel Attitude Estimation Algorithm Based on the Non-Orthogonal Magnetic Sensors. Sensors. 2016; 16(5):730. https://doi.org/10.3390/s16050730

Chicago/Turabian Style

Zhu, Jianliang, Panlong Wu, and Yuming Bo. 2016. "A Novel Attitude Estimation Algorithm Based on the Non-Orthogonal Magnetic Sensors" Sensors 16, no. 5: 730. https://doi.org/10.3390/s16050730

APA Style

Zhu, J., Wu, P., & Bo, Y. (2016). A Novel Attitude Estimation Algorithm Based on the Non-Orthogonal Magnetic Sensors. Sensors, 16(5), 730. https://doi.org/10.3390/s16050730

Note that from the first issue of 2016, this journal uses article numbers instead of page numbers. See further details here.

Article Metrics

Back to TopTop