Reconstruction of Hourly Gap-Free Sea Surface Skin Temperature from Multi-Sensors

Abstract

The sea surface skin temperature (SST_skin) is of critical importance with regard to air–sea interactions and marine carbon circulation. At present, no single remote sensor is capable of providing a gap-free SST_skin. The use of data fusion techniques is therefore essential for the purpose of filling these gaps. The extant fusion methodologies frequently fail to account for the influence of depth disparities and the diurnal variability of sea surface temperatures (SSTs) retrieved from multi-sensors. We have developed a novel approach that integrates depth and diurnal corrections and employs advanced data fusion techniques to generate hourly gap-free SST datasets. The General Ocean Turbulence Model (GOTM) is employed to model the diurnal variability of the SST profile, incorporating depth and diurnal corrections. Subsequently, the corrected SSTs at the same observed time and depth are blended using the Markov method and the remaining data gaps are filled with optimal interpolation. The overall precision of the hourly gap-free SST_skin generated demonstrates a mean bias of −0.14 °C and a root mean square error of 0.57 °C, which is comparable to the precision of satellite observations. The hourly gap-free SST_skin is vital for improving our comprehension of air–sea interactions and monitoring critical oceanographic processes with high-frequency variability.

1. Introduction

The sea surface temperature (SST) plays a crucial role in numerous physical, biological, and chemical processes within the Earth system, serving as a fundamental variable in weather studies and climate change assessments [1]. The SST exhibits a complex vertical structure characterized by diurnal variations (DVs) and can be categorized into various layers of temperature: interface temperature (SST_int), skin temperature (SST_skin), subskin temperature (SST_subskin), depth temperature (SST_depth), and foundation temperature (SST_fnd) [2]. Numerical weather prediction and climate simulations frequently rely on weekly or daily mean SST_fnd to calculate the air–sea fluxes, which introduces two primary sources of error. It is important to note that SST_fnd does not correspond to the temperature at the air–sea interface. A discrepancy of approximately 1 °C exists between SST_fnd and SST_skin [3]. Secondly, the use of mean SST values has the effect of obscuring DVs in the SST. The direct measurement of the SST at the air–sea interface is a challenging process; therefore, the measurement of SST_skin below the interface represents an ideal alternative. The DV in SST_skin has a significant impact on a number of processes, including near-surface ocean mixing, energy and material exchange at the air–sea interface, cloud formation, Madden–Julian Oscillation, and inter-annual atmospheric variations [4]. Indeed, the DV of SST_skin contributes to an estimated error of approximately 4.5 W/m² in global air–sea fluxes, which rises to approximately 10 W/m² in tropical regions [5] and 30–50 W/m² for extreme days in the northwestern Pacific Ocean [6]. Moreover, SST_skin is a fundamental parameter for satellite remote sensing of sea surface salinity [7]. It is therefore of the utmost importance to obtain accurate spatially and temporally gap-free SST_skin data for use in marine, meteorological, and remote sensing applications.

No single sensor can directly provide gap-free SST_skin nowadays [8]. SSTs are primarily obtained in three ways: in situ measurement, satellite remote sensing, and numerical simulation. The traditional in situ SSTs are mainly measured by thermometers equipped on vessels, buoys, and monitoring stations, which are deployed at a certain depth rather than at the actual skin layer. For example, the drifting buoy measures the SST at approximately 0.2 m, the anchor buoy measures at approximately 1 m, and the ship is even deeper. Conversely, the spatial–temporal distributions of in situ SST are sparse and limited. Despite the potential of infrared radiometers to measure the SST_skin directly, their high cost limits their widespread applicability [9]. Furthermore, long-term buoy observations indicated the existence of DVs of the SST at a certain depth, with a mean climatic state of approximately 0.2–0.6 °C [10]. The depth at which different satellite sensors measure the SST varies. Infrared sensors capture SST_skin at a depth of 10–20 μm, while microwave sensors measure SST_subskin at approximately 1 mm depths. This results in discrepancies in SST observations. Furthermore, the SST can also be obtained by numerical modeling. The general ocean model is typically employed to simulate SST_fnd at the surface layer, which is typically assumed to be a uniform layer between 5 and 10 m in depth. It would be very time-consuming to model such a thin layer of the skin layer or subskin layer.

The process of multi-sensor data fusion represents a significant approach to addressing the issue of data missing and creating a gap-free SST product. Fusion methods have evolved from subjective analysis to objective analysis [11], including polynomial fitting, stepwise correction, optimal interpolation (OI), empirical orthogonal functions, and so forth. Other forms of data fusion, such as those based on data assimilation [12,13] (e.g., the variational method and Kalman filtering), are employed to integrate observational and modeling data. The aforementioned fusion methods were typically based on the assumption that the SST remains constant over the course of a day and that the multi-source data used for reconstruction have the same spatial and temporal reference. However, discrepancies may exist between the depth and time of observation of the SST by multi-sensors. In order to eliminate or reduce the influence of these differences, the conventional algorithm generally reconstructs the daily SST_fnd or the weekly average SST [14,15], utilizing only night-time data or data with a large wind speed. In order to obtain a higher time resolution SST, such as the hourly SST_skin, it is necessary to implement depth and diurnal variation correction (or normalization) on multi-sensor SSTs before fusion. This processing not only effectively improves the accessibility of data but also greatly reduces the uncertainty of reconstruction results.

