GNSS-R sea surface wind speed inversion method

Abstract

The invention relates to the technical field of global navigation satellite system reflection, in particular to a GNSS-R sea surface wind speed inversion method, which comprises the following steps of: firstly, comprehensively considering multiple constraint conditions such as inversion error distribution characteristics, interval quantity, sample distribution uniformity and interval minimum sample requirements by utilizing the global optimal search characteristic of an SA algorithm; a multi-objective optimization function is constructed to preliminarily divide a wind speed interval, and then a GD algorithm is adopted to accurately adjust the boundary of the preliminarily divided interval. And constructing a special XGBoost prediction model for each wind speed interval, and finally, fusing wind speed prediction results of a plurality of sub-interval models through a stacked integrated learning architecture to obtain a wind speed inversion value. According to the method provided by the invention, the problem that the inversion error of a traditional single model in different wind speed ranges is extremely unstable can be effectively relieved, the inversion precision of the whole wind speed range is improved, and meanwhile, the problem that a real-time observation sample cannot correspond to a sub-model under an interval model in practical application is solved; and the application conversion of the GNSS-R sea surface wind speed inversion technology is promoted.

Description

Technical Field

[0001] This invention relates to the field of global navigation satellite system reflection technology, specifically to a GNSS-R sea surface wind speed inversion method. Background Technology

[0002] Global Navigation Satellite System Reflection (GNSS-R) technology, as an emerging marine remote sensing technology, retrieves parameters such as sea surface wind speed by receiving satellite navigation signals reflected from the sea surface. This technology boasts significant advantages such as all-weather, all-time availability, high coverage, and low cost, making it valuable for applications in marine environmental monitoring, climate change research, and extreme weather early warning.

[0003] However, existing GNSS-R sea surface wind speed inversion methods mainly suffer from the following technical challenges: Traditional single inversion models use a method of directly fitting all observations with wind speeds, resulting in significant error differences in inversion results for samples with different wind speed ranges, especially with a significant decrease in inversion accuracy for high wind speed samples; Existing wind speed interval modeling methods mostly rely on fixed empirical threshold divisions, lacking a comprehensive consideration of the distribution characteristics of sample data and the characteristics of inversion errors, and thus failing to achieve optimal interval division; Traditional interval division methods only establish independent sub-models, ignoring the correspondence and matching problem between observations and different interval models in practical applications, and the sample inversion results at interval boundaries will exhibit discontinuities during the inversion process. Summary of the Invention

[0004] The purpose of this invention is to provide a GNSS-R sea surface wind speed inversion method, which aims to solve the technical problems existing in the current method, such as unreasonable wind speed interval division, significant differences in inversion error between different wind speed ranges, and discontinuous boundary of inversion results, so as to achieve high-precision and robust inversion of sea surface wind speed.

[0005] To achieve the above objectives, this invention provides a GNSS-R sea surface wind speed inversion method, comprising the following steps:

[0006] Step 1: Acquire spaceborne GNSS-R observation data and ERA5 reanalysis data, and perform data processing on the acquired CYGNSS satellite observations;

[0007] Step 2: Perform quality control on the solved CYGNSS data, interpolate the ERA5 reanalysis data, and perform spatiotemporal matching between the quality-controlled CYGNSS data and the interpolated ERA5 data to construct the dataset;

[0008] Step 3: Use the improved SA algorithm to perform preliminary division of wind speed ranges for multi-objective optimization;

[0009] Step 4: Use the GD algorithm to further refine the boundaries of the initially defined wind speed intervals;

[0010] Step 5: Construct a dedicated XGB prediction model for each of the defined wind speed intervals;

[0011] Step 6: Use a stacked ensemble learning architecture to fuse the wind speed prediction results of the models in each interval and output the final sea surface wind speed inversion value.

[0012] Optionally, the CYGNSS satellite observations include the normalized bistatic radar cross section and the slope of the leading edge of the integral delay waveform, which need to be calculated from the satellite observation data;

[0013] Step 1 also considers six other feature variables, namely the longitude of the specular reflection point, the latitude of the specular reflection point, the incident angle of the specular reflection point, the gain of the receiving antenna at the specular reflection point, the distance from the transmitting point to the specular reflection point, and the distance from the receiving point to the specular reflection point, in order to improve the prediction accuracy of the model. The above six feature variables can be directly obtained from satellite observation data.

