A strapdown inertial navigation heave measurement method based on improved singular spectrum analysis

Abstract

The present invention provides a strapdown inertial navigation heave measurement method based on improved singular spectrum analysis. Conventional singular spectrum analysis is used to perform singular value decomposition on the celestial acceleration time series to obtain an acceleration initialization singular value matrix and an initialization eigenvector matrix. Grouping and reconstruction are performed based on the singular value characteristics to obtain the noise-filtered acceleration, and the velocity is obtained by integration. SVD decomposition is performed on the velocity time series to obtain a velocity initialization singular value matrix and an initialization eigenvector matrix, and initialization is completed. Recursive updating is performed using celestial acceleration drive data for real-time updating to obtain the acceleration singular value matrix and eigenvector matrix at the current moment. Grouping and reconstruction are performed to obtain the noise-filtered acceleration, and the acceleration is integrated to complete the velocity update. The velocity update is recursively used to update the velocity singular value matrix and velocity eigenvector matrix at the current moment, and grouping and reconstruction are performed to obtain the noise-filtered velocity. The velocity is integrated to complete the heave displacement update.

Description

Technical Field

[0001] The invention relates to a strapdown inertial navigation heave measurement method based on improved singular spectrum analysis, and belongs to the technical field of instantaneous linear motion measurement of strapdown inertial navigation equipment. Background Art

[0002] When ships sail at sea, they are affected by environmental factors such as waves and winds, passively generating heave motion along the celestial axis. Real-time and accurate heave measurements play an important role in providing a reference for equipment such as the replenishment of supplies between ocean-going ships, aircraft takeoff and landing, sonar error correction, and heave compensation devices for drilling platforms. Currently, the heave motion of marine vessels is measured primarily using strapdown inertial navigation equipment (SIN). The inertial measurement unit (IMU) measures the vertical acceleration in real time and obtains heave information through quadratic integration. Due to the divergent altitude channel of the inertial navigation system, a digital high-pass filter is required to filter out accelerometer bias and low-frequency noise signals such as Schuler oscillations, Earth oscillations, and Foucault oscillations. However, high-pass filters introduce phase errors, making it difficult to achieve the required accuracy and real-time performance of heave measurements. Therefore, it is necessary to develop a filtering algorithm suitable for heave measurement using SIN.

[0003] Most existing technologies compensate for the amplitude attenuation and phase advance caused by digital high-pass filters. For example, patent application number CN202011358976.8, entitled "A Strapdown Inertial Navigation System Heave Measurement Method with Phase Compensation," measures the dominant wave frequency and designs an all-pass filter to compensate for the phase advance caused by the high-pass filter. However, this does not address the amplitude attenuation issue and is only effective for single-frequency heave. When a maritime vessel is subject to environmental influences such as waves and experiences heave motion, the primary wave energy is concentrated in a single frequency. Accurate phase compensation for a single frequency is insufficient. Another example is patent application number CN201710202159.5, entitled "A Ship Heave Measurement Method Based on Band-Limited Fourier Linear Combinations," which uses a band-limited Fourier linear combination algorithm to model the heave displacement output by the heave filter and performs phase and amplitude compensation at each frequency point in the model. A heave filter is essentially a series connection between a high-pass filter and an integrator, and still suffers from the errors introduced by digital filtering techniques. The Fourier linear combination algorithm and its improved versions provide high-precision frequency estimation of signals. However, they suffer from frequency resolution and picket fence effects, and can only compensate for fundamental frequencies, making them unsuitable for complex sea conditions. For example, the article "Ship Heave Motion Measurement Method Based on IMU and UKF" in Volume 47, Issue 7 of the Journal of Beijing University of Aeronautics and Astronautics analyzes the ship heave motion model, establishes state equations and measurement equations for heave motion measurement, and uses an unscented Kalman filter for filtering. Compared to digital filter algorithms, this algorithm eliminates phase issues. However, the algorithm's accuracy relies on accurately estimating the primary cosine component of the heave waveform, requiring spectral analysis of the heave acceleration signal. Furthermore, using fast Fourier transforms (FFTs) to estimate frequency and phase is not real-time, thus failing to obtain high-precision heave information. Summary of the Invention