DV models can be categorized into two principal categories: empirical and physical. Empirical models are relatively straightforward to implement but have limited adaptability and lack a robust physical foundation [16]. In contrast, physical models offer greater applicability across a wider range of contexts and are less constrained by temporal and geographical factors. These models can be further subdivided into two categories: bulk models and diffusion models. Bulk models simulate the heat and momentum balance of the entire mixed layer, adjusting for momentum and buoyancy effects. However, these models assume vertical uniformity within the layer, which restricts their accuracy in modeling the vertical structure near the ocean surface [17]. To address this limitation, some models incorporate specific vertical decay structures in the near-surface layer. Nevertheless, the actual vertical structures in this region can vary due to changing weather conditions, indicating that specific vertical structure models may not fully capture the complexities of real-world scenarios. In contrast, diffusion models are capable of reasonably parameterizing the turbulent mixing and diffusion processes within the mixed layer, thereby enabling a more accurate simulation of DV near the surface in comparison to empirical and bulk models. The General Ocean Turbulence Model (GOTM), which integrates numerous validated diffusion models, has been widely utilized to simulate the diurnal variations in the vertical distribution of the SST [18,19].

The purpose of this study is to make full use of the SST retrieved from multi-sensors onboard various platforms to generate hourly gap-free SST_skin at a spatial resolution of 0.05°. The GOTM will be employed to address discrepancies in depth and DVs across multi-sensor SST data. Subsequently, the corrected SSTs will be merged using Markov methods, with any remaining gaps filled through the use of OI.

2. Materials and Methods

This section outlines the satellite data used, the GOTM employed to simulate diurnal variability, and the methodologies employed for SST correction, fusion, and gap filling.

2.1. Study Aera

The study area encompasses the East China Seas (ECSs), including the Bohai Sea, Yellow Sea, and the East China Sea, within the geographical bounds of 117–130°E and 22–42°N. The region is characterized by a complex natural geography, hydrology, and meteorology, each of which exerts a considerable influence on the variability of SST. From a geographical perspective, the ECSs are bordered by the eastern coastline of China and extend toward the Ryukyu Islands to the east. The Bohai Sea, a semi-enclosed bay, is characterized by shallow depths and a high degree of sensitivity to atmospheric changes. In contrast, the Yellow Sea and East China Sea are relatively deeper, yet still exhibit substantial seasonal and diurnal variations due to their proximity to continental land masses and exposure to monsoonal winds.

In terms of hydrology, the ECS is subject to a combination of warm and cold currents, including the Kuroshio Current and coastal currents, which give rise to pronounced SST gradients. This is further compounded by the input of freshwater from rivers, particularly the Yangtze River, which affects salinity and thermal structure, thereby intensifying stratification and SST variability. The ECS region is subject to the influence of the East Asian Monsoon system, which gives rise to pronounced seasonal shifts in SST across summer and winter. The SST DV is notably complex in this region, driven by a number of factors including high solar insolation, variable cloud cover, and fluctuations in wind speed. These conditions contribute to the formation of a dynamic SST pattern, characterized by significant diurnal heating and cooling cycles that can result in temperature differences of up to several degrees Celsius within a day. The intricate nature of this DV presents a significant challenge for the accurate monitoring of SST. In order to achieve reliable SST_skin products, it is necessary to employ sophisticated modeling and data fusion methods, as will be discussed in this study.

2.2. Satellite Data

The primary data sources for this study are infrared and microwave SST observations from geostationary and polar-orbiting satellites, with the year 2007 serving as a case example. The objective of this study is to develop a method for the reconstruction of an hourly gap-free SST_skin product derived from multi-sensors. Accordingly, the selection of the specific year is not a critical factor. Furthermore, 2007 represents an optimal experimental year as it provides access to a diverse range of available remote sensing data types. A summary of the characteristics of each dataset is provided in

Table 1. Summary of the satellite data used.

Sensor/Type	Platform/Orbit	Resolution (km)	Temporal Resolution (Local Time)	Accuracy (RMS, °C)
JAMI/infrared	MTSAT-1R/geostationary	5	1 h	0.8–1.0
AATSR/infrared	ENVISAT/polar	1	twice per day (10:30 a.m./p.m.)	0.2–0.3
AVHRR/infrared	NOAA-18/polar	4	twice per day (02:00 a.m./p.m.)	0.4–0.7
MODIS/infrared	Aqua, Terra/polar	4	twice per day (Aqua: 01:30 a.m./p.m. Terra: 10:30 a.m./p.m.)	0.3–0.5
TMI/microwave	TRMM/near-equatorial	25	twice per day (varying times)	0.6
AMSRE/microwave	Aqua/polar	25	twice per day (01:30 a.m./p.m.)	0.5–0.8
WindSAT/microwave	Coriolis/polar	25	twice per day (06:00 a.m./p.m.)	0.5–0.8