[0014] During the execution of step 1, the easterly and northerly wind components of the sea surface wind speed at 10 meters are extracted from the ERA5 reanalysis dataset, and the vector summation is calculated to obtain the wind speed value at the corresponding location.

[0015] Optionally, in step 2, the following quality control measures were taken for the two calculated feature variables and the six directly obtained feature variables: observation samples with missing values ​​in the original data were removed; observation samples with negative slope values ​​of the leading edge of the integral delay waveform were removed; observation samples located in land areas were removed; observation samples with a gain value of less than 0 for the receiving antenna at the specular reflection point were removed; observation samples with an incident angle greater than 60° at the specular reflection point were removed; the receiver gain of each sample was calculated, and observation samples with a receiver gain less than 10 were removed.

[0016] The interpolation process for ERA5 reanalysis data specifically involves performing temporal and spatial interpolation on the ERA5 data. Bilinear interpolation is used to adjust the temporal and spatial resolutions of the ERA5 data to be consistent with those of the CYGNSS data. Then, the CYGNSS data with consistent temporal and spatial information is spatiotemporally matched with the ERA5 data to form data pairs. All matched data pairs constitute the dataset.

[0017] Optionally, in step 3, the SA algorithm sets the initial temperature, cooling rate, and initial weights. For each iteration, it generates a perturbation and calculates the objective function value corresponding to the new boundary. It decides whether to accept the newly divided interval based on the acceptance probability, while updating the weight parameters and temperature, and returns the optimal interval boundary that minimizes the objective function value.

[0018] Construct a multi-objective optimization function that comprehensively considers the characteristics of inversion error distribution, number of intervals, uniformity of sample distribution, and minimum sample requirement of intervals. It includes the sum of inversion errors of each interval, interval number penalty term, sample distribution variance penalty term, and minimum sample penalty term.

[0019] A multi-objective optimization strategy based on a fixed normalization range is adopted to ensure that each factor can play a normal regulatory role in the optimization process.

[0020] Optionally, in step 4, the GD algorithm uses numerical differentiation to calculate the gradient. Specifically, the gradient is estimated by perturbing the boundary values ​​with small values ​​and observing the changes in the objective function. A fixed learning rate is used to control the step size of each update, and the boundary values ​​are updated according to the corresponding rules.

[0021] Optionally, in step 5, the objective function of the XGB model corresponding to each interval includes the loss between the wind speed prediction value and the actual value, and a complexity regularization term. The regularization term is used to control the model complexity and prevent overfitting.

[0022] Optionally, in step 6, the prediction performance is improved by using a two-level stacked integrated model structure. The first-level model is the base model, which contains N heterogeneous models, each trained using samples from different wind speed ranges. The second-level model is the meta-model, which also uses the XGB model to integrate the output of the base model. Its input features are a combination of the prediction results of the base model and the original features.

[0023] For training the meta-model, a meta-feature matrix is ​​generated through K-fold cross-validation. The training set is divided into K subsets, and one subset is retained as the validation set each time. The remaining subsets are used to train the base model and generate meta-features. The meta-features of all folds are integrated to train the meta-model to avoid data leakage.

[0024] The XGB model is used as the meta-model, and the model complexity is controlled by the regularization coefficient. Finally, the wind speed inversion value is output.

[0025] This invention provides a GNSS-R sea surface wind speed inversion method. First, the SA algorithm is used to efficiently search for global wind speed interval division schemes. Taking into account multiple constraints such as inversion error distribution characteristics, number of intervals, sample distribution uniformity, and minimum sample requirement for each interval, a multi-objective optimization function is constructed for initial wind speed interval division. Then, the GD algorithm is used to further refine the boundaries of the initially divided intervals. For each divided wind speed interval, a dedicated XGBoost (XGB) prediction model is constructed. Finally, the wind speed prediction results from multiple interval models are fused using a stacked ensemble learning architecture to obtain the wind speed inversion value. The method proposed in this invention effectively addresses the problem of extremely unstable inversion errors in traditional single models across different wind speed ranges, improving the inversion accuracy across the entire wind speed range. It also solves the problem in practical applications where real-time observation samples under interval models cannot correspond to sub-models, thus promoting the application and transformation of GNSS-R sea surface wind speed inversion technology. Attached Figure Description

