Method for estimating wave direction and wave period from X-band radar image

Abstract

The invention discloses a method for estimating a wave direction and a wave period from an X-band radar image, and the method comprises the steps: dividing each X-band radar image into a plurality of sub-images, extracting the texture features of the sub-images through two-dimensional empirical curve transformation (2D-ECT), and decomposing the sub-images into a plurality of two-dimensional curve components; selecting a curvelet component with the highest relevancy with the sub-image, and estimating a wave direction by using empirical ridge transform (ERT); meanwhile, 180-degree ambiguity of the wave direction is eliminated by using a cross correlation coefficient; in order to estimate a wave period, calculating a periodogram power spectrum by the sub-images after the texture features are extracted; then, obtaining the frequency corresponding to the maximum power spectral density in the azimuth angle pixels of the sub-images, and estimating the wave period through the reciprocal of the frequency; according to the estimation method, the sea wave texture features can be effectively obtained, and an effective means and method are provided for obtaining high-precision sea wave directions and wave periods.

Description

Technical Field

[0001] The invention belongs to the technical field of ocean remote sensing and relates to a method for estimating wave direction and wave period from an X-band radar image. Background Art

[0002] With the rapid development of ocean development and maritime transportation, accurately obtaining the characteristics of waves (including wave direction and wave period) has become crucial in the fields of marine meteorological forecasting, ship navigation, and offshore structure design. At present, X-band ocean radar has been widely used to obtain wave parameters because of its advantages such as being unaffected by weather conditions and being able to monitor a large area of ​​the ocean surface around the clock. However, traditional wave parameter inversion methods usually face problems such as low accuracy and complex data processing when dealing with complex marine environments. Two-dimensional empirical wavelet transform is an adaptive time-frequency analysis method that can effectively decompose different frequency components in radar images, thereby more accurately extracting the direction and period information of waves. Compared with traditional methods, two-dimensional empirical wavelet transform can better capture the local characteristics of waves, especially when the marine environment is complex and the wave parameters vary widely. Summary of the invention

[0003] The purpose of the present invention is to propose an improved wave direction and wave period estimation method to address the problem of low accuracy in the above-mentioned prior art.

[0004] In order to solve the above technical problems, the present invention proposes a method for estimating wave direction and wave period from an X-band radar image, comprising the following steps:

[0005] Step 1: Radar original wavefront image acquisition. Use the navigation radar to collect a set of original radar wavefront image sequences and store them according to the radar protocol to obtain the radar wavefront image I (r, θ) in polar coordinates, where r is the distance from the sea surface point to the radar and θ is the azimuth;

[0006] Step 2: Select the wave parameter inversion area. In the effective wave texture image screened from the X-band radar wave surface image, select a sub-area I with a distance range of d1*d2 t (r',θ'), t=1,2,…,S as the wave parameter inversion area. t (r',θ') is transformed into Cartesian coordinates to obtain S wave inversion regions I in Cartesian coordinates t (x,y);

[0007] Step 3: Radar sub-image texture feature extraction. For the selected S wave inversion areas I t(x, y) is subjected to a 2D empirical curvelet transform (2D-ECT). Each sub-inversion region divided in step 2 is decomposed into multiple curvelet component images (curvelet coefficients). Image processing includes co-frequency interference, target interference, wave texture feature extraction, etc., to obtain a processed radar wave surface image.

[0008] Step 4: Estimate the wave crest direction. Perform a correlation test on the decomposed curvelet component image in each sub-region in step 3 and the original sub-region, retain the curvelet component image with the greatest correlation with each sub-region, use the empirical ridgelet transform (ERT) to act on the maximum curvelet component image in each sub-region, and then solve the standard deviation of the pixel intensity in each projection direction in each sub-region image. Take the direction with the largest standard deviation as the texture dominant direction β. For all sub-region images of each radar image, take the average of the standard deviations of all sub-images after the ERT method for each projection direction, and take the direction with the largest average standard deviation as the roughly estimated direction α. ​​Select the sub-image whose β is within the predetermined deviation range of α, and take the correlation median of the selected sub-image as the final texture dominant direction γ of a radar image.

[0009] Preferably, a sub-image whose β is within the range of α±10° is selected.

