Preliminary analysis of wave retrieval from Chinese Gaofen-3 SAR imagery in the Arctic Ocean

ABSTRACT Synthetic aperture radar (SAR) is a powerful tool for sea surface monitoring. The 2017–2020 mission of the Chinese Gaofen-3 (GF-3) over the Arctic Ocean provided a unique opportunity to investigate the feasibility of wave retrieval in the Arctic Ocean. In this study, several GF-3 SAR images acquired in quad-polarization (QPS) mode were applied for wave retrieval using two theoretical-based algorithms, i.e. the parameterized first-guess spectrum method (PFSM) (with improved tilt modulation) and quad-polarization (Q-P). These images were collocated with wave fields simulated using the third-generation numeric wave model WAVEWATCH-III (WW3) and measurements from the Haiyang-2B (HY-2B) altimeter. A comparison of the HY-2B winds with the wind speeds generated by the geophysical model function (GMF) CMOD5.N showed a 3.36 m/s root mean square error (RMSE) with a correlation (COR) of 0.66. The retrieved significant wave height (SWH) validated against the WW3-simulated results had a 0.50 m RMSE using PFSM and 0.91 m RMSE with the Q-P algorithm. These results show that PFSM is suitable for wave retrieval from GF-3 SAR images. We also used PFSM to investigate the characteristics of the wave spectrum in the presence of sea ice. Although the analysis concludes that GF-3 SAR has the capability for wave monitoring in Arctic Ocean due to the high spatial resolution of SAR-derived wave spectra, an optimal wave retrieval algorithm needs to be developed for improving the retrieval accuracy.


Introduction
A gradual reduction of perennial sea ice in the Arctic Ocean has been occurring since the late 1990s due to global climate change (Shimada et al., 2006). The reduction in sea ice is conducive to an increase in the quantity of waves (Thomson & Rogers, 2014) caused by an increase in the fetch (Dobson et al., 1989). In past decades, the major method for dynamic processes monitoring in the Arctic Ocean has been ship-based observation, which requires a significant amount of manpower and other resources. The third-generation numerical wave model WAVEWATCH-III (WW3; Tolman, 1991) could be used for wave simulation in Arctic Ocean (Shao et al., 2022;Wojtysiak et al., 2018). With the development of remote sensing, satellites operating at a microwave frequency can be used to detect sea surface phenomena such as wind and waves (Gourrion et al., 2002). Currently, the operational products derived from a scatterometer, e.g. Advanced Scatterometer (ASCAT; Chou et al., 2013), and an altimeter, e.g. Haiyang-2B (HY-2B; H. F. Zhang et al., 2015), are the primary sources of remote-sensing information for oceanography research. However, the spatial resolution of altimetermeasured waves (10 km) is relatively coarse.
Two categories of SAR wave retrieval include theoretical-based algorithms (Alpers et al., 1981) and empirical models. Three types of modulations, e.g. tilt modulation (Lyzenga, 1986), hydrodynamic modulation (Feindt et al., 1986), and velocity bunching (Alpers et al., 1981) have been well studied in the last decades. However, these modulations have not yet been applied at the polar regions because the tilt and hydrodynamic modulation (MTF) needs to be adapted to consider the non-Bragg backscattering caused by sea ice (Schulz-Stellenfleth & Lehner, 2002). The first co-polarization wave retrieval algorithm (VV and HH) is called the Max-Planck Institute Algorithm (MPI; Hasselmann & Hasselmann, 1991), which uses simulated data from numeric models as the inputs and solves three modulations by minimizing a cost function. Continuously, a semi-parametric retrieval algorithm (SPRA) has been proposed (Mastenbroek & Valk, 2000), in which winds from a scatterometer are used to obtain a first-guess wind-sea spectrum using an empirical wave model. A partition rescaling and shift algorithm (PARSA; X. M. Li et al., 2010;Schulz-Stellenfleth et al., 2005) has also been developed based on the SAR look-cross spectrum. The advantages of MPI and SPRA are included in the parameterized firstguess spectrum method (PFSM; Sun & Guan, 2006;Zhu et al., 2018). For instance, a wind-sea and swell SAR intensity spectrum is empirically separated. He et al. (2006) proposed an algorithm for quadpolarization (Q-P) using information from both coand cross-polarization channels. Validation against moored buoys has proved the applicability of wave retrieval for quad-polarization RADARSAT-2 (R-2; B. Zhang et al., 2010). In the meantime, empirical functions have been proposed for significant wave height (SWH) retrieval of co-polarized (Schulz-Stellenfleth et al., 2007) and cross-polarized SAR (Shao et al., 2020a) waves, e.g. CWAVE (X. M. Li et al., 2010;Stopa & Mouche, 2017) and XWAVE (Bruck & Lehner, 2013;Shao et al., 2020b). Although these empirical models are usually applied without resolving the MTF of non-linear velocity bunching, the specific formulations are various defined (Grieco et al., 2016;Zhu et al., 2019).
Empirical wave retrieval algorithms have been adopted for specific SAR sensors. However, the effect of sea ice has not been considered until now and the empirical algorithms could not be applied for SAR wave retrieval in the background of sea ices. Moreover, the wave spectrum cannot be retrieved using the empirical algorithms, which is essential for research on wave dynamics. Therefore, in this work, we focused on wave retrieval from quad-polarization GF-3 SAR images of the Arctic Ocean using PFSM and Q-P algorithms. As demonstrated in Schulz-Stellenfleth and Lehner (2002), tilt modulation is stronger within sea ice. We updated this parameter in PFSM using the incidence angle so as to include the ice-induced effect on tilt MTF. We investigated the performance the algorithms by comparing retrieval results with WW3-simulated wave fields and measurements from the HY-2B altimeter. The main purpose of our work was to not only confirm the feasibility of wave monitoring in Arctic Ocean using the GF-3 SAR image but also to select a suitable algorithm for wave retrieval. We also performed error analysis for various sea states to provide a perspective for improving the adaptability of the algorithm. This paper is organized in six sections. The collocated dataset is described in Section 2, including the GF-3 SAR images, the measurements from the HY-2B satellite, the European Centre for Medium-Range Weather Forecasts (ECMWF) reanalysis data, and the wave fields simulated using the WW3 model. The methodology of the two types of wave retrieval algorithms is presented in Section 3. The validation of winds and waves against auxiliary data is presented in Section 4. Our conclusions are summarized in Section 5.

