Adaptive time-frequency supported frequency modulation signal decomposition method and system

Abstract

The application provides a self-adapting time-frequency supported frequency modulation signal decomposition method and system, which comprises the following steps: calculating a parameterized time-frequency transform of a signal; extracting a rough instantaneous frequency ridge curve of a sub-signal with maximum energy from a time-frequency representation of the parameterized time-frequency transform of the signal; taking the extracted instantaneous frequency curve of the sub-signal as an initial instantaneous frequency of signal decomposition; calculating an instantaneous frequency kernel function matrix of the sub-signal according to the instantaneous frequency of the sub-signal; calculating a demodulation signal of the sub-signal; calculating an initial instantaneous frequency increment of the sub-signal; calculating an optimized initial instantaneous frequency increment of the sub-signal; and updating the instantaneous frequency of the sub-signal. The application realizes the decomposition of a complex multi-component signal with strong noise, can accurately estimate the instantaneous frequency (IF) of each sub-signal and reconstruct each sub-signal, effectively reveals the change rule of the original complex signal, can be applied to the actual engineering signal processing in many fields, and has strong engineering practical value.

Description

Technical Field

[0001] This invention relates to the technical field of engineering signal processing, and more specifically, to an adaptive time-frequency supported frequency modulation signal decomposition method and system. Background Technology

[0002] In practical industrial applications, there are numerous amplitude-frequency modulation (AM-FM) signals whose instantaneous frequency and amplitude vary over time. Examples include vibration signals from rotating machinery operating at varying speeds, spindle vibration signals from cutting equipment, sound signals, and radar detection signals. Many practically measured AM-FM signals are inevitably affected by noise, especially under harsh operating conditions, where they can be severely disrupted by strong noise. Strong noise significantly interferes with the estimation of the signal's instantaneous frequency (IF) and instantaneous amplitude (IA). Existing signal decomposition methods cannot accurately decompose complex AM-FM signals with strong noise, and are prone to mode aliasing problems under strong noise interference. Decomposing complex, multi-component AM-FM signals with strong noise can effectively estimate the instantaneous frequency of each sub-signal and reconstruct each sub-signal, thereby effectively filtering out strong noise interference, revealing the signal's time-frequency variation patterns, and effectively understanding the key information contained in the signal.

[0003] Regarding the aforementioned related technologies, the inventors believe that the signal decomposition methods have shortcomings. Therefore, a new technical solution is needed to improve upon these technical problems. Summary of the Invention

[0004] In view of the shortcomings of the prior art, the purpose of this invention is to provide an adaptive time-frequency supported frequency modulation signal decomposition method and system.

[0005] According to the present invention, an adaptive time-frequency supported frequency modulation signal decomposition method is provided, the method comprising the following steps:

[0006] Step S1: Calculate the parameterized time-frequency transform (PTFT) of signal y;

[0007] Step S2: Extract the sub-signal y with the highest energy from the time-frequency representation (TFR) of the parameterized time-frequency transform of signal y. m The rough instantaneous frequency ridge curve IF m ;

[0008] Step S3: Extract the sub-signal y m Instantaneous frequency curve IF m The initial instantaneous frequency f for signal decomposition m ;

[0009] Step S4: Based on the sub-signal y m instantaneous frequency f mcalculating the instantaneous frequency kernel function matrix Φ of the sub-signal m ;

[0010] Step S5: calculating the demodulation signal u m and v m of the sub-signal y m ;

[0011] Step S6: calculating the initial instantaneous frequency increment Δf m of the sub-signal y

[0012] Step S7: calculating the optimized initial instantaneous frequency increment Δf m of the sub-signal y m ;

[0013] Step S8: updating the instantaneous frequency f m of the sub-signal y m ;

[0014] Step S9: reconstructing the sub-signal y m ;

[0015] Step S10: repeating the steps S3-S9 until the reconstructed sub-signal y m satisfies the condition of terminating iteration;

[0016] Step S11: fitting the sub-signal y m with the estimated instantaneous frequency f m using the K-order Fourier series;

[0017] Step S12: updating the kernel function parameters of the parametric time-frequency transform of the signal y and performing the parametric time-frequency transform of the signal y again;

[0018] Step S13: repeating the steps S1-S12 until the instantaneous frequency f m of the sub-signal y m satisfies the condition of terminating iteration;

[0019] Step S14: reconstructing the sub-signal y m ;

[0020] Step S15: removing the sub-signal y m from the current residual signal r m , and then continuing to decompose the next sub-signal y m+1 and repeating the steps S1-S14 until each sub-signal in the original signal is completely decomposed.

[0021] Preferably, the kernel function parameters of the PTFT in the step S1 are Fourier series, and the calculation formula is:

[0022]

[0023] in The calculation formula is:

[0024]

[0025] In the formula, K is the order of the Fourier series, f0 is the fundamental frequency, f0 = 1 / (2T), and T is the signal duration. The coefficients of the Fourier series; w h The calculation formula is:

[0026]

[0027] In step S2, E(IF) is solved. m The local maximum value of ) yields the sub-signal y. m Instantaneous frequency ridge curve IF m E(IF) m The formula for calculating ) is:

[0028]

[0029] In the formula S(t,IF m (t) represents the PTFT of the signal, and λ and β are penalty coefficients;

[0030] The instantaneous frequency kernel function matrix Φ in step S4 m The calculation formula is:

[0031] Φ m =[C m ,S m ]

[0032] In the formula C m Let S represent the sine function matrix. m Represents the cosine function matrix;

[0033] C m The calculation formula is:

[0034] C m =diag[cos(θ) m (t1))…cos(θ m (t N ))]

[0035] In the formula θ m represents the sub-signal y m The instantaneous phase;

[0036] Instantaneous phase θ m The calculation formula is:

[0037]

[0038] wherein denotes the demodulation signal of the sub-signal y m the estimated instantaneous frequency;

[0039] S m is calculated by:

[0040] S m = diag[sin(θ m (t1))…sin(θ m (t N ))]

[0041] wherein θ m denotes the instantaneous phase of the sub-signal y m .

[0042] Preferably, the demodulation signal u m and v m of the sub-signal y m in the step S5 are calculated by:

[0043]