[0010] Step 5: Estimate the peak period. Taking the center point of the sub-region image in the dominant direction of the final texture as the center, expand the sub-image to p×q pixels, solve the periodogram power spectrum density for each column of pixel data in the expanded sub-image, and obtain the frequency corresponding to the maximum value of the power spectrum density of each column. Then, find the maximum value of each column frequency as the final frequency of the radar image, and then estimate the peak period;

[0011] Furthermore, the step 3 may include:

[0012] Step 31: 2D-ECT method is used to analyze the sub-region I t (x, y) for feature extraction, in I t (x, y) is transformed into pseudo polar coordinates. When the frequency adopted on the pseudo grid is (ω 1 ,ω 2 ), assuming that the image input pixels are N*N, the pseudo-polar Fourier transform is:

[0013]

[0014] Where x, y are the coordinate values ​​of the sub-region, N is the image pixel size of the sub-region, and F P It is a pseudo-polar coordinate transformation.

[0015] Step 32: Perform Fourier edge detection on the pseudo-polar Fourier transform, assuming that the number of scales (N s ) and the number of angular sectors (N q ) is known, and the scale boundary is Ω ω ={ω n},n=0,...,N s , the angle boundary is Ω θ ={θ m},η=1,...,N θ , generate an adaptive filter to filter each subband coefficient of the empirical curvelet decomposition to facilitate the necessary processing of each high and low frequency coefficient. Low pass filter φ 1 It is expressed as:

[0016]

[0017] Where β is an arbitrary function in [0,1] that satisfies the following relationship:

[0018]

[0019] In the formula, γ ensures that two consecutive transition intervals do not overlap and must satisfy:

[0020]

[0021] Among them, min() is the minimum function.

[0022] In the Fourier domain, a polar wedge ψ nm (n and m are the scale and angle index respectively) is a radial window W n and a polar window V m The product of s -1, W n It is expressed as:

[0023]

[0024] If N=N s -1 hour:

[0025]

[0026] Polar Window V m It is expressed as:

[0027]

[0028] in, (ω,θ) are polar coordinates in the Fourier plane. Furthermore, the present invention only considers real transforms, ie, the filters are symmetric about the origin.

[0029] Step 33: Establish an empirical curve wave tight support framework, using the first generation empirical curve wave framework, and the Fourier boundary detection method is: the scale and angle are detected independently, first detect the scale Ω ω , and then detect different angle sets for each scale The detection method is the local extreme value method. The first generation of empirical curve wave tight support framework is shown as follows:

[0030]

[0031] Input image I t The detail coefficients of the 2D-ECT at (x,y) are first radially windowed W n , rear pole window V m Definition, the definition of detail coefficient at this time:

[0032]

[0033] The approximate coefficients are expressed as:

[0034]

[0035] Where F 2 is the Fourier transform, F 2 * is the inverse Fourier transform.

[0036] Furthermore, the step 4 may further include:

[0037] Step 41: The detail coefficients decomposed in step 33 are combined with the atomic region I t (x, y) is used for correlation detection, and the detail coefficient with the greatest correlation with the sub-region is retained, which is expressed as U t (x,y).

[0038] Step 42: U t (x, y) is transformed empirically to estimate the wave direction. First, U t (x, y) is transformed by Radon to obtain the projections of the image at different angles. These projections are one-dimensional functions that represent the integral of the image in all directions. The Radon transform is shown in the following formula:

[0039]

[0040] Among them, ρ is the vertical distance from the origin, that is, the shortest distance from the projection line to the origin. θ is the projection angle, the angle with the positive direction of the X-axis, ranging from 0-180°, and the integration is performed along the line perpendicular to the θ direction.

[0041] Step 43: Perform a one-dimensional empirical wavelet transform on R(ρ,θ) in the Radon domain to obtain detail coefficients and approximate coefficients, detect the Fourier boundary by taking the average spectrum obtained by performing a pseudo-polar Fourier transform on the image, and then establish a wavelet framework. Among them, φ 1 As described in step 32 F 2 (φ 1 ), ψ n As shown below:

[0042]

[0043] If n = N-1,

[0044]

[0045] The empirical ridgelet transform is defined as:

[0046]

[0047] Among them, W t E It is a one-dimensional empirical wavelet transform based on B with respect to t. Because the empirical wavelet is defined in the Fourier domain, the above formula is transformed into:

[0048]

[0049] In the formula is the detail coefficient, is the approximate coefficient, F 1,τ is the one-dimensional Fourier transform, is the one-dimensional inverse Fourier transform, and τ is the time.

[0050] Step 44: Calculate the standard deviation of the pixel intensity in each projection direction for the detail coefficient after the empirical ridgelet transform, and the direction with the largest standard deviation is the texture dominant direction β. For all sub-region images of each radar image, take the average of the standard deviations of all sub-images after the ERT method for each projection direction, and take the direction with the largest average standard deviation as the roughly estimated direction α, select the sub-image with β within the range of α±10°, and take the relevant median of the selected sub-image as the final texture dominant direction γ of a radar image;

[0051] Step 4.5: Eliminate the 180° wave direction blur by continuously calculating the correlation coefficient of the final texture dominant direction γ of multiple sub-regions. Assume that two consecutive radar sub-region images are I t and I t+1 , calculate the cross-correlation function R(Δx,Δy):

[0052]

[0053] in and I t and I t+1 The average intensity value of image I t and I t+1 Correlation at offset (Δx, Δy).

[0054] Furthermore, the step 5 may include:

[0055] Step 51: Taking the center point of the sub-image of the final texture dominant direction γ as the center, expand the regional image to a size of p×p pixels. The pixel data of each column of the expanded sub-image pixel matrix is ​​expressed as: S(:, j), j = 1, 2, ..., p. The power spectrum density of the periodogram is solved for each column of pixel data in the sub-image, which is expressed as:

[0056]

[0057] Among them, f s is the sampling frequency, which can be set according to the specific situation. S(n,j) represents the nth pixel value of the jth column of the expanded sub-image.

[0058] Step 52: Find the frequency (f max ), and then from f max The maximum frequency is selected as the frequency of the radar sub-image (f M ), the peak period of this image is:

[0059] T=1 / f M

[0060] Where M is the number of radar sub-images after expansion.

[0061] The beneficial effects of the present invention include: using two-dimensional empirical curvelet transform to extract features from X-band radar images, which can eliminate co-frequency interference and target interference and extract texture features in radar images; using empirical ridgelet transform to estimate the wave crest direction from the image after feature extraction, and using the mutual correlation coefficient to eliminate the 180° blur of the wave direction. Image expansion is performed on the sub-region of the final texture dominant direction, and the peak period is estimated by the periodogram power spectrum method. The advantages of two-dimensional empirical wavelet transform over traditional three-dimensional Fourier transform are processing non-stationary signals, capturing local transformations of wave texture images, performing multi-scale analysis, and robustness in noisy environments. BRIEF DESCRIPTION OF THE DRAWINGS

[0062] Figure 1 4 is a flow chart of a method for estimating a wave peak direction and a wave peak period according to an embodiment of the present invention.

[0063] Figure 2 Shown is the selection of sub-regions of X-band radar images in polar coordinates.

[0064] Figures 3a to 3f Six selected sub-region images are shown respectively.

[0065] Figures 4a to 4j Shown are the curvelet coefficients, Fourier-partitioned images and original sub-images after decomposition using the 2D-ECT method.

[0066] Figure 5a to Figure 5f Shown are the curvelet coefficients with the greatest image correlation in each sub-region.

[0067] Figure 6a Shown are the results after empirical ridgelet transform.

[0068] Figure 6b Shown is the standard deviation of pixel intensity in each direction.

[0069] Figure 7 Shown is a 180° blur of the wave direction.

[0070] Figure 8a Shown is a sub-region image used to estimate the peak period.

[0071] Figure 8b Shown is a sub-region image expansion for estimating the peak period.

[0072] Figure 9a It is the power spectrum density of the periodogram of one column of pixel data in the sub-image.

[0073] Figure 9b It is the frequency corresponding to the maximum power spectral density of all azimuth pixels in a sub-image. DETAILED DESCRIPTION

[0074] The embodiments of the present invention are further described in detail below in conjunction with the accompanying drawings.

[0075] Figure 1 The present invention is a flow chart of a method for estimating wave direction and wave period from an X-band radar image.

[0076] The present invention proposes a method for estimating wave direction and wave period from an X-band radar image, comprising the following steps:

[0077] Step 1: Collect the original radar wavefront image to obtain the X-band radar wavefront image.

[0078] Step 2: Select (divide) an ocean wave parameter inversion region for the X-band radar wave surface image.

[0079] Step 3: Use two-dimensional empirical curvelet transform to extract radar sub-image texture features in the wave parameter inversion area.

[0080] Step 4: Empirical ridgelet transform is performed on the wave parameter inversion area to estimate the wave crest direction.

[0081] Step 5: Estimate the peak period.

[0082] Specifically, step 1 of collecting radar original wavefront images may include collecting a set of original radar wavefront image sequences, namely, X-band radar wavefront images, by using a navigation radar, and storing them according to a radar protocol to obtain a radar wavefront image I(r,θ) in polar coordinates, wherein r is the distance from the sea surface point to the radar, and θ is the azimuth;

[0083] Step 2: Selecting the wave parameter inversion area may include selecting a sub-area I with a distance range of d1*d2 from the effective wave texture image screened from the X-band radar wave surface image. t (r',θ'), t=1,2,…,S as the wave parameter inversion area. t (r',θ') is transformed into Cartesian coordinates to obtain S wave inversion regions I in Cartesian coordinates t (x,y);

[0084] Step 3: Radar sub-image texture feature extraction may include the following steps: t (x, y) is subjected to a 2D empirical curvelet transform (2D-ECT). Each sub-inversion region divided in step 2 is decomposed into multiple curvelet component images (curvelet coefficients). Image processing includes co-frequency interference, target interference, wave texture feature extraction, etc., to obtain a processed radar wave surface image.

[0085] Step 4, estimating the wave crest direction, may include performing a correlation test on the decomposed curvelet component image in each sub-region in step 3 and the original sub-region, retaining the curvelet component image with the greatest correlation with each sub-region, using the empirical ridgelet transform (ERT) to act on the maximum curvelet component image in each sub-region, and then solving the standard deviation of the pixel intensity in each projection direction in each sub-region image, taking the direction with the largest standard deviation as the texture dominant direction β, taking the average of the standard deviations of all sub-images after the ERT method for each projection direction for all sub-region images of each radar image, taking the direction with the largest average standard deviation as the roughly estimated direction (α), selecting the sub-image whose β is within the predetermined deviation range of α, and taking the correlation median of the selected sub-image as the final texture dominant direction (γ) of the radar image;

[0086] Preferably, a sub-image whose β is within the range of α±10° is selected.

[0087] Step 5 estimating the peak period may include taking the center point of the sub-region image in the dominant direction of the final texture as the center, expanding the sub-image to a size of p×q pixels, solving the power spectrum density of the periodogram for each column of pixel data in the expanded sub-image, obtaining the frequency corresponding to the maximum value of the power spectrum density of each column, and then obtaining the maximum value of each column frequency as the final frequency of the radar image, thereby estimating the peak period;

[0088] Furthermore, the step 3 may include:

[0089] Step 31: 2D-ECT method is used to analyze the sub-region I t (x, y) for feature extraction, in I t (x, y) is transformed into pseudo polar coordinates. When the frequency adopted on the pseudo grid is (ω 1 ,ω 2 ), assuming that the image input pixels are N*N, the pseudo-polar Fourier transform is:

[0090]

[0091] Where x, y are the coordinate values ​​of the sub-region, N is the image pixel size of the sub-region, and F P is the pseudo-polar coordinate transformation;

[0092] Step 32: Perform Fourier edge detection on the pseudo-polar Fourier transform, assuming that the number of scales (N s ) and the number of angular sectors (N q ) is known, and the scale boundary is Ω ω ={ω n},n=0,...,N s , the angle boundary is Ω θ ={θ m},η=1,...,N θ , generate an adaptive filter to filter each subband coefficient of the empirical curvelet decomposition to facilitate the necessary processing of each high and low frequency coefficient. Low pass filter φ 1 It is expressed as:

[0093]

[0094] Where β is an arbitrary function in [0,1] that satisfies the following relationship:

[0095]

[0096] In the formula, γ ensures that two consecutive transition intervals do not overlap and must satisfy:

[0097]

[0098] Among them, min() is the minimum function.

