MEMS navigation attitude correction method and device based on motion acceleration estimation

Abstract

The invention provides an MEMS navigation attitude correction method and device based on motion acceleration estimation, and belongs to the field of MEMS navigation algorithms. The method comprises the following steps: carrying out filtering and moving average processing on real-time angular velocity and acceleration data, adding Markov modeling on a motion acceleration state on the basis of a current navigation system model, and carrying out one-step prediction through state transition. And when the heading attitude filtering admission condition is met, separating the motion acceleration in the measurement matrix, and estimating the current navigation state of the carrier by using Kalman filtering. And finally, east and north mathematical platform deflection angles with motion acceleration influence deducted are used during filtering correction. According to the method, accurate modeling can be carried out on the motion acceleration under the large maneuvering condition, the attitude error correction probability of the MEMS attitude and heading system is reduced, the convergence time after a carrier is converted into stable flight from maneuvering is shortened, and the product performance and the use experience are improved.

Description

Technical Field

[0001] This invention relates to the field of MEMS navigation algorithms, and in particular to a MEMS attitude correction method and apparatus based on motion acceleration estimation. Background Technology

[0002] The heading and attitude system is a flight instrument that measures, displays, and provides aircraft heading and attitude angle signals. This system can instruct the pilot on the heading angle required for navigation and the roll and pitch angles required for piloting. Strapdown inertial navigation (SINS) is an inertial navigation system where inertial sensors (mainly gyroscopes and accelerometers) are directly fixed to the vehicle. Strapdown attitude and heading systems are based on the principles of strapdown inertial navigation, using gyroscopes and accelerometers to measure the angular velocity and acceleration of the vehicle to obtain its attitude and heading information. This has become a very important direction in the development of inertial technology.

[0003] In a strapdown attitude and heading system, the main task of the MEMS attitude and heading system is to collect signals from various sensors, primarily gyroscopes and accelerometers. It then performs navigation calculations on the collected angular velocity and acceleration information and finally sends the calculated attitude information to a display terminal. In a strapdown inertial attitude and heading system (attitude and heading system), MEMS gyroscopes are used to sense the angular velocity of the carrier. After signal preprocessing, this provides the raw information for carrier attitude calculation. However, due to the low accuracy and large drift of gyroscopes, the required attitude accuracy cannot be guaranteed over long periods of operation. By using accelerometers to measure the components of gravity along the three axes of the carrier, the tilt angle of the carrier relative to the local horizontal plane can be calculated.

[0004] However, when the attitude control system's carrier accelerates or decelerates, the specific force sensed by the accelerometer includes the acceleration along the corresponding axis. Directly using the accelerometer's sensed value as the filter's measurement would lead to extremely serious attitude miscalibration. To avoid this, the attitude filtering thresholds generally include a rough determination of whether the carrier is in motion, depending on the carrier's motion characteristics. However, even small accelerations or decelerations within the aforementioned thresholds can cause significant attitude miscalibration. Furthermore, once miscalibrated, the time required to return to a normal state increases, significantly reducing system performance and user experience. Summary of the Invention

[0005] The purpose of this invention is to design a MEMS attitude correction method and device based on motion acceleration estimation. By performing Markov modeling on motion acceleration, motion acceleration is separated in the measurement matrix to achieve real-time estimation of the carrier's motion acceleration, and finally obtains the angle correction value that eliminates motion acceleration interference.

[0006] The first aspect of this invention provides a MEMS attitude correction method based on motion acceleration estimation, comprising the following steps:

[0007] S1. Acquire angular velocity and acceleration according to a preset period, and perform filtering and moving average processing on the acquired angular velocity and acceleration to obtain filtered and smoothed angular velocity and acceleration.

[0008] S2, based on the state and variance of the attitude system at the previous moment, as well as the preset noise variance matrix and one-step transition matrix, obtain the predicted state and variance of the attitude system at the current moment; the one-step transition matrix is ​​obtained based on first-order Markov modeling;