[0044] wherein n denotes the iteration number, Φ m denotes the instantaneous frequency kernel function matrix, ρ denotes the penalty coefficient, and D denotes the diagonal matrix of the second-order difference matrix Ω r m denotes the current residual signal;

[0045] The expression of the second-order difference matrix Ω is:

[0046]

[0047] Ω is an N-row and N-column matrix, wherein N denotes the number of discrete sampling points of the signal;

[0048] The calculation formula of the initial instantaneous frequency increment Δf of the sub-signal y m in the step S6 is:

[0049]

[0050] wherein and denote the demodulation signal of the sub-signal y m at the n-th iteration;

[0051] The calculation formula of the optimized initial instantaneous frequency increment Δf m of the sub-signal y m in the step S7 is:

[0052]

[0053] In the formula, Ω represents the second-order difference matrix, and I represents the identity matrix.

[0054] Preferably, in step S8, the sub-signal y m instantaneous frequency f m The update formula is:

[0055]

[0056] In the formula, n represents the number of iterations;

[0057] In step S9, the sub-signal y m The reconstruction calculation formula is:

[0058] y m =Φ m x m

[0059] In the formula Φ m represents the sub-signal y m The instantaneous frequency kernel function matrix, x m represents the sub-signal y m The demodulated signal;

[0060] In step S10, the initial iteration count n is set to 1. After one iteration is completed, the iteration count n is incremented by 1; the reconstructed sub-signal y m The termination condition for the iteration is:

[0061]

[0062] In the formula, ε1 represents the termination iteration threshold; the reconstructed sub-signal y m The termination iteration threshold ε1 is set to 10. -8 .

[0063] Preferably, in step S11, the sub-signal y is used for fitting. m instantaneous frequency f m The formula for the Kth order Fourier series is:

[0064]

[0065] In the formula The coefficients of the Fourier series;

[0066] In step S13, the sub-signal y m instantaneous frequency f m The termination condition for the iteration is:

[0067]

[0068] Wherein h is the iteration number, and ε2 represents the termination iteration threshold; the termination iteration threshold ε2 is set to 10 -8 ;

[0069] The reconstruction calculation formula of the sub-signal y m in the step S14 is:

[0070] y m = Φ m x m

[0071] Wherein Φ m represents the instantaneous frequency kernel function matrix of the sub-signal y m , and x m represents the demodulation signal of the sub-signal y m ;

[0072] In the step S15, m is set to 1 at the beginning, and the current residual signal r1 is equal to the original signal y; after the decomposition of the mth sub-signal is completed, the sub-signal y m is removed, at this time, the current residual signal r m+1 =r m -y m , and m is set to m+1.

[0073] The application further provides an adaptive time-frequency support frequency modulation signal decomposition system, which comprises the following modules:

[0074] Module M1: calculating the parameterized time-frequency transform PTFT of the signal y;

[0075] Module M2: extracting the rough instantaneous frequency ridge curve IF m of the sub-signal y m with the maximum energy from the time-frequency representation TFR of the parameterized time-frequency transform of the signal y;

[0076] Module M3: taking the instantaneous frequency curve IF m of the extracted sub-signal y m as the initial instantaneous frequency f m of the signal decomposition;

[0077] Module M4: calculating the instantaneous frequency kernel function matrix Φ m of the sub-signal according to the instantaneous frequency f m of the sub-signal y m ;

[0078] Module M5: calculating the demodulation signal u m and v m of the sub-signal y m ;

[0079] Module M6: calculating the sub-signal y mthe initial instantaneous frequency increment of the sub-signal y

[0080] Module M7: calculating the optimized initial instantaneous frequency increment of the sub-signal y m

[0081] Module M8: updating the instantaneous frequency f of the sub-signal y m m

[0082] Module M9: reconstructing the sub-signal y m

[0083] Module M10: repeatedly calling the modules M3-M9 until the reconstructed sub-signal y m satisfies the condition of terminating iteration

[0084] Module M11: fitting the sub-signal y with a K-order Fourier series to estimate the instantaneous frequency f m m

[0085] Module M12: updating the kernel function parameters of the parametric time-frequency transform of the signal y and performing the parametric time-frequency transform of the signal y again

[0086] Module M13: repeatedly calling the modules M1-M12 until the instantaneous frequency f of the sub-signal y m m satisfies the condition of terminating iteration

[0087] Module M14: reconstructing the sub-signal y m

[0088] Module M15: removing the sub-signal y from the current residual signal r and then continuing to decompose the next sub-signal y and repeatedly calling the modules M1-M14 until each sub-signal in the original signal is completely decomposed m m m+1

[0089] Preferably, the kernel function parameters of the PTFT in the module M1 are a Fourier series, and the calculation formula is:

[0090]

[0091] wherein the calculation formula of is:

[0092]

[0093] wherein K is the order of the Fourier series, f0 is the fundamental frequency, f0 = 1 / (2T), T is the signal duration, and ​​​​​​​​​​​​ The coefficients of the Fourier series; w h The calculation formula is:

[0094]

[0095] In module M2, E(IF) is solved. m The local maximum value of ) yields the sub-signal y. m Instantaneous frequency ridge curve IF m E(IF) m The formula for calculating ) is:

[0096]

[0097] In the formula S(t,IF m (t) represents the PTFT of the signal, and λ and β are penalty coefficients;

[0098] The instantaneous frequency kernel function matrix Φ in module M4 m The calculation formula is:

[0099] Φ m =[C m ,S m ]

[0100] In the formula C m Let S represent the sine function matrix. m Represents the cosine function matrix;

[0101] C m The calculation formula is:

[0102] C m =diag[cos(θ) m (t1))…cos(θ m (t N ))]

[0103] In the formula θ m represents the sub-signal y m The instantaneous phase;

[0104] Instantaneous phase θ m The calculation formula is:

[0105]

[0106] In the formula represents the sub-signal y m Estimated instantaneous frequency;

[0107] S m The calculation formula is:

[0108] S m= diag [sin (θ m (t1))…sin(θ m (t N ))

[0109] where θ m represents the instantaneous phase of the sub-signal y m .

[0110] Preferably, the calculation formula of the demodulation signals u m and v m of the sub-signal y m in the module M5 is as follows:

[0111]

[0112] where n represents the iteration number, Φ m represents the instantaneous frequency kernel function matrix, ρ represents the penalty coefficient, and D represents the diagonal matrix of the second-order difference matrix Ω r m represents the current residual signal;

[0113] The expression of the second-order difference matrix Ω is as follows:

[0114]

[0115] Ω is an N-row and N-column matrix, and N represents the number of discrete sampling points of the signal;

[0116] The calculation formula of the initial instantaneous frequency increment Δf m of the sub-signal y in the module M6 is as follows:

[0117]

[0118] where u and v represent the demodulation signals of the sub-signal y m at the n-th iteration;

[0119] The calculation formula of the optimized initial instantaneous frequency increment Δf m of the sub-signal y m in the module M7 is as follows:

[0120]

[0121] where Ω represents the second-order difference matrix, and I represents the unit matrix.

[0122] Preferably, the update formula of the instantaneous frequency f m of the sub-signal y m in the module M8 is as follows:

[0123]

[0124] wherein n represents the iteration number;

[0125] The reconstruction calculation formula of the sub-signal y m in the module M9 is:

[0126] y m = Φ m x m

[0127] wherein Φ m represents the instantaneous frequency kernel function matrix of the sub-signal y m , and x m represents the demodulation signal of the sub-signal y m ;

[0128] In the module M10, the iteration number n is initially set as 1, and after completing one iteration, the iteration number n is increased by 1; and the termination iteration condition of the reconstructed sub-signal y m is:

[0129]

[0130] wherein ε1 represents the termination iteration threshold; and the termination iteration threshold ε1 of the reconstructed sub-signal y m is set as 10 -8 .

[0131] Preferably, the formula for fitting the K-order Fourier series of the instantaneous frequency f m of the sub-signal y m in the module M11 is:

[0132]

[0133] wherein a is the coefficient of the Fourier series;

[0134] In the module M13, the termination iteration condition of the instantaneous frequency f m of the sub-signal y m is:

[0135]

[0136] wherein h is the iteration number, and ε2 represents the termination iteration threshold; and the termination iteration threshold ε2 is set as 10 -8 .

[0137] The reconstruction calculation formula of the sub-signal y m in the module M14 is:

[0138] y m = Φ m xm

[0139] where Φ m represents the instantaneous frequency kernel matrix of the sub-signal y m ; x m represents the demodulation signal of the sub-signal y m .

[0140] The module M15 is initially set as m=1, and the current residual signal r1 is equal to the original signal y; after the decomposition of the mth sub-signal is completed, the sub-signal y m is removed, and the current residual signal r m+1 = r m -y m is set as m=m+1.

[0141] Compared with the prior art, the adaptive time-frequency supported FM signal decomposition method and system has the following beneficial effects:

[0142] 1. The adaptive time-frequency supported FM signal decomposition method and system combines the IF arctangent demodulation optimization and the iterative optimization of the time-frequency representation of the signal, uses the time-frequency representation information of the signal to assist the iteration of the IF of each sub-signal, can effectively suppress the interference of strong noise on the signal decomposition, and further accurately estimates the IF of the signal.

[0143] 2. The adaptive time-frequency supported FM signal decomposition method and system can effectively decompose the AM-FM signal with strong noise, accurately estimate the IF of the signal and reconstruct the sub-signal. BRIEF DESCRIPTION OF DRAWINGS

[0144] Other characteristics, objects and advantages of the present application will become more apparent from the following detailed description of non-restrictive embodiments, made with reference to the attached drawings:

[0145] Figure 1 is a flow chart of the adaptive time-frequency supported FM signal decomposition method of the present application;

[0146] Figure 2 is a time domain curve diagram of the simulation signal with strong noise of the present application;

[0147] Figure 3 is a time-frequency representation diagram of the short-time Fourier transform of the simulation signal with strong noise of the present application;

[0148] Figure 4 is an actual instantaneous frequency diagram of the simulation signal of the present application;

[0149] Figure 5 is an actual time domain curve diagram of each sub-signal of the simulation signal of the present application;

[0150] Figure 6The IF curve of the simulated signal estimated by the adaptive time-frequency supported frequency modulation signal decomposition method of the present invention is shown in the figure.

[0151] Figure 7 The time-domain curves of the two sub-signals of the simulated signal reconstructed by the adaptive time-frequency supported frequency modulation signal decomposition method of the present invention are shown.

[0152] Figure 8 This is a time-domain curve of the bearing housing vibration signal when the bearing with an outer ring fault due to noise interference is running at variable speed.

[0153] Figure 9 This is a time-domain curve of the short-time Fourier transform of the bearing housing vibration signal of the present invention;

[0154] Figure 10 The time-frequency representation diagram is obtained by performing signal decomposition on the bearing housing vibration signal using the adaptive time-frequency supported frequency modulation signal decomposition method of the present invention. Detailed Implementation

[0155] The present invention will now be described in detail with reference to specific embodiments. These embodiments will help those skilled in the art to further understand the present invention, but do not limit the invention in any way. It should be noted that those skilled in the art can make several changes and improvements without departing from the concept of the present invention. These all fall within the scope of protection of the present invention.

[0156] Example 1

[0157] According to the present invention, an adaptive time-frequency supported frequency modulation signal decomposition method is provided, the method comprising the following steps:

[0158] Step S1: Calculate the parameterized time-frequency transform (PTFT) of signal y; the kernel function parameter of PTFT is a Fourier series, and its calculation formula is:

[0159]

[0160] in The calculation formula is:

[0161]

[0162] In the formula, K is the order of the Fourier series, f0 is the fundamental frequency, f0 = 1 / (2T), and T is the signal duration. The coefficients of the Fourier series; w h The calculation formula is:

[0163]

[0164] Step S2: extracting a sub-signal y from the time-frequency representation TFR of the parametric time-frequency transform of the signal y with the largest energy m The rough instantaneous frequency ridge curve IF of the sub-signal y m The rough instantaneous frequency ridge curve IF of the sub-signal y m is obtained by solving the local maximum of E(IF m ) m The calculation formula of E(IF m ) is as follows:

[0165]

[0166] Where S(t, IF m (t)) represents the PTFT of the signal, and λ and β are penalty coefficients.

[0167] Step S3: taking the instantaneous frequency curve IF m of the sub-signal y m as the initial instantaneous frequency f m of the signal decomposition;

[0168] Step S4: calculating the instantaneous frequency kernel function matrix Φ m of the sub-signal y m according to the instantaneous frequency f m of the sub-signal; the calculation formula of the instantaneous frequency kernel function matrix Φ m is as follows:

[0169] Φ m = [C m , S m ]

[0170] Where C m represents a sine function matrix, and S m represents a cosine function matrix;

[0171] The calculation formula of C m is as follows:

[0172] C m = diag [cos(θ m (t1))…cos(θ m (t N ))]

[0173] Where θ m represents the instantaneous phase of the sub-signal y m ;

[0174] The calculation formula of the instantaneous phase θ m is as follows:

[0175]

[0176] In the formula The demodulation signal u m The estimated instantaneous frequency

[0177] S m The calculation formula of the demodulation signal u

[0178] S m = diag[sin(θ m (t1))…sin(θ m (t N ))]

[0179] In the formula m The instantaneous phase of the sub-signal y m .

[0180] Step S5: Calculate the demodulation signal u m and v m of the sub-signal y m ; The calculation formula of the demodulation signal u m and v m of the sub-signal y m is as follows:

[0181]

[0182] In the formula, n represents the number of iterations, Φ m represents the instantaneous frequency kernel function matrix, ρ represents the penalty coefficient, and D represents the diagonal matrix of the second-order difference matrix Ω r m represents the current residual signal;

[0183] The expression of the second-order difference matrix Ω is as follows:

[0184]

[0185] Ω is an N-row and N-column matrix, and N represents the number of discrete sampling points of the signal.

[0186] Step S6: Calculate the initial instantaneous frequency increment m of the sub-signal y The calculation formula of the initial instantaneous frequency increment m of the sub-signal y is as follows:

[0187]

[0188] In the formula and represent the demodulation signal of the sub-signal y m at the n-th iteration.

[0189] Step S7: Calculate the demodulation signal u mThe optimized initial instantaneous frequency increment Δf m ;Sub signal y m The optimized initial instantaneous frequency increment Δf m The calculation formula is:

[0190]

[0191] In the formula, Ω represents the second-order difference matrix, and I represents the identity matrix.

[0192] Step S8: Update the sub-signal y m instantaneous frequency f m ;Sub signal y m instantaneous frequency f m The update formula is:

[0193]

[0194] In the formula, n represents the number of iterations.

[0195] Step S9: Reconstruct the sub-signal y m ;Sub signal y m The reconstruction calculation formula is:

[0196] y m =Φ m x m

[0197] In the formula Φ m represents the sub-signal y m The instantaneous frequency kernel function matrix, x m represents the sub-signal y m The demodulated signal.

[0198] Step S10: Repeat steps S3-S9 until the reconstructed sub-signal y is obtained. m The condition for terminating the iteration is met; initially, the iteration count n is set to 1, and after each iteration, the iteration count n is incremented by 1; the reconstructed sub-signal y m The termination condition for the iteration is:

[0199]

[0200] In the formula, ε1 represents the termination iteration threshold; the reconstructed sub-signal y m The termination iteration threshold ε1 is set to 10. -8 .

[0201] Step S11: Fit the sub-signal y using a K-order Fourier series. m Estimated instantaneous frequency f m Used for fitting sub-signal y m instantaneous frequency f mThe formula of the K-order Fourier series of y is:

[0202]

[0203] wherein is the coefficient of the Fourier series.

[0204] Step S12: updating the kernel function parameter of the parametric time-frequency transform of the signal y, and performing the parametric time-frequency transform of the signal y again;

[0205] Step S13: repeating the step S1-Step S12 until the instantaneous frequency f m of the sub-signal y m satisfies the termination iteration condition; the instantaneous frequency f m of the sub-signal y m satisfies the termination iteration condition; the termination iteration condition of the instantaneous frequency f

[0206]

[0207] wherein h is the iteration number, and ε2 represents the termination iteration threshold; the termination iteration threshold ε2 is set to 10 -8 .

[0208] Step S14: reconstructing the sub-signal y m ; the reconstruction calculation formula of the sub-signal y m is:

[0209] y m = Φ m x m

[0210] wherein Φ m represents the instantaneous frequency kernel function matrix of the sub-signal y m , and x m represents the demodulation signal of the sub-signal y m .

[0211] Step S15: removing the sub-signal y m from the current residual signal r m , and then continuing to decompose the next sub-signal y m+1 , and repeating the step S1-Step S14 until each sub-signal in the original signal is completely decomposed; the step S15 is first set as m=1, and the current residual signal r1 is equal to the original signal y; after the decomposition of the mth sub-signal is completed, the sub-signal y m is removed, at this time, the current residual signal r m+1 =r m -y m , and m is set as m=m+1.

[0212] The present invention also provides an adaptive time-frequency supported frequency modulation signal decomposition system, which can be implemented by executing the process steps of the adaptive time-frequency supported frequency modulation signal decomposition method. That is, those skilled in the art can understand the adaptive time-frequency supported frequency modulation signal decomposition method as a preferred embodiment of the adaptive time-frequency supported frequency modulation signal decomposition system.

[0213] Example 2

[0214] The present invention also provides an adaptive time-frequency supported frequency modulation signal decomposition system, the system comprising the following modules:

[0215] Module M1: Calculates the parameterized time-frequency transform (PTFT) of signal y; the kernel function parameter of the PTFT is a Fourier series, and its calculation formula is:

[0216]

[0217] in The calculation formula is:

[0218]

[0219] In the formula, K is the order of the Fourier series, f0 is the fundamental frequency, f0 = 1 / (2T), and T is the signal duration. The coefficients of the Fourier series; w h The calculation formula is:

[0220]

[0221] Module M2: Extracts the sub-signal y with the highest energy from the time-frequency representation (TFR) of the parameterized time-frequency transform of the signal y. m The rough instantaneous frequency ridge curve IF m By solving E(IF) m The local maximum value of ) yields the sub-signal y. m Instantaneous frequency ridge curve IF m E(IF) m The formula for calculating ) is:

