Sea surface wind speed inversion method and system based on CYGNSS and Willoughby models

Abstract

The invention discloses a sea surface wind speed inversion method and system based on a CYGNSS and Willoughby model, and belongs to the technical field of ocean remote sensing. The method comprises the following steps: firstly, extracting a cyclone parameter from IBTrACS data through a Willoughby model, and encrypting the cyclone parameter to 1 second for sampling by utilizing segmented three times of Hermite interpolation; screening CYGNSS mirror reflection points and constructing a data set in combination with a space matching algorithm; processing CYGNSS L1-level data and extracting observation parameters; establishing a basic GMF (Gaussian Mixture Function) according to incidence angle grouping, and optimizing the weight through nonlinear fitting; the Willoughby model wind speed and the basic GMF are fused to generate an enhanced GMF, the DDM signal-to-noise ratio is dynamically adjusted in combination with minimum variance estimation, and the precision is improved. According to the method, the GNSS-R technology is utilized to overcome the defects of high cost, limited coverage and the like of a traditional method, the sea surface wind speed monitoring reliability under the extreme weather is enhanced through collaborative optimization of a physical model and actually measured data, and the method has the advantages of low cost, wide coverage and interference resistance and is suitable for refined marine meteorological monitoring.

Description

Technical Field

[0001] The present invention belongs to the field of marine remote sensing technology, and in particular relates to a sea surface wind speed inversion method and system based on CYGNSS and the Willoughby model. Background Art

[0002] Ocean wind speed plays a very important role in weather and climate research, ocean currents, sea-air interactions, and marine shipping safety. GNSS-R ocean remote sensing technology uses GNSS navigation and positioning signals that are already ubiquitous in space, so there is no need to design a separate remote sensing equipment transmitter system. It has the advantages of low cost, wide coverage, and all-weather operation. The inversion of ocean surface wind speed is a core mission objective of many Global Navigation Satellite System Reflectometry (GNSS-R) missions. With the continuous development of technology, a variety of inversion methods have emerged, but all methods have different degrees of accuracy reduction when inverting high sea surface wind speeds. Therefore, a new method is needed to improve the inversion accuracy at high wind speeds. To this end, improving its inversion accuracy in the high wind speed range is the focus of current research. Summary of the Invention

[0003] The purpose of the present invention is to provide a sea surface wind speed inversion method and system based on CYGNSS and Willoughby model, so as to overcome the problems of high cost and insufficient accuracy in the prior art.

[0004] The purpose of the present invention is achieved through the following technical solutions:

[0005] A sea surface wind speed inversion method based on CYGNSS and the Willoughby model, the specific steps are as follows:

[0006] Step 1: Based on the synthetic cyclone model Willoughby, cyclone parameters are obtained from the IBTrACS dataset. The sampling interval is encrypted to 1 second by using the time interpolation method of piecewise cubic Hermite interpolation.

[0007] Step 2: Filter the CYGNSS specular reflection points, compile the corresponding IBTrACS / Willoughby model wind speed at the specular reflection points, and calculate the Haversine distance based on the spatial matching algorithm;

[0008] Step 3: Obtain CYGNSS L1 GNSS-R data, ERA5 wind speed and direction, process the CYGNSS data, extract observation variable parameters, perform data quality control and data filtering, build a model dataset, and complete spatiotemporal matching;

[0009] Step 4: Construct the basic GMF, establish the mapping relationship according to the incident angle group, perform weighted processing on the training data, and use nonlinear function fitting;

[0010] Step 5: Convert the wind speed generated by the Willoughby model into NBRCS to form GMF willoughby and basic GMF base The enhanced GMF is generated by mixing weight functions, and finally the minimum variance estimation method is used to dynamically adjust the DDM signal-to-noise ratio changes caused by observation geometry changes to improve the accuracy of sea surface wind speed monitoring data.