[0099] In the Fourier domain, a polar wedge ψ nm (n and m are the scale and angle index respectively) is a radial window W n and a polar window V m The product of s -1, W n It is expressed as:

[0100]

[0101] If N=N s -1 hour:

[0102]

[0103] Polar Window V m It is expressed as:

[0104]

[0105] in, (ω,θ) are polar coordinates in the Fourier plane. In addition, the present invention only considers real number transforms, i.e., these filters are symmetric about the origin; and

[0106] Step 33: Establish an empirical curve wave tight support framework, using the first generation empirical curve wave framework, and the Fourier boundary detection method is: the scale and angle are detected independently, first detect the scale Ω ω , and then detect different angle sets for each scale The detection method is the local extreme value method. The first generation of empirical curve wave tight support framework is shown as follows:

[0107]

[0108] Input image I t The detail coefficients of the 2D-ECT at (x,y) are first radially windowed W n , rear pole window V m Definition, the definition of detail coefficient at this time:

[0109]

[0110] The approximate coefficients are expressed as:

[0111]

[0112] Where F 2 is the Fourier transform, F 2* is the inverse Fourier transform.

[0113] Furthermore, the step 4 may further include:

[0114] Step 41: The detail coefficients decomposed in step 33 are combined with the atomic region I t (x, y) is used for correlation detection, and the detail coefficient with the greatest correlation with the sub-region is retained, which is expressed as U t (x,y);

[0115] Step 42: U t (x, y) is transformed empirically to estimate the wave direction. First, U t (x, y) is transformed by Radon to obtain the projections of the image at different angles. These projections are one-dimensional functions that represent the integral of the image in all directions. The Radon transform is shown in the following formula:

[0116]

[0117] Where ρ is the vertical distance from the origin, that is, the shortest distance from the projection line to the origin. θ is the projection angle, the angle with the positive direction of the X-axis, ranging from 0-180°, and the integration is performed along the line perpendicular to the θ direction;

[0118] Step 43: Perform a one-dimensional empirical wavelet transform on R(ρ,θ) in the Radon domain to obtain detail coefficients and approximate coefficients, detect the Fourier boundary by taking the average spectrum obtained by performing a pseudo-polar Fourier transform on the image, and then establish a wavelet framework. Among them, φ 1 As described in step 32 F 2 (φ 1 ), ψ n As shown below:

[0119]

[0120] If n = N-1,

[0121]

[0122] The empirical ridgelet transform is defined as:

[0123]

[0124] Among them, W t E It is a one-dimensional empirical wavelet transform based on B with respect to t. Because the empirical wavelet is defined in the Fourier domain, the above formula is transformed into:

[0125]

[0126] In the formula is the detail coefficient, is the approximate coefficient, F 1,τ is the one-dimensional Fourier transform, is the one-dimensional inverse Fourier transform, τ is the time;

[0127] Step 44: Calculate the standard deviation of the pixel intensity in each projection direction for the detail coefficient after the empirical ridgelet transform, and the direction with the largest standard deviation is the texture dominant direction β. For all sub-region images of each radar image, take the average of the standard deviations of all sub-images after the ERT method for each projection direction, and take the direction with the largest average standard deviation as the roughly estimated direction α, select the sub-images with β within the range of α±10°, and take the relevant median of the selected sub-images as the final texture dominant direction γ of a radar image; and

[0128] Step 45: Eliminate the 180° wave direction blur by continuously calculating the correlation coefficient of the final texture dominant direction γ of multiple sub-regions. Assume that two consecutive radar sub-region images are I t and I t+1 , calculate the cross-correlation function R(Δx,Δy):

[0129]

[0130] in and I t and I t+1 The average intensity value of image I t and I t+1 Correlation at offset (Δx, Δy).

[0131] Furthermore, the step 5 may include:

[0132] Step 51: Taking the center point of the sub-image of the final texture dominant direction γ as the center, expand the regional image to a size of p×p pixels. The pixel data of each column of the expanded sub-image pixel matrix is ​​expressed as: S(:, j), j = 1, 2, ..., p. The power spectrum density of the periodogram is solved for each column of pixel data in the sub-image, which is expressed as:

[0133]

[0134] Among them, f s is the sampling frequency, which can be set according to the specific situation. S(n,j) represents the nth pixel value of the jth column of the expanded sub-image.