[0222]

[0223] In the formula S(t,IF m (t) represents the PTFT of the signal, and λ and β are penalty coefficients.

[0224] Module M3: Extracts the sub-signal y m Instantaneous frequency curve IF m The initial instantaneous frequency f for signal decompositionm ;

[0225] Module M4: Based on the sub-signal y m instantaneous frequency f m Calculate the instantaneous frequency kernel function matrix Φ of the sub-signal. m Instantaneous frequency kernel function matrix Φ m The calculation formula is:

[0226] Φ m =[C m ,S m ]

[0227] In the formula C m Let S represent the sine function matrix. m Represents the cosine function matrix;

[0228] C m The calculation formula is:

[0229] C m =diag[cos(θ) m (t1))…cos(θ m (t N ))]

[0230] In the formula θ m represents the sub-signal y m The instantaneous phase;

[0231] Instantaneous phase θ m The calculation formula is:

[0232]

[0233] In the formula represents the sub-signal y m Estimated instantaneous frequency;

[0234] S m The calculation formula is:

[0235] S m =diag[sin(θ)] m (t1))…sin(θ m (t N ))]

[0236] In the formula θ m represents the sub-signal y m The instantaneous phase.

[0237] Module M5: Calculates the sub-signal y m demodulated signal u m and v m ;Sub signal y m demodulated signal um and v m The calculation formula of is:

[0238]

[0239] where n represents the iteration number, Φ m represents the instantaneous frequency kernel function matrix, p represents the penalty coefficient, and D represents the diagonal matrix of the second-order difference matrix Ω r m represents the current residual signal;

[0240] The expression of the second-order difference matrix Ω is:

[0241]

[0242] Ω is an N-row and N-column matrix, and N represents the number of discrete sampling points of the signal.

[0243] Module M6: calculating the initial instantaneous frequency increment Δf m of the sub-signal y The initial instantaneous frequency increment Δf m of the sub-signal y The calculation formula of the initial instantaneous frequency increment Δf

[0244]

[0245] where and represent the demodulation signal of the sub-signal y m at the n-th iteration.

[0246] Module M7: calculating the optimized initial instantaneous frequency increment Δf m of the sub-signal y m ; and the optimized initial instantaneous frequency increment Δf m of the sub-signal y m The calculation formula of the optimized initial instantaneous frequency increment Δf

[0247]

[0248] where Ω represents the second-order difference matrix, and I represents the unit matrix.

[0249] Module M8: updating the instantaneous frequency f m of the sub-signal y m ; and the updating formula of the instantaneous frequency f m of the sub-signal y m

[0250]

[0251] where n represents the iteration number.

[0252] ​Module M9: reconstructing the sub-signal y m ; the reconstruction formula of the sub-signal y m is:

[0253] y m = Φ m x m

[0254] where Φ m represents the instantaneous frequency kernel function matrix of the sub-signal y m , and x m represents the demodulated signal of the sub-signal y m .

[0255] Module M10: repeatedly calling the modules M3-M9 until the reconstructed sub-signal y m satisfies the termination iteration condition; the iteration number n is initially set as 1, and the iteration number n is increased by 1 after one iteration is completed; the termination iteration condition of the reconstructed sub-signal y m is:

[0256]

[0257] where ε1represents the termination iteration threshold; the termination iteration threshold ε1of the reconstructed sub-signal y m is set as 10 -8 .

[0258] Module M11: estimating the instantaneous frequency f m of the sub-signal y m by using the K-order Fourier series fitting; the formula of the K-order Fourier series fitting the instantaneous frequency f m of the sub-signal y m is:

[0259]

[0260] where a is the coefficient of the Fourier series.

[0261] Module M12: updating the kernel function parameters of the parameterized time-frequency transform of the signal y, and performing the parameterized time-frequency transform on the signal y again;

[0262] Module M13: repeatedly calling the modules M1-M12 until the instantaneous frequency f m of the sub-signal y m satisfies the termination iteration condition; the termination iteration condition of the instantaneous frequency f m of the sub-signal y m is:

[0263]

[0264] where h is the iteration number, and ε2 represents the termination iteration threshold; the termination iteration threshold ε2 is set to 10 -8 .

[0265] Module M14: reconstructing the sub-signal y m ; the reconstruction calculation formula of the sub-signal y m is:

[0266] y m = Φ m x m

[0267] where Φ m represents the instantaneous frequency kernel function matrix of the sub-signal y m , and x m represents the demodulation signal of the sub-signal y m .

[0268] Module M15: removing the sub-signal y m from the current residual signal r m , and then continuing to decompose the next sub-signal y m+1 , and repeatedly calling the modules M1-M14 until each sub-signal in the original signal is completely decomposed; the module M15 is initially set as m = 1, and the current residual signal r1 is equal to the original signal y; after the decomposition of the mth sub-signal is completed, the sub-signal y m is removed, at this time the current residual signal r m+1 = r m - y m , and m = m + 1 is set.

[0269] Example 3

[0270] In order to overcome the shortcomings of the existing signal decomposition methods, more effectively decompose complex AM-FM signals with strong noise, and accurately estimate the instantaneous frequency of each sub-signal and reconstruct the sub-signal, we invented an adaptive time-frequency supported frequency modulation signal decomposition method and system. The invented signal decomposition method optimizes the iteration process of the instantaneous frequency of each sub-signal, uses the time-frequency representation signal of the signal to assist the iteration of the instantaneous frequency of the signal, and can effectively filter out the interference of strong noise.

[0271] An adaptive time-frequency supported frequency modulation signal decomposition method, the steps of which comprise:

[0272] Step S1: calculating the parameterized time-frequency transform (PTFT) of the signal y.

[0273] Step S2: extracting the rough instantaneous frequency ridge curve IFm of the sub-signal ym with the maximum energy from the time-frequency representation (TFR) of the parameterized time-frequency transform of the signal y.

[0274] Step S3: The instantaneous frequency curve IFm of the extracted sub-signal ym is taken as the initial instantaneous frequency fm of the signal decomposition.

[0275] Step S4: The instantaneous frequency kernel function matrix Φm of the sub-signal ym is calculated according to the instantaneous frequency fm of the sub-signal ym.