[0004] The purpose of the present invention is to provide a strapdown inertial navigation heave measurement method based on improved singular spectrum analysis in order to solve the problems of low accuracy, poor real-time performance and being only effective for single-frequency heave motion in actual offshore heave measurement.

[0005] The object of the present invention is achieved in that the steps are as follows:

[0006] Step 1: Initialize the strapdown inertial navigation device, complete the initial alignment and enter the navigation working mode, collect the output signals of the gyroscope and accelerometer in the three axes of the carrier coordinate system in real time, and use the strapdown navigation solution algorithm to obtain the real-time input celestial acceleration a u .

[0007] Step 2: Initialize the recursive singular spectrum analysis (RSSA), including the acceleration recursive singular spectrum analysis initialization and velocity recursive singular spectrum analysis initialization. Use the output N moments of the celestial acceleration sequence AN =(a u，1 ,a u，2 ,…,a u,N ), construct the trajectory matrix X, perform singular value decomposition (SVD) on X, and obtain the acceleration initialization singular value matrix Σ N And initialize the eigenvector matrix U N .

[0008] Step 3: Group the trajectory matrix X and reconstruct the noise-filtered acceleration sequence A heave =(a heave，1 ,a heave，2 ,…,a heave,N ), for a heave Integrate to obtain the heave velocity sequence B N =(v u,1 ,v u,2 ,…,v u,N ), use this sequence to construct the trajectory matrix Y, perform singular value decomposition (SVD) on Y, and obtain the velocity initialization singular value matrix Σ′ N And initialize the eigenvector matrix U′ N .

[0009] Step 4: Initialize the singular value matrix Σ of the acceleration in step 2 N , initialize the eigenvector matrix U N As the initial state of acceleration, take the celestial acceleration as the driving data, and start to progressively update to obtain the acceleration singular value matrix Σ at time k k and the acceleration eigenvector matrix U k According to the contribution rate of the singular value, the corresponding component is selected, the reconstructed component is calculated, and the filtered acceleration is obtained.

[0010] Step 5: Integrate the filtered acceleration in step 4 to obtain the heave velocity.

[0011] Step 6: Initialize the singular value matrix Σ′ of the velocity in step 3 N , initialize the eigenvector matrix U′ N As the initial state of velocity, start progressive updating to obtain the velocity singular value matrix Σ′ at time k k and the velocity eigenvector matrix U k According to the contribution rate of the singular value, the corresponding component is selected, the reconstructed component is calculated, and the filtered velocity is obtained.

[0012] Step 7: Integrate the filtered velocity in step 6 to obtain the heave displacement.

[0013] Step 8: Transform the acceleration singular value matrix Σ at time k into k and the acceleration eigenvector matrix U k, velocity singular value matrix Σ′ k and the velocity eigenvector matrix U k ' is used as the initial value of the state quantity at the next moment, and steps 4-7 are repeated until the navigation working state ends.

[0014] The present invention also includes:

[0015] 1. The recursive methods for the singular value matrix and eigenvector matrix in steps 4 and 6 are the same, as follows (taking step 4 as an example):

[0016] The acceleration singular value matrix Σ at time k k and the eigenvector matrix U k It can be expressed as Σ k =(Σ k-1 +P), U k =U k-1 (I+Q), Σ k-1 is the singular value matrix at time k-1, U k-1 is the eigenvector matrix at time k-1, I is the identity matrix, and P and Q are small perturbation matrices.

[0017]

[0018]

[0019] Where, ρ i，(k-1) Σ k-1 The i-th element of the diagonal, α i is an intermediate variable The i-th element in a u,k is the celestial acceleration at time k.

[0020] Compared with existing technologies, the present invention offers the following advantages: It utilizes a signal stratification and denoising method based on improved singular spectrum analysis, treating acceleration bias, random errors, and low-frequency noise signals as trend terms, and heave-related acceleration signals as periodic terms. Trend terms are filtered out, while periodic terms are retained. This method is a completely data-driven analysis technique that requires no basis functions and is an independent, adaptive time-domain filtering technique. High-precision heave velocity information is obtained by integrating the filtered acceleration, and recursive singular spectrum analysis is used to output heave information in real time. BRIEF DESCRIPTION OF THE DRAWINGS