[0026] To more clearly illustrate the technical solutions in the embodiments of the present invention or the prior art, the drawings used in the description of the embodiments or the prior art will be briefly introduced below. Obviously, the drawings described below are only some embodiments of the present invention. For those skilled in the art, other drawings can be obtained based on these drawings without creative effort.

[0027] Figure 1 This is a schematic diagram illustrating the process principle of a GNSS-R sea surface wind speed inversion method according to the present invention.

[0028] Figure 2 This is a graph showing the change of the objective function during the optimization process of the SA algorithm in a specific embodiment of the present invention.

[0029] Figure 3 This is a graph showing the change of the objective function during the fine-tuning process of the GD algorithm in a specific embodiment of the present invention.

[0030] Figure 4 This is a comparison chart of scatter density maps of different prediction methods in a specific embodiment of the present invention.

[0031] Figure 5 This is a comparison chart of the RMSE performance of different methods in various wind speed ranges in specific embodiments of the present invention. Detailed Implementation

[0032] Embodiments of the present invention are described in detail below, examples of which are illustrated in the accompanying drawings, wherein the same or similar reference numerals denote the same or similar elements or elements having the same or similar functions throughout. The embodiments described below with reference to the accompanying drawings are exemplary and intended to explain the present invention, and should not be construed as limiting the present invention.

[0033] The following explains some of the abbreviations of English terms used in this invention:

[0034] GNSS-R: Global Navigation Satellite Systems - Reflectometry.

[0035] CYGNSS: Cyclone Global Navigation Satellite System;

[0036] ECMWF: European Centre for Medium-Range Weather Forecasts;

[0037] ERA5: The fifth-generation global climate and atmosphere reanalysis dataset;

[0038] SA: Simulated Annealing Algorithm;

[0039] GD: Gradient Descent Algorithm;

[0040] NBRCS: Normalized Bistatic Radar Cross-Section;

[0041] LES: Slope of the leading edge of the integral delay waveform;

[0042] RCG: Receiver Channel Gain;

[0043] XGB: Extreme gradient hints.

[0044] Please see Figure 1 This invention provides a GNSS-R method for inverting sea surface wind speed, comprising the following steps:

[0045] S1: Acquire spaceborne GNSS-R observation data and ERA5 reanalysis data, and perform data processing on the acquired CYGNSS satellite observations;

[0046] S2: Perform quality control on the solved CYGNSS data, interpolate the ERA5 reanalysis data, and perform spatiotemporal matching between the quality-controlled CYGNSS data and the interpolated ERA5 data to construct a dataset;

[0047] S3: Preliminary division of wind speed ranges for multi-objective optimization using an improved SA algorithm;

[0048] S4: The GD algorithm is used to further refine the boundaries of the initially defined wind speed intervals;

[0049] S5: Construct a dedicated XGB prediction model for each defined wind speed range;

[0050] S6: Employs a stacked ensemble learning architecture to fuse the wind speed prediction results from models in different regions, and outputs the final sea surface wind speed inversion value.

[0051] The following provides further explanation in conjunction with the specific implementation steps:

[0052] In step S1, the extracted observations include NBRCS and LES, which are obtained by solving the satellite observation data. The NBRCS calculation formula is as follows:

[0053]

[0054] in, For the calibrated signal power, , For atmospheric loss correction. For instrument wear and tear correction, This indicates the satellite's transmitted signal power at the mirror point; Indicates the wavelength of the L-band; This represents the transmitter antenna gain at the reflection point. This represents the receiver antenna gain at the reflection point. For distance loss from the transmitter to the sea surface and from the sea surface to the receiver, , These are the code correlation function and the Doppler frequency shift function, respectively.

[0055] The LES calculation formula is as follows:

[0056]

[0057] in, This represents the time delay value of the integral delay waveform. The amplitude value of the integral delay waveform. The number of sampling points used to calculate the slope of the front edge.

