Blending Sea Surface Winds from the HY-2 Satellite Scatterometers Based on a 2D-Var Method

Abstract

The launch of the Haiyang-2 (HY-2) satellite constellation fosters the quick acquisition of global sea surface vector winds from the perspective of remote sensing. This study intends to develop a six-hourly mesoscale analysis of sea surface winds based on the microwave scatterometers onboard the HY-2 satellite series, with the objective of meeting the considerable demand for accurate and gap-free ocean wind forcing products. First, the accuracy of HY-2 scatterometers (HSCATs) in measuring wind is evaluated. In particular, the standard deviation (SD) errors of HSCATs data are assessed using the collocated buoy measurements with different temporal windows in order to account for the temporal representativeness errors in the blending analysis. Afterwards, a two-dimensional variational (2D-Var) method is implemented to blend the HSCATs measured winds and the European Centre for Medium-Range Weather Forecasts (ECMWF) reanalysis winds over the global ocean surface. This approach is different from existing methods in that it takes both wind error and background error correlation into account. The results show that the blended wind product is of a promising quality compared with independent wind references. Interestingly, the blended winds derived from 2D-Var in combination with an empirical background error correlation show smaller SD errors than those using a Gaussian error correlation function. Overall, the blended wind product should be valuable for forcing global ocean models or describing air-sea interaction processes on a scale close to the scatterometer measurements.

1. Introduction

Sea surface vector wind is an essential climate variable that influences many physical processes over the air-sea interface [1,2,3]. There is a growing demand for accurate sea surface wind products with higher resolution and, notably, more frequent coverage, for tracking transient features and resolving air-sea interaction and its diurnal cycle, etc. For spaceborne scatterometers, which are the most relevant sensors in sea surface wind observation, the World Meteorological Organization (WMO) specifies a 25 km spatial resolution and six-hourly temporal sampling as the operational goals [4]. However, a single scatterometer can only acquire wind data with a temporal resolution of 12 h, covering less than 90% of the global sea surface per day. Consequently, numerous methodologies have been developed to blend the satellite wind in order to enhance the temporal coverage of ocean vector winds without reducing the short-scale true variance resolved by the remote sensing technique.

In the early days, the number of operational scatterometers was too few for developing a gap-free and high-temporal resolution (e.g., six-hourly) wind product only based on the scatterometer-derived wind data. Consequently, the scatterometer data are usually blended with ancillary wind data, such as buoy measurements, radiometer wind speeds, and numerical weather prediction (NWP) winds. For instance, the Cross-Calibrated Multi-Platform (CCMP) ocean vector wind analysis used a two-dimensional variational analysis (2D-Var) method to combine satellite scatterometer and radiometer winds, in-situ buoy observations, and NWP background field to fill the sampling gaps [5,6,7]. A similar product was described in ref. [8], which used the kriging technique instead of the 2D-Var. Moreover, there were two other blended wind products with lower temporal resolution (daily), which were developed in the context of the Objectively Analyzed air-sea Fluxes (OAFlux) project and the merging of Advanced Scatterometer (ASCAT) winds with the NWP reanalysis [9,10], respectively. These wind products extend more than 20 years and have made significant contributions to the research on mesoscale air-sea interaction.

More recently, international collaboration is highlighted with the objective of developing a global virtual constellation of wind scatterometers for high-resolution ocean forcing. Although the number of operational scatterometers is still insufficient to fully meet the above WMO requirements [11], the attempt in ref. [12] shows that the spatial scale of global winds can be effectively improved by blending the ASCATs data with the European Centre for Medium-range Weather Forecasts (ECMWF) reanalysis. Later in 2021, more scatterometer systems measuring at the different local times have been available. In particular, the Haiyang-2 (HY-2) constellation, which consists of three satellites, has been successfully developed by China. Notably, HY-2B is operated in a sun-synchronous orbit, while HY-2C and -2D are operated in drifting orbits, such that the three HY-2 satellite scatterometers (HSCATs) are able to provide unprecedented short revisit intervals over the global sea surface [13,14,15,16].

This study aims at developing a six-hourly mesoscale wind product over the global sea surface only based on the HSCATs wind data and the ECMWF reanalysis winds. Section 2 describes the datasets and the methods used in this paper. Section 3 presents a brief overview of the wind quality for the three HSCATs. In Section 4, the blended wind product is verified comprehensively using buoy data as a reference. Finally, conclusions are summarized in Section 5.

2. Data and Methods

2.1. Data

The main parameters of HY-2B, -2C and -2D satellites are summarized in

Table 1. Main parameters of the HY-2 satellites.

