Minimizing heterogeneity of VOS observations for improving estimates of short-wave radiation over the oceans

The paper describes the development of climatologies of short-wave radiation over the ocean using visually observed cloud cover from the International Comprehensive Ocean-Atmosphere Data Set (ICOADS). For producing climatology of short-wave radiation, we applied bulk parameterization (known as SAIL/IOARAS). For a proper computation of short-wave radiation values special effort was taken to reduce sampling errors in cloud cover from ICOADS. Specifically, we consider the uncertainties associated with temporal inhomogeneity of observations and potentially implying biases in short-wave radiation estimates related to the astronomy factors critical for short-wave radiation. For this purpose, we propose simple and effective algorithms.


Introduction
Heat fluxes between the ocean and the atmosphere play a key role in the understanding of global climate changes and the ocean's role in climate variability. Currently the main source of long-term information about ocean-atmosphere heat fluxes are the re-analyses, satellite data, and estimates based on the use of bulk parameterizations applied to voluntary observing ship (VOS) data. Heat fluxes from re-analyses are characterized by global coverage and rather high resolution, but they are strongly dependent on model settings and parameterizations embedded into models. This is especially critical for cloud cover, which is a key parameter for computing short wave radiation. Despite the progress in representation of cloud cover in reanalyses, the quality of reanalysis cloud cover products is still low [1]. Satellite data are quite accurate and relatively homogenous in time and in space, however they are still short in time being very limited prior 1980s and, thus limiting capabilities of development climate quality long term time series [2]. In this respect computations based on VOS data using simplified bulk parameterizations have an advantage of covering multidecadal time periods, thus providing information about long-term changes in surface radiation.
Recent decades were characterized by substantial progress in the development of methods for computations of heat fluxes at the ocean surface. However, the main efforts were assigned to the improvement of the parameterization of sensible and latent turbulent heat fluxes and improvement of the accuracy of their computations which was considered to be low compared to radiative fluxes. As a result, developed during 2000s COARE-3 parameterization and its updates now provides computations of turbulent fluxes with an accuracy of 5 W/m 2 for most conditions [3]. At the same time, the development of parameterizations of radiative fluxes and specifically of shortwave radiation was not very active after the period 1980s and 1990s when some advanced parameterizations based on the use  [4]. However, the accuracy of radiative fluxes computed using these parameterizations does not meet already the requirements of climate science.
Global and basin-scale estimates of the ocean -atmosphere fluxes [5] are characterized by systematic and random errors, locally exceeding 30 W/m 2 . Given the relatively high accuracy of the modern estimation of turbulent fluxes and the progress in minimization of the impact of sampling errors [6], radiative fluxes become the most likely source of inaccuracies.
Empirical parameterizations of short-wave radiation are based on the solar altitude and the cloud cover. Solar altitude is based on astronomical parameters and cloud cover enters parameterizations to account for atmospheric transmission factor. The main source of cloud data over the oceans is the International Comprehensive Ocean-Atmosphere Data Set (ICOADS) [7]. ICOADS is a set of historical visual observations made by Voluntary Observing Ships (VOS). This data set contains a wide range of meteorological variables. Cloud cover in ICOADS is given as a visually observed total and low cloud cover and morphological types of cloud observed by marine officers. The main problem for using ICOADS cloud data is spatially and temporally inhomogeneous sampling. For a proper computation of radiative fluxes it is necessary to develop a scheme minimizing sampling errors. In this paper we focus on reducing sampling errors in cloud data. Computations are performed for the boreal summer (JAS). This choice is justified by our focus on the Northern Hemisphere oceans, given that the number of ICOADS reports in the Southern Ocean (south of 40S) is strongly limited. All calculations were carried out for the period 1950-2017 for individual summer seasons. Gridding is performed using averaging of reports for 5-degree grid cells.