. The geostationary satellite Multifunction Transport Satellite-1R (MTSAT-1R) provides hourly temporal resolution SSTs with a spatial resolution of 5 km. The root mean square error (RMSE) of MTSAT-1R when compared to measurements taken from buoys has been found to range from approximately 0.8 to 1.0 °C [20]. The Advanced Along-Track Scanning Radiometer (AATSR) on the EnviSat satellite, which operates in polar orbit, was designed to deliver SST retrievals with an accuracy better than 0.3 °C and long-term stability exceeding 0.1 °C per decade. The Moderate-Resolution Imaging Spectroradiometer (MODIS) sensor, which is aboard the Terra/Aqua satellites, is positioned in a sun-synchronous near-polar orbit at an altitude of 705 km with a descending node at 10:30 a.m. It can achieve a high spatial resolution of approximately 1 km at nadir for both Level 1 and Level 2 datasets. However, this study employs remapped Level 3 data with a resolution of approximately 4 km, which reduce pixel noise and are comparable to MTSAT data. The accuracy of the Level 3 MODIS SST is approximately 0.3 to 0.5 °C, while the AVHRR SST ranges from 0.4 to 0.7 °C, as reported by the SST Quality Monitor (SQUAM) [21]. Despite the availability of numerous infrared SST datasets, those exhibiting significant cloud cover were excluded from the merging process. In instances where infrared SST data were unavailable, SST measurements were derived from microwave radiometers onboard the Tropical Rainfall Measuring Mission’s (TRMM) Microwave Imager (TMI), and the Advanced Microwave Scanning Radiometer for EOS (AMSR-E) and WindSat were used for fusion, despite their lower spatial resolution of 25 km. The use of microwave data is less susceptible to the effects of cloud cover, thereby facilitating the acquisition of more comprehensive datasets. However, the absence of SST data may still occur due to rain, proximity to land or ice, and high wind speeds. The RMSE for AMSR-E and WindSat SST data ranges from 0.5 to 0.8 °C [22].

The data retrieved from these sensors are accessible via the Group for High-Resolution Sea Surface Temperature (GHRSST) Data Assembly Center (ftp.nodc.noaa.gov/, accessed on 1 March 2024). All satellite data, including infrared and microwave observations, have undergone pre-processing to remove contamination from clouds, rain, or land/ice proximity. To guarantee the dependability of the dataset, quality flags for SST data, in accordance with the GHRSST Data Processing Specification (GDS), were employed. Only data flagged with quality levels 4 (“good”) and 5 (“excellent”) were included in the final analysis. The data were then resampled to the target resolution of 5 km using bilinear interpolation.

2.3. General Ocean Turbulence Model (GOTM)

The GOTM is a one-dimensional turbulence model that describes the basic thermodynamic and hydrodynamic processes. It is able to simulate SST_skin with continuous time variation in the water column [18]. The equations are solved based on the transport of heat, mass, and momentum associated with the vertical mixing processes [23]. For the purposes of this study, a non-uniform grid has been constructed for GOTM, which consists of 50 vertical levels resolving a depth down to 50 m. The grid has a much finer resolution near the surface in order to resolve the diurnal thermocline. This grid zooming is applied according to the formula in which the GOTM uses a hyperbolic tangent function to non-uniformly scale the mesh and encrypt it at the top to simulate the temperature of the skin and subskin layers. The source code for the GOTM and all related experimental data are freely available at http://www.gotm.net/ (accessed on 15 May 2023). This allows for full replication of the simulation setup and results.

The GOTM is driven by prescribed or dynamic air–sea fluxes. We chose to calculate the dynamic air–sea fluxes from the forecast weather variables by the GOTM. The internal extinction of solar radiation in the water column was modeled using 9-band ocean radiant heat parameterization based on Jerlov water types. The GOTM also uses parameterization for surface fluxes, which accounts for the cool-skin effect and diurnal warming [24]. The upwelling shortwave radiation is calculated based on sea surface albedo, while the upwelling longwave radiation is derived from the model output of SST_skin. The integration time step for the model is set to 60 s, with the simulation period spanning the entirety of the year 2007, and the results are outputted on an hourly basis. The hourly GOTM results are employed as the background field for optimal interpolation.

2.4. Ancillary Data

The meteorological variables, including the 2 m air temperature, dew point temperature, and 10 m wind speed, are accessible via the fifth-generation European Centre for Medium-Range Weather Forecasts reanalysis (ERA5) public datasets, available at https://cds.climate.copernicus.eu/datasets?q=ERA5+single+levels+reanalysis&limit=30 (accessed on 1 March 2024).

In situ SST measurements obtained from the International Quality Monitor (iQuam) are used to validate the fusion SST_skin product. The in situ data underwent rigorous quality control procedures, and only the highest quality data were retained. Drifter buoys, which measure the SST at approximately 0.2 m depth, were used for comparison with the modeled and fused SST_skin. The depth correction procedure (described in Section 3.2) was applied to align the in situ SST data with the skin temperature, thereby ensuring an accurate comparison. The estimated error for drifter buoy SST is estimated to be approximately 0.24 °C [25].

2.5. Methods

The overall flowchart is presented in Figure 1. The multi-source SST data used in this study were primarily derived from polar-orbiting infrared and microwave remote sensing, as well as geostationary satellite infrared observations. Weather variables from the ECMWF reanalysis were employed as inputs for the GOTM to simulate diurnal variations. After the preprocessing steps, including quality control, diurnal and depth corrections were applied. Then, the Markov method was employed to blend the normalized SST data, with OI utilized to fill any remaining gaps. The principal processing workflow is outlined as follows:

(1)