[0021] Figure 1 Flowchart of the present invention.

[0022] Figure 2 This is the flow chart of the recursive singular spectrum analysis algorithm. DETAILED DESCRIPTION

[0023] The present invention will be described in further detail below with reference to the accompanying drawings and specific embodiments.

[0024] The present invention is implemented as follows: the inertial measurement element of the strapdown inertial navigation system is fully preheated and the initial alignment of the strapdown inertial navigation system is completed, and then the navigation working state can be entered, at which time the celestial acceleration a can be normally output. u . Collect ten minutes of celestial acceleration as the initialization data for recursive singular spectrum analysis, calculate the acceleration initialization singular value matrix and initialization eigenvector matrix, velocity initialization singular value matrix and initialization eigenvector matrix as the initial state of recursive singular spectrum analysis. After initialization, use the celestial acceleration as the driving data for real-time recursive update, treat the acceleration zero bias, random error and low-frequency noise signal as trend items, and the acceleration signal related to heave as the periodic item, select the periodic component, calculate the reconstructed component, obtain the filtered acceleration, and integrate to output the heave velocity. Perform the same processing on the heave velocity to output the heave displacement. The specific steps are as follows:

[0025] Step 1: Fully preheat the strapdown inertial navigation system and complete the initialization settings. Complete the initial alignment and enter the navigation working mode. Real-time acquisition of the output signals of the gyroscope and accelerometer in the three axes of the carrier coordinate system. Use the strapdown navigation solution algorithm to obtain the real-time input celestial acceleration a. u .

[0026] Step 2: Collect the celestial acceleration a at N moments u As the initialization data for recursive singular spectrum analysis, the one-dimensional acceleration sequence A is used. N =(a u，1 ,a u，2 ,…,a u,N ), initialize the recursive singular spectrum analysis (RSSA), including acceleration recursive singular spectrum analysis initialization and velocity recursive singular spectrum analysis initialization. Use the output celestial acceleration sequence a u , construct the trajectory matrix X, perform singular value decomposition (SVD) on X, and obtain the acceleration initialization singular value matrix Σ N And initialize the eigenvector matrix U N Specifically: using the measured one-dimensional acceleration time series A N =(a u，1 ,a u，2 ,…,a u,N ), construct a multidimensional sequence, namely the trajectory matrix X;

[0027]

[0028] Where L is the selected window length, generally 2≤L≤N / 2, and K=N-L+1.

[0029] Perform singular value decomposition on X and get X = U N Σ N V N T , where U N =(U1 U2…U L ), U i It's XX T Eigenvectors, called left singular vectors of X; is a diagonal matrix, λ i It's XX T The eigenvalues of is called the X singular value; V N =(V1V2…V L ), V i is a column vector, called the right singular vector of X; then The signal composition of X, X=X1+X2+…+X L , i=1,…,L.

[0030] Step 3: Replace the X in step 2 i The purpose of grouping is to separate the target signal component from other signal components, and to convert X=X1+X2+…+X L The signal composition in is divided into trend terms and periodic terms, and the matrix contained in the periodic terms is added to obtain the reconstructed trajectory matrix The matrix Converted into one-dimensional time series data, the noise-filtered acceleration sequence A is obtained. heave =(a heave，1 ,a heave，2 ,…,a heave,N ), for a heave Integrate to obtain the heave velocity sequence B N =(v u,1 ,v u,2 ,…,v u,N ), use this sequence to construct the trajectory matrix Y, perform singular value decomposition (SVD) on Y, and obtain the velocity initialization singular value matrix Σ′ N And initialize the eigenvector matrix U′ N .

[0031]

[0032] Where L, N, and K have the same values as in step 2.

[0033] Perform singular value decomposition on Y and get Y = U′ N Σ′ N V N ' T , where U′ N =(U1′ U2′ … UL ′), U i ' is YY T Eigenvectors are called left singular vectors of Y; is a diagonal matrix, λ i ' is YY T The eigenvalues of is called the Y singular value; V N ′=(V1′V2′…V L ′), V i ′ is a column vector, called the right singular vector of Y.

