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
, 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
.
where
is the pitch angle, which includes magnetic dip,
;
is the relationship that can be described by the direction cosine matrix.
Vector
Hc in carrier coordinates
Oxcyczc is
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]
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
,
and
h:
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
Furthermore, we can obtain
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
When
, because
, then
. Thus, Equation (11) can be described as
Substituting Equation (12) into
gives
When
,
, Equation (13) has real solution
where
,
, and
.
When , Hs2 has two zero points.
When , Hs2 has one zero point.
When , 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
. 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
. 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.
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.
- (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.
where
is the sampling value at time
and
is the least squares filtering value at time
.
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 .
Step 2: Calculating the integral value of
and
.
- 1)
Determine the integral expression .
- 2)
Determine the integral expression .
Step 3: Calculating the integral ratio ..
One 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 .
During one projectile rotation cycle, assuming that both the included angle and the heading angle are invariable, the pitch angle 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
, the included angle between magnetic sensor S2 and axis
is
, and the data range of pitch angle is
. The noise of the magnetic sensor is zero mean white noise.
Figure 5 and
Figure 6 show a comparison of the pitch angle error
for different noise intensities (
and
, 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
. 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
where
x,
y and
z are the position components and
,
and
are the velocity components.
The ballistic coefficient , D is the diameter of projectile, i is the elastic coefficient, and is the mass of projectile.
, is the standard virtual temperature. is the virtual temperature, is the function of air density, kg/m3, , and is the flight height.
is the resistance function, and .
is the gravitational acceleration, m/s2, and m.
The equations of the projectile attitude is denoted as
The initial conditions of the flight trajectory are set as follows: D = 152 mm, L = 1300 mm, m = 52.8 kg, rad/s, and J = 7.5 kg·m2. The muzzle velocity is 550m/s, the initial heading angle is , the pitch angle is , and the rolling angle is . The included angle between the pitch angle of magnetic sensor S2 and axis is . 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.