Dataset
Twenty GF-3 SAR images in QPS mode of the Arctic Ocean from the last two years were acquired from the National Satellite Ocean Application Service (NSOAS). Specifically, information on sea ice is included in most those images which proceeded to be a Level-1A (L-1A) product. The geographic locations of the collected images are shown in Figure 1. The calibration method is applied using the following equation : in which σ 0 is the polarimetric SAR Normalized Radar-Cross Section (NRCS); DN is the intensity; and variables M and N are external constants accompanied with raw data. It has been found that NSOAS has recently released a document so as to remove the uncertainties in GF-3 calibration process since 2020 and the calibration document is available on the server vai ftp://osddsftp.nsoas.org.cn. Figure 1 shows an example of GF-3 SAR images at the VV-polarization channel (a), HHpolarization channel (b), HV-polarization channel, (c) and VH-polarization channel (d), which was taken at 08:25 UTC on 25 May 2017. In this example, the sea ice coverage is identified from the image. Numeric wave models, e.g. Simulating Waves Nearshore (SWAN) and WW3 models, are useful for global wave distribution research (Zheng et al., 2016), even in tropical cyclones (Hu et al., 2020). Sea ice has been considered in the WW3 model (version 6.07; WAVEWATCH III Development Group, 2019). Studies have shown that the WW3 model has the capability to simulate waves in the Arctic Ocean (Pingree-Shippee et al., 2016;Shao et al., 2022). In order to study the distribution of the SAR-derived spectra across the sea ice, the waves were simulated using the WW3 model. ECMWF reanalysis winds on a 0.25° grid with 1-hour intervals were used as the forcing fields. These winds were used as the forcing field in previous studies concerning wave simulation. The simulated region is 40°N, 85°N latitude 0°E, 360°E longitude. The water depths are derived from the General Bathymetric Chart of the Oceans (GEBCO) bathymetric data with ~1-km horizontal resolution). The sea ice concentration and ice thickness data from the Copernicus Marine Environment Monitoring Service (CMEMS) are used as forcing data. Maps of sea ice concentration and thickness in September 2017 are shown in Figure 2(a and b), respectively. The black spot in these figures represents the spatial coverage of the GF-3 SAR image in Figure 1(b). The WW3 model provides eight switches for ice-wave interaction and wave simulation in the presence of sea ice. A recent study (Shao et al., 2022) concluded that switch IC4_M1 (J. Li et al., 2019;Rogers et al., 2016;Wadhams et al., 1988) is most suitable for wave simulation. Figure 3(a and b) show the ECWMF and WW3-simulated SWH maps at 08:00 UTC on 25 May 2017. The black rectangle in these figures represents the geographic location of the image in Figure 1.
Wave observations are reliably available for most oceanography research. However, it is difficult to obtain real-time wave data in the Arctic Ocean. Remotely-sensed measurements (Campbell et al., 1987), e.g. winds from scatterometer and waves from   Figure 1 altimeter, are one source for polar sea studies. The civilian marine dynamics satellite HY-2B (Xu et al., 2016) is designed for operational sea surface winds and waves monitoring. It has been proven that the measurements from HY-2B altimeters (Jia et al., 2020) and scatterometer (Wang et al., 2020) have good quality from global wind and wave monitoring. The maximum latitude of the available data from the HY-2B satellite can extend as far as 80°N. Therefore, in this work, available HY-2B measured wind speeds are used to validate the SAR-derived wind speeds, and HY-2B measured SWH data are used to validate the WW3-simulated SWH due to there are few footprints of HY-2B altimeter passing the spatial coverage of collected GF-3 SAR images. Specifically, the matchups are selected following the criterion: spatial distance between WW3 grids and HY-2B footprints is within 3 km and the time difference is less than 30-min. Figure 4 shows the SWH maps from the HY-2B altimeter following the satellite's orbit in May of 2017. The validation of 0.1° gridded WW3-simulated SWH against 0.25° gridded ECMWF and HY-2B SWH is presented in Figure 5(a and b). The analysis yielded 0.42 m RMSE of SWH with 0.84 correlation (COR) as compared with ECMWF data and 0.37 m RMSE of SWH with 0.82 COR as compared with HY-2B data. Although the WW3 model has a relatively poor performance at low states (SWH < 0.3 m), we think that the simulations from the WW3 model are reliable enough for this work.

