A Doppler center frequency estimation method for SAR inversion of two-dimensional sea surface flow field

Abstract

The present invention relates to a Doppler center frequency estimation method for two-dimensional sea surface flow field SAR inversion, comprising: performing block processing on SAR single-view complex image data; obtaining Doppler center frequency estimation values ​​by using an amplitude method; obtaining Doppler center frequency estimation values ​​by using a phase method; screening out Doppler center frequency estimation values ​​with the best imaging quality as initial estimation results; performing iterative operations on the initial estimation results; completing radial current SAR inversion; and realizing two-dimensional real-time ocean current vector inversion based on the radial current SAR inversion results. The present invention has the following beneficial effects: the present invention can improve the imaging quality and accuracy of Doppler center frequency estimation results; it is an effective supplement to the Doppler center frequency estimation method in the existing SAR ocean current inversion, and can improve the accuracy of SAR ocean current inversion.

Description

Technical Field

[0001] The invention relates to the field of remote sensing image processing, and more specifically, to a Doppler center frequency estimation method for SAR inversion of two-dimensional sea surface flow fields. Background Art

[0002] The Doppler frequency shift of the scattering echo of Synthetic Aperture Radar (SAR) has the ability to reflect the dynamic characteristics of the sea surface and has become a key sea surface dynamic remote sensing parameter. In the monitoring and mapping of small and medium-scale ocean currents, the Doppler frequency anomaly (Doppler Centroid Anomaly) algorithm is one of the mainstream methods. The essential requirement of this algorithm is to accurately separate the Doppler frequency shift contribution caused by non-ocean current motion from the SAR Doppler center frequency. Therefore, accurately estimating the SAR Doppler center frequency is crucial to improving the accuracy of ocean current SAR inversion.

[0003] At present, only the SAR data of Sentinel-1 satellite and ENVISAT satellite bandwidth mode provide Doppler center parameters. Most other SAR data do not provide this key parameter, such as Gaofen-3, RADARSAT and TerraSAR-X. For SAR data that do not provide the Doppler center frequency, the Doppler center must be estimated before it can be used for SAR ocean current inversion. There are two major types of Doppler center frequency estimation methods: one is calculated based on satellite orbit parameters, and the other is estimated from SAR echo signals. In practical applications, the orbit parameters provided by satellites are usually not accurate enough. Therefore, it is necessary to estimate the Doppler center frequency from SAR echo data, including amplitude and phase methods.

[0004] The amplitude method refers to the use of the Fourier transform duality between the signal azimuth power spectrum and the correlation function to calculate the mutual correlation coefficient of the signal in the time domain, thereby realizing the target Doppler center frequency estimation, such as the Correlation Doppler Estimation (CDE) and the Sign Doppler Estimation (SDE). Unlike the amplitude method, the core of the phase method is based on the symmetry principle of the Doppler spectrum of the radar beam in the azimuth direction, and uses different weighting functions to perform convolution operations with the Doppler power spectrum to obtain the Doppler center frequency from the convolution result. Such as the Energy Balance (EB), Match Correlation (MC) and Optimal estimation (OP). It should be noted that the average value and imaging quality of the estimation results of different Doppler center frequency estimation methods are different. In addition, different Doppler center frequency estimation algorithms also have different imaging qualities for SAR images with different scanning modes and spatial resolutions. Summary of the invention

[0005] The purpose of the present invention is to address the deficiencies of the prior art and to propose a Doppler center frequency estimation method for SAR inversion of two-dimensional sea surface flow fields.

[0006] In a first aspect, a Doppler center frequency estimation method for two-dimensional sea surface flow SAR inversion is provided, comprising:

[0007] Step 1: Acquire SAR single-view complex image data, and perform block processing on the SAR single-view complex image data;

[0008] Step 2: Obtain a Doppler center frequency estimation value using an amplitude method; the amplitude method includes a correlation Doppler method and a sign Doppler method;

[0009] Step 3, using a phase method to obtain a Doppler center frequency estimation value; the phase method includes an energy balance method, a matching correlation method, and an optimal estimation method;

[0010] Step 4: Calculate the signal-to-noise ratio of the Doppler center frequency estimation values ​​obtained by different methods, and select the Doppler center frequency estimation value with the best imaging quality as the initial estimation result, which is recorded as

[0011] Step 5: Based on the maximum likelihood principle, iteratively calculate the initial estimation results;

[0012] Step 6: Remove the predicted Doppler frequency shift f caused by satellite attitude geo , get the Doppler frequency shift anomaly f Dca, then subtract the Doppler center frequency contributed by wind and wave deviation to complete the radial current SAR inversion;

