A data analysis method for inertial sensors based on frequency domain processing

Abstract

The present application relates to the technical field of inertial sensor state data analysis, and discloses a data analysis method for an inertial sensor based on frequency domain processing, which realizes frequency domain analysis on the actual error composition of the inertial sensor, thereby improving the accuracy of effective coefficient performance analysis of the inertial sensor. Specifically, S1, a MEMS inertial measurement unit of the sensor is used to collect original output signals; S2, a time domain signal in the original output signal is decomposed into a solid mode equation and a residual term through a modal decomposition algorithm; the solid mode equation is extracted through empirical mode decomposition using a screening algorithm, and the residual is updated; the center frequency is selected using signal characteristics and prior knowledge through variational modal decomposition, the solid mode equation at the center frequency is extracted, and modeling is performed as a function of amplitude, phase and frequency; S3, Hilbert spectrum of time-frequency domain joint distribution is generated by performing Hilbert-Huang transformation on the solid mode equation component; and S4, a sensor performance evaluation model is constructed based on the Hilbert spectrum.

Description

Technical Field

[0001] This invention relates to the field of inertial sensor state data analysis technology, and more specifically to a data analysis method for inertial sensors based on frequency domain processing. Background Technology

[0002] An Inertial Measurement Unit (IMU) is an instrument that measures the motion of a vehicle based on Newton's law of inertia. Common IMUs include gyroscopes and accelerometers, which measure the angular velocity and linear acceleration of the vehicle, respectively. MEMSMU, on the other hand, refers to the integration of inertial sensors with micromechanical structures, characterized by low cost and small size.

[0003] However, compared to traditional measuring instruments, the readings of MEMSIMUs often contain a large number of different types of noise, which is particularly noticeable in low-cost products. Therefore, in actual production and applications, algorithm software is used to evaluate the sensor status and model the error. Common analysis methods in the industry include Allen's analysis of variance and standard deviation analysis. Common error models include polynomial interpolation or linear interpolation.

[0004] Common inertial measurement units (IMUs) analysis methods struggle to address the fact that error components such as output signal, white noise, and cumulative walk noise are superimposed in the time domain, making direct separation and identification difficult. This results in an inability to effectively analyze the actual error composition and the contribution of each error. Summary of the Invention

[0005] The purpose of this invention is to provide a data analysis method for inertial sensors based on frequency domain processing, and to solve the following technical problems: How to accurately distinguish the output time signal and perform frequency domain analysis on the actual error composition of the inertial sensor, thereby improving the effective performance analysis of the inertial sensor.

[0006] The objective of this invention can be achieved through the following technical solutions: A data analysis method for inertial sensors based on frequency domain processing, the method comprising: S1. The original output signal is acquired by the MEMS inertial measurement unit of the sensor; S2. The time-domain signal in the original output signal is decomposed into solid-state mode equations and residual terms through mode decomposition algorithms. Mode decomposition algorithms include empirical mode decomposition and variational mode decomposition. The mode decomposition algorithms obtain IMFs through the sieve factorial free factorial (SIFT) method. The sieve factorial free factorial (SIFT) method first selects the IMFs that meet two properties, including: each consecutive local extremum must be sandwiched between a zero point and the local mean of each IMF must be zero. Empirical Mode Decomposition (EMD) uses numerical iterative methods to extract the frequency-banded portion of the original signal, thereby extracting solid-state mode equations. Variational mode decomposition uses signal characteristics and prior knowledge to select the center frequency, extracts the solid-state mode equations at the center frequency, and models them as functions of amplitude, phase, and frequency. S3. Perform Hilbert transform on the solid-state mode equation components to generate a Hilbert spectrum with joint time-frequency domain distribution; S4. Construct a sensor performance evaluation model based on the Hilbert spectrum.

[0007] Preferably, the implementation steps of empirical mode decomposition are as follows: S21. Find the maximum and minimum values ​​of each interval using cubic spline interpolation, classify the maximum and minimum values, construct upper and lower envelopes, and obtain the cubic spline interpolation coefficients of all maximum and minimum values. S22. Perform LU decomposition using the Doolittle algorithm to obtain an upper triangular matrix U and a lower triangular matrix L, satisfying A=LU. Initialize L as an identity matrix of dimension n, and U as an empty matrix of dimension n. The calculation formula is: , ; where the subscripts i, j, and k all represent iteration variables and satisfy: Perform a forward replacement; Perform backward calculation; substitute the values ​​that satisfy the condition. A system of linear equations of the form x is used to obtain the final result x. S23. Calculate the average value of the upper and lower envelopes of a continuous interval. The algorithm determines whether the obtained function is an IMF. If the relative energy ratio is less than a preset threshold or the number of iterations exceeds the preset maximum number of iterations, the iteration stops, and the average signal is subtracted from the residual signal to obtain an IMF. S24. Continue to repeat step S23 iterative sieving until the external loop iteration stop condition is met; the external loop stop condition is that the number of solid-state model equations reaches the preset maximum number, or the energy ratio is less than the set ratio, or the number of poles is less than the preset number.

