Global Mean Sea Surface Height Estimated from Spaceborne Cyclone-GNSS Reflectometry

Mean sea surface height (MSSH) is an important parameter, which plays an important role in the analysis of the geoid gap and the prediction of ocean dynamics. Traditional measurement methods, such as the buoy and ship survey, have a small cover area, sparse data, and high cost. Recently, the Global Navigation Satellite System-Reflectometry (GNSS-R) and the spaceborne Cyclone Global Navigation Satellite System (CYGNSS) mission, which were launched on 15 December 2016, have provided a new opportunity to estimate MSSH with all-weather, global coverage, high spatial-temporal resolution, rich signal sources, and strong concealability. In this paper, the global MSSH was estimated by using the relationship between the waveform characteristics of the delay waveform (DM) obtained by the delay Doppler map (DDM) of CYGNSS data, which was validated by satellite altimetry. Compared with the altimetry CNES_CLS2015 product provided by AVISO, the mean absolute error was 1.33 m, the root mean square error was 2.26 m, and the correlation coefficient was 0.97. Compared with the sea surface height model DTU10, the mean absolute error was 1.20 m, the root mean square error was 2.15 m, and the correlation coefficient was 0.97. Furthermore, the sea surface height obtained from CYGNSS was consistent with Jason-2′s results by the average absolute error of 2.63 m, a root mean square error (RMSE) of 3.56 m and, a correlation coefficient (R) of 0.95.


Introduction
About three-quarters of the Earth's surface is covered by the sea, and the rise of sea level affects the survival of the land and tropical coastal areas. These changes will bring huge economic losses and human casualties [1]. Meanwhile, the mean sea surface height (MSSH) is one of the concerns of marine science and environmental science, which is widely used by geodesists and geophysicists to analyze the geoid gap, Earth's internal dynamic mechanism, and the prediction of ocean dynamics. Oceanographers use the mean sea surface as a reference to study sea surface height, ocean circulation and other issues [2,3]. Traditional ocean observation methods, such as ship survey, buoy and coastal monitoring station, have high cost, sparse data and poor repeatability, while altimetry satellites only measure subsatellite points with low spatial and temporal resolution, and cannot obtain ocean information with large spatial scale in a short time range [4]. As early as 1968, Barrick et al. (1968) studied and obtained the scattering characteristics of rough surfaces [5].
nearly elliptic) and Doppler line (iso-Doppler, close to hyperbolic), and maps the light areas of each point to Doppler and Delay space and can get a Delay Doppler Map. GNSS signal is a kind of direct sequence spread spectrum signal, the satellite signal distribution within a wide frequency band, and, as a result of the limitation of the satellite launch power the space of long-distance transmission caused by the free space attenuation, and the existence of reflection, scattering the signal will leave high gain buried in noise, only through correlation processing to complete after capturing and tracking the measurement, as shown in Figure 1. Due to the reflection surface roughness, the reflection signals show the signal amplitude attenuation and different time Delay and Doppler signal superposition. The different delay and Doppler are opposite to the different reflection units of the reflection surface. The reflection signal correlation value, therefore, needs, from two aspects of time delay and the Doppler frequency shift, the delay-Doppler correlation function Y , , , see formula (1), and it reflects that the reflex zones at each time delay line and Doppler wire cross-area reflected signal correlation value are the most comprehensive description of the reflected signal.
where, is the integration time, is the output signal receiving antenna, is the local pseudo random noise (PRN) code copy code, is the center of the received signal frequency, is the reference point estimate Doppler. Each CYGNSS satellite had four DDMIs to observe four different highlights, so ideally CYGNSS could obtain 32 DDMs per second. The development of CYGNSS provided us with the first opportunity to estimate the global mean sea surface height (MSSH) based on the DDM obtained from the CYGNSS.