[0013] Step 7: Based on the radial current SAR inversion results, use the two-dimensional Ekman flow inversion method to achieve two-dimensional real-time ocean current vector inversion.

[0014] Preferably, step 2 comprises:

[0015] Step 2.1: Use the correlation Doppler method to estimate the Doppler center frequency, expressed as:

[0016]

[0017] Among them, x(m)=x′(mT) is defined as two random processes, that is, two adjacent sub-image blocks with an interval of η, and their correlation coefficient is recorded as N is the azimuth length of each image block, x * (m) is the conjugate complex number of x(m), η is the delay; arg is the phase angle function, T is the time interval between two adjacent sampling signals, which is equal to 1 / PRF; PRF is the pulse frequency;

[0018] Step 2.2, using the symbol Doppler method to estimate the Doppler center frequency;

[0019] Step 2.3: Use the sliding window method to traverse each image block from left to right and from top to bottom.

[0020] Preferably, step 2.2 comprises:

[0021] Step 2.2.1, define the sign function sv, and obtain the sign correlation expression of the single-view complex image, which is expressed as:

[0022]

[0023] Among them, N y and N x are the azimuth and distance lengths of each image block, respectively. sx(i+k,j) and sx(i,j) are two local windows / samples with an interval of k.

[0024] Step 2.2.2: Derive the normalized correlation coefficient ρ xy (η) and the complex correlation coefficient of the single-view complex image It is expressed as:

[0025]

[0026] in, is the correlation coefficient between the real part I and the imaginary part Q, and R sx,sy (k) have the same argument;

[0027] Step 2.2.3: Substitute the multiple correlation coefficient into The Doppler center frequency is obtained

[0028] Preferably, step 3 comprises:

[0029] Step 3.1, perform Fourier transform on each image, convert the time domain image into the frequency domain image, and calculate the azimuth Doppler spectrum;

[0030] Step 3.2, select different weighting functions and perform convolution operation with the azimuth Doppler spectrum;

[0031] Step 3.3: Search the energy peak or zero-crossing point of the convolution operation result from a pulse transmission frequency by integration. The corresponding Doppler frequency is the Doppler center frequency. Each sub-image returns a Doppler center frequency estimate f Dc_block ; The estimation results of energy balance method, matching correlation method and optimal estimation method are recorded as and

[0032] Preferably, step 3.2 comprises:

[0033] Step 3.2.1, select three weighting functions: energy balance, correlation matching and optimal estimation, denoted as B 1 (f), B 2 (f) and B 3 (f), the formula is;

[0034]

[0035] Among them, p(f) is the azimuth Doppler spectrum, and E(p(f)) is the power spectrum density of the signal;

[0036] Step 3.2.2: Use a sliding window to traverse each image block from left to right and from top to bottom, and perform convolution operation on the weighted function and the azimuthal spectrum.

[0037] Preferably, in step 3.3, the integral equation for searching for the zero-crossing point of the convolution operation result is:

[0038]

[0039] Among them, F(φ) is the integral function for searching the zero-crossing point, and B(f) corresponds to the three weighting functions in step 3.2.1. Thus, the estimation results of the three phase methods of energy balance method, matching correlation method and optimal estimation method are obtained, which are recorded as and

[0040] Preferably, in step 4, the signal-to-noise ratio calculation formula is:

[0041]

[0042] Among them, S 0 and S 1 are the main peak and side lobe of the SAR signal in azimuth direction respectively; by comparing the signal-to-noise ratios of the five Doppler center frequency estimation results, the estimation result with the best signal-to-noise ratio is selected as the initial estimation value before iterative calculation, denoted as

[0043] Preferably, in step 5, the iterative calculation formula is:

[0044]

[0045] -PRF / 2≤Fa≤PRF / 2

[0046] Where Fa is the real Doppler center, F ′ (ψ) is the integral function used to search for the zero-crossing point of the azimuth spectrum in the iterative calculation, f Dc is the Doppler center frequency after iterative calculation, s2 / σ2 is the variance; the Doppler center frequency corresponding to the zero crossing point is equal to

[0047] Preferably, step 6 comprises:

[0048] Step 6.1: Calculate the Doppler frequency shift f contributed by the satellite attitude geo :

[0049] f geo =d 0 +d 1 (t s -t 0 )+d 2 (t s -t 0 ) 2 +d 3 (t s -t 0 ) 3 +d 4 (t s -t 0 ) 4

