A Signal Instantaneous Frequency Search Method Based on Centroid Dynamic Path Planning

Abstract

The present invention discloses a signal instantaneous frequency search method based on centroid dynamic path planning, belonging to the technical field of signal processing, and comprising the following steps: S1. Based on the time-frequency distribution of non-stationary signals, use the CFSA method to preliminarily search for a possible set of instantaneous frequency ridge lines; S2. Based on the set of instantaneous frequency ridge lines, window all elements and calculate the centroid of the ridge segments to obtain the centroid sparse matrix of the set of instantaneous frequency ridge lines; S3. Construct a dynamic path optimization function to predict the position of the instantaneous frequency points to be searched, and determine the optimal centroid points in the centroid sparse matrix; S4. Use the points with the maximum mass in each column of the centroid sparse matrix V as the starting points, and repeat S2-S3 to obtain a ridge line set based on the centroid sparse matrix V; S5. According to the energy of each element of the ridge line set, make a decision from the ridge line set to obtain the optimal instantaneous frequency ridge line. The signal instantaneous frequency search method based on centroid dynamic path planning provided by the present invention improves the accuracy and robustness of ridge line search.

Description

Technical Field

[0001] The present invention belongs to the technical field of signal processing, and in particular, relates to a signal instantaneous frequency search method based on centroid dynamic path planning. Background Art

[0002] The gearbox will generate non-stationary vibration signals related to the rotational speed under variable rotational speed conditions. It is difficult to diagnose gearbox faults using traditional stationary signal analysis methods. The instantaneous rotational speed is not only an important parameter reflecting the operating state of the gearbox, but also a necessary parameter for non-stationary signal analysis methods such as order tracking algorithms and generalized demodulation analysis. Therefore, accurately obtaining the instantaneous rotational speed of a variable-speed gearbox is of great significance for gearbox health condition monitoring and fault diagnosis.

[0003] The rotational speed corresponding to the vibration signal can be obtained by using a high-precision rotational speed measurement device. However, due to the constraints of the actual working environment, most gearboxes do not have the conditions to install a rotational speed measurement device. After obtaining the time-frequency distribution of the signal through time-frequency analysis methods such as short-time Fourier transform (STFT), wavelet transform (WT), and synchrosqueezing transform (SST), the instantaneous frequency ridge line of the signal can be obtained through an instantaneous frequency estimation algorithm. However, existing instantaneous ridge line estimation methods, such as peak search algorithm (PSA), Viterbi algorithm (VA), cost function-based search algorithm (CFSA), etc., all use discrete time-frequency points on the time-frequency diagram as the search object. The local amplitude extrema generated by strong noise will cause strong interference to the ridge line search. At the same time, these algorithms use the set of local optimal points as the final ridge line, so their results are not the global optimal solution.

[0004] Therefore, we urgently need to propose a more reliable instantaneous frequency search algorithm to improve the accuracy and robustness of the time-frequency ridge line. Summary of the Invention

[0005] The purpose of the present invention is to provide a signal instantaneous frequency search method based on centroid dynamic path planning to solve the problem that the local amplitude extrema generated by strong noise interfere with the ridge line search, resulting in low accuracy and robustness of the ridge line in the above-mentioned technology.

[0006] To achieve the above purpose, the present invention provides a signal instantaneous frequency search method based on centroid dynamic path planning, including the following steps:

[0007] S1. Based on the time-frequency distribution of the non-stationary signal, use the CFSA method to preliminarily search for a possible set of instantaneous frequency ridge lines {L k};

