Wignerhoff Transform Algorithm Based on SG Filtering and Switched Median Filtering

Abstract

This invention proposes a Wigner-Hough transform algorithm based on Savitzky-Golay (S-G) filtering and switched median filtering, suitable for parameter estimation of linear frequency modulated (LFM) signals under low signal-to-noise ratio (SNR) conditions. The algorithm estimates two parameters of the LFM signal: the starting frequency and the modulation slope. First, the LFM signal containing Gaussian white noise is subjected to Savitzky-Golay filtering (S-G filtering). Then, a Wigner-Hough transform is applied to the signal to obtain its Wigner-Hough time-frequency distribution. This distribution is then subjected to a second filtering using switched median filtering. Following this, a Wough transform is performed, and the peak point coordinates are searched. Substituting these coordinates into a specific formula yields the estimated values ​​of the two parameters. This invention improves the algorithm's noise robustness without reducing estimation accuracy. The main contribution of this algorithm is that, compared to the traditional WHT algorithm, it improves the noise robustness of both the modulation slope and starting frequency estimations by 2 dB without reducing estimation accuracy.

Description

Technical Field

[0001] This invention belongs to the field of radar reconnaissance signal processing in radar countermeasures, and in particular relates to a Wigner Hough Transform (WHT) algorithm based on SG filtering and switched median filtering. Background Technology

[0002] Radar reconnaissance, jamming, and electronic defense are the main components of radar countermeasures. Radar reconnaissance aims to use radar reconnaissance equipment to intercept, modulate, identify, measure parameters, analyze location, and sort radar signals emitted by non-cooperative parties, thereby obtaining key information such as radar-related technical parameters, location deployment, and system type. Effective and accurate parameter estimation is crucial for radar reconnaissance signal processing, playing a vital role in radar jamming, radar electronic defense, and even combat strategies. However, with the increasing complexity of the electromagnetic environment, the signal-to-noise ratio of intercepted signals is low, severely affecting the accuracy of subsequent parameter estimation. Therefore, it is essential to research algorithms that can still effectively estimate signals under low signal-to-noise ratio conditions. Linear Frequency Modulation (LFM) signals, as a typical non-stationary signal, are widely used in communications, radar, sonar, and seismic exploration. Parameter estimation of LFM signals has always been a hot research topic and is also the subject of this invention.

[0003] There are many existing parameter estimation algorithms for LFM signals, among which the Wigner-Ville Transform (WHT)-based algorithm is an important class. The traditional WHT combines the Wigner-Ville Distribution (WVD) from the bilinear time-frequency distribution with the Hough Transform from image feature extraction methods. WVD exhibits good energy concentration for LFM signals, while the Hough Transform transforms the parameter estimation problem into a problem of searching for local maxima and their corresponding coordinates in the parameter space. The traditional WHT algorithm is a classic signal parameter estimation algorithm, offering good estimation accuracy for LFM signals and exhibiting strong robustness. However, it is susceptible to noise interference at low signal-to-noise ratios, leading to a sharp performance degradation. Summary of the Invention

[0004] To address the aforementioned technical problems, this invention proposes a Wigner-Ville transform algorithm based on SG filtering and switched median filtering. First, the LFM signal containing Gaussian white noise is subjected to Savitzky-Golay filtering (SG filtering). Then, the signal undergoes a Wigner-Ville transform to obtain the Wigner-Ville distribution (WVD) time-frequency graph. This distribution is then subjected to a second filtering using switched median filtering. Following the WVD, the peak point coordinates are searched, and by substituting them into a specific formula, the estimated values ​​of the two parameters can be obtained, further improving the algorithm's noise resistance.

[0005] To achieve the above objectives, the technical solution adopted by the present invention is as follows:

[0006] A Wignerhoff transform algorithm based on SG filtering and switched median filtering includes the following steps:

[0007] Step 1: Perform SG filtering on the noisy, low signal-to-noise ratio LFM signal;

[0008] Step 2: Perform a Wigner-Nawell transform on the low signal-to-noise ratio LFM signal after SG filtering to obtain the Wigner-Nawell time-frequency distribution diagram of the signal;