[0008] Preferably, the cubic spline interpolation method employs natural boundary conditions, where the second derivative is zero at the boundaries of the interpolation interval, and the interpolation function... Because the interpolation intervals overlap, and The linear system of interpolation coefficients includes: constructing equations based on second derivatives. , ,Bundle Substitute and get Construct an equation about the second derivative for the four parameters a, b, c, and d. , ; According to the definition Substituting a, c, and d into the equation, we can obtain... ; According to the first derivative equation Substituting the forms of a, b, c, and d with respect to their second derivatives, we get... ; Based on this equation, a system of linear equations can be constructed to solve for all second derivatives: .

[0009] Preferably, the implementation steps of variational mode decomposition are as follows: select the center frequency using signal characteristics and prior knowledge, extract the solid-state mode equation at the center frequency, and model it as a function of amplitude, phase, and frequency; SS1. Construct the objective function of the constrained variational problem, and apply it to the IMF components. and specific center frequency Solution: ; Among them, It is the k-th IMF component. It is the specific center frequency of the k-th IMF. For the Dirac function, It is the imaginary unit, and ; For time differential operators, for Norm, This represents the change in the frequency domain obtained by modulating the modulo signal u onto a carrier wave; SS2. Iteratively solve for the optimal solution set using the alternating direction multiplier method. .

[0010] Preferably, SS2 includes: Update the IMF frequency domain expression: ;in, It is relative to continuous frequency Discrete frequencies; It is relative to continuous frequency The discrete frequency of the k-th iteration; n is the iteration number; Discretized Fourier transform; Update center frequency: ;in, It is a continuous frequency, and , The frequency domain of the k-th IMF component serves as the input source for generating the Hilbert spectrum; The center frequency of the k-th IMF component; Update the Lagrange multipliers: ;in, Update the step size for the multiplier; The conjugate of the frequency-domain Lagrange multiplier serves as a constraint force adjustment factor. To augment the Lagrange function For conjugate Lagrange multipliers The gradient; The frequency domain of the original output signal is used as the Fourier transform of the sensor output. These are the Lagrange multipliers before the update; For the updated Lagrange multipliers; and Discretization and subsequent input yield the following results: .

[0011] Preferably, the Hilbert transform includes: The continuous Hilbert transform is defined as: ;in, To update the step size for the multiplier, For the first IMF frequency components of a step size; This represents the result of a 90° phase shift in all frequency components of the original signal after the Hilbert transform. Constructing analytical signals: ,in, The original signal, To analyze the signal, The imaginary unit; Transform polar coordinates: ,in It is the instantaneous amplitude. It is the instantaneous frequency.

[0012] Preferably, the evaluation conditions for constructing the sensor performance evaluation model based on the Hilbert spectrum are as follows: When the Hilbert spectrum of EMD-HHT has no spikes and a uniform amplitude distribution at the same amplitude level in the full time domain and the full frequency domain, it is judged to be in normal working condition; otherwise, it is judged to be at risk of failure in working condition. When the Hilbert spectrum of VMD-HHT has a uniform amplitude distribution at a specific frequency, it is considered to be in normal working condition; otherwise, it is considered to be in a faulty working condition.

[0013] Preferably, S4 further includes compressing the Hilbert spectrum using numerical methods: Rational number interpolation is constructed based on Chebyshev's alternation theorem: ; According to the original function: Error equations obtained through rational number interpolation ;in, , , Representing elements respectively The natural number to the power of, , , They represent , , polynomials; The extreme points of the error are calculated by Remez iteration and interpolation is continued to satisfy the Shchev alternation theorem. After convergence, the optimal solution is obtained.

[0014] Preferably, the method for solving for rational number interpolation coefficients includes: The error bit is obtained by simplifying the expression using the properties of Chebyshev polynomials. Given k = 1, 2, 3...N, construct an N+1 dimensional system of linear equations. ; Map the obtained Chebyshev nodes to the frequency domain: Let i = 1, 2, 3...m; Construct the matrix: ; Using normal equations Solve the system of linear equations to obtain the coefficients that satisfy the definition of Chebyshev rational number interpolation. and ;in, This indicates solving the system of equations. express The transpose of .

