[0033] In order to make the technical means, creative features, goals and effects achieved by the present invention easy to understand, the present invention will be further described below in conjunction with specific embodiments.
[0034] refer to figure 1, this specific example adopts the following technical scheme: a kind of adaptive denoising method based on wavelet analysis, aiming at molecular liquid ring type angular accelerometer signal, its steps are: 1, carry out wavelet layered energy entropy calculation to angular accelerometer signal , to obtain the optimal number of layers M of wavelet decomposition, and perform M-layer discrete wavelet decomposition on the signal to obtain the wavelet coefficients on each decomposition scale; 2. Calculate the wavelet decomposition threshold threshold by coupling with the median filter algorithm, and use the threshold to threshold the wavelet coefficients Quantization processing; 3. Finally, use the processed wavelet coefficients to perform one-dimensional wavelet inverse transformation to complete the denoising processing.
[0035] 1. Adaptive determination of wavelet decomposition layers
[0036] In traditional engineering applications, the following method is usually used to determine the number of discrete wavelet decomposition layers: by performing m=1 layer discrete wavelet decomposition on the noisy signal, and whitening the detail coefficients obtained after decomposition, if the details of the layer If the coefficient sequence is a white noise sequence, then the number of decomposition layers is m+1. When the detail coefficients fail the whitening test, the number of discrete wavelet decomposition layers is M=m-1. However, when the angular accelerometer is applied in a multi-disturbance environment, its real signal will also show high dynamic characteristics that are correlated with white noise. At this time, it is impossible to distinguish the real signal simply by whitening the decomposed detail coefficients. with noise signal.
[0037] Therefore, according to the maximum entropy principle proposed by E.T.Jaynes in 1957, the present invention proposes a wavelet energy entropy algorithm to realize the adaptive determination of wavelet decomposition layers. This method is based on the difference between the noise signal and the real signal energy distribution to determine Decompose layers. When the angular accelerometer is used in a multi-disturbance environment, the detail coefficients of the output signal decomposed will contain the high dynamic part of the real signal, and their energy distribution will be significantly different from the detail coefficients of the white noise signal decomposed. Through Setting the threshold of this difference can determine the number of decomposition layers more reliably and accurately.
[0038] The wavelet energy entropy algorithm starts from the energy of the signal and is related to the complexity of the signal and the number of decomposition layers. The specific algorithm is as follows:
[0039] Define Y(t) to represent the angular accelerometer signal sequence with length n, then its energy is:
[0040]
[0041] Then the detail coefficient of each layer after performing M-level discrete wavelet decomposition on the Y(t) signal is Y i (t), i=1,2,...,M, for the i-th layer wavelet energy entropy WEE(i) and the ratio coefficient P of the detail coefficient energy of each layer to the total detail coefficient energy sum j is given by:
[0042]
[0043]
[0044] The specific operation steps are given below:
[0045] 1) if figure 2 As shown, the angular accelerometer signal needs to be obtained first. What the present invention adopts is the molecular type liquid ring type angular accelerometer developed by the 33rd Research Institute of the Third Research Institute of China Aerospace Science and Industry Corporation, which is installed on an adjustable-speed turntable. Voltage signals, so as to obtain multiple sets of voltage signals output by the angular accelerometer at different frequencies.
[0046] Define the required variables as follows:
[0047] XN={XN i |i=1,2,...,n} represents the original noisy signal sequence whose length is n;
[0048] XN sam ={XN i |i=1,2,...,m,m=min{n,1000}} indicates the sample noisy signal, where min{n,1000} indicates that if the original noisy signal length n is greater than 1000 , then take the first 1000 data of the original noisy signal as the sample sequence, otherwise take all the data of the original noisy signal as the sample sequence;
[0049] X={X i |i=1,2,...,m} represents the signal sequence obtained after normalizing the sample noise signal;
[0050] N sam ={N i |i=1,2,...,m} means to generate a white noise sequence with the same length as the sample noise signal;
[0051] N={N i |i=1,2,...,m} represents the signal sequence obtained after normalizing the white noise sequence;
[0052] S={S i |i=1,2,...,n} represents the real signal sequence;
[0053] Indicates the signal sequence after denoising;
[0054] m max Indicates the maximum number of decomposable layers;
[0055] M={M i |i=1,2,...,M max} represents the actual number of layers for discrete wavelet decomposition of the signal;
[0056] For the discrete wavelet decomposition method, for a signal of length n, the maximum number of decomposition layers is given by the following formula:
[0057] m max = log 2 (length(XN sam ))(4)
[0058] In order to improve the calculation speed and reduce the amount of calculation, it is stipulated that the maximum number of decomposition layers does not exceed 9 layers, then the length of the noisy signal length (XN) < 1024 can be calculated from (4). For the obtained angular accelerometer signal XN, if its length is greater than 1000, then take the first 1000 data as the sample data for calculating the optimal number of decomposition layers; if its length is less than 1000, then take all the data for calculation. From this, the signal sample X used to determine the number of decomposition layers is obtained sam , and its maximum number of decomposition layers M max It is 9 floors.
[0059] 2) Generate a signal sample X with sam White noise sequence samples N of the same length sam , and use the following formula for X sam and N sam Simultaneously normalize:
[0060]
[0061]
[0062] 3) Carry out M=1 layer wavelet decomposition on the signals X and N, and the wavelet energy entropy WEE of X can be calculated by formulas (2) and (3) X and the wavelet energy entropy WEE of N N;
[0063] 4) Calculate the difference quotient of the layered energy entropy of X and N obtained by M layer discrete wavelet decomposition, the formula is as follows:
[0064]
[0065] 5) We set the threshold T to detect the obvious difference between the energy entropy of the signal X and the noise N layers to judge whether the optimal decomposition layer M is obtained. If T<0.1, then the number of decomposition layers M=M+1, skip to step 3) to continue discrete wavelet decomposition; if T≥0.1, stop the operation at this time, and obtain the final optimal number of decomposition layers M=M-1.
[0066] 2. Determination of wavelet decomposition threshold
[0067] In order to better separate the real signal of the angular accelerometer from the noise signal and achieve a better denoising effect, in addition to obtaining a reasonable number of wavelet decomposition layers, the selection of the wavelet decomposition threshold is also very important.
[0068] When separating the noise from the detail coefficients decomposed by wavelet, it is necessary to reasonably estimate the strength of the noise, so as to ensure that the real signal part in the detail coefficients is preserved to the greatest extent. Median filtering is a common nonlinear signal smoothing technology. It is based on the theory of sorting statistics and has a good filtering effect on impulse noise. It can protect the contour of the signal and remove large burrs while effectively suppressing noise. By using the median filter, the contour of the real signal can be quickly and accurately estimated to a certain extent, and then subtracted from the original noisy signal to realize the estimation of the noise intensity.
[0069] Therefore, the present invention defines the concept of the average half-width W of noise, which is the average value of the absolute value of the noise intensity, and the calculation formula is as follows:
[0070]
[0071] in, Represents a sequence of estimated noise strengths of length n.
[0072] Furthermore, we can calculate the statistical characteristic of the standard deviation σ of the noise, and combine the 3σ rule to set the wavelet decomposition threshold threshold. When the angular accelerometer is in a working environment with many disturbances and heavy loads, its output signal will be mixed with some high-dynamic and high-intensity real signals in the detail coefficients after discrete wavelet decomposition. In order to preserve the real signals to the greatest extent, To filter out the noise, after calculating the noise half-width, apply the 3σ rule on this baseline as the threshold threshold, and keep the detail coefficients of "abnormal" intensity beyond the range of white noise intensity, and the rest of the detail coefficients are in the probability statistics. The range basically belongs to the noise signal, and it can be set to zero.
[0073] The specific operation steps for wavelet threshold determination are given below:
[0074] 1) As attached image 3 As shown, define the required variables as follows:
[0075] XN={XN i |i=1,2,...,n} represents a noisy signal sequence with a length of n;
[0076] XM={XM i |i=1,2,...,n} represents the signal sequence obtained after the noisy signal is filtered by the median;
[0077] represents the estimated noise intensity sequence;
[0078] W N Indicates the average half-width of the estimated noise intensity sequence;
[0079] gate L Indicates the upper bound of the threshold threshold, gate H Represents the lower bound of the threshold threshold;
[0080] First, the window size is selected to be 5% of the sequence length n, and the median filter is performed on the noisy signal XN to obtain the filtered sequence XM;
[0081] 2) Calculate the estimated noise intensity sequence, the calculation formula is as follows:
[0082]
[0083] 3) calculate The half-width W N and standard deviation Calculated as follows:
[0084]
[0085]
[0086] 4) Finally, the lower bound and upper bound of the threshold threshold are obtained, and the calculation formula is as follows:
[0087]
[0088]
[0089] So far, the determination of the threshold threshold for discrete wavelet decomposition of the angular accelerometer signal has been completed.
[0090] 3. Threshold processing and one-dimensional wavelet reconstruction
[0091] The present invention can adaptively determine the threshold threshold of discrete wavelet decomposition layers and detail coefficients through the above two-step algorithm. To complete the denoising of the diagonal accelerometer signal, two steps of threshold processing and one-dimensional wavelet inverse transformation are required. Specifically The calculation steps are as follows:
[0092] 1) First, use the db3 wavelet basis function to perform M-level discrete wavelet decomposition on the original noisy signal XN to obtain the approximation coefficient cAM and the detail coefficient cD1, cD2,...cDM;
[0093] 2) According to the upper and lower bounds of the threshold threshold calculated in 2, perform threshold processing on the detail coefficients of each layer, and perform threshold processing on the coefficients between (gate L , gate H ) between the detail coefficients are all set to zero, and the part exceeding the threshold value is reserved, so as to obtain the processed detail coefficients cD1′, cD2′...cDM′;
[0094]
[0095] 3) Use the approximation coefficient cAM and the detail coefficients cD1', cD2'...cDM' of each layer after thresholding to perform one-dimensional wavelet inverse transformation to obtain the final denoising signal.
[0096] 4. Example analysis
[0097] 1) Non-stationary random signal simulation experiment
[0098] The present invention adopts Matlab to compile the angular acceleration signal denoising method program based on wavelet analysis, at first select Blocks, Bumps, Heavysine, Doopler in Matlab to simulate the non-stationary random signal similar to the angular accelerometer signal feature to test, as Figure 4 shown. The sampling points of the four signals are all 2048, and Gaussian white noise is superimposed on the signals to obtain noisy signals with signal-to-noise ratios of 10dB and 15dB, such as Figure 5 , Image 6 shown. Apply the wavelet energy entropy algorithm proposed by the present invention to determine adaptively the decomposition layers of each non-stationary random signal as shown in Table 1-Table 2:
[0099] Table 1 The best adaptive decomposition layers for each noisy non-stationary random signal when the SNR is 10dB
[0100]
[0101] Table 2 The best adaptive decomposition layers for each noisy non-stationary random signal when the signal-to-noise ratio is 15db
[0102]
[0103] From Daubechies wavelet function and Symlets wavelet function in MATLAB wavelet library, db3 wavelet is selected for processing.
[0104] To measure the signal noise reduction effect, the signal-to-noise ratio (SNR) and mean square error (MSE) are usually used. The signal XN(i) after the noise-free signal S(i) is superimposed on the Gaussian white noise signal N(i), XN(i)=S(i)+N(i), i=1, 2, 3... ..n, where n is the number of sampling points. is the signal after denoising, then the calculation formula of signal-to-noise ratio is formula (15):
[0105]
[0106] The formula for calculating the signal-to-noise ratio after noise reduction is as follows:
[0107]
[0108] The mean square error is defined as:
[0109]
[0110] The smaller the root mean square error between the original signal and the denoised signal, and the higher the signal-to-noise ratio of the signal, the closer the denoised signal is to the real original signal, and the better the denoising effect is. Using traditional thresholding methods: heuristic threshold selection methods, generic thresholding The maximum and minimum value threshold method and the threshold processing method proposed by the present invention carry out noise reduction processing to four kinds of non-stationary random random signals respectively, and the results are shown in Table 3 and Table 4. It can be known that the threshold processing method proposed by the present invention can obtain Higher signal-to-noise ratio and smaller mean square error.
[0111] Table 3 Simulation signal noise reduction results (the original signal-to-noise ratio is 10dB)
[0112]
[0113] Table 4 Simulation signal noise reduction results (the original signal-to-noise ratio is 15dB)
[0114]
[0115]
[0116] Using 4 kinds of threshold value processing methods to process the 4 kinds of signals respectively, the results are as follows: Image 6 — Figure 14 As shown, it can be intuitively seen that the wavelet decomposition threshold determination method determined by the present invention has significantly better noise reduction effects on the four types of signals than other traditional threshold functions.
[0117] 2) Molecular type liquid ring angular accelerometer real signal denoising effect
[0118] The present invention adopts the molecular type liquid ring angular accelerometer developed by the 33rd Research Institute of the Third Research Institute of China Aerospace Science and Industry Corporation, and installs it on an adjustable-speed turntable. signal, so as to obtain multiple sets of voltage signals output by the angular accelerometer at different frequencies. Select one of the signal sequences with an output voltage frequency of 5HZ and a signal strength of 940.5mV to denoise it. The real original signal is as follows: Figure 15 As shown, using the method proposed by the present invention, the optimal number of wavelet decomposition layers is calculated to be 4 layers, and compared with the traditional threshold selection method, the effect after denoising is as follows Figure 16 shown.
[0119] To sum up, the above are only preferred embodiments of the present invention, and are not intended to limit the protection scope of the present invention. Any modifications, equivalent replacements, improvements, etc. made within the spirit and principles of the present invention shall be included within the protection scope of the present invention.
