A high-precision calibration method for residual drift of rotationally modulated inertial navigation

Abstract

This invention belongs to the field of inertial navigation technology and relates to a high-precision calibration method for residual drift in a rotating modulation inertial navigation system, suitable for long-endurance flight environments. This invention establishes a mathematical model of the latitude and longitude navigation error of the rotating inertial navigation system's geographic system, along with the geographic system's azimuth residual drift, northward residual drift, and azimuth mathematical platform error angle. Through the design of long-endurance navigation experiments at eight positions with azimuth angles of 45° intervals from 0 to 360°, accurate identification of azimuth and northward residual drift, as well as alignment heading errors, is achieved at different azimuths. This invention offers high accuracy in residual drift identification, significantly improves long-endurance navigation performance, is highly autonomous, and is easy to implement.

Description

Technical Field

[0001] This invention belongs to the field of inertial navigation technology and relates to a high-precision calibration method for residual drift of rotationally modulated inertial navigation systems. Background Technology

[0002] Rotationally modulated inertial navigation (INS) utilizes the alternating, symmetrical, and orderly rotation of a rotating mechanism to drive the sequential rotation of the IMU components. This enables automatic compensation of the constant errors of the three-axis gyroscope and accelerometer in the geographic frame, significantly improving the long-endurance navigation performance of the INS. However, due to the combined effects of various factors such as gyroscope drift, carrier maneuvering, magnetic fields, and temperature fields, engineering practice has revealed that gyroscope drift cannot be completely modulated, leaving a certain magnitude of equivalent geographic frame residual drift. This residual drift restricts further improvements in the performance of rotationally modulated INS and is one of the key challenges that urgently needs to be addressed and resolved. Summary of the Invention

[0003] The purpose of this invention is to propose a high-precision calibration method for residual drift of rotating inertial navigation systems, which can further improve the long-endurance navigation performance of rotating inertial navigation systems.

[0004] The technical solution of the present invention: In order to achieve the above-mentioned objective, according to the first aspect of the present invention, a high-precision calibration method for residual drift of rotating modulation inertial navigation is proposed. A simplified nine-dimensional Kalman filter model is designed for the residual drift error unique to rotating inertial navigation, and the filtering process is performed at eight different positions. Finally, the filtering results at each position are fitted with Fourier series to obtain the fitting coefficients of the three drifts.

[0005] Specifically, the steps include the following:

[0006] Step 1: After the rotation modulation inertial navigation system completes alignment, it enters navigation mode, collects real-time updated values ​​of longitude and latitude during the navigation process, and records the heading angle ψ at the end of alignment.

[0007] Step 2: Based on the real-time updated longitude and latitude values, a simplified nine-dimensional Kalman filter model is used to estimate the azimuth residual drift, northward residual drift, and azimuth mathematical platform error angle. The inertial navigation attitude angle is set to (0,0,0), and the attitude matrix is ​​always an identity matrix. By removing the three-axis gyroscope drift and three-axis accelerator zero-bias state variables, and adding two-dimensional state variables for geographic system northward drift and geographic system azimuth drift, accurate separation of geographic system drift and alignment heading error is achieved.

[0008] Step 3: Rotate the turntable 45° and perform steps 1 to 2;

[0009] Step 4: Repeat steps 1-3 to make the system perform single-position drift estimation at heading angles ψ+K×45°, K=1,2,...,7 respectively;

[0010] Step 5: Using the 8 sets of azimuth angles and the estimated values ​​of the skyward residual drift and the northward residual drift, Fourier series is used for fitting to obtain the final calibration values ​​of the skyward residual drift and the northward residual drift;

[0011] Step 6: During navigation, the remaining drift calibration values ​​of the azimuth and north directions and the real-time heading angle of the aircraft system are used to calculate the angular velocity compensation and superimpose it onto the original gyroscope data. Then, strapdown calculation is performed to update the velocity, position and attitude to obtain long-endurance high-precision navigation results.