[0015] The beneficial effects of this invention are: (1) This invention can obtain different results that satisfy the basic formula by considering different numerical methods. The time domain signal is decomposed into solid mode equations and residual terms by mode decomposition algorithm. The combination analysis of empirical mode decomposition (EMD) and variational mode decomposition (VMD) can provide accurate analysis of the signal from different angles.

[0016] (2) This invention generates a Hilbert spectrum with joint time-frequency distribution by performing Hilbert transform on the IMF components. Based on the Hilbert spectrum, a sensor performance evaluation model is constructed. By combining time information, frequency information, and amplitude information through the Hilbert spectrum, the changes in instantaneous frequency and instantaneous amplitude can be observed, which provides significant assistance for the analysis of sensor data. The modulus separation methods of EMD and VMD bring different evaluation perspectives to the interpretation of data. Firstly, EMD-HHT decomposes the DC residual, which can be used to decompose the fixed deviation and wandering information. In the Hilbert spectrum of EMD-HHT, if there are no peaks and the amplitude distribution is uniform at the same amplitude level in the full time and full frequency domains, it is considered to have good performance and normal working status. Secondly, in the Hilbert spectrum of VMD-HHT, if the amplitude distribution is at a specific frequency, such as the working frequency, and the amplitude is relatively uniform in the time domain, it is considered to have good performance and normal working status. HHT can provide significant assistance for sensor performance analysis and fault diagnosis.

[0017] Of course, any product implementing this invention does not necessarily need to achieve all the advantages described above at the same time. Attached Figure Description

[0018] To more clearly illustrate the technical solutions of the embodiments of the present invention, the accompanying drawings used in the description of the embodiments will be briefly introduced below. Obviously, the drawings described below are only some embodiments of the present invention. For those skilled in the art, other drawings can be obtained based on these drawings without creative effort.

[0019] Figure 1 This is a flowchart illustrating the steps of a data analysis method for an inertial sensor based on frequency domain processing according to the present invention. Figure 2 This diagram illustrates the implementation steps of the empirical mode decomposition method of this invention. Detailed Implementation

[0020] 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 embodiments of the present invention, and not all embodiments. 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.

[0021] Please see Figure 1 As shown, this invention provides a data analysis method for inertial sensors based on frequency domain processing. The method steps are as follows: S1. The original output signal is acquired by the MEMS inertial measurement unit of the sensor; S2. The time-domain signal in the original output signal is decomposed into solid-state mode equations and residual terms through mode decomposition algorithms. Mode decomposition algorithms include empirical mode decomposition and variational mode decomposition. The mode decomposition algorithms obtain IMFs through the sieve factorial free factorial (SIFT) method. The sieve factorial free factorial (SIFT) method first selects the IMFs that meet two properties, including: each consecutive local extremum must be sandwiched between a zero point and the local mean of each IMF must be zero. Empirical Mode Decomposition (EMD) uses numerical iterative methods to extract the frequency-banded portion of the original signal, thereby extracting solid-state mode equations. Variational mode decomposition uses signal characteristics and prior knowledge to select the center frequency, extracts the solid-state mode equations at the center frequency, and models them as functions of amplitude, phase, and frequency. S3. Perform Hilbert transform on the solid-state mode equation components to generate a Hilbert spectrum with joint time-frequency domain distribution; S4. Construct a sensor performance evaluation model based on the Hilbert spectrum.

[0022] In the above technical solutions, considering that the output objects of the sensors (accelerometers and gyroscopes) can be subdivided into the following different areas, namely output signal, fixed deviation, white noise, accumulated wander noise, etc., it is difficult to directly distinguish them in the output time signal, so it is difficult to analyze the actual error composition of the sensor.

[0023] Against this technical backdrop, the designed data analysis method for inertial sensors based on frequency domain processing first acquires the raw output signal using a MEMS inertial measurement unit in step S1. Then, in step S2, the output signal is analyzed in the frequency domain, dividing the frequency domain information into several different single-frequency modes, called Intrinsic Mode Functions (IMFs). The basic formula is as follows: ,in The output time-domain signal, It's the IMFs. For the number of IMFs, and ∈ , It is the residual; in the formula The result indicates that the time-domain output signal is the sum of all IMFs and residuals. Since IMFs are artificially created through numerical methods and do not strictly correspond to actual physical quantities, the choice of IMFs is not unique. Different numerical methods can yield different results that satisfy the basic formula, providing different perspectives for analysis and interpretation. The time-domain signal is decomposed into solid-state mode equations and residual terms through mode decomposition algorithms. Based on this, different methods include Empirical Mode Decomposition (EMD) and Variational Mode Decomposition (VMD).