[0009] Step 3: Perform a switching median filter on the Wignerville time-frequency distribution plot of the signal to obtain the filtered time-frequency distribution plot;

[0010] Step 4: Perform Hough transform on the time-frequency distribution graph after switching median filtering;

[0011] Step 5: Perform peak search on the peak values ​​formed by the signal energy after Hough transform, and calculate the estimated values ​​of the parameters using the coordinates of the peak points.

[0012] Furthermore, step 1 includes the following steps:

[0013] Step 1.1, set the fitting polynomial y = a0 + a1x + a2x 2 +...+a l-1 x l-1 The fitting order; where x is the data to be fitted at a certain measurement point within the filter window, y is the output data after fitting at that measurement point, and a p (p = 0, 1, ..., l-1) are the weighting coefficients to be solved;

[0014] Step 1.2: Select the window width;

[0015] Step 1.3: Determine the weighting coefficients using the least squares method;

[0016] Step 1.4: Move the window to the right sequentially and repeat step 1.3 until the entire signal has been traversed.

[0017] Furthermore, step 3 includes the following steps:

[0018] Step 3.1: Search for the maximum amplitude of energy in the Wigner-Ville time-frequency distribution obtained in Step 2. And set the threshold for determining whether the amplitude of an image point corresponds to the magnitude of noise energy to... in, Let f be a Wignerwell distribution, where t is the time variable and f is the frequency.

[0019] Step 3.2: For each image point in the obtained Wignerville time-frequency distribution map, sort the amplitude values ​​of all image points within a 3×3 window centered on that image point in ascending order;

[0020] Step 3.3: Compare the absolute value of the amplitude of the image point at the center of the window with a threshold. If the absolute value of the amplitude of the image point is less than or equal to the set threshold, it is determined that the amplitude of the image point corresponds to the amount of noise generated, and the median of the amplitudes of the sorted image points is taken as the amplitude of the current image point. Otherwise, the amplitude of the image point remains unchanged.

[0021] Step 3.4: Slide the window from the top to the bottom of the Wigner-Williams time-frequency distribution plot, moving from left to right. Repeat step 3.3 until the entire Wigner-Williams time-frequency distribution plot has been traversed.

[0022] Furthermore, step 4 includes the following steps:

[0023] Let the Wigner-Weil distribution of the signal after the switching median filter be: The Wignerhoff transform of the signal is:

[0024]

[0025] Where f0 is the starting frequency of the LFM signal, k is the frequency modulation slope; s1(t) represents the noisy signal after SG filtering, τ represents the time delay, * represents taking the conjugate, and j represents the imaginary number;

[0026] In polar coordinates, the Wignerhoff transform of the signal is:

[0027]

[0028] Where ρ is the length of the perpendicular line passing through the origin, and θ is the angle between ρ and the X-axis.

[0029] Furthermore, step 5 includes the following steps:

[0030] Step 5.1: Perform peak search on the Wignerhoff transform result obtained in Step 4 to obtain the coordinates corresponding to the peak points;

[0031] Step 5.2: Estimate the parameters of the coordinates corresponding to the peak points; for a given LFM signal, the sampling rate is f. s Let the number of sampling points be n³, the time-resolution unit on the Wigner-Hough distribution time-frequency graph be Δt, the frequency-resolution unit be Δf, the length and width of the Hough transform image be n³, and the coordinates of the peak point after the Wigner-Hough transform be (ρ, θ). Then the following relationship holds:

[0032]

[0033]

[0034] The estimation formulas for the starting frequency f0 and the frequency modulation slope k of the LFM signal are as follows:

[0035]

[0036]

[0037] Beneficial effects:

[0038] This invention is applicable to parameter estimation of linear frequency modulated (LFM) signals under low signal-to-noise ratio (SNR) conditions, primarily focusing on estimating two parameters: the starting frequency and the modulation slope. This invention improves the algorithm's noise robustness without compromising estimation accuracy. The main contribution of this invention's algorithm is that, compared to the traditional WHT algorithm, it improves the noise robustness of both the modulation slope and starting frequency estimations by 2 dB without reducing estimation accuracy. Attached Figure Description