[0050] f Dca =f Dc -f geo

[0051] Among them, d i(i=0,1,2,3,4) is the Doppler coefficient, t s and t 0 They are slope range time and standard slope range time;

[0052] Step 6.2: Use the CDOP geophysical model function to estimate the Doppler frequency shift f contributed by wind and wave deviation ww and change it from f Dca Subtract it from the current to obtain the Doppler frequency shift f contributed by the ocean current osc :

[0053]

[0054] Among them, u 10 and are wind speed and relative wind direction, θ and pol are radar incident angle and polarization mode, respectively;

[0055] Step 6.3, calculate the radial flow velocity U:

[0056]

[0057] In the formula, k r is the wave number of radar electromagnetic wave.

[0058] In a second aspect, a Doppler center frequency estimation system for two-dimensional sea surface flow field SAR inversion is provided, which is used to execute any method described in the first aspect, including:

[0059] A first acquisition module is used to acquire SAR single-view complex image data and perform block processing on the SAR single-view complex image data;

[0060] A second acquisition module is used to acquire a Doppler center frequency estimation value by using an amplitude method; the amplitude method includes a correlation Doppler method and a symbol Doppler method;

[0061] A third acquisition module is used to obtain a Doppler center frequency estimation value by using a phase method; the phase method includes an energy balance method, a matching correlation method and an optimal estimation method;

[0062] The calculation module is used to calculate the signal-to-noise ratio of the Doppler center frequency estimation values ​​obtained by different methods, and select the Doppler center frequency estimation value with the best imaging quality as the initial estimation result, which is recorded as

[0063] Iteration module, used to perform iterative calculations on the initial estimation results based on the maximum likelihood principle;

[0064] The first inversion module is used to remove the predicted Doppler frequency shift f caused by the satellite attitude. geo , get the Doppler frequency shift anomaly f Dca, then subtract the Doppler center frequency contributed by wind and wave deviation to complete the radial current SAR inversion;

[0065] The second inversion module is used to realize two-dimensional real-time ocean current vector inversion based on the radial ocean current SAR inversion result and using the two-dimensional Ekman flow inversion method.

[0066] The beneficial effects of the present invention are as follows: the present invention proposes a Doppler center frequency estimation method for SAR inversion of two-dimensional sea surface flow fields, which can improve the imaging quality and accuracy of Doppler center frequency estimation results; it is an effective supplement to the existing Doppler center frequency estimation method in SAR ocean current inversion, can improve the accuracy of SAR ocean current inversion, and has important practical application significance. BRIEF DESCRIPTION OF THE DRAWINGS

[0067] Figure 1 It is a technical flow chart of the Doppler center frequency estimation method for two-dimensional sea surface flow field SAR inversion of the present invention;

[0068] Figure 2 A schematic diagram of block size and step length in an embodiment of the present invention;

[0069] Figure 3 Schematic diagram of the results of the correlation Doppler method and the sign Doppler method in an embodiment of the present invention;

[0070] Figure 4 A schematic diagram of an azimuth Doppler spectrum and a weighting function in an embodiment of the present invention;

[0071] Figure 5 Schematic diagram of estimation results of three phase methods and optimization method in an embodiment of the present invention;

[0072] Figure 6 Schematic diagram of radial ocean current inversion results in an embodiment of the present invention;

[0073] Figure 7 It is a schematic diagram of the two-dimensional real-time ocean current inversion result in an embodiment of the present invention. DETAILED DESCRIPTION

[0074] The present invention is further described below in conjunction with embodiments. The description of the following embodiments is only used to help understand the present invention. It should be noted that for ordinary persons in the art, without departing from the principle of the present invention, the present invention can also be modified in some ways, and these improvements and modifications also fall within the scope of protection of the claims of the present invention.

[0075] Embodiment 1:

[0076] In order to solve the problems of the prior art, Example 1 of the present application provides a Doppler center frequency estimation method for two-dimensional sea surface flow field SAR inversion. The initial estimation result with the best imaging quality is obtained from a variety of Doppler center frequency estimation methods through signal-to-noise ratio, and then iterative calculation is performed to optimize it to obtain the Doppler center frequency estimation result with the smallest variance. This can improve the SAR inversion accuracy of radial currents and two-dimensional currents, and all calculations are implemented in the MATLAB language.

[0077] Specifically, the technical process of this method is as follows Figure 1 As shown, including:

[0078] Step 1: Acquire SAR single-view complex image data, and perform block processing on the SAR single-view complex image data.