Specular Reflection Point Computation
Sea surface height (SSH) is the height from sea level to the WGS84 ellipsoidal surface. As shown in Figure 2, the CYGNSS needs to find its reference point on the WGS84 ellipsoidal surface after receiving signals from the effective reflection area to further obtain the average sea surface height of the effective reflection area. In Figure 2, is the receiver and is the GNSS satellites, point is the mirror reflection point, point is the center of the Earth, and is the radius of the Earth. The specular reflection point is on the Earth's surface. Of all the points on the Earth's surface, the distance from the transmitter to the point and then to the receiver is the shortest. The specular reflection point satisfies three conditions: (1) the specular reflection point is on the Earth's surface; (2) the shortest The reflection signal correlation value, therefore, needs, from two aspects of time delay and the Doppler frequency shift, the delay-Doppler correlation function Y R−Delay−Dop (t 0 , τ, f ), see formula (1), and it reflects that the reflex zones at each time delay line and Doppler wire cross-area reflected signal correlation value are the most comprehensive description of the reflected signal.
where, T i is the integration time, u R is the output signal receiving antenna, a is the local pseudo random noise (PRN) code copy code, f L is the center of the received signal frequency,f R is the reference point estimate Doppler. Each CYGNSS satellite had four DDMIs to observe four different highlights, so ideally CYGNSS could obtain 32 DDMs per second. The development of CYGNSS provided us with the first opportunity to estimate the global mean sea surface height (MSSH) based on the DDM obtained from the CYGNSS.

Specular Reflection Point Computation
Sea surface height (SSH) is the height from sea level to the WGS84 ellipsoidal surface. As shown in Figure 2, the CYGNSS needs to find its reference point on the WGS84 ellipsoidal surface after receiving signals from the effective reflection area to further obtain the average sea surface height of the effective reflection area. In Figure 2, RX is the receiver and TX is the GNSS satellites, point SP is the mirror reflection point, point O is the center of the Earth, and r is the radius of the Earth. The specular reflection point is on the Earth's surface. Of all the points on the Earth's surface, the distance from the transmitter to the point and then to the receiver is the shortest. The specular reflection point satisfies three conditions: (1) the specular reflection point is on the Earth's surface; (2) the shortest path Remote Sens. 2020, 12, 356 4 of 16 between the transmitter, receiver, and specular reflection point; (3) the incident angle at the specular reflection point is equal to the reflected angle. According to these conditions, the mirror reflection point position can be roughly estimated, and then based on the angle bisector variable step length correction of specular reflection point position, through iterative calculation of the precise positioning of the specular reflection, compared to Gleason's [17] and Wu's [18] method, it can be the results of high precision by the faster iteration [19].
Remote Sens. 2020, 12, 356 4 of 17 path between the transmitter, receiver, and specular reflection point; (3) the incident angle at the specular reflection point is equal to the reflected angle. According to these conditions, the mirror reflection point position can be roughly estimated, and then based on the angle bisector variable step length correction of specular reflection point position, through iterative calculation of the precise positioning of the specular reflection, compared to Gleason's [17] and Wu's [18] method, it can be the results of high precision by the faster iteration [19]. As shown in Figure 3, is the receiver, is the GNSS, point is the initial position of the specular reflection point, point is the center of the Earth, is the Earth's radius, ̂ is the correct direction, is correction step, is a new mirror point location, based on the geometric relation in the diagram. We calculated the specular reflection point steps, which are as follows: (1) Enter the location of the receiver and transmitter obtained from the navigation information.
where, and are the heights of receiver and GNSS satellite, respectively. In the WGS84 coordinate system, the long radius is 6,378,137 m and the eccentricity is 0.081819190 84264.
(3) Correct direction ̂ calculated using Formulas (6) and (7). As shown in Figure 3, R is the receiver, T is the GNSS, point S is the initial position of the specular reflection point, point O is the center of the Earth, r is the Earth's radius,ŝ is the correct direction, K is correction step, S temp is a new mirror point location, based on the geometric relation in the diagram. We calculated the specular reflection point S steps, which are as follows: (1) Enter the location of the receiver and transmitter obtained from the navigation information.
(2) Using the formulas (2-5) roughly use S as the initial position. where, H R and H T are the heights of receiver and GNSS satellite, respectively. In the WGS84 coordinate system, the long radius a WGS84 is 6,378,137 m and the eccentricity e WGS84 is 0.081819190 84264. (3) Correct directionŝ calculated using Formulas (6) and (7).
(5) Formulas (9) and (10) are used to calculate the new mirror point position S temp .
(5) Formulas (9) and (10) are used to calculate the new mirror point position .  After calculating the position of the specular reflection point on the WGS84 ellipsoid, the time delay and Doppler frequency shift of the point can be calculated, which is the first step to calculate the mean sea surface height (MSSH).