[0009] S3: Determine the preset attitude filtering state based on the filtered and smoothed angular velocity and acceleration; establish the Kalman filter gain matrix based on the predicted state of the attitude system at the current moment, the measurement noise matrix, and the observation matrix; obtain the attitude filtering measurement matrix based on the filtered acceleration; obtain the state of the attitude system at the current moment based on the attitude filtering measurement matrix, the observation matrix, the Kalman filter gain matrix, the predicted state and variance of the attitude system at the current moment; the observation matrix is ​​a matrix that separates motion acceleration.

[0010] S4 corrects the attitude of the attitude system based on the current state of the attitude system.

[0011] Optionally, the method further includes:

[0012] Based on the attitude and attitude filter measurement matrix, observation matrix, Kalman filter gain matrix, and the predicted state and variance of the attitude and attitude system at the current moment, the variance of the attitude and attitude system at the current moment is obtained.

[0013] Optionally, the one-step transition matrix can be represented by the following formula:

[0014] Φ K / K-1 =I+F·T F +(F·T F ) 2 / 2;

[0015] Where, Φ K / K-1 This represents the state transition from the previous time k-1 to the current time k. I is the unit matrix, and F is the continuous state transition matrix. The elements on the diagonal of the matrix take the value (-α-1) / T. F T F The state transition period is α = 1 / τ, where τ is the correlation time of the Markov process.

[0016] Optionally, a preset attitude filtering state is determined based on the filtered and smoothed angular velocity and acceleration, including:

[0017] Based on the filtered and smoothed acceleration, obtain the composite acceleration, and then obtain the absolute value of the difference between the composite acceleration and the gravitational acceleration.

[0018] Determine whether the absolute value of the difference is less than a given force threshold and whether the filtered and smoothed angular velocity is less than a given angular rate threshold. If so, the attitude filtering state is satisfied.

[0019] Optionally, in the observation matrix H[f E ,mE]=H[f N ,mN]=H[f U [mU] = 1;

[0020] Where H is the observation matrix, mE, mN, and mU represent the east, north, and celestial accelerations, respectively; fE, fN, and fU represent the east, north, and celestial accelerations, respectively.

[0021] Optionally, the attitude of the attitude system can be corrected based on its current state, including:

[0022] The eastward and northward platform yaw angles in the current attitude system state are subjected to amplitude limiting processing to obtain the geographic system correction quantity φ. E ,φ N ;

[0023] According to the geographical system, the quantity to be corrected is φ. E ,φ N Obtain the attitude angles to be corrected, Δθ and Δγ;

[0024] The attitude of the attitude system is corrected based on the attitude angles to be corrected, Δθ and Δγ.

[0025] Optionally, the values ​​of certain elements on the diagonal of the measurement noise matrix, system noise matrix, and system error variance matrix are determined as follows;

[0026] The diagonal elements of the measurement noise array R are:

[0027] R = diag(R) 11 ,R 22 ,R 33 )

[0028]

[0029] in, The equivalent zero-point mean value is determined by the selected device, where g is the gravitational acceleration.

[0030] The partial diagonal elements of the systematic error variance matrix P take the following values:

[0031]

[0032] Where P[x,y] represents the element in the x-row, y-column of matrix P. The standard deviation of acceleration is estimated from the system's trajectory.

[0033] The partial diagonal elements of the system noise variance matrix Q are:

[0034]

[0035] A second aspect of the present invention provides a MEMS attitude correction device based on motion acceleration estimation for performing the method as described in any one of the first aspects.

[0036] The advantages and effects of this invention can be:

[0037] (1) The motion acceleration prediction model based on the Markov process is simple to model, has strong robustness, and is suitable for time-varying parameter prediction.

[0038] (2) Kalman filtering is used to calculate the attitude after deducting motion acceleration. The algorithm is relatively simple and the estimation effect is good. It is very suitable for the scenario of attitude estimation under large vehicle maneuver.

[0039] (3) The algorithm is universal and suitable for attitude estimation of attitude systems with different working modes. Attached Figure Description

[0040] To more clearly illustrate the technical solutions implemented in this invention, a simple explanation of the accompanying drawings used in the description of this invention will be provided below. Obviously, the drawings described below are merely some embodiments of this invention. For those skilled in the art, other drawings can be obtained based on these drawings without any creative effort.

[0041] Figure 1 This is a schematic diagram illustrating the principle of a MEMS attitude algorithm based on motion acceleration estimation according to an embodiment of the present invention.