[0079] Specifically, according to the spatial resolution of the SAR complex image, that is, the number of pixels in the range and azimuth directions, the image is divided into M×N sub-images. In addition, the image block size principle in step 1 is: the resolution of each image block is as close as possible to the homogeneous sea surface of 1 km×1 km.

[0080] The code for image segmentation is:

[0081]

[0082] For example, assuming that the spatial resolution of the SAR complex image is 4.3–4.9 m / s×1.7–3.6 m / s (azimuth×range), the number of pixels in the azimuth and range directions are 36889 and 19098 respectively. The block size and step size are set to 256×128 and 8×4 respectively, a total of 306×154 sub-images are divided, and the spatial resolution of the generated Doppler center frequency image is close to 1 km×1 km.

[0083] Step 2: Obtain the Doppler center frequency estimation value using the amplitude method; the amplitude method includes the correlation Doppler method and the sign Doppler method.

[0084] Step 2 includes:

[0085] Step 2.1: Use the correlation Doppler method to estimate the Doppler center frequency, expressed as:

[0086]

[0087] Among them, x(m)=x′(mT) is defined as two random processes, that is, two adjacent sub-image blocks with an interval of η, and their correlation coefficient is recorded as N is the azimuth length of each image block, x *(m) is the conjugate complex number of x(m), η is the delay; arg is the phase angle function, T is the time interval between two adjacent sampling signals, which is equal to 1 / PRF; PRF is the pulse frequency, which can be obtained from the metadata of the SAR image.

[0088] Since most SAR sensors are sampled at a rate slightly higher than the Nyquist sampling rate, as η increases, the correlation coefficient will quickly approach 0. Therefore, in actual calculations, η is generally taken as 1. The Doppler center frequency estimation is achieved by continuously estimating the single cross-correlation coefficient:

[0089]

[0090] Step 2.2, using the symbol Doppler method to estimate the Doppler center frequency;

[0091] Step 2.2 includes:

[0092] Step 2.2.1. Since the SAR echo signal is a quasi-Gaussian process, the "arcsine law" of the Gaussian process can be used to calculate the Doppler center frequency shift of the SAR data. Based on the Gaussian property, assuming that x(t) and y(t) are two real Gaussian processes with zero mean, a symbol function can be defined:

[0093]

[0094] Taking the single-view complex image (SLC) as an example, the complex calculation process can be expressed as h(k)=I(k)+jQ(k). Where I and Q are both real functions, and the symbol correlation of the SLC image can be expressed as:

[0095]

[0096] Among them, N y and N x They are the azimuth and distance lengths of each image block, sx(i+k,j) and sx(i,j) are two local windows / samples with an interval of k; usually equal to 1.

[0097] Step 2.2.2: The normalized correlation coefficient ρ can be quantitatively derived using the arcsine xy (η) and the complex correlation coefficient of the single-view complex image It is expressed as:

[0098]

[0099] in, is the correlation coefficient between the real part I and the imaginary part Q, and R sx,sy(k) have the same argument;

[0100] Step 2.2.3: Substitute the multiple correlation coefficient into The Doppler center frequency is obtained

[0101] Step 2.3: Use the sliding window method to traverse each image block from left to right and from top to bottom, and store the results in step 1. variable.

[0102] Step 3: Obtain the Doppler center frequency estimation value using a phase method; the phase method includes an energy balance method, a matching correlation method, and an optimal estimation method.

[0103] Step 3 includes:

[0104] Step 3.1: Perform Fourier transform on each image, convert the time domain image into the frequency domain image, and calculate the azimuth Doppler spectrum.

[0105] For example, the fft function in the matlab language is used to perform Fourier transform on the sub-image, from the time domain to the frequency domain; and the fftshift function is used to transfer the Fourier transform result to the center of the image, and after calculating the absolute value and square calculation, the azimuth Doppler frequency shift can be obtained, which is a left-right symmetrical azimuth energy spectrum with Gaussian distribution.

[0106] Step 3.2: Select different weighting functions and perform convolution operation with the azimuth Doppler spectrum.

[0107] Step 3.2 includes:

[0108] Step 3.2.1, select three weighting functions: energy balance, correlation matching and optimal estimation, denoted as B 1 (f), B 2 (f) and B 3 (f), the formula is;

[0109]

[0110] Among them, p(f) is the azimuth Doppler spectrum, and E(p(f)) is the power spectrum density of the signal;

[0111] Step 3.2.2: Use a sliding window to traverse each image block from left to right and from top to bottom, and perform convolution operation on the weighted function and the azimuthal spectrum.