Satellite	Altitude (km)	Inclination (°)	LTDN *	Repeat Cycle
				The Early Stage	The End Stage
HY-2B	971	99.3	6:00 am	14 days	168 days
HY-2C	957	66	Drifting	10 days	400 days
HY-2D	957	66	Drifting	10 days	400 days

* LTDN: Local time of the descending node.

. These three satellites were launched in October 2018, September 2020, and May 2021, respectively. Consequently, the overlapped period is not adequate for developing a long-term time series of global sea surface vector wind products. The main objective of this paper is to verify the proposed 2D-Var method for blending HY-2 satellite winds, and, thus, only six months (December 2021–May 2022) of the HSCATs 25-km Level 2 (L2) wind data are used in the study. Such data are provided by the National Satellite Ocean Application Service (NSOAS) through the Ocean Satellite Data Distribution System (OSDDS). Note that, OSDDS provides two types of HSCAT winds, namely “pwp” and “dps” separately. This paper uses the HSCAT “pwp” data.

Figure 1 illustrates the six-hourly coverage of the three HSCATs on 7 May 2022, at UTC 00 h, 06 h, 12 h, and 18 h, respectively. Each panel shows the observations that are less than 3 h distance in time to the above-mentioned hours. It is clear that the three HSCATs are able to observe more than 70% of sea surface area in six hours, e.g., the coverage ratio over the sea surface is about 70%, 75.2%, 74.6%, and 75.6% for Figure 1a–d, respectively.

The fifth generation of ECMWF atmospheric reanalysis (ERA5) winds is used as the background data when blending the HSCATs data. This study uses the ERA5 neutral wind vectors at a 10 m height with approximately 0.25° × 0.25° spatial resolution and six-hourly temporal resolution. Before the data analysis, both ERA5 background and HSCATs observations are interpolated bi-linearly to a regular latitude/longitude grid of 0.25° resolution. Consequently, one can collocate the HSCATs data with the ERA5 background winds easily by setting a temporal criterion of 3 h distance.

Finally, the moored buoy winds are used for verification purposes. Such data are disseminated through the Global Telecommunications System (GTS) stream and, subsequently, quality controlled and archived at the ECMWF Meteorological Archive and Retrieval System (MARS). Note that the buoy winds are reported hourly and segregated into 1 m/s speed bins and 10° direction bins. The measured wind speed at a given anemometer height is converted to 10 m equivalent neutral winds using the Liu-Katsaros-Businger (LKB) model to make them comparable with the HSCATs and the ERA5 winds [17]. Eventually, the buoy measurements that are 5 min distance in time and 0.125° distance in space from the ERA5 grids are used to verify both the HSCATs and the ERA5 winds.

2.2. 2D-Var Method

Following the basic assumptions in meteorological variational data assimilation, a 2D-Var method is applied to blend the HSCATs observations with the ERA5 background wind field. Conceptually, given a set of scatterometer observations $x_o^k$ (x denotes a vector, k denotes the kth satellite scatterometer) and a vector of background winds $x_b$, the analyzed wind vector x is obtained by minimizing the following objective function:

$$\min \left[J\left(x_o^k,x,x_b\right)\right]=J_b\left(x_b,x\right)+J_o\left(x_o^k,x\right)=\frac{1}{2}{\left(x-x_b\right)}^T{R}_b^{-1}\left(x-x_b\right)+\frac{1}{2}{\left(x-x_o^k\right)}^T{R}_o^{-1}\left(x-x_o^k\right)$$

(1)

where J_b and J_o denote the background (B) and the observational (O) terms, respectively. R_b is the background error covariance matrix, and R_o is the observation error covariance matrix. For computational convenience, the state vector x is replaced by the analysis increment $\Delta x=x-x_b$. In practice, an alternative formulation of the observational term, which was developed for the ambiguity removal of SeaWinds on QuikSCAT [18], is used in this study. That is, the observation term J_o is written as:

$$J_o={\displaystyle \sum _{i,j}^{N_1,N_2}{\left\{{\displaystyle \sum _k^{M_{ij}}{\left[\frac{{\left(\Delta u_{ij}-\Delta v_{ij,k}^{(o)}\right)}^2}{{\sigma }_u^2}+\frac{{\left(\Delta u_{ij}-\Delta v_{ij,k}^{(o)}\right)}^2}{{\sigma }_v^2}-2\ln p_k\right]}^{-\lambda }}\right\}}^{-1/\lambda }}$$

(2)