[0024] In one embodiment, Method 1, IMF is determined through EMD: IMFs are extracted using a screening algorithm, and the residuals are updated; Empirical Mode Decomposition (EMD) iteratively extracts the frequency-carrying components from the original signal step by step. Since the residuals are non-frequency DC signals, while VMD artificially selects key center frequencies based on signal characteristics and prior knowledge, it extracts the modulus at the center frequency. The residuals represent the remaining, less important frequency information. Through these characteristics, different information can be extracted, thus obtaining IMFs. The first step of EMD is an algorithm called Hilbert Huang Transform (HHT). EMD obtains IMFs through an algorithm called Sieve-Induced Frames (SIFT). Therefore, the IMFs selected by the SIFT algorithm must meet two properties: First, each consecutive local extremum must be separated by a zero point; Second, the local mean of each IMF must be zero; For details, please refer to Figure 2 As shown, the implementation steps of empirical mode decomposition are as follows: S21. Find the maximum and minimum values ​​of each interval using cubic spline interpolation, classify the maximum and minimum values, construct upper and lower envelopes, and obtain the cubic spline interpolation coefficients of all maximum and minimum values. S22. Perform LU decomposition using the Doolittle algorithm to obtain an upper triangular matrix U and a lower triangular matrix L, satisfying A=LU. Initialize L as an identity matrix of dimension n, and U as an empty matrix of dimension n. The calculation formula is: , ; where the subscripts i, j, and k all represent iteration variables and satisfy: Perform a forward replacement; Perform backward calculation; substitute the values ​​that satisfy the condition. A system of linear equations of the form x is used to obtain the final result x. S23. Calculate the average value of the upper and lower envelopes of a continuous interval. The algorithm determines whether the obtained function is an IMF. If the relative energy ratio is less than a preset threshold or the number of iterations exceeds the preset maximum number of iterations, the iteration stops, and the average signal is subtracted from the residual signal to obtain an IMF. S24. Continue to repeat step S23 iterative sieving until the external loop iteration stop condition is met; the external loop stop condition is that the number of solid-state model equations reaches the preset maximum number, or the energy ratio is less than the set ratio, or the number of poles is less than the preset number.

[0025] In the above technical solution, firstly, the local maximum and minimum values ​​are found. The maximum and minimum values ​​are then classified and upper and lower envelopes are constructed, thus obtaining the cubic spline interpolation coefficients for all maximum and minimum values. Specifically, the cubic spline interpolation method calculates an independent cubic function in each small interval, ensuring that the curve is smooth over the entire interval—that is, the original function is smooth, the first derivative is smooth, and the second derivative is smooth.

[0026] To satisfy this condition, a natural boundary condition needs to be introduced, that is, the second derivative is zero at the boundary of the interpolation interval. In this case, the interpolation function... Because the interpolation intervals overlap, and Next, let the interpolation function be... Therefore, each cubic function segment requires four parameters a, b, c, and d; according to the Power Rule, we can deduce... , Let the step size be... Let the second derivative Therefore, a first-order derivative equation can be constructed at the points of overlap in the interpolation interval. Substituting the interpolation definition into the equation, we get... , ,so Construct equations based on second derivatives using the same method. ,so ,Bundle By substituting, we can rewrite it and obtain... It is possible to construct equations about the second derivative for all four parameters a, b, c, and d. , Then, according to the definition, Substituting a, c, and d into the equation, we can obtain... According to the first derivative equation Substituting the forms of a, b, c, and d with respect to their second derivatives, we obtain the equation. Based on this equation, a system of linear equations can be constructed to solve for all second derivatives: This linear system has a unique solution. There are many methods for solving linear systems. In this technique, the Doolittle algorithm is adopted by comprehensively considering numerical stability, algorithm complexity, and the form of the system. It performs LU decomposition to obtain an upper triangular matrix U and a lower triangular matrix L, satisfying A = LU, and then solves the linear system. Since the Doolittle algorithm is an iterative algorithm in sequence, it requires that the matrix being decomposed be a square matrix, i.e., the number of rows equals the number of columns. First, L is initialized as an identity matrix of dimension n, and U is an empty matrix of dimension n. The calculation formula is as follows: , Where i, j, k represent iteration variables, and the calculated results satisfy... Perform a forward replacement. Perform backward calculations, substituting the original values ​​that satisfy the condition. The linear equations of this form yield the final result x; in cubic spline interpolation, this is the required second-order guiding quantity J; based on this calculation, the complete expressions for the upper and lower envelopes are obtained. , This is the first step.