[0058] In addition, six other characteristic variables are also considered: the longitude of the specular reflection point, the latitude of the specular reflection point, the angle of incidence of the specular reflection point, the gain of the receiving antenna at the specular reflection point, the distance from the transmitting point to the specular reflection point, and the distance from the receiving point to the specular reflection point, to improve the prediction accuracy of the model. These six characteristic variables can be directly obtained from satellite observation data. The easterly and northerly wind components of the sea surface wind speed at 10 meters are extracted from the ERA5 reanalysis dataset, and the vector summation is used to obtain the wind speed value at the corresponding location.

[0059] In step S2, the following quality control measures were taken for the two calculated feature variables and the six directly obtained feature variables: observation samples with missing values ​​in the original data were removed; observation samples with negative slope values ​​of the leading edge of the integral delay waveform were removed; observation samples located in land areas were removed; observation samples with a gain value of less than 0 for the receiving antenna at the specular reflection point were removed; observation samples with an incident angle greater than 60° at the specular reflection point were removed; the receiver gain of each sample was calculated, and observation samples with a receiver gain less than 10 were removed. The formula for calculating the receiver gain (RCG) is as follows:

[0060]

[0061] in, This represents the distance from the transmitter to the point of reflection on the mirror. This represents the distance from the receiver to the point of reflection on the mirror. This represents the receiver antenna gain at the point of specular reflection.

[0062] The interpolation process for ERA5 reanalysis data specifically involves performing temporal and spatial interpolation on the ERA5 data. Bilinear interpolation is used to adjust the temporal and spatial resolutions of the ERA5 data to be consistent with those of the CYGNSS data. Then, the CYGNSS data with consistent temporal and spatial information is spatiotemporally matched with the ERA5 data to form data pairs. All matched data pairs constitute the dataset.

[0063] In step S3, the specific process of the SA algorithm is as follows:

[0064] (1) Set the initial temperature Cooling rate Initial weights , ;

[0065] (2) Determine the normalization range of each term of the objective function by sampling. ;

[0066] (3) For each iteration With probability Perform interval number adjustment to generate disturbance:

[0067]

[0068] in, Indicates the new interval boundary. Indicates the current interval boundary. Indicates the variance of the disturbance. This indicates the current iteration number of the simulated annealing algorithm.

[0069] Calculate the objective function value corresponding to the new boundary, i.e. Calculate the reception probability according to the formula:

[0070]

[0071] in, Let be the temperature of the nth iteration, when When generating a random number a in the range of 0-1, if a ≥ If the condition is met, the newly defined interval will be rejected, and the currently defined interval will be maintained; otherwise, the newly defined interval will be accepted.

[0072] (4) Update the weight parameters and temperature, as shown in the following formula:

[0073]

[0074] in, The temperature at the (n+1)th iteration;

[0075] (5) Return the optimal interval boundary that minimizes the objective function value. .

[0076] During the execution of step (3), a multi-objective optimization function is constructed that comprehensively considers the characteristics of the prediction error distribution, the number of intervals, the uniformity of sample distribution, and the minimum sample requirement of the interval:

[0077]

[0078] in, , , and These are the sum of prediction errors for each interval after normalization, the interval quantity penalty term, the sample distribution variance penalty term, and the minimum sample penalty term, respectively. and These are the corresponding weighting coefficients.

[0079] The total error term represents the sum of prediction errors for each interval, reflecting the prediction accuracy of each interval. When the intervals are appropriately divided, the errors in each interval decrease, and the overall error decreases. Its calculation formula is:

[0080]

[0081] For the first The mean square error of each interval For the first The error weight for each interval is based on the interval average wind speed. Dynamic calculation:

[0082]

[0083] This is the base value for calculating the error weight, i.e., the weight value at the lowest wind speed. This is the error weighting growth coefficient, which controls the magnitude of the error weighting increase with wind speed. and These represent the maximum and minimum wind speeds within the interval, respectively.

[0084] For the interval quantity penalty term, an improved penalty function is used, which is characterized by producing slight changes even within the target interval. The formula is:

[0085]

[0086] in, Indicates the number of intervals. and These represent the target lower and upper limits for the number of intervals, respectively, and the penalty coefficients of 2.0, 0.5, and 0.1 are optimization parameters determined through experiments.