[0039] Figure 1 This is a graph showing the variation of the LFM signal parameter NRMSE with the window length.

[0040] Figure 2 This is a flowchart of the Wignerhoff transform algorithm based on SG filtering and switched median filtering of the present invention.

[0041] Figure 3a A comparison of the time-domain waveforms of the LFM signal before and after noise addition.

[0042] Figure 3b This is a comparison of the time-domain waveforms of the LFM signal after noise addition.

[0043] Figure 4 The time-domain waveform of the noisy LFM signal after SG filtering.

[0044] Figure 5 WVD time-frequency plot of the noisy LFM signal.

[0045] Figure 6 The WVD time-frequency diagram is the WVD of the noisy LFM signal after SG filtering.

[0046] Figure 7 This is the WVD time-frequency diagram of the signal after switching median filtering.

[0047] Figure 8 The WHT 3D plot of the noisy LFM signal with SNR = -13dB.

[0048] Figure 9 The WHT 3D plot is the WHT of the LFM signal with added noise at SNR = -13dB after two filtering processes.

[0049] Figure 10 The figure shows the results of 100 Monte Carlo experiments for estimating LFM signal parameters using traditional WHT.

[0050] Figure 11 The figure shows the results of 100 Monte Carlo experiments on the parameter estimation of LFM signals using the novel WHT. Detailed Implementation

[0051] To make the objectives, technical solutions, and advantages of this invention clearer, the invention will be further described in detail below with reference to embodiments. It should be understood that the specific embodiments described herein are merely illustrative and not intended to limit the invention.

[0052] like Figure 2 As shown, the Wignerhoff transform algorithm based on SG filtering and switched median filtering of the present invention includes the following steps:

[0053] Step 1: Perform SG filtering on the noisy, low signal-to-noise ratio LFM signal;

[0054] Step 2: Perform Wigner-Vell transform on the filtered signal to obtain the Wigner-Vell (WVD) time-frequency distribution diagram of the signal;

[0055] Step 3: Perform on-off median filtering on the Wignerville time-frequency distribution plot of the signal;

[0056] Step 4: Perform Hough transform on the time-frequency distribution map after switching median filtering;

[0057] Step 5: Search for the coordinates of the formed peak.

[0058] Furthermore, step 1 includes the following steps:

[0059] The linear frequency modulated signal model containing noise is defined as follows:

[0060]

[0061] Where A0 is the amplitude of the signal, t is the time variable, T is the pulse width, f0 is the starting frequency of the signal, k is the frequency modulation slope, w(t) is zero-mean Gaussian white noise uncorrelated with the signal, and exp[] is an exponential function with the natural constant e as the base.

[0062] The key to SG filtering lies in solving the matrix operator. Let the width of the filtering window be n1 = 2m + 1, where n1 is a positive odd number, m is half the width of the filtering window, and x is the measurement point within the window. q (q=-m,-m+1,...,0,1,...,m-1,m), a polynomial of degree l-1 is used to fit the data of a certain measurement point within the window:

[0063] y = a0 + a1x + a2x 2 +...+a l-1 x l-1 (2)

[0064] Where x is the data to be fitted at a certain measurement point, y is the fitted output data, and a p (p = 0, 1, ..., l-1) are the weighting coefficients to be solved. There are a total of n1 measurement points in the filtering window, so n1 = 2m + 1 equations similar to equation (2) can be obtained, forming a system of l linear equations. In order for the system of equations to have a solution, n1 should be greater than or equal to l. Usually, n1 > l is chosen, and the weighting coefficients a are determined by least squares fitting. p b q (q = -m, -m+1, ..., 0, ..., m-1, m) represents the fitting difference between the output data and the original data. Therefore, we can obtain:

[0065]

[0066] Step 1.1, setting the order of the polynomial series fitting. A higher order results in a smoother filtered signal, but also takes longer. Furthermore, when the order is large, the fitting process becomes problematic due to window length limitations, and high-frequency curves may become straight lines. To avoid excessive computation and prevent the curve from becoming a straight line, this invention uses a polynomial series fitting order of l = 3.