[0027] At this point, summing the expressions for the upper and lower envelopes over a continuous interval and dividing by two yields the average value of the upper and lower envelopes. After obtaining the expression for the average value, analyze whether the obtained function is an IMF. If the relative energy ratio (the ratio of the energy of the average signal to the energy of the original signal) is less than the set threshold, or the number of iterations exceeds the set maximum number of iterations, then the current iteration stops. Subtract the average signal from the residual signal (if it is the first iteration, the residual signal is the original signal) to obtain an IMF. This is the second step.

[0028] Repeat the above process until the outer loop iteration stops: the number of IMFs reaches the preset maximum number, or the energy percentage is less than the set percentage, or the number of poles is less than the preset number. Then the EMD algorithm is completed; this is the third step.

[0029] In one embodiment, Method 2 involves determining the IMF via VMD: A center frequency is selected based on signal characteristics and prior knowledge; the IMF at the center frequency is extracted and modeled as a function of amplitude, phase, and frequency; Variational Mode Decomposition (VMD): VMD selects a series of IMFs satisfying the form... Where A represents the envelope amplitude, The mode represents the phase and satisfies the following two properties: first, the envelope energy is positive and changes slowly; second, the instantaneous frequency of the mode is non-decreasing and concentrated at a certain center frequency. As can be seen above, unlike EMD which models the IMF using numerical methods, VMD models the IMF as a function of amplitude, phase, and frequency. This is the core difference between the two methods.

[0030] Specifically, the methods for determining the IMF through VMD include the following: SS1. Construct the objective function for the constrained variational problem: ; SS2. Iteratively solve for the optimal solution set using the alternating direction multiplier method. .

[0031] In the above technical solution, the definition of VMD requires that the decomposed IMFs satisfy the minimum bandwidth property. Therefore, the process of finding IMFs can be transformed into a constrained variational problem. The transformation method is to first list the objective function of the optimization problem. This modulates the analog signal to the center frequency, converting the original analog signal to the baseband signal that needs processing. It is the k-th IMF component. It is the specific center frequency of the k-th IMF. It is the imaginary unit, and , For time differential operators, for Norm, This represents the frequency domain variation obtained by modulating the modulo signal u onto the carrier wave; its physical meaning is the baseband bandwidth. VMD requires the obtained IMFs to have the minimum bandwidth, so the optimization problem can be understood as needing to optimize the IMFs to the minimum bandwidth. and The solution is needed to satisfy the condition that the sum of the bandwidths of all moduli reaches the theoretical minimum. In practical applications, to prevent redundancy, the Hilbert transform is required to obtain the analytic signal, and the original expression becomes... ,in The objective function is the Dirac function, and a 90-degree phase shift is added, thus constructing the objective equation for the optimization problem.

[0032] Constraint Problem: According to the definition of VMD, the sum of all IMFs should equal the original signal; therefore, this property can be used as a constraint condition to construct a constraint function expression. ,in, The original signal represents the output time-domain signal. Let k be the k-th IMF component, and This term evaluates the squared error of the sum of the original signal and all the decomposed IMFs, and the optimal solution needs to be found through an optimization process.

[0033] Lagrange multiplier term: Another Lagrange multiplier term is needed to limit errors in the signal. This is because the definition of decomposition introduces a residual term *r*. If the residual is zero, this term can be omitted. However, in actual calculations, undecomposed clutter inevitably exists in the signal. The expression for the Lagrange multiplier term is: ,in For Lagrange multipliers, another term It is a residual.

[0034] Construct the augmented Lagrange equation: Finally, add a penalty factor to the above three expressions. α The final function expression is:

[0035] The process of finding IMFs is essentially the process of solving this function: specifically, firstly, due to the constraints of VMD on frequency domain characteristics, in order to reduce marginal effects, it is necessary to first extend the function in both directions, using the midpoint of the signal as the boundary; then, the augmented Lagrangian function needs to be transformed into the frequency domain before further processing, which yields... , , ; The output time-domain signal, For the k-th IMF component, Let be the Lagrange multipliers; then, the augmented Lagrange function can be written as: By unifying the integration interval through symmetry, the above equation can be rewritten as: Next, the Alternating Direction Multiplier Method (ADMM) is used to iteratively solve the equation. The core idea of ​​this algorithm is to fix one variable, differentiate the expression for the other variable, set the result to 0, and then solve the equation. First, fixed ,right Differentiate: Solving the equation yields:

[0036] After discretization, the following is obtained during the iterative process:

[0037] in, For time differential operators, α As a penalty factor, It is the specific center frequency of the k-th iteration. For continuous frequency, ; Discrete frequency; It is relative to continuous frequency The discrete frequency of the k-th iteration; n is the iteration number; for Discrete values ​​of the Fourier transform; For the conjugate of the frequency domain Lagrange multipliers, Indicates the optimal solution; Second, fixed right Differentiate:

[0038] Solving the equation, we get:

[0039] Converted to a discrete approximate expression:

[0040] Substituting the discretization into the iterative process yields:

[0041] Third, the Lagrange multiplier update:

[0042] in , Update the step size for the multiplier; The conjugate of the frequency domain Lagrange multiplier serves as a constraint force adjustment factor; The frequency domain of the original output signal is used as the Fourier transform of the sensor output. For continuous frequency, Discrete frequency; After discretization and substitution, we get:

[0043] When relative energy threshold and absolute energy threshold satisfy:

[0044]

[0045] In summary, the VMD algorithm is complete; where u represents the decomposed modulus. This indicates that the optimization process has reached a plateau and the iteration is complete, represented by a relative threshold. This represents the absolute threshold.

[0046] Step S3 involves performing a Hilbert transform on the IMF components to generate a Hilbert spectrum of the joint time-frequency distribution; the Hilbert transform includes: Define the continuous Hilbert transform as , The transformed signal, in its physical sense, represents a 90° phase shift across all frequency components of the signal; where, This represents the convolution stride; in the discrete implementation, a discrete Hilbert transform filter needs to be constructed first, when 0 < k, N / 2. When N / 2 < k < N, When k = 0, N / 2, That is, if the positive frequency components are phase-shifted by -90° and the negative frequencies by 90°, then DC and Nyquist are not processed. Compared to performing convolution in the time domain, processing in the frequency domain is more convenient; simply multiply the frequency domain signal by the Hilbert transform filter. Then you can get the result by performing IDFT. ; Let n be the discrete digital signal, and n represents the digital signal.

[0047] In practice, the Hilbert transform is used in HHT to construct analytic signals: ,in, To analyze the signal, we use a complex signal representation. The imaginary unit is used; the orthogonal component is labeled; the analytic signal is obtained by transforming it from Cartesian coordinates to polar coordinates. ,in It is the instantaneous amplitude, representing the signal envelope; It is the instantaneous frequency, representing the frequency change. The complex exponent is used as a phase rotation factor; finally, the Hilbert spectrum is derived: ,in, This represents the Hilbert spectrum, thus completing the HHT algorithm.

[0048] In step S4, a sensor performance evaluation model is constructed based on the Hilbert spectrum. By combining time, frequency, and amplitude information through the Hilbert spectrum, changes in instantaneous frequency and amplitude can be observed, significantly aiding in the analysis of sensor data. The EMD and VMD modeling methods offer different evaluation perspectives for data interpretation. Firstly, EMD-HHT decomposes the DC residual, which can be used to decompose fixed bias and wandering information. In the Hilbert spectrum of EMD-HHT, if there are no spikes and the amplitude distribution is uniform at the same amplitude level across the entire time and frequency domains, the performance is considered good and the operating condition is normal. Secondly, in the Hilbert spectrum of VMD-HHT, if the amplitude distribution is relatively uniform in the time domain at a specific frequency, such as the operating frequency, the performance is considered good and the operating condition is normal. In summary, HHT provides significant assistance in sensor performance analysis and fault diagnosis.

[0049] Specifically, the numerical methods for implementing the above technical process are as follows: Because the dimension of the Hilbert spectrum is m×n, where m is the number of frequency bins and n is the length of the time-domain signal, in practical engineering, even though the Hilbert spectrum is generally sparse, it still contains a lot of information. Therefore, in order to compress the information, numerical methods are used to compress the spectrum. By using Chebyshev's alternation theorem: If there is a function There are polynomials Cross The set of powers of natural numbers. Then the error. ,if: ,satisfy The smallest is called Yes The optimal estimate; if there are n points on the function's definition interval C that satisfy the function If the signs of the functions (i=1,2,3,...,n) are opposite at adjacent points and their absolute values ​​are equal to the maximum value, then the function satisfies the alternating extremum condition. This set of points is called the Chebyshev alternating point set; Chebyshev's alternation theorem specifies that... in the case of The optimal estimate requires that the polynomial has at least n+2 points forming an alternating set of points; and the rational function has at least m+n+2 points forming an alternating set of points. Satisfying this condition guarantees, by definition, that the error curve has an equiripple structure and that the maximum error is minimized.