with (i,j) the indices of the 2D-Var grid cell, N₁ and N₂ the number of 2D-Var grids along the zonal and the meridional direction, M_ij the number of HSCATs observations at point (i,j). Moreover, Δu_ij and Δv_ij stand for the analysis wind component increments at grid (i,j), $\Delta u_{ij,k}^{(o)}$ and $\Delta v_{ij,k}^{(o)}$ are the difference between the kth observation and the background wind components, ${\sigma }_u$ and ${\sigma }_v$ denote the errors in the observed wind fields.

Note that, Equation (2) is similar to the cost function used in the scatterometer ambiguity removal scheme, such that λ = 4 is used without any adaption [19,20]. Furthermore, note that the index k indicates the satellite number rather than the ambiguous solution. Consequently, the a priori probability of the kth observation being the “true” wind is the inverse of the number of observations within the grid (i,j), i.e., p_k = 1/M_ij, following the assumption that the three HSCATs contribute equally to the cost function. As shown in Figure 1, there is usually only one HSCAT observation in a grid cell (except for the grids at high latitudes), such that the term 2lnp_k is usually equal to zero. The background wind term of Equation (1) is written as:

$$J_b(\Delta x)={(\Delta x)}^{\mathrm{T}}{B}_{u,v}^{-1}(\Delta x)$$

(3)

where $B_{u,v}^{-1}$ denotes the background error covariance matrix. In the numerical implementation, Equation (3) is transformed to the frequency domain and is expressed in terms of velocity potential (χ) and stream function (ψ) increments, and the background error covariance is factorized into standard deviation (SD) errors and error correlations, which both can be expressed as diagonal matrices.

As described in ref. [21], the 2D-Var analysis field can be optimized by adjusting the background error correlation function, and the background or observation error variances. Consequently, on the one hand, the SD errors of HSCATs winds are estimated for different collocation windows and then used to construct the observational error covariance matrix. On the other hand, two different background error correlation functions (as shown in Figure 2), namely the Gaussian structure function (GSF) and the empirical background error correlation (EBEC), are evaluated in order to figure out an appropriate 2D-Var implementation. The exact formulas of GSF and EBEC can be found in refs. [21,22], such that they are not reproduced here in order to keep the presentation concise. Note that the GSF is with different length scales for the Tropics (solid thin curve, L_χ = L_ψ = 212 km) and the extratropics (solid bold curve, L_χ = L_ψ = 424 km), while the EBEC is with the same length scales (L_χ = 494 km, L_ψ = 552 km) over the global region.

3. Data Verification

Before blending the HSCATs winds with the ERA5 background data based on the 2D-Var method, the differences among HSCATs, ERA5, and buoy winds are evaluated in order to verify the consistency of the used data sources. In particular, the SD difference between HSCATs and buoy winds is assessed for different collocation windows, with the objective of accounting for the temporal representativeness errors in the blending analysis.

Figure 3 shows the contour plots of HSCATs versus ERA5 winds. Here, the time difference between HSCATs and ERA5 is less than 30 min, and only the data accepted by the quality control scheme are used in the analysis. The upper left corner of each panel illustrates the correlation coefficient (cc), the bias, and the SD values of HSCATs with respect to (with regards to) ERA5 winds. Compared with the ERA5 background field, the bias of HSCATs winds is negligible, and the SD differences of wind speed and direction are generally less than 0.95 m/s and 10°, indicating that the observational data are of excellent consistency with the background fields. In addition, winds derived from the three HSCATs are with similar statistical scores with regards to. ERA5, such that it is not necessary to perform a bias correction before the 2D-Var analysis.

Figure 4 verifies the consistency between the ERA5 wind and the collocated buoy winds. Note that the moored buoy data have been already assimilated into the ERA5 model, such that the verification is indeed not independent in terms of quality assessment. Nonetheless, statistical scores of the ERA5 winds with respect to buoy measurements are shown in the upper-right corner of each panel, which indicates that the two wind data sets are consistent with each other.

Eventually, the collocated HSCAT-buoy data sets are segregated into six categories according to the time difference between the two data sources, i.e., 0–0.5 h, 0.5–1.0 h, …, and 2.5–3.0 h. The statistical scores between HSCATs and buoy winds are evaluated for each category and summarized in

Table 2. Biases of the HSCATs winds with regards to. buoy measurements for different time windows.

Scatterometer		Wind Speed (m/s)			Wind Direction (°)
		HSCAT-B	HSCAT-C	HSCAT-D	HSCAT-B	HSCAT-C	HSCAT-D
Time difference	0.0–0.5 h	−0.28	−0.19	−0.18	0.5	1.6	0.3
	0.5–1.0 h	−0.26	−0.19	−0.16	0.5	0.9	1.3
	1.0–1.5 h	−0.28	−0.23	−0.15	0.8	1.5	1.2
	1.5–2.0 h	−0.31	−0.24	−0.20	0.3	0.1	0.2
	2.0–2.5 h	−0.30	−0.20	−0.18	1.0	0.9	0.8
	2.5–3.0 h	−0.25	−0.25	−0.18	1.3	1.0	1.1