Tracking Point Delay Estimation Method
The previous section introduced the estimation method of the WGS84 specular reflection point. In the actual measurement, to obtain the delay difference, it is necessary to estimate the delay of the actual reflection, which is the delay of the shortest path reflected to the sea surface. Cardellach et al. proposed three methods to obtain the actual delay [20], all of which were obtained based on the delayrelated waveform, as follows: (1) First derivative spectrophotometry This estimation method is to put the actual delay as the solution of Equation (11). After calculating the position of the specular reflection point on the WGS84 ellipsoid, the time delay and Doppler frequency shift of the point can be calculated, which is the first step to calculate the mean sea surface height (MSSH).

Tracking Point Delay Estimation Method
The previous section introduced the estimation method of the WGS84 specular reflection point. In the actual measurement, to obtain the delay difference, it is necessary to estimate the delay of the Remote Sens. 2020, 12, 356 6 of 16 actual reflection, which is the delay of the shortest path reflected to the sea surface. Cardellach et al. proposed three methods to obtain the actual delay [20], all of which were obtained based on the delay-related waveform, as follows: (1) First derivative spectrophotometry This estimation method is to put the actual delay τ spec as the solution of Equation (11).
In the equation, W is delay relative waveforms (DM), which can be caused by the time delay Doppler figure (DDM). We can know that this point from the formula delays is the peaks of the wave's first derivative.
(2) Rule of Maximum Area Another tracking point delay estimation method is based on the GNSS navigation time delay estimation method. In GNSS navigation, the time delay of the peak power is treated as the time delay of the navigation signal. When the surface is smooth, the algorithm only involves mirror delay, in which case, the tiny surface roughness does not affect the shape of the waveform. However, in Marine GNSS reflection measurement, the power waveform and peak are essentially determined by the roughness. Therefore, this method can be used in Global Navigation Satellite System-Reflectometry (GNSS-R) altimetry if the surface is smooth enough [21].

(3) Half Peak Method
This algorithm is a simplified version of the algorithm used in single base radar altimetry. At the front edge of the waveform, the time delay corresponding to a part of the peak power is regarded as the actual time delay. In the altimetry measurement of the radar, the coefficient is very close to 0.5, which is the reason why the algorithm is named. In GNSS-R, this coefficient is generally selected as 0.75 because its result is very close to that of the first-order derivative method. It is important to note that when the coefficient is very close to 1, like the peak itself, it is strongly affected by the roughness of the ocean surface.
In order to reduce the inversion error of the actual Delay, the method by Mashburn et al. can be used to perform preliminary processing to improve the results. In other words, the interpolation processing of the DM can be performed using the Whitker-Shannon theorem [22], as shown in Formula (12).
where the x[n] is related to the power sequence of the original time delay, δ is the extended time delay sequence, and d is the sampling period. By comparing the three methods of actual time delay estimation, the influence of different inversion methods on the measurement accuracy is analyzed, and finally the first method is selected for the subsequent sea surface height calculation.

Sea Surface Height Calculation
After calculating the delay of the WGS84 specular reflection point and the actual delay, according to the geometric relationship in Figure 2, the sea surface height (SSH) of the effective reflection area can be calculated, according to Formulas (13)(14)(15)(16) [23] proposed by Clarizia et al.
Remote Sens. 2020, 12, 356 Because the effective scattering area is very small, the sea level change in this very small area is also very small. This inversion method takes advantage of this feature and assumes that the longitude and latitude of the tracking points are the same as those of the WGS84 specular reflection points, thus reducing unknown variables to sea surface height (SSH) and achieving the goal of inversion sea surface height.