[0087] The variance of the sample distribution quantifies the degree of unevenness in the distribution of samples across intervals, and its calculation formula is as follows:

[0088]

[0089] in, and These represent the number of samples in the interval and the average number of samples in the interval, respectively.

[0090] As a minimum sample penalty term, this invention addresses the problems of model overfitting and biased interval data distribution estimation caused by insufficient samples by dynamically constraining intervals with insufficient sample size. The penalty term is designed based on wind speed characteristics, and its calculation formula is as follows:

[0091]

[0092] in, For the penalty function, Indicates the first The dynamic minimum sample requirement for each interval.

[0093] To avoid excessively large differences in the variation amplitudes of different objective functions during multi-objective optimization, which could mask objectives with smaller variations and thus affect the effectiveness of the optimization process, this invention proposes a multi-objective optimization strategy based on a fixed normalization range. Before optimization begins, this strategy determines the fixed normalization range of each objective function through initial sampling, and its calculation formula is as follows:

[0094]

[0095] in, The range of the total error term. The range of the interval penalty term. For the range of the sample variance term, This represents the range of the minimum sample penalty term. The fixed range is determined by collecting statistical information on each factor through multiple random perturbations of the initial interval.

[0096]

[0097] in, and The first The minimum and maximum values ​​of each factor during the sampling process.

[0098] After normalization, the components of the objective function are mapped to the same numerical range, ensuring that each factor can play a normal regulatory role in the optimization process.

[0099] In step S4, the gradient descent algorithm uses a numerical approximate differentiation method to calculate the gradient:

[0100]

[0101] in, Take 0.1, Indicates the first For each interval boundary value, the gradient is estimated by perturbing the boundary values ​​with small perturbations and observing changes in the objective function, and the boundary values ​​are updated according to the following rules:

[0102]

[0103] in, To represent a fixed learning rate, the step size for each update is controlled; and These represent the current interval boundary value and the interval boundary value after gradient descent adjustment, respectively.

[0104] In step S5, the objective function of the XGB model for each interval is:

[0105]

[0106] in, The loss between the predicted and actual wind speed values ​​(e.g., mean square error) For the first The complexity regularization term for each tree is used to control model complexity and prevent overfitting; The number of samples; The total number of trees; For the first The complexity of the tree is a regularization term.

[0107] In step S6, a two-level stacked ensemble model structure is used to improve prediction performance. The first-level model (base model) contains N heterogeneous models, each trained in different wind speed ranges. The second-level model (meta-model) integrates the output of the base models using the XGB model. The input features are as follows:

[0108]

[0109] in, For the nth base model, Original features These are variables extracted from the onboard dataset. For the nth interval, the objective function of the base model is as follows:

[0110]

[0111] The first item represents the error between the actual value and the predicted value. The number of leaf nodes. Leaf weight, and These are the hyperparameters that control the complexity of the tree and the weight regularization, respectively.

[0112] For training the meta-model, the meta-feature matrix is ​​generated through K-fold cross-validation as shown in the following equation:

[0113]

[0114] The meta-feature matrix is ​​generated by K-fold cross-validation. The training set is divided into K subsets. One subset is kept as the validation set each time, and the remaining subsets are used to train the base model and generate meta-features. The meta-features of all folds are integrated to train the meta-model and avoid data leakage.

[0115] Using XGB as the meta-model Its optimization objective is as follows:

[0116]

[0117] in, The regularization coefficient is . Represents the weight parameter vector of the meta-model; Indicates the first Meta-features of each sample. The final wind speed prediction is as follows:

[0118]

[0119] Furthermore, to verify the effectiveness of the proposed method, specific embodiments are proposed, and multiple experiments are conducted to verify the model inversion performance of the proposed method. The experiments used CYGNSS L1 data and ERA5 reanalysis data from May and June 2024 to construct datasets, with 80% of the data used as the training set and 20% as the test set. Several comparative experimental schemes were designed, including: a global single model (without interval division), a traditional experience-based wind speed numerical equal-four interval modeling method, a wind speed numerical equal-three interval modeling method, a sample size equal-three interval modeling method, and the SA-GD-XGB-SE modeling method proposed in this invention. The performance of each method was evaluated and compared using three metrics: mean squared error (MSE), root mean square error (RMSE), and the square of the correlation coefficient (R²).