Depth correction: SST observations from infrared and microwave sensors, which measure different depths, were normalized using the GOTM-modeled vertical profile. The temperature differential between SST_skin (10–20 μm) and SST_subskin (1 mm) was determined for each observation point, and the microwave SST data were adjusted to align with the skin temperature by subtracting this differential. The SST_fnd is taken at an approximate depth of 10 m. The satellite-derived skin or subskin layer temperature is adjusted by subtracting the difference between the GOTM-simulated temperature at the corresponding depth and the temperature at 10 m.

(2)

Diurnal Variation Correction: SST observations at varying local times were normalized to the same observation time frame by using the GOTM-modeled diurnal variability. This correction was applied to both infrared and microwave satellite data to eliminate biases that may have been introduced by temporal differences.

(3)

Markov Estimation for Fusion: Following the depth and diurnal variation correction, the SST data from multiple sensors were fused using the Markov estimation method, which assigns weights to each observation based on its accuracy. This step combines SST_skin data from multiple sensors to create a more comprehensive and accurate representation of the spatial coverage, while also reducing the impact of random errors.

(4)

Optimal Interpolation (OI) for Gap Filling: Notwithstanding the enhancements in data coverage subsequent to fusion, some gaps persisted due to cloud contamination or the absence of sensor data. The OI technique was then applied to fill these remaining gaps. A covariance function was used to estimate the SST values for missing points based on the surrounding data, considering both spatial and temporal correlations. The parameters for OI, including spatial correlation scales, were determined empirically based on the satellite data.

3. Results

3.1. Depth and Diurnal Correction for Multi-Sources SST

The discrepancies in depth and time between SST observations from different sensors were corrected using the GOTM. It is imperative that these corrections are applied in order to produce a consistent, high-resolution SST_skin product from infrared and microwave satellite data, which measures different depths and times.

3.1.1. DV Modelling

The accurate simulation of the DV of an SST profile provides a solid foundation for the correction of multi-source SSTs acquired at different depths and times. Since this process involves complicated nonlinear interactions at the air–sea interface, several key issues are of concern in this study, i.e., the grid resolution, choice of forcing data, penetration of solar radiation, and choice of turbulence parameters (for details, see Section 2.3). The meteorological variables are used to calculate the surface fluxes and force in the model. The penetration of solar radiation follows the exponential law, and the absorption coefficients depend on the water type and are prescribed by means of Jerlov classes [26]. The scheme devised by Fairall et al. [24] is implemented into the turbulence parameterizations.

The DV of SST_skin is defined as the difference between the SST_skin and the SST_fnd. Figure 2 illustrates the typical diurnal variability of SST_skin on 8 May 2007 modeled by the GOTM. The DV is very small during nighttime and close to zero in the early morning, the coldest period during the day. The temperature increases gradually from 9:00 a.m. onwards and approaches the maximum of up to 5 °C around 3:00 p.m. Thereafter, warming begins to decrease, finally disappearing at midnight. The corresponding DV of SST_skin derived from MTSAT observations is used to validate the DV simulated by the GOTM.

Figure 3 illustrates that the distributions of the differences between MTSAT observations and the GOTM are quite close to Gaussian. The empirical histograms display a slight asymmetry and peakedness, as indicated by the skewness (−0.317) and kurtosis (2.052). The global differences (mean and median) are close to zero. The maximum and minimum differences are 3.367 °C and −2.16 °C, respectively. Most of the data range from −1 °C to 1 °C with a variance of 0.116 °C. It can therefore be concluded that the GOTM is an effective tool for modeling the DV of the SST, and the resulting data are valuable for the normalization of multi-sensor SSTs. The modeled diurnal signals were validated against MTSAT satellite observations, and the results were closely matched, with a mean difference of less than 0.1 °C and a variance of 0.116 °C. This provides empirical verification that the DV simulated by the GOTM is applicable to our study.

3.1.2. Depth Correction

Infrared and microwave SST retrievals are complementary and can be combined to obtain the gap-free SST. The infrared sensor measures the temperature at the skin layer (10–20 μm, referred to as SST_skin), whereas microwave instruments are able to measure the subskin layer (approximately 1 mm, referred to as SST_subskin). Typically, SST_skin is observed to be slightly cooler than SST_subskin by a few tenths of a degree, which is thought to be mainly due to the vertical heat flux at the air–sea interface. It is essential to take this discrepancy into account when fusing infrared and microwave SST data in order to guarantee consistency between the two sources. The GOTM is capable of capturing the diurnal thermocline, including the skin and subskin layers, which correspond to the infrared and microwave observations, respectively.

Figure 4a shows an example of the vertical variation profile of the SST at 2 p.m. on 8 May 2007, located at 122°E and 35.25°N, as simulated by GOTM. The temperature profile near the surface initially exhibits an increase with depth, reaching its peak in the subskin layer before its subsequent decrease and stabilization below approximately 5 m, which is known as SST_fnd. The temperature difference between SST_skin and SST_subskin is evident in this profile, with the skin layer being slightly cooler than the subskin layer due to molecular diffusion processes at the surface. Figure 4b illustrates the spatial distributions of the difference between SST_skin and SST_subskin, as simulated by the GOTM. In the low-wind-speed areas, such as the Yellow Sea and the northern part of the East China Sea, the SST_skin is typically approximately 0.1–0.2 °C cooler than the SST_subskin. In contrast, this difference diminishes or disappears in regions with high wind speeds, such as the southern part of the East China Sea. This finding highlights the potential inaccuracies of using a simple empirical value or formula for depth correction as it may not fully account for regional variations in wind speed and other environmental factors [27].