MSSH from CYGNSS
The mean sea surface height (MSSH) is the focus of current Earth science and environmental science. Compared with the reference ellipsoid, it contains the information of geoid and sea surface topography, so it is widely used to study geoids, instantaneous sea surface height, crustal deformation, ocean circulation, sea level change, and other issues. The global mean sea surface height model is obtained by using abundant satellite altimetry-measuring data (ERS-1, ERS-2, Jason series, etc.) and corresponding data processing. Generally, the model is based on grid data, which can be interpolated to obtain the sea surface height of a certain location in the world. The global mean sea surface height (MSSH) models used in this paper are CNES_CLS2015, provided by AVISO and DTU-10 models.
Since the resolution of the MSSH model is different from that of CYGNSS inversion, data matching is required for comparison. The bilinear interpolation method can be used to match two points. The lon 0 and lat 0 are latitude and longitude, respectively, and the inversion results are H measured . The position of this point on the MSSH grid is shown in Figure 4. According to Formulas (17)(18)(19), the H measured of the sea surface height of the SP point in the MSSH grid can be calculated.
Among them, H Q 11 , H Q 21 , H Q 12 , and H Q 22 can be found in the MSSH model grid, according to lon 1 , lon 2 , lat 1 , and lat 2 . After matching, a series of corresponding lon 0 , lat 0 , H measured , and H model can be evaluated. In this paper, the bias average absolute error (MAE), root mean square error (RMSE), and correlation coefficient (R) are calculated by using Formulas (20)(21)(22)(23). interpolated to obtain the sea surface height of a certain location in the world. The global mean sea surface height (MSSH) models used in this paper are CNES_CLS2015, provided by AVISO and DTU-10 models.
Since the resolution of the MSSH model is different from that of CYGNSS inversion, data matching is required for comparison. The bilinear interpolation method can be used to match two points. The and are latitude and longitude, respectively, and the inversion results are . The position of this point on the MSSH grid is shown in Figure 4. According to Formulas (17)(18)(19), the of the sea surface height of the point in the MSSH grid can be calculated. , , , Among them, , , , and can be found in the MSSH model grid, according to , , , and . After matching, a series of corresponding , , , and can be evaluated. In this paper, the bias average absolute error (MAE), root mean square error (RMSE), and correlation coefficient (R) are calculated by using Formulas (20)(21)(22)(23). ,