[0008] S2. Based on the set of instantaneous frequency ridge lines {L k}, window all elements and calculate the centroid of the ridge line segment to obtain a set of instantaneous frequency ridge lines {Lk The centroid sparse matrix V of {

[0009] S4. Take the point with the maximum mass in each column of the centroid sparse matrix V as the starting point, and repeat S2 - S3 to obtain the ridge line set based on the centroid sparse matrix V

[0010] S5. According to the ridge line set For the energy of each element in Decide the optimal instantaneous frequency ridge line

[0011] Preferably, the instantaneous frequency ridge line set {L k}} contains all possible ridge lines in the time - frequency distribution. For the instantaneous frequency ridge line set {L k}, there exists k0 ∈ {1, 2…n} and a set such that:

[0012] |f k (t l ) - f * (t l )| ≤ f w (1)

[0013] There is a certain section on the instantaneous frequency ridge line that exactly falls within [f - f , f + f w , f + f w . The set of all ridge line points that satisfy equation (1) constitutes the optimal instantaneous frequency ridge line

[0014] Preferably, the specific calculation process of the instantaneous frequency ridge line set {L k}} in S1 is as follows:

[0015] S101. For the signal z(t) ∈ L 2 (R), the window function g(t), obtain the frequency distribution of the signal Z(t) in different time periods by sliding the window, and obtain the time - frequency spectrum distribution T(t, f) of the signal

[0016]

[0017] where t, u are time, f is frequency, and the discrete form is expressed as:

[0018]

[0019] where m and n are discrete time and frequency, and N0 is the length of the signal;

[0020] S102. Define the ridge line search cost function C at the m - th time point mis:

[0021]

[0022] Wherein, t m , f m are respectively the time and frequency at the m-th moment, e is the weight coefficient, and the amplitude T and the frequency change rate |f m -f m-1 | 2 are the contribution rates to the cost function C, and f w is the ridge search radius;

[0023] S103. At T(t, f), take K - 1 points as the initial points for ridge search, and divide the signal into K segments on average in the time domain. The time interval of each segment is t s = |N0 / K|, and initialize k = 0;

[0024] S104. Let k = k + 1, m = (k - 1)·t s + 1, and the ridge point at the initial moment is which is the cost function at the m-th moment point on the k-th ridge line;

[0025] S105. Forward search: Let m = m + 1, and search for the ridge point at t w according to formula (4) with a radius of f m {where f m-1 -f w ≤ f m ≤ f m-1 + f w};

[0026] S106. Repeat S105 until m > N0 to stop the iteration, and m = (k - 1)·t s + 1;

[0027] S107. Backward search: Let m = m - 1, and search for the position of the ridge line at t w according to formula (4) with a radius of f m {where f m+1 -f w ≤ f m ≤ f m+1 + f w};

[0028] S108. Repeat S107 until m < 1 to stop the iteration, and obtain the time-frequency ridge line L k (t, f k (t), C k (t));

[0029] S109. Repeat steps S104 - S108 until k > K to stop iteration, and obtain the instantaneous frequency ridge line set {L k}

[0030] Preferably, the calculation process of the centroid sparse matrix V in S2 is as follows:

[0031] S201. Regard C k as the mass of the ridge line L k . Let the window function be g(m), the time window radius be Δ g , the window translation step be λ, and the window centroid of the k-th ridge line L k be expressed as:

[0032]

[0033] In the formula, and are the time and frequency corresponding to the centroid, M k is the centroid mass, is the centroid sequence, represents the ceiling operation;

[0034] S202. Calculate the window centroid of {L k}, and obtain the centroid set of the ridge lines as:

[0035]

[0036] S203. Store all the centroid points in in chronological order using the COO method to obtain the centroid sparse matrix V of T(t, f):

[0037]

[0038] In the formula, V row , V col and V value are the row matrix, column matrix and value matrix of the centroid set respectively.

[0039] Preferably, the specific process of S3 is as follows:

[0040] S301. Define the dynamic order {o1, o2} and the judgment threshold ξ in the dynamic path planning;

[0041] S302. When searching for the i L -th ridge line i L = 1, 2,..., B, let b = i L , and search for the point with the maximum energy in V(b, k) as the search starting point Specifically as follows

[0042]

[0043] In the formula, and are respectively the position and cost function value of the ridge point in V(b,k);

[0044] S303. Define the dynamic path planning function as:

[0045]

[0046] In the formula, are respectively the time set and frequency set of the prediction point , i = 1, 2; S is the sequence of ridge points participating in the fitting. Assume are the b ridge points that have been searched. When b < ξ, i = 1, S = {1, 2,... b}; when b > ξ, i = 2, S = {b - ξ + 1, b - ξ + 2,... b};

[0047] S304. Construct the centroid search cost function as:

[0048]

[0049] In the formula is the predicted frequency of the k-th point in, d(·) represents the distance. The centroid search cost function is a function related to the distance and the mass V value (b + 1, k);

[0050] S305. The limit point at b + 1 is expressed as:

[0051]

[0052] In the formula, V row , V col are respectively the row matrix and column matrix of the centroid set ;

[0053] S306. Repeat S302 - S305 to obtain the ridge set composed of B limits

[0054] Preferably, take the sum of the cost function values of all points on each ridge in the ridge set as the energy of the ridge, and select the ridge with the maximum energy as the optimal instantaneous frequency ridge The specific expression is as follows:

[0055]

[0056] i L represents the number of ridge lines, represents the ridge line set, represents the value of the cost function.

[0057] Therefore, the present invention adopts the above-mentioned method for searching the instantaneous frequency of a signal based on centroid dynamic path planning. First, by performing windowed centroid operation on the time-frequency distribution of the signal, the search object is converted from time-frequency discrete points to centroid points, reducing the interference of noise on the ridge line search and reducing the dimension of the time-frequency distribution. Secondly, the dynamic path planning algorithm is used to search for the time-frequency ridge line, fully considering the influence of historical data on the ridge line search, and improving the accuracy and robustness of the ridge line search.

[0058] The following will further describe the technical solution of the present invention in detail through the accompanying drawings and embodiments. Description of the Drawings

[0059] Figure 1 is the overall flowchart of a method for searching the instantaneous frequency of a signal based on centroid dynamic path planning according to the present invention;

[0060] Figure 2 is a partial view of the window centroid operation result of the present invention;

[0061] Figure 3 is the time-frequency distribution diagram obtained by performing short-time Fourier transform on the real signal Z(t) of the present invention;

[0062] Figure 4 is the sparse time-frequency distribution diagram of the signal Z(t) of the present invention. Detailed Embodiments

[0063] Embodiment

[0064] The following detailed description of the embodiments of the present invention provided in the accompanying drawings is not intended to limit the scope of the present invention claimed, but merely represents selected embodiments of the present invention. All other embodiments obtained by those of ordinary skill in the art based on the embodiments of the present invention without creative efforts fall within the scope of protection of the present invention.

[0065] Please refer to Figure 1 , a method for searching the instantaneous frequency of a signal based on centroid dynamic path planning, comprising the following steps:

[0066] S1. Based on the time-frequency distribution of the non-stationary signal, use the CFSA method to preliminarily search for the possible instantaneous frequency ridge line set {L k};

[0067] S2. Based on the instantaneous frequency ridge line set {Lk}, window all the elements and calculate the centroid of the ridge segments to obtain the instantaneous frequency ridge set {L k} of the centroid sparse matrix V;

[0068] S3. Construct a dynamic path optimization function to predict the position of the instantaneous frequency points to be searched, and determine the optimal centroid points in the centroid sparse matrix;

[0069] S4. Take the points with the maximum mass in each column of the centroid sparse matrix V as the starting points, and repeat S2 - S3 to obtain the ridge set based on the centroid sparse matrix V

[0070] S5. According to the energy of each element in the ridge set in, make a decision from to obtain the optimal instantaneous frequency ridge

[0071] The instantaneous frequency ridge set {L k} contains all possible ridges in the time - frequency distribution. Although noise may cause that there may not be a ridge L in the ridge set {L k} that can well describe the true time - frequency distribution of the signal, but for the instantaneous frequency ridge set {L k}, there exist k0 ∈ {1, 2…n} and a set such that:

[0072] |f k (t l ) - f * (t l )| ≤ f w (1)

[0073] There is a certain section on the instantaneous frequency ridge that exactly falls within [f - f w , f + f w . The set of all ridge points that satisfy equation (1) constitutes the optimal instantaneous frequency ridge

[0074] The specific calculation process of the instantaneous frequency ridge set {L k} in S1 is as follows:

[0075]

[0075] S101. For the signal z(t) ∈ L 2 (R), window function g(t), obtain the frequency distribution of the signal Z(t) in different time periods by sliding the window, and obtain the time - frequency spectrum distribution T(t, f) of the signal

[0076]

[0077] where \(t\) and \(u\) are time, \(f\) is frequency, and the discrete form is expressed as:

[0078]

[0079] where \(m\) and \(n\) are discrete time and frequency, and \(N_0\) is the length of the signal;

[0080] S102. Define the ridge search cost function \(C\) at the \(m\)-th time point m as:

[0081]

[0082] where \(t\) m , \(f\) m are the time and frequency at the \(m\)-th moment respectively, \(e\) is the weight coefficient, and the contribution rate of the amplitude \(T\) and the frequency change rate \(|f\) m - \(f\) m-1 |\) 2 to the cost function \(C\) is determined, and \(f\) w is the ridge search radius;

[0083] S103. At \(T(t, f)\), take \(K - 1\) points as the initial points for ridge search, and divide the signal into \(K\) segments evenly in the time domain. The time interval of each segment is \(t\) s = \(|N_0 / K|\), and initialize \(k = 0\);

[0084] S104. Let \(k = k + 1\), \(m=(k - 1)\cdot t\) s + 1, and the ridge point at the initial moment is the cost function at the \(m\)-th moment on the \(k\)-th ridge line;

[0085] S105. Forward search: Let \(m = m + 1\), and search for the ridge point at \(t\) w according to formula (4) with a radius of \(f\) m {where \(f\) \(f\) m-1 - \(f\) w \(\leq f\) m \(\leq f\) m-1 + \(f\) w};

[0086] S106. Repeat S105 until \(m > N_0\) to stop the iteration, and \(m=(k - 1)\cdot t\) s + 1;

[0087] S107. Backward search: Let \(m = m - 1\), and search for the position of the ridge at \(t\) w according to formula (4) with a radius of \(f\) m {where \(f\) \(f\) m+1 - \(f\) w \(\leq f\) m \(\leq f\)m+1 +f w};

[0088] S108. Repeat S107 until m < 1 to stop the iteration, and obtain the time-frequency ridge line L k (t, f k (t), C k (t));

[0089] S109. Repeat steps S104 - S108 until k > K to stop the iteration, and obtain the set of instantaneous frequency ridge lines {L k}}.