[0050] Chebyshev rational number interpolation: The rational number interpolation equation constructed based on the above conditions can achieve higher accuracy with fewer solutions compared to polynomial interpolation, and compared to the Pad approximation, it can obtain better global consistency and equiripple structure. It is particularly advantageous in sensor error analysis because the spectrum is very sensitive to additional introduced errors; isolating numerical errors allows the compressed information to better match the original signal. The interpolation method that satisfies the above conditions is Chebyshev rational number interpolation, specifically: First, write out the polynomial rational form. The original function f has a Maclaurin series form: The error equation is ,in, , , Representing elements respectively The natural number to the power of, , , They represent , , The polynomial; using the properties of polynomials to simplify the expression and obtain the error bit. (k = 1, 2, 3...N), thus we can construct an N+1 dimensional system of linear equations. .

[0051] The extreme points of the error are calculated by Remez iteration and interpolation is continued to satisfy the Shchev alternation theorem. After convergence, the optimal solution is obtained.

[0052] Chebyshev polynomials: Chebyshev polynomials Within the interval (-1, 1), relative to the weight function Orthogonal; the Chebyshev polynomial is defined as follows: According to the definition, the values ​​of the first two Chebyshev polynomials can be obtained. With the values ​​of the first two polynomials, we can use the iterative relation. Determine the value of the polynomial where n is any positive integer (n>1); Basis transformation: First, convert the monomial basis into a Chebyshev polynomial basis. , ( ), The orthogonality of Chebyshev polynomials can be used to derive... , Then, the x-coordinates of the original data are mapped to Chebyshev nodes, and... Replace with The integral can be rewritten. In order to calculate the efficiency coefficient Multiply by 2 In this way, the data points are mapped to Chebyshev nodes; further, the Gauss-Chebyshev quadrature formula is used. After discarding the error term, substituting the original expression, we can obtain... , After obtaining the coefficient 'a', a matrix can be constructed:

[0053] Finally, the normal equation can be used. Solving this linear system yields coefficients that satisfy the definition of Chebyshev rational number interpolation. and ;in, This indicates solving the system of equations. express The transpose of .

[0054] Furthermore, Chebyshev's rational function interpolation at this point cannot provide the optimal solution that minimizes the error; Remez iteration is still required to obtain the optimal solution. After obtaining Chebyshev's rational number interpolation, the error is calculated. The extreme points of the error are calculated, and the original data is interpolated again at these extreme points until Chebyshev's alternation theorem is satisfied. According to the theorem, the converged result is theoretically the global optimal solution.

[0055] In implementation, it is necessary to select energy feature points within the Hilbert spectrum to extract amplitude and phase to establish interpolation, ensuring good noise immunity and good signal fidelity.

[0056] This concludes the methodology. This method enables the analysis of sensor data from frequency domain analysis and data reconstruction to model building, providing a superior approach and improving the accuracy of inertial sensor analysis.

[0057] The various embodiments in this specification are described in a progressive manner. Similar or identical parts between embodiments can be referred to mutually. Each embodiment focuses on describing the differences from other embodiments. In particular, the embodiments of apparatus, devices, and non-volatile computer storage media are basically similar to the method embodiments, so the descriptions are relatively simple; relevant parts can be referred to the descriptions of the method embodiments.

[0058] The foregoing has described specific embodiments of this specification. Other embodiments are within the scope of the appended documents. In some cases, the actions or steps described in this application may be performed in a different order than that shown in the embodiments and still achieve the desired results. Furthermore, the processes depicted in the drawings do not necessarily require the specific or sequential order shown to achieve the desired results. In some embodiments, multitasking and parallel processing are also possible or may be advantageous.

[0059] The above content is merely an example and illustration of the concept of the present invention. Those skilled in the art can make various modifications or additions to the specific embodiments described or use similar methods to replace them, as long as they do not deviate from the concept of the invention or exceed the scope defined in this application, they should all fall within the protection scope of the present invention.

Claims