Comparison with AVISO
The first global mean sea surface height model used in this paper is CNES_CLS2015, provided by AVISO. The model uses 20 years (1993-2012) of altimetry data, mainly including Topex/Poseidon, ERS-2, GFO, Jason-1, Jason-2, Envisat, ERS-1, Jason-1, and Cryosat-2. The model currently offers a 1 min resolution global grid product available on the AVISO website. In addition to CNES_CLS2015, AVISO has previously released several related products, including CLS_SHOM98.2, CNES_CLS01, CNES_CLS10, and CNES_CLS2011. In this study, CYGNSS L1 was used in August 2017 to obtain the global mean sea surface height distribution, mainly including bistatic radar cross section (BRCS) data, spacecraft position, specular point position, GPS spacecraft position, and the incidence angle using the method in Section 2. We firstly used the position data named 'tx_pos' and 'sc_pos' to compute the position of the specular point on WGS84, according to the method in Section 2.2, and then we used the data named 'sp_pos', 'sp_lat', and 'sp_lon' to check our results. After the tested position was computed, we got the specular delay simply. In this experience, the spatial resolution was 0.1 • × 0.1 • , and each grid had about 30 measurements. The sea surface high monitoring area was the maximum range that the CYGNSS constellation could reach, including the full longitude range, and the latitude range was about 40 • N~40 • S. The average sea surface height distribution of the monitored sea areas obtained by the CNES_CLS2015 model is shown in Figure 5.
The inversion results were generally consistent with the CNES_CLS2015 model. The average sea surface height in the western Pacific is high and positive, which means that the average sea surface height here was higher than the WGS84 ellipsoid, and the average sea surface height in the Northern Indian Ocean is low and negative, which means that the average sea surface here was lower than the WGS84 ellipsoid. Considering Section 3.2, about the comparison of evaluation methods, as shown in Figure 6, you can see CYGNSS estimation results and CNES_CLS2015 model mainly concentrated in the deviation plus or minus 10 m and a small amount of distribution in the plus or minus 10 to 20 m. The graphics are very close to a Gaussian probability density distribution, and, at the same time, two set of values of the average absolute error was 1.33 m, root mean square error (RMSE) was 2.26 m, and the correlation coefficient (R) was 0.97, which is very similar in both groups. The results of the inversion results within a certain precision are credible. The inversion results were generally consistent with the CNES_CLS2015 model. The average sea surface height in the western Pacific is high and positive, which means that the average sea surface height here was higher than the WGS84 ellipsoid, and the average sea surface height in the Northern Indian Ocean is low and negative, which means that the average sea surface here was lower than the WGS84 ellipsoid. Considering Section 3.2, about the comparison of evaluation methods, as shown in Figure 6, you can see CYGNSS estimation results and CNES_CLS2015 model mainly concentrated in the deviation plus or minus 10 m and a small amount of distribution in the plus or minus 10 to 20 m. The graphics are very close to a Gaussian probability density distribution, and, at the same time, two set of values of the average absolute error was 1.33 m, root mean square error (RMSE) was 2.26 m, and the correlation coefficient (R) was 0.97, which is very similar in both groups. The results of the inversion results within a certain precision are credible.

Comparison with DTU-10
Another global mean sea surface height model used in this paper is DTU-10, which is the global mean sea surface height model proposed by Andersen et al., from Denmark Technical University in 2010, which is based on the DNSC08 model and improved. The DNSC08 model was also proposed by Andersen et al., which mainly used nine altimetry satellites, including Jason-1, T/P, T/P interleaved mission, ERS-1 GM, ERS-2 Exact Repeat Mission (ERM), Geosat GM, Geosat follow-on (GFO)-ERM, Envisat ERM, and ICESat, and 12 years of observation data were averaged. Compared with DNSC08, DTU-10 extended data coverage, averaged 17 years of data, and took ENVISAT data from ERS-2 and ENVISAT again. Dtu-10 data mainly include 1 min and 2 min resolution global grid products, and

Comparison with DTU-10
Another global mean sea surface height model used in this paper is DTU-10, which is the global mean sea surface height model proposed by Andersen et al., from Denmark Technical University in 2010, which is based on the DNSC08 model and improved. The DNSC08 model was also proposed by Andersen et al., which mainly used nine altimetry satellites, including Jason-1, T/P, T/P interleaved mission, ERS-1 GM, ERS-2 Exact Repeat Mission (ERM), Geosat GM, Geosat follow-on (GFO)-ERM, Envisat ERM, and ICESat, and 12 years of observation data were averaged. Compared with DNSC08, DTU-10 extended data coverage, averaged 17 years of data, and took ENVISAT data from ERS-2 and ENVISAT again. Dtu-10 data mainly include 1 min and 2 min resolution global grid products, and the specific products can be obtained from the official website of the DTU-10 server. Similar to Section 3.2, Figures 7 and 8 can be obtained. The sea surface height obtained through CYGNSS inversion is consistent with the value of the DTU10 model, because DTU10 model and CNSE_CLS2015 model have strong consistencies. It can be seen that the differences between CYGNSS results and the DTU10 model are mainly concentrated within 10 m and the distributions are very close to a Gaussian probability density distribution. The average absolute error was 1.20 m, the root mean square error (RMSE) was 2.15 m, and the correlation coefficient (R) was 0.97.