[0034] Step 4: Use the acceleration initialization singular value matrix and the initialization eigenvector matrix in step 2 as the initial state of acceleration, and start progressive updating to obtain the acceleration singular value matrix Σ at time k k and the acceleration eigenvector matrix U k , k>N. According to the contribution rate of the singular value, select the component corresponding to the periodic term, calculate the reconstructed component, and obtain the filtered acceleration. The specific recursive method is as follows:

[0035] The acceleration singular value matrix Σ at time k k and the eigenvector matrix U k It can be expressed as Σ k =(Σ k-1 +P), U k =U k-1 (I+Q), Σ k-1 is the singular value matrix at time k-1, U k-1 is the eigenvector matrix at time k-1, I is the identity matrix, and P and Q are small perturbation matrices.

[0036]

[0037]

[0038] ρ i，(k-1) Σ k-1 The i-th element of the diagonal, α i is an intermediate variable The i-th element in a u,k is the celestial acceleration at time k.

[0039] Calculate the right singular vector at time k and the k-time trajectory matrix X k Signal composition Complete signal grouping and reconstruction The filtered acceleration

[0040] Step 5: Filter the acceleration a in step 4 heave,kIntegrate to get the heave velocity v u,k .

[0041] Step 6: Initialize the singular value matrix Σ′ of the velocity in step 3 N , initialize the eigenvector matrix U′ N As the initial state of velocity, start progressive updating to obtain the velocity singular value matrix Σ′ at time k k and the velocity eigenvector matrix U k ′. According to the contribution rate of the singular value, select the component corresponding to the periodic term, calculate the reconstructed component, and obtain the filtered velocity. The specific recursive method is similar to step 4:

[0042] Velocity singular value matrix Σ′ at time k k and the eigenvector matrix U k ′ can be expressed as Σ′ k =(Σ′ k-1 +P′), U k ′=U k ' -1 (I+Q′),Σ′ k-1 is the singular value matrix at time k-1, U k ' -1 is the eigenvector matrix at time k-1, I is the identity matrix, and P′ and Q′ are small perturbation matrices.

[0043]

[0044]

[0045] ρ i ' (k-1) is Σ′ k-1 The i-th element of the diagonal, υ i is an intermediate variable The i-th element in v u,k is the heave speed at time k.

[0046] Calculate the right singular vector at time k and the k-time trajectory matrix Y k Signal composition Complete signal grouping and reconstruction The filtered speed

[0047] Step 7: Filter the velocity in step 6 The integral gives the heave displacement.

[0048] Step 8: Transform the acceleration singular value matrix Σ at time k into k and the acceleration eigenvector matrix U k , velocity singular value matrix Σ′k and the velocity eigenvector matrix U k ' is used as the initial value of the state quantity at the next moment, and steps 4-7 are repeated until the navigation working state ends.

[0049] This completes the strapdown inertial navigation heave measurement method based on improved singular spectrum analysis.

[0050] In order to illustrate the effectiveness of the algorithm, the algorithm was simulated. The simulation conditions were set as follows: the frequency band of the ship's heave motion is 0.05-0.25 Hz, the amplitude is 1 meter, the heave acceleration is the superposition of the second-order derivative of the true value of the heave displacement and the low-frequency noise, the data sampling frequency is 100 Hz, and the simulation time is 7200 seconds. Compared with the traditional quaternion fourth-order Runge-Kutta method and the equivalent rotating vector method of angular rate fitting, the results are shown in Table 1: (Group A is the simulation result of the traditional high-pass filter method, Group B is the simulation result of the digital filter method with phase compensation, and Group C is the simulation result of the present invention)

[0051] Table 1 Algorithm error comparison

[0052]