SAIL parameterization of incoming short-wave radiation over the Ocean
In the late 1990s, a significant update of the ICOADS data set [7] was carried out. As a result of this update, not only the total cloud cover amount but also the low cloud cover amount and morphological cloud types were included in VOS reports assimilated in ICOADS. This allows for engaging new cloud characteristics in the empirical parameterization of the short-wave radiation. Existing parameterizations use only a total cloud cover amount [4] or total and low cloud cover amount. In order to use the new information about clouds from the ICOADS we developed new parameterization of short-wave radiation over the ocean (SAIL parameterization herein) [8].
This parameterization is based on 4 years of highly accurate in-situ ship observations in different regions of the Atlantic Ocean. Incoming short-wave radiation Q at the ocean surface can be expressed as: where Q0 is the solar constant set to 1368×sinh, W/m 2 , h is solar altitude, ai and bi are the empirical coefficients dependent on total cloud amount or morphological cloud types under the overcast or close to overcast conditions (7-8 octas of total cloud cover). The use of the logarithmic function instead of the linear in (1) (which is different from e.g. [4]) allows for highly accurate approximation of the dependence of the transmission factor on solar altitude under small solar declinations. As the other parameterizations SAIL scheme provides maximum accuracy under clear skies and a small cloud amounts, and the lowest accuracy under overcast conditions. To improve the accuracy for 7-8 octas, SAIL scheme engages cloud categories based on morphological types of cloudiness.
3. Minimizing the errors associated with temporal heterogeneity SAIL parameterization of short-wave radiation was applied to the cloud data from ICOADS data set on individual sampling basis. In doing this, we suggest a method for reducing errors associated with temporal data inhomogeneity.
The temporal inhomogeneity results from the intra-daily and intra-seasonal inhomogeneity. To solve the problem of inhomogeneous distribution of data during the day we apply the "virtual clockwise rotation" implying computation of daily mean radiation for each ICOADS report: where Qcc is the daily mean radiation, hi is solar altitude, computed for each of 48 steps using latitude and longitude reported in ICOADS with a time step of 30 minutes, a and b are empirical coefficients corresponding to cloud conditions. This approach can lead to a slight error associated with the diurnal variation of cloud cover, but this error is much less compared to the error related to the inhomogeneity of data.   Figures 1a and 1b present the climatologies of short-wave radiation, calculated with and without the "virtual clockwise rotation". The climatology of short-wave radiation calculated without applying the "virtual clockwise rotation" (Figure 1a), does not reflect the realistic spatial pattern of short-wave radiation being largely affected by biases in astronomy associated with inhomogeneous distribution of the reports over day. When we use the "virtual clockwise rotation", climatology of short-wave radiation becomes very reasonable ( figure 1b). The maximum values of short-wave radiation (250-280 W/m 2 ) in JAS are observed in the Northern Hemisphere tropics and gradually decrease towards the north and the south. Figure 1c shows the differences in the climatology of short-wave radiation between the short-wave radiation calculated directly from the ICOADS reports and climatology calculated using "virtual clockwise rotation". These differences are positive over most of the Atlantic Ocean and the eastern Pacific. In the Indian Ocean, especially in the eastern part, the differences are negative.    as reported by ICOADS observations of cloud cover for 5-degree grid cells. Since the solar altitude depends exclusively on astronomical factors its distribution should be strictly zonal, which is not the case in Figure 2a. Figure  2b shows the averaged solar altitude as provided by the calculations using "virtual clockwise rotation". Here, the spatial pattern takes the almost zonal distribution. A small deviation from the latitudinal distribution is explained by the fact that the "virtual clockwise rotation" excludes only heterogeneous distribution during the day, but there will remain heterogeneity on the intra-seasonal scale.
The maximum differences between short-wave radiation, calculated with and without "virtual clockwise rotation" are observed in the regions where the strongest biases in the computed averaged   There are several ways to solve the problem of heterogeneous sampling within a season. The easiest way is to attribute all the calculations to the central date of the season (August 16 for JAS). The climatology of short-wave radiation calculated in this way is shown in Figure 3a. If calculations are carried out for the central date of the season for each report, a realistic zonally distribution of short-wave radiation is preserved. Differences between the climatology calculated using dates from each ICOADS report and climatology attributed to the central day of the season are shown in Figure 3b. The differences are negative in the northern parts of the Atlantic and Pacific Oceans (with the strongest deviations being 7-15 W/m 2 ) and are close to zero over the remaining part of the ocean.
Differences in figure 3a are close to zero in the tropics that can be explained by very minor impact of biases in solar altitude near the Equator (Figure 4). In the Southern Hemisphere mid-latitudes due to the low solar altitude in JAS, the values of the incoming short-wave radiation are relatively small. This partly explains also small biases in the computed short-wave radiation here. Negative differences in the mid-latitudes of the Northern Hemisphere are due to the mismatch between the central season date and the date of the fall equinox. In the first half of the season, the decrease in the average solar altitude and the duration of the day is slower than at the end of the season. Thus, when using approach based on the central season date for calculations (August 16), the contribution from the first half of the season is smaller compared to the second half. Positive differences near 40S should be considered with caution, due to poor sampling in this region. More effective way to solve the problem of heterogeneity of data distribution within the season is to apply the computation for all days of the individual season for a given report, i.e. using actual cloud cover information and changing the dates. Conceptually, this is similar to the "virtual clockwise rotation" approach for minimizing biases associated with uneven distribution of the reports over the day. Resulting this case, the mean seasonal short-wave radiation is computed as the average over all days during season (92 values for JAS). This may introduce errors due to the intra-seasonal variation of cloud cover. However, these errors are smaller than those associated with the heterogeneous data distribution or non-linearity of solar altitude within the season. Climatology of short-wave radiation calculated using this approach is shown in Figure 3c. Differences between the climatologies of short-wave radiation computed for the central date of the season and with "running" over all days approach (Figure 3d) are the positive between 30N and 45N (6-8 W/m 2 ). Negative differences reach their strongest values south of 40S (also 6-8 W/m 2 ), while in the tropics of the Southern Hemisphere (between the Equator and 20S), these differences are close to zero. The nature of such differences can be explained based on figure 4, which shows the changes in the averaged daily solar altitude over the season for different latitudes. In the tropics of the Southern Hemisphere, the end of the season is associated with a rapid increase in the averaged solar altitude and the duration of the day, while in the Northern Hemisphere the solar altitude and the duration of the day decrease. The fastest changes in both the solar altitude and the duration of the day are observed on dates close to the fall equinox. Thus, it is the second half of the season (JAS) that makes the greatest contribution to the averaged seasonal flux of the incoming short-wave radiation. This leads to decreasing values computed with "running" for all days approach. Near the equator, where the changes in the solar altitude are small, these differences are close to zero.
Summarizing, we note that the use of "virtual clockwise rotation" avoids the temporal inhomogeneity associated with the heterogeneous distribution of VOS reports during the day. Approach based upon "running" overall days of the individual season computations reduces the errors associated with the heterogeneous distribution of reports within the season and takes into account the nonlinearity of changes in the solar altitude during the season. As a result, calculations this approach result in increase in the mean seasonal value of short-wave radiation by 7-12 W/m 2 in the North Atlantic and the North Pacific and 3-8 W/m 2 south of 20S in JAS compared to climatology that takes into account only "virtual clockwise rotation". These differences are quite large, given that the accuracy of computing turbulent fluxes using the COARE-3 parameterization being around 5 W/m 2 .