[0135] Step 52: Find the frequency (f max ), and then from f maxThe maximum frequency is selected as the frequency of the radar sub-image (f M ), the peak period of this image is:

[0136] T=1 / f M

[0137] Where M is the number of radar sub-images after expansion.

[0138] like Figures 2 to 9b A specific experimental example shown may include the following steps:

[0139] Step 1: Radar original wavefront image acquisition. Use the navigation radar to collect a set of original radar wavefront image sequences as X-band radar images, and store them according to the radar protocol to obtain the radar wavefront image I (r, θ) in polar coordinates, where r is the distance from the sea surface point to the radar, and θ is the azimuth;

[0140] Step 2: Figure 2 The wave parameter inversion area is selected. In the effective wave texture image screened out by the X-band radar image, a sub-area I with a distance range of 1125m*1125m is selected. t (r',θ'), t=1,2,…,S, in this case S=6, as the wave parameter inversion area. t (r',θ') is transformed into Cartesian coordinates to obtain 6 wave inversion regions I in Cartesian coordinates t (x, y), the pixel size of each sub-region is 150×150, such as Figure 3a , Figure 3b , Figure 3c , Figure 3d , Figure 3e , Figure 3f As shown;

[0141] Step 3: Right Figure 3a , Figure 3b , Figure 3c , Figure 3d , Figure 3e , Figure 3f I 1 -I 6 The feature extraction is performed on the 6 radar sub-area images. The filter used in this experiment is a first-generation empirical curve support framework. First, scale Ω ω Directional decomposition, the number of decompositions is set to 2, and the decomposition method is the local maximum and minimum values ​​on the spectrum; then the angle Ω θ Directional decomposition, the number of decompositions is 8, and the decomposition method is also the local maximum and minimum values ​​on the spectrum. Figure 4a to Figure 4h Among them Figure 3a 2D-EWTC decomposition of neutron image I1, Figure 4a to Figure 4h The sub-image is decomposed into 8 curvelet coefficient components, and Figure 4i is the low-pass filtering of the sub-region image, Figure 4i The white dividing line in the middle is the dividing line in the angle direction. Figure 4j is the radar sub-image I1. Figure 4a As can be seen from the sub-image, Figure 3a The red circle in the figure is the interference noise of the radar sub-image. After 2D-ECT decomposition, the noise interference in the radar image can be separated from the wave texture features. Figure 4c , Figure 4e is the noise interference part in the sub-image, Figure 4a It is the wave texture part in the decomposed radar image;

[0142] Step 4: Estimate the wave crest direction. Perform correlation test on the decomposed curvelet component image in each sub-region in step 3 and the original sub-region, and retain the curvelet component image with the greatest correlation with each sub-region. After 2D-ECT of the sub-region in Figure 3, the curvelet coefficient with the greatest correlation with the atomic region is as follows: Figure 5a to Figure 5f The ERT method is used for Figure 5a to Figure 5f The curvelet coefficient with the greatest correlation with the atomic region is used to estimate the wave crest direction using the transformed third sub-band image. The sub-band image after the ERT method is as follows: Figure 6a As shown in the figure, the horizontal axis is the projection direction θ, and the vertical axis is the distance ρ from the edge image to the origin. The image pixels in the corresponding texture line projection direction will show the largest intensity value change. The standard deviation of the pixel intensity in each projection direction of the transformation result is calculated. The direction with the largest standard deviation is the texture dominant direction (β). The standard deviation of the pixel intensity in each direction of the sub-band image after ERT transformation is shown in Figure 6b As shown. Furthermore, for each sub-image of the radar image, the average standard deviation of all sub-images after ERT is taken for each projection direction, and the direction with the largest average standard deviation is selected as the roughly estimated direction α. ​​Next, the texture dominant direction β is compared with the roughly estimated direction α, and the sub-images with β within the range of α±10° are selected, and the relevant median of the selected sub-images is taken as the final texture dominant direction γ of a radar image. Figure 6b The estimated peak direction is 114°.