[0011] Furthermore, in step 1, the cyclone center position, maximum wind speed, maximum wind speed radius and wind field profile index are obtained by IBTrACS data processing and used as inputs of the Willoughby model;

[0012]

[0013] Among them, V f Defined as the output wind speed, V m Defined as V in , V out or V tr is the model wind speed, V m Depends on the distance from a point to the center of the cyclone separator; V t is the typhoon translation speed, θ t is the direction angle of motion, R max is the maximum wind speed radius, V in is the radius of the cyclone center R max Wind speed inside;

[0014] in, r is the radial distance from the storm center to the ground point;

[0015] V out Describes distances greater than R max +25km surface wind, A is the attenuation coefficient, X1, X2 are characteristic length parameters, R max Internal and R max The wind speed within +25km is ultimately a mixture of the inner and outer zone winds weighted by ω;

[0016] V tr (r)=V in (1-ω)+ωV out .

[0017] Furthermore, the step 1 is a piecewise cubic Hermite interpolation method in the time interpolation method:

[0018] P(t)=h 00 (t)p0+h 10 (t)m0+h 01 (t)p1+h11 (t)m1

[0019] Among them, the specific form of the Hermite basis function is:

[0020] h 00 (t) = 2t 3 -3t 2 +1

[0021] h 10 (t) = t 3 -2t 2 +t

[0022] h 01 (t)=-2t 3 +3t 2

[0023] h 11 (t) = t 3 -t 2

[0024] Calculation of the normalized parameter t:

[0025]

[0026] Among them, x0, x1 are position parameters, p0, p1 are the positions of the interpolation points, and m0, m1 are the tangent values of the interpolation points.

[0027] Furthermore, the step 2 screens CYGNSS specular reflection points to identify CYGNSS specular reflection points within 250 km of the cyclone center and within a time window of ±15 minutes;

[0028] The Haversine distance calculation formula is:

[0029]

[0030] Among them, R = 6371km is the average radius of the earth, is latitude, λ is longitude, Δλ=λ2-λ1, the distance threshold is set to 250km, the time window is ±15 minutes, and the spatial resolution is 0.1°.

[0031] Furthermore, when extracting the new observation variable parameters in step 3, the entire process ensures that each CYGNSS observation point can obtain the corresponding ERA5 wind speed through time matching and spatial interpolation; eliminates low-quality measurements in the data, and selects valid data by calculating the distance correction gain (RCG). The specific formula is:

[0032]

[0033] in, Indicates the receiver antenna gain at the mirror point, in dB. They represent the distance from the transmitter to the mirror and the distance from the mirror point to the receiver respectively; when the spatial resolution of the surface inverted wind speed reaches 25km, the samples outside the 95% confidence interval are eliminated.

[0034] Furthermore, the elimination of low-quality measurement values in the data includes the following conditions:

[0035] The uncertainty of BRCS is less than 1 nm and the attitude status of the star tracker is “OK”;

[0036] The absolute value of the spacecraft roll is between 1 and 30 degrees, the pitch is between 1 and 10, or the yaw angle is between 1 and 5;

[0037] The receiving antenna gain in the direction of the mirror reflection point indicated by quality_flags is greater than 0dBi;

[0038] The range correction gain RCG figure of merit FOM of DDM is greater than 0;

[0039] LES is greater than 0;

[0040] Zenith signal-to-noise ratio is greater than 0dB.

[0041] Furthermore, in step 4, the basic GMF is constructed based on the observation data, which is:

[0042] obs=a0+a1u -1 +a2u -2

[0043] In high wind speed area (u>20m / s):

[0044] obs=b0+b1u+b2u 2

[0045] Through the weight function fusion method: the coefficients a0, a1, a2, b0, b1, b2 are determined by the least squares method, where obs represents the DDM observation and u represents the true wind speed.

[0046] Furthermore, in step 5, GMF willoughby for:

[0047] GMF willoughby =A(θ)[1+B(θ)tan 2 (θ)]V γ(θ) exp(-C(θ)V)