Minimizing the errors associated with spatial heterogeneity
Another important problem in the ICOADS data is a spatial inhomogeneity of sampling. Most of the ICOADS reports are concentrated along the main ship routes in the North Atlantic and the North Pacific, while in the Southern Ocean there are regions with almost no data. So, it is necessary to find an approach to minimizing sampling errors.
To minimize errors associated with spatial inhomogeneity of ICOADS reports we used different grid cells in different regions of the ocean. In most regions, we used 5×5˚ cells, however south of 40S the size of grid cells was extended to 10×20˚. As a justification for this we can note large spatial correlation of short-wave radiation, especially in latitudinal direction, as for the other atmospheric parameters in this area [9]. Also, due to the proximity to the pole, grid cells in high latitudes are considerably smaller than those in mid-and low latitudes. Nevertheless, even in this case, there are cells without data. To provide gap-filling we used objective interpolation scheme [10] and developed the seasonally mean fields of short-wave radiation fluxes over the ocean.

Discussion and conclusions
We proposed algorithms for reducing errors in computed incoming short-wave radiation at the sea surface using ICOADS archive and bulk parameterizations. It is shown that to reduce errors associated with temporal heterogeneity of data distribution over the day "virtual clockwise rotation" method should be used. To reduce the errors associated with the heterogeneity of the data distribution within the season, we used the "running' computations over all days of the season, which allows us account for the nonlinearity of changes in the solar altitude. It is shown that the use of the "running" over all days computations leads to an increase in the flux of incoming short-wave radiation by 7-12 W/m 2 (in the North Atlantic and the North Pacific).
For minimizing biases associated with spatial heterogeneity of the ICOADS data, we suggested to use variable grid cell size in conjunction with objective interpolation [12]. Further, to reduce errors associated with spatial heterogeneity a useful extension of our approach would be application of the virtual move of the latitudes of the observations within 5-degree grid cell with spatial step of 0.1˚. This may help to further improve the accuracy of computations.