and

Table 3. SD differences of the HSCATs winds with regards to. buoy measurements for different time windows.

Scatterometer		Wind Speed (m/s)			Wind Direction (°)
		HSCAT-B	HSCAT-C	HSCAT-D	HSCAT-B	HSCAT-C	HSCAT-D
Time difference	0.0–0.5 h	1.11	1.05	1.20	17.7	17.9	18.1
	0.5–1.0 h	1.20	1.06	1.21	16.6	18.0	17.4
	1.0–1.5 h	1.18	1.16	1.13	17.1	16.6	17.0
	1.5–2.0 h	1.24	1.28	1.38	17.0	18.1	17.6
	2.0–2.5 h	1.40	1.45	1.42	18.7	19.9	19.3
	2.5–3.0 h	1.39	1.53	1.55	20.3	21.1	22.1
		u component (m/s)			v component (m/s)
		HSCAT-B	HSCAT-C	HSCAT-D	HSCAT-B	HSCAT-C	HSCAT-D
Time difference	0.0–0.5 h	1.65	1.69	1.75	1.85	2.22	2.19
	0.5–1.0 h	1.70	1.65	1.71	1.76	2.05	2.22
	1.0–1.5 h	1.62	1.64	1.60	2.01	1.96	1.77
	1.5–2.0 h	1.75	1.82	1.67	1.77	2.10	1.89
	2.0–2.5 h	1.95	1.93	1.96	2.08	2.21	2.49
	2.5–3.0 h	2.00	2.06	2.23	2.28	2.33	2.45

. Note that, the bias of HSCATs winds (with regards to. buoy) does not show a clear dependency on the temporal window (see Table 2), while the SD difference does increase as the time difference increases (see Table 3). Consequently, it is not necessary to do a time-dependent bias correction for the HSCATs observations. However, the temporal representativeness errors, which are induced by the time difference between O and B, should be considered in the 2D-Var cost function. In practice, the ratio of observation to background error, rather than their absolute values, impacts the 2D-Var cost function. Following the results of Figure 4 and Table 3, a time-dependent observation-to-background error ratio is used in the following study, as shown in

Table 4. Time-dependent observation-to-background error ratio.

Scatterometer		u Component			v Component
		HSCAT-B	HSCAT-C	HSCAT-D	HSCAT-B	HSCAT-C	HSCAT-D
Time difference	0.0–0.5 h	1.00	1.02	1.06	0.97	1.17	1.15
	0.5–1.0 h	1.03	1.00	1.03	0.93	1.08	1.17
	1.0–1.5 h	0.98	0.99	0.97	1.06	1.03	0.93
	1.5–2.0 h	1.05	1.10	1.00	0.93	1.11	0.99
	2.0–2.5 h	1.18	1.16	1.18	1.10	1.16	1.31
	2.5–3.0 h	1.20	1.24	1.34	1.20	1.23	1.29

. Regarding that the buoy observations are not independent for the verification of ERA5 wind quality, the actual error ratio in Table 4 may be even smaller, which could be tested further once an independent data source is available.

4. Results

In this section, the 2D-Var approaches with different background error correlation functions are tested in order to look for an optimal configuration for blending the HSCATs winds. Hereinafter, the results derived from different 2D-Var schemes are denoted by the used error correlation functions, namely GSF and EBEC, respectively. Section 4.1 illustrates a typical 2D-Var experiment with different settings, and then a general verification is presented in Section 4.2.

4.1. Test Case

Figure 5a shows the HSCATs winds collocated with the ERA5 background wind field (Figure 5b) on January 15th, 2022, at UTC 00:00. The 2D-Var analysis wind fields (i.e., the blended winds), based on the configurations of GSF and EBEC, are shown in Figure 5c,d, respectively. For the sake of comparison, the wind speed difference between the blended winds and the background field is shown in Figure 5e,f. It is clear that the differences are observed globally, notably over the tropics (e.g., the intertropical convergence zone, ITCZ), the western boundary ocean current systems (WBC), and the high latitudes (e.g., lat < 55°S). Such differences, on the one hand, are due to the misrepresented physical processes in the ERA5 model (e.g, ITCZ and WBC). On the other hand, they may be attributed to the persistent wind speed bias in the HSCATs data (e.g., high latitudes), because the current HSCATs wind inversion does not take the effects of sea surface temperature into account [20]. Furthermore, Figure 5f contains more small-scale variance than Figure 5e, probably due to that the used EBEC is with larger length scales than the GSF. Having different correlation length scales means they are filtering different-scale analysis increments. One could change the parameters of GSF, letting it have a similar length scale as EBEC in order to produce a finer structure than Figure 5e, but this could be tested further.

