1961 – 1990 monthly high-resolution solar radiation climatologies for Italy

We present a methodology for estimating solar radiation climatologies from a sparse network of global radiation and /or sunshine duration records: it allows to obtain high-resolution grids of monthly normal values for global radiation (and for the direct and di ffuse components), atmospheric turbidity, and surface absorbed radiation. We discuss the application of the methodology to a preliminary version of an Italian global radiation and sunshine duration data set, which completion is still in progress and present the resulting 1961–1990 monthly radiation climatologies.


Introduction
High-resolution datasets of monthly climatological normals (i.e.high-resolution climatologies) have proved to be increasingly important in the recent past, and they are likely to become even more important in the future.They are used in a variety of models and decision support tools in a wide spectrum of fields such as, just to cite a few, energy, agriculture, engineering, hydrology, ecology and natural resource conservation (Daly et al., 2002;Daly, 2006).
One of the most important variables for a lot of possible applications is solar radiation.Even though some examples of solar radiation climatologies are already available for Italy (see e.g.Lavagnini and Jibril, 1991;Petrarca et al., 2000;Suri and Hofierka, 2004), they suffer the lack of a dense network of long records of observational data.
In this context we set up a research program with the aim of (i) setting up an extensive data base of Italian global radiation and sunshine duration records and (ii) developing a methodology for estimating high resolution solar radiation climatologies from these records.Sunshine duration records have the great advantage, with respect to global radiation records, of a much larger data availability, especially when long-term records are considered.
The methodology for estimating solar radiation climatologies from global radiation and/or sunshine duration records has been developed within a Ph.D. thesis recently concluded at Milan University (Spinoni, 2010).It consists in (i) con-verting sunshine duration data into global radiation (if global radiation data are not available for the site), (ii) decomposing global radiation into a direct and a diffuse component, (iii) gridding direct and diffuse components of global radiation, (iv) evaluating atmospheric turbidity over the same grid by means of the direct component of global radiation, (v) calculating direct, diffuse and reflected components of global radiation for any cell of the used grid, taking into account its slope and aspect and considering shading, (vi) calculating the corresponding absorbed radiation by means of landuse-based albedo estimations.The first application of the methodology consisted in the estimation of global radiation climatologies (Spinoni, 2010) that have been used as proxies to support the construction of 1961-1990 temperature climatologies for Italy.For this application we were mainly interested to the result produced at the last point.Other applications, however, might require the results produced at the other points or modified versions of them (e.g. for solar energy production it might be interesting to produce the result of point v, substituting the slope and the aspect of the gridcell, with values corresponding to an hypothetical panel in a solar plant).
The paper aims at presenting this methodology and showing a preliminary version of radiation climatologies obtained applying it to the Italian global radiation and sunshine duration records that are already available in digital form.The climatologies are represented on the USGS GTOPO30 Digital Elevation Model grid (USGS, 1996), i.e. with a Published by Copernicus Publications.the W-E direction and about 900 m in the S-N direction.This DEM has also been used to estimate the slope and the aspect of the surface and the rate of shading due to the surrounding areas.

