Evaluation of the J-OFURO3 Sea Surface Net Radiation and Inconsistency Correction

: A new satellite-based product containing daily sea surface net radiation ( R n ) values at a spatial resolution of 0.25 ◦ from 1988 to 2013, named the Japanese Ocean Flux Data Sets with Use of Remote Sensing Observations, version 3 (J-OFURO3), was recently generated and released. In this letter, the performance of the J-OFURO3 sea-surface R n product was fully evaluated by using observations from 55 global moored buoy sites. The overall accuracy was satisfactory, with root-mean-square difference (RMSD) of 24.05 and 10.76 Wm − 2 at daily and monthly scales, respectively. However, an inconsistency issue was found in the long-term variations in the J-OFURO3 sea-surface R n values in approximately 2000; this inconsistency may be due to the replacement of the input dataset. To address this issue, a simple but effective inconsistency correction method was developed and conducted in this study. The analysis results demonstrated that the variations in the corrected J-OFURO3 sea-surface R n data were more reasonable and that its daily validation accuracy was signiﬁcantly improved by decreasing the bias from 4.67 to 0.27 Wm − 2 before the year 2000. Thereby, it is suggested that the inconsistency correction method should be applied before using the J-OFURO3 sea-surface R n data. However, the data users still should be cautious about another discontinuity issues caused by the quality of the input dataset itself.


