Evaluation of the Reanalysis Surface Incident Shortwave Radiation Products from Ncep, Ecmwf, Gsfc, and Jma Using Satellite and Surface Observations

Solar radiation incident at the Earth's surface (R s) is an essential component of the total energy exchange between the atmosphere and the surface. Reanalysis data have been widely used, but a comprehensive validation using surface measurements is still highly needed. In this study, we evaluated the R s estimates from six current representative global reanalyses (NCEP–NCAR, NCEP-DOE; CFSR; ERA-Interim; MERRA; and JRA-55) using surface measurements from different observation networks [GEBA; BSRN; GC-NET; Buoy; and CMA] (674 sites in total) and the Earth's Radiant Energy System (CERES) EBAF product from 2001 to 2009. The global mean biases between the reanalysis R s and surface measurements at all sites ranged from 11.25 W/m 2 to 49.80 W/m 2. Comparing with the CERES-EBAF R s product, all the reanalyses overestimate R s , except for ERA-Interim, with the biases ranging from´2.98 W/m 2 to 21.97 W/m 2 over the globe. It was also found that the biases of cloud fraction (CF) in the reanalyses caused the overestimation of R s. After removing the averaged bias of CERES-EBAF, weighted by the area of the latitudinal band, a global annual mean R s values of 184.6 W/m 2 , 180.0 W/m 2 , and 182.9 W/m 2 were obtained over land, ocean, and the globe, respectively.