Methodology of wave retrieval algorithm
In this section, two types of theory-based wave retrieval algorithms (PFSM for co-polarized SAR and Q-P for quad-polarized SAR) are discussed in detail. Because winds are the prior information in the process of SAR wave retrieval, the GMF for C-band SAR wind retrieval is also briefly introduced.

SAR wind retrieval algorithm
In general, CMOD family is popular for wind retrieval from C-band SAR image. At present, the GMF CMOD5.N has good performance for GF-3 SAR wind retrieval with an accuracy of 2 m s −1 RMSE as validated against the measurements from scatterometer , which takes the basic function of form (Hersbach, 2010): in which σ 0 is the linear NRCS, φ is the wind direction related to the radar look direction, and matrices B 0 to B 2 are empirical functions of sea surface wind speed at 10 m height U 10 and radar incidence angle. After employing 0.25° gridded ECMWF wind directions, the wind speeds are inverted from SAR NRCS images.

Algorithm PFSM
The wave retrieval scheme MPI based on SAR mapping mechanism was initially proposed in (Hasselmann & Hasselmann, 1991) by solving nonlinear mapping mechanism velocity bunching. Using the wave spectrum simulated from a numeric model, the principle of MPI is used to minimize the cost function J in Equation. (3): where k is the wavenumber, I k is the SAR-derived wave spectrum at wave number k, Ī k is prior information on wave spectrum from WAM model, Ē k is the mapped SAR spectrum using the three MTFs, E k is the calculated SAR intensity spectrum by twodimensional Fast Fourier Transform (FFT-2), and μ is the weight coefficient. The small positive number B is assumed to be 0.001 in order to ensure convergence of variation in the solution. The SAR MTFs T S k (Schulz-Stellenfleth & Lehner, 2002) includes three portions: the tilt modulation T t k , hydrodynamic modulation T h k and velocity bunching T v k : in which k x and k y are the wave number at azimuth and range directions, V is the satellite flight velocity, R is the satellite slant range at each pixel, θ is the radar incidence angle, and ω is the circular wave frequency. In order to improve the efficiency of wave retrieval algorithm, algorithm SPRA was constructed (Mastenbroek & Valk, 2000), in which the "firstguess" wave spectrum is produced by an empirical wind-sea function, such as JONSWAP and Elfouhaily model (E-spectrum; Elfouhaily et al., 1997). The windsea retrieval follows the MPI scheme using a SAR intensity spectrum. However, the systematic error induced by the difference between the wind-sea and total wave spectrum is inherent in the swell retrieval process. Algorithm PFSM was developed recently to reduce the error in algorithm PFSM. It divides the SAR intensity spectrum into two portions: non-linear wind-sea mapping and linear swell mapping. The separation threshold k s is calculated using prior information on wind speed at 10 m height U 10 and other imaging parameters, as described by Equation. (7): 10 cos 2 ϕ sin 2 ϕsin 2 θ þ cos 2 ϕ ð Þ � � 0:33 (8) in which g is gravitational acceleration at 9.8 m/s 2 , and ϕ is the angle of the wave propagation direction relative to the radar look direction. It has been proven that algorithm PFSM works for C-band and X-band SAR (Shao et al., 2015), yielding about 0.6 m root mean square error (RMSE) of SWH for homogeneous SAR images with look-like wave pattern. In this study, wind speeds were inverted from VV-polarized SAR images using a GMF CMOD5.N (Hersbach, 2010) and ECMWF wind directions. SAR wind-sea I and swell II were retrieved from corresponding portions. I. When using the E-spectrum model, three variables (wind speed, dominant wave propagation velocity, and direction) are necessary. Therefore, wind-sea retrieval from the portion of SAR spectrum at wave-number k≥ k s requires 5 steps: (1) calculating velocity c p and direction ϕ at peaks of SAR intensity; (2) parameterizing dominate wave propagation velocity from 0.4c p to 1.6c p at an interval of 0.2 m/s and parameterizing dominant wave propagation direction from 0.8ϕ to 1.2ϕ at an interval of 2°; (3) producing wave spectra corresponding to each parameterized group using an E-spectrum function and SAR mapping spectra simulated by MTFs; (4) obtaining the best "first-guess" wave spectrum by minimizing the series of simulated SAR spectra and SAR wind-sea portion; and (5) retrieving wind-sea l by solving Equation.
(2), in which "first-guess" wave spectrum and SAR wind-sea portion are the inputs.
II. For the SAR swell portion, velocity bunching becomes weak at wave-number k< k s , indicating that the MTF of velocity bunching can be ignored. Under this circumstance, the SAR mapping mechanism includes only tilt modulation and hydrodynamic modulation. The swell spectrum is obtained by directly solving the linear portion of the SAR spectrum.
The SAR-derived wave spectrum is a composite of the wind-sea and swell retrieval. The SWH H s is calculated by integrating a wave spectrum W k at wavenumber k: ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ð W k dk s (9)