Introduction
Sea areas cover more than 70% of Earth's surface, and they store large amounts of energy, significantly affecting the general atmospheric and oceanic circulations [1]; therefore, global and regional climates are changed through frequent air-sea flux exchanges. Hence, accurate estimations of the oceanic heat fluxes (e.g., net radiation flux and turbulent heat flux) are of great significance for evaluating the energy balance of the Earth system. The sea surface net all-wave radiation (R n ) is the difference between the total downward and upward shortwave and longwave radiation on the sea surface and describes the radiative energy balance at the air-sea interface. As an essential parameter, the sea-surface R n can directly affect the heat content of the ocean and the sea-surface temperature (SST) [2]; moreover, the sea-surface R n is also one of the most important inputs to most ocean surface physical models [3]. However, the sea-surface R n has received less attention than the turbulent flux components (i.e., latent heat flux (LHF) and sensible heat flux (SHF)), and it has been either ignored or removed from simple climatology models or parameterizations in most previous studies [4].
A new satellite-based product containing radiative fluxes on the sea surface, named the Japanese Ocean Flux Datasets with Use of Remote Sensing Observations, version 3 (J-OFURO3, hereinafter), and including data from 1988 to 2013 at a 0.25 • spatial resolution, was developed and recently released by Nagoya University, Japan (https://j-ofuro.scc. u-tokai.ac.jp/en/, accessed on 1 December 2019) [5]. The J-OFURO project focused on developing datasets of surface heat, momentum, and freshwater fluxes over global icefree sea areas. J-OFURO3 is the latest version of this project and contains significant improvements over J-OFURO [6] and J-OFURO2 [7], including accuracy improvements and a new validation scheme [5]. Various parameters in J-OFURO3 have been validated with satisfactory results, such as the sea-surface LHF and SHF [5]; hence, J-OFURO3 is popular and is widely used in various applications [8,9]. J-OFURO3 provides daily net shortwave radiation (R ns ) and net longwave radiation (R nl ) data, and these values can be directly combined to obtain daily sea-surface R n values. The sea-surface R ns and R nl values in J-OFURO3 are calculated by subtracting the upward shortwave radiation (R su ) and longwave radiation (R lu ) from the downward shortwave radiation (R si ) and longwave radiation (R li ), respectively. Specifically, R si , R li , and R su are all obtained from the International Satellite Cloud Climatology Project radiative flux D-series product (ISCCP-FD) [10] comprising the period from January 1988 to February 2000 and from the Clouds and the Earth's Energy System Synoptic Radiative Fluxes and Clouds edition 3A (CERES-3A) product at a 1 • spatial resolution [11], covering the period from March 2000 to December 2013. These products were re-gridded to a spatial resolution of 0.25 • over open waters, and their near-coastal areas were further processed with the creeping sea fill (CSF) extrapolation method [12], while the R lu values were calculated according to the Plank function, in which the SST values were obtained from J-OFURO3 as the ensemble median values based on multiple global sea-surface temperature products [5] and the sea-surface emissivity was defined as 0.984 [13]. Notably, J-OFURO3 only covers ice-free global oceans, and the sea-ice mask used was obtained from the operational sea-surface temperature and sea-ice analysis (OSTIA, [5]). To date, the performance of the J-OFURO3 sea-surface R n data is still unknown.
In this paper, observations from 55 moored buoy sites were collected to evaluate the performance of the J-OFURO3 sea-surface R n data from the 1988-2013 period at both daily and monthly scales over global ice-free oceans. Then the consistency of the long-term J-OFURO3 sea-surface R n data was examined and analyzed. To address the inconsistency issue, a simple correction scheme was proposed and conducted thereafter, and finally, the performance of the corrected J-OFURO3 sea-surface R n data was fully evaluated and discussed, including the spatiotemporal analysis. The objective of this paper is to inform the data users of the comprehensive performance of the J-OFURO3 sea-surface R n data and present an effective inconsistency correction scheme for this product.

Data Processing
Sea-surface R n is not a routine measurement on the sea surface; instead, R si , R li , and SST (unit: K) are often measured. Hence, the sea-surface R n at each moored buoy site can be calculated as follows: where α ocean is the daily sea-surface shortwave broadband albedo, ε ocean is the daily sea-surface broadband emissivity, and σ is the Stefan-Boltzmann constant (5.67 × 10 −8 W·m −2 ·K −4 ). The terms α ocean and ε ocean were usually defined as constants in many previous studies without consideration of their variations [13,14]; therefore, the use of the α ocean dataset (2000-2013) from Feng [15] and the ε ocean dataset (1988-2013) from Cheng [16], both of which have 0.05 • spatial resolutions at the daily scale used in this study, would be more reasonable than the use of data from other studies [17,18]. Note that the α ocean dataset starts in 2000; thus, the multiannual average daily α ocean values from 2000 to 2013 were used for the period before 2000 (1988 to 1999) by considering the relatively lower variations in the values of this parameter from year to year. The units for all radiative components are Wm −2 , the time format is GMT, and the downward direction is positive in this study. The R si , R li , and SST measurements were collected from 55 global moored buoy sites in 5 measuring networks (Table 1) and have been commonly used in previous studies [19,20]. Figure 1 shows the distribution of these buoy sites, which were mainly located within the sea areas between 50 • S and 50 • N. After strict quality control and a manual inspection, only the measurements with the highest quality were kept [21]. For measurements with sampling frequencies less than one day (i.e., the Ocean Sustained Inter-Disciplinary Time-Series Environment Observation System (OS) and Upper Ocean Processes Group (UOP) data), the daily mean R si , R li , and SST values were calculated by averaging all hourly observations as long as no missing data appeared in one day. Afterwards, the daily R n value was obtained according to Equations (1)-(3), and then the monthly R n value was calculated only when more than 60% of the daily R n samples in the corresponding month were available. Finally, totals of 56,460 daily and 1773 monthly in situ sea-surface R n samples were obtained and used for the accuracy validation. According to the locations of the 55 moored buoy sites, the J-OFURO3 daily seasurface R n data, which were calculated by combining the J-OFURO3 daily R ns and R nl values directly, were extracted from 1988 to 2013 and processed into a monthly scale. The sea-surface R n data in J-OFURO3 were obtained with a spatial resolution of 0.25 • from 1988 to 2013 at daily and monthly scales and then extracted according to the locations of the 55 buoy sites. Moreover, sea-surface R n data were derived from four popular reanalysis products and a ship-based product, the Japanese 55-Year Reanalysis (JRA55, 0.56 • × 0.56 • ) from the Japanese Meteorology Administration (JMA) [22], ERA5 (0.25 • × 0.25 • ) from the European Centre for Medium-Range Weather Forecasts (ECWMF) [23], the Modern Era Retrospective Analysis for Research and Applications, version 2 (MERRA2, 0.625 • × 0.5 • ) [24] from NASA's Global Modeling Assimilation Office (GMAO), the third-generation, highresolution version of the ISCCP radiative flux profile product named ISCCP-FH (1 • × 1 • , 1983-2017) [25], and the ship-based product NOCS v2.0 (1 • × 1 • , 1973-2014) by in situ observations from voluntary observing ships contained in the International Comprehensive Ocean-Atmosphere Data Set (ICOADS) [26], extracted at daily and monthly scales from 1988 to 2013, and used for comparison in this study.

Methodology
First, the accuracy of the J-OFURO3 sea-surface R n was validated against the moored samples at both daily and monthly scales. Three commonly used statistical indices (coefficient of determination (R 2 ), root-mean-square difference (RMSD), and bias) were applied to measure the uncertainties. Note that the in situ sea-surface R n used in this study was not the direct measurement, so it is more reasonable to use RMSD to represent the difference between J-OFURO3 and in situ sea-surface R n . Then the consistency performance of the J-OFURO3 sea-surface R n dataset was assessed by examining the variations in its averaged annual mean values over global ice-free seas. After that, a correction method was proposed for addressing the inconsistency issues in the J-OFURO3 sea-surface R n , as Figure 2 shows.
Three major steps were contained in this method. In Step 1, the thresholds were defined to determine the pixels with inconsistency issues. By considering the relatively stable performances and available time lengths, three products (ERA5, JRA5, and ISCCP-FH) were taken as the references. The differences between the average annual mean sea-surface R ns (∆R ns ) and R li (∆R li ) values before and after 2000 (the averages for 2000-2013, minus the averages for 1993-1999) for the three products were computed and then interpolated into a 0.25 • grid, respectively; then the thresholds for each 0.25 • pixel were defined as the sum of the median and the standard deviation of the ∆R ns and ∆R li values of the three products. In Step 2, all the J-OFURO3 pixels that needed to be corrected were determined as long as their ∆R ns or ∆R li values exceeded the corresponding thresholds. Then, in Step 3, because of the superiority of the performance in the CERES-3A radiation product both for land and sea surfaces, as indicated by previous studies [32,33], a set of regression models were built pixel by pixel by regressing the sea-surface R ns or R li value from the ISCCP-FD dataset to that from the CERES-3A dataset during their overlapping period from 2001 to 2009 to correct these inconsistent pixels, and then these models were applied on the R ns and/or R li value of these selected pixels from 1988 to 1999 by assuming that the relationships of the sea-surface R ns or R li values between ISCCP-FD and CERES-3A remained constant. The correction model is shown in Equation (4): where X is the daily sea-surface R ns or R li value from the ISCCP-FD dataset from 2001 to 2009, Y is the daily sea-surface R ns or R li value from the corresponding CERES-3A dataset from 2001 to 2009, and a and b are the regression coefficients. Finally, the annual trends observed in the sea-surface R n data from 1988 to 2013 were calculated by using the corrected J-OFURO3 sea-surface R n values with a linear regression model [34], and the results were compared to those calculated from ERA5. Figure 3 presents the explicit spatial distribution of the J-OFURO3 multiannual mean sea-surface R n data from 1988 to 2013. The sea-surface R n values decrease from the equatorial to the high-latitude sea areas, and nearly no data are available when the latitude is higher than 60 • .

Evaluation against the In Situ Observations
The scatter plots of the J-OFURO3 daily/monthly sea-surface R n data from 1988 to 2013 against the in situ sea-surface R n data are shown in    Figure 5 presents the spatial distribution of the J-OFURO3 sea-surface R n uncertainties at each individual site presented with the RMSD values. The increased RMSD values from the sites in open-sea areas to those in coastal areas at both the daily and monthly scales indicate that the performance of the J-OFURO3 sea-surface R n gradually worsened; it was speculated that the worse accuracy of this dataset in coastal areas than in open-sea areas was possibly caused by pixels being extrapolated from nearby pixels, even though no satellite observations were available [12]; further quality checks should be performed in coastal areas, as pointed out by Tomita [5].

Inconsistency Examination
The annual mean J-OFURO3 sea-surface R n values from 1988 to 2013 were calculated and compared with the values contained in three reanalysis products (JRA55, ERA5 and MERRA2) over the global ice-free oceans, as defined by J-OFURO3 and shown in Figure 6. Generally, the variations in the annual mean sea-surface R n data contained in the three reanalysis products were similar to each other, with gentle but slightly decreasing trends over the past 26 years, compared to the more dissimilar trend observed in the J-OFURO3 dataset (the black line in Figure 6), especially after 2000. Moreover, two abrupt turning points were observed in the J-OFURO3, with a sudden decrease in approximately 1992 and a sharp increase during 1999/2000. The sudden drop in 1992 also appeared in the ERA5 and MERRA2 products; this drop was thought to be caused by the satellite observations taken as the inputs of the three products being influenced by the eruption of Mount Pinatubo in 1991 [36,37]. However, the other sharp increase from 1999 to 2000 with a value larger than 3 Wm −2 only appeared in the J-OFURO3 dataset. To explore the reasonability of this sudden increase, the differences in the annual mean sea-surface R n values before 1999 and after 2000 (calculated by the means of 2000/2001 minus the means of 1998/1999) in the J-OFURO3 dataset were calculated and compared with those from the other three reference datasets of three different types, including the satellite-based product ISCCP-FD, the reanalysis product ERA5, and the ship-based product NOCS 2.0. Results shown in Figure 7 indicate that only the annual mean sea-surface R n values of the J-OFURO3 dataset changed significantly both in magnitude and spatial distribution between the periods before 1999 and after 2000, with the values decreasing even more than 20 Wm −2 , mainly in the equatorial seas, and with increasing values over nearly all other oceans (Figure 7a); the results obtained from the other three reference datasets were more similar to each other than to the J-OFURO3 dataset (Figure 7b-d). Therefore, it was assumed that the abrupt increase observed during 1999/2000 in the J-OFURO3 sea-surface R n time series was unreasonable. To determine the reason for this discrepancy, the consistency of each radiative component in the J-OFURO3 dataset was further examined. As mentioned above, the major radiative components (R ns and R li ) in the J-OFURO3 dataset were interpolated from the ISCCP-FD (before 2000) and CERES-3A datasets; hence, the average annual mean seasurface R ns and R li values over the global ice-free oceans obtained from the two data sources (ISCCP-FD from 1988-2009 and CERES-3A from 2001-2013) were examined, and then, the SST values from J-OFURO3 during 1988-2013 were also examined, as they determined the R lu values. All results are shown in Figure 7. Obvious discrepancies were seen between ISCCP-FD and CERES-3A during their overlapping period (2000-2009) in both the sea-surface R ns (Figure 8a) and R li (Figure 8b) values, especially for R li , while the relatively stable variations observed in the SST data (Figure 8c) indicated that the sharp increase in the J-OFURO3 sea-surface R n values during 1999/2000 was most likely caused by the changes in the input data sources of R ns and R li in approximately 2000.

Inconsistency Correction
As described in Section 2.2, the pixels with the inconsistency observed in the J-OFURO3 sea-surface R n time series in approximately 2000 were determined and then were addressed by correcting their R ns and R li values from the ISCCP-FD product before 2000 (1988-1999) to those from the CERES-3A product according to Equation (4). Figure 9 presents all pixels that needed to be corrected in J-OFURO3.    Figure 10a,b shows the interannual variations in the annual mean sea-surface R ns or R li over the global ice-free oceans before and after the inconsistency correction. Compared to the original values, the values of the corrected sea-surface R ns and R li and, thereby, the corrected R n values before 2000 (the black line in Figure 10) were more consistent with the values after 2000 but without much change in their variation trends; specifically, the sharp increase in approximately 2000 caused by the changes in the input datasets was eliminated in the corrected time series. After correction and comparison with the original dataset, it was found that the inconsistency issues in the J-OFURO3 sea-surface R n were caused by different radiative component at different regions (Figure 11). Overall speaking, the inconsistency in R n over the seas within the tropical region (20 • N-20 • S) was mainly caused by the inconsistent R ns (Figure 11a); while the inconsistency in R li lead to the inconsistent R n over the seas for other regions (Figure 11b,c). Moreover, it shows that the corrected results before 2000 for all regions consistent with the sea-surface R n values after 2000 very well with little changes in their variations.
The direct validation results in Figure 12 indicated that the uncertainties in the J-OFURO3 sea-surface R n during 1988-2000 decreased significantly by reducing RMSD from 39.03 to 38.59 Wm −2 and bias from 4.67 to 0.27 Wm −2 respectively after correction. Taking one moored buoy site UOP-SUB-SE (18.0 • N, 22.0 • W) as an example, Figure 13 shows that the variations in the monthly sea-surface R n time series from the corrected J-OFURO3 in 1992 were closer to the in situ R n than the original J-OFURO3. Therefore, the inconsistency correction method was so effective that the corrected J-OFURO3 sea-surface R n values performed more reasonably.
However, note that the correction method proposed in this study only addressing the discontinuity issues occurring around 2000 resulted from the inputs replacement, while the other inconsistency appearing at mid-high latitude seas in approximately 2007 ( Figure 11) caused by CERES-3A [38] was not addressed in this study.    inconsistency correction method was so effective that the corrected J-OFURO3 sea-surface Rn values performed more reasonably.  However, note that the correction method proposed in this study only addressing the discontinuity issues occurring around 2000 resulted from the inputs replacement, while the other inconsistency appearing at mid-high latitude seas in approximately 2007 ( Figure 11) caused by CERES-3A [38] was not addressed in this study.

Spatiotemporal Trend Analysis Based on the Corrected J-OFURO3 Sea-Surface R n Values
The spatiotemporal variations in the annual mean sea-surface R n over the tropical seas (20°N-20°S) from 1988 to 2013 were analyzed based on the corrected J-OFURO3, using the linear regression model. For comparison, the spatiotemporal analysis was also conducted on the sea-surface R n from the original J-OFURO3 and ERA5, and all results are given in Figure 14.
Comparing with the other two datasets, the results from the corrected J-OFURO3 dataset (Figure 14a) were more consistent with the results from the ERA5 dataset ( Figure  14c), which was similar to previous studies [39,40], but it had a large discrepancy with the results from the original J-OFURO3 dataset (Figure 14b), especially in the regions with increasing trends shown in Figure 14a,c. Therefore, this result again emphasized the importance and necessity of correcting the inconsistency in the original J-OFURO3 seasurface R n time series. Specifically, based on the results obtained from the corrected J-

Spatiotemporal Trend Analysis Based on the Corrected J-OFURO3 Sea-Surface R n Values
The spatiotemporal variations in the annual mean sea-surface R n over the tropical seas (20 • N-20 • S) from 1988 to 2013 were analyzed based on the corrected J-OFURO3, using the linear regression model. For comparison, the spatiotemporal analysis was also conducted on the sea-surface R n from the original J-OFURO3 and ERA5, and all results are given in Figure 14.
Comparing with the other two datasets, the results from the corrected J-OFURO3 dataset (Figure 14a) were more consistent with the results from the ERA5 dataset (Figure 14c), which was similar to previous studies [39,40], but it had a large discrepancy with the results from the original J-OFURO3 dataset (Figure 14b), especially in the regions with increasing trends shown in Figure 14a,c. Therefore, this result again emphasized the importance and necessity of correcting the inconsistency in the original J-OFURO3 sea-surface R n time series. Specifically, based on the results obtained from the corrected J-OFURO3, the seasurface R n increased at a mean rate of 0.257 Wm −2 over most regions of the Western and Central Tropical Pacific, which was consistent with the studies of Liu [39] and Cook [40]. According to Cook [40], the increased sea-surface R n in this region was possibly caused by the increased R ns , which were resulted from the decreased cloud amounts and the enhanced aerosol distributions or reduced water vapor droplet size. It was also seen that the sea-surface R n over the Eastern Tropical Pacific mostly decreased at a mean rate of 0.316 Wm −2 per year, which was possibly caused by the decreased top-of-atmosphere (TOA) R n within the same region as discussed by Allan [41]. The decreased sea-surface R n would enhance the sea-surface cooling and then enhance TOA albedo, therefore further enhancing the reflected shortwave radiation and reducing the net downward energy flux radiation [42], which could explain the recent slowing rate of global surface temperature rise in this region [43,44]. Moreover, another area worth noting is the Northern Indian Ocean, where the mean trend in the J-OFURO3 sea-surface R n changed from decreasing at a rate of 0.395 Wm −2 per year to increasing at a rate of 0.288 Wm −2 per year after correction and this was similar to ERA5 and more reasonable. According to the studies by Arora [45], the sea-surface R n over the Tropical Indian Ocean has shown a rise trend in the recent decade because of the increased SST in this region caused by the anomalous cyclonic circulations on either side of equator and anomalous easterlies along the Tropical Pacific Ocean resulted from the weakened relation between the Indian Ocean and the Pacific Ocean [46]. Therefore, the corrected J-OFURO3 could not only provide correct tempo-spatial variations in the sea-surface R n , but also more information on the regions different from the one from ERA5.

Conclusions
The sea-surface R n values obtained from a new satellite-based dataset, J-OFURO3, with a 0.25 • spatial resolution from 1988 to 2013, were objectively assessed in this study. Observations collected from 55 global moored buoys in five measuring networks from 1988 to 2013 were applied to evaluate the performances of the daily and monthly J-OFURO3derived sea-surfaces R n values; then the inconsistency problems in the J-OFURO3 seasurface R n time series in approximately 2000 were deeply explored. Finally, a simple but effective correction method was developed and conducted. Based on all the results, some major conclusions were drawn.

1.
The uncertainties in the J-OFURO3 sea-surface R n were accepted, with overall R 2 values of 0.86 and 0.96, RMSD values of 24.05 and 10.76 Wm −2 , and biases of 0.16 and 0.22 Wm −2 at daily and monthly scales, respectively.

2.
An abrupt increase appearing in approximately 2000 in the J-OFURO3 sea-surface R n time series was very possibly caused by the replacement of the input data sources from the ISCCP-FD dataset to the CERES-3A dataset.

3.
A simple correction method is proposed by regressing the radiative components (R ns and R li ) from the ISCCP-FD dataset to those from the CERES-3A dataset separately, pixel by pixel, and the uncertainties in the sea-surface R n were decreased remarkably by reducing the bias from 4.67 to 0.27 Wm −2 after correction. 4.
The tempo-spatial variations in the J-OFURO3 sea-surface R n were more reasonable after correction especially over tropical seas.
In summary, the performance of the J-OFURO3 sea-surface R n time series was generally satisfactory, with a relative low uncertainty, and has significant potential for wide uses in various applications. However, the limitations in this product, for example, its limited coverage (only for ice-free seas), its poor performance in the near-coastal areas, and especially the inconsistency issues cannot be ignored. It is also suggested that the correction method proposed in this study should be conducted before this product is used. More efforts are needed to improve the quality of J-OFURO3. Moreover, reliable and long-time-series sea-surface R n data with high accuracy and spatiotemporal resolution are urgently required.