[0048] Where A(θ) is the overall scaling factor, B(θ) is the modulation coefficient related to the incident angle, γ(θ) is the wind speed exponent (usually between 2.0 and 2.5), and C(θ) is the attenuation coefficient; θ is the incident angle; and V is the wind speed generated by the Willoughby model.

[0049] The enhanced GMF first needs to determine the weight function:

[0050] ω=exp(-(uu t ) 2 / 2σ 2 )

[0051] Then perform model fusion according to the following formula:

[0052] GMF enhanced (u,θ)=[1-ω]GMF base +ωGMF willoughby

[0053] Among them, θ is the incident angle, a0, a1, a2, b0, b2, b3 are empirical coefficients, u is the current wind speed, u t is the threshold wind speed, which is set to 30 m / s, and σ is the standard deviation parameter, which controls the width of the transition region.

[0054] Furthermore, the minimum variance estimation method is used to dynamically adjust the DDM signal-to-noise ratio change caused by the observation geometry change, specifically:

[0055] The complete form of the observation equation is:

[0056] y=H(x)+ε

[0057] Where y is the observation vector, which contains NBRCS observations; H is the nonlinear observation operator, which represents the enhanced GMF; χ is the wind speed vector to be estimated; ε is the observation noise; it is assumed to obey a Gaussian distribution with mean 0;

[0058] The minimum variance estimate is:

[0059]

[0060] The estimated error variance is:

[0061]

[0062] in, is the estimated wind speed, n is the number of observation samples, C is the observation error covariance matrix, represents the inverse matrix element of the covariance matrix, is the inverse of the i-th row of the covariance matrix C, y i is the i-th observation, is the variance of the wind speed estimate; the minimum variance estimator is the optimal weight coefficient in the nonlinear observation operator that determines the mapping of different observations to wind speed.

[0063] A computer device / equipment / system comprising a memory, a processor, and a computer program stored in the memory, wherein the processor executes the computer program to implement the steps of a sea surface wind speed inversion method based on CYGNSS and the Willoughby model.

[0064] The beneficial effects of the present invention are:

[0065] Compared with existing technologies, traditional sea surface wind speed measurement methods such as buoy observations and scatterometers have the disadvantages of expensive equipment, limited coverage and greater weather influence. GNSS-R technology has the advantages of low cost, wide coverage and strong anti-interference ability. The retrieval of sea surface wind speed is the core mission objective of GNSS-R missions, such as NASA's Cyclone Global Navigation Satellite System (CYGNSS) mission. However, CYGNSS measurements are affected by multiple sources of uncertainty. As one of the main sources of data sets used in existing research inversion models, the introduction of additional models to provide data on relevant real sea wind speed information provides valuable help for further improving the accuracy and reliability of wind speed inversion. This invention proposes a new wind speed inversion method that cleverly combines the Willoughby physical model and CYGNSS observation data. BRIEF DESCRIPTION OF THE DRAWINGS

[0066] Figure 1 Schematic diagram of the process of the present invention;

[0067] Figure 2 A graph showing the relationship between wind speed and observed values at different incident angles during the implementation of the method of the present invention;

[0068] Figure 3 The relationship between the SNR and distance at different wind speeds according to the method of the present invention is shown below:

[0069] Figure 4 This is an exemplary wind field diagram of the Willoughby model at the CYGNSS location of the present invention. DETAILED DESCRIPTION

[0070] The present invention will be further described below with reference to the accompanying drawings.

[0071] The present invention provides a sea surface wind speed inversion method and system based on CYGNSS and Willoughby model, which can improve the performance of global wind speed estimation of CYGNSS data. From the perspective of the model itself, the inversion model is integrated into the inversion model to make the inversion result more accurate. The method integrates the Willoughby storm model into the CYGNSS wind speed inversion program to improve the high wind speed inversion accuracy. Figure 1 The specific steps are as follows:

[0072] Step 1: Generate high wind speed data that partially matches the CYGNSS data in time and space through the Willoughby storm model, and fuse it with other matching datasets such as ERA5 to form a dataset.

[0073] The Willoughby storm model is used to reproduce the wind speed within a few hundred kilometers of the cyclone center while retaining a simple analytical form. The model uses only the latitude and longitude of the cyclone center at time t, the maximum sustained surface wind (SSW) V max 、Compass angle θ of the cyclone center translation t This equation represents the max ,n,X1,X2,A and ω determine the quantity V max and The equation for a translational storm outputs the wind speed V f Determined as:

[0074]

[0075] Among them, V m Defined as V in , V out or V tr is the model wind speed, V m Depends on the distance from a point to the center of the cyclone. V t is the typhoon translation speed, θ t is the direction angle of motion, R max is the maximum wind speed radius, V in is the radius of the cyclone center R max Wind speed inside.

[0076]

[0077] Where r is the radial distance from the storm center to the ground point. Then, V out Describes distances greater than R max +25km surface wind:

[0078]

[0079] Among them, A is the attenuation coefficient, X1, X2 are characteristic length parameters, R max Internal and R maxThe wind speed within +25km is ultimately a mixture of the inner and outer winds weighted by ω:

[0080] V tr (r)=V in (1-ω)+ωV out

[0081] Storm center location and V required by the Willoughby storm model max The time course information is obtained from the International Climate Management Track Archive (IBTrACS) dataset, which unifies global storm best track data from multiple sources, including Hurricane DATabase (HURDAT), the Automated Tropical Cyclone Forecasting System (ATCF), and the Japan Meteorological Agency (JMA), into a single, self-contained data release.

[0082] Storm position and maximum wind speed information is generated between reporting intervals using linear interpolation. CYGNSS mirror reflection points are identified within 250 km of the cyclone center and within a time window of ±15 minutes. Corresponding IBTrACS / Willoughby model wind speeds are compiled at the mirror reflection points. Interpolation is performed using the piecewise cubic Hermite interpolation method. The time interpolation method is piecewise cubic Hermite interpolation:

[0083] P(t)=h 00 (t)p0+h 10 (t)m0+h 01 (t)p1+h 11 (t)m1

[0084] Among them, the specific form of the Hermite basis function is:

[0085] h 00 (t) = 2t 3 -3t 2 +1

[0086] h 10 (t) = t 3 -2t 2 +t

[0087] h 01 (t)=-2t 3 +3t 2

[0088] h 11 (t) = t 3 -t 2

[0089] Calculation of the normalized parameter t:

[0090]

[0091] Where x0 and x1 are position parameters, p0 and p1 are the positions of the interpolation points, and m0 and m1 are the tangent values at the interpolation points. This interpolation method has the following advantages: maintaining data smoothness, avoiding interpolation overshoot, and maintaining the continuity of physical quantities.

[0092] The spatial matching algorithm includes the Haversine distance calculation:

[0093]

[0094] Among them, R = 6371km is the average radius of the earth, is latitude, λ is longitude, Δλ=λ2-λ1, set the distance threshold to 250km, the time window to ±15 minutes, and the spatial resolution to 0.1°. See the implementation example. Figure 4 Demonstration wind field diagram of the Willoughby model at the CYGNSS location.

[0095] During the generation of the CYGNSS matching dataset, data quality control is also required to ensure high-quality DDM observations: in the new variables of the CYGNSS data file, the entire process ensures that each CYGNSS observation point has a corresponding ERA5 wind speed through time matching and spatial interpolation. High-quality DDM observations are essential for the model to correctly and intuitively understand how to map different ocean wind speeds to their corresponding observations. Therefore, it is essential to eliminate low-quality measurements in the data as follows:

[0096] The uncertainty of BRCS (ddm_brcs_uncert) is less than 1 nanometer. The attitude status of the star tracker is "OK", that is, nst_att_status=0.