[0090] The calculation process of the centroid sparse matrix V in S2 is as follows:

[0091] S201. Regard C k as the mass of the ridge line L k . Let the window function be g(m), the time window radius be Δ g , the window translation step be λ, and the window centroid of the k-th ridge line L k be expressed as:

[0092]

[0093] In the formula, and are the time and frequency corresponding to the centroid, M k is the centroid mass, is the centroid sequence, represents the ceiling operation;

[0094] Through the window centroid operation, the influence of abnormal ridge line points on the time-frequency distribution is weakened, so as to achieve the purpose of eliminating noise interference. At the same time, a centroid point reflects the overall distribution of ridge line points within its window, so the dimension of the time-frequency distribution is reduced. As Figure 2 shown, the gray and black alternating line segments represent the windowed ridge line segments, "*" is the centroid of this section of the ridge line, and the light gray dotted line represents the theoretical curve.

[0095] S202. Calculate the window centroids of {L k}, and obtain the centroid set of the ridge lines as:

[0096]

[0097] S203. Store all the centroid points in in chronological order using the COO method to obtain the centroid sparse matrix V of T(t, f):

[0098]

[0099] In the formula, Vrow , V col and V value are the row matrix, column matrix, and numerical matrix of the centroid set respectively. By constructing a centroid sparse matrix, the time-frequency ridge search object is transformed from the time-frequency distribution T(t,f) to a sparse time-frequency distribution V, retaining the key information in the time-frequency diagram and greatly compressing the data volume. When searching for the ridge points at a certain moment, only the target elements need to be searched in the COO matrix at that moment, thus avoiding traversing all the data of the sparse time-frequency distribution and achieving the purpose of reducing the computational amount, as Figure 3 and 4 shown.