In this study, the microwave SST data are corrected to the depth of the infrared SST (skin layer) by subtracting the modeled temperature difference between SST_subskin and SST_skin, based on the GOTM simulation. This methodology provided a more accurate depth correction compared to using constant correction values. The same approach was applied to correct the in situ SST data, ensuring consistency across all SST observations.

3.1.3. Diurnal Variation Correction

Since the GOTM models the DV of the SST well, it allows the satellite-derived SST from different local times on the same day to be normalized to a common reference time. The DV correction is applied as follows:

$${\mathrm{SST}}_{\mathrm{SAT}}\left(\mathrm{j}\right)={\displaystyle \frac{{{\displaystyle \sum }}_{\mathrm{i}=1}^{\mathrm{N}}\left({\mathrm{SST}}_{\mathrm{SAT}}\left(\mathrm{i}\right)+\left({\mathrm{SST}}_{\mathrm{GOTM}}\left(\mathrm{j}\right)-{\mathrm{SST}}_{\mathrm{GOTM}}\left(\mathrm{i}\right)\right)\right)}{\mathrm{N}}}$$

(1)

where ${\mathrm{SST}}_{\mathrm{SAT}}$ represents the satellite observations, j is the target time for normalization, i represents the effective observation times during the day, excluding the target time, and N is the total number of effective observations for each sensor per day (if there are no other effective observations, N = 0, and time normalization is not required).${ \mathrm{SST}}_{\mathrm{GOTM}}$ is the SST value modeled by the GOTM for the specified time interval between i and j. For example, polar-orbiting satellites provide two observations per day. In the event that there is one effective observation aside from the target time, then N = 1. In the case of geostationary satellites, there are 24 potential observations per day. In the event that the objective is to normalize to a specific target time, for example, 15:00, there will be 23 effective observations for other times, thus N = 23. The normalized SST at time j is calculated as the average of the satellite observations, adjusted by the GOTM DV.

To demonstrate this procedure, consider the geostationary satellite MTSAT SST observations on 8 May 2007. Figure 5a shows the original SST from MTSAT before diurnal correction. Due to sensor issues, there are missing observations at 13:00 and 21:00. Additionally, other gaps exist at various times, caused by cloud cover, heavy aerosols, and solar glint. The high temporal resolution of MTSAT (hourly data) helps to reduce these issues, as the contamination tends to move or dissipate over time. After applying the diurnal variation correction, the SST data from different times are normalized to the same reference time, and any missing data within the day are filled when possible. This significantly increases the data availability, as demonstrated in Figure 5b, where the normalized SST data fill in many of the gaps present in the original data. For instance, before normalization, large portions of the north Yellow Sea and Bohai Sea lacked SST observations between 09:00 and 24:00. After normalization, the coverage rate was greatly improved. In this way, time normalization can also be employed as an effective technique for filling gaps in satellite remote sensing data, particularly for the purpose of time interpolation.

Compared to traditional mathematical interpolation methods, the physical basis of this study is clearer as it relies on the thermal and dynamic processes modeled by the GOTM rather than purely statistical techniques. In Figure 5b, the diurnal variation of the Taiwan Warm Current (TWC) is clearly visible. As solar radiation intensifies throughout the day, the temperature of the TWC increases, making the temperature difference between the current and the surrounding waters more pronounced. The maximum temperature is reached at approximately 15:00, after which a gradual decline in temperature is observed. Similar patterns are observed in other smaller branches of the Kuroshio Current.

This finding emphasizes the necessity for high-temporal-resolution SST data, which are essential for capturing dynamic oceanographic phenomena with high-frequency changes. The use of daily average SST data alone would result in the obscuration of important details pertaining to diurnal processes and their interactions with the surrounding environment.

3.2. Generation of Hourly Gap-Gree SSTskin

3.2.1. Markov Estimation

After normalizing the SST observations to the same time, the spatial coverage of satellite-derived SST data was significantly improved. As an example, we consider the SST observations from 10:30 a.m. on 8 May 2007. At this time, SST data were available from eight different sensors across the central East China Sea, as illustrated in Figure 6a (for a comprehensive list of the sensors, please refer to Table 1). In the majority of regions, data were provided by more than six satellite sensors, while in coastal areas, up to three sensors contributed, though some microwave data were contaminated near land due to the nature of the signal. Overall, the availability of SST data increased dramatically to 87%, largely due to the high temporal resolution of geostationary satellites like MTSAT, which helped mitigate the effects of cloud cover. When considering all sensors, the availability reached 94%, except in certain parts of the northern Bohai Sea where data gaps remained.

In regions where multiple observations from different sensors are available, the challenge is to combine these data in a way that maximizes accuracy while minimizing errors. A common approach is to use a least-squares method, where the observations from each sensor are averaged to estimate the true SST value. However, this method does not account for the varying levels of uncertainty and error in the observations from different sensors. Consequently, the application of equal weighting to all sensor data may result in the introduction of inaccuracies.