[0097] The absolute value of the spacecraft roll is between 1 and 30 degrees, the pitch is between 1 and 10, or the yaw angle is between 1 and 5, and the receive antenna gain in the direction of the specular reflection point (sp_rx_gain) indicated by quality_flags is greater than 0dBi.

[0098] The range correction gain (RCG) figure of merit (FOM) of the DDM (prn_fig_of_merit) is greater than 0.

[0099] LES(ddm_les) is greater than 0.

[0100] The Zenith (direct) signal-to-noise ratio (direct_signal_snr) is greater than 0 dB.

[0101] And the valid data is screened by calculating the range correction gain (RCG). The specific formula is:

[0102]

[0103] in, Indicates the receiver antenna gain at the mirror point, in dB. where represents the distance from the transmitter to the mirror and the distance from the mirror to the receiver, respectively. The spatial resolution of the DDM is affected by the geometry of the bistatic mine. When the spatial resolution of the surface wind speed inversion reaches 25 km, the LES can be further time-averaged to improve the signal-to-noise ratio of the observation.

[0104] In addition, because CYGNSS measurements or ERA5 estimates may also be contaminated by unexpected errors that usually do not affect the main population, samples outside the 95% confidence interval are discarded.

[0105] Step 2: Establish a complete GMF inversion model to invert the wind speeds in the high, medium and low data segments.

[0106] The direct observation quantity of spaceborne GNSS-R technology is DDM. The actual wind speed inversion algorithm is to obtain the empirical mapping relationship by extracting the wind speed-sensitive observation quantity in DDM and matching it with the corresponding wind speed regression of other wind speed observation systems. Therefore, it is necessary to understand the properties of different observation quantities before building the model. Figure 2 It shows the relationship between wind speed and observed value at different incident angles. Figure 3 The relationship between SNR and distance at different wind speeds is shown. After understanding the relationship between the observations, the core innovation lies in the construction method of the enhanced GMF. Before establishing the GMF, the mirror point incident angle is used as the independent variable, and the sampling of consecutive epochs before and after the observation is averaged in different incidence angle intervals. The time averaging rule is as follows:

[0107] Table 1 Time averaging rules

[0108] Time average sample size Incident angle range (°) 5 0°<θ<17° 4 17°<θ<31° 3 31°<θ<41° 2 41°<θ<48° 1 ＞48°

[0109] The GMF established by the inversion model is actually a function cluster, that is, a set of mapping relationships between NBRCS, LES and wind speed is established at different incident angles, and is represented by scatter points of wind speed and corresponding observations with a certain step size. When establishing the GMF, the incident angle range of the training data is 0° to 68°. In the incident angle dimension, the GMF is established at different incident angles starting from 0.5° with a step size of 1°. At a certain incident angle, the weighted average observations at different wind speeds are calculated with a step size of 0.1m / s starting from a wind speed of 0.05m / s, and a discrete empirical mapping relationship is established. At a certain incident angle, all data within the two step size intervals (step size is 1°) on the left and right are taken as the model training data at the incident angle. At a certain wind speed, the data within the two step size intervals on the left and right are also taken to calculate the weighted DDM observations. It should be noted that the wind speed interval step size needs to be re-determined according to the wind speed probability density. The wind speed step size selection strategy is shown in Table 2:

[0110] Table 2 Wind speed step selection strategy

[0111] Step length (m / s) Wind speed range (m / s) 0.4 U10<5 0.3 5<U10<8 0.2 8<U10<12 0.4 12<U10<16 0.6 16<U10<20 0.8 20<U10<24 1.0 24<U10

[0112] To ensure the accuracy of the empirical model, we first calculated the maximum probability density interval of the training data wind speed distribution at each incident angle. To eliminate model jitter caused by insufficient training data and noise, we calculated the corresponding weighted DDM observations starting from the wind speed corresponding to the maximum wind speed probability density interval. We then calculated the observations for wind speeds above and below this maximum probability density interval in steps of 0.1 m / s, forcing the discrete GMF to be a monotonic function. The empirical GMF was then smoothed using a piecewise nonlinear function.