Algorithm Q-P
The abundant information available from the quadpolarization channel should be extremely useful for wave retrieval from polarimetric SAR (Schuler et al., 2004). The empirical function CWAVE has been recently adopted for cross-polarized GF-3 SAR (Shao et al., 2020a). The basic principle of algorithm Q-P is that the wave slope spectrum S rms is correlated with a linearly polarized image σ p in the SAR azimuth and range direction, in which the winds and complex MTFs of mapping mechanism are skipped. SWH H s can be quantitatively calculated from the rms wave elevation S rms calculated from wave slope, as described in the following equation (He et al., 2006): As proposed in B. Zhang et al. (2010), the rms wave elevation S rms is empirically calculated by: S rms ¼ ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffiffi @� @x � sinϕ in which ϕ is the dominant wave propagation direction and @�=@x and @�=@y satisfies following equations (He et al., 2004): Based on the SAR images at the HH-, VV-, and HV-polarization channels, the linearly polarized image σ 0 pp is estimated when the polarization orientation angle η induced by target rotation about radar look beam is set to 45°, which changes the target scattering characteristics:  in which σ 0 VV is VV-polarized NRCS, σ 0 HH is HHpolarized NRCS, σ 0 HV is HV-polarized NRCS, σ 0 HHVV is correlation coefficient between the HH-and VVpolarizations, and <½� represents the real parts of the indicated quantities. Using the above equations, algorithm Q-P is applied for GF-3 SAR. However, implementation of this algorithm still needs to be further studies when applied to data from the Arctic Ocean.