To address this issue, we employed Markov estimation [28], a more sophisticated approach that assigns weights to each sensor’s observations based on their accuracy. In essence, Markov estimation adjusts the contribution of each sensor’s data according to the reliability of the data, as determined by the error characteristics of the sensor in question. The application of this weighting ensures that sensors with greater accuracy exert a greater influence on the final SST estimate, while sensors with lesser accuracy contribute less.

This method minimizes the overall error in the final SST estimate by incorporating both the observed data and the uncertainty associated with each sensor. The application of Markov estimation to the fusion process yields a more accurate representation of the true SST at each location. This approach leverages the strengths of each sensor while accounting for their respective limitations, thereby enhancing the overall accuracy of the estimation. The fusion results from all available sensors at 10:30 a.m. on 8 May 2007 are shown in Figure 6b. The use of Markov estimation allowed for more accurate and seamless integration of multi-sensor SST data, enabling the capture of significant oceanographic characteristics with high fidelity. This method not only improved the spatial coverage of the dataset but also enhanced the overall quality of the SST dataset by effectively managing the uncertainties inherent in the observations.

3.2.2. Gap Filling via Optimal Interpolation

Optimal interpolation (OI) is applied to address the remaining gaps after Markov estimation. In this approach, both the distance between observations and gaps, as well as the influence of changes in the relative positions of observations, are used to estimate the interpolation weights, resulting in more reasonable estimations. For further details regarding the implementation of OI, please refer to [10], which produced daily and weekly fused SST data. In this study, the focus is on the parameterization of the background error covariance function $\mathrm{B}\left(∆\mathrm{x},∆\mathrm{t}\right)$, which distinguishes hourly from daily OI schemes. This covariance function is expressed as the product of spatial and temporal components

$$\mathrm{B}\left(∆\mathrm{x},∆\mathrm{t}\right)={\mathrm{b}}_{\mathrm{s}}\left(∆\mathrm{x}\right){\mathrm{b}}_{\mathrm{t}}\left(∆\mathrm{t}\right)$$

(2)

where $∆\mathrm{x}$ and $∆\mathrm{t}$ represent the spatial and temporal distances, respectively. A negative-squared exponential function is adopted to describe the spatial correlation between SST observations

$${\mathrm{b}}_{\mathrm{s}}\left(∆\mathrm{x}\right)=\mathrm{Aexp}\left[{\displaystyle \frac{-{\left({\mathrm{x}}_{\mathrm{i}}-{\mathrm{x}}_{\mathrm{j}}\right)}^2}{{\lambda }_{\mathrm{x}}^2}}+{\displaystyle \frac{-{\left({\mathrm{y}}_{\mathrm{i}}-{\mathrm{y}}_{\mathrm{j}}\right)}^2}{{\lambda }_{\mathrm{y}}^2}}\right]+\mathrm{c}$$

(3)

where x and y are coordinates in the zonal and meridional directions, respectively, i and j are different observations, and $\lambda$ is the spatial correlation scale. Given that SST is predominantly influenced by solar radiation, which varies more significantly in the zonal direction, the spatial components are separated into zonal and meridional directions, simplifying to:

$${\mathrm{b}}_{\mathrm{s}}\left(\mathrm{x}\right)={\mathrm{A}}_{\mathrm{x}}\exp \left[{\displaystyle \frac{-{\left({\mathrm{x}}_{\mathrm{i}}-{\mathrm{x}}_{\mathrm{j}}\right)}^2}{{\mathsf{\lambda }}_{\mathrm{x}}^2}}\right]+{\mathrm{c}}_{\mathrm{x}}$$

(4)

and

$${\mathrm{b}}_{\mathrm{s}}\left(\mathrm{y}\right)={\mathrm{A}}_{\mathrm{y}}\exp \left[{\displaystyle \frac{-{\left({\mathrm{y}}_{\mathrm{i}}-{\mathrm{y}}_{\mathrm{j}}\right)}^2}{{\mathsf{\lambda }}_{\mathrm{y}}^2}}\right]+{\mathrm{c}}_{\mathrm{y}}$$

(5)

Parameters ${\mathrm{A}}_{\mathrm{x}}$, ${ \mathrm{A}}_{\mathrm{y}}$, λ_x, λ_y,${ \mathrm{c}}_{\mathrm{x}}$, and ${\mathrm{c}}_{\mathrm{y}}$ are determined through least-squares fitting to MTSAT observations (Figure 7a,b) and then the final empirical function where A is 0.410936 (the average of ${\mathrm{A}}_{\mathrm{x}}$ and ${ \mathrm{A}}_{\mathrm{y}}$), C 0 is.503 (average of ${\mathrm{c}}_{\mathrm{x}}$ and ${ \mathrm{c}}_{\mathrm{y}}$), λ_x is 85 km, and λ_y is 100 km:

$${\mathrm{b}}_{\mathrm{s}}\left(∆\mathrm{x}\right)=0.410936 \exp \left[-0.0100627\mathrm{dx}-0.01168\mathrm{dy}\right]+0.503$$

(6)

Similarly, the temporal covariance component ${\mathrm{b}}_{\mathrm{t}}\left(∆\mathrm{t}\right)$ is derived from SST data based on a time lag of Δt (the black line in Figure 7c), which reveals a decreasing correlation over the first 12 h, followed by periodic oscillations every 24 h. Multi-sensor data within ±3 h of the interpolation time over a ±2-day window are used in the OI input fields. The temporal behavior is fitted by an analytical exponential function, as shown in Figure 7c.