[0067] Step 1.2: Select the window width. With an order l = 3, both n1 > l and n1 being a positive odd number must be satisfied simultaneously. Therefore, window lengths of 5, 7, 9, and 11 are chosen. With a signal-to-noise ratio of -12dB, 100 Monte Carlo experiments are conducted to calculate the changes in the signal modulation slope and the normalization root mean square error (NRMSE) of the initial frequency as a function of the window length. The experimental results are as follows: Figure 1 As shown in the figure, as the window length increases, the NRMSE also increases. Therefore, in this invention, the window length is selected as n1 = 5.

[0068] Step 1.3: Apply the least squares method to determine the weighting coefficients. Based on steps 1.1 and 1.2, this invention uses the five-point cubic smoothing formula, i.e., m = 2, n1 = 2 × 2 + 1 = 5, l = 3. Substituting these into equation (3), we get:

[0069]

[0070] The matrix representation is as follows:

[0071] Y 5×1 =X 5×3 ·A 3×1 +B 5×1 (5)

[0072] Among them, Y 5×1 X represents the output value to be solved. 5×3 For the observations within the selected window, A 3×1 Least square solution for:

[0073]

[0074] Then Y 5×1 Filter value for:

[0075]

[0076] The estimated value of the center point within the window is obtained using the fitted polynomial.

[0077] Step 1.4: Move the window to the right sequentially and repeat step 1.3 until the entire signal has been traversed.

[0078] Furthermore, step 2 includes the following steps:

[0079] Step 2.1: Perform a Wigner-Raphson transform on the noisy signal s1(t) after SG filtering to obtain the Wigner-Raphson distribution of the signal. Its expression is:

[0080]

[0081] Where f is the frequency, s1(t) represents the noisy signal after SG filtering, τ represents the time delay, * represents taking the conjugate, and j represents the imaginary number.

[0082] The WVD of LFM signals exhibits excellent time-frequency convergence; after the WVD, the energy of the LFM signal is concentrated on a straight line in the time-frequency domain. Since the amplitudes of signal energy and noise energy concentrated in the time-frequency domain differ significantly, the concept of switched median filtering is introduced. This involves comparing the absolute amplitude of each image point in the time-frequency graph with a preset threshold to determine whether the image point is a noise point.

[0083] Furthermore, step 3 includes the following steps:

[0084] Step 3.1: Search for the maximum amplitude of energy in the time-frequency distribution map obtained in Step 2. A threshold is set based on the search results. If the threshold is set too high, it can easily filter out image points corresponding to the signal energy. Therefore, the threshold set in this invention is [value missing].

[0085] Step 3.2: Select a window of appropriate size and sort the pixels within the window. Let x... ij This represents the amplitude of image point (i,j) in the WVD time-frequency plot obtained in step 2. This represents the magnitude x of the image point centered at point (i,j). ij Perform n²×n²=2N+1 window operations (where n² is a positive odd number and N is half the window width, a positive integer). Since a smaller window size preserves more details of the time-frequency distribution map, this invention uses a 3×3 window size. The amplitude values ​​corresponding to each image point within the window are then sorted in ascending order. Let the magnitude of the k-th image point after sorting within the window be:

[0086]

[0087] Step 3.3: Compare the absolute value of the amplitude of the image point at the center of the window with the threshold. If the absolute value of the amplitude of the image point is less than the threshold... Then, determine the magnitude of the noise energy corresponding to the amplitude of the image point, and take the median of the amplitudes of the sorted image points as the amplitude of the current image point. Let... Indicates the window The median magnitude of all image points within the range is taken. If y ij x is the amplitude of the image point. ij The output value after switching median filtering is:

[0088]

[0089] Here, med(·) indicates taking the median value.

[0090] Step 3.4: Slide the window from the top to the bottom of the time-frequency plot, moving from left to right. Repeat step 3.3 until the entire Wignerville time-frequency distribution plot has been traversed.

[0091] Furthermore, step 4 includes the following steps:

[0092] Step 4.1: Perform a Hough transform on the time-frequency distribution obtained from the switched median filter in Step 3. Since the Wigner-Weil distribution of the LFM signal is a straight line, a coordinate transformation is performed using the Hough transform to convert the points on the straight line in the original image into curves in polar coordinates. Let the Wigner-Weil distribution of the signal after the switched median filter be... The WHT of the signal is:

[0093]

[0094] Where f0 is the initial frequency of the LFM signal, and k is the frequency modulation slope. In polar coordinates, the WHT of the signal is:

[0095]

[0096] Where ρ is the length of the perpendicular line passing through the origin, and θ is the angle between ρ and the X-axis.

[0097] Furthermore, step 5 includes the following steps:

[0098] Step 5.1: Perform peak search on the WHT result of the signal obtained in Step 4 to obtain the coordinates of the peak point.

[0099] Step 5.2: Substitute the peak point coordinates into a specific formula to calculate the estimated value. For a given LFM signal, the sampling rate is f. s Let the number of sampling points be n3, the time resolution unit on the Wigner-Wilwell distribution time-frequency plot be Δt, the frequency resolution unit be Δf, and the length and width of the Hough transform image both be n3. Given the coordinates (ρ, θ) of the peak point after the WHT, the following relationship holds:

[0100]

[0101]

[0102] The estimation formulas for the starting frequency f0 and the frequency modulation slope k of the LFM signal are as follows:

[0103]

[0104]

[0105] To verify the effectiveness of this invention, we performed the following simulation. The LFM signal parameters were set as follows: amplitude A0 = 1, starting frequency f0 = 50Hz, frequency modulation slope k = 40Hz / s, and signal-to-noise ratio SNR = -13dB. The time-domain waveforms of the signal before and after noise addition are shown in the following figure. Figure 3a , Figure 3b As shown. An LFM signal with added Gaussian white noise and a signal-to-noise ratio of -13dB is subjected to SG filtering, resulting in the time-domain waveform of the filtered signal as shown. Figure 4 As shown, the WVD time-frequency diagram of the signal is obtained by performing a Wigner-Vell transform on the SG-filtered signal. To more clearly compare the effect of SG filtering, Figure 5 and Figure 6 The WVD time-frequency plots of the noisy signal and the noisy signal after SG filtering are presented respectively. Clearly, the time-frequency plot of the LFM signal without SG filtering is completely covered by noise, making it impossible to extract signal information. However, the WVD time-frequency plot of the SG-filtered signal still shows a straight line composed of energy containing signal information, demonstrating the effectiveness of SG filtering in suppressing Gaussian white noise, which is beneficial for subsequent steps.

[0106] The WVD of the SG-filtered signal is then subjected to a switched median filter to obtain the time-frequency plot of the switched median filter, as shown below. Figure 7 As shown, a Hough transform is further performed on the time-frequency graph after the switching median filter, resulting in the WHT three-dimensional graph as shown. Figure 9 As shown. For comparison, Figure 8 The WHT three-dimensional plot of the noisy LFM signal without two filtering operations is given. Obviously, the algorithm proposed in this invention has a significant noise suppression effect, and the peak formed by the accumulation of signal energy is more prominent.

[0107] To further verify the effectiveness of the present invention and the improvement in noise resistance, the normalized root mean square error (NRMSE) was selected to evaluate the performance of the algorithm. NRMSE is defined as follows:

[0108]

[0109] Where L represents the number of Monte Carlo experiments, x real This represents the actual value of the parameter. This represents the estimated value of the parameter. This represents the difference between the true value and the estimated value of the parameter in the z-th experiment.

[0110] The traditional WHT algorithm and the proposed novel WHT algorithm based on SG filtering and switched median filtering were used to estimate the parameters of the LFM signal. One hundred Monte Carlo experiments were conducted to calculate the normalized root mean square error (NRMSE) of the two parameters, the starting frequency and the frequency modulation slope. The curves of NRMSE versus signal-to-noise ratio were plotted, and the results are shown below. Figure 10 (Results of 100 Monte Carlo experiments for estimating LFM signal parameters using traditional WHT: NRMSE of both parameters < 0.05 when the signal-to-noise ratio is greater than -10dB) Figure 11 (The results of 100 Monte Carlo experiments on LFM signal parameter estimation using the novel WHT show that the NRMSE of both parameters is <0.05 when the signal-to-noise ratio is greater than -12dB). It can be seen that the algorithm proposed in this invention, compared to the traditional WHT, improves noise immunity while maintaining the estimation accuracy of the two parameters, the starting frequency and the modulation slope of the LFM signal. Specifically, the noise immunity of this algorithm is improved by 2dB for both the modulation slope and the starting frequency estimation.