[0112] Step 3.3: Search the energy peak or zero-crossing point of the convolution operation result from a pulse transmission frequency by integration. The corresponding Doppler frequency is the Doppler center frequency. Each sub-image returns a Doppler center frequency estimate fDc_block ; The estimation results of energy balance method, matching correlation method and optimal estimation method are recorded as and

[0113] In step 3.3, the integral equation used to search for the zero-crossing point of the convolution operation result is:

[0114]

[0115] Among them, F(φ) is the integral function for searching the zero-crossing point, B(f) corresponds to the three weighting functions in step 3.2.1; the Doppler center frequency f of each sub-image is adjusted by sliding window. Dc_block Stored in step 1 Variable. The estimation results of the three phase methods, energy balance method, matching correlation method and optimal estimation method, are obtained and recorded as and

[0116] Step 4: Calculate the signal-to-noise ratio of the Doppler center frequency estimation values ​​obtained by different methods, and select the Doppler center frequency estimation value with the best imaging quality as the initial estimation result, which is recorded as

[0117] In step 4, the signal-to-noise ratio calculation formula is:

[0118]

[0119] Among them, S 0 and S 1 are the main peak and side lobe of the SAR signal in azimuth direction respectively; by comparing the signal-to-noise ratios of the five Doppler center frequency estimation results, the estimation result with the best signal-to-noise ratio is selected as the initial estimation value before iterative calculation, denoted as

[0120] Step 5: Based on the maximum likelihood principle, perform iterative operations on the initial estimation results. When the variance between the result and the initial estimation result is the smallest, that is, the Cramér–Rao bound is satisfied, the iteration is stopped, and the estimation result with the smallest variance with the simulated true Doppler center (-PRF / 2≤Fa≤PRF / 2) is obtained, which is recorded as f Dc .

[0121] In step 5, the iterative calculation formula is:

[0122]

[0123] -PRF / 2≤Fa≤PRF / 2

[0124] Where Fa is the true Doppler center, F′(φ) is the integral function used to search for the zero-crossing point of the azimuth spectrum in the iterative calculation, and fDc is the Doppler center frequency after iterative calculation, s2 / σ2 is the variance; it should be noted that the Doppler center frequency corresponding to the zero crossing point here is equal to

[0125] Step 6: Remove the predicted Doppler frequency shift f caused by satellite attitude geo , get the Doppler frequency shift anomaly f Dca , then subtract the Doppler center frequency contributed by wind and wave deviation to complete the radial current SAR inversion.

[0126] Step 6 includes:

[0127] Step 6.1: Calculate the Doppler frequency shift f contributed by the satellite attitude geo :

[0128] f geo =d 0 +d 1 (t s -t 0 )+d 2 (t s -t 0 ) 2 +d 3 (t s -t 0 ) 3 +d 4 (t s -t 0 ) 4

[0129] f Dca =f Dc -f geo

[0130] Among them, d i (i=0,1,2,3,4) is the Doppler coefficient, t s and t 0 They are slant range time and standard slant range time respectively. The above parameters can be obtained from SAR metadata;

[0131] Step 6.2: Use the CDOP geophysical model function to estimate the Doppler frequency shift f contributed by wind and wave deviation ww and change it from f Dca Subtract it from the current to obtain the Doppler frequency shift f contributed by the ocean current osc :

[0132]

[0133] Among them, u 10 and are wind speed and relative wind direction, θ and pol are radar incident angle and polarization mode, respectively;

[0134] Step 6.3, calculate the radial flow velocity U:

[0135]

[0136] In the formula, k r is the wave number of radar electromagnetic wave.

[0137] Step 7: Based on the radial current SAR inversion results, the radial velocity is decomposed into longitudinal and latitudinal components using the two-dimensional Ekman flow method to achieve two-dimensional real-time current inversion. The formula is:

[0138]

[0139] in, θ is equal to 30.75°, f ρ0 is the density of seawater, equal to 1.02×10 3 kg / m 3 It can be obtained that the solution of the velocity component is composed of the Ekman depth D, the component of the wind stress τ x and τ y Control, while τ x and τ y It can be obtained by calculating the wind field inverted by SAR. Therefore, it is only necessary to solve the value of Ekman depth D to invert the two-dimensional sea surface vector flow field. The relationship between vector modulus mag and u, v can be expressed as:

[0140]

[0141] Assuming that mag is equal to the radial flow velocity U, combining the above three formulas, we can get:

[0142]

