Sea Ice Detection Using GNSS ‐ R Data From TechDemoSat ‐ 1

A new method for the detection of sea ice using GNSS ‐ R (Global Navigation Satellite Systems Re ﬂ ectometry) is presented and applied to 33 months of data from the U.K. TechDemoSat ‐ 1 mission. This method of sea ice detection shows the potential for a future GNSS ‐ R polar mission, attaining an agreement of over 98% and 96% in the Antarctic and Arctic, respectively, when compared to the European Space Agency's Climate Change Initiative sea ice concentration product. The algorithm uses a combination of two parameters derived from the delay ‐ Doppler Maps to quantify the spread of power in delay and Doppler. Application of thresholds then allows sea ice to be distinguished from open water. Differences between the TechDemoSat ‐ 1 sea ice detection and comparison data sets are explored. The results provide information on the seasonal and multiyear changes in sea ice distribution of the Arctic and Antarctic. Future potential and applications of this technique are discussed.


Introduction
Monitoring sea ice is important to observe and understand global climate processes and their changes, due to the dynamic role that sea ice plays in the physics of the Earth's climate, not least in altering the Earth's albedo (Stroeve et al., 2011). As a key element of the climate system, sea ice also affects the ocean circulation and weather (Comiso et al., 2003). From a practical perspective, there is also interest in detecting and mapping sea ice for maritime operations such as exploitation of oil and gas reserves and to determine possible shipping routes, especially in the Arctic (Bird et al., 2008).
Global Navigation Satellite Systems Reflectometry (GNSS-R) uses reflected microwave navigation signals (L-Band from GPS or similar systems) to characterize the reflective surface and was first proposed in 1988 (Hall & Cordey, 1988). This technique is at present largely used to study winds over the ocean (e.g., Foti et al., 2015;Foti et al., 2017;Ruf & Balasubramaniam, 2018). However, other Earth observation applications, such as ice sheet and ocean altimetry (e.g., Cardellach et al., 2004;Cartwright et al., 2018;, the monitoring of soil moisture (Chew et al., 2016), and the cryosphere (e.g., Belmonte Rivas et al., 2010;Gleason, 2006), are also being investigated. Cryospheric applications include the detection and altimetry of sea ice (e.g., Alonso-Arroyo et al., 2017;Hu et al., 2017;Li et al., 2017;Yan & Huang, 2016), the altimetry of glacial ice (e.g., Cardellach et al., 2004;Cartwright et al., 2018;Rius et al., 2017), and permittivity of reflective surfaces (Belmonte Rivas et al., 2010;Fabra et al., 2012). The use of these signals of opportunity and the need only for a receiver in space mean that GNSS-R is an extremely low-cost and, being a passive system, low-power method of remote sensing. Combined with the low mass of the sensor itself, this makes the approach ideal for applications in constellation missions, as in the case of CYGNSS (Cyclone GNSS), launched in 2016, for monitoring of tropical cyclone winds (Ruf

10.1029/2019JC015327
Key Points: • Sea ice detection is presented from satellite GNSS-R over the Arctic and Antarctic • Sea ice presence is predicted based on multiple parameters across the first 33 months of data from U.K. TechDemoSat-1 • Results show excellent agreement with European Space Agency Climate Change Initiative sea ice product as well as other extent products Supporting Information: • Supporting Information S1 • Movie S1 S5 et al., 2013). The eight CYGNSS microsatellites achieve an average revisit time of 4 hr (Ruf et al., 2013). Such a high temporal resolution, combined with the spatial coverage enabled by this quantity of measurements, would allow the monitoring of the dynamic sea ice system. However, the low inclination (35°) orbit of the CYGNSS satellites means that they cannot be used to monitor sea ice. A similar system for cryosphere applications (G-TERN) has been proposed (Cardellach et al., 2018). In orbit and operational until recently, TechDemoSat-1 (TDS-1) carries a GNSS-R receiver and is in a high inclination orbit (98.4°). Data from it can be used to detect sea ice, and thus, these data are the focus of this paper.
Previous studies (Alonso-Arroyo et al., 2017;Rivas et al., 2010;Yan et al., 2017;Yan & Huang, 2016;Zhu et al., 2017) have been shown to correctly distinguish between sea ice and open ocean in up to 98.4% of cases in the detection of sea ice compared to colocated passive microwave data. However, these studies were limited in scope, examining small subsets of TDS-1 observations only. This study builds upon and extends these previous studies by developing and applying a new method of sea ice detection to the complete data set collected during the initial mission of TDS-1. In what follows, section 2 discusses the satellite, the sensor, and the data used and then section 3 describes the sea ice detection methods tested. Section 4 validates the results of the algorithm against the most recent ESA CCI (European Space Agency Climate Change Initiative) sea ice concentration product (Toudal Pedersen et al., 2017) as well as the National Snow and Ice Data Centre (NSIDC) sea ice concentration product (Cavalieri et al., 1996), Ocean and Sea Ice Satellite Application Facility (OSI SAF) ice edge product (Breivik et al., 2012;OSI SAF, n.d.), and Multisensor Analyzed Sea Ice Extent (MASIE) operational ice edge product (Fetterer et al., 2010). The likely sources of differences are then discussed. Conclusions of the study and a discussion of the future potential of GNSS-R for sea ice monitoring are given in section 5.

TechDemoSat-1 and GNSS-R
Data from the GNSS-R sensor on TDS-1 are provided in the form of delay-Doppler maps (DDMs). These represent the scattered power in the delay and Doppler domains, with a smooth surface yielding a coherent reflection with the majority of the scattered power coming from the specular point, and very little from the glistening zone around it (Zavorotny & Voronovich, 2000). This is visible in the DDM as a limited spread in the power distribution and a strong peak. Roughness is considered at the scale of the wavelength (19 cm), with sea ice or flat leads producing a more coherent DDM as compared to one from the open ocean, as the latter results in more scattering in the delay and Doppler domains due to the presence of waves on the sea surface. This distinct difference in the distributions of scattering from these two types of surface allows information on their properties to be determined. Where reflections are coherent, such as over sea ice, the footprint of the sensor is believed to be significantly smaller than passive microwave at 6 km by 0.4 km along track and across track, respectively (Alonso-Arroyo et al., 2017;Hu et al., 2017). The elongation in the alongtrack dimension is primarily due to the onboard incoherent averaging (1 Hz) and, were this not applied, the footprint would be on the order of 400 m by 400 m.
The asynchronism of the orbit of TechDemoSat-1 with those of the GPS satellites, in combination with the ability of the sensor to receive up to four reflections at once, allows a varying web of specular points to be built up over time. This has the potential to provide high-resolution data in both time and space and therefore capture the dynamics of sea ice, even over the geographical North Pole.
TDS-1 was designed as a satellite technology demonstration platform by Surrey Satellite Technology Ltd. and carries eight separate payloads operated in an 8-day cycle. The data for this study were obtained from the SGR-ReSI, a low-power and low-mass sensor built from commercial components (full details can be found in  and in operation two days in every cycle. The satellite was launched in 2014 and placed into a quasi-Sun-synchronous orbit (with a local time ascending node drift of 1.42 hr/year; Foti et al., 2015) with an inclination of 98.4°at an altitude of 635 km. In 2017, it came to the end of its demonstration mission and has since been recommissioned with the SGR-ReSI in operation continuously. Only data acquired before the mission extension will be considered in this study as data of antenna gain less than 0 dB has not been downlinked during the 24-hr operational phase, eliminating much of the data at high latitudes. The SGR-ReSI creates 1-ms DDMs, which are averaged over 1-s intervals prior to downlinking. These are then provided along with the metadata for the same period in 6-hr collection windows and are available at merrbys.co.uk website. The DDMs consist of 128 delay bins and 20 Doppler bins, with resolutions of 0.252 chips (~0.246 μs) and 500 Hz, respectively.
The data considered in this study have been collected over almost 3 years of the initial TDS-1 mission (September 2014 to May 2017), using observations from both polar regions (north of 50°N and south of 50°S) and applying the land mask given in the ESA CCI data (Toudal Pedersen et al., 2017). The TDS-1 data were then filtered following the approach of Cartwright et al. (2018), which ensures the reduction of noise, as well as the removal of data where the reflection has left the tracking window. Data were also removed where flagged as containing the direct (unreflected) component of the signal. In addition, data from Collection Period 12 (September 2016) were discarded due to inconsistencies in the DDMs resulting from changes in sensor settings during that period. Contrary to previous studies, data were not filtered according to gain control mode, antenna gain nor incidence angle to retain the maximum amount of data, particularly over the North Pole. It should be noted that the antenna gain of data points in the Southern Hemisphere is naturally higher than that of the Northern Hemisphere due to the distance of the sea ice region from the South Pole. This therefore leads to lower incidence angles and concomitant higher antenna gains.
The recently released CCI sea ice concentration product (version 2.1) was used as a comparison to ensure the highest availability of high-quality reprocessed data for evaluation. This product is provided on a 25-km grid and was linearly interpolated to the location of the TDS-1 specular points. Any concentration above 0% was deemed to be sea ice and below this water. This threshold was used, rather than the conventional 15% concentration for sea ice coverage due to the sensitivity of the SGR-ReSI to the presence of sea ice. This choice was found to give the best agreement with derived parameters from the TDS-1 data and the largest separation of the two populations of DDMs associated with the different surface types.
Other sea ice products such as that of the NSIDC Special Sensor Microwave Imager Sounder (Cavalieri et al., 1996) were also tested and found to agree with the CCI product in 96% of cases using the 0% ice concentration threshold. This means that the CCI and Special Sensor Microwave Imager Sounder products agree 96% of the time regarding their classification of the surface as sea ice or open ocean. As such, any degree of agreement between the GNSS-R sea ice detection method and either data set that exceeds 96% is an excellent result, given that the best available comparable data disagree at that level.

Methods and Results
A number of parameters were derived from the DDMs and compared to the presence of sea ice at each observation point. A simple method (Figure 1) was applied to assess the extent to which each of the derived parameters could be used to detect sea ice through the use of a straightforward threshold. This was done using histograms to assess the point at which the threshold should be placed to maximize the probability of correct assignment (the point at which the histograms of the seawater and sea ice cross). The DDM-derived parameters tested have been defined below, many of which quantify the Trailing Edge Slope (TES). Although the DDMs from TDS-1 are not calibrated with regard to direct signal power, observables from the maximum-normalized DDMs were considered, such as the ratio of power in the signal box to that of the noise and the percentage of total power contained in signal boxes of varying sizes. Unfortunately, the lack of phase information prohibits the use of some previously used GNSS-R parameters such as DDM variance and Allan DDM variance (Clarizia et al., 2014). Journal of Geophysical Research: Oceans 1. OCOG (Offset Centre of Gravity). Similar to metrics used by Yan and Huang (2016). This parameter is based on that used in altimetry for the retracking of waveforms and identification of topographical changes (Wingham et al., 1986). The waveform is selected and recentered on the maximum power (in order to account for shifts of the maximum amplitude away from the zero-Doppler and zero-delay point).
The center of mass of the waveform is then identified and the difference in delay between this and the maximum of the waveform calculated. This has been calculated for the waveform taken through the maximum amplitudes of the DDM (varying in delay space) and the Doppler-and delay-integrated waveforms. 2. dx and dy. The difference between the maximum amplitude and the location on the DDM where the scattered power has decreased to a given fraction of its maximum. The position of the maximum power of the DDM is determined in the Doppler domain and the corresponding waveform selected (i.e., a slice through the DDM corresponding to the maximum power). The distance is then calculated between the maximum power and the point at which it has fallen to a nominal level. See Figure 2 for a visual depiction of this calculation. In this case 85% of the maximum was chosen moving down the DDM. Also investigated were 65%, 75%, and 95%. This was calculated in all four directions on the DDM (up and down: dy; left and right: dx); dy in the down direction can be considered another method of calculation of the TES. 3. DDMA (DDM average). As used by Alonso-Arroyo et al. (2017). The average value in a box around the maximum amplitude for the box sizes below. This average was computed both as a mean and a median for the maximum-normalized DDM as well as a median of the DDMs of proportional power (percentage of total DDM power contained in each of the Doppler and delay bins).
a 3 (Doppler) × 3 (delay) bins b 5 (Doppler) × 3 (delay) bins c 7 (Doppler) × 3 (delay) bins d 20 (Doppler) × 3 (delay) bins e κ (kurtosis). The kurtosis has been calculated for the DDM of scattered power, for the waveform in the delay domain through the maximum amplitude and for both the Doppler-and delay-integrated waveforms. f DDMS (DDM sum). Sum of values in a box around the maximum amplitude for the same box sizes considered above (DDMA). This is calculated from the normalized DDM. g TES. The gradient calculated from the maximum amplitude of the waveform (taken in the delay domain through the maximum amplitude) to the point on the waveform three, six, or nine bins later in delay. h Maximum power. This is the maximum amplitude of the DDM in scattered power. This measurement is taken directly from the supplied DDMs and is uncalibrated. i Noise floor. Average (median) value of a box in the signal-free region (first five delay bins, and all Doppler bins) of the normalized DDM (scattered power). j Signal-to-noise ratio. Ratio of the signal (DDMA median 5 × 3, as above) to the noise floor calculated in scattered power.
The parameters of OCOG (using the waveform through the maximum of the DDM) and dy (down the DDM to 85% of the maximum power) performed best in terms of agreement with the passive microwave comparison data set and data retention, in both the Arctic and Antarctic and are thus applied in this study. Full results of all parameters' combined performances can be found in the supporting information. Essentially, these parameters (OCOG and dy) quantify the peakiness of the waveform (sea ice waveforms are much peakier than those for open ocean) and will be used to distinguish whether the DDM is due to a reflection from sea ice or open ocean.
A training data set was extracted from the observations-in this case the first whole year of data in the Antarctic (2015), which is approximately 20% of the entire data set. The training data set was taken from the Antarctic rather than the Arctic, as the data are of a higher quality. This is due to the geographical difference in sea ice location, with a lower-latitude data set yielding lower incidence angles and higher antenna gains. The process was also performed with Arctic training data, but with lower agreement, so this approach was not pursued.
Thresholds were calculated on a month-by-month basis using the Antarctic training data in order to account for any seasonal changes in the ice surface. The median of these monthly thresholds was then calculated to give one threshold for each observable (OCOG = 0.2537 chips, dy = 0.4772 chips). Thresholds were then combined from the different derived parameters, with both criteria having to be met to determine the correct classification. This requirement results in the loss of some data (a further 7% after the filtering process) but improves the agreement with the comparison data set overall (94% or 84% when the respective OCOG or dy threshold is applied alone, rising to 98% where combined).

Validation Against ESA CCI Sea Ice Concentration
Where the combined thresholds are used, a high level of agreement is seen between the TDS-1 sea ice detection approach and that of the CCI data, at 98.4% in the Antarctic and 96.6% in the Arctic over the entire data set (98.2% in the Antarctic where the training data set is excluded). Assignment of points can be seen in Figure 3 for months as close as possible to the ice extent maximums and minimums for 2016 in the Arctic and Antarctic. Due to inconsistent DDMs, it is not possible to show the true minimum and maximum in September 2016 (this is Collection Period 12 and was subject to changes to the sensor settings) while 2015, as the training year would not have independently confirmed the success of the technique. Gaps in Figure 3 result from the varying amounts of data received in different months, due to the shared and experimental nature of TDS-1. Although the algorithm was trained using CCI data, similarly promising results are obtained when the results are compared to NSIDC detection (at >0% concentration; Cavalieri et al., 1996), OSI SAF ice edge product (derived from both active and passive sources; Breivik et al., 2012), and the operational product MASIE for Arctic data (Fetterer et al., 2010). Mapped comparisons of the ESA CCI, OSI SAF, and MASIE data sets can be found in the supporting information.
Figures 4 and 5 display the results of the validation as compared to data from the ESA CCI product. The results are broken down into categories as seen in Table 1. The seasonal cycles are well captured, with the sea ice maxima (largest proportion of correct positives) visible in the expected months, February and September for the Arctic and Antarctic respectively. This 6-month offset is also seen in the sea ice minima, with the largest proportion of correct negatives seen in September in the Arctic and in February in the Antarctic. The double peak and increased false negative seen in September 2016 (Collection Period 12) of Figure 5 can be attributed to the low sample numbers for this month, with only 3,000 useful data points in each hemisphere. It must be noted that the false negatives, although present in both hemispheres are such a small proportion of the total that their presence is not clearly visible in the figures. Other points to note are the poorer performance of the algorithm in the summer months in both hemispheres, with a large increase in the false positives in these months, at times to over 10% of the data. A possible cause of these differences in detection is the various sensitivities of the different sensors to the sea ice and the differing behavior of microwave emissions from the surface captured by passive microwave radiometers and used in the CCI product as

10.1029/2019JC015327
Journal of Geophysical Research: Oceans compared to the L-Band reflections used by GNSS-R. It is believed that the GNSS-R sensors may be more sensitive to the presence of sea ice (Alonso-Arroyo et al., 2017), and it is for this reason that the 0% concentration threshold for sea ice assignment was used in training the algorithm. The use of the more widely used 15% concentration threshold with these parameters yielded agreement of 97.7% and 95.4% with the CCI product in the Antarctic and Arctic respectively. This effect along with that of snow on passive radiometer measurements (Rostosky et al., 2018) may explain the lower agreement in the marginal ice zone and the increase in false positives (water categorized as sea ice) in the summer months, when higher melt and thin ice may lead to categorization of these points as water by microwave radiometers (operating at higher microwave frequencies than L-band), as used in the CCI algorithm. The highly coherent reflections resulting from leads or meltwater ponds in the melt season may also cause the increase in false positives in the summer time.
In addition to the above discussion, there are some known issues with the TDS-1 data set (covered in detail in Foti et al., 2015 and, including attitude and orbit uncertainties as well as experimental changes in settings due to its use as a technology demonstration platform. These uncertainties related to the September cannot be shown as the majority of the data is in filtered out due to inconsistent delay-Doppler maps (Collection Period 12). Red displays a sea ice assignment and blue an open water assignment. Black shaded area represents the data in the European Space Agency Climate Change Initiative product (Toudal Pedersen et al., 2017) that transitioned between sea ice and water during the month and thus may contain correct ice or water assignments. Poleward of this shaded ring is ice, and equatorward, water.
location of the specular point are believed to be less than 3 km (SSTL, 2017) although this is thought to be an upper limit with likely errors less than a few hundred meters. As this is significantly less than the ESA CCI data resolution, this makes little difference to the comparison. As can be seen in Figure 3, many of the "wrongly" categorized points (ice outside the shaded ring, or water inside) follow the ground tracks of TDS-1, which may suggest that the errors are related to the track as a whole, which may reflect a miscalculation of specular point. These are often in areas which could be easily screened out (e.g., in the Baltic Sea), as performed with the ESA CCI comparison data and other similar products.
It has been found that where the algorithm was trained using Antarctic data, higher agreement is obtained for both hemispheres. This is believed to be due to the geographic differences in the regions of interest, with the Antarctic sea ice zone being further from the pole, and therefore having, on average, data of a lower

10.1029/2019JC015327
Journal of Geophysical Research: Oceans incidence angle and higher antenna gain. This results in higher-quality DDMs, with lower noise, aiding both metrics in the discrimination between sea ice and open ocean. It must also be noted that the higher presence of landmasses in the sea ice margin in the Arctic has an effect on agreement here, with 52% of false positives in the Arctic lying within 50 km of a landmass. This is within the footprint of the sensor where a rough surface (such as land) is concerned. Regardless, the application of the same threshold in both hemispheres and the resulting good discrimination between sea ice and open ocean displays the robustness of this technique.
As can be seen in Figures 4 and 5, there is at times a lack of data within the month. This is partly due to the shared nature of the platform and partly due to the combined threshold approach. The use of the filters described in section 2 reduces the data to approximately 46% of the original observations, and the use of a single observable would allow the use of the entirety of this data set. Therefore, a compromise must be reached between agreement and quantity of data. Where OCOG alone is used to detect sea ice, this agreement falls to from 98.4% to 94.4% in the Antarctic but retains 7% more of the original data.
It should be noted that while previous studies (Alonso-Arroyo et al., 2017;Yan et al., 2017;Yan & Huang, 2016;Zhu et al., 2017) have attained similar or higher accuracies than those seen here, the detection techniques described in those studies are often used across only one season or a single track and have not been applied across the data set as a whole. The method developed and applied here is trained on 20% of the filtered data and shown to produce excellent agreement with CCI across the whole of the remaining data.

Conclusions and Future Work
This study has shown that sea ice detection using GNSS-R is possible with an agreement of 98.4% in the Antarctic and 96.6% in the Arctic through the use of two combined parameters and appropriate thresholds (where compared with state-of-the-art passive microwave data). This is the first time that sea ice detection has been performed using GNSS-R across such a large data set, spanning 33 months of the TDS-1 mission. Maps of sea ice have been created on a monthly basis, including the area over the North Pole that is unavailable for most active remote sensing methods.
This technique shows great potential for future missions, with the extremely low-power, low-cost, and lightweight nature of the sensor lending itself to a polar constellation mission, in a format similar to that of CYGNSS. Such harnessing of these signals of opportunity by a platform and sensor optimized for such measurements would allow high spatial and temporal resolution monitoring of the dynamic sea ice system (average of 4-hr revisit time for CYGNSS).
Improved modeling of the interaction of L-band with the sea ice itself will allow further increases in accuracy and applications of GNSS-R to sea ice classification and snow cover (not considered here). Comparisons of GNSS-R data with results from passive radiometers will facilitate this understanding, especially data from Lband missions such as SMOS (e.g., Kaleschke et al., 2012), which operate at similar frequencies to GNSS-R.