[0053] In summary, the present invention provides a method for heave measurement of a strapdown inertial navigation device. Based on an improved singular spectrum analysis method, the celestial acceleration output by the strapdown inertial navigation device is used as the driving data to output the heave velocity and heave displacement in real time. The entire algorithm includes recursive singular spectrum analysis (RSSA) initialization and recursive update. First, the celestial acceleration a is analyzed using conventional singular spectrum analysis (SSA). u Perform singular value decomposition (SVD) on the time series to obtain the acceleration initialization singular value matrix and initialization eigenvector matrix. Reconstruct the noise-filtered acceleration by grouping based on the singular value characteristics, and integrate to obtain the velocity. Perform SVD decomposition on the velocity time series to obtain the velocity initialization singular value matrix and initialization eigenvector matrix, completing initialization. Recursively update the initialization singular value matrix and eigenvector matrix as the initial state, and use the celestial acceleration drive data for real-time updates to obtain the acceleration singular value matrix and eigenvector matrix at the current moment. Reconstruct the noise-filtered acceleration by grouping, and complete the velocity update by integrating the acceleration. Recursively update the velocity singular value matrix and velocity eigenvector matrix at the current moment using the velocity update, and reconstruct the noise-filtered velocity by grouping. Update the heave displacement by integrating the velocity.

Claims

1. A strapdown inertial navigation heave measurement method based on improved singular spectrum analysis, characterized in that: Here are the steps: Step 1: Initialize the strapdown inertial navigation device, complete the initial alignment and enter the navigation working mode, collect the output signals of the gyroscope and accelerometer in the three axes of the carrier coordinate system in real time, and use the strapdown navigation solution algorithm to obtain the real-time input celestial acceleration a u ; Step 2: Initialize the recursive singular spectrum analysis (RSSA), including the initialization of acceleration recursive singular spectrum analysis and velocity recursive singular spectrum analysis; use the output N moments of the celestial acceleration sequence A N =(a u，1 ,a u，2 ,…,a u,N ), construct the trajectory matrix X, perform singular value decomposition (SVD) on X, and obtain the acceleration initialization singular value matrix Σ N And initialize the eigenvector matrix U N ; Step 3: Group the trajectory matrix X and reconstruct the noise-filtered acceleration sequence A heave =(a heave，1 ,a heave，2 ,…,a heave,N ), for a heave Integrate to obtain the heave velocity sequence B N =(v u,1 ,v u,2 ,…,v u,N ), use this sequence to construct the trajectory matrix Y, perform singular value decomposition (SVD) on Y, and obtain the velocity initialization singular value matrix Σ′ N And initialize the eigenvector matrix U′ N ; Step 4: Initialize the singular value matrix Σ with the acceleration in step 2 N , initialize the eigenvector matrix U N As the initial state of acceleration, take the celestial acceleration as the driving data, and start to progressively update to obtain the acceleration singular value matrix Σ at time k k and the acceleration eigenvector matrix U k ; According to the contribution rate of the singular value, select the corresponding component, calculate the reconstructed component, and obtain the filtered acceleration; Step 5: Integrate the filtered acceleration in step 4 to obtain the heave velocity; Step 6: Initialize the singular value matrix Σ′ of the velocity in step 3 N , initialize the eigenvector matrix U′ N As the initial state of velocity, start progressive updating to obtain the velocity singular value matrix Σ′ at time k k and the velocity eigenvector matrix U k '; According to the contribution rate of the singular value, select the corresponding component, calculate the reconstructed component, and obtain the filtered velocity; Step 7: Integrate the filtered velocity in step 6 to obtain the heave displacement; Step 8: Transform the acceleration singular value matrix Σ at time k into k and the acceleration eigenvector matrix U k , velocity singular value matrix Σ′ k and the velocity eigenvector matrix U k ' is used as the initial value of the state quantity at the next moment, and steps 4-7 are repeated until the navigation working state ends.

2. The strapdown inertial navigation heave measurement method based on improved singular spectrum analysis according to claim 1, characterized in that: The recursive method for singular values and eigenvectors in step 4 is as follows: The acceleration singular value matrix Σ at time k k and the eigenvector matrix U k It can be expressed as Σ k =(Σ k-1 +P), U k =U k-1 (I+Q), Σ k-1 is the singular value matrix at time k-1, U k-1 is the eigenvector matrix at time k-1, I is the identity matrix, P and Q are small perturbation matrices; Where, ρ i，(k-1) Σ k-1 The i-th element of the diagonal, α i is an intermediate variable The i-th element in a u,k is the celestial acceleration at time k.

Images