Results
Because wind speed is necessary for the PFSM algorithm, geophysical model function (GMF) CMOD5.N (Hersbach, 2010) was used for inverting winds from VV-polarized SAR images in 2019-2020. Figure 6 shows the SAR-derived wind map corresponding to the image in Figure 1(a). Note that the float ices are manually identified where brightness on the SAR image is due to reflection rather than Bragg scattering and the inverted winds are masked at such regions. The comparison of available match-ups between SARderived winds and measurements from the HY-2B altimeter and scatterometer shows a 3.36 m/s RMSE with a 0.66 COR (Figure 7). SAR-derived wind speeds are distorted in the presence of sea ice due to non-Bragg backscattering, resulting in a larger error than the typical accuracy of 2 m/s of SAR-derived wind speeds for GF-3 SAR ).
An image of the VV-polarized sub-scene in which the wave patterns are not apparent is presented in Figure 8(a) so as to examine the applicability of two wave retrieval algorithms. The corresponding VVpolarized SAR intensity spectrum and the quadpolarized SAR slope spectrum are shown in Figure 8 (b and c), respectively. The symmetrical structure of the quad-polarized SAR slope spectrum is not apparent, in contrast with the VV-polarized SAR intensity spectrum. The wave direction could not be directly measured by SAR, that is, the ambiguity is included in the signal at the radar look direction. Therefore, the wave propagation direction from WW3 model is employed in order to remove that ambiguity of inversion result. The two-dimensional SAR-derived wave spectrum of the sub-scene by PFSM is exhibited in Figure 8(d). In presence of sea ices, small change in the slope angles perpendicular to the radar incidence plane synchronously induces change in the polarization orientation angle due to the volume scattering at HV-polarization channel, however, the polarization orientation angle is set as 45° in Q-P. The asymmetric pattern from quad-polarized SAR slope spectrum is probably caused by the missing information on tilt and polarization orientation modulations by sea ices. The SAR-derived SWHs from the sub-scenes (extracted from GF-3 SAR images using two algorithms) are validated against the SWHs from the HY-2B altimeter. We found that the RMSE of SWH is 0.5 m using PFSM with a 0.77 COR (Figure 9(a)), which is better than the 0.9 m with a 0.42 COR obtained using Q-P (Figure 9(b)). The overestimation is apparently observed at low state SWH < 1.2 m. Because of the non-Bragg backscattering caused by sea ice and wave breaking mixed with the Bragg backscattering, the SAR imaging mechanism is truly complicated. Because the development of Q-P excluded the modulation of non-Bragg backscattering, the Q-P has bad performance for SWH retrieval in Arctic Ocean.
We selected several spots where floating sea ice exists for further study of the distribution of wave spectra ( Figure 10). Specifically, there is floating sea ice at selected spots, e.g. A and B, C and D. The SARderived wave spectra using PFSM and WW3simulated wave spectra at spots A to D are presented in Figure 11. Generally, the patterns of SAR-derived wave spectra are consistent with WW3-simulated wave spectra, especially at spot C, which is away from the floating sea ice. In this sense, SAR is a promising technique for wave analysis in the Arctic Ocean due to the high spatial resolution of SAR-derived wave spectra. However, our work only gives preliminary analysis of wave spectrum retrieval from GF-3 SAR image, more images are necessary to make a further study.