[0143] Among them, all the parameters in the above formulas are known variables, which can be used to solve the Ekman depth D and then realize the inversion of two-dimensional real-time ocean currents.

[0144] Embodiment 2:

[0145] Based on Example 1, Example 2 of the present application aims at the average value and imaging quality difference of different Doppler center frequency estimation methods, and selects a Sentinel-1SM mode HH polarization single-view complex image of the east coast of the United States on September 29, 2016 to perform Doppler center frequency estimation, error correction, radial current and two-dimensional current inversion. Figure 1 As shown, the method of this embodiment includes the following steps:

[0146] Step 1: The selected Sentinel-1SM mode VV polarization single-view complex image is processed in blocks. Its spatial resolution is 4.3–4.9m / s×1.7–3.6m / s (azimuth×range), and the number of pixels in azimuth and range are 36889 and 19098 respectively. The ocean area of ​​1×1km can usually be regarded as a homogeneous target. The block size and step size can be set to 256×128 and 8×4 respectively, and a total of 306×154 sub-images are divided. The spatial resolution of the generated Doppler center frequency image is close to 1×1km. The schematic diagram of the block size and step size can be seen in Figure 2 ;

[0147] Step 2: Use the two classic amplitude methods, correlation Doppler method and symbol Doppler method, to calculate the Doppler center frequency of each sub-image. The results are as follows: Figure 3 Compared with the correlation Doppler method, the advantage of the symbol Doppler method is that it does not give a greater weight to bright targets and is insensitive to non-uniform scenes with large changes;

[0148] Step 3: The difference between the phase method and the amplitude method is that the former belongs to the frequency domain method, while the latter belongs to the time domain method. Therefore, it is necessary to perform Fourier transform on each single-view complex sub-image, convert the time domain image to the frequency domain image, and calculate the azimuth Doppler spectrum, such as Figure 4 (a)

[0149] Step 4: There are three commonly used phase methods, namely, energy balance method, matching correlation method and optimal estimation method. The difference lies in the use of three different weighting functions to perform convolution operations with the azimuth Doppler spectrum. The weighting function is as follows: Figure 4 As shown in (b)-(d);

[0150] Step 5. Since the azimuth Doppler spectrum is even symmetrical and Gaussian distributed, the Doppler frequency corresponding to the energy peak (slope crossing zero point) is the Doppler center frequency. Here, we can search for the slope crossing zero point, which corresponds to the Doppler center frequency, from the convolution operation result of the weighting function and the azimuth Doppler spectrum. The above calculation can be implemented using the interp1 function in matlab. The estimation results of the three phase methods are as follows: Figure 5 In and As shown;

[0151] Step 6: From steps 1 to 6, the results of 5 commonly used Doppler center frequency methods can be obtained. By calculating the signal-to-noise ratio of the estimation results of different methods, the Doppler center frequency estimation method with the best imaging quality is selected as the initial estimation result. The Doppler center frequency estimation method with the highest signal-to-noise ratio here is the optimal estimation method

[0152] Step 7: Based on the maximum likelihood principle, the initial estimation result is optimized by iterative calculation, so that the estimation result is closer to the simulated real Doppler center, corresponding to Figure 5 f Dc ;

[0153] Step 8: Read the Doppler coefficient, slant range time, and standard slant range time from the SAR metadata, calculate the predicted Doppler shift contributed by the satellite attitude, and calculate the Doppler frequency shift from f Dc Then, the ECWMF reanalysis wind field is used as the wind vector constant input to the CDOP geophysical model to estimate the Doppler shift contributed by wind and waves. After removing the wind and wave deviation, the Doppler shift contributed by the current can be obtained. According to the linear relationship between the Doppler shift and the radial current velocity, the radial current velocity can be inverted. The result is as follows: Figure 6 As shown;

[0154] Step 9: Based on the radial current inversion results, the recently proposed two-dimensional Ekman flow inversion method can be used to achieve two-dimensional current inversion. The results are as follows: Figure 7 shown.

[0155] After calculation using the above-mentioned optimized Doppler center frequency estimation method, the Doppler center frequency, a key current monitoring parameter, can be obtained from the single-view complex image, which can improve the accuracy of radial current and two-dimensional current inversion and is very practical.

[0156] It should be noted that the parts in this embodiment that are the same or similar to those in Embodiment 1 can be referenced to each other and will not be described in detail in this application.

[0157] Embodiment 3:

[0158] Based on Example 1, Example 3 of the present application provides a Doppler center frequency estimation system for two-dimensional sea surface flow field SAR inversion, including:

[0159] A first acquisition module is used to acquire SAR single-view complex image data and perform block processing on the SAR single-view complex image data;

[0160] A second acquisition module is used to acquire a Doppler center frequency estimation value by using an amplitude method; the amplitude method includes a correlation Doppler method and a symbol Doppler method;

[0161] A third acquisition module is used to obtain a Doppler center frequency estimation value by using a phase method; the phase method includes an energy balance method, a matching correlation method and an optimal estimation method;

[0162] The calculation module is used to calculate the signal-to-noise ratio of the Doppler center frequency estimation values ​​obtained by different methods, and select the Doppler center frequency estimation value with the best imaging quality as the initial estimation result, which is recorded as

[0163] Iteration module, used to perform iterative calculations on the initial estimation results based on the maximum likelihood principle;

[0164] The first inversion module is used to remove the predicted Doppler frequency shift f caused by the satellite attitude. geo , get the Doppler frequency shift anomaly f Dca , then subtract the Doppler center frequency contributed by wind and wave deviation to complete the radial current SAR inversion;

[0165] The second inversion module is used to realize two-dimensional real-time ocean current vector inversion based on the radial ocean current SAR inversion result and using the two-dimensional Ekman flow inversion method.

[0166] Specifically, the system provided in this embodiment is a system corresponding to the method provided in Embodiment 1. Therefore, the parts in this embodiment that are the same or similar to those in Embodiment 1 can be referenced to each other and will not be repeated in this application.

Claims

1. A Doppler center frequency estimation method for SAR inversion of two-dimensional sea surface flow field, characterized in that: include: Step 1: Acquire SAR single-view complex image data, and perform block processing on the SAR single-view complex image data; Step 2: Obtain a Doppler center frequency estimation value using an amplitude method; the amplitude method includes a correlation Doppler method and a sign Doppler method; Step 2 includes: Step 2.1: Use the correlation Doppler method to estimate the Doppler center frequency, expressed as: Among them, x(m)=x′(mT) is defined as two random processes, that is, two adjacent sub-image blocks with an interval of η, and their correlation coefficient is recorded as N is the azimuth length of each image block, x * (m) is the conjugate complex number of x(m), η is the delay; arg is the phase angle function, T is the time interval between two adjacent sampling signals, which is equal to 1 / PRF; PRF is the pulse frequency; Step 2.2, using the symbol Doppler method to estimate the Doppler center frequency; Step 2.3: Use the sliding window method to traverse each image block from left to right and from top to bottom. Step 3, using a phase method to obtain a Doppler center frequency estimation value; the phase method includes an energy balance method, a matching correlation method, and an optimal estimation method; Step 4: Calculate the signal-to-noise ratio of the Doppler center frequency estimation values ​​obtained by different methods, and select the Doppler center frequency estimation value with the best imaging quality as the initial estimation result, which is recorded as Step 5: Based on the maximum likelihood principle, iteratively calculate the initial estimation results; Step 6: Remove the predicted Doppler frequency shift f caused by satellite attitude geo , get the Doppler frequency shift anomaly f Dca , then subtract the Doppler center frequency contributed by wind and wave deviation to complete the radial current SAR inversion; Step 7: Based on the radial current SAR inversion results, use the two-dimensional Ekman flow inversion method to achieve two-dimensional real-time ocean current vector inversion.

2. The Doppler center frequency estimation method for two-dimensional sea surface flow field SAR inversion according to claim 1 is characterized in that: Step 2.2 includes: Step 2.2.1, define the sign function sv, and obtain the sign correlation expression of the single-view complex image, which is expressed as: Among them, N y and N x are the azimuth and distance lengths of each image block, respectively. sx(i+k,j) and sx(i,j) are two local windows / samples with an interval of k. Step 2.2.2: Derive the normalized correlation coefficient ρ xy (η) and the complex correlation coefficient of the single-view complex image It is expressed as: in, is the correlation coefficient between the real part I and the imaginary part Q, and R sx,sy (k) have the same argument; Step 2.2.3: Substitute the multiple correlation coefficient into The Doppler center frequency is obtained 3. The Doppler center frequency estimation method for two-dimensional sea surface flow field SAR inversion according to claim 2 is characterized in that: Step 3 includes: Step 3.1, perform Fourier transform on each image, convert the time domain image into the frequency domain image, and calculate the azimuth Doppler spectrum; Step 3.2, select different weighting functions and perform convolution operation with the azimuth Doppler spectrum; Step 3.3: Search the energy peak or zero-crossing point of the convolution operation result from a pulse transmission frequency by integration. The corresponding Doppler frequency is the Doppler center frequency. Each sub-image returns a Doppler center frequency estimate f Dc_block ; The estimation results of energy balance method, matching correlation method and optimal estimation method are recorded as and 4. The Doppler center frequency estimation method for two-dimensional sea surface flow field SAR inversion according to claim 3 is characterized in that: Step 3.2 includes: Step 3.2.1, select three weighting functions: energy balance, correlation matching and optimal estimation, which are represented by B1(f), B2(f) and B3(f) respectively. The formula is: Among them, p(f) is the azimuth Doppler spectrum, and E(p(f)) is the power spectrum density of the signal; Step 3.2.2: Use a sliding window to traverse each image block from left to right and from top to bottom, and perform convolution operation on the weighted function and the azimuthal spectrum.