$${\mathrm{b}}_{\mathrm{t}}\left(∆\mathrm{t}\right)=\left\{\begin{array}{c}0.8277\exp \left(-0.28607∆\mathrm{t}\right)+0.285978,\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}0\le ∆\mathrm{t}<24\\ 0.4966\exp \left(-0.28607\left(∆\mathrm{t}-24\right)\right)+0.285978,\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}24\le ∆\mathrm{t}<48\\ 0.3599\exp \left(-0.28607\left(∆\mathrm{t}-48\right)\right)+0.285978,\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}48\le ∆\mathrm{t}<72\end{array}\right.$$

(7)

After normalization, correction, fusion, and interpolation, gap-free hourly SST_skin data are generated, as illustrated in Figure 8. It is noteworthy that a significant diurnal warming event is observed west of the Kuroshio, with a diurnal cycle demonstrated at 124°E and 28°N (Figure 9). The minimum SST_skin occurs at approximately 06:00 a.m. and then increases until 14:00 p.m. with a remarkable magnitude of approximately 2.5 °C. Then, SST_skin decreases in the late afternoon and tapers off at night. These results suggest that multi-sensor SST_skin data effectively capture diurnal variations, warranting further validation against in situ data.

3.3. Validation with In-Situ Data

The in situ SST measurements from iQuam were utilized to validate the fusion results for the entire year of 2007. The iQuam data, known for undergoing critical quality control, are widely used for satellite calibration and validation purposes. Since the in situ SST measurements are taken at a certain depth, it is necessary to apply a depth correction to adjust the measurements to the skin layer, as described in Section 3.1. The time collocation is set to ±1 h and the spatial collocation window is taken as 0.1° × 0.1°. Figure 10a illustrates that the majority of in situ SST_skin and fusion SST_skin data points are evenly distributed around the 1:1 reference line. The overall statistics are based on 40,959 matchups, indicating a mean bias of −0.14 °C and a root mean square error of 0.57 °C. These results suggest that the accuracy of the fusion SST_skin is comparable to the satellite observations. A systematic negative bias is observed when temperatures are below 9 °C, and the reasons for this bias require further investigation.

The hourly mean bias and standard deviation are illustrated in Figure 10b. A small negative bias of approximately −0.1 °C is observed during the nighttime, which is likely attributable to two underlying factors. Firstly, in situ sensors, such as drifting buoys, fluctuate in response to wave action, resulting in alterations to the depth of observation. As drifting buoys lack recorded depth data, a standard depth of 0.2 m is assumed for depth correction, which may introduce errors in the resulting data. Under calm conditions, the depth of a drifting buoy typically exhibits fluctuations of approximately 5 cm. It can be reasonably assumed that the DV between SST_0.2m and SST_fnd will reach a value of 1–5 °C and decrease linearly with depth, while the temperature variation resulting from buoy depth fluctuations would be approximately 0.01–0.05 °C. In rougher sea conditions, buoy depth variations may increase; however, under such conditions, the stratification of the surface layer is disrupted, leading to more homogeneous mixing. As a result, the impact of buoy depth variations on temperature measurements becomes negligible. Secondly, flaws in the nighttime cloud detection algorithm may cause cloud contamination of the satellite-retrieved SST, which could introduce negative bias if these contaminated pixels are incorporated into the fusion. A positive bias, ranging from 0.1 °C to 0.2 °C, is observed during the daytime, potentially resulting from errors in in situ data depth correction, especially during significant diurnal warming events in the skin layer.

4. Discussion

4.1. Uncertainty of the Procedure

The uncertainty of the procedure can be attributed to three main sources: the remote sensing SSTs, the results of the GOTM, and the interpolation algorithm.

4.1.1. Uncertainty of Satellite-Derived SST

It is generally accepted that the retrieval of SST from infrared sensors, with the exception of the AATSR, does not provide an accurate representation of the skin temperature. The AATSR uses a “dual-view” algorithm that enhances atmospheric correction by employing two distinct atmospheric path lengths, thereby facilitating more accurate SST_skin observations. However, other infrared sensors receive radiance from approximately 10 µm depths at spectral channels centered around 3.7 µm (nighttime), 10.8 µm, and 12.0 µm. These data are processed using regression against drifting buoy SST measurements, typically taken at a depth of between 0.2 and 0.3 m to estimate SST_subskin. To account for the cool-skin effect, a constant value of 0.17 °C is subtracted from SST_subskin to obtain SST_skin [29]. The propagation of uncertainties in these corrections has the potential to introduce biases into the fused SST_skin products, as illustrated in Figure 10. The accuracy of satellite-derived SSTs is generally robust; however, it is susceptible to variations in atmospheric conditions, calibration techniques, and sensor types.

4.1.2. Uncertainty of GOTM Simulation