[0042] Figure 2 This is a schematic diagram of the MEMS attitude algorithm based on motion acceleration estimation according to an embodiment of the present invention. Detailed Implementation

[0043] To make the objectives, technical solutions, and advantages of the embodiments of the present invention clearer, the technical solutions of the embodiments of the present invention will be clearly and completely described below with reference to the accompanying drawings. Obviously, the described embodiments are only some, not all, of the embodiments of the present invention. Based on the embodiments of the present invention, all other embodiments obtained by those skilled in the art without creative effort are within the scope of protection of the present invention.

[0044] The features of various aspects of the embodiments of the present invention will now be described in detail. In the following detailed description, numerous specific details are set forth to provide a thorough understanding of the invention. However, it will be apparent to those skilled in the art that the invention can also be practiced without these specific details. The following description of the embodiments is merely intended to provide a better understanding of the invention by illustrating examples. The invention is not limited to any specific setups and methods provided below, but covers all improvements, substitutions, etc., to product structures and methods without departing from the spirit of the invention. In the various drawings and the following description, well-known structures and techniques are not shown to avoid unnecessarily obscuring the invention.

[0045] It should be noted that, unless otherwise specified, the embodiments of the present invention and the features thereof can be combined with each other, and the various embodiments can be referenced and cited in each other. The present invention will now be described in detail with reference to the accompanying drawings and embodiments.

[0046] Figure 1 This is a schematic diagram illustrating the principle of a MEMS attitude algorithm based on motion acceleration estimation according to an embodiment of the present invention.

[0047] like Figure 1 As shown, the method may include the following steps:

[0048] Step 1: Filter and move average the real-time angular velocity and acceleration data;

[0049] Step 2: Perform Markov modeling on the motion acceleration state and make a one-step prediction through state transition;

[0050] Step 3: Separate the motion acceleration from the measurement matrix and use Kalman filtering to estimate the current navigation state of the vehicle;

[0051] Step 4: Correct the attitude using the platform deflection angle after deducting motion acceleration.

[0052] Optionally, S1 performs filtering and moving average processing on the real-time angular velocity and acceleration data, including:

[0053] Collect real-time angular velocity and acceleration data ω i ,f i i = x, y, z;

[0054] The collected angular velocity and acceleration data ω i ,f i Filtering and moving average preprocessing are performed.

[0055] Optionally, the collected angular velocity and acceleration data ω i ,f iFiltering and moving average preprocessing are performed, including:

[0056] S11, for the collected angular velocity and acceleration data ω i ,f i Perform low-pass filtering and calculate Lω i (t),Lf i (t);

[0057] S12, with Lω i (t),Lf i Using (t) as the time base, calculate the 5-step moving average to obtain...

[0058] S13, Calculate the composite acceleration using the processed acceleration.

[0059]

[0060] And calculate |g′-g|, where g is the given gravitational acceleration;

[0061] Optionally, S2, performs Markov modeling of the motion acceleration state and makes a one-step prediction through state transition, including:

[0062] S21, perform a first-order Markov model of the motion acceleration, i.e.

[0063] X K [mj]=-αX K-1 [mj]+ε k j = E, N, U;

[0064] Where α = 1 / τ, τ is the correlation time of the Markov process, and ε is the white noise excitation of the motion acceleration;

[0065] S22, let F[mi,mi]=(-α-1) / T F T F This is the state transition period;

[0066] Calculate the continuous state transition matrix (F matrix) using the current system navigation parameters;

[0067] S23, Calculate the one-step transition matrix

[0068] Φ K / K-1 =I+F·T F +(F·T F ) 2 / 2