[0111] The above description is only a preferred embodiment of the present invention. It should be noted that for those skilled in the art, several improvements and modifications can be made without departing from the principle of the present invention, and these improvements and modifications should also be considered within the scope of protection of the present invention.

Claims

1. A Wigner-Hough transform algorithm based on S-G filtering and switch median filtering, characterized by, The method comprises the following steps: Step 1, S-G filtering of a low signal-to-noise ratio LFM signal containing noise, comprising the following steps: Step 1.1, setting the fitting polynomial wherein, is the data to be fitted for a measurement point within the filter window, is the output data after fitting for the measurement point, a p is the weighted coefficient to be solved, p = 0, 1, …, l - 1; Step 1.2, selection of window width; Step 1.3, application of a least square method to determine a weighting coefficient; Step 1.4, sequentially moving the window to the right and repeating step 1.3 until the entire signal is traversed; Step 2, Wigner-Ville transformation of the S-G filtered low signal-to-noise ratio LFM signal to obtain a Wigner-Ville time-frequency distribution of the signal; Step 3, switch median filtering of the Wigner-Ville time-frequency distribution of the signal to obtain a filtered time-frequency distribution, comprising the following steps: Step 3.1, search for the maximum amplitude of energy in the Wigner-Ville time-frequency distribution obtained in step 2 ; wherein, is the Wigner-Ville distribution, is the time argument, is the frequency; Step 3.2, for each image point in the obtained Wigner-Ville time-frequency distribution, the amplitudes of all image points in the window centered at the image point are sorted in ascending order; the amplitudes of all image points in the window centered at the image point are sorted in ascending order; Step 3.3, compare the absolute value of the amplitude of the window center image point with the threshold value, if the absolute value of the amplitude of the image point is less than then determine that the amplitude of the image point corresponds to the energy size of the noise, take the median of the sorted image point amplitudes as the amplitude of the current image point; otherwise, the amplitude of the image point remains unchanged; Step 3.4, sliding the window from the top to the bottom of the Wigner-Ville time-frequency distribution, from left to right, repeating step 3.3 until the entire Wigner-Ville time-frequency distribution is traversed; Step 4, Hough transformation of the switch median filtered time-frequency distribution; Step 5, peak searching of a peak value formed by signal energy after Hough transformation, and calculation of an estimated value of a parameter using a peak point coordinate.

2. The S-G filter and switch median filter based Wigner-Hough transform algorithm according to claim 1, characterized in that, The step 4 comprises the following steps: Let the Wigner-Ville distribution of the signal filtered by the switch median filter be The Wigner-Hough transform of the signal is （12） wherein, is the start frequency of the LFM signal, is the frequency modulation slope; denotes the S-G filtered noisy signal, denotes the time delay, denotes the conjugate, denotes the imaginary number; In polar form, the Wigner Hough transformation of the signal is: （13） wherein, is the length of the perpendicular through the origin, is the angle with the X-axis.

3. The S-G filter and switch median filter based Wigner-Hough transform algorithm according to claim 2, characterized in that, The step 5 comprises the following steps: Step 5.1, peak searching of a result of the Wigner Hough transformation obtained in step 4 to obtain a coordinate corresponding to a peak point; Step 5.2, parameter estimation using the coordinates of the peak point; for a given LFM signal, the sampling rate is , let the number of sampling points be , the time resolution unit on the time-frequency diagram of the Wigner-Ville distribution is , the frequency resolution unit is , the length and width of the image of the Hough transform are both , and the coordinates of the peak point after the Wigner-Hough transform are , then the following relationship holds: （14） （15） The starting frequency of the LFM signal And the estimation formula of the frequency modulation slope k is as follows: （16） （17）。

Images