[0120] Figure 2 The results show the trend of the objective function value changing with the number of iterations during the SA algorithm optimization process. It can be observed that as the iteration proceeds, the interval boundary is continuously updated and adjusted, and the objective function value gradually decreases and tends to stabilize, verifying the effectiveness of the SA algorithm in adaptive wind speed interval division.

[0121] Figure 3 The results demonstrate the changes in the objective function during the fine-tuning process of the GD algorithm. The results show that the adjustment range of parameters such as boundary values ​​and objective function values ​​gradually decreases with the increase of the number of iterations and eventually tends to stabilize. The whole process gradually converges, achieving precise adjustment of the interval boundary.

[0122] Figure 4 (a)-(e) show a comparison of scatter plots for different prediction methods, where Figure 4 (a) Scatter density plot for the global single model method Figure 4 (b) Scatter density plot for the method of modeling wind speed values ​​in equal quarter intervals. Figure 4 (c) Scatter density plot of the method for modeling wind speed values ​​in three equal intervals. Figure 4 (d) Scatter density plot of the three-interval modeling method with equal sample size and Figure 4 (e) Scatter density plot of the SA-GD-XGB-SE method. As can be seen from the figure, the method of the present invention effectively improves the overall inversion accuracy of the model, the data points are more densely distributed, and the performance is more stable, especially in areas with wind speeds greater than 15 m / s. The correlation between the predicted values ​​and the actual values ​​is significantly improved.

[0123] Figure 5 (a)-(d) show the comparison of RMSE performance of different methods in various wind speed ranges, where Figure 5 (a) is a comparison chart of the effects of the equal four-interval modeling method for wind speed values ​​and the single-model method. Figure 5 (b) is a comparison chart of the effects of the three-interval modeling method and the single-model method for wind speed numerical equalization. Figure 5 (c) is a comparison chart of the performance of the three-interval modeling method with equal sample sizes and the single-model method. Figure 5 (d) is a comparison of the performance of the SA-GD-XGB-SE method and the single-model method. The results show that the improved method of this invention exhibits significant performance advantages in all wind speed ranges. Compared with the global single-model method, the RMSE of the method of this invention is reduced from 1.683 m / s to 1.366 m / s, an improvement of 18.8%, and the R² is increased from 0.6016 to 0.7336, an improvement of 21.9%. Compared with the traditional four-part interval modeling method, the RMSE is improved by 5.6%; compared with the three-part interval modeling method with equal wind speed values, the RMSE is improved by 11.0%; and compared with the three-part interval modeling method with equal sample sizes, the RMSE is improved by 12.8%.

[0124] This invention achieves adaptive optimal division of wind speed ranges through the collaborative optimization of SA and GD algorithms. By combining XGB model construction and stacked ensemble learning multi-model fusion, it significantly reduces the significant inversion error differences of traditional single models in different wind speed ranges. At the same time, it solves the problem of difficulty in matching observations with different range models in practical applications, effectively improving the accuracy and robustness of GNSS-R sea surface wind speed inversion.

[0125] The above description discloses only one or more preferred embodiments of the present invention, and should not be construed as limiting the scope of the present invention. Those skilled in the art will understand that all or part of the processes of the above embodiments can be implemented, and equivalent changes made in accordance with the claims of the present invention are still within the scope of the invention.

Claims

1. A GNSS-R sea surface wind speed retrieval method, characterized in that, It comprises the following steps: Step 1: Obtain satellite-borne GNSS-R observation data and ERA5 reanalysis data, and perform data solving on the obtained CYGNSS satellite observation data; Step 2: Perform quality control on the solved CYGNSS data, perform interpolation processing on the ERA5 reanalysis data, and perform spatio-temporal matching on the quality-controlled CYGNSS data and the interpolated ERA5 data to construct a data set; Step 3: Perform preliminary division of the wind speed interval by using the improved SA algorithm for multi-objective optimization; Step 4: Further fine-tune the boundaries of the preliminary divided wind speed interval by using the GD algorithm; Step 5: Construct a special XGB prediction model for each divided wind speed interval; Step 6: Fuse the wind speed prediction results of each interval model by using a stacked ensemble learning architecture to output the final sea surface wind speed inversion value.