Comparison with Jason-2
Jason-2/Ocean Surface Topography Mission (OSTM) took over and continued Topex/Poseidon and Jason-1 missions in 2008, in the frame of a cooperation between Center National d'Etudes Spatiales (CNES), Eumetsat, NASA and National Oceanic and Atmospheric Administration (NOAA). It carried two predecessors for a high-accuracy altimetry mission: a Poseidon-class altimeter, a radiometer, and three location systems. Its orbital inclination was 66° and the altitude was 1336 km. From 2008 to October 2016, Jason-2 was located on its nominal orbit. From October 2016 (at the end of cycle 303 and until cycle 327), after more than 8 years of service on this nominal ground track, Jason-2 swifted to the interleaved orbit. From July 2017 (beginning cycle 500), Jason-2 operated on a new long-repeat orbit (LRO) at roughly 1309.5 km altitude. From July 2018 (beginning cycle 600), Jason-2 operated on an interleaved long repeat orbit (i-LRO) for the second geodetic cycle, on a ground track in the middle of the grid defined by the first geodetic cycle. The sea surface height anomalies (SSHA) data can be download from AVISO, including the information of mean sea surface. The reference of Jason-2 is CLS01 MSS. Using Formula (24), SSH data can be calculated. , where, is the reference that can be got from the data. Data for August 2017 were collected to calculate the mean sea surface and then compared with the CYGNSS result. Similar to Sections 3.2 and 3.3, Figures 9 and 10 can be obtained. The sea surface height obtained through CYGNSS inversion

Comparison with Jason-2
Jason-2/Ocean Surface Topography Mission (OSTM) took over and continued Topex/Poseidon and Jason-1 missions in 2008, in the frame of a cooperation between Center National d'Etudes Spatiales (CNES), Eumetsat, NASA and National Oceanic and Atmospheric Administration (NOAA). It carried two predecessors for a high-accuracy altimetry mission: a Poseidon-class altimeter, a radiometer, and three location systems. Its orbital inclination was 66 • and the altitude was 1336 km. From 2008 to October 2016, Jason-2 was located on its nominal orbit. From October 2016 (at the end of cycle 303 and until cycle 327), after more than 8 years of service on this nominal ground track, Jason-2 swifted to the interleaved orbit. From July 2017 (beginning cycle 500), Jason-2 operated on a new long-repeat orbit (LRO) at roughly 1309.5 km altitude. From July 2018 (beginning cycle 600), Jason-2 operated on an interleaved long repeat orbit (i-LRO) for the second geodetic cycle, on a ground track in the middle of the grid defined by the first geodetic cycle. The sea surface height anomalies (SSHA) data can be download from AVISO, including the information of mean sea surface. The reference of Jason-2 is CLS01 MSS. Using Formula (24), SSH data can be calculated.
where, SSH is the reference that can be got from the data. Data for August 2017 were collected to calculate the mean sea surface and then compared with the CYGNSS result. Similar to Sections 3.2 and 3.3, Figures 9 and 10 can be obtained. The sea surface height obtained through CYGNSS inversion is consistent with the value of Jason-2 s results. It can be seen that the deviation of CYGNSS estimation results and Jason-2 s results were mainly concentrated in the plus or minus 10 m and a small amount of distribution in the plus or minus 10 to 20 m. The graphics were very close to a Gaussian probability density distribution, and, at the same time, two set of values of the average absolute error were 2.63 m, the root mean square error (RMSE) was 3.56 m, the correlation coefficient (R) was 0.95, that value was very similar in both groups. The results of inversion results within a certain precision are credible.
Remote Sens. 2020, 12, 356 13 of 17 is consistent with the value of Jason-2′s results. It can be seen that the deviation of CYGNSS estimation results and Jason-2′s results were mainly concentrated in the plus or minus 10 m and a small amount of distribution in the plus or minus 10 to 20 m. The graphics were very close to a Gaussian probability density distribution, and, at the same time, two set of values of the average absolute error were 2.63 m, the root mean square error (RMSE) was 3.56 m, the correlation coefficient (R) was 0.95, that value was very similar in both groups. The results of inversion results within a certain precision are credible.

