Erythemal ultraviolet irradiation trends in the Iberian Peninsula from 1950 to 2011

Erythemal ultraviolet (UVER) irradiation was reconstructed at nine Spanish locations, with series starting around 1950 in five of these places. Each series was checked by applying homogeneity tests in order to discard nonhomogeneous series. Available series were used to create an averaged Iberian Peninsula UVER series. Results indicate that annual UVER irradiation in the Iberian Peninsula increased by 155 Jm (6.5 %) between 1950 and 2011 due to the decrease observed in atmospheric ozone rather than changes in aerosol and clouds. By contrast, annual UVER irradiation increased by 135 Jm (5.6 %) between 1985 and 2011, mainly due to changes in aerosol and clouds.


Introduction
Among other effects, ultraviolet (UV) radiation, which is a part of total solar shortwave (SW) radiation, produces harmful effects on human skin, such as erythema (sunburn) induction (UNEP, 2003).Nevertheless, UV radiation exposure is also beneficial for the health, for example, by contributing towards human Vitamin D synthesis (Webb, 2006).The effectiveness of UV radiation in producing erythema on human skin is usually quantified by the erythemal action spectrum (McKinlay and Diffey, 1987), and the UV radiation weighted by this spectrum is erythemal ultraviolet (UVER) radiation.
The damage caused to human skin by UVER radiation is cumulative and is proportional to exposure time (WHO, 1995).It is therefore important to know both present-day as well as past UVER radiation levels in order to estimate future epidemiological data related to diseases caused by sun exposure.However, it was not until the 1970s and 1980s that the first long-term UVER measurement databases appeared (Scotto et al., 1988;Blumthaler and Ambach, 1990;Correl et al., 1992;Borkowski, 2000).In fact, the longest UVER series in Spain started in 1995.In order to obtain UVER data previous to the use of instrumentation, several authors have reconstructed UVER series using other available records like SW radiation, sunshine duration, and cloud cover (Lindfors et al., 2003(Lindfors et al., , 2007;;den Outer et al., 2005den Outer et al., , 2010;;Rieder et al., 2008;Walker, 2010;Antón et al., 2011a;Bilbao et al., 2011, among others).
UVER radiation is sensitive to factors such as ozone, clouds, and aerosol particles in the atmosphere.Over the last few decades, the presence of these factors in the atmosphere has changed, and might have affected past UVER levels.
Solar radiation reaching the Earth's surface ranges mainly from 305 nm to 2800 nm, its wavelength integration being known as SW radiation, and it plays a key role in the Earth's energy balance and radiation budget.SW radiation decreased between 1950 and the mid-1980s in the Northern Hemisphere, a phenomenon known as "global dimming" (Stanhill and Cohen, 2001).SW radiation began to increase in the mid-1980s in the Northern Hemisphere, a phenomenon known as "global brightening" (Wild et al., 2005).Dimming was caused because aerosol loads increased in the Northern Hemisphere between 1950 and the mid-1980s, absorbing and scattering (aerosol direct effect) more radiation, although after the mid-1980s the aerosol load started to decrease (Wild, 2009(Wild, , 2012)).Changes in aerosols can modify the microphysical properties of clouds (aerosol indirect effect), since aerosols act as condensation nuclei, which contributes to enhancing the dimming and brightening phenomena.The mentioned aerosol and cloud changes might cause variations in the amount of UVER radiation reaching the Earth.Dimming and brightening phenomena were also observed in the Iberian Published by Copernicus Publications on behalf of the European Geosciences Union.
The atmospheric total ozone column (TOC) evidenced major depletion after the late 1970s up to the mid-1990s due to strong atmospheric emission of halogen gases between the 1960s and 1980s (WMO, 2011).TOC evolution in the Iberian Peninsula was studied in depth by Román et al. (2014d) in the dimming and brightening periods, with usually negative but not statistically significant trends being reported in both periods, and a statistically significant trend of −0.7 % dc −1 (percent per decade units) in the annual TOC between 1950 and 2011.Antón and Mateos (2013) showed that these TOC changes have no significant influence on SW radiation (e.g. the SW radiative effects due to TOC variations between the late 1970s and 2012 were smaller than 0.2 Wm −2 over the Iberian Peninsula).Nevertheless, these TOC changes might prove extremely relevant for UVER evolution over the last few decades, marking the difference between the SW and UVER trends in the past.
Many authors have found statistically significant positive trends over the last decades for UVER radiation at different places in the following European countries: Austria, Czech Republic, Finland, Greece, Germany, the Netherlands, Norway, Sweden, and Switzerland (Lindfors et al., 2003(Lindfors et al., , 2007;;Rieder et al., 2008;Walker, 2010;den Outer et al., 2010;Krzyscin et al., 2011).From these papers, the increase in UVER radiation values was mainly attributed to ozone depletion since the late 1970s; some of them also suggested that the UVER increase was partially caused by a diminishing of cloudiness or aerosols.
UVER radiation presents high levels in the Iberian Peninsula due to the great values of the solar elevation, making it an interesting area to study the evolution of UVER radiation.However, in the Iberian Peninsula, UVER radiation has only been reconstructed at Valladolid since 1991 (Bilbao et al., 2011), and at Badajoz and Cáceres (only in the summer months) since 1950 (Antón et al., 2011a), a significant rise in UVER radiation levels having been reported at the three sites.
Recently, Román et al. (2014d) obtained the TOC trends from 1950 to 2011 on the Iberian Peninsula, and Román et al. (2014a) the SW and temperature trends for the same period and region, but not the UVER trends due to the lack of UVER data.As a first step, Román et al. (2014b) simulated the UVER irradiance under cloudless conditions on the Iberian Peninsula with a radiative transfer model, using as inputs monthly climatological values of aerosols, water vapour, etc.They also characterised the uncertainty of these simulations caused by the monthly variability (standard deviation of the used monthly climatological values) and uncertainty of the inputs.These cloudless simulations are useful for reconstructing UVER radiation under all conditions (e.g.Lindfors et al., 2003Lindfors et al., , 2007;;den Outer et al., 2005;Walker, 2010;Bilbao et al., 2011).Therefore, in this framework, the main objectives of the present paper are (1) to develop a new UVER reconstruction model and to optimise another one (both models based on cloudless UVER simulations), (2) to reconstruct UVER irradiation series since 1950 using the developed models with data at certain Spanish locations over the Iberian Peninsula, and (3) to analyse their evolution and trends of the reconstructed data set.A novel issue of this paper is the analysis of the annual and seasonal (not only summer) UVER trends on the Iberian Peninsula from 1950, since these trends were not yet calculated in previous studies.Another novel issue is the quantification of the role of aerosols and clouds (both together), and ozone in the changes in UVER during the last six decades.
The paper is structured as follows: Sect. 2 shows the relevant information concerning the locations, the instrumentation, and the explanation of all the data used.The methods used to obtain the reconstructed UVER series are developed and explained in detail in Sect.3. Section 4 presents the main results for the evolution and trends of UVER.The factors not taken into account in the work are mentioned in Sect. 5. Finally, Sect.6 summarises the main results and conclusions.

Places and instrumentation
All data used in this paper were taken at nine Spanish radiometric stations located in the Iberian Peninsula.These locations are marked in Fig. 1 and their coordinates are also shown in Table 1.The Iberian Peninsula is well covered by these locations.One of these stations is controlled by the University of Valladolid and is located in the village of Villalba de los Alcores (de Miguel et al., 2012).The rest are controlled by the Spanish Meteorological Agency (Moreta et al., 2013).Hourly UVER and SW irradiance were measured at all these stations, although sunshine duration and temperature were not measured at the Villalba station.Therefore, the sunshine duration and meteorological variables, as with temperature, measured at the Valladolid Airport AEMet station were considered the same as at the Villalba station, since the two stations are located just a few kilometres away from each other.
Hourly UVER irradiance was recorded at the nine locations using UVB-1 pyranometers (Yankee Environmental Systems Inc.).These pyranometers were periodically calibrated by a two-step method (Vilaplana et al., 2009), which provides a combined uncertainty (68 % confidence) of between 5.4 and 8.0 % for the measured hourly UVER data (Hülsen and Gröbner, 2007).This uncertainty was considered the maximum (8.0 %) in this work.The oldest UVER data recorded in Spain date from 1 November 1995 and continue up to the present day at the Madrid station.Hourly SW irradiance was initially measured at each location using a CM6B pyranometer (Kipp & Zonen), whose spectral response ranges from 305 to 2800 nm.The expanded uncertainty (95 % confidence) of the hourly SW recorded by this pyranometer is 8 %, and was the expanded uncertainty assumed by all available hourly SW measurements even when more recent records were taken using improved pyranometers (displaying less uncertainty).The longest SW data in Spain date from 11 July 1973 and continue up to the present day at the Madrid station.
In general, it was not possible to obtain information on the instruments used for sunshine duration records at the various stations, although most were probably Campbell-Stokes heliographs (Sánchez-Lorenzo et al., 2007).This heliograph comprises a spherical lens which concentrates direct radiation from the sun onto a dark paper card, which is burned when direct radiation exceeds a certain threshold.The combined uncertainty of the sunshine duration records was assumed to be 15 min (0.25 h).The daily sunshine fraction (F ) is the ratio of the measured sunshine duration to the same sunshine duration under cloudless conditions (SunDu cl ).This variable was calculated by the following equation (Iqbal, 1983): where SunDu cl is in hours, δ is the solar declination, is the location latitude, and θ S is the solar zenith angle (SZA) at sunset and sunrise (equal to 87 • in this case).SunDu cl was calculated between the solar zenith angle of 87 • near sunrise and the SZA equal to 87 • near sunset, since direct solar radiation usually is not enough to burn the dark paper card even under cloud-free conditions for a SZA below 87 • .The oldest available F data in Spain date from 1 January 1920 and continue up to the present day at the Madrid station.F data have been available at certain locations (A Coruña, Madrid, San Sebastián, Tortosa, and Villalba) since the 1950s.
The instruments used to take all the mentioned measurements were well calibrated on a regular basis, following World Meteorological Organization (WMO) recommendations (Webb et al., 2006;WMO, 2008) for instrument maintenance, and involved bubble levelling of the instruments, cleaning domes, monitoring and replacing desiccant, etc.In addition, a high-quality control of UVER, SW, and F data was applied to all available data in order to reject possible spurious and outlier data.Daily UVER and SW irradiation data were obtained by integrating the hourly values each day.

Additional data
Other atmospheric variables were also used in this work: aerosol optical depth (AOD), Angström exponent, aerosol single scattering albedo (SSA), water vapour column (w), TOC and surface albedo.The daily data of these variables (except TOC) were used to calculate a monthly climatology in order to be used as input in a radiative transfer model to calculate UVER and SW irradiance under cloudless conditions (see Sect. 3.1).These daily data are described in this section, and the monthly climatology values were obtained from them.
Daily AOD at 433 and 670 nm (MISR-Terra Prod.ver.31:MIL3DAE.004product) from the MISR (Multi-angle Imaging SpectroRadiometer) instrument were obtained at each location from 2000 to 2012.The Angström exponent was calculated directly using the AOD values at the two wavelengths.The daily w (MODIS-Terra Ver.5.1: MOD08_D3.051product) from MODIS (MODerate resolution Imaging Spectroradiometer) was also obtained at each location between 2000 and 2012.AOD and w data were downloaded from the GIOVANNI (GES-DISC Interactive Online Visualization ANd aNalysis Infrastructure; http: //disc.sci.gsfc.nasa.gov/giovanni)application as an averaged 0.2 • × 0.2 • square centred at each location (Acker and Leptoukh, 2007).The SSA at 354 and 500 nm retrieved from the OMI (Ozone Monitoring Instrument) instrument between 2004 and 2011 were also obtained for all locations as overpass files available at AVDC (Aura Validation Data Center).These data are the same as those used by Román et al. (2014c), who calculated the uncertainty of some of these products in the Iberian Peninsula.The combined uncertainty of AOD at 433 nm and 670 nm is 0.074 and 0.054, respectively.The combined uncertainty of the Angström exponent is below 0.5 when AOD at 433 nm is above 0.25, except for high Angström exponent values.The combined uncertainty of the water vapour column is between 0.38 cm (w = 0.5 cm) and 0.52 cm (w = 3 cm).
A daily TOC series for 1950 to 2011 was available for each location.These series comprised different databases: groundbased ozone; ozone retrieved from the TOMS (Total Ozone Mapping Spectrometer) instrument on board the Nimbus-7, Meteor-3, and Earth-probe satellites; TOC from OMI; retrieved TOC from the GOME (Global Ozone Monitoring Experiment) and GOME-2 instruments on board the ERS-2 and MetOp-A satellites; and reconstructed ozone data from the COST-726 project (Krzyscin, 2008;www.cost726.org).The construction of these TOC series was explained by Román et al. (2014d), who, by means of an intercomparison with ground measurements, calculated that the combined uncertainty of the daily TOC values of these series was around 10 Dobson units (DU).

Simulations under cloudless conditions
Global, diffuse and direct horizontal UVER and SW irradiances were simulated under cloudless conditions using a ra-diative transfer model (UVSPEC/libRadtran) from 1950 to 2011 each hour for all locations shown in Fig. 1.UVSPEC is the main tool of the libRadtran (version 1.7 in this work) software package developed by Mayer and Kylling (2005).For UVER simulations, irradiance was calculated every 1 nm from 280 to 400 nm under cloud-free conditions using the "cdisort" solver with six streams (Buras et al., 2011) and the "SUSIM SL2" extraterrestrial spectrum (Van Hoosier et al., 1988), these obtained values then being weighted with the erythemal action spectrum.For SW simulations, the model was run under cloudless conditions using the "twostr" solver (Kylling et al., 1995), the extraterrestrial spectrum from Kurucz (1992), and the "SBDART" pseudo-spectral k distribution from Ricchiazzi et al. (1998).Irradiance was calculated from 305 to 800 nm in 2 nm intervals, from 800 to 1600 nm in 5 nm intervals, and from 1600 to 2800 nm in 10 nm intervals, these spectral values then being spectrally integrated to obtain SW irradiance.An hourly (UVER or SW) irradiance value was simulated at a fixed SZA given by the averaged cosine of the SZA over the hour.
The UVSPEC model was run using standard vertical profiles.A mid-latitude summer atmosphere with springsummer aerosol profiles was used as input for the months from May to October, with a mid-latitude winter atmosphere with autumn-winter aerosol profiles being selected for the other months (Anderson et al., 1987;Shettle, 1989).These vertical profiles were rescaled with monthly climatological tables of water vapour, AOD at 443 nm, Angström exponent, and SSA (at 354 nm for UVER simulations and at 500 nm for SW ones).These climatology tables (one per location and variable) comprised 12 monthly averaged (using all available data) values for each variable, said climatological tables being available in Román et al. (2014c).The daily TOC at each location was included in the inputs, changing the value for each location on each day.Finally, the surface albedo values were also monthly averaged for each wavelength (twelve monthly values at each location for each available wavelength), these monthly values being linearly interpolated to obtain the monthly surface albedo at all the other wavelengths to be then used as input in the radiative transfer code (it is well detailed in Román et al., 2014d).
Both the combined and expanded uncertainties of all simulations were calculated using the results obtained by Román et al. (2014b), who calculated the maximum variations in simulated UVER and SW irradiance caused by the uncertainty of the inputs.Simulated hourly SW and UVER values were also compared with global SW and UVER irradiance measurements under cloudless conditions by Román et al. (2014b).It was found there was better agreement for low SZA values, and that the differences between simulations and measurements were in agreement within the uncertainty.Daily UVER and SW irradiation were calculated for each day at each location, and the uncertainties of these daily values were also calculated.Román et al. (2014b) compared these simulations with measured irradiation under cloudless conditions, with better agreement being found for the spring and summer months.For all months and locations together, a mean bias error (MBE) of −0.1 % and a root mean square error (RMSE) of 3.6 % for the SW case, and a MBE of 2.9 % and a RMSE of 7.7 % for the UVER case, were also reported.

Cloud modification factor
The cloud modification factor (CMF) is defined as the ratio between measured radiation and simulated radiation under cloudless conditions: where R can be SW or UVER, the "me" index indicating the measured R value, and the "cl" index standing for simulated R radiation under cloudless conditions.CMF can be calculated for hourly irradiance or daily irradiation values.CMF is a useful variable to quantify cloudiness, CMF near 1 indicating cloudless conditions, CMF close to 0 being for high overcast cloudiness, and CMF above 1 indicating enhancement effects (Sabburg and Parisi, 2006;Sabburg and Calbó, 2009;Piedehierro et al., 2014).CMF for SW (CMF SW ) is usually different for UVER (CMF UVER ).The relationships between CMF UVER and other variables, like CMF SW , are an important research topic, since they prove useful for reconstructing UVER data.If a relationship between CMF UVER and CMF SW is given by a function f , then UVER radiation can be obtained when the measured SW radiation, UVER and SW simulations under cloudless conditions, and the f function are known: where UVER is the calculated UVER radiation, and UVER cl is the simulated UVER under cloudless conditions.The main goal of this section is to obtain relationships such as the one given by Eq. ( 3) in order to reconstruct UVER data as a function of other variables as in Eq. ( 4).

Model based on hourly SW irradiance measurements
Some UVER reconstruction models are based on different measured variables (Calbó et al., 2005), with several UVER reconstruction models being based on SW radiation measurements, which were normally measured before UVER (Bodeker and McKenzie, 1996;Kaurola et al., 2000;Matthijsen et al., 2000;den Outer et al., 2000den Outer et al., , 2005;;Lindfors et al., 2007;Walker, 2010;Antón et al., 2011b;Bilbao et al., 2011).The data used to obtain a model based on SW were the measured values when both SW and UVER irradiance measurements were available.A total of 294 047 pairs of data (SW and UVER) were available for SZA below 85 where c and d are two parameters which can be calculated by a least square fit.Bilbao et al. (2011) calculated c and d using all data measured at one location and did not take uncertainty into account in the variables of Eq. ( 5).It motivates this paper to improve the reconstruction method; hence, a novelty in this paper is that c and d are calculated considering a similar number of data for different cloudiness and SZA conditions, and taking into account the uncertainty in these data, which optimised the mentioned method.To this end, a data set of CMF UVER -CMF SW pairs of data was selected in this paper.
All available data pairs of CMF UVER -CMF SW from all locations were divided into data groups, taking into account the cosine of SZA.Eighty-seven data groups were obtained separated by intervals of the cosine of SZA, from 0.095 to 0.965 in 0.01 steps.Each data group was divided into a further 15 subgroups considering the CMF SW value in the intervals from 0 to 1.5 in 0.1 bins.A total of 1305 subgroups were available.Fifty pairs of CMF UVER -CMF SW data were chosen randomly for each subgroup.Some groups with high CMF SW values had fewer than 50 data.Finally, 49 777 data pairs were randomly selected.This method of choosing data ensures that the number of data selected is similar and is representative of any cloudiness and SZA condition.
The combined uncertainty of each selected datum was known; hence, the c and d parameters of Eq. ( 5) were calculated by the weighted least square method using the 49 777 chosen data.The weight in the fit was the inverse of the square combined uncertainty.Uncertainty in SZA was considered null.The result of the fit gave a c value equal to 0.6106 with a combined uncertainty of 0.0014, and a d value equal to 0.358 with a combined uncertainty of 0.002.These calculated values are different to those obtained by Bilbao et al. (2011).Once the c and d values have been calculated, the f function of Eq. ( 4) can be expressed as f (CMF SW ) = (CMF SW ) (0.6106±0.0014)+(0.358±0.002)cos(SZA) .

Model based on daily sunshine fraction measurements
Another variable used to reconstruct UVER radiation is sunshine fraction (Lindfors et al., 2003), due to the availability of longer series than SW; it has the advantage of reconstructing older data.The sunshine fraction represents the fraction of the day when the sun was not blocked by clouds.Therefore, by way of an initial approach, daily direct (on horizontal) UVER irradiation (UVER dir ) can be expressed as where UVER dir cl is the daily direct UVER irradiation simulated under cloudless conditions on a horizontal surface.Assuming Eq. ( 7), daily global UVER irradiation can thus be expressed as the sum of the direct (on a horizontal surface) and diffuse (UVER dif ) components: As a second approximation, daily diffuse UVER irradiation can be assumed to be a function (g), which depends on the daily diffuse UVER under cloudless conditions (UVER dif cl ) and on F : Eq. ( 9) in Eq. ( 8) gives Global UVER irradiation can be obtained through F measurements and cloudless simulations using Eq. ( 10) once the g function is known.The g function can be calculated by F and UVER measurements: In order to calculate the g values, a new data set made up of daily UVER irradiation and F measurements was chosen.Days without SW, UVER or F measurements were discarded in the new data set.Available data (all locations) were then divided by the season: winter (January, February, and December), spring (March, April, and May), summer (June, July, and August) and autumn (September, October, and November).10 % of the available data for each season were randomly selected and used to obtain the g function.
Figure 2 shows the values of g calculated by Eq. ( 11) as a function of UVER dif cl for four F intervals.The g values increase linearly with cloudless diffuse UVER irradiation, showing a correlation coefficient (r) above 0.87 for all F intervals.Moreover, the dependence of g on UVER dif cl varies with F , and shows a higher slope when F increases.The y intercept of the fits is close to zero; then, assuming that UVER dif should be equal to zero when UVER dif cl was zero, the g function could be expressed as where B is a parameter which depends on F .B was calculated through a least square fit using the same data as in Fig. 2, but for the following F intervals: F = 0, 0 < F < 0.2, 0.2 ≤ F < 0.4, 0.4 ≤ F < 0.6, 0.6 ≤ F < 0.8, 0.8 ≤ F < 1, F ≥ 1.These intervals were chosen in order to have more resolution of F for the model development, while the intervals of Fig. 2 were wider because they were enough to find the kind of dependence on F of the g function vs.
UVER dif cl .The values obtained are shown in Fig. 3 with their combined uncertainty as a function of F .This B parameter increases with F , and is close to 1 when F tends to 1. Uncertainty is around 0.2, except for F equal to 1, when it strongly increases up to 0.55 (caused in part by the high uncertainty in UVER dir cl , which is not reduced by multiplying by F ).The B value can be obtained by interpolating the values of Fig. 3 for an F measurement.This interpolated value of B can be used to retrieve global UVER irradiation, taking into account Eq. ( 12) in Eq. ( 10): The method for retrieving global UVER irradiation by finding a B value by an F measurement interpolation in Fig. 3 and using this B value with the measured F and cloudless simulations in Eq. ( 12) was called "M-F" in this paper.

Models vs. measurements
Hourly UVER irradiance was reconstructed with model M-SW using the available data which were not used to develop model M-SW.Twenty-five UVER data (3 days) were rejected by applying a new quality control.A total of 220 105  In order to quantify the agreement between the models and the daily UVER measurements, all data not used to develop model M-F (90 % of data per season) were selected.The reconstructed UVER irradiance these days was integrated to obtain the reconstructed daily UVER irradiation by model M-SW.Daily irradiance was also reconstructed by model M-F.The available number of pairs of data for reconstructed and measured daily UVER irradiation was 21 349, as can be seen in Table 2. Table 2 also shows that model M-SW fits better (with a RMSE of 8.4 %) within measurements than model M-F (RMSE higher than 20 %); moreover, model M-F overestimates the daily measurements, showing a MBE of 5.1 %.RMSE and MBE were also calculated for different F intervals (not shown) and model M-SW presents RMSE values lower than 10 % and MBE closest to 0 % for all F intervals; model M-F shows good agreement with measurements for F higher than 0.5, but RMSE and MBE increases when F decreases.
All available data were used to obtain daily UVER irradiation, and these values were monthly averaged (at least 25 days per month) for each location and year.The same monthly averages were reconstructed with the models.The MBE and RMSE obtained using the monthly averages of reconstructed and measured UVER irradiation are shown in Table 2.These values are similar for both models, the MBE being around 2 % and the RMSE around below 6.5 %.The number of data for the models is different due to the lack of F data in some periods at Madrid in the 2000s decade.Monthly averages were used to obtain the annual averaged values (always with all 12 months available) of UVER irradiation at each location.Annual reconstructed and measured data were compared also in Table 2 showing the statistical estimators of the comparison.Both models present a good agreement, with measurements of RMSE being lower than 3 % and MBE below 1 % in absolute values.95 % of the reconstructed annual UVER data with model M-SW shows a difference below 5 % with the measurements.
Finally, the differences between measured and reconstructed values are in agreement within uncertainty for both models and for all temporal resolutions (not shown).Model M-F measurements are not the best for reconstructing daily UVER irradiation, although it evidences great agreement with monthly and annual averaged measurements.A depth comparison not shown in this paper of these models can be found in Román (2014).

Final reconstructed series
Hourly UVER irradiance was reconstructed at each location, when hourly SW records were available, using the method referred to as M-SW.Daily reconstructed UVER irradiation was obtained by integrating the hourly reconstructed values.Table 1 (M-SW column) shows the number of daily UVER irradiation data reconstructed by the model based on SW records.This number is around 10 000 (27 years) for A Coruña, Cáceres, and Murcia, and is over 13 000 (36 years) for Madrid.Villalba shows the lowest number of reconstructed UVER data by this model due to the scant number of SW records at this location.When SW records were not available, UVER irradiation was reconstructed using method M-F and based on F measurements.Table 1 (M-F column) shows the number of daily UVER irradiation data reconstructed by the model based on F records.This number varies with the location, and is higher for Villalba, San Sebatián, and Tortosa.The reconstructed daily UVER irradiation series were completed with daily measured UVER irradiation when SW and F records were not available.The number of measured data used to form the UVER series is shown in Table 1 ("Measured data" column).This number of data is low compared to the other reconstructed data, is less than 20 at all locations, and is zero in Murcia and Villalba.Finally, a long-term daily UVER irradiation series was obtained at each location using models and measurements.The total number of data of these series and the year of the first UVER irradiation value are shown in Table 1.A Coruña, Madrid, San Sebastián, and Villalba show UVER series with more than 20 000 data, all commencing in the 1950s.The lowest number of daily UVER data is around 10 000 for Cáceres and Murcia, whose UVER series started in the mid-1980s.The combined uncertainty of all daily UVER irradiation values of the obtained series was also calculated taking into account the uncertainty in the cloudless simulations and measured values.

Anomalies and homogeneity
Monthly averages of daily UVER irradiation were calculated using the available reconstructed series taking into account at least 20 daily UVER data per month, year, and location.The monthly series obtained were deseasonalised calculating their monthly anomalies considering the reference period from 1985 to 2011.UVER anomaly (A) in month "m" and year "y" is thus calculated as where N is the number of data used in the sum of Eq. ( 14).Monthly UVER irradiation anomalies were calculated for all months at all locations.The monthly anomalies of UVER irradiation at the nine locations were averaged, obtaining a new monthly series of anomalies representative of the Iberian Peninsula.This was then called the "Iberian Peninsula" series.Annual UVER anomalies were calculated by averaging the monthly anomalies when at least six monthly data were available for each year.Seasonal anomalies were also calculated when at least two monthly data were available in winter (December, January, and February), spring (March, April, and May), summer (June, July, and August), and autumn (September, October, and November).Winter anomalies were calculated with the January and February anomalies for a specific year, together with the December anomaly of the previous year.Homogeneity of all these averaged daily UVER irradiation anomaly series was tested in a similar way to the TOC and SW series analysed by Román et al. (2014a, d).The homogeneity tests are mainly applied in order to know whether the series are valid for trend studies or, on the contrary, whether they are not valid due to instrumentation/measurement problems; additionally, the homogeneity tests can provide information about break points in some climatic trends.In this case, the null hypothesis assumes that a temporal series is homogenous, and this hypothesis was verified using the standard normal homogeneity test (SNHT), the Pettitt test, the Buishand test, and the Von Neumann ratio (Wijngaard et al., 2003).Hakuba et al. (2013) assumed that a temporal series is non-homogeneous when the null hypothesis is rejected with a confidence of 99 % by at least three of the four mentioned tests.The four tests were directly carried out on the annual UVER series.Eight (all except Madrid and Murcia) of the ten (nine locations plus the averaged Iberian Peninsula) annual series show non-homogeneity around the mid-1980s, which could be caused by a climate change in UVER from the dimming to brightening periods.Homogeneity analysis was thus performed for the same annual UVER series for the 1950-1984 and 1985-2011 periods.The first period evidences inhomogeneities in San Sebastián and A Coruña, and the second period is free of inhomogeneities.The same result was obtained by Román et al. (2014c) for the annual SW irradiation series.
The monthly series was also subjected to homogeneity analysis by applying the four tests to the synthetic reference series generated with UVER data from the other locations (Alexandersson and Moberg, 1997;Sánchez-Lorenzo et al., 2013c).No monthly UVER series shows inhomogeneities for the 1985-2011 period, and only one month shows inhomogeneity in Madrid and San Sebastián in the 1950-1984 period.These results thus indicate that all the UVER anomaly series can be considered homogeneous, or at least not inhomogeneous enough to change the UVER series values.

UVER temporal evolution
Figure 4 shows the annual UVER daily irradiation anomaly series for the nine locations and the averaged Iberian Peninsula series.Anomalies are shown with their combined uncertainty.Figure 4 panels present a 21-year Gaussian low-pass filter to reduce noise in the evolution.Moreover, the linear trends calculated by the least square method are plotted for the 1950-2011, 1950-1984, and 1985-2011 periods.Madrid shows a strong peak increasing at the beginning of 1970s, while A Coruña exhibits a negative peak around 1980.These peaks look strange, but the homogeneity tests do not consider them to be inhomogeneities, and they also appeared in the sunshine fraction and SW irradiation series in the same places (Román et al., 2014a).In addition, a relative peak in the daily maximum temperature series appeared in the date periods in the same places (not shown), which shows that the UVER peaks were real, but which does not explain the cause of them.All annual UVER anomalies of Fig. 4 display an in- crease for the periods 1950-2011and 1985-2011, except in A Coruña. By contrast, in the period 1950-1984, the annual UVER series, which contains most data, shows a slight decrease.San Sebastián exhibits the clearest change in UVER in the 1980s, which explains the break point detected in this location with the homogeneity tests.Román et al. (2014a) and Román et al. (2014d) found that SW irradiation and TOC decreased until the mid-1980s (the dimming period) at the same locations.These results suggest that UVER in the dimming period experienced two opposite effects: UVER tended to decrease, due mainly to the increase in aerosols (as in SW irradiation), and, at the same time, UVER tended to increase due to ozone depletion.The two effects offset one another, leading to a not statistically significant change in UVER irradiation over the dimming period.
An example of the balanced effects between the impact of aerosol increase and ozone depletion can be seen after a major volcanic eruption during which vast amounts of gaseous compounds can be shot into the stratosphere.These are precursors to the atmospheric formation of sulfate aerosol particles, which in turn provide surfaces for heterogeneous processes on polar stratospheric clouds in the lower stratosphere, enhancing ozone depletion (Peter, 1997;Solomon, 1999;Rieder et al., 2013).In sum, aerosol load increases and ozone decreases after a violent volcanic eruption.The volcanic eruptions at El Chichón (México) and Pinatubo (Philippines) in 1982 and 1991, respectively, are highlighted in Fig. 4. Román et al. (2014a) found that in the years after these eruptions, a reduction in SW irradiation due to the in-crease in sulfate aerosol particles is apparent.However, by contrast, UVER irradiation shows little increase in most of the panels in Fig. 4. The UVER increase is caused by the strong decrease in TOC after volcanic eruptions detected at these locations by Román et al. (2014d), which leads to an increase in UVER more than a UVER decrease caused by aerosols.

UVER trends: dimming and brightening periods
Figure 4 shows the annual UVER evolution and trends in qualitative but not quantitative terms.In order to quantify them, the temporal trends of the monthly, seasonal, and annual UVER anomaly series were assumed to be the trends obtained by the Theil-Sen (TS) estimator and its 95 % confidence interval (95CI).The statistical significance of each calculated trend was evaluated by the non-parametric Mann-Kendall test, considering three types of trends: with a confidence of 99 % (p < 0.01), with a confidence of 95 % but not 99 % (p < 0.05), or non-significant at least at 95 % confidence (p ≥ 0.05).All these estimators were calculated following the methods of Gilbert (1987).If the Mann-Kendall test considered a trend to be statistically significant with at least 95 % confidence, this trend was then assumed to be just significant.A trend was only calculated when a series comprised more than 10 data.TS o3 is the TS trend calculated with the same UVER irradiation series but simulated under cloudless conditions.The TS o3 value gives the UVER trend caused by changes in TOC because aerosol and clouds changes are not included in cloudless simulations.TS ac is defined as TS minus TS o3 , and represents the UVER trend brought about by changes in aerosol and clouds (both together).
The trends (and their significance and 95CI) of monthly, seasonal, and annual UVER irradiation series were calculated for all locations at three periods: 1950-1984 (considered the dimming period), 1985-2011 (considered the brightening period), and 1950-2011 (whole period).Fig. 5 shows all these values following the methodology used by Walker (2010), who plotted the trends as points whose type depends on the trend significance, and he also added the 95 % confidence intervals as error bars.The significant seasonal and annual trends are also shown in Table 3, which also shows TS o3 and TS ac trends.The most representative trends of the dimming period are those obtained for the San Sebastián, Madrid, Villalba, and Iberian Peninsula series, since they have fewer missing data.The mentioned series show no significant trends in the dimming period, except San Sebastián, which evidences a significant trend of −211 Jm −2 dc −1 (−7.7 % dc −1 ) in May, and an annual trend of −48 Jm −2 dc −1 (−2.8 % dc −1 ).The negative annual trend found in San Sebastián is mainly caused by changes in aerosol and clouds radiative effects rather than ozone effects, since TS ac is much higher than TS o3 .Román et al. (2014a) found many more significant trends in SW irradiation for the same locations during the dimming period, underpinning the key role played by ozone decrease in UVER trends during the dimming period.
The brightening period  has the advantage that all UVER series are complete.UVER trends are mainly significant in summer and in the annual series.All series, except for A Coruña and Madrid, present significant trends in summer, and are 109 Jm −2 dc −1 (2.5 % dc −1 ) for the Iberian Peninsula series.The annual trend of the Iberian Peninsula series is 50 Jm −2 dc −1 (2.1 % dc −1 ).The TS ac and TS o3 values in Table 3 reveal that UVER irradiation increased in the brightening period due to a reduction in aerosol and cloud radiative effects and in ozone effects.The ozone changes represent approximately two thirds of the trend caused by aerosol and cloud changes.
As regards the 1950-2011 period, the most interesting series are San Sebastián, Madrid, Villalba, and the Iberian Peninsula series for the same reason as during the dimming period.All annual UVER series trends, except Madrid, are significant at 99 % (95 % for Murcia), this trend being 25 Jm −2 dc −1 (1.1 % dc −1 ) in the Iberian Peninsula, which indicates an increase of 155 Jm −2 (6.5 %) in annual UVER irradiation over the last 62 years in the Iberian Peninsula.March presents positive and significant trends for the four most complete series, ranging from 31 Jm −2 dc −1 (San Sebastián: 2.2 % dc −1 ) to 74 Jm −2 dc −1 (Villalba: 4.1 % dc −1 ).June and July exhibit the highest UVER trends, all of them proving significant, except Madrid in July.The UVER trend in the Iberian Peninsula is 83 Jm −2 dc −1 (1.9 % dc −1 ) in June and 47 Jm −2 dc −1 (1.0 % dc −1 ) in July.As regards seasonal trends, all are significant in spring and summer, except for spring in Madrid.The trend in the Iberian Peninsula series is 32 Jm −2 dc −1 (1.2 % dc −1 ) in spring and 54 Jm −2 dc −1 (1.2 % dc −1 ) in summer.UVER trends in the 1950-2011 period are mainly caused by changes in TOC, because these trends are similar to the obtained values of TS o3 , TS ac being around zero.This result indicates that aerosol and cloud presence increasing during dimming and decreasing during brightening are well compensated for over the 1950-2011 period.

UVER trends: comparison with literature
In the previous section, UVER trends were calculated for three specific periods.Other authors, however, have calculated UVER trends for other periods in the literature.UVER trends were thus recalculated with the same UVER anomaly series but for other periods in this section in order to compare them with the results of other works.Lindfors et al. (2003) reconstructed UVER irradiation at Sodankylä (Finland) between 1950 and 1999 also using sunshine duration records, and reported two significant trends for this period: 3.9 % dc −1 in March and −3.3 % dc −1 in July.Slightly higher significant trends were found in March for Madrid (5.3 % dc −1 ), Villalba (6.6 % dc −1 ), and the Iberian Peninsula (4.4 % dc −1 ) in the same period.These series do not present statistically significant trends in July for the 1950-1999 period.Josefsson (2006) analysed the measured UVER at Norrköping (Sweden) between 1983 and 2003, finding significant trends in UVER irradiation in spring (7.8 % dc −1 ), autumn (8.2 % dc −1 ), and in the annual series (5.2 % dc −1 ).Significant although lower trends were also detected in the same period in the series analysed in this paper.The Iberian Peninsula series showed a significant trend of 4.4 % dc −1 in spring and 2.7 % dc −1 for the annual series.Lindfors et al. (2007) reconstructed UVER irradiation using SW records from 1983 to 2005 in Bergen (Norway), Jokionen (Finland), Norrköping and Sodankylä, and reported a significant increase in annual UVER at Sodankylä (4.1 % dc −1 ).For the same period, trends in the annual UVER series at Cáceres, Murcia, and the Iberian Peninsula were similar and also significant, being 3.2 % dc −1 , 2.7 % dc −1 , and 2.7 % dc −1 , respectively.
den Outer et al. (2010) obtained UVER irradiation between 1980 and 2006 using different reconstruction models at eight European locations: Sodankylä, Jokionen, Norrköping, Potsdam (Germany), Lindenberg (Germany), Bilthoven (the Netherlands), Hradec Kralove (Czech Republic), and Thessaloniki (Greece).The annual UVER trends obtained by den Outer et al. (2010) range from 2.8 to 5.8 % dc −1 .The trends in the annual UVER series of this paper are all significant for the same period, except for Murcia, and range from 1.8 % dc −1 (Madrid) to 5.3 % dc −1 (San Sebastián), with the Iberian Peninsula trend being 3.2 % dc −1 .These results are similar to those obtained in the literature.Walker (2010) analysed the trends of reconstructed UVER irradiation between 1981 and 2007 at four Swiss locations: Davos, Payerne, Locarno, and Jungfraujoch, with the UVER irradiation trend proving significant in March and June for all four locations during this period.Spring and summer months present the highest number of significant trends in the UVER series of this work for the 1981-2007 period.Annual trends in Switzerland were similar to Spain (eight of the ten are significant), with values varying between 2 and 4 % dc −1 from 1981 to 2007 in both countries, with 3.0 % dc −1 being the trend in the annual Iberian Peninsula series.Krzyscin et al. (2011) studied the UVER radiation observed at Belsk (Poland) between 1976 and 2008, and found an annual trend of 5.6 % dc −1 .In the same period, the trend in the annual Iberian Peninsula series was 2.8 % dc −1 , half that of Belsk.
As regards UVER trends in the Iberian Peninsula obtained by other authors, three papers are well known.Bilbao et al. (2011) reconstructed UVER at Valladolid from 1991 to 2010 and found significant trends in summer (3.5 % dc −1 ) and autumn (4.1 % dc −1 ), similar to those obtained for the 1991-2010 period at Valladolid with the reconstructed series used in this paper: 3.6 % dc −1 (summer) and 5.9 % dc −1 (autumn).Antón et al. (2011a) reconstructed UVER irradiance at solar noon in summer from 1950 to 2000 at Badajoz (Spain) and Cáceres, and obtained a trend of 4.9 % dc −1 for the 1979-2000 period at Cáceres.The significant UVER irradiation trend at Cáceres during the same period was 5.2 % dc −1 using data from the present work, a similar value to that obtained by Antón et al. (2011a).This means that the UVER irradiance trend at solar noon was similar to the daily UVER irradiation trend at Cáceres.Finally, Ialongo et al. (2011) calculated UVER irradiation from 1979 to 2010 over the whole world using satellite images and found a UVER trend in the Iberian Peninsula of around 2.5 % dc −1 in March and October.The annual Iberian Peninsula series of the present work shows a similar trend (2.8 % dc −1 ) for the same period.
A comparison between the results obtained and those reported in the literature reveals that UVER changes in Europe have been similar, at least over the last few decades.

Effect of UVER uncertainty on trends
The trends obtained were calculated without considering uncertainty in the UVER irradiation values, although uncertainty might influence the value and significance of trends.The effect of uncertainty on trends was studied following the method used by Román et al. (2014a, d).
For one specific series of UVER anomalies A, with N values A i , and σ (A i ) being the uncertainty of the i value of A, a normally distributed (centred on zero with a standard deviation of σ (A i )) random number (b i ) is generated for each i Table 3. Statistically significant UVER irradiation trends with a confidence of 99 % (95 % marked with an asterisk) and their 95 % confidence interval, at different seasons and locations for the 1950-2011, 1950-1984, and 1985-2011 periods. TS o3 periods. TS o3 and TS ac are also included.N is the number of data used.

Location
Period Season N TS TS TS o3 TS ac 95CI (Jm −2 dc −1 ) (% dc −1 ) (Jm −2 dc −1 ) (Jm −2 dc −1 ) (Jm −2 dc −1 )   1950-2011, 1950-1984, and 1985-2011 is about 68 %, a probability that increases to 95 % when the interval is [A i − 2σ (A i ), A i + 2σ (A i )].Once the N values of b i for all values of A are obtained, a new synthetic series (SS) is formed as the sum of the original A series and the random numbers, the i value of the SS series being equal to A i +b i .The SS series is physically valid since it is inside the uncertainty of the A series. 10 000 SS series were randomly obtained for each series analysed in the previous section following the methodology described, and their trends and significance were cal-culated.The percentage of the 10 000 series whose trend is statistically significant at 95 % (P (p < 0.05)), and at 99 % (P (p < 0.01)), was calculated.Only the series with values of P (p < 0.05) and P (p < 0.01) above 95 and 99 % were considered statistically significant at 95 and 99 %, respectively.55 % (65 out of 119, considering the 10 locations and the 12 monthly + 4 seasonal + 1 annual) of the significant trends (in the 1950-2011, 1950-1984, and 1985-2011 periods) with at least 95 % confidence are considered significant (95 % confidence) using the criterion based on uncertainties.All series not considered significant in the previous section are not significant with the uncertainty-based criterion.
The P (p < 0.05) and P (p < 0.01) values of the seasonal and annual UVER series considered statistically significant with at least 95 % confidence taking into account the uncertainty are shown in Table 4. High P (p < 0.05) values do not imply high P (p < 0.01) values.The median and standard deviation of the 10 000 values of TS calculated from all synthetic series are also included in Table 4.These trends are similar to the trends in Table 3, the differences between them being below 5 Jm −2 dc −1 (0.1 % dc −1 ) in most cases, although this difference reaches 18 Jm −2 dc −1 (0.4 % dc −1 ) for Ciudad Real in summer in the 1985-2011 period.These differences are always below the standard deviation of the trend given in Table 4 (except for Ciudad Real in summer in the 1985-2011 period).

Factors not taken into account
The results of this paper were obtained using reconstructed data series by models.However, the paper is not without certain limitations.Changes in aerosol optical depth, surface albedo or water vapour column were not considered, these variables being used in the radiative transfer model as monthly climatology.
The lack of AOD data before 2000 led to the use of a climatological table which does not contain the aerosol changes in the 1950-2011 period.However, the aerosol effect is included in the SW and F measurements (like clouds) and, as a first approximation, the reconstruction models transfer this effect to the UVER radiation.This approximation is not valid for the case of water vapour, because it affects SW and F , but not UVER.
If the water vapour column had increased, cloudless SW irradiation would have decreased and, therefore, reconstructed UVER should be higher.Daily water vapour column trends were calculated between 1957 and 2002 in the Iberian Peninsula using the daily ERA-40 data (Uppala et al., 2005;Lindfors et al., 2007).These trends (not shown) indicate a slight water vapour decrease in recent decades, which did not always prove significant, but which might contribute to reducing the UVER trends obtained.
Trends in surface albedo at 360 nm (Sect.2.2) from 1958 to 2002 were calculated (not shown), revealing that, apart from a slight decrease in winter months, surface albedo has suffered no significant changes in recent decades.This result indicates that the UVER trend obtained in winter might be slightly lower due to changes in albedo, but should not affect the remaining months.
A further factor to take into account should be the uncertainty of the data used (caused in part by the monthly variability of the radiative transfer inputs), since certain trends cannot be considered significant when uncertainty is taken into account, as can be seen in Sect.4.5.Moreover, the averaged Iberian Peninsula anomaly series was calculated using nine locations, with only four or five locations having data available for the years prior to 1970.This number of locations might not be sufficient to obtain a representative averaged result for the Iberian Peninsula.However, by way of an initial approximation, the number of locations was considered representative since, when using the same locations, Román et al. (2014a) obtained similar results in SW irradiation to Sánchez-Lorenzo et al. (2013a), who used more locations to obtain an averaged Iberian Peninsula series.

Conclusions
UVER irradiation in the Iberian Peninsula increased by 2.1 % dc −1 (annual) and 2.5 % dc −1 (summer) for the 1985-2011 period, and by 1.1 % dc −1 (annual) and 1.3 % dc −1 (summer) for the 1950-2011 period.The amount of ozone in the atmosphere is returning to pre-1980 levels due to the reduction in halogen gases subsequent to the Montreal Protocol.This decrease supports the belief that increased UVER radiation over the past 27 years has in large part been due to a reduction in the release of aerosols into the atmosphere, also reducing the cloud presence (brightening).However, increased UVER was mainly caused by ozone depletion during the 1950-2011 period.A significant UVER trend can become non-significant when uncertainty in the reconstructed data is taken into account.

Figure 2 .
Figure 2. The g function against the simulated diffuse UVER under cloudless conditions for four sunshine fraction (F ) intervals.A linear fit calculated with the N data, and its equation and correlation coefficient (r), are included in each panel.

Figure 3 .
Figure 3. Parameters a (left), b (middle), and B (right) calculated for different sunshine fraction (F ) intervals.The error bar represents the combined uncertainty.

Figure 4 .
Figure 4. Evolution of the annual UVER irradiation anomalies and their uncertainties for ten series.The red line corresponds to a linear fit between 1950 and 2011, and green lines to linear fits in the 1950-1984 and 1985-2011 periods.The solid black line is a 21-year Gaussian low-pass filter, and the years 1982 and 1991 are marked with a dashed black line.

Figure 5 .
Figure 5. UVER trends for different months and seasons, and for annual and three periods.The error bars are the 95 % confidence interval and their colour represents the location of the legend.The green circles represent statistically significant trends with 99 % confidence (p < 0.01), yellow squares represent statistically significant trends with 95 % confidence (p < 0.05), and red triangles represent non-statistically significant trends with at least 95 % confidence.

Table 1 .
Coordinates of the AEMet stations used, and number of data of the final reconstructed series which were filled with the M-SW model, the M-F model or with measurements.The total number of data of each final reconstructed series and the year when it began are included.
Figure 1.Distribution of selected Spanish stations located in the Iberian Peninsula.
Bilbao et al. (2011)ions into account.Hourly CMF was calculated for hourly SW and UVER irradiance with the mentioned available data.The reconstruction model proposed in this paper was based on the model developed byBilbao et al. (2011), suggesting the next relationship:

Table 2 .
Statistical estimators calculated with N pairs of measured and reconstructed UVER radiation for M-SW and M-F models and for different temporal resolutions.The hourly values were obtained using data with SZA below 80 • .

Table 4 .
Statistically significant UVER irradiation trends considered as the median of 10 000 trends (standard deviation in parentheses), with a confidence of 99 % (95 % marked with an asterisk), and the P (p < 0.05) and P (p < 0.01) values at various locations and seasons for the periods.