[0012] In one possible embodiment, in step 2, the state variables of the simplified nine-dimensional Kalman filter model are as follows: X K =[δφ,δλ,δV E ,δV N ,φ E ,φ N ,φ U D N D U ] T

[0013] Where δφ is the latitude error, δλ is the longitude error, and δV is the longitude error. E For eastward velocity error, δV N For northbound velocity error, φ E Eastward deflection angle, φ N Northward deflection angle, φ U The celestial deflection angle is D. N Geography Department Northward Residual Drift, D U For the remaining drift of the celestial sphere in the geography system;

[0014] The state matrix F is a 9×9 matrix, where the element F(i,j) in the i-th row and j-th column of matrix F is defined as follows: any undefined element is 0.

[0015]

[0016] F(3,6)=-f U ;

[0017]

[0018]

[0019] Where R is the Earth's radius, h is the local altitude, and ω ie ω is the angular rate of Earth's rotation. n ω u V represents the northward and celestial components of the Earth's rotation angular rate; φ represents the local latitude; E V represents the eastward velocity of the geographic system. N For the northward velocity of the geographic system; f EFor the eastward direction of the geography system; f N For the northward comparison of the geography system; f U For the celestial comparator of the geography system; The element in the first row and first column of the attitude transformation matrix of the machine system relative to the geographic system;

[0020] The measurement matrix H has a dimension of 4×9, and all elements are 0 except for the following:

[0021] H(1,1)=1,H(2,2)=1,H(3,3)=1,H(4,4)=1

[0022] The Kalman filter P-array is initialized as a 9×9 dimensional matrix:

[0023] P=diag(6e-8,6e-8,0.25,0.25,3e-4,3e-4,6e-2,1e-12,1e-12)

[0024] The Kalman filter Q-matrix is ​​initialized as a 25×25 matrix:

[0025] Q=diag(0,0,5e-5,5e-5,5e-13,5e-13,5e-13,0,0);

[0026] Kalman filter R array:

[0027] R=diag(2.47e-2,2.47e-2,1,1).

[0028] In one possible embodiment, in step 2, the attitude matrix is ​​set only for its diagonal elements. All other elements of the attitude matrix are set to By setting the rotating inertial navigation attitude matrix to a constant value, the influence of attitude matrix changes caused by rotation on geographic system drift estimation can be avoided.

[0029] In one possible embodiment, in step 2, the Kalman filter estimate of the astronomical deflection at the first hour is taken. As the result of the celestial drift estimation, the Kalman filter estimates at the 6th hour of the northward drift and celestial drift are taken. As estimates of northward drift and celestial drift Simultaneously, the equivalent eastward drift is calculated using the estimated yaw angle: in, This is an estimate of the eastward drift.

[0030] In one possible embodiment, in step 5, after completing the drift measurement test at 8 different azimuth angles, Fourier series fitting is performed using each set of heading angles and drift estimates according to the following formula;

[0031]

[0032] in, These are the estimated residual drift values ​​for the i-th test location in the east, north, and sky directions, respectively; A E B E These are the first eastward residual drift Fourier fitting coefficients and the second eastward residual drift Fourier fitting coefficients, respectively; A N B N These are the first northward residual drift Fourier fitting coefficients and the second northward residual drift Fourier fitting coefficients, respectively; A U B U These are the Fourier fit coefficients for the first day's residual drift and the Fourier fit coefficients for the second day's residual drift, respectively.

[0033] In one possible embodiment, in step 6, the gyroscope angle increment correction value is calculated using the following formula.

[0034] in, This is the gyroscope angle increment correction value. ψ is the attitude transformation matrix relative to the geographic frame; b DT represents the real-time heading angle of the aircraft system, and DT represents the calculation period.

[0035] In one possible embodiment, in step 6, the gyroscope angle increment output value is:

[0036] in, This is the output after gyroscope angle increment correction. This is the raw output of the gyroscope angle increment.

[0037] In one possible embodiment, the calculation period is 0.0025s.

[0038] Beneficial technical effects of the present invention:

[0039] This invention addresses the problem of residual drift severely limiting the long-endurance navigation performance of rotating inertial navigation systems (INS). Based on the equivalent drift mechanism, a simplified 9-dimensional Kalman state filter is established. Through an eight-position turntable experiment, the residual drift of the INS is calibrated and calculated with high precision, offering advantages such as high accuracy and strong practicality. It effectively solves the problem of residual drift in INS and can be widely applied to various long-endurance INS navigation fields, including airborne, land-based, and vehicle-mounted systems. Attached Figure Description

[0040] Figure 1 A schematic diagram illustrating the use of a single-line lidar to collect environmental information by three robots in a preferred embodiment of the present invention;

[0041] Figure 2 A schematic diagram comparing the northward residual drift of the preferred embodiment of the present invention with that of a comparative example;

[0042] Figure 3 A schematic diagram comparing the celestial residual drift of the preferred embodiment of the present invention with that of a comparative example. Detailed Implementation

[0043] To make the objectives, technical solutions, and advantages of this invention clearer, the embodiments of this invention are described in detail below. It should be noted that, unless otherwise specified, the embodiments and features described in this application can be arbitrarily combined with each other.

[0044] Example 1:

[0045] A high-precision calibration method for residual drift in rotationally modulated inertial navigation systems, such as Figure 1 As shown, it includes the following steps:

[0046] Step 1: The rotating inertial navigation system (INS) is fixedly installed on the single-axis azimuth turntable. After power-on initialization, the rotating frame returns to zero. After alignment, it enters navigation mode. Longitude and latitude update values ​​are collected during navigation, and the heading angle ψ at the end of alignment is recorded.

[0047] Step 2: Calculate the latitude and longitude navigation error using real-time updated latitude and longitude values. Estimate the azimuth residual drift, northward residual drift, and azimuth mathematical platform error angle using the established Kalman filter model. The Kalman filter state variables are reduced to 9 dimensions: X K =[δφ,δλ,δV E ,δV N ,φ E ,φ N ,φ U D N D U ] T

[0048] Where δφ is the latitude error, δλ is the longitude error, and δV is the longitude error. E For eastward velocity error, δV N For northbound velocity error, φ E Eastward deflection angle, φ N Northward deflection angle, φ U The celestial deflection angle is D. N Geography Department Northward Residual Drift, D U This refers to the remaining drift of the celestial sphere in the geography system.

[0049] The state matrix F is a 9×9 matrix, where the element F(i,j) in the i-th row and j-th column of matrix F is defined as follows, and any undefined element is 0.

[0050]

[0051] F(3,6)= - f U ;

[0052]

[0053]

[0054] Where R is the Earth's radius, h is the local altitude, and ω ie ω is the angular rate of Earth's rotation. n ω u V represents the northward and celestial components of the Earth's rotation angular rate. φ represents the local latitude. E V represents the eastward velocity of the geographic system. N For the northward velocity of the geographic system; f E For the eastward direction of the geography system; f N For the northward comparison of the geography system; f U It is a comparison of the celestial directions in the geography system. The element in the first row and first column of the attitude transformation matrix of the machine system relative to the geographic system.

[0055] The measurement matrix H has a dimension of 4×9, and all elements are 0 except for the following:

[0056] H(1,1)=1,H(2,2)=1,H(3,3)=1,H(4,4)=1

[0057] The Kalman filter P-array is initialized as a 9×9 dimensional matrix:

[0058] P=diag(6e-8,6e-8,0.25,0.25,3e-4,3e-4,6e-2,1e-12,1e-12)

[0059] The Kalman filter Q-matrix is ​​initialized as a 25×25 matrix:

[0060] Q=diag(0,0,5e-5,5e-5,5e-13,5e-13,5e-13,0,0);

[0061] Kalman filter R array:

[0062] R = diag(2.47e-2, 2.47e-2, 1, 1)

[0063] Step 3: Rotate the turntable 45° and perform steps 1 to 2.

[0064] Step 4: Repeat steps (1) to (3). Make the system perform single-position drift estimation at each heading angle ψ+K×45° (K=1,2,...,7).