2. The GNSS-R sea surface wind speed inversion method of claim 1, wherein the CYGNSS satellite observation data includes a normalized bistatic radar cross section and an integrated delay waveform front slope, and the satellite observation data needs to be solved; In step 1, six other feature variables are also considered, namely the longitude of the specular reflection point, the latitude of the specular reflection point, the incident angle of the specular reflection point, the gain of the receiving antenna at the specular reflection point, the distance from the transmitting point to the specular reflection point, and the distance from the receiving point to the specular reflection point, to improve the prediction accuracy of the model. The above six feature variables can be directly obtained from the satellite observation data; During the execution of step 1, the east wind and north wind component values of the 10-meter sea surface wind speed are extracted from the ERA5 reanalysis data set, and the vector sum is calculated to obtain the wind speed value at the corresponding position.

3. The GNSS-R sea surface wind speed inversion method of claim 1, wherein in step 2, the following quality control measures are taken for the two feature variables obtained by solving and the six feature variables directly obtained: removing observation samples containing missing values in the original data; removing observation samples with negative integrated delay waveform front slope values; removing observation samples located in land areas; removing observation samples with receiving antenna gain values at the specular reflection point less than 0; removing observation samples with incident angles of the specular reflection point greater than 60°; calculating the receiver gain of each sample, and removing observation samples with receiver gain less than 10; The interpolation processing process of the ERA5 reanalysis data is as follows: the ERA5 data is subjected to time and space interpolation processing, the time resolution and spatial resolution of the ERA5 data are adjusted to be consistent with the time and spatial resolution of the CYGNSS data by bilinear interpolation, and then the CYGNSS data and the ERA5 data with consistent time and space information are matched in time and space to form a data pair, and all the matched data pairs constitute a data set.

4. The GNSS-R sea surface wind speed inversion method of claim 1, wherein ​ ​ In step 3, the SA algorithm generates perturbations and calculates the objective function value corresponding to the new boundary for each iteration by setting the initial temperature, cooling rate, and initial weight, decides whether to accept the newly divided interval according to the acceptance probability, and updates the weight parameter and temperature, returning the optimal interval boundary that minimizes the objective function value; A multi-objective optimization function is constructed by considering the characteristics of inversion error distribution, the number of intervals, the uniformity of sample distribution, and the minimum sample requirement of the interval, including the total sum of inversion error of each interval, the penalty term of the number of intervals, the penalty term of sample distribution variance, and the minimum sample penalty term. A multi-objective optimization strategy based on fixed normalization range is adopted to ensure that each factor can play a normal adjusting role in the optimization process.

5. The GNSS-R sea surface wind speed inversion method of claim 1, wherein, In step 4, the GD algorithm uses numerical differentiation to calculate the gradient, specifically by adding a small perturbation to the boundary value and observing the change in the objective function to estimate the gradient, using a fixed learning rate to control the step size of each update, and updating the boundary value according to the corresponding rules.

6. The GNSS-R sea surface wind speed inversion method of claim 1, wherein, In step 5, the objective function of the XGB model corresponding to each interval includes the loss of wind speed prediction value and true value and the complexity regular term, which is used to control the model complexity and prevent overfitting.

7. The GNSS-R sea surface wind speed inversion method of claim 1, wherein, In step 6, the two-level stacked ensemble model structure is used to improve the prediction performance, where the primary model is the base model, which includes N heterogeneous models trained using different wind speed interval samples; the secondary model is the meta-model, which also uses the XGB model to integrate the base model output, with the input features being the combination of base model prediction results and original features; For the training of the meta-model, a meta-feature matrix is generated through K-fold cross-validation, the training set is divided into K subsets, one subset is reserved as the validation set each time, and the remaining subsets are used to train the base model and generate meta-features, the meta-features of all folds are integrated to train the meta-model, avoiding data leakage; The XGB model is used as the meta-model, and the model complexity is controlled by the regularization coefficient, and the final wind speed inversion value is output.

Images