[0100] The specific process of S3 is as follows:

[0101] S301. Define the dynamic order {o1, o2} and the judgment threshold ξ in the dynamic path planning;

[0102] S302. When searching for the i L th ridge i L = 1, 2,..., B, let b = i L , and search for the point with the maximum energy in V(b,k) as the search starting point Specifically as follows:

[0103]

[0104] In the formula, and are the position and cost function value of the ridge point in V(b,k) respectively;

[0105] S303. Define the dynamic path planning function as:

[0106]

[0107] In the formula, are the time set and frequency set of the prediction point respectively, i = 1, 2; S is the sequence of ridge points participating in the fitting. Assume are the b ridge points that have been searched. When b < ξ, i = 1, S = {1, 2,... b}; when b > ξ, i = 2, S = {b - ξ + 1, b - ξ + 1,... b};

[0108] S304. Construct the centroid search cost function as:

[0109]

[0110] In the formula is The predicted frequency of the k-th point in, where d(·) represents the distance, and the centroid search cost function is a function related to the distance and the mass V value (b + 1, k);

[0111] S305. The limit point at b + 1 is expressed as:

[0112]

[0113] In the formula, V row , V col are respectively the row matrix and column matrix of the centroid set ;

[0114] S306. Repeat S302 - S305 to obtain a ridge line set composed of B limits