[0069] S24, let Q k [mi,mi]=1+E[ε k The noise variance matrix (Q) of the discrete system is calculated using set values.k (array);

[0070] S25, perform a one-step prediction of the state.

[0071] X K / K-1 =Φ K / K-1 X K-1 ;

[0072] S26, Calculate the mean square error of the one-step prediction.

[0073]

[0074] Optionally, S3, the motion acceleration is separated in the measurement matrix, and the current navigation state of the vehicle is estimated using Kalman filtering, including:

[0075] S31, if |g′-g| is less than the given specific force threshold, and If the angular rate is less than the given threshold, the attitude filtering state is satisfied, and calculations S33-S36 and S41-S43 are performed; otherwise, the attitude filtering state is not satisfied, and the process is terminated.

[0076] S32, Calculate attitude filter measurement The coordinate transformation matrix from the machine system to the navigation system is derived from the inertial calculation.

[0077] S33, separate the motion acceleration in the measurement matrix, that is

[0078] H[f E ,mE]=H[f N ,mN]=H[f U [mU] = 1;

[0079] S34, Calculate the measurement noise array R

[0080] S35, Calculate the Kalman filter gain matrix K

[0081] K = P K / K-1 ·H T ·(H·P K / K-1 ·H T +R) -1

[0082] S36, Estimating the new state X K

[0083] X K =X K / K-1 +K·(ZH·X K / K-1 )

[0084] S37, Update the error variance matrix P K

[0085] P K =(IK·H)·P K / K-1

[0086] Optionally, S4, attitude correction is performed using a platform yaw angle that has been subtracted from motion acceleration, including:

[0087] S41, calculated X[φ E ],X[φ N After being limited in width, the value φ is used as the correction value for the geography system. E ,φ N ;

[0088] S42, using φ E ,φ N Calculate the attitude angles to be corrected, Δθ and Δγ;

[0089] Δθ=-φ E cosψ+φ N sinψ

[0090] Δγ=-φ E sinψ / cosθ-φ N cosψ / cosθ

[0091] Where θ and ψ are the pitch and yaw angles of the current attitude system, calculated by inertial analysis. Δθ and Δγ are the values ​​to be corrected for the pitch and roll angles.

[0092] S43, use Δθ, Δγ to correct the system attitude angle, i.e. θ = θ - Δθ, γ = γ - Δγ.

[0093] like Figure 2 As shown, the method may include the following steps:

[0094] Step 1: Execute 101 to begin and set the initial state of the parameters;

[0095] Step 2: Execute 102 to perform low-pass filtering and moving average processing on the collected acceleration and angular velocity;

[0096] Step 3: Execute 103, and calculate the composite acceleration when the filtering period ends.

[0097] Otherwise, execute 125;

[0098] Step 4: Execute 104, use the pre-processed angular velocity and acceleration to determine whether absolute stillness is satisfied. If it is satisfied, it is determined to be absolutely still; otherwise, it is not satisfied.

[0099] Step 5: Execute 105. If the system changes from absolute rest to non-absolute rest, then set P. K[mi,mi] are the initial values;

[0100] Step 6: Execute 106, calculate F[mi,mi]=(-α-1) / T F i = E, N, U;

[0101] Step 7: Execute 107 to calculate other parameters of matrix F;

[0102] Step 8: Execute step 108 to calculate the one-step transition matrix.

[0103] Φ K / K-1 =I+F·T F +(F·T F ) 2 / 2;

[0104] Step 9: Execute 109 to calculate the discrete system noise variance matrix (Q) using the initial settings. k (array), among which;

[0105] Step 10: Execute 110 to make a one-step prediction of the state.

[0106] X K / K-1 =Φ K / K-1 X K-1 ;

[0107] Step 11: Execute step 111 to calculate the one-step prediction mean square error.

[0108]

[0109] Step 12: Execute 112, use the preprocessed angular velocity and acceleration to determine whether the attitude filtering state is met. If it is met, it is determined that the attitude filtering condition is met, and execute S113-S124. Otherwise, the attitude filtering condition is not met, and execute S125.

[0110] Step 13: Execute 113 to calculate attitude filtering measurements. The coordinate transformation matrix from the machine system to the navigation system is derived from the inertial calculation.

[0111] Step 14: Execute 114 to separate the motion acceleration in the measurement matrix, i.e.

[0112] H[f E ,mE]=H[f N ,mN]=H[f U [mU] = 1

[0113] Step 15: Execute step 115 to calculate the other values ​​in the H matrix, i.e., H[f E ,φ N ]=-g,H[fN ,φ E ] = g;

[0114] Step 16: Execute 116 to calculate the measurement noise array R.

[0115] Step 17: Execute step 117 to calculate the Kalman filter gain matrix.

[0116] K = P K / K-1 ·H T ·(H·P K / K-1 ·H T +R) -1

[0117] Step 18: Execute 118 to estimate the new system state.

[0118] X K =X K / K-1 +K·(ZH·X K / K-1 )

[0119] Step 19: Execute step 119 to update the error variance matrix.

[0120] P K =(IK·H)·P K / K-1

[0121] Step 20: Execute step 120 to perform symmetry processing on matrix P.

[0122]

[0123] Step 21: Execute 121 to save the previous data P. K-1 =P K ,X K-1 =X K ;

[0124] Step 22: Execute 122 to calculate X[φ] E ],X[φ N After being limited in width, the value φ is used as the correction value for the geography system. E ,φ N ;

[0125] Step 23: Execute 123, using φ E ,φ N Calculate the attitude angles to be corrected, Δθ and Δγ.

[0126] Δθ=-φ E cosψ+φ N sinψ

[0127] Δγ=-φ E sinψ / cosθ-φ N cosψ / cosθ

[0128] Where θ and ψ are the pitch and yaw angles of the current attitude system, calculated by inertial analysis. Δθ and Δγ are the values ​​to be corrected for the pitch and roll angles.

[0129] Step 24: Execute 124, using Δθ and Δγ to correct the system attitude angles, i.e.

[0130] θ = θ - Δθ, γ = γ - Δγ;

[0131] Step 25: Execute 125 to complete;

[0132] In some embodiments, the MEMS attitude estimation algorithm based on motion acceleration is described in detail below.

[0133] First, initialize the filter parameters after the system is powered on.

[0134] The initial value of the X matrix is ​​0, and its dimension is 15×1.

[0135]

[0136] The initial value of matrix P is:

[0137]

[0138] in:

[0139]

[0140] P δVE =P δVN = (10.0 m / s) 2 ;

[0141] P φE =P φN = (0.5 / 57.3) 2 ;

[0142]

[0143] P εx =P εy =P εz = (30 × RADS) 2 ;

[0144] P ax =P ay =P az = (0.001g) 2 ;

[0145] P δH =1e4;

[0146] P δVU=100

[0147] P fmE =P fmN =P fmU =(0.5) 2

[0148] The dimension of matrix P(0) is 15×15, and the remainder of matrix P(0) is 0. The initial value of matrix Q is:

[0149]

[0150] σ 2 amE =σ 2 amN =σ 2 amU =(0.05) 2

[0151] The initial value of array R is...

[0152] R(0) 3×3 =diag(R0,R0,R0)

[0153] R0 = 0.0292

[0154] After the interrupt is enabled, the real-time acquired angular velocity and acceleration are preprocessed. In this invention, data preprocessing involves performing moving average and low-pass filtering on the angular velocity and acceleration.

[0155] Calculate the angular velocity and acceleration using the moving average method:

[0156] Taking the X-axis as an example, the instantaneous angular velocity of the machine system relative to the geographic system calculated every Ts is ω. x Calculate the current (5×Ts) moving average angular velocity.

[0157] m1=t / T SA

[0158]

[0159] Where: t is the moving average time period, which is (5×Ts) here;

[0160] T SA The period for calculating the moving average is Ts in this case.

[0161] It is the earliest one among the m1 numbers in the previous group;

[0162] It is the latest one among the current m1 numbers.

[0163] n represents the number of cumulative calculations, with an initial value of 1;

[0164] This is the average of the earliest m1 numbers.

[0165] Calculate angular velocity and acceleration using a 1Hz low-pass filter (sampling at 400Hz):

[0166] The result of the moving average calculation per Ts is used as input for filtering. Let... Let Ts be the instantaneous moving average angular velocity. Let Lω be the sliding average angular velocity of the previous Ts. i (n) represents the current filtered angular velocity of Ts, Lω i The result of the previous Ts filtering (n-1) is then the discretized filter is:

[0167]

[0168] Where i = x, y, z;

[0169] a=0.9844141274160970.007792936291951550.00779293629195155

[0170] Angular velocity result Lω calculated using low-pass filtering i (n) is the input, and the recursive average angular velocity over 1 second is calculated.

[0171]

[0172] in: The initial value is 0; the initial value of k is 1. After each 1-second recursion, the value is set to 0. k = 1.

[0173] Next, the composite acceleration is calculated using the preprocessed acceleration.

[0174]

[0175] These are the low-pass filtered recursive average values ​​of the mechanical system's horizontal, vertical, and normal axes within the current 1 second.

[0176] The difference between the resultant acceleration and the local gravitational acceleration, along with the pre-processed angular velocity, is used to determine whether the absolute rest condition is met. The threshold is: |g′-g|≤0.01g.

[0177] These are the low-pass filtered recursive average values ​​of the machine system's angular velocities relative to the geographic system's horizontal, vertical, and normal axes within the current 1 second.

[0178] When the state transition period expires (1 second), a state transition occurs:

[0179] Calculate the F matrix

[0180]

[0181] If the current state transition period is determined to be absolutely static, then

[0182] F(13,13)=F(14,14)=F(15,15)=-0.9

[0183] otherwise

[0184] F(13,13)=F(14,14)=F(15,15)=-0.05

[0185] If the system transitions from absolute rest to non-absolute rest, it is necessary to restart the motion acceleration estimation and set P. fmE =P fmN =P fmU =(0.5) 2 .

[0186] Calculate the one-step transfer matrix

[0187] Φ K / K-1 =I+F+F 2 / 2

[0188] Computational driving noise array

[0189] Q K-1 =Q + (FQ + (FQ) T ) / 2

[0190] Calculate one-step prediction

[0191] X K / K-1 =Φ K / K-1 ·X K-1

[0192] Calculate the error variance matrix

[0193]

[0194] Save P after this calculation is complete. K-1 =P K / K-1 X K-1 =X K / K-1

[0195] When the filtering period (1s) expires, the difference between the synthetic acceleration and the local gravitational acceleration, along with the preprocessed angular velocity, is used to determine whether the filtering conditions are met. The threshold is: If the condition is met, attitude filtering and attitude correction calculations will be performed; otherwise, they will not be performed.

[0196] The attitude filtering process is as follows:

[0197] Z-array

[0198]

[0199] H array

[0200]

[0201] Calculate the filter gain matrix

[0202] K = P K / K-1 ·H T ·(H·P K / K-1 ·H T +R) -1

[0203] Calculate state estimation

[0204] X K =X K / K-1 +K·(ZH·X K / K-1 )

[0205] Update the error variance matrix

[0206] P K =(IK·H)·P K / K-1

[0207] Save P after this calculation is complete. K-1 =P K

[0208] Calculate attitude error parameters:

[0209] φ E =X K (5), φ N =X K (6)

[0210] If |X K (i)|≥L ACCi , then δ i =L ACCi ·sign(X K (i)), X K-1 (i)=X K (i)-δ i (i = 5, 6)

[0211] Where L ACC5 =L ACC6 =3°

[0212] Otherwise: δ i =X K (i), XK-1 (i) = 0 (i = 5, 6)

[0213] Use φ E ,φ N Calculate the attitude angles to be corrected, Δθ and Δγ.

[0214] Δθ=-φ E cosψ+φ N sinψ

[0215] Δγ=-φ E sinψ / cosθ-φ N cosψ / cosθ

[0216] Finally, Δθ and Δγ are used to correct the system attitude angles, i.e.

[0217] θ K =θ K-1 -Δθ,γ K =γ K-1 -Δγ

[0218] It should be noted that the above process operations can be combined to varying degrees. For the sake of simplicity, the implementation methods of various combinations will not be elaborated further. Those skilled in the art can flexibly adjust or combine the order of the steps of the above method (or the position of the components of the product) according to the actual situation.

[0219] It should be noted that the functional components shown in the above embodiments can be implemented in hardware, software, or a combination of both. When implemented in hardware, they can be electronic circuits, application-specific integrated circuits (ASICs), plug-ins, function cards, etc. When implemented in software, they can be programs or code segments used to perform the required tasks. The programs or code segments can be stored in a machine or readable medium, or they can be transmitted over a transmission medium or communication link via data signals carried in a carrier wave.

[0220] The above detailed embodiments are a description of the present invention. It should not be considered that the specific embodiments of the present invention are limited to these descriptions. For those skilled in the art, several simple deductions and substitutions can be made without departing from the concept of the present invention, and all of these should be considered to fall within the protection scope of the present invention.

Claims

1. A MEMS attitude correction method based on motion acceleration estimation, characterized in that, Includes the following steps: S1. Acquire angular velocity and acceleration according to a preset period, and perform filtering and moving average processing on the acquired angular velocity and acceleration to obtain filtered and smoothed angular velocity and acceleration. S2, based on the state and variance of the attitude system at the previous moment, as well as the preset noise variance matrix and one-step transition matrix, obtain the predicted state and variance of the attitude system at the current moment; the one-step transition matrix is ​​obtained based on first-order Markov modeling; S3, determine the preset attitude filtering state based on the filtered and smoothed angular velocity and acceleration; establish the Kalman filter gain matrix based on the predicted state of the attitude system at the current moment, the measurement noise matrix and the observation matrix; obtain the attitude filtering measurement matrix based on the filtered acceleration; Based on the attitude and attitude filter measurement matrix, observation matrix, Kalman filter gain matrix, and the predicted state and variance of the attitude and attitude system at the current moment, the state of the attitude and attitude system at the current moment is obtained. The observation matrix is ​​a matrix that separates the motion acceleration; S4 corrects the attitude of the attitude system based on the current state of the attitude system.

2. The method according to claim 1, characterized in that, The method further includes: Based on the attitude and attitude filter measurement matrix, observation matrix, Kalman filter gain matrix, and the predicted state and variance of the attitude and attitude system at the current moment, the variance of the attitude and attitude system at the current moment is obtained.

3. The method according to claim 1, characterized in that, The one-step transition matrix is ​​represented by the following formula: Φ K / K-1 =I+F·T F +(F·T F ) 2 / 2 Where, Φ K / K-1 This represents the state transition from the previous time k-1 to the current time k. I is the unit matrix, and F is the continuous state transition matrix. The elements on the diagonal of the matrix take the value (-α-1) / T. F T F The state transition period is α = 1 / τ, where τ is the correlation time of the Markov process.

4. The method according to claim 1, characterized in that, Based on the filtered and smoothed angular velocity and acceleration, a preset attitude filtering state is determined, including: Based on the filtered and smoothed acceleration, obtain the composite acceleration, and then obtain the absolute value of the difference between the composite acceleration and the gravitational acceleration. Determine whether the absolute value of the difference is less than a given force threshold and whether the filtered and smoothed angular velocity is less than a given angular rate threshold. If so, the attitude filtering state is satisfied.

5. The method according to claim 1, characterized in that, In the observation matrix H[f E ,mE]=H[f N ,mN]=H[f U ,mU]=1; Where H is the observation matrix, mE, mN, and mU represent the east, north, and celestial accelerations, respectively; fE, fN, and fU represent the east, north, and celestial accelerations, respectively.

6. The method according to claim 1, characterized in that, Based on the current state of the attitude system, correct the attitude of the attitude system, including: The eastward and northward platform yaw angles in the current attitude system state are subjected to amplitude limiting processing to obtain the geographic system correction quantity φ. E ,φ N ; According to the geographical system, the quantity to be corrected is φ. E ,φ N Obtain the attitude angles to be corrected, Δθ and Δγ; The attitude of the attitude system is corrected based on the attitude angles to be corrected, Δθ and Δγ.

7. The method according to claim 1, characterized in that, The values ​​of certain elements on the diagonal of the measurement noise matrix, system noise matrix, and system error variance matrix are determined as follows; The diagonal elements of the measurement noise array R are: R=diag(R 11 ,R 22 ,R 33 ) in, The equivalent zero-point mean value is determined by the selected device, where g is the gravitational acceleration. The partial diagonal elements of the systematic error variance matrix P take the following values: Where P[x,y] represents the element in the x-row, y-column of matrix P. The standard deviation of acceleration is estimated from the system's trajectory. The partial diagonal elements of the system noise variance matrix Q are:

8. A MEMS attitude correction device based on motion acceleration estimation, characterized in that, Used to perform the method as described in any one of claims 1-7.

Images