Errors Discussion
In the whole sea surface height inversion process, the biggest error comes from the inversion accuracy of the actual time delay, as shown in Figure 4. Formula (25) can be obtained from the geometric relationship.
Take the partial derivative and get Formula (26).
• cos cos , Due to SSH ≪ , , ≈ , and ≈ , Equation (26) can be simplified to Equation (27). The error of actual measurement delay directly affects the inversion of mean sea level height.

Errors Discussion
In the whole sea surface height inversion process, the biggest error comes from the inversion accuracy of the actual time delay, as shown in Figure 4. Formula (25) can be obtained from the geometric relationship.
Take the partial derivative and get Formula (26).
Due to SSH R TS , R RS , R TS ≈ R TS , and R RS ≈ R RS , Equation (26) can be simplified to Equation (27). The error of actual measurement delay directly affects the inversion of mean sea level height.
Due to the delay of the CYGNSS resolution, which was about 0.255/1023000 s, according to Formula (27), it can be calculated for the lowest height inversion effect, which was about 37 m. When the incident angle was small, the time Delay error of height error was bigger, based on one dimensional time delay waveform figure (delay map) inversion when actual Delay. In order to reduce the actual inversion error of time delay, the Whitaker-Shannon theorem of DM interpolation processing can improve the result.
Except for the time delay estimation error and the receiver position, the GNSS position error was also bigger, because the location precision of the CYGNSS providing for level meters. This will bring the ideal mirror reflection point estimation error and can eventually bring the mean sea level height inversion error, which was about 3 m. Using the precise ephemeris and using appropriate methods to remove the ionosphere, troposphere, tidal, and pressure on the influence of orbit determination, this error was greatly reduced. Finally, the inherent errors of the instrument should also be considered. The noise of the instrument also affected the time delay, orbit determination, and ultimately the inversion of the mean sea level. With improvement of GNSS-R receivers and more spaceborne GNSS-R missions in the future, precise ocean environmental remote sensing will be expected from spaceborne GNSS-R observations [24][25][26][27][28].

Conclusions
In this paper, DDM data, provided by CYGNSS space-borne GNSS-R observations was used to invert the time delay of reflected signals according to the relation between peak position, half-peak position, and first-order derivative peak position of DM and the propagation time of reflected signals. In August 2017, the MSSH was obtained from CYGNSS data, which are validated by satellite altimetry CNES_CLS2015, DTU10 and Janson-2 results. The mean absolute error between the MSSH from the CYGNSS and the MSSH from the CNES_CLS2015 model was 1.33 m, the root mean square error was 2.26 m, and the correlation coefficient was 0.97. Compared with the DTU10 sea surface height model, the mean absolute error was 1.20 m, the root mean square error was 2.15 m, and the correlation coefficient was 0.97. Compared with Jason-2 results in August 2017, the mean absolute error was 2.63 m, the root mean square error was 3.56 m, and the correlation coefficient was 0.95. Therefore, the MSSH estimation from the CYGNSS can provide important support for marine shipping development, marine environmental protection, marine disaster warnings, and forecasting, etc. Here, the precision was global, and the performance in each grid will be further investigated in the future.
The CYGNSS MSSH estimation method requires high spatial resolution, and the time required for the grid to produce high-precision data is also long, so the time resolution of the obtained estimation value is low. Furthermore, the influence of other satellite parameters, such as the satellite orbital error, should be analyzed in the future in order to reduce these effects. Due to the covering limitation of the CYGNSS observation area, the CYGNSS can only estimate the MSSH of the area within 40 degrees north and south latitude. In addition, the estimation of MSSH based on the CYGNSS just considered one single parameter. In the future, multi-parameters, such as sea surface wind speed and significant wave height, are brought into the MSSH estimation, which may improve the accuracy of the MSSH estimation.