[0113] The basic GMF is established based on observation data and adopts the following form in low wind speed areas (u<20m / s):

[0114] obs=a0+a1u -1 +a2u -2

[0115] Use in high wind speed areas (u>20m / s):

[0116] obs=b0+b1u+b2u 2

[0117] Through weight function fusion: the coefficients a0, a1, a2, b0, b1, and b2 are determined by the least squares method. Where obs represents the DDM observation and u represents the true wind speed.

[0118] After generating the wind speed through the previous Willoughby model, it is converted into NBRCS. The specific formula is:

[0119] GMFwilloughby =A(θ)[1+B(θ)tan 2 (θ)]V γ(θ) exp(-C(θ)V)

[0120] Where A(θ) is the overall scale factor, B(θ) is the modulation coefficient related to the incident angle, γ(θ) is the wind speed index (usually between 2.0 and 2.5), and C(θ) is the attenuation coefficient. θ is the incident angle, and V is the wind speed generated by the Willoughby model. The calculation formula for the specific parameters is:

[0121] A(θ)=a0+a1θ+a2θ 2

[0122] B(θ)=b0+b1exp(-b2θ)

[0123] The enhanced GMF first needs to determine the weight function:

[0124] ω=exp(-(uu t ) 2 / 2σ 2 )

[0125] Then perform model fusion according to the following formula:

[0126] GMF enhanced (u,θ)=[1-ω]GMF base +ωGMF willoughby

[0127] Among them, θ is the incident angle, a0, a1, a2, b0, b2, b3 are empirical coefficients, u is the current wind speed, u t is the threshold wind speed (taken as 30 m / s), and σ is the standard deviation parameter, which controls the width of the transition region.

[0128] The main challenge is to determine the appropriate segmentation points to ensure optimal model fit while smoothing at the transition points. First, smoothing is performed at different discrete points to find the discrete point with the smallest fitting residual. Then, a search interval is set around the discrete point, and the point with the smallest fitting residual is searched again with a smaller step size as the final transition point.

[0129] The minimum variance estimation method dynamically adjusts the DDM signal-to-noise ratio changes caused by changes in observation geometry. That is, in different RCG intervals, the wind speeds mapped by different NBRCS and LES are weighted to obtain the optimal wind speed estimate. The RCG can first select a smaller interval and iteratively adjust it based on the estimation results.

[0130] The complete form of the observation equation is:

[0131] y=H(x)+ε

[0132] Where y is the observation vector, which contains the NBRCS observations; H is the nonlinear observation operator, which represents the enhanced GMF; χ is the wind speed vector to be estimated, and ε is the observation noise, which is assumed to obey a Gaussian distribution with mean 0.

[0133] The complete analytical solution for the minimum variance estimate is:

[0134]

[0135] The estimated error variance is:

[0136]

[0137] in, is the estimated wind speed, n is the number of observation samples, C is the observation error covariance matrix, represents the inverse matrix element of the covariance matrix, is the inverse of the i-th row of the covariance matrix C, y i is the i-th observation, is the variance of the wind speed estimate. The minimum variance estimator is the optimal weight coefficient in the nonlinear observation operator that determines the mapping of different observations to wind speed.

[0138] The inversion model was evaluated during a test period from May 1 to July 2017. Because the focus of the study was on the effects of introducing physical constraints, the model performance was quantified using the following formula:

[0139]

[0140] in, is the predicted wind speed obtained from the i-th inversion, v i is the i-th reference wind speed from ERA5. The overall RMSE and bias are determined based on a total of n samples across all models. The root mean square error applies to all inversion architectures as well as the MVE traditional model.

[0141] The present invention also provides a sea surface wind speed inversion system based on CYGNSS and the Willoughby model, and the system is used for the above-mentioned sea surface wind speed inversion method based on CYGNSS and the Willoughby model.