Summary
Due to rapid climate change, the sea fetch in the Arctic Ocean has grown as ice has melted over the past decades (Duan et al., 2019). Till now, numeric wave models (i.e. WW3 and SWAN) have been widely employed for wave monitoring in Arctic Ocean. However, the output from the numerical model relies on accurate winds and considering the non-linear icewave interactions is still a challenge in the modeling community (Liu et al., 2020). The Chinese satellite GF-3 carrying C-band SAR has been operating over global seas, including the Arctic Ocean, since August 2016. The purpose of this study was to investigate the applicability of wave monitoring using GF-3 SAR in the Arctic Ocean, with a particular focus on wave propagation across sea ice. Twenty images from the period 2017-2020 were available for our work.
Those images were collocated with HY-2B measurements and wave spectra simulated from the WW3 model. The SAR-derived wind speeds using GMF CMOD5.N were validated against measurements from the HY-2B scatterometer, yielding a 3.36 m/s RMSE with a 0.66 COR. Two kinds of wave retrieval algorithms, e.g. PFSM and Q-P, were implemented. The wind-sea wave spectrum function JONSWAP in PFSM was replaced by a full wavenumber spectrum function, denoted as E-spectrum (Elfouhaily et al., 1997). The MTF of tilt modulation considering the term of incidence angle proposed in Schulz-Stellenfleth and Lehner (2002) was initially adopted for GF-3 SAR wave retrieval in the Arctic Ocean. The retrieval SWHs were compared with WW3-simulated results. The RMSE of SWH using PFSM was 0.5 m with a 0.77 COR, which is better than the 0.9 m with a 0.42 COR obtained using Q-P. It is necessary to figure out that non-Bragg backscattering is ignored in the principle basis of Q-P, which is why Q-P is not suitable for SWH retrieval in the Arctic Ocean. In this sense, the accuracy of wave retrieval could be improved when including ice-induced effects in the methodology of Q-P, i.e. tilt and polarization orientation modulations.
SAR retrieval using PFSM and WW3-simulated wave spectra at several locations was selected to study the wave spectrum characteristics. We found that a comparison of SAR-derived wave spectra with the model simulations illustrated the SARderived wave spectra have more details than WW3simulated wave spectra as encountering the floating sea ices, due to the high spatial resolution of SARderived wave spectra. Therefore, the SAR observations are useful for wave monitoring around small float ices in near real-time. It is concluded that GF-3 SAR has the capability for wave monitoring in Figure 11. The SAR-derived wave spectra using PFSM and WW3-simulated wave spectra: (a) spot A; (b) spot B; (c) spot C; and (d) spot D.
Arctic Ocean, although an optimal wave retrieval algorithm is necessary for improving the retrieval accuracy. In the near future, a collaboration between GF-3 SAR, HY-2C/2D, and the China-France Oceanography Satellite (CFOSAT) in Arctic Ocean will be proposed to further investigate the characteristics of waves, especially with respect to air-sea and wave-ice interactions.

Acknowledgments
We appreciate the National Ocean Satellite Application Center (NSOAS) providing the Gaofen-3 (GF-3) SAR images and Haiyang-2B (HY-2B) product and the National Center for Environmental Prediction (NCEP) of National Oceanic and Atmospheric Administration (NOAA) providing the source code for the WAVEWATCH-III (WW3) model free of charge. The European Centre for Medium-Range Weather Forecasts (ECMWF) wind data were accessed via http://www. ecmwf.int. General Bathymetry Chart of the Oceans (GEBCO) data were downloaded via: ftp.edcftp.cr.usgs.gov. Additionally, we truly thank the Copernicus Marine Environment Monitoring Service (CMEMS) for offering sea ice concentration and thickness data via http://marine.copernicus.eu.

Disclosure statement
No potential conflict of interest was reported by the author(s).

Data availability statement
Due to the nature of this research, participants of this study did not agree for their data to be shared publicly, so supporting data are not available.