[0276] Step S5: The demodulation signals um and vm of the sub-signal ym are calculated.

[0277] Step S6: The initial instantaneous frequency increment Δf of the sub-signal ym is calculated.

[0278] Step S7: The optimized initial instantaneous frequency increment Δf of the sub-signal ym is calculated. m .

[0279] Step S8: The instantaneous frequency fm of the sub-signal ym is updated.

[0280] Step S9: The sub-signal ym is reconstructed.

[0281] Step S10: Steps S3 to S9 are repeated until the reconstructed sub-signal ym meets the condition for terminating iteration.

[0282] Step S11: The estimated instantaneous frequency fm of the sub-signal ym is fitted with a K-order Fourier series.

[0283] Step S12: The kernel function parameters of the parametric time-frequency transform of the signal y are updated, and the parametric time-frequency transform of the signal y is performed again.

[0284] Step S13: Steps S1 to S12 are repeated until the instantaneous frequency fm of the sub-signal ym meets the condition for terminating iteration.

[0285] Step S14: The sub-signal ym is reconstructed.

[0286] Step S15: The sub-signal ym is removed from the current residual signal rm, and then the next sub-signal y is decomposed. m+1 and Steps S1 to S14 are repeated until each sub-signal in the original signal is completely decomposed.

[0287] The kernel function parameters of the PTFT in Step S1 are a Fourier series, and the calculation formula is: wherein The calculation formula of is:

[0288]

[0289] In the formula, K is the order of the Fourier series, f0 is the base frequency, f0 = 1 / (2T), T is the signal duration, and are coefficients of the Fourier series. The formula for calculating whis:

[0290] The instantaneous frequency ridge curve IFmof the sub-signal ymis obtained by solving the local maximum value of E(IFm) in step S2. The formula for calculating E(IFm) is:

[0291]

[0292] , where S(t, IF m (t)) represents the PTFT of the signal, and λ and β are penalty coefficients.

[0293] The formula for calculating the instantaneous frequency kernel function matrix Φm in step S4 is: m = [C m , S m ], where Cmrepresents a sine function matrix, and Smrepresents a cosine function matrix.

[0294] The formula for calculating Cmin step S4 is: m = diag [cos(θ m (t1))…cos(θ m (t N )))], where θmrepresents the instantaneous phase of the sub-signal ym.

[0295] The formula for calculating the instantaneous phase θmin step S4 is: , where represents the estimated instantaneous frequency of the sub-signal ym.

[0296] The formula for calculating Sm in step S4 is: m = diag [sin(θ m (t1))…sin(θ m (t N )))], where θmrepresents the instantaneous phase of the sub-signal ym.

[0297] The formula for calculating the demodulation signals umand vmin step S5 of the sub-signal ymis: , where n represents the number of iterations, Φmrepresents the instantaneous frequency kernel function matrix, ρ represents a penalty coefficient, D represents the diagonal matrix of the second-order difference matrix Ω rmrepresents the current residual signal.

[0298] The expression of the second-order difference matrix Ω in step S5 is: Ω is an N-row and N-column matrix, and N represents the number of discrete sampling points of the signal.

[0299] The initial instantaneous frequency increment of the sub-signal ymin step S6 is: The calculation formula of the initial instantaneous frequency increment Δf of the sub-signal ym in step S7 is: In the formula, Ω represents the second-order difference matrix, and I represents the unit matrix. And represents the demodulation signal of the sub-signal ym at the n th iteration.

[0300] The calculation formula of the initial instantaneous frequency increment Δf of the sub-signal ym in step S7 is: m The calculation formula of the initial instantaneous frequency increment Δf of the sub-signal ym in step S7 is: In the formula, Ω represents the second-order difference matrix, and I represents the unit matrix.

[0301] The update formula of the instantaneous frequency fm of the sub-signal ym in step S8 is: In the formula, n represents the iteration number.

[0302] The reconstruction calculation formula of the sub-signal ym in step S9 is: m = Φ m x m , in which Φm represents the instantaneous frequency kernel function matrix of the sub-signal ym, and xm represents the demodulation signal of the sub-signal ym.

[0303] In step S10, the iteration number n is set to 1 at the beginning, and after completing one iteration, the iteration number n is increased by 1.

[0304] The termination iteration condition of the reconstructed sub-signal ym in step S10 is: In the formula, ε1 represents the termination iteration threshold.

[0305] The termination iteration threshold ε1 of the reconstructed sub-signal ym in step S10 is generally set to 10 -8 .

[0306] The formula for fitting the K-order Fourier series of the instantaneous frequency fm of the sub-signal ym in step S11 is: In the formula, is the coefficient of the Fourier series.

[0307] In step S13, the termination iteration condition of the instantaneous frequency fm of the sub-signal ym is: In the formula, h is the iteration number, and ε2 represents the termination iteration threshold.

[0308] In step S13, the termination iteration threshold ε2 is generally set to 10 -8 .

[0309] The reconstruction calculation formula of the sub-signal ym in step S14 is: m = Φ m x m , in which Φm represents the instantaneous frequency kernel function matrix of the sub-signal ym, and xm represents the demodulation signal of the sub-signal ym.

[0310] In the step S15, m is set to 1 at the beginning, and the current residual signal r1 is equal to the original signal y. After the decomposition of the mth sub-signal is completed, the sub-signal ym is removed, and the current residual signal r m+1 m m m is set to m+1.

[0311] The method provided by the application is used for the simulation signal and the actual signal with complex variation law and strong noise interference. Figure 2 and Figure 3 is the time domain curve and the TFR of the short-time Fourier transform of the simulation signal. Figures 4 to 5 is the actual IF of the simulation signal and the time domain curve of each sub-signal. Figures 6 to 7 is the result of the decomposition of the simulation signal by the method of the application, and the result shows that the application can accurately estimate the IF of each sub-signal and reconstruct the sub-signal. Figure 8 is the time domain curve of the actual bearing vibration signal with variable rotating speed and noise interference. Figure 9 is the TFR of the short-time Fourier transform of the bearing vibration signal. Figure 10 is the TFR obtained by the decomposition of the bearing vibration signal by the method of the application, and the result shows that the application can effectively filter the noise of the original vibration signal, accurately estimate the IF of each sub-signal and reconstruct the sub-signal, and further obtain the TFR with strong energy concentration which accurately reflects the time-frequency variation law of the original signal.