1. A data analysis method for inertial sensors based on frequency domain processing, characterized in that, The method includes: S1. The original output signal is acquired by the MEMS inertial measurement unit of the sensor; S2. The time-domain signal in the original output signal is decomposed into solid-state mode equations and residual terms using a mode decomposition algorithm. The mode decomposition algorithm includes empirical mode decomposition and variational mode decomposition. The mode decomposition algorithm obtains IMFs through a sieving method. The sieving method first selects the IMFs that meet two properties, including: each consecutive local extremum must be sandwiched between a zero point and the local mean of each IMF must be zero. The empirical mode decomposition uses a numerical iterative method to extract the frequency-banded portion of the original signal, thereby extracting the solid-state mode equations. The variational mode decomposition utilizes signal characteristics and prior knowledge to select the center frequency, extracts the solid-state mode equation at the center frequency, and models it as a function of amplitude, phase, and frequency. S3. Perform Hilbert transform on the solid-state mode equation components to generate a Hilbert spectrum with joint time-frequency domain distribution; S4. Construct a sensor performance evaluation model based on the Hilbert spectrum; The implementation steps of the empirical mode decomposition are as follows: S21. Find the maximum and minimum values ​​of each interval using cubic spline interpolation, classify the maximum and minimum values, construct upper and lower envelopes, and obtain the cubic spline interpolation coefficients of all maximum and minimum values. S22. Perform LU decomposition using the Doolittle algorithm to obtain an upper triangular matrix U and a lower triangular matrix L, satisfying A=LU. Initialize L as an identity matrix of dimension n, and U as an empty matrix of dimension n. The calculation formula is: , ; where the subscripts i, j, and k all represent iteration variables and satisfy: Perform a forward replacement; Perform backward calculation; substitute the values ​​that satisfy the condition. A system of linear equations of the form x is used to obtain the final result x. S23. Calculate the average value of the upper and lower envelopes of a continuous interval. The algorithm determines whether the obtained function is an IMF. If the relative energy ratio is less than a preset threshold or the number of iterations exceeds the preset maximum number of iterations, the iteration stops, and the average signal is subtracted from the residual signal to obtain an IMF. S24. Continue to repeat step S23 iterative screening until the external loop iteration stop condition is met; the external loop stop condition is that the number of solid mode equations reaches the preset maximum number, or the energy ratio is less than the set ratio, or the number of poles is less than the preset number. The implementation steps of the variational mode decomposition are as follows: select the center frequency using signal characteristics and prior knowledge, extract the solid-state mode equation at the center frequency, and model it as a function of amplitude, phase, and frequency. SS1. Construct the objective function of the constrained variational problem, and apply it to the IMF components. and specific center frequency Solution: ； Among them, It is the k-th IMF component. It is the specific center frequency of the k-th IMF. For the Dirac function, It is the imaginary unit, and ; For time differential operators, for Norm, This represents the change in the frequency domain obtained by modulating the modulo signal u onto a carrier wave; SS2. Iteratively solve for the optimal solution set using the alternating direction multiplier method. ; The evaluation conditions for the sensor performance evaluation model based on the Hilbert spectrum are as follows: When the Hilbert spectrum of EMD-HHT has no spikes and a uniform amplitude distribution at the same amplitude level in the full time domain and the full frequency domain, it is judged to be in normal working condition; otherwise, it is judged to be at risk of failure in working condition. When the Hilbert spectrum of VMD-HHT has a uniform amplitude distribution at a specific frequency, it is considered to be in normal working condition; otherwise, it is considered to be in a faulty working condition.

2. The data analysis method for an inertial sensor based on frequency domain processing according to claim 1, characterized in that, The SS2 includes: Update the IMF frequency domain expression: ;in, It is relative to continuous frequency Discrete frequencies; It is relative to continuous frequency The discrete frequency of the k-th iteration; n is the iteration number; Discretized Fourier transform; Update center frequency: ;in, It is a continuous frequency, and , The frequency domain of the k-th IMF component serves as the input source for generating the Hilbert spectrum; The center frequency of the k-th IMF component; Update the Lagrange multipliers: ;in, Update the step size for the multiplier; The conjugate of the frequency-domain Lagrange multiplier serves as a constraint force adjustment factor. To augment the Lagrange function For conjugate Lagrange multipliers The gradient; The frequency domain of the original output signal is used as the Fourier transform of the sensor output. These are the Lagrange multipliers before the update; For the updated Lagrange multipliers; and Discretization and subsequent input yield the following results: .

3. The data analysis method for an inertial sensor based on frequency domain processing according to claim 1, characterized in that, The Hilbert transform includes: The continuous Hilbert transform is defined as: ;in, To update the step size for the multiplier, For the first IMF frequency components of a step size; This represents the result of a 90° phase shift in all frequency components of the original signal after the Hilbert transform. Constructing analytical signals: ,in, The original signal, To analyze the signal, The imaginary unit; Transform polar coordinates: ,in It is the instantaneous amplitude. It is the instantaneous frequency.

4. The data analysis method for an inertial sensor based on frequency domain processing according to claim 1, characterized in that, S4 also includes using numerical methods to compress the Hilbert spectrum: Rational number interpolation is constructed based on Chebyshev's alternation theorem: ; According to the original function: Error equations obtained through rational number interpolation ;in, , , Representing elements respectively The natural number to the power of, , , They represent , , polynomials; The extreme points of the error are calculated by Remez iteration and interpolation is continued to satisfy the Shchev alternation theorem. After convergence, the optimal solution is obtained.

Images