5. The Doppler center frequency estimation method for two-dimensional sea surface flow field SAR inversion according to claim 4 is characterized in that: In step 3.3, the integral equation used to search for the zero-crossing point of the convolution operation result is: Among them, F(φ) is the integral function for searching the zero-crossing point, and B(f) corresponds to the three weighting functions in step 3.2.

1. Thus, the estimation results of the three phase methods, namely, the energy balance method, the matching correlation method and the optimal estimation method, are obtained, which are respectively recorded as and 6. The Doppler center frequency estimation method for two-dimensional sea surface flow field SAR inversion according to claim 5 is characterized in that: In step 4, the signal-to-noise ratio calculation formula is: Among them, S0 and S1 are the main peak and side lobe of the SAR signal in azimuth, respectively. By comparing the signal-to-noise ratios of the five Doppler center frequency estimation results, the estimation result with the best signal-to-noise ratio is selected as the initial estimation value before iterative calculation, which is recorded as 7. The Doppler center frequency estimation method for two-dimensional sea surface flow field SAR inversion according to claim 6 is characterized in that: In step 5, the iterative calculation formula is: Where Fa is the real Doppler center, F ′ (φ) is the integral function used to search for the zero-crossing point of the azimuth spectrum in the iterative calculation, f Dc is the Doppler center frequency after iterative calculation, s2 / σ2 is the variance; the Doppler center frequency corresponding to the zero crossing point is equal to Fa- 8. The Doppler center frequency estimation method for two-dimensional sea surface flow field SAR inversion according to claim 7, characterized in that: Step 6 includes: Step 6.1: Calculate the Doppler frequency shift f contributed by the satellite attitude geo : f geo =d0+d1(t s -t0)+d2(t s -t0) 2 +d3(t s -t0) 3 +d4(t s -t0) 4 f Dca =f Dc -f geo Among them, d i (i=0,1,2,3,4) is the Doppler coefficient, t s and t0 are the slant range time and the standard slant range time respectively; Step 6.2: Use the CDOP geophysical model function to estimate the Doppler frequency shift f contributed by wind and wave deviation ww and change it from f Dca Subtract it from the current to obtain the Doppler frequency shift f contributed by the ocean current osc : Among them, u 10 and are wind speed and relative wind direction, θ and pol are radar incident angle and polarization mode, respectively; Step 6.3, calculate the radial flow velocity U: In the formula, k r is the radar electromagnetic wave number.

9. A Doppler center frequency estimation system for two-dimensional sea surface flow field SAR inversion, characterized in that: Used to perform the method according to any one of claims 1 to 8, comprising: A first acquisition module is used to acquire SAR single-view complex image data and perform block processing on the SAR single-view complex image data; A second acquisition module is used to acquire a Doppler center frequency estimation value by using an amplitude method; the amplitude method includes a correlation Doppler method and a symbol Doppler method; A third acquisition module is used to obtain a Doppler center frequency estimation value by using a phase method; the phase method includes an energy balance method, a matching correlation method and an optimal estimation method; The calculation module is used to calculate the signal-to-noise ratio of the Doppler center frequency estimation values ​​obtained by different methods, and select the Doppler center frequency estimation value with the best imaging quality as the initial estimation result, which is recorded as Iteration module, used to perform iterative calculations on the initial estimation results based on the maximum likelihood principle; The first inversion module is used to remove the predicted Doppler frequency shift f caused by the satellite attitude. geo , get the Doppler frequency shift anomaly f Dca , then subtract the Doppler center frequency contributed by wind and wave deviation to complete the radial current SAR inversion; The second inversion module is used to realize two-dimensional real-time ocean current vector inversion based on the radial ocean current SAR inversion result and using the two-dimensional Ekman flow inversion method.

Images