Introduction
Solar radiation incident at the Earth's surface (R s ) is one essential component of the total energy exchange between the atmosphere and the Earth's surface.Worldwide monitoring of R s from ground-based stations began in the late 1950s [1].It is still insufficient to derive the global radiation distribution from ground measurements alone due to the sparsity and heterogeneity of stations [2][3][4].However, observational measurements provide the best estimate of the state of R s where they are taken.Currently, a number of gridded global R s products exist from remote sensing [2,[5][6][7] and reanalysis [8][9][10][11][12].Satellite remote sensing is one of the most practical ways to derive R s with relatively higher spatial resolution and accuracy, but the temporal coverage is always limited.
Besides satellite remote sensing methods, estimates from reanalyses are another feasible way to produce long-term global R s products.A reanalysis generates a dynamically consistent global analysis of the state of the atmosphere over an extended period of time and is a combination of model (geophysical fluid-dynamical model of the atmosphere containing important physical processes like radiative transfer and convection) and measurements, using observations to constrain the dynamical model to optimize the properties of complete coverage and accuracy [13].The insertion of observational information into model integration is called data assimilation.Although reanalyses estimates are comprehensive in nature, the land surface products obtained from reanalyses do have some problems resulting from the assimilation of data obtained primarily from atmospheric profiles, forcing near-surface meteorology with the estimated precipitation.These along with other error sources (e.g., clouds and aerosols,) tend to cause systematic biases in reanalysis R s products.Therefore, it is necessary to verify reanalysis R s products with available ground measurements before they are applied to climate studies.
The National Centers for Environmental Prediction (NCEP), the European Centre for Medium-Range Weather Forecasts (ECMWF), the National Aeronautics and Space Administration (NASA) Goddard Space Flight Center (GSFC)'s Global Modeling and Assimilation Office (GMAO), and the Japan Meteorological Agency (JMA) are the four representative organizations that produce the "second generation" reanalysis datasets.These reanalyses have been used in various applications, for example, to detect the long-term climate trends [14][15][16][17], to generate the atmospheric forcing data [18,19], and to investigate the water and energy budgets between the Earth's surface and the atmosphere [13,[20][21][22][23][24][25].
Although the R s products from these reanalyses are the best approximation of R s based on both data and models, the various reanalysis R s products have been proven to have various deficiencies.Comparisons and verifications of reanalysis R s products against ground measurements and satellite estimates of R s have been performed over individual sites, specific regions, and the globe.Comparison with the ground measurements have revealed larger positive biases in current reanalysis R s products [21, [23][24][25][26][27][28][29].For instance, Brotzge [23] found that the NCEP-National Center for Atmospheric Research (NCEP-NCAR) global reanalysis data overestimated R s by 17%-27% using two years' ground measurements at two oasis sites.Wang and Zeng [25] also found that six reanalysis R s products (MERRA, NCEP-NCAR, ERA-40, ERA-Interim, CFSR, and GLDAS) slightly overestimated R s by 1.56 W/m 2 -5.00 W/m 2 at a daily time scale at nine stations from the Coordinated Enhanced Observing Period (CEOP)-Asia-Australia Monsoon Project (CAMP/Tibet).NCEP reanalysis R s data exceeded ground measurements by 40 W/m 2 to more than 100 W/m 2 at 46 sites provided by the Climate Data Center of the Chinese Meteorological Administration (CDC/CMA) [26].Based on the satellite R s products, large uncertainties were also found in the current reanalysis R s products [13,27,[30][31][32].Positive biases in the NCEP-NCAR data ranging from 40 W/m 2 to 80 W/m 2 were also found in Europe [30].Although these R s products have been widely evaluated in previous studies, the number of observations used in those studies was limited, and the selected sites were located in different regions.In addition, fewer studies have shown whether the reanalysis R s products are consistent with regard to seasonal variations, spatial distributions, and global annual mean values.Moreover, comparisons of the reanalysis R s product consistency in trend analysis with satellite-derived R s products are rare, especially in the case of new satellite missions that have improved our knowledge of the energy exchange between the sun, earth, and space [33].Additionally, critical verifications of these reanalysis R s products are helpful to researchers who will use them for land-atmosphere interaction and climate change studies.
Therefore, satellite R s estimates from the Earth's Radiant Energy System, Energy Balanced And Filled (CERES-EBAF) and quality-controlled ground measurements from the Global Energy Balance Archive (GEBA), the Baseline Surface Radiation Network (BSRN), the Greenland Climate Network (GC-NET), Buoy observations over tropical oceans, and the Climate Data Center of Chinese Meteorological Administration (CDC/CMA) were used to assess the estimated R s of the latest versions of six representative global reanalysis R s products in this study.Furthermore, seasonal variations, spatial distributions, global annual mean values, and long-term trends of R s from these six gridded reanalysis R s products were also compared.This paper is organized as follows.The data used in this study are introduced in Section 2; the results and analysis are introduced in Section 3; a short summary and conclusions are provided at the end of the paper.

Global Reanalyses
Six reanalysis products [NCEP-NCAR reanalysis, NCEP-DOE reanalysis, Climate Forecast System Reanalysis (CFSR) from NCEP; ECMWF Interim Reanalysis (ERA-Interim) from ECMWF; Modern-Era Retrospective Analysis for Research and Applications (MERRA) reanalysis from GSFC; and JRA-55 from JMA] were used in this study.These products are different in many aspects, such as the numerical schemes, the physical parameterizations in their numerical models, the qualities and quantities of observational data used in the assimilation processes, and the assimilation schemes [24,25].Monthly mean R s data from these six reanalysis datasets were evaluated in this research.The details of these datasets are described in the following subsections.

ERA-Interim
ERA-Interim is a global atmospheric reanalysis produced by ECMWF [40] with data available from ECMWF since 1979.The ERA-Interim project was conducted in part to prepare for a new atmospheric reanalysis to replace ERA-40.The data assimilation method used to produce ERA-Interim is based on an updated version of the ECMWF forecasting model [version cycle 31r1 (Cy31r2)].Besides the increased resolution compared to ERA-40, ERA-Interim includes a four-dimensional variational analysis (4D-Var).The 4D-Var is a temporal extension of 3D-Var [42].The cost function was used to minimize the difference between the model and the observation with a 12 h assimilation window.The spatial resolution of the data set is 0.75 ˝(approximately 80 km) on 60 vertical levels from the surface up to 0.1 hPa.Radiation is calculated based on the Rapid Radiation Transfer Model (RRTM) [43].

MERRA
The MERRA product [11] is a second reanalysis project from NASA for the satellite era, from 1979 to present, using an updated new version of the Goddard Earth Observing System Data Assimilation System Version 5 (GEOS-5).The GEOS-5 data assimilation system implements an Incremental Analysis Update (IAU) to slowly adjust the model estimates to the observations [11].MERRA achieves significant advances in the representation of the water cycle by using a catchment-based hydrological land surface model [11].A major difference between MERRA and other global reanalysis products is the spatial and temporal resolution of the archived data.The two-dimensional diagnostics (surface fluxes, single level meteorology, vertical integrals, and land states) are generated at 1-h intervals.These data and the 6-h 3-Dimensional atmospheric analyses are available at the full resolution of 0.5 ˝ˆ0.667 ˝.Extensive 3-h 3-Dimensional atmospheric diagnostics on 42 pressure levels are also provided by MERRA, but with coarser (1.25 ˝) resolution.The shortwave scheme used in MERRA is based on the scheme proposed by Chou and Suarez [44].

NCEP-NCAR
The NCEP-NCAR reanalysis [10,45] is a global reanalysis from 1948 almost to present, with a temporal resolution of 6-h, a spatial resolution of T62 Gaussian grids (1.9 ˝ˆ1.9 ˝), and 28 vertical layers.The land surface model is the Oregon State University (OSU) land surface model with two vertical layers.The analysis scheme is a 3-Dimensional variational (3D-Var) scheme cast in spectral space, which is used to assimilate the observations into the model.The observations include rawinsonde observations, satellite-derived temperature and winds from cloud-tracking, airborne temperature and wind measurements, and land and oceanic surface pressure measurements [45].The shortwave radiation parameterization method is based on that proposed by Lacis and Hansen [46].

NCEP-DOE
The NCEP-DOE reanalysis project [41] was established to provide an improved version of the original NCEP-NCAR reanalysis [10,45].The resolution of the NCEP-DOE model is the same as in the NCEP-NCAR Reanalysis project, T62 Gaussian grids (1.9 ˝ˆ1.9 ˝) with 28 vertical sigma levels.The temporal resolution is 6-h as in NCEP-NCAR.The main purpose of NCEP-DOE is to correct known errors (e.g., humidity diffusion and oceanic albedo) in the NCEP-NCAR project and also to update the parameterizations of the physical processes (e.g., shortwave radiation, radiation calculation, and cloud cover analysis).The NCEP-DOE is the "second generation" reanalysis, which focuses more strongly on accuracy, resolution, and long-term trends using an improved assimilation procedure based on 4D-Var.Improvements have been reported in hydrological cycles, shortwave radiation flux, etc. [41].The new shortwave radiation parameterization method by Chou [47] and Chou and Lee [48] replaced the Lacis and Hansen [46] parameterization method to mitigate the overestimation of surface solar radiation at the surface in NCEP-NCAR [42].

CFSR
The CFSR reanalysis [42] is the newest reanalysis dataset from NCEP, which was initially provided over the 31-year period from 1979 to 2009, and has been currently extended to March 2011.The atmospheric model used is the Global Forecast System (GFS) at a horizontal resolution of T382 Gaussian grids (0.3 ˝) with 64 vertical layers.Atmospheric analysis is performed at 6-h intervals and uses a grid point statistical interpolation (GSI) technique like MERRA [24].Shortwave radiation is parameterized according to the NASA approach [49,50].The parameterizations use random cloud overlap with shortwave radiation being called every hour.

JRA-55
The JMA conducted JRA-55 [39], the second Japanese global atmospheric reanalysis project.The JRA-55 covers 55 years, spanning from 1958 to 2013, coinciding with the establishment of the global radiosonde observation system.Compared to its predecessor, JRA-25, JRA-55 is based on a new data assimilation and prediction system that mitigates many deficiencies in the JRA-25 reanalysis [39].These improvements include a higher spatial resolution of T319 Gaussian grids (0.5625 ˝), a new radiation scheme, 4D-Var with Variational Bias Correction (VarBC) for satellite radiances, and introduction of greenhouse gases with time varying concentrations.The resulting products are considerably better than the corresponding JRA-25 products [39].Shortwave radiation fluxes are calculated by a two-stream method with delta-Eddington approximation [51].Computations for shortwave radiation are performed at 1 hour intervals rather than at every time step to cut down on the computing costs [39].

Ground Measurements
The surface observations used to evaluate R s estimates from the reanalyses were obtained from five databases: the Global Energy Balance Archive (GEBA, [34]), the Baseline Surface Radiation Network (BSRN, [35]), the Greenland Climate Network (GC-NET) [36], Buoy observation provided by the Tropical Atmosphere Ocean (TAO) project office of the National Oceanic and Atmospheric Administration (NOAA), and also a set of quality-controlled measurements [4] provided by CMA.R s measurements were originally collected from 1677 stations from GEBA (1551 without China), 58 stations from BSRN, 22 stations from GC-NET, 294 stations from NOAA, and 122 stations from CMA.
GEBA is a fundamental database, developed and maintained at ETH Zurich, of worldwide instrumentally measured energy fluxes at the surface [34].To date, GEBA provides 2500 stations with 450,000 monthly mean values of various surface energy balance components.The relative random error (root mean square errors) of R s at the surface in GEBA is 5% for monthly means and 2% for annual means [34].BSRN was established to provide radiation measurements of high quality and high temporal resolution at a limited number of sites around the world.At present, the BSRN project has archived more than 60 sites covering a wide latitudinal range from ´89.98 ˝to 82.49˝, a longitudinal range from ´156.61 ˝to 169.69 ˝, since January 1992.The GC-NET consists of more than 20 automatic weather stations (AWS) distributed primarily in the accumulation zones of ice sheet [36].GC-NET provides hourly radiation observations with instruments located between 0.1 and 5 m in height above the surface depending on local accumulation rates and tower stations.GC-NET instruments are factory calibrated.Because BSRN and GC-NET provide only the instantaneous values of R s , the daily integrated R s was obtained from the instantaneous values through a sinusoidal interpolation method [52].Critical quality control procedures were performed to estimate the monthly R s data from the instantaneous values of the surface measurements from BSRN and GC-NET.If the instantaneous values were missing for more than 20% of total observations during the daytime in one day, the daily integrated values were excluded.If daily R s data were missing for more than five days in one month, the monthly R s data were excluded from the comparison.Buoy stations from the TAO project office of NOAA are located over the tropical oceans covering a latitudinal range from ´10 ˝to 10 ˝.Unlike other networks, the Buoy observations are usually influenced by their cosine response, the thermal offset, and also wind and water waves [53].Please note that the monthly mean R s from Buoy observations was directly obtained from the TAO project office of NOAA without critical quality control.Daily and Monthly meteorological data at 122 CMA routine weather stations released by the CMA Meteorological Information Center were used in this study.The CMA surface observations were dubious due to technical failures and operating problems (e.g., instrument replacement) [54,55].The quality of the CMA surface observations used in this study was controlled based on the surface reconstructed data using routine meteorological data including air temperature, air pressure, relative humidity, sunshine duration, and precipitation [4,56].
Most of the reanalysis R s products are near real-time products, but the CFSR reanalysis product ended in 2009.Satellite retrieval from CERES began in March 2000.Therefore, evaluations in this study were conducted for the overlapping period from 2001 to 2009.For this time period, more than one year of R s data were required to be collected at each station.A total of 674 stations matching these data-availability requirements were used in this study, including 458 stations from GEBA, 44 stations from BSRN, 16 stations from GC-NET, 60 stations from Buoy observations, and 96 stations from CMA. Figure 1 shows the geographical distributions of the observation sites used in this study.

CERES_EBAF Product
Although several global surface radiation observation networks have been used in this study, many surfaces are under-represented, especially for the oceanic surfaces.The evaluation of the reanalysis R s datasets was performed at a limited number of stations.Therefore, it is still deficient on a global scale.Besides the ground measurements, satellite-derived R s data provide another effective way to validate the reanalysis products.The latest satellite retrieval by the CERES-EBAF has been reported to show higher accuracy than other gridded R s products [4, 33,53,57], because it utilized more accurate input data on clouds for calculating R s [5].Section 3.2 will show that the CERES-EBAF R s data compare more favorably with ground measurements than the R s reanalysis data.
The CERES-EBAF incorporates more accurate cloud information than that in reanalysis.The CERES-EBAF R s values were obtained using cloud and aerosol properties derived from instruments on the A-train Constellation (the Cloud Aerosol Lidar with Orthogonal Polarization (CALIOP) on the Cloud Aerosol Lidar and Infrared Pathfinder Satellite Observation (CALIPSO) satellite, the CloudSat Cloud Profiling Radar (CPR), and the Aqua Moderate Resolution Imaging Spectrometer (MODIS)) based on a radiative transfer model of CERES with k-distribution and correlated-k for radiation (FLCKKR) with a two stream approximation using the independent column approximation considering observational constraints of TOA irradiance from CERES [5].Besides the cloud and aerosol properties data from the A-train satellite, the other inputs included temperature and humidity profiles, ozone amounts, ocean spectral surface albedo, TOA albedo, and emissivities.The ocean heat uptake and CERES derived TOA irradiances from EBAF were used to constrain surface irradiance computations so that computed TOA irradiances would be consistent with CERES observations.The CERES-EBAF R s data set is a gap filled product with a monthly temporal resolution and 1 ˝spatial resolution from March 2000 to February 2015.

Validation Using Ground Measurements
In this study, the estimated R s values on the original grid scale were compared with ground-based observations within the grid cells.To perform the comparison, all of the selected global reanalysis R s products were projected onto a 1 ˝spatial resolution using a bilinear interpolation of a weighted average of pixels in the nearest 2-by-2 neighborhood to match that of the CERES-EBAF (Ed2.7)R s data.Please note that the spatial representativeness of ground measurements is one potential error source for R s evaluation.Monthly representation errors at the surface sites with respect to their 1 surroundings were on average 3.7% (4 W/m 2 ) [58].Four measures were used to characterize model performance: the correlation coefficients (R), root mean squared errors (RMSE), bias, and the ratio of the standard deviation of the reanalysis or the satellite product to that of the surface measurements (σ).This last is also a good indicator to test whether the reanalysis matches the surface measurements well.
Table 2 shows the validation results of the R s monthly values of the six reanalysis products by comparison with surface measurements at the selected 674 stations from 2001 to 2009.The statistical measures calculated using the various surface observation network sources as the evaluation data were shown in Table 2.It was apparent that all selected reanalysis R s products overestimated monthly R s when compared to the surface measurements from the five networks.The ERA-Interim monthly R s product had an overall R of 0.98, a positive bias of 4.15, an RMSE of 19.56 W/m 2 , and a ratio (σ) of 1.02 when the BSRN observations were used as validation data.When GEBA, Buoy, GC-NET, and CMA were used as validation data, the biases were 10.17 W/m 2 , ´2.86 W/m 2 , ´10.84 W/m 2 , and 23.15 W/m 2 , respectively.The ERA-Interim monthly R s product had an overall R of 0.95, a bias of 11.25 W/m 2 , an RMSE of 27.70 W/m 2 , and a ratio (σ) of 1.03.Because high correlation coefficients profit from the strong latitudinal dependencies of both surface-observed and reanalysis estimated R s , the correlation coefficients after removal of the R s seasonal cycle are also provided in Table 2.The latter coefficients were calculated to reduce the impacts of seasonal cycles on the quality of R s agreement between surface observations and the reanalysis product, except that both include a common seasonal cycle.The MERRA R s product overestimated R s compared to most of the surface observational networks, except for GC-NET.When BSRN, GEBA, Buoy, GC-NET, and CMA were used as validation data, the biases of MERRA monthly R s product were 13.80 W/m 2 , 22.17 W/m 2 , 3.60 W/m 2 , ´6.09 W/m 2 , and 34.38 W/m 2 , respectively.The MERRA monthly R s product had an overall R of 0.95, a bias of 22.36 W/m 2 , an RMSE of 34.13 W/m 2 , and a ratio (σ) of 1.05.
One of the main purposes of NCEP-DOE was to update the parameterizations of the physical processes in NCEP-NCAR including the shortwave radiation parameterization.Clearly, the NCEP-DOE R s product had a better accuracy than NCEP-NCAR compared to most of the surface observation networks, except for the Buoy observations.The biases of the NCEP-NCAR monthly R s product were 40.12 W/m 2 , 48.95 W/m 2 , 5.77 W/m 2 , 35.73 W/m 2 , and 71.95 W/m 2 for the BSRN, GEBA, Buoy, GC-NET, and CMA observations, whereas these values were 20.48 W/m 2 , 33.67 W/m 2 , ´20.14 W/m 2 , 12.20 W/m 2 , and 45.75 W/m 2 for NCEP-DOE.The NCEP-DOE monthly R s product had an overall R of 0.90, a bias of 30.83W/m 2 , an RMSE of 46.82 W/m 2 , and a ratio (σ) of 1.06.The CFSR showed biases of 8.97 W/m 2 , 13.99 W/m 2 , 30.91 W/m 2 , 5.84 W/m 2 , and 30.69W/m 2 compared to the BSRN, GEBA, Buoy, GC-NET, and CMA observations, respectively.The overall R, bias, RMSE, and ratio (σ) of the CFSR monthly R s product was 0.95, 18.57W/m 2 , 32.23 W/m 2 , and 1.10.Similarly to NCEP-DOE, the JRA-55 monthly R s product overestimated R s compared to most of the surface observational networks, except for the Buoy observations.The biases were 12.46 W/m 2 , 20.06 W/m 2 , -15.58 W/m 2 , 4.75 W/m 2 , and 32.23 W/m 2 for the BSRN, GEBA, Buoy, GC-NET, and CMA observations, respectively.The overall R, bias, RMSE, and ratio (σ) of the JRA-55 monthly R s product were 0.95, 19.15 W/m 2 , 31.71W/m 2 , and 1.00, respectively.
The reanalysis monthly mean R s products were divided into summer (JJA) and winter (DJF) seasons to assess the seasonal dependency of their accuracy.Table 3 displays the seasonal validation results for the monthly R s of the six reanalysis data sets using surface measurements collected from 674 stations from 2001 to 2009.It was apparent that almost all the reanalysis Rs products showed better accuracy in DJF than that in JJA, except for the ERA-Interim Rs product at the CMA and BSRN stations, for the MERRA, JAR-55, and NCEP-DOE Rs products at the Buoy stations, and for the JAR-55 at the GC-net stations.They showed higher bias in DJF, but lower bias in JJA.
Furthermore, the evaluations were also conducted over different regions to address the spatial representativeness issue because the measuring instruments and measurement methods change by country more than by climate regime [53].The whole world was divided into nine sub-regions, which include North American, South American, Asia, Europe, Africa, Australia, Oceania, Antarctic, and the oceans.Table 4 summarizes the statistical measures for the evaluation results of the monthly R s products from the reanalysis data using surface measurements at the selected 674 stations from 2001 to 2009.Relatively small discrepancies were found for the ERA-Interim R s product.The biases of the ERA-Interim monthly R s were 2.32 W/m 2 , 4.64 W/m 2 , 2.06 W/m 2 , 4.71 W/m 2 , and ´2.20 W/m 2 for North America, Europe, Oceania, and Antarctica, whereas these values were 20.62 W/m 2 , 12.64 W/m 2 , 13.79 W/m 2 , and 10.92 W/m 2 for Asia, Africa, South America, and the whole of the oceans, respectively.According to the validation results shown in Tables 2-4 we concluded that the limited number of surface observations cannot entirely represent the spatial variability of R s globally.
To investigate the latitudinal dependencies of the monthly R s from the reanalysis products, the biases for the monthly R s products were calculated as a function of latitude and are shown in Figure 2. In this figure, the calculated biases were averages over the biases of the monthly reanalysis R s at sites located within latitudinal bands of 5 ˝.Most of the reanalysis R s products showed a pronounced latitudinal dependence in their biases, except for the ERA-Interim and JRA-55 products.The averaged bias weighted by the area of the latitudinal band was also calculated for each reanalysis R s product.The averaged bias values weighted by the area of each latitudinal band were 5.32 W/m 2 , 11.25 W/m 2 , 21.62 W/m 2 , 9.43 W/m 2 , 9.95 W/m 2 , and 7.02 W/m 2 for ERA-Interim, MERRA, NCEP-NCAR, NCEP-DOE, CFSR, and JRA-55, respectively.The biases obtained in this way were smaller than those of the reanalysis R s products, as shown in Table 2, which indicates that the mean model biases at the 674 sites are effected by the which they are determined and which further confirms that a limited number of surface observations cannot represent the spatial variability in the R s evaluation results.

Comparison with CERES-EBAF
Although the R s values from the reanalysis products have been evaluated using 674 sites from various observation networks around the world, many land and oceanic surfaces are still under-represented.Moreover, the reanalysis systems differ in their model structure, physical parameterization, spatial resolution and vertical resolution, and in the data assimilation methods.It is impossible to discuss in detail all the differences among the reanalysis methods in the context of this paper.Rather, to address this issue, a broad view of the seasonal time scales of the differences will be provided by comparing the reanalysis R s data with independent data from satellite retrieval of R s .
The satellite retrieval of R s from the CERES-EBAF product was used in this study.Firstly, the CERES-EBAF monthly R s data were validated against the surface observed R s at the 674 sites, as shown in Figure 3.The CERES-EBAF R s product correlates well with the surface measurements by the small biases, and RMSE was found at most sites for the CERES-EBAF R s product.The absolute error was less than 10 W/m 2 at 484 out of 674 sites and the RMSE was less than 20 W/m 2 at 546 out of 674 sites.Table 5 shows that the CERES-EBAF R s product agreed better with surface observations than that of the reanalysis products, as shown in Tables 2 and 3.The biases of the CERES-EBAF R s product were ´1.28 W/m 2 , 2.10 W/m 2 , 4.39 W/m 2 , 8.77 W/m 2 , and 7.31 W/m 2 for the BSRN, GEBA, Buoy, GC-NET, and CMA observations, respectively.The seasonal validation results for the CERES-EBAF R s product are also shown in Table 5.Clearly, the biases and the RMSE have been greatly reduced compared to the values calculated from the reanalysis products, whereas the correlation coefficients increased, especially for the Buoy, GC-NET, and CMA observations.For example, the CERES-EBAF had a bias of 7.31 W/m 2 at the CMA sites, whereas the biases of the six reanalysis R s products ranged from 23.15 W/m 2 to 71.95 W/m 2 .Previous studies have pointed out that the reliability of data obtained in China has often been dubious because of the technical failures and operating problems [54,55].However, this study found that the bias and RMSE of the CERES-EBAF R s product were 11.36 W/m 2 and 19.38 W/m 2 , respectively, in DJF using the quality-controlled CMA surface observations as the validation data, whereas these values were 1.30 W/m 2 and 18.04 W/m 2 in JJA.Generally, the validation results should be worse than those derived in JJA when heavy clouds occur.Clouds and aerosols are two of the most important factors in R s estimation.The only data available at the CMA sites were the climatologies of clouds (cloud fraction (CF)) and aerosols (<10 km visibility frequencies).It was found that the frequencies of visibility lower than 10 km in DJF were 20% higher than in JJA.Therefore, it can be concluded that the quality-controlled CMA R s observations are good enough to be used to validate the R s products, either from reanalysis or from satellite observations.It can also be concluded that the differences between R s products and surface measurements might have been related to aerosols, clouds and their interactions as represented in the algorithm, especially in DJF when heavy air pollution occurs.The evaluation results for the CERES-EBAF R s product in the various sub-regions are shown in Table 6 and further confirm that the CERES-EBAF R s product achieves high precision and is therefore qualified to represent surface observations on a large scale.Therefore, CERES-EBAF can be used to validate the R s obtained from the six reanalysis products.Before comparison using CERES-EBAF R s data, the obtained six global reanalysis R s products under study were projected onto a 1 ˝spatial resolution using the bilinear interpolation of a weighted average of pixels in the nearest 2-by-2 neighborhood to match that of the CERES-EBAF R s data.Figures 4 and 5 display the CERES-EBAF monthly R s climatologies for DJF (December, January, and February) and JJA (June, July, and August) for the nine-year period (2001)(2002)(2003)(2004)(2005)(2006)(2007)(2008)(2009), and also the biases between the reanalyzed R s product and the CERES-EBAF R s product (reanalysis minus CERES MODIS) for DJF and JJA.As shown in Figure 4, the CERES-EBAF R s climatology for DJF indicates that greater R s was received in the latitudinal bands between 0 ˝S and 30 ˝S, and also the Antarctic in DJF, whereas smaller R s was received in the Northern hemisphere.It was also shown that the greater R s values were received in the Northern Pacific, the Northern Atlantic, the tropical oceans, and also in Greenland in JJA.Solar transmission is determined primarily by clouds and aerosols (which were not available for comparison in this study).Previous studies have reported that the discrepancies in the reanalysis CF products might be reflected in the surface radiation fields [59,60].In this study, the CERES CF climatologies for DJF and JJA, and the biases between the reanalyzed CF and the CERES CF (reanalysis minus CERES) from 2001 to 2009 were also examined in this study, as shown in Figures 6 and 7.The observed CF from CERES indicated that clouds occurred more frequently in the Southern Ocean, the Northern Pacific, and the Northern Atlantic Oceans, and the Intertropical Convergence Zone (ITCZ) in both JJA and DJF.However, clouds occurred less frequently over the central areas of the Pacific and Atlantic Oceans, the Saharan Desert, Western Asia, Australia, Southwestern America, and in DJF.Clouds occurred more frequently in the Arctic and Antarctic in DJF than that in JJA.From the combination of Figures 4-7 it is clear that the high values of R s were accompanied by lower CF in most of the regions within the same latitudinal bands (regardless of the seasonal cycles).According to Figures 4 and 6 the six reanalyses overestimated R s in the Southern Oceans and Eastern China in DJF, which were related to the underestimation of CF in DJF.The overestimation of NCEP-NCAR on most of the continents and the underestimation of NCEP-DOE in tropical oceans in DJF were related to the underestimation and overestimation of CF in the corresponding regions, respectively.Similar results were also found in JJA, with the exception of China.Although almost all of the six reanalysis products overestimated R s in China, the CF values from the reanalysis products were slightly overestimated, except for NCEP-DOE and JRA-55.Previous studies have shown that aerosols might be another potential source of the larger positive biases in China [4, 26,53], where heavy pollution is occurring due to rapid economic development and high population density.Similarly, the interactions between clouds and aerosols may amplify this overestimation in China.However, because quantitative clouds and aerosols measurements are not available, it is impossible to determine the quantitative contributions of clouds and aerosols to the larger positive biases of R s in China.In any case, efforts are needed to improve the cloud schemes for different regions and seasons in reanalysis estimates.The gridded comparison results for the monthly R s from the six global reanalysis products using the CERES-EBAF R s are shown in Table 7.The biases of the reanalysis R s products ranged from ´2.98 W/m 2 to 21.97 W/m 2 over the globe, ´3.26 W/m 2 to 32.13 W/m 2 over land, and ´2.85 W/m 2 to 16.99 W/m 2 over the oceans.Almost all of the reanalysis overestimated the monthly R s over land, the oceans, and the globe, except for ERA-Interim which slightly underestimated R s .These results were different from the validation results presented in Section 3.1.The biases and RMSE values of the global reanalyses compared with the CERES-EBAF R s data were smaller than those between the global reanalyses and surface observations.For instance, the ERA-Interim had a bias of 11.25 W/m 2 using surface observations at 674 sites as the validation data, but the bias was only ´2.85 W/m 2 when compared with CERES-EBAF R s data.Overall, most of the reanalyzed R s were overestimated, although in some regions, R s was underestimated, and the differences between the surface observations and the reanalyses were sometimes greater than 100 W/m 2 .The gridded comparison results for the monthly mean CF estimated by the six global reanalysis products using the CERES CF as the validation data are also summarized in Table 8.The biases of the reanalysis CF products ranged from ´14.17 (%) to 2.42 (%) over the globe, ´15.61 (%) to 3.94 (%) over land, and ´13.47 (%) to 1.93 (%) over the oceans, as shown in Table 8. Figure 8 shows the multi-year annual mean CF and R s obtained from the global reanalyses and CERES data from 2001 to 2009 over land, oceans, and the globe.From the combination of Tables 7 and 8 and Figure 8, it is apparent that most of the reanalysis CF data underestimated over land, oceans, and the globe compared with the CERES MODIS observed CF data, whereas most of the reanalysis R s data overestimated over land, oceans, and the globe compared with the CERES-EBAF R s data.These results further proved that the discrepancies in the reanalysis CF products were reflected in the surface radiation fields, and that the clouds scheme should be improved in the reanalysis estimates.The long-term global annual variability of R s from the global reanalyses were also examined using CERES-EBAF R s data.The global reanalysis R s data showed quite different long-term global annual variability pattern from those derived from the CERES-EBAF data.As shown in Figure 10, the CF anomalies of global reanalyses showed a highly negative correlation with that of R s from the global reanalyses.This means that the long-term global annual variability of R s was primarily determined by the CF.Since the cloud scheme may have introduced large uncertainties to the estimates of R s in reanalysis, the reanalysis R s cannot be used to determine long-term variability in R s .

Global Annual Mean
As is well known, R s is an essential component of the total energy exchange between the atmosphere and the Earth's surface.over oceans, respectively.The annual mean R s calculated from the nine-year CERES-EBAF data were 186.7 W/m 2 , 188.2 W/m 2 , and 186.1 W/m 2 over the globe, land, and the oceans, respectively.It is obvious that the multi-year annual mean R s calculated from the nine-year global reanalysis data are greater than those derived from CERES-EBAF satellite retrievals except for ERA-Interim.Previous studies have always calculated the annual mean after removing the bias between the obtained R s from gridded products and surface observations.However, the R s obtained from gridded products may show strong latitudinal dependencies in their biases compared with the surface observations, which may introduce large uncertainties into the annual mean R s estimates.In this study, a bias, which is the average of biases in the latitude bands (5o latitudinal bands) weighted by the area of the latitudinal band, was also calculated.The averaged biases weighted by the area of the latitudinal band were 3.6 W/m 2 , 6.1 W/m 2 , and 3.8 W/m 2 over land, the oceans, and the globe, whereas these values were 3.0 W/m 2 , 4.2 W/m 2 , and 3.3 W/m 2 (Tables 5 and 6) using the surface observations as the validation data directly.It is clear that the biases determined in this way are different from those calculated directly using the surface observations.The global annual mean R s values of 184.6 W/m 2 , 180.0 W/m 2 , and 182.9 W/m 2 were obtained after removing the averaged bias weighted by the area of the latitudinal band over land, the oceans, and the globe.Note that the averaged bias weighted by the area of the latitudinal band over the oceans was greater than those over land and the globe, because the surface observations over the oceans are limited and most of the oceanic surfaces are under-represented in the high latitude regions such as the Northern Pacific Ocean and the Southern Ocean.Therefore, the annual mean R s over the oceans may contain larger uncertainties than those derived for land.The global land mean R s of 184.6 W/m 2 obtained in this study is very close to the land mean R s of 184 W/m 2 obtained by Trenberth and Fasullo [64] and Wild [57]

Discussion
The biases of current global reanalyses at the CMA sites ranged from 23.15 W/m 2 to 71.95 W/m 2 , whereas the reanalyses significantly underestimated CF compared to the CERES MODIS observed CF over China, as shown in Figures 6 and 7.Besides clouds, an underestimation of aerosols may also have introduced an overestimation of R s in China.Previous studies have pointed out that the reliability of data obtained in China has often been dubious because of technical failures, instrument replacement, and other factors.The quality of surface measurements from CMA was controlled in this study based on reconstruction using routine meteorological data maintained by CMA.CERES-EBAF had a smaller bias, especially in JJA, with a bias of 1.30 W/m 2 and an RMSE of 18.04 W/m 2 compared to those values derived from global reanalyses.The bias and RMSE of CERES-EBAF were 11.36 W/m 2 and 19.38 W/m 2 in DJF, when heavy air pollution occurred due to intensive use of coal.Therefore, it can be concluded that clouds and aerosol representation schemes over China, especially for the winter season (DJF), need to be further improved for both satellite retrievals and reanalysis estimates.
It is clear that aerosols and cloud properties represented by the global reanalyses R s models are two important factors in regulating the estimated R s .Besides cloud and aerosols, measurement errors (instrument sensitivity, drift, and urbanization effects) and spatial representativeness of surface measurements may cause potential errors in R s evaluation.Measurement errors may result in relatively low values of surface measurements and amplify overestimation.Monthly representation errors at the surface sites with respect to their 1 ˝surroundings were on average 3.7% (4 W/m 2 ) [58].
The global reanalysis R s and CF products were also evaluated by comparison with independent data from the satellite retrieval of CERES-EBAF in this study.It was found that the overestimation and the underestimation of the reanalysis R s were related to the underestimation and overestimation of CF in the reanalysis products in the corresponding regions, respectively.The gridded comparison results for the monthly R s from the six global reanalysis products using the CERES-EBAF R s as the validation data are also shown in this study (Table 7).The biases of the reanalysis R s products ranged from ´2.98 W/m 2 to 21.97 W/m 2 over the globe, ´3.26 W/m 2 to 32.13 W/m 2 over land, and ´2.85 W/m 2 to 16.99 W/m 2 over the oceans.Almost all of the monthly reanalysis R s data were overestimated over land, the oceans, and the globe except for ERA-Interim, which slightly underestimated R s .The reanalysis R s data showed similar seasonal variability to the CERES-EBAF R s data, but with a relatively high value for global monthly mean compared to the CERES-EBAF R s data.The reanalysis R s data exhibited quite different global long-term variation trends compared to CERES-EBAF.The global long-term variability in reanalysis R s data was primarily determined by the CF estimates in global reanalyses.
The global annual mean R s values calculated directly from the nine-year CERES-EBAF data were 186.7 W/m 2 , 188.2 W/m 2 , and 186.1 W/m 2 over the globe, land, and the oceans, respectively.After removing the averaged bias weighted by the area of the latitudinal band, the global annual mean estimates of R s were 184.6 W/m 2 , 180.0 W/m 2 , and 182.9 W/m 2 over land, the oceans, and the globe over the nine-year period from 2001 to 2009, respectively.However, the annual global annual mean R s for the oceans and the globe may contain larger uncertainties than those derived for land because the surface observations over the oceans are limited and most of the oceanic surfaces are under-represented in the high latitude regions.The global land mean R s of 184.6 W/m 2 obtained in this study is very close to the land mean R s of 184 W/m 2 obtained by Trenberth and Fasullo [64] and Wild [57].

Conclusions
The current study presents the validation and the inter-comparison of R s estimates provided by the current representative from six current representative global reanalyses (NCEP-NCAR, NCEP-DOE, CFSR; ERA-Interim; MERRA; and JRA-55) using the quality-controlled surface measurements at 674 sites from five different observational networks including GEBA, BSRN, GC-NET, Buoy, and CMA of 2001-2009.The reanalysis products of R s were also compared with the satellite retrievals of CERES-EBAF, including accuracy, spatial distribution, seasonal variation, and global annual mean values of R s .
Overall, the biases between the reanalysis products of R s and surface measurements at all sites ranged from 11.25 W/m 2 to 49.80 W/m 2 .Current global reanalyses overestimated R s compared to most of the ground measurements, especially for the CMA sites in China.It was also found that the biases of cloud fraction (CF) in the reanalyses caused the overestimation of R s .After removing the averaged bias of CERES-EBAF, weighted by the area of the latitudinal band, a global annual mean R R s values of 184.6 W/m 2 , 180.0 W/m 2 , and 182.9 W/m 2 were obtained over land, ocean, and the globe, respectively.Generally, further improvements and efforts are required for global energy budget studies, in terms of collecting high accuracy surface measurements, improving R s estimates both from satellite retrievals or reanalysis products, and refining clouds and aerosols schemes in their physical models.

Figure 1 .
Figure 1.Geographical distribution of observation sites (674 sites in total) used in this study from GEBA (458 sites in blue asterisks), from BSRN (44 sites in red triangles), from Buoy observations (60 sites in green circles), from GC-NET (16 sites in green squares), and from CMA (96 sites in black asterisks).

Figure 2 .
Figure 2. Monthly R s biases at the 674 sites as a function of latitude, for six different reanalysis and one satellite-derived (CERES-EBAF) R s products.Biases are averaged over sites within 5 ˝latitudinal bands from 2001 to 2009.Units are W/m 2 .

Figure 3 .Table 5 .
Figure 3. Spatial distributions of model biases (a) and RMSE (b) between surface observed R s at 674 sites and R s obtained from CERES-EBAF R s products.Units are W/m 2 .

Figure
Figure9a-cshow the seasonal variations in calculated surface R s from CERES-EBAF and the reanalysis R s data from 2001 to 2009 for land, oceans, and the globe.As shown in the figure, although the reanalysis R s data showed similar seasonal variability to the CERES-EBAF R s data, the reanalysis R s data, with the exception of the ERA-Interim, exhibited a higher global monthly mean value compared to the CERES-EBAF R s data, especially in the case of the NCEP-NCAR data, with the maximum values occurring in JJA.The ERA-Interim showed a slightly lower value compared to the CERES-EBAF R s data at a global scale.

Figure 8 .
Figure 8. Multi-year annual mean cloud fraction (CF) and R s obtained from the global reanalyses and CERES data from 2001 to 2009: (a) CF for land; (b) R s for land; (c) CF for oceans; (d) R s for oceans; (e) CF for the global; (f) R s for the globe.The numbers in the brackets and the red lines indicate the multi-year annual mean CF and R s calculated from CERES data.

Figure 9 .
Figure 9. Seasonal variability in R s calculated from the global reanalyses and CERES-EBAF R s data from 2001 to 2009: (a) land; (b) oceans; (c) the globe.

Figure 10 .
Figure 10.Global long-term variability of cloud fraction (CF) and R s obtained from the global reanalyses and CERES: (a) ERA-Interim; (b) MERRA; (c) NCEP-NCAR; (d) NCEP-DOE; (e) CFSR; (f) JRA-55; (g) CERES.3.3.Global Annual Mean Much effort has been expended to estimate the global annual mean values from various data sources including surface observations, satellite retrievals, reanalysis products, and global climate simulations [2,4,5,53,57,61-64].Figures 8 and 11 display the annual mean R s calculated directly from the CERES-EBAF and the reanalysis R s data for 2001-2009.The multi-year annual mean R s calculated from the nine-year global reanalysis data ranged from 184.8 W/m 2 to 205.9 W/m 2 over the globe, 185.8 W/m 2 to 223.2 W/m 2 over land, and 184.4 W/m 2 to 198.9 W/m 2

Figure 11 .
Figure 11.Box plots of the annual mean R s from the global reanalyses and CERES-EBAF data from 2001 to 2009.The red and black lines within the box indicate the mean and median annual mean R s values, respectively: (a) land; (b) oceans; (c) the globe.

Table 1 .
Details of the data products used in this comparison.

Table 2 .
Evaluation of the monthly R s from six reanalysis data sets using surface measurements collected from 674 stations from 2001 to 2009.The numbers in the brackets are the correlation coefficients after removing of the seasonal cycle.Units are W/m 2 for bias and RMSE.

Table 3 .
Seasonal validation results for the monthly R s of six reanalysis data sets using surface measurements collected from 674 stations from 2001 to 2009.DJF represents December, January, and February, whereas JJA represents June, July, and August.Units are W/m 2 for bias and RMSE.

Table 4 .
Evaluation of the monthly R s from six reanalysis data sets using surface measurements in various sub-regions from 2001 to 2009.The numbers of observational stations are shown in brackets in the first column.Units are W/m 2 for Bias and RMSE.

Table 6 .
Evaluation of the monthly R s for CERES-EBAF data using surface measurements in various sub-regions from 2001 to 2009.The numbers of observational stations are shown in brackets in the first column.The numbers in the brackets in the second column are the correlation coefficients after removing the seasonal cycle.Units are W/m 2 for bias and RMSE.

Table 7 .
Evaluation of monthly R s global reanalyses using CERES-EBAF data from 2001 to 2009 at 1 ˝ˆ1 ˝grid scale over land, the oceans and the globe.Units are W/m 2 for bias and RMSE.

Table 8 .
Evaluation of monthly cloud fraction (CF) from global reanalyses using CERES MODIS observed CF data from 2001 to 2009 at 1 ˝ˆ1 ˝grid scale over land, the oceans, and the globe.