[0312] Those skilled in the art can understand the embodiment as a more specific description of the embodiment 1 and the embodiment 2.

[0313] Those skilled in the art know that, in addition to implementing the system and each device, module and unit thereof provided by the application in the form of pure computer readable program code, the system and each device, module and unit thereof provided by the application can also be implemented in the form of logic gate, switch, special integrated circuit, programmable logic controller and embedded microcontroller by logically programming the method steps to achieve the same function. Therefore, the system and each device, module and unit thereof provided by the application can be considered as a hardware component, and the devices, modules and units included therein for achieving various functions can also be considered as structures in the hardware component; the devices, modules and units for achieving various functions can also be considered as both software modules for implementing the method and structures in the hardware component.

[0314] The specific embodiments of the application are described above. It should be understood that the application is not limited to the above specific embodiments, and those skilled in the art can make various changes or modifications within the scope of the claims, which does not affect the essential content of the application. The embodiments of the application and the features in the embodiments can be arbitrarily combined with each other without conflict.​​

Claims

1. A method of adaptive time-frequency supported chirplet decomposition, characterized in that, The method comprises the following steps: Step S1: calculating a parametric time-frequency transform PTFT of the signal y; Step S2: extracting the sub-signal y with the highest energy from the time- frequency representation TFR of the parametric time-frequency transform of the signal y m the coarse instantaneous frequency ridge curve IF m ; Step S3: extracting the sub-signal y m of the instantaneous frequency curve IF m as the initial instantaneous frequency f m of the signal decomposition; Step S4: Calculate the instantaneous frequency f m of the sub-signal y m Step S5: Calculate the instantaneous frequency kernel matrix Φ m of the sub-signal y Step S5: calculating the demodulation signal u m of the sub-signal y m and v m ; Step S6: Calculate the sub-signal y m of the initial instantaneous frequency increment Step S7: calculating the sub-signal y m the optimized initial instantaneous frequency increment Δf m ; Step S8: updating the sub-signal y m of the instantaneous frequency f m ; Step S9: Reconstructing the sub-signal y m ; Step S10: Steps S3-S9 are repeated until the reconstructed sub-signal y m the condition for terminating the iteration is met. Step S11: fitting the sub-signal y with a Kth order Fourier series m estimated instantaneous frequency f m ; Step S12: updating the kernel function parameters of the parametric time-frequency transform of the signal y and performing the parametric time-frequency transform on the signal y again; Step S13: Steps S1-S12 are repeated until the instantaneous frequency f m of the sub-signal y m satisfies the condition for terminating the iteration; Step S14: reconstructing the sub-signal y m ; Step S15: remove this sub-signal y m from the current residual signal r m and continue to decompose the next sub-signal y m+1 and repeat steps S1-S14 until every sub-signal in the original signal is fully decomposed.