[0142] The foregoing description is merely a preferred embodiment of the present invention and is not intended to limit the present invention. Those skilled in the art will readily appreciate that various modifications and variations of the present invention are possible. Any modifications, equivalent substitutions, or improvements made within the spirit and principles of the present invention shall be included within the scope of protection of the present invention.

Claims

1. A sea surface wind speed inversion method based on CYGNSS and the Willoughby model, characterized by: The specific steps are as follows: Step 1: Based on the synthetic cyclone model Willoughby, cyclone parameters are obtained from the IBTrACS dataset. The sampling interval is encrypted to 1 second by using the time interpolation method of piecewise cubic Hermite interpolation. Step 2: Filter the CYGNSS specular reflection points, compile the corresponding IBTrACS / Willoughby model wind speed at the specular reflection points, and calculate the Haversine distance based on the spatial matching algorithm; Step 3: Obtain CYGNSS L1 GNSS-R data, ERA5 wind speed and direction, process the CYGNSS data, extract observation variable parameters, perform data quality control and data filtering, build a model dataset, and complete spatiotemporal matching; Step 4: Construct the basic GMF, establish the mapping relationship according to the incident angle group, perform weighted processing on the training data, and use nonlinear function fitting; Step 5: Convert the wind speed generated by the Willoughby model into NBRCS to form GMF willoughby and basic GMF base The enhanced GMF is generated by mixing weight functions, and finally the minimum variance estimation method is used to dynamically adjust the DDM signal-to-noise ratio changes caused by observation geometry changes to improve the accuracy of sea surface wind speed monitoring data.

2. The sea surface wind speed inversion method based on CYGNSS and the Willoughby model according to claim 1, characterized in that: The step 1 obtains the cyclone center position, maximum wind speed, maximum wind speed radius and wind field profile index by IBTrACS data processing, and uses them as inputs of the Willoughby model; Among them, V f Defined as the output wind speed, V m Defined as V in , V out or V tr is the model wind speed, V m Depends on the distance from a point to the center of the cyclone separator; V t is the typhoon translation speed, θ t is the direction angle of motion, R max is the maximum wind speed radius, V in is the radius of the cyclone center R max Wind speed inside; in, r is the radial distance from the storm center to the ground point; V out Describes distances greater than R max +25km surface wind, A is the attenuation coefficient, X1, X2 are characteristic length parameters, R max Internal and R max The wind speed within +25km is ultimately a mixture of the inner and outer zone winds weighted by ω; V tr (r)=V in (1-ω)+ωV out 。 3. The sea surface wind speed inversion method based on CYGNSS and the Willoughby model according to claim 2, characterized in that: The step 1 is the time interpolation method piecewise cubic Hermite interpolation: P(t)=h 00 (t)p0+h 10 (t)m0+h 01 (t)p1+h 11 (t)m1 Among them, the specific form of the Hermite basis function is: h 00 (t)=2t 3 -3t 2 +1 h 10 (t)=t 3 -2t 2 +t h 01 (t)=-2t 3 +3t 2 h 11 (t)=t 3 -t 2 Calculation of the normalized parameter t: Among them, x0, x1 are position parameters, p0, p1 are the positions of the interpolation points, and m0, m1 are the tangent values of the interpolation points.

4. The sea surface wind speed inversion method based on CYGNSS and the Willoughby model according to claim 1, characterized in that: Step 2 screens CYGNSS specular reflection points and identifies CYGNSS specular reflection points within 250 kilometers of the cyclone center and within a time window of ±15 minutes; The Haversine distance calculation formula is: Among them, R = 6371km is the average radius of the earth, is latitude, λ is longitude, Δλ=λ2-λ1, the distance threshold is set to 250km, the time window is ±15 minutes, and the spatial resolution is 0.1°.