[0143] At this time, the final texture dominant direction calculated has a 180° ambiguity problem. The present invention eliminates the 180° ambiguity of the wave direction by calculating the cross-correlation function between multiple consecutive sub-images of the final texture dominant direction γ. The normalized cross-correlation coefficient matrix of two consecutive radar sub-images is as follows: Figure 7As shown in the figure, the waves move in a specific direction in the image, and this movement will appear as a significant peak in the cross-correlation function. The red arrow in the figure indicates the direction of wave propagation, eliminating the problem of solving the 180° ambiguity of the wave direction from the radar image.

[0144] Step 5: Estimate the peak period. Centered on the center point of the sub-region image in the dominant direction of the final texture, expand the sub-image to 256×256 pixels, as shown in Figure 8a , Figure 8b shown. Figure 9a , Figure 9b The process of peak period estimation is shown, where Figure 9(a) is the power spectrum density of the periodogram of a column of pixel data in the sub-image. Figure 9b is the frequency corresponding to the maximum power spectral density of all azimuth pixels in a sub-image, Figure 9a , Figure 9b For example, the maximum frequency of this sub-image is 0.125 Hz and the peak period is 8 s.

[0145] The above-described embodiments merely express the implementation methods of the present invention, but they cannot be understood as limiting the scope of the patent of the present invention. It should be pointed out that for those skilled in the art, several modifications and improvements can be made without departing from the concept of the present invention, which all belong to the protection scope of the present invention.

Claims

1. A method for estimating wave direction and wave period from X-band radar images, characterized in that: Includes steps: Step 1: Use the navigation radar to collect a set of X-band radar wave surface image sequences, and store them according to the radar protocol to obtain the X-band radar wave surface image I (r, θ) in polar coordinates, where r is the distance from the sea surface point to the radar, and θ is the azimuth; Step 2: Filter out the effective wave texture image of the X-band radar wave surface image, and select the sub-region image I with a distance range of d1*d2 from the effective wave texture image. t (r′,θ′), t=1,2,…,S as the wave parameter inversion area image; the S wave parameter inversion area images I in polar coordinates t (r′,θ′) is transformed into Cartesian coordinates to obtain S wave parameter inversion area images I in Cartesian coordinates t (x,y); Step 3: Invert the regional image I for the S wave parameters t (x, y) is subjected to a two-dimensional empirical curvelet transform, wherein each of the divided wave parameter inversion area images is decomposed into a plurality of curvelet component images, and image processing is performed on the curvelet component images to obtain processed radar wave surface images, wherein the image processing includes wave texture feature extraction; Step 4: perform correlation detection on each of the curvelet component images and the corresponding wave parameter inversion area image, retain the curvelet component image with the greatest correlation with each of the wave parameter inversion area images, use the empirical ridgelet transform (ERT) to act on each of the largest curvelet component images, and then solve the standard deviation of the pixel intensity in each projection direction in each of the wave parameter inversion area images, and take the direction with the largest standard deviation as the texture dominant direction (β); for all the wave parameter inversion area images of each radar image, take the average value of the standard deviation of each of the wave parameter inversion area images after ERT for each projection direction, take the direction with the largest average value of the standard deviation as the rough estimation direction (α), select the wave parameter inversion area image whose β is within the predetermined deviation range of α, and take the correlation median of the selected wave parameter inversion area image as the final texture dominant direction (γ) of the effective wave texture image; Step 5: estimating the peak period, including expanding the wave parameter inversion area image to p×q pixels with the center point of the wave parameter inversion area image of the final texture dominant direction (γ) as the center, solving the power spectrum density of each column of pixel data in the expanded wave parameter inversion area image, obtaining the frequency corresponding to the maximum value of the power spectrum density of each column of pixel data, and then calculating the maximum value of each column frequency as the final frequency of the effective wave texture image, and then estimating the peak period based on the final frequency.

2. The method according to claim 1, characterized in that: The predetermined deviation is ±10°.

3. The method according to claim 1, characterized in that: The image processing also includes same-frequency interference extraction and target object interference extraction.