[0115] Take the sum of the cost function values of all points on each ridge line in the ridge line set as the energy of the ridge line, and select the ridge line with the maximum energy as the optimal instantaneous frequency ridge line The specific expression is as follows:

[0116]

[0117] i L represents the number of ridge lines, represents the ridge line set, represents the cost function value.

[0118] Therefore, the present invention adopts the above signal instantaneous frequency search method based on centroid dynamic path planning to solve the problem that the local amplitude extreme values generated by strong noise in the traditional technology will interfere with the ridge line search, resulting in low accuracy and robustness of the ridge line. The present invention first performs windowed centroid operation on the time-frequency distribution of the signal, converts the search object from time-frequency discrete points to centroid points, reduces the interference of noise on the ridge line search, and reduces the dimension of the time-frequency distribution. Secondly, the dynamic path planning algorithm is used to search for the time-frequency ridge line, fully considering the influence of historical data on the ridge line search, and improving the accuracy and robustness of the ridge line search.

[0119] Finally, it should be noted that the above embodiments are only used to illustrate the technical solutions of the present invention and are not intended to limit them. Although the present invention has been described in detail with reference to the preferred embodiments, those of ordinary skill in the art should understand that they can still modify or equivalently replace the technical solutions of the present invention, and these modifications or equivalent replacements cannot make the modified technical solutions deviate from the spirit and scope of the technical solutions of the present invention.

Claims