[0065] Step 5: In the single-location experiment, take the Kalman filter estimate of the celestial deflection angle at the first hour. As the result of the celestial drift estimation, the Kalman filter estimates of the northward drift and celestial drift at the 6th hour are taken as the results of the northward drift and celestial drift estimation. Simultaneously, the equivalent eastward drift is calculated using the celestial drift estimates. in, This is an estimate of the eastward drift.

[0066] After completing the drift measurement test at 8 different azimuth angles, Fourier series fitting was performed using the heading angles and drift estimates of each group.

[0067]

[0068] in, These are the estimated residual drift values ​​for the i-th test location in the east, north, and sky directions, respectively. A E B E A N B N A U B U These are the residual drift Fourier fitting coefficients for the east, north, and sky directions, respectively.

[0069] Step 6: During navigation, using the remaining drift calibration values ​​for the azimuth and north directions and the real-time heading angle of the aircraft system, calculate the angular velocity compensation and superimpose it onto the original gyro data. Then, perform strapdown calculations to update the velocity, position, and attitude to obtain high-precision navigation results for long-endurance navigation. Calculate the gyro angle increment correction value for 2.5ms using the fitting coefficients.

[0070]

[0071] in, This is the gyroscope angle increment correction value. Let ψ be the attitude transformation matrix relative to the geographic frame of reference. b The real-time heading angle of the aircraft system is given by DT, which is the calculation period of 0.0025s.

[0072] The gyroscope angle increment output value is:

[0073] in, This is the output after gyroscope angle increment correction. This is the raw output of the gyroscope angle increment.

[0074] Comparative Example

[0075] Step 1: The rotating inertial navigation system (INS) is fixedly installed on the single-axis azimuth turntable. After power-on initialization, the rotating frame returns to zero. After alignment, it enters navigation mode. Longitude and latitude update values ​​are collected during navigation, and the heading angle ψ at the end of alignment is recorded.

[0076] Step 2: Calculate the latitude and longitude navigation error using real-time updated values, and estimate the azimuth residual drift and northward residual drift using a Kalman filter model. The Kalman filter state variables are 25-dimensional.

[0077]

[0078] Note: δφ is latitude error, δλ is longitude error, δh is altitude error, and δV is latitude error. E For eastward velocity error, δV N For northbound velocity error, δV U For the upward velocity error, φ E Eastward deflection angle, φ N Northward deflection angle, φ U The celestial deflection angle is D. X X-axis gyroscope drift, D Y For Y-axis gyroscope drift, D Z For Z-axis gyroscope drift, Add zero position to the X-axis, Add zero position to the Y-axis, Add zero position to Z-axis and D-axis N For the northward drift of the geography system, D U For the celestial drift of the geography system, α XZ Install the deflection angles XZ and α on the gyroscope. XY Install deflection angles XY and α on the gyroscope YX Install deflection angles YX and α on the gyroscope ZX Install the deflection angle ZX, ε on the gyroscope XY To increase the installation deflection angle XY, δK GZ For the Z-axis gyroscope scale coefficient error, δK AX Add scale coefficient error and δK to X. AZ Add a scale coefficient error to the Z-axis.

[0079] The state matrix F is a 25×25 matrix, where the element F(i,j) in the i-th row and j-th column of matrix F is defined as follows, and any undefined element is 0.

[0080]

[0081] F(3,6) = 1;

[0082]

[0083] F(4,9)=f N F(4:6,13:15)=C;

[0084]

[0085] F(5,7)=f U F(5,9)=-f E ;

[0086]

[0087] F(7:9,10:12) = C;

[0088]

[0089]

[0090] Note: The average angular velocity within the filtering period of the IMU along the Y and Z axes; This represents the average specific force within the filtering cycle along the X and Z axes of the IMU. ω ie f is the Earth's rotational angular rate. E For the eastward direction of the geography system; f N For the northward comparison of the geography system; f U It is a comparison of the celestial directions in the geography system.

[0091] The measurement matrix H has a dimension of 3×25, and all elements are 0 except for the following:

[0092] H(1,4)=1,H(2,5)=1,H(3,6)=1

[0093] The Kalman filter P-matrix is ​​initialized as a 25×25 dimensional matrix:

[0094] P=diag(6e-8,6e-8,0,0.25,0.25,0,3e-4,3e-4,6e-2,1e-12,1e-12,1e-12,

[0095] 2.5e-5, 2.5e-5, 2.5e-5, 2.35e-15, 2.35e-15,0,…,0)

[0096] The Kalman filter Q-matrix is ​​initialized as a 25×25 matrix:

[0097] Q=diag(0,0,0,5e-5,5e-5,5e-5,5e-13,5e-13,5e-13,0,...,0,0,0);

[0098] Kalman filter R array:

[0099] R = diag(4e-2, 4e-2, 4e-2)

[0100] Step 3: In the single-location experiment, take the Kalman filter estimate at the 3rd hour of the celestial deflection angle. As the result of the celestial drift estimation, the Kalman filter estimates of the northward drift and celestial drift at the 3rd hour are taken as the results of the northward drift and celestial drift estimation. Simultaneously, the equivalent eastward drift is calculated using the celestial drift estimates. in, This is an estimate of the eastward drift.

[0101] Step 4: During navigation, using the estimated residual drift values ​​for the azimuth and north directions and the real-time heading angle of the aircraft system, calculate the angular velocity compensation and superimpose it onto the raw gyro data. Then, perform strapdown calculations to update the velocity, position, and attitude to obtain high-precision navigation results for long-endurance navigation. Calculate the gyro angle increment correction value for 2.5ms using the fitting coefficients.

[0102]

[0103] in, This is the gyroscope angle increment correction value. Let ψ be the attitude transformation matrix relative to the geographic frame of reference. b The real-time heading angle of the aircraft system is given by DT, which is the calculation period of 0.0025s.

[0104] The gyroscope angle increment output value is:

[0105] in, This is the output after gyroscope angle increment correction. This is the raw output of the gyroscope angle increment.

[0106] like Figure 2 As shown, the embodiment and comparative examples show a -0.001° / h northward residual drift, as Figure 3 As shown in the comparison of the estimation results of the 0° / h celestial residual drift, it can be seen that the model used in the embodiment has a fast convergence speed and high estimation accuracy in estimating the residual drift.

[0107] For illustrative purposes, the specific embodiments described above are merely exemplary and are intended to enable those skilled in the art to better understand this patent. They should not be construed as limiting the scope of this patent. All technical solutions obtained by means of equivalent substitution or equivalent transformation fall within the protection scope of this invention.

Claims

1. A high-precision calibration method for residual drift of a rotationally modulated inertial navigation system, characterized in that, Includes the following steps: Step 1: After the rotation modulation inertial navigation system completes alignment, it enters navigation mode, collects real-time updated values ​​of longitude and latitude during the navigation process, and records the heading angle ψ at the end of alignment. Step 2: Based on the real-time updated longitude and latitude values, use a simplified nine-dimensional Kalman filter model to estimate the directional residual drift, northward residual drift, and directional mathematical platform error angle. Set the inertial navigation attitude angle to (0,0,0), and the attitude matrix is ​​always an identity matrix. Step 3: Rotate the turntable 45° and perform steps 1 to 2; Step 4: Repeat steps 1-3 to make the system perform single-position drift estimation at heading angles ψ+K×45°, K=1,2,...,7 respectively; Step 5: Using the 8 sets of azimuth angles and the estimated values ​​of the skyward residual drift and the northward residual drift, Fourier series fitting is performed to obtain the final calibration values ​​of the skyward residual drift and the northward residual drift; In step 5, after completing the drift measurement test at 8 different azimuth angles, Fourier series fitting is performed using each set of heading angles and drift estimates according to the following formula. in, These are the estimated residual drift values ​​for the i-th test location in the east, north, and sky directions, respectively; A E B E These are the first eastward residual drift Fourier fitting coefficients and the second eastward residual drift Fourier fitting coefficients, respectively; A N B N These are the first northward residual drift Fourier fitting coefficients and the second northward residual drift Fourier fitting coefficients, respectively; A U B U These are the Fourier fit coefficients for the first day's residual drift and the Fourier fit coefficients for the second day's residual drift, respectively. Step 6: During navigation, the remaining drift calibration values ​​of the azimuth and north directions and the real-time heading angle of the aircraft system are used to calculate the angular velocity compensation and superimpose it onto the original gyroscope data. Then, strapdown calculation is performed to update the velocity, position and attitude to obtain long-endurance high-precision navigation results.