4. The method according to claim 1, characterized in that: The step 3 comprises: Step 31: Use 2D-ECT method to identify sub-region I t (x, y) for feature extraction, in I t (x, y) is transformed into pseudo-polar coordinates. When the frequency used on the pseudo grid is (ω1, ω2), assuming that the image input pixels are N*N, the pseudo-polar Fourier transform is: Where x, y are the coordinate values ​​of the sub-region, N is the image pixel size of the sub-region, and F P is the pseudo-polar coordinate transformation; Step 32: Perform Fourier edge detection on the pseudo-polar Fourier transform, assuming that the number of scales (N s ) and the number of angular sectors (N q ) is known, and the scale boundary is Ω ω ={ω n },n=0,...,N s , the angle boundary is Ω θ ={θ m },η=1,...,N θ , generate an adaptive filter to filter each sub-band coefficient of the empirical curvelet decomposition to facilitate the necessary processing of each high and low frequency coefficient; The adaptive filter is a low-pass filter φ1, expressed as: Where β is an arbitrary function in [0,1] that satisfies the following relationship: In the formula, γ ensures that two consecutive transition intervals do not overlap and must satisfy: Among them, min() is the minimum function; In the Fourier domain, a polar wedge ψ nm (n and m are the scale and angle index respectively) is a radial window W n and a polar window V m The product of s -1, W n It is expressed as: If N=N s -1 hour: Polar Window V m It is expressed as: in, (ω,θ) are polar coordinates in the Fourier plane, where the adaptive filter is symmetric about the origin; and Step 33: Establish an empirical curve wave tight support framework, using the first generation empirical curve wave framework, and the Fourier boundary detection method is: the scale and angle are detected independently, first detect the scale Ω ω , and then detect different angle sets for each scale The detection method is the local extremum method.

5. The method according to claim 4, characterized in that: The step 4 further comprises: Step 41: The detail coefficients decomposed in step 33 are combined with the atomic region I t (x, y) is used for correlation detection, and the detail coefficient with the greatest correlation with the sub-region is retained, which is expressed as U t (x,y); Step 42: U t (x, y) is transformed empirically to estimate the wave direction. First, U t (x, y) performs Radon transform to obtain the projections of the image at different angles. These projections are one-dimensional functions, which represent the integral of the image in all directions. The Radon transform is shown below: Where ρ is the vertical distance from the origin, that is, the shortest distance from the projection line to the origin, θ is the projection angle, the angle with the positive direction of the X-axis, ranging from 0-180°, and the integration is performed along the line perpendicular to the θ direction; Step 43: Perform a one-dimensional empirical wavelet transform on R(ρ,θ) in the Radon domain to obtain detail coefficients and approximate coefficients, detect the Fourier boundary by taking the average spectrum obtained by performing a pseudo-polar Fourier transform on the image, and then establish a wavelet framework. Wherein, φ1 is as shown in F2(φ1) in step 32, ψ n As shown below: If n = N-1, The empirical ridgelet transform is defined as: Among them, W t E It is based on the one-dimensional empirical wavelet transform of B with respect to t; because the empirical wavelet is defined in the Fourier domain, the above formula is transformed into: In the formula is the detail coefficient, is the approximate coefficient, F 1,τ is the one-dimensional Fourier transform, is the one-dimensional inverse Fourier transform, τ is the time; Step 44: Calculate the standard deviation of the pixel intensity in each projection direction for the detail coefficient after the empirical ridgelet transform, and the direction with the largest standard deviation is the texture dominant direction β; for all sub-region images of each radar image, take the average value of the standard deviation of all sub-images after the ERT method for each projection direction, and take the direction with the largest average standard deviation as the roughly estimated direction α, select the sub-image whose β is within the preset deviation range of α, and take the relevant median of the selected sub-image as the final texture dominant direction γ of a radar image; and Step 45: Eliminate the wave direction 180° blur by continuously calculating the mutual correlation coefficient of the final texture dominant direction γ of multiple sub-regions.

6. The method according to claim 5, characterized in that: The step 5 comprises: Step 51: Taking the center point of the sub-image of the final texture dominant direction γ as the center, expand the regional image to a size of p×p pixels. The pixel data of each column of the expanded sub-image pixel matrix is ​​expressed as: S(:, j), j = 1, 2, ..., p. The power spectrum density of the periodogram is solved for each column of pixel data in the sub-image, which is expressed as: Among them, f s is the sampling frequency, which can be set according to the specific situation. S(n,j) represents the nth pixel value of the jth column of the expanded sub-image; Step 52: Find the frequency (f max ), and then from f max The maximum frequency is selected as the frequency of the radar sub-image (f M ), the peak period of this image is: T = 1 / f M , where M is the number of radar sub-images after expansion.

Images