1. A signal instantaneous frequency search method based on centroid dynamic path planning, characterized in that Including the following steps: S1. Based on the time-frequency distribution of non-stationary signals, use the CFSA method to preliminarily search for a set of possible instantaneous frequency ridge lines {L k}; S2. Based on the instantaneous frequency ridge line set {L k}, window all elements and calculate the centroid of the ridge line segments to obtain the centroid sparse matrix V of the instantaneous frequency ridge line set {L k}; S3. Construct a dynamic path optimization function to predict the position of the instantaneous frequency point to be searched, and determine the optimal centroid point in the centroid sparse matrix; S4. Take the point with the maximum mass in each column of the centroid sparse matrix V as the starting point, and repeat S2 - S3 to obtain the ridge line set based on the centroid sparse matrix V S5. According to the energy of each element in the ridge line set decide the optimal instantaneous frequency ridge line from ​ 2. The signal instantaneous frequency search method based on centroid dynamic path planning according to claim 1, wherein The set of instantaneous frequency ridge lines {L k} contains all possible ridge lines in the time-frequency distribution. For the set of instantaneous frequency ridge lines {L k}, there exists k0 ∈ {1, 2…n} and a set such that: Instantaneous frequency ridge There is a certain section that exactly falls within [[f * - f w , f * + f w . The set of all ridge points that satisfy Equation (1) constitutes the optimal instantaneous frequency ridge 3. A method for searching the instantaneous frequency of a signal based on centroid dynamic path planning according to claim 2, wherein The instantaneous frequency ridge line set {L k} in S1 is calculated as follows: S101. For the signal z(t) ∈ L 2 (R), using the window function g(t), obtain the frequency distribution of the signal Z(t) in different time periods by sliding the window, and obtain the time-frequency spectrum distribution T(t, f) of the signal In the formula, t and u are time, f is frequency, and the discrete form is expressed as: In the formula, m and n are discrete time and frequency, and N0 is the length of the signal; S102. Define the ridge line search cost function \(C\) at the \(m\)-th moment point m as follows: where t m , f m are the time and frequency at time m respectively, e is the weight coefficient, and the amplitude T and the frequency change rate |f m - f m-1 | 2 are the contribution rates to the cost function C m , and f w is the ridge search radius; S103. At T(t, f), take K - 1 points as the initial points for ridge line search, and evenly divide the signal into K segments in the time domain, with the time interval for each segment being t s = |N0 / K|, initialize k = 0; S104. Let \(k = k + 1\) and \(m=(k - 1)\cdot t+1\). The ridge point at the initial moment is s The cost function for the \(m\)-th moment point on the \(k\)-th ridge line is ​ S105. Forward search: Let m = m + 1, and search for the ridge point at t according to Equation (4) with a radius of f w Search for t m at the ridge point S106. Repeat S105 until m > N0 to stop iteration, where m = (k - 1)·t s + 1; S107. Negative search: Let m = m - 1, and search for the position of the ridge line at t according to Equation (4) with a radius of f w Search for t m at the position of the ridge line S108. Repeat S107 until m < 1 to stop the iteration, and obtain the time-frequency ridge line L k (t, f k (t), C k (t)); S109. Repeat steps S104 - S108 until k > K to stop the iteration and obtain the instantaneous frequency ridge set {L k}.

4. A method for searching the instantaneous frequency of a signal based on centroid dynamic path planning according to claim 3, characterized in that, The calculation process of the centroid sparse matrix V in S2 is as follows: S201. Regard C k as the quality of the ridge line L k . Assume the window function is g(m), the time window radius is Δ g , the window translation step size is λ, and the window centroid of the k-th ridge line L k is expressed as: where r t k and are the time and frequency corresponding to the centroid, M k is the centroid mass, is the centroid sequence, represents the ceiling operation; S202. Calculate the window centroid of {L k}, and obtain the centroid set of the ridge line as follows: S203. Store all the centroid points in using the COO method in chronological order to obtain the centroid sparse matrix V of T(t ,f) : wherein, V row , V col and V value are respectively the row matrix, column matrix and numerical matrix of the centroid set .

5. A method for searching the instantaneous frequency of a signal based on centroid dynamic path planning according to claim 4, characterized in that, The specific process of S3 is as follows: S301. Define the dynamic order {o1, o2} and the judgment threshold ξ in the dynamic path planning; S302. Search for the i L th ridge line i L = 1, 2,..., B, let b = i L , and search for the point with the maximum energy in V(b, k) as the search starting point Specifically as follows Wherein, and are the position and cost function value of the ridge point in V(b,k), respectively; S303. Define the dynamic path planning function as: In the formula, are respectively the time set and frequency set of the prediction point , where i = 1, 2; S is the sequence of ridge points participating in fitting. Assume are the b ridge points that have been searched. When b < ξ, i = 1, S = {1, 2, … b}; when b > ξ, i = 2, S = {b - ξ + 1, b - ξ + 2, … b}; S304. Construct the centroid search cost function It is as follows: In the formula is the predicted frequency of the k-th point in , d(·) represents the distance, and the centroid search cost function is a function related to the distance and the mass V value (b + 1, k); S305. Limit point at b + 1 Expressed as: Where, V row , V col are respectively the row matrix and column matrix of the centroid set ; S306. Repeat S302 to S305 to obtain a set of ridge lines composed of B limits 6. The signal instantaneous frequency search method based on centroid dynamic path planning according to claim 5, wherein: Take the sum of the cost functions of all points on each ridge line in the ridge line set as the energy of the ridge line, and select the ridge line with the maximum energy as the optimal instantaneous frequency ridge line The specific expression is as follows: i L represents the number of ridge lines represents the ridge line set represents the cost function value

Images