Sunshine and solar radiation data
The activities aiming at setting up an extensive data base of Italian sunshine duration records are in progress.They will include the digitisation of a great amount of data that are available only on paper forms, allowing to greatly extend the records that are available at present time.The result will not only be used for describing the spatial distribution of solar radiation normal values, but will also allow the temporal trends over different climatic regions of Italy to be studied.
The new data set is however not available yet.So, we started to develop our methodology with a more limited data availability.In particular we considered monthly records which are already available in digital form from Italian Air Force, ENAV (www.enav.it),CRA-CMA (http://www.cra-cma.it),some regional Environmental Agencies and the European Solar Radiation Atlas (Kasten et al., 1984).The standard method for sunshine duration observation was the Campbell-Stokes instrument; we have however not yet full information on the measures as the data and metadata collection is still in progress.
A significant part of the data used in this paper is available from the web (see http://www.scia.sinanet.apat.it/and http://www.cra-cma.it/).The main deficit of these data is that many records cover rather short time periods and only a minor fraction of them has no missing data in the 1961-1990 period.So the monthly station climatic normals which can be obtained from these records are not completely representative of this period.This problem may introduce a significant bias as it is well demonstrated that sunshine duration and solar radiation are not constant through decades: records show a "global early brightening" period approximately between 1940 and 1950, a "global dimming" period between 1950 and 1980 and a "global late brightening" period after 1980 (Wild et al., 2005;Ohmura, 2006;Wild, 2009).Furthermore, solar radiation and sunshine duration are influenced by great volcanic explosions such as, e.g.Eyjafjallajökull (Iceland, 2010).Due to these and other facts, which cause variability over a wide range of time scales (see e.g.Brunetti et al., 2009), sunshine duration data that are not related to the reference period should be handled with care.A complete solution to this problem will be possible only when the full data set will become available, as the conversion of all station normal values to 1961-1990 normals requires the knowledge of the spatial distribution of the time evolution of sunshine duration over Italy.At present time we only fixed a minimum length of the records (10 yr) and performed some preliminary comparisons among neighbouring stations in order to exclude from the analyses the records which seem to exhibit the largest bias due to missing data in the 1961-1990 period.
Besides the records from the listed sources, we considered also monthly 1961-1990 normals from a monographic book (Cat Berro et al., 2005) and from the web for a few more Italian stations and for about 30 stations of the surrounding countries (France, Switzerland, Austria and Slovenia).
The final data set used in this paper consists of 158 stations (see Fig. 1).For all stations, being they declared at WMO standard, we assume they are far away from surrounding shading obstacles.For 31 of the stations, besides the sunshine duration records, we also use global radiation 1961-1990 monthly normals.These stations are from the Italian Air Force network (see www.meteoam.it).
3 Estimating flat surface solar radiation grids from station data

Converting sunshine duration to global radiation
The first step of our procedure consists in estimating the monthly clearness index (K t ) values for all the stations of our data set that do not have global radiation data.The clearness index is the ratio between the global radiation received by a surface (H T ) and the exo-atmospheric radiation received by the same surface (H 0 ) (Gueymard, 2001).K t includes cloudiness and turbidity.Following the Angstrom-Black's equation (Black et al., 1954), it can be estimated from relative sunshine duration (i.e. from the ratio between the number of sun hours measured by a sunshine recorders (S ) and the solar day length from sunrise to sunset (S 0 )): In Eq. ( 1) we consider sunshine duration as representative of the 15th day of each month and we use corresponding values for S 0 and H 0 : they only depend on astronomical and geographical factors and they can be calculated according to standard procedures (see Spinoni, 2010 for details).Actually, before using Eq. ( 1) for estimating the monthly clearness indexes of all our 158 sites, we use the sites with global radiation data to estimate the coefficients a and b.The results we obtain with present time data availability (31 stations) are shown in Table 1, together with their standard errors.Beside these coefficients, 8 other alternative sets of coefficients found in the literature (Rietveld, 1978;Landsberg, 1981;Andretta et al., 1982;Iqbal, 1983;Newland, 1989;Gopinathan and Soler, 1995;Akinoglu and Ecevit, 1990;Coppolino, 1994) were tested.The errors turned however out to be larger than the ones obtained with the coefficients tuned on Italian data and so we prefer using these coefficients.The Mean Absolute Error (MAE) over all months and stations gives a synthetic information on the ability of Eq. (1) to get the observed clearness index: it results 0.021 with our coefficients whereas it ranges between 0.025 and 0.108 with the 8 other alternative sets of coefficients.When the full dataset will be available, the estimation of the clearness index from sunshine duration data has to be studied more in detail by investigating the use of local coefficients (see e.g.Scharmer et al., 2000) and by trying to consider other variables beside sunshine duration.
Once the monthly clearness index values are available, the corresponding global radiation normals can simply be calculated as:

Decomposition models: from global radiation to direct and diffuse radiation
After global radiation normals are available for all stations, we estimate the direct and diffuse components of solar radia-tion by means of the so called decomposition models (we assume that our sunshine duration data, being measured under WMO standard conditions, are not influenced by reflected radiation).In particular, we use a decomposition model based on Eq. (3) (Iqbal, 1983).
K dif is the diffuse radiation fraction of the global radiation received by a surface, K dir is the corresponding direct radiation fraction, H dif and H dir are the diffuse and direct components of global radiation.Decomposition models should be based on local data.However, when local data are not available, models which are reasonably valid worldwide can be used, e.g. the third order polynomial model used in the European solar Radiation Atlas (Erbs et al., 1982;Scharmer et al., 2000) or the models proposed by Page (1964); Iqbal (1983); Reindl et al. (1990) and Gopinathan and Soler (1995).In our methodology we use (Spinoni, 2010) the following relation (Gopinathan and Soler, 1995): Therefore, in the second step of our procedure, we calculate, by means of Eqs. ( 3)-( 4), the monthly diffuse and direct fractions for all our stations.

Gridding of direct, diffuse, and global radiation for flat surfaces
The third step of our procedure consists in using, for each month, direct and diffuse components of station global radiation to construct high-resolution grids covering all the Italian territory.This gridding procedure is performed, on each node of the USGS GTOPO30 DEM, by means of an Inverse Distance Gaussian Weighting (IDGW) spatialisation model, using Gaussian radial weights (w rad i ) for the contribution of each station: where d i (x,y) is the distance between the i-th station and the considered grid-cell and c d is a coefficient regulating the decrease of the weighting factor with distance: it is chosen in order to have weight equal to 0.5 at distance d: On the basis of present time data availability, we choose d = 50 km.When the full dataset will be available, we will also use more complex spatialisation techniques, trying to take into account the effect of geographic variables.The gridding procedure allows obtaining monthly highresolution fields of direct, diffuse and global radiation that are representative of flat and non-shaded surfaces.Figure 2 shows, as an example, the yearly average global radiation that we obtain for such surfaces with present time data availability.

Evaluation of the turbidity of the atmosphere
The fourth step of our procedure consists in evaluating the spatial distribution of atmospheric turbidity.This evaluation is based on the following relation (Iqbal, 1983): where H dir is the direct component of global radiation calculated for the 15 day of each month as described in the previous sections, E 0 is the eccentricity factor (i.e. the correction due to the elliptical orbit of the Earth), I 0 is the solar constant, θ inc is the solar angle of incidence and the exponential part explains the attenuation due to the atmosphere: T F is the turbidity factor, m A is the optical air mass, δ R is the Rayleigh's depth of the atmosphere.T F represents the turbidity of the vertical column of the atmosphere over the grid cell: clouds, water vapor, pollution, fog, ozone, and many other factors are included in T F .
For each point and each month we search for the T F best matching the H dir in Eq. ( 7).We consider for the integration 5-min time intervals (dh) and calculate the time dependent variables (θ inc , m A and δ R ) over each interval.In the calculation of m A and δ R we also consider elevation and take into account the refraction of the atmosphere (see Kasten and Young, 1989;Rigollier et al., 2000).Details on the calculations can be found in Spinoni (2010).This step of the proce-dure allows obtaining monthly atmospheric turbidity 1961-1990 normals over the same grid used to spatialise H dir .

Solar radiation model for inclined surfaces
Once we know the turbidity of the atmosphere, we can calculate the solar radiation received by inclined surfaces.This calculation requires the knowledge of the slope and the aspect of each grid-cell, as well as the evaluation of the shading due to the surrounding grid-cells: this information is obtained by means of the GTOPO30 DEM.
Direct radiation for inclined surfaces (H incl dir ) is calculated with a slightly modified version of Eq. ( 7), i.e. introducing a binary factor J that represents shading: it is obtained by exploring the grid-cells surrounding each node of the GTOPO30 DEM and checking, with a 5-min time resolution, if the path from the node to the sun does or does not intercept the DEM surface.If the grid cell is shadowed in the 5-min interval that we use in the integration, J is set to 0, otherwise it is set to 1.In this case θ inc is naturally calculated taking into account the slope and the aspect of the surface.
Actually, in spite of the analogies of Eqs. ( 7) and ( 8), they are used in a completely different way: in fact in Eq. ( 7), we know H dir and use it to get T F ; on the contrary, in Eq. ( 8) we know T F and use it to get H incl dir .Diffuse radiation for inclined surfaces (H incl dif ) is calculated considering diffuse radiation as isotropic (as it is usual in solar radiation models made for climate-related purposes).In order to obtain the grids, we just multiply the diffuse radiation received by a flat surface by the sky view factor (V F ), i.e. the visible fraction of the sky from the grid-cell.In our procedure this factor is assumed to be dependent only on the slope (s) of the grid-cell itself.More details on this assumption can be found in Chung and Yun (2004).
Reflected radiation for inclined surfaces (H incl ref ) is calculated as: where O SF is the obstructed sky factor (see Chung and Yun, 2004), i.e. the obstructed portion of the sky and α is the ground albedo.In our procedure we assume for the albedo the value which we attribute to the grid-cell itself, even though the reflection is due to the surrounding cells.This approach is justified as the very limited contribution of reflected radiation that we have with a DEM resolution of 30 arc-seconds, does not justify the much greater complexity which would be necessary in order to take into account the slope, aspect and albedo of the surrounding grid-cells.Albedo is estimated by means of the GLC2000 land cover grid provided by Joint Research Center (see the website: http://bioval.jrc.ec.europa.eu/products/glc2000/glc2000.php)and on the basis of literature albedocloud cover relations (see e.g.Hummel and Reck, 1979;Henderson-Sellers and Wilson, 1983;Wilson and Henderson-Sellers, 1985).Albedo is not corrected for fresh snow cover in winter, because snow climatologies for Italy are not available.
Once direct, diffuse and reflected radiation are available, we simply calculate global radiation summing them.All the procedure is naturally performed for all the grid-points of Italy.The absorbed radiation is then obtained simply by considering the albedo factor.

Validation
We performed a validation of the preliminary climatologies that we present in this paper.This validation is based on a subset of 28 of the Air Force stations with 1961-1990 global radiation normals: they were selected focusing only on stations located in non inclined grid-cells.In other terms we require, not only that a station is located on flat surface, but also that the grid-cell in which it is located is completely flat.The agreement between the climatologies and the station data was evaluated by means of the mean absolute error (MAE) and the relative mean absolute error (MAER).Due to the fact that we used the same data set to evaluate the clearness index, we performed a leave-one-out validation, removing from the input data the station that is, in turn, evaluated.The results are shown in Table 2: MAERs are smaller in late winter and spring than in autumn and early winter, but no systematic over or underestimation was found.The average MAER is under the threshold of 5 %, but such a validation is based on a very small data set and has to be considered as preliminary.

Conclusions and area for further work
We described a methodology for the construction of solar radiation and atmospheric turbidity normal value grids.It requires solar radiation or sunshine duration data.All the points of the methodology described in the paper have been encoded in an unique program which allows the user to handle the different steps of the calculations.Several points require future improvements.The first point on which we are already working consists in a significant enlargement of the data set: it will allow both to enlarge the number of records and to extend the time coverage of each record; moreover it will allow to make available, beside longterm sunshine records, also shorter radiation records that will be used to better study the relation between global radiation and sunshine duration.Another important point that will be considered in the future consists in the use of a DEM with higher resolution: this will improve the evaluation of the influence of the grid-cells surrounding each grid-cell: they regulate shading and influence also direct radiation (by means of the sky view factor) and reflected radiation.

Figure 1 .
Figure 1.Stations with sunshine duration records.Red dots: also global radiation available.

Table 1 .
Monthly a, b coefficients in Eq. (1) and corresponding standard errors based on the 31 Italian Air Force station with both global radiation and sunshine duration data.