2. The high-precision calibration method for residual drift of a rotationally modulated inertial navigation system according to claim 1, characterized in that, In step 2, the state variables of the simplified nine-dimensional Kalman filter model are as follows: X K =[δφ,δλ,δV E ,δV N ,φ E ,φ N ,φ U D N D U ] T Where δφ is the latitude error, δλ is the longitude error, and δV is the longitude error. E For eastward velocity error, δV N For northbound velocity error, φ E Eastward deflection angle, φ N Northward deflection angle, φ U The celestial deflection angle is D. N Geography Department Northward Residual Drift, D U For the remaining drift of the celestial sphere in the geography system; The state matrix F is a 9×9 matrix, where the element F(i,j) in the i-th row and j-th column of matrix F is defined as follows: any undefined element is 0. F(3,6)=-f U ; F(4,5)=f U ; Where R is the Earth's radius, h is the local altitude, and ω ie ω is the angular rate of Earth's rotation. n ω u V represents the northward and celestial components of the Earth's rotation angular rate; φ represents the local latitude; E V represents the eastward velocity of the geographic system. N For the northward velocity of the geographic system; f E For the eastward direction of the geography system; f N For the northward comparison of the geography system; f U For the celestial comparator of the geography system; The element in the first row and first column of the attitude transformation matrix of the machine system relative to the geographic system; The measurement matrix H has a dimension of 4×9, and all elements are 0 except for the following: H(1,1)=1,H(2,2)=1,H(3,3)=1,H(4,4)=1 The Kalman filter P-array is initialized as a 9×9 dimensional matrix: P=diag(6e-8,6e-8,0.25,0.25,3e-4,3e-4,6e-2,1e-12,1e-12) The Kalman filter Q-matrix is ​​initialized as a 25×25 matrix: Q=diag(0,0,5e-5,5e-5,5e-13,5e-13,5e-13,0,0); Kalman filter R array: R=diag(2.47e-2,2.47e-2,1,1).

3. The high-precision calibration method for residual drift of a rotationally modulated inertial navigation system according to claim 2, characterized in that, In step 2, the attitude matrix is ​​set to 1 for only the diagonal elements: All other elements of the attitude matrix are set to 0:

4. The high-precision calibration method for residual drift of a rotationally modulated inertial navigation system according to claim 3, characterized in that, In step 2, the Kalman filter estimate of the celestial deviation angle at the first hour is taken. As the result of the celestial drift estimation, the Kalman filter estimates at the 6th hour of the northward drift and celestial drift are taken. As estimates of northward drift and celestial drift Simultaneously, the equivalent eastward drift is calculated using the estimated yaw angle: in, This is an estimate of the eastward drift.

5. The high-precision calibration method for residual drift of a rotationally modulated inertial navigation system according to claim 1, characterized in that, In step 6, the gyroscope angle increment correction value is calculated using the following formula. in, This is the gyroscope angle increment correction value. ψ is the attitude transformation matrix relative to the geographic frame; b DT represents the real-time heading angle of the aircraft system, and DT represents the calculation period.

6. The high-precision calibration method for residual drift of a rotationally modulated inertial navigation system according to claim 5, characterized in that, In step 6, the gyroscope angle increment output value is: in, This is the output after gyroscope angle increment correction. This is the raw output of the gyroscope angle increment.

7. The high-precision calibration method for residual drift of a rotationally modulated inertial navigation system according to claim 5, characterized in that, The calculation period is 0.0025s.

Images