5. The sea surface wind speed inversion method based on CYGNSS and the Willoughby model according to claim 1, characterized in that: When extracting the new observation variable parameters in step 3, the entire process ensures that each CYGNSS observation point can obtain the corresponding ERA5 wind speed through time matching and spatial interpolation; eliminates low-quality measurements in the data, and selects valid data by calculating the range correction gain (RCG). The specific formula is: in, Indicates the receiver antenna gain at the mirror point, in dB. They represent the distance from the transmitter to the mirror and the distance from the mirror point to the receiver respectively; when the spatial resolution of the surface inverted wind speed reaches 25km, the samples outside the 95% confidence interval are eliminated.

6. The sea surface wind speed inversion method based on CYGNSS and the Willoughby model according to claim 5, characterized in that: The low-quality measurement values in the data are eliminated, and the judgment conditions for eliminating the data include: The uncertainty of BRCS is less than 1 nanometer. The attitude status of the star tracker is "OK". The absolute value of the spacecraft roll is between 1 and 30 degrees, the pitch is between 1 and 10, or the yaw angle is between 1 and 5; The receiving antenna gain in the direction of the mirror reflection point indicated by quality_flags is greater than 0dBi; The range correction gain RCG figure of merit FOM of DDM is greater than 0; LES is greater than 0; Zenith signal-to-noise ratio is greater than 0dB.

7. The sea surface wind speed inversion method based on CYGNSS and the Willoughby model according to claim 1, characterized in that: In step 4, the basic GMF is constructed based on the observation data. In the low wind speed area (u<20m / s), it is: obs=a0+a1u -1 +a2u -2 In high wind speed area (u>20m / s): obs=b0+b1u+b2u 2 Through the weight function fusion method: the coefficients a0, a1, a2, b0, b1, b2 are determined by the least squares method, where obs represents the DDM observation and u represents the true wind speed.

8. The sea surface wind speed inversion method based on CYGNSS and the Willoughby model according to claim 1, characterized in that: GMF in step 5 willoughby for: GMF willoughby =A(θ)[1+B(θ)tan 2 (i)]V γ(θ) exp(-C(θ)V) Where A(θ) is the overall scaling factor, B(θ) is the modulation coefficient related to the incident angle, γ(θ) is the wind speed exponent (usually between 2.0 and 2.5), and C(θ) is the attenuation coefficient; θ is the incident angle; and V is the wind speed generated by the Willoughby model. The enhanced GMF first needs to determine the weight function: ω=exp(-(uu t ) 2 / 2σ 2 ) Then perform model fusion according to the following formula: GMF enhanced (u,θ)=[1-ω]GMF base +ωGMF willoughby Among them, θ is the incident angle, a0, a1, a2, b0, b2, b3 are empirical coefficients, u is the current wind speed, u t is the threshold wind speed, which is set to 30 m / s, and σ is the standard deviation parameter, which controls the width of the transition region.

9. The sea surface wind speed inversion method based on CYGNSS and the Willoughby model according to claim 8, characterized in that: The method of using the minimum variance estimation method to dynamically adjust the DDM signal-to-noise ratio change caused by the observation geometry change is specifically as follows: The complete form of the observation equation is: y=H(x)+ε Where y is the observation vector, which contains NBRCS observations; H is the nonlinear observation operator, which represents the enhanced GMF; χ is the wind speed vector to be estimated; ε is the observation noise; it is assumed to obey a Gaussian distribution with mean 0; The minimum variance estimate is: The estimated error variance is: in, is the estimated wind speed, n is the number of observation samples, C is the observation error covariance matrix, represents the inverse matrix element of the covariance matrix, is the inverse of the i-th row of the covariance matrix C, y i is the i-th observation, is the variance of the wind speed estimate; the minimum variance estimator is the optimal weight coefficient in the nonlinear observation operator that determines the mapping of different observations to wind speed.

10. A computer device / apparatus / system comprising a memory, a processor, and a computer program stored in the memory, characterized in that: The processor executes the computer program to implement the steps of the method according to any one of claims 1 to 8.

Images