The accuracy of the GOTM simulations is highly sensitive to input weather parameters, particularly wind speed and solar radiation. The wind speed, relative humidity, and cloud cover used in this study are sourced from ERA5 reanalysis datasets with spatial resolutions of 0.25° × 0.25°. However, these resolutions are insufficient to capture the finer-scale atmospheric dynamics, especially in coastal regions where variability is higher. Additionally, the GOTM is a one-dimensional model, which does not account for horizontal advection, meaning it cannot simulate the diurnal variations caused by horizontal ocean currents. It has been demonstrated that factors such as turbidity and phytoplankton blooms can exert an influence on the diurnal variation of the SST [30]. However, the GOTM only provides a rough approximation of these effects by considering six general water types. In future studies, the incorporation of ocean color remote sensing data, including the transparency, total suspended matter, and chlorophyll-a concentration, could improve the accuracy of modeled solar radiation transmission through the water column.

4.1.3. Limitation of the OI Method

In this study, microwave SST data with a lower spatial resolution (0.25° × 0.25°) were interpolated onto a finer grid (0.05° × 0.05°) using bilinear interpolation. Although this process is computationally efficient, it has the potential to smooth out smaller oceanic features and create artifacts, such as jagged edges, in the SST field. Furthermore, the OI method relies on a fixed-scale error covariance model, which may not be adequate for representing the complex variability in coastal areas where environmental conditions are prone to rapid change. A more adaptive, multi-scale error covariance approach could better capture these complexities. Furthermore, due to computational constraints, only a limited set of data points surrounding the target location are used for interpolation, which may lead to suboptimal estimations on a global scale. The implementation of a three-dimensional variational data assimilation technique could potentially improve the accuracy of SST estimates by incorporating a more comprehensive dataset.

4.2. Implications

The hourly gap-free SST_skin products generated in this study offer valuable insights into sea–air interactions. It has been demonstrated that improving SST_skin accuracy by 1 °C can reduce heat flux uncertainty by approximately 10 W/m² [5]. Moreover, the rate of gas exchange, particularly for carbon dioxide, is closely correlated with the temperature gradient between the ocean surface and the atmosphere [31]. It is therefore imperative that SST_skin data be corrected at high temporal resolution in order to capture high-frequency physical processes, such as the ocean carbon cycle. Furthermore, these data have the potential to enhance the monitoring of dynamic oceanic processes, including suspended matter concentration [32], frontal systems [33], and red-tide events [34], especially during the daytime. Finally, the hourly SST_skin data can be used as input for atmospheric, oceanic, or biogeochemical models, thereby enhancing the simulations of upper-ocean circulation and thermodynamic structures across different time scales.

5. Conclusions

In light of the fact that no single sensor is currently capable of providing continuous, gap-free SST_skin data, data fusion emerges as a promising solution to fill these gaps and create a comprehensive SST_skin product. Infrared and microwave sensors aboard geostationary and polar-orbiting satellites complement each other in terms of their respective spatial and temporal coverage. Significant progress has been made in the generation of daily and weekly SST_fnd products, which typically rely on data collected during nighttime hours or in conditions of high wind speeds. The aim of this study is to investigate a procedure for the production of an hourly gap-free SST_skin dataset. The key conclusions are as follows:

(1)

Correction of observation depth differences between multi-sensors based on the GOTM: Infrared sensors measure SST_skin at a depth of 10–20 µm, while microwave instruments capture SST_subskin at a depth of 1 mm. The GOTM effectively simulates the diurnal thermocline, including both the skin and subskin layers. This allows for the adjustment of microwave SST to be adjusted to match the infrared-measured SST at the same depth (skin layer). This is achieved by subtracting the modeled difference, rather than applying a constant correction. Similarly, in situ SST measurements at varying depths can be adjusted to represent the skin layer for the purposes of validation.

(2)

Correction of observation time differences between multi-sensors using the GOTM: The diurnal variation in SST_skin is a well-known natural fluctuation. Accordingly, it is essential to account for discrepancies in the observation times of satellite observations on the same day, particularly during periods of significant diurnal warming. The GOTM can accurately simulate diurnal SST signals, thereby enabling the normalization of satellite-derived SST data collected at different local times to a common time. This time normalization considerably enhances spatial coverage and serves as an effective technique for filling temporal gaps in satellite remote sensing.

(3)

The accuracy of the hourly gap-free SST_skin data is comparable to that of satellite observations. Markov estimation is employed to fuse the normalized SST data, and OI is used to fill the remaining gaps. The covariance function for background error in the OI process is parameterized with data acquired within ±3 h of the target time, spanning a ±2 day window. The spatial correlation scales have been set at 100 km in the zonal direction and 85 km in the meridional direction. The hourly gap-free SST_skin product was validated against in situ data for the entire year of 2007, demonstrating an overall bias of −0.14 °C and a root mean square error of 0.57 °C, which is comparable to satellite-derived SST measurements.

The hourly gap-free SST_skin produced in this study captures the continuous spatio-temporal distribution of the diurnal cycle, offering valuable insights into air–sea interactions and the physical pump mechanisms involved in the ocean carbon cycle. It is also well-suited for monitoring high-frequency oceanic processes. Additionally, the hourly gap-free SST_skin product can be used to force atmospheric and marine models, enhancing predictions of upper-ocean circulation and thermal structures across different temporal scales.

Although this study has demonstrated the feasibility of an operational diurnal gap-free SST_skin product, further work is required. The ability to resolve diurnal variations is heavily dependent on the availability of comprehensive input data. Moreover, additional validation and refinement of the GOTM is necessary to enhance its ability to reproduce diurnal signals across multiple sensors.