2. The adaptive time-frequency supported FM signal decomposition method of claim 1, wherein, The kernel function parameters of the PTFT in the step S1 are Fourier series, and the calculation formula is: whereinThe calculation formula of is: where K is the order of the Fourier series, f0is the fundamental frequency, f0= 1 / (2T), T is the signal duration, are the coefficients of the Fourier series; w h The calculation formula is: In step S2, E(IF) is solved. m The local maximum value of ) yields the sub-signal y. m Instantaneous frequency ridge curve IF m E(IF) m The formula for calculating ) is: In the formula S(t,IF m (t) represents the PTFT of the signal, and λ and β are penalty coefficients; The instantaneous frequency kernel function matrix Φ in the step S4 m The calculation formula is: Φ m = [C m , S m ] where C m denotes a sine function matrix, S m denotes a cosine function matrix; C The calculation formula is: m The calculation formula is: C m = diag [cos(θ m (t1))…cos(θ m (t N ))] where θ m represents the instantaneous phase of the sub-signal y m ; Instantaneous phase θ m The formula for calculating θ is: In the formula sub-signal y m estimated instantaneous frequency; S The calculation formula is: m The calculation formula is: S m = diag[sin(θ m (t1))…sin(θ m (t N ))] where θ m represents the instantaneous phase of the sub-signal y m .

3. The adaptive time-frequency supported FM signal decomposition method of claim 1, wherein, The sub-signal y in the step S5 m The demodulation signal u m and v m The calculation formula is: where n represents the number of iterations, Φ m represents the instantaneous frequency kernel function matrix, p represents the penalty coefficient, and D represents the diagonal matrix of the second-order difference matrix Ω r m represents the current residual signal; The expression of the second-order difference matrix Ω is: Ω is an N-row and N-column matrix, and N represents the number of discrete sampling points of the signal; The step S6 sub-signal y m The initial instantaneous frequency increment The calculation formula is: wherein and denotes the sub-signal y m demodulated signal at the n-th iteration; The step S7 of calculating the sub-signal y m The initial instantaneous frequency increment Δf m The calculation formula is: In the formula, Ω represents the second-order difference matrix, and I represents the unit matrix.

4. The adaptive time-frequency supported FM signal decomposition method of claim 1, wherein, The step S8 sub-signal y m The instantaneous frequency f m The update formula is: In the formula, n represents the number of iterations. The reconstruction calculation formula of the sub-signal y m in the step S9 is: y m = Φ m x m where Φ m represents the instantaneous frequency kernel matrix of the sub-signals y m x m represents the demodulated signals of the sub-signals y m ​ The iteration number n is set to 1 in the step S10 at first, and the iteration number n is added by 1 after one iteration is completed; the reconstructed sub-signal y m The termination iteration condition of the reconstructed sub-signal y where ε1represents a termination iteration threshold; the reconstructed sub-signal y m The termination iteration threshold ε1of the reconstruction sub-signal y -8 is set to 10.

5. The adaptive time-frequency supported FM signal decomposition method of claim 1, wherein, The formula for the Kth Fourier series of the instantaneous frequency f m of the sub-signal y m in the step S11 for fitting the sub-signal y wherein are coefficients of the Fourier series; In the step S13, the instantaneous frequency f m of the sub-signal y m The termination iteration condition is: where h is the iteration number, and ε2represents the termination iteration threshold; the termination iteration threshold ε2is set to 10 -8 ; The reconstruction calculation formula of the sub-signal y in the step S14 is: m y = y1 + y2 y m = Φ m x m where Φ m represents the instantaneous frequency kernel matrix of the sub-signals y m x m represents the demodulated signals of the sub-signals y m ​ The step S15 is initially set m = 1, the current residual signal r1 is equal to the original signal y; after completing the decomposition of the mth sub-signal, the sub-signal y m is removed, and the current residual signal r m+1 = r m - y m is set m = m + 1.

6. A self-adapting time-frequency supported FM signal decomposition system, characterized in that, The system comprises the following modules: Module M1: calculating a parametric time-frequency transform PTFT of the signal y; Module M2: extracting the sub-signal y of maximum energy from the time- frequency representation TFR of the parametric time-frequency transform of the signal y m the coarse instantaneous frequency ridge curve IF m ; Module M3: the instantaneous frequency curve IF of the extracted sub-signal y m m as the initial instantaneous frequency f of the signal decomposition m ;​ Module M4: from this sub-signal y m the instantaneous frequency f m calculates the instantaneous frequency kernel matrix Φ m of this sub-signal Module M5: computing the demodulation signal u m of the sub-signal y m and v m ; Module M6: computing the initial instantaneous frequency increment of the sub-signal y m ​ Module M7: computing the optimized initial instantaneous frequency increment Δf of the sub-signal y m m ;​ Module M8: updating the sub-signal y m of the instantaneous frequency f m ; Module M9: reconstructing the sub-signal y m ; Module M10: Repeatedly call modules M3 - M9 until the reconstructed sub-signal y m The condition for terminating the iteration is met; Module M11: The sub-signal y is fitted with a Kth order Fourier series m estimated instantaneous frequency f m ; Module M12: updating the kernel function parameters of the parametric time-frequency transform of the signal y and performing the parametric time-frequency transform on the signal y again; Module M13: Repeatedly call modules M1 - M12 until the instantaneous frequency f m of the sub-signal y m satisfies the condition for terminating the iteration; Module M14: reconstructing the sub-signal y m ; Module M15: remove this sub-signal y m from the current residual signal r m and continue to decompose the next sub-signal y m+1 and repeat the call of modules M1 - M14 until every sub-signal in the original signal is fully decomposed.

7. The adaptive time-frequency supported FM signal decomposition system of claim 6, wherein, The kernel function parameters of the PTFT in the module M1 are Fourier series, and the calculation formula is: whereinThe calculation formula of is: where K is the order of the Fourier series, f0is the fundamental frequency, f0= 1 / (2T), T is the signal duration, are the coefficients of the Fourier series; w h The calculation formula is: The instantaneous frequency ridge curve IF m of the sub-signal y m is obtained by solving the local maximum of E(IF m ); the calculation formula of E(IF m ) is: In the formula S(t,IF m (t) represents the PTFT of the signal, and λ and β are penalty coefficients; The instantaneous frequency kernel function matrix Φ in the module M4 m The calculation formula is: Φ m = [C m , S m ] where C m denotes a sine function matrix, S m denotes a cosine function matrix; C The calculation formula is: m The calculation formula is: C m = diag [cos(θ m (t1))…cos(θ m (t N ))] where θ m represents the instantaneous phase of the sub-signal y m ; Instantaneous phase θ m The formula for calculating θ is: In the formula sub-signal y m estimated instantaneous frequency; S The calculation formula is: m The calculation formula is: S m = diag[sin(θ m (t1))…sin(θ m (t N ))] where θ m represents the instantaneous phase of the sub-signal y m .

8. The adaptive time-frequency supported FM signal decomposition system of claim 6, wherein, The module M5 demodulates the signal u m The calculation formula of the demodulated signal u m and v m is where n represents the number of iterations, Φ m represents the instantaneous frequency kernel function matrix, p represents the penalty coefficient, and D represents the diagonal matrix of the second-order difference matrix Ω r m represents the current residual signal; The expression of the second-order difference matrix Ω is: Ω is an N-row and N-column matrix, and N represents the number of discrete sampling points of the signal; The module M6 sub-signals y m The initial instantaneous frequency increment The calculation formula is: wherein and denotes the sub-signal y m demodulated signal at the n-th iteration; The module M7 outputs the sub-signal y m The optimized initial instantaneous frequency increment Δf m The calculation formula is: In the formula, Ω represents the second-order difference matrix, and I represents the unit matrix.

9. The adaptive time-frequency supported FM signal decomposition system of claim 6, wherein, The module M8 outputs the sub-signal y m The instantaneous frequency f m The update formula is: In the formula, n represents the number of iterations. The reconstruction calculation formula of the sub-module M9 sub-signal y m is: y m = Φ m x m where Φ m represents the instantaneous frequency kernel matrix of the sub-signals y m m represents the demodulated signals of the sub-signals y m ​​ The iteration number n is initially set to 1 in the module M10, and after one iteration is completed, the iteration number n is increased by 1; the reconstructed sub-signal y m The termination iteration condition is: where ε1represents a termination iteration threshold; the reconstructed sub-signal y m The termination iteration threshold ε1of the reconstruction sub-signal y -8 is set to 10.

10. The adaptive time-frequency supported FM signal decomposition system of claim 6, wherein, The formula for the Kth Fourier series of the instantaneous frequency f m of the sub-signal y m is wherein are coefficients of the Fourier series; In said module M13, the sub-signals y m the instantaneous frequency f m of the sub-signals y the termination iteration condition is: where h is the iteration number, and ε2represents the termination iteration threshold; the termination iteration threshold ε2is set to 10 -8 ; The reconstruction calculation formula of the sub-signal y m in the module M14 is: y m = Φ m x m where Φ m represents the instantaneous frequency kernel matrix of the sub-signal y m m represents the demodulated signal of the sub-signal y m ;​ The module M15 is initially set with m = 1, and the current residual signal r1 is equal to the original signal y; after the decomposition of the mth sub-signal is completed, the sub-signal y m is removed, and the current residual signal r m+1 = r m -y m is set with m = m + 1.

Images