4.2. General Validation

The 2D-Var analysis results are further verified using the collocated buoy winds. The total amount of collocations is about 50 k over the above-mentioned period. Figure 6 illustrates the blended winds versus the buoy measurements for the GSF and EBEC experiments, respectively. Both results show good consistency with respect to the buoy references. Compared with the GSF experiment, the blended winds derived from EBEC are closer to both buoy measurements and ERA5 winds (see below Figure 7), which again demonstrates the advantage of the proposed EBEC in the 2D-Var wind blending.

To further validate the proposed method, Figure 7 shows the normalized histograms of wind u (a) and v (b) differences for ERA5 versus buoys (black), GSF-derived results versus buoys (red), and EBEC-derived results versus buoys (blue), respectively. The spread of u_ERA5−u_buoy (v_ERA5−v_buoy) shows a smaller deviation than that of u_GSF−u_buoy or u_EBEC-u_buoy (v_GSF−v_buoy or v_EBEC−v_buoy), which indicates assimilating HSCAT retrieved winds actually leads to “degraded” winds with respect to the background winds. Figure 7c,d are the same as the upper panels, but they used HSCAT winds as a reference. The spread of u_GSF−u_HSCAT or u_EBEC−u_HSCAT is much smaller than that of u_ERA5−u_HSCAT. This implies that bias may persist in the input data sets (e.g., the HSCAT wind retrievals), and, as such, bias correction should be carried out to further improve the blending of HSCAT winds with ERA5 background data [23].

As shown in Figure 1, about 30% of sea surface area is actually not observed by the HSCATs during the six-hourly 2D-Var analysis window. Therefore, a more comprehensive study is needed to better understand the characteristics of the blended wind fields. Practically, the blended winds are separated into two categories according to the corresponding inputs, i.e., the one with HSCATs observation(s) (denoted as ‘Obs’) and the other without any observation (denoted as ‘no Obs’), and then verified against buoy measurements separately. The statistical scores are summarized in

Table 5. SD scores of the blended winds with regards to. buoy measurements.

Statistical Scores		Speed (m/s)	Direction (°)	u (m/s)	v (m/s)
GSF	Obs	1.20	15.6	1.61	1.89
	No Obs	1.77	19.3	2.16	2.31
	Overall	1.51	17.5	1.90	2.11
EBEC	Obs	1.17	15.8	1.62	1.89
	No Obs	1.60	18.0	1.99	2.15
	Overall	1.40	16.9	1.81	2.02

.

5. Conclusions

This paper examines the feasibility of developing a six-hourly global sea surface wind product based on the scatterometers onboard the HY-2 satellite constellation. Benefiting from the orbit configuration of HY-2 satellites, i.e., HY-2B operated in a sun-synchronous orbit and HY-2C & -D operated in an inclined orbit; the three HSCATs together are able to observe more than 70% of sea surface area within six hours, providing unprecedented temporal coverage in the perspective of remote sensing. Moreover, the HSCATs vector winds are verified using both buoy measurements and ERA5 background fields. The results show that the winds derived from the three HSCATs are consistent with each other, and they are in good agreement with both ERA5 and buoy winds. Consequently, it does not need a bias correction for the wind blending procedure.

A 2D-Var approach, which is adapted from the ambiguity removal procedure for spaceborne scatterometers, is proposed to blend HSCATs winds with the ERA5 background fields. Two relevant aspects of the 2D-Var cost function are considered in this study. On the one hand, the observation error is estimated as a function of the time difference between O and B, and then a time-dependent observation-to-background error ratio is implemented in the 2D-Var procedure. On the other hand, two different background error correlation functions, namely GSF and EBEC, are tested in order to look for an appropriate 2D-Var setting for wind blending. A comprehensive verification against the buoy winds demonstrates that the use of EBEC generally leads to a better-blended wind product than the GSF, probably due to the larger length scales used for EBEC than for GSF. As such, the used EBEC adds more small-scale variance than GSF to the blended winds.

Furthermore, note that the quality of blended wind data degrades remarkably in case where there is no HSCATs observation presented (see Table 5). Such cases occur over ~30% of the sea surface. Consequently, more scatterometer observations are demanded to fill the HSCATs gaps in order to further improve the blended wind quality. Future research activities will include the blending of HSCATs winds with ASCAT wind data or the radiometer wind speed data.