Using Shortwave Radiation to Evaluate the HARMONIE-AROME Weather Model

Abstract

Evaluation of global shortwave irradiance forecasts from the HARMONIE-AROME weather prediction model is presented in this paper. We give examples of how such an evaluation can be used when testing a weather model or reanalysis product. We specifically use the non-dimensional clear sky and variability indices. We have tested seven months of HARMONIE-AROME 40h1.1 output against Danish global irradiance stations and 35 years of the Irish Met Éireann reanalysis (MÉRA) simulations. MÉRA, which is run with HARMONIE-AROME 38h1.2, is shown to have a significantly lower bias than the previously available global horizontal irradiance (GHI) reanalysis data from the ERA-Interim dataset. The Danish HARMONIE-AROME 40h1.1 has a negative bias during the summer months that is not seen in the Irish HARMONIE-AROME 38h1.2. For both model runs, we find a negative bias in the shortwave irradiance forecasts on days with thick clouds. This suggest that the model has too much cloud water in thick clouds.

1. Introduction

Forecasting clouds poses the biggest challenge in everyday weather model operations. This is a problem particularly for aviation and solar energy purposes. Reanalysis datasets have also been tested and found to be lacking in quality compared to satellite-derived cloud datasets [1]. Clearly, cloud forecasting needs to be improved. However, this is not easy because clouds are affected by several atmospheric processes, the model surface physics and assumed climatologies.

In order to improve cloud forecasts and reanalyses, it is necessary to perform a detailed quantitative evaluation. Cloud cover is often used to evaluate cloud forecasts. However, this is not a good metric due to inconsistent methods for assessing the cloud cover. Some ceilometers only measure clouds up to a height of 7 km (e.g., [2]), and thus, cirrus clouds are not included as part of the cloud cover. In some countries, the cloud cover is assessed by a human observer who can see the full sky unlike a ceilometer, which normally measures only at zenith. Tests have shown that this can cause markedly different cloud cover evaluation results across country boundaries [3].

Global horizontal irradiance (GHI), also referred to as “global radiation”, provides an objective and quantitative measure for evaluating cloud forecasts during the daytime. GHI is the shortwave irradiance integrated over all downward directions and wavelengths in the solar spectrum received by a horizontal surface of unit area. Furthermore, non-dimensional indices for solar energy resource assessment have been developed in recent years and decades that are very useful. One such index is the clear sky index (CSI), which is the GHI divided by the theoretical GHI during clear sky conditions, e.g., [4,5]. In this paper, we use the theoretical GHI clear sky model of [6,7], which includes coefficients that account for variable integrated atmospheric water vapour, aerosols and ozone. Another non-dimensional index that we will use here is the daily variability index (VI) as introduced by [8]. This VI is essentially the relative length of the daily curve of actual measured GHI to the daily curve of GHI under clear sky conditions. In this manner the relative variability of morning and evening periods with low GHI are not weighted disproportionally to the periods with high GHI during the middle of the day. In our analyses, we use the daily VI and a daily averaged GHI weighted CSI to make bins in which we can test the HARMONIE-AROME weather model GHI forecast under specific cloud conditions. Our hypothesis is that such a binned GHI evaluation will help in identifying weaknesses in the model.

The paper is outlined as follows. The observation networks and numerical weather prediction data are described in Section 2. Results and the discussion are included in Section 3 and Section 4 with conclusions in Section 5.

2. Datasets/Models

2.1. Observation Networks

2.1.1. The Irish Network

A GHI network was established in Ireland in 1954 when the first pyranometers were installed by the Irish Meteorological Service, Met Éireann, at Valentia Observatory in the southwest of the country. Over the next few decades, a further six stations were added to the network, the locations of which are shown in Figure 1a [9]. From 2008 onwards, the network underwent a process of automation with some of the original stations moving to more appropriate nearby sites. This automated network (referred to as the TUCSON network) is shown in Figure 1b. The stations are equipped with Kipp and Zonen CMP11 pyranometers (2628 XH Delft, Netherlands) which are “Secondary Standard pyranometers”, but are more accurate than “First Class Pyranometers”.

2.1.2. The Danish Network

Since the beginning of 2001, the Danish Meteorological Institute, DMI, has measured GHI at 13 stations across Denmark. In the years up to 2011, this network was expanded to include 28 stations that are all still operational. The distribution of these stations across the country can be seen in Figure 2. Before 2001, the sunshine duration has been monitored since 1920 [10]. The GHI stations are all equipped with Star 8101 pyranometers (Schenk GmbH, A-1210 Wien, Austria), which are “First Class Pyranometers” according to the ISO 9060 standard.

2.2. Models

2.2.1. MÉRA

Met Éireann has carried out a 35-year very high resolution (2.5-km horizontal grid) regional climate reanalysis for Ireland using the ALADIN-HIRLAM numerical weather prediction system and in particular the HARMONIE-AROME canonical configuration of this system (see [11] for details on the model, acronym definitions and an overview of the reanalysis). This Met Éireann reanalysis, abbreviated MÉRA, spans the years 1981 to 2015. An analysis of surface parameters including screen-level and soil temperatures, 10-m wind speeds, mean sea-level pressure and precipitation is presented in [11,12]. This paper extends the analysis of MÉRA to include an evaluation of GHI with particular emphasis on its use in the validation of clouds.

Version 38h1.2 of the ALADIN-HIRLAM system was used to carry out the MÉRA simulations. This canonical configuration of HARMONIE-AROME was run on a horizontal grid of 2.5-km spacing, using 65 vertical levels and a model top of 10 hPa. The data assimilation component of the model is described in [13,14]; the forecast component is described in [15] with more recent updates included in [16]. The SURFEX externalized surface scheme Version 7.2 [17] is used for surface data assimilation and the modelling of surface processes. A summary of the model details is presented in

Table 1. HARMONIE-AROME configuration used for MÉRA. SYNOP, synoptic station .

Model Version	HARMONIE-AROME 38h1.2
Domain	540 × 500 grid points ($\Delta x=$ 2.5 km)
Vertical Levels	65 levels up to 10 hPa, first level at approximately 12 m
Forecast Cycle	3 h
Data Assimilation	Optimal interpolation for surface parameters
	3DVAR assimilation for upper air parameters
Observations	Pressure from SYNOP, SHIP and DRIBU
	Temperature and winds from AIREP and AMDAR
	Winds from PILOT
	Temperature, winds and humidity from TEMP
Forecast	3-h forecasts, but a 33-h forecast at 00 Z

.

Conventional observations (i.e., from synoptic stations (SYNOP), ships, drifting buoys (DRIBU), radiosonde ascents (TEMP) and aircraft (AIREP, AMDAR) were assimilated. These observations are the same as those used by the European Centre for Medium Range Weather Forecasts’ (ECMWF) ERA-Interim reanalysis [18]. Locally available SYNOP observations were used to fill gaps in the ERA-Interim SYNOP observation archive to ensure cycle continuity and successful data assimilation each cycle. MÉRA uses a three-hour forecast cycle with surface and upper-air data assimilation. Three-hour forecasts were produced for each cycle except the midnight (00 UTC) cycle when a 33-h forecast was produced.

2.2.2. Danish Model Runs

Operational weather model runs from DMI are included for the seven-month period from June to December 2017. Specifically, output from the HARMONIE-AROME 40h1.1 model run covering the Northern Europe A (NEA) domain of DMI is included. The physics in this model version has been described in detail by [16] with the specific NEA setup described by [19]. In

Table 2. HARMONIE-AROME configuration used for the Danish Meteorological Institute (DMI) Northern Europe A (NEA) domain.

Model Version	HARMONIE-AROME 40h1.1
Domain	1280 × 1080 grid points ($\Delta x=$ 2.5 km)
Vertical Levels	65 levels up to 10 hPa, first level at approximately 12 m
Forecast Cycle	3 h
Data Assimilation	Optimal interpolation for surface parameters
	3DVAR assimilation for upper air parameters
Observations	Pressure from SYNOP, SHIP and DRIBU
	Temperature and winds from AIREP, AMDAR and MODE-S-EHS
	Radiances from AMSU-A, MHS and ATMS
	Atmospheric Motion Vectors (AMV) from MSG
	GNSS-RO bending angles from Metop-A, Metop-B and COSMIC
	Radar volume data from OPERA
	Winds from PILOT
	Temperature, winds and humidity from TEMP
Forecast	60-h forecasts

, an overview of the NEA setup is given. As can be seen, this setup has the same horizontal and vertical resolution as MÉRA, but includes significantly more observations for data assimilation. In particular, the DMI setup includes additional aircraft data (MODE-S, EHS), satellite data (AMSU-A, MHS, ATMS, MSG, Metop-A, Metop-B, COSMIC) and radar data from the European OPERA radar network.

HARMONIE-AROME 40h1.1 is the model version release that followed immediately after HARMONIE-AROME 38h1.2 used for MÉRA. The differences between the model versions as run operationally at DMI (NEA; 40h1.1) and run for the long-term reanalysis at Met Éireann (MÉRA; 38h1.2) are too extensive to list here. The following major changes, which affect the cloud cover and cloud transmittance, we would, however, like to mention. Firstly, the externalized surface scheme has been updated to Version 7.3 [20]. Secondly, the so-called cloud-inhomogeneity factor has been changed from 0.7 to 1.0 as suggested by [21], which means that the cloud optical thickness is no longer reduced by 30% before the cloud transmittance is computed. Thirdly, the Cuxart, Bougeault and Redelsperger (CBR) turbulence scheme [22] has been replaced by the turbulence scheme from the RACMO model [23].

Note that all of the results for Irish stations are compared to all of the MÉRA reanalysis dataset, which spans the years 1981 to 2015. The MÉRA domain does not include Denmark. Hence, the Danish analysis has been carried out using using HARMONIE-AROME operational NWP simulations for May to December 2017.

3. Results

3.1. GHI Statistics

The seven original stations in the Irish network were used to create Figure 3, which shows the daily mean GHI as a function of the day of the year for observations and the MÉRA dataset. 24-h GHI MÉRA forecasts for the years 1981 to 2015 are included. These forecasts were initiated at 00 UTC. Output from the ERA-Interim reanalysis is also included as a comparison and highlights the benefits of the high resolution MÉRA dataset over the coarse resolution ERA-Interim (79 km grid spacing approximately). Twenty four-hour forecasts of accumulated GHI were extracted from the MÉRA and ERA-Interim datasets, interpolated to each Irish station using bilinear interpolation and averaged over the full day. The biases in ERA-Interim are small and mostly negative in winter (a few Wm⁻²), but are otherwise positive and up to 40 Wm⁻² in summer. On the other hand, the biases in MÉRA are smaller (±10 Wm⁻²) (positive except for a period in July and August, where the bias is slightly negative overall). Because the MÉRA dataset is clearly more suitable than ERA-Interim as regards GHI because of its higher resolution, the remainder of the reanalysis comparisons are restricted to the MÉRA dataset. In addition, hourly outputs of GHI are available from MÉRA, whereas only daily accumulations of GHI are available from ERA-Interim.

In Figure 4, the observed and modelled time series of the average of all Danish GHI stations are shown. The averages are hourly averages over the full day. As only seven months of data were available, these results are less clear than the MÉRA results, but an average negative bias of −10 to −20 Wm⁻² can be seen in the period from June to September, while the last three months have little bias.

To analyse the biases in MÉRA relative to observations in more detail, the mean daily biases in GHI by month averaged over all available years (2008 to 2014) were calculated for the TUCSON network of 20 stations. The biases in MÉRA relative to observations are positive between February and June and in September, the exceptions being some of the western stations. In July, the biases are positive in the east, but negative in the west. In August, the biases are negative except for a few eastern stations. From October to January, the biases are mostly negative with some exceptions for stations near the coasts (see Figure 5 for biases for a selection of months). RMSE (root mean square error) plots are also included for completeness in Figure 6.

Likewise for the Danish stations, monthly statistics for four selected months (June, July, August, September) are shown in Figure 7, Figure 8 and Figure 9. This includes maps of bias (Figure 7) and RMSE (Figure 8) for each station and a plot of the same data shown graphically in order to make the variations easier to quantify. From the Danish model runs, hourly GHI forecasts up to 12 hours ahead initiated at 00 and 12 UTC are used in these statistics.

3.2. Clear Sky Index

The definition of clear sky index (CSI) used here involves the ratio of GHI divided by the theoretical GHI during clear sky conditions, e.g., [4,5], which is dependent on the location, date and time. In particular, we used the theoretical GHI clear sky model of [6,7], which includes coefficients that account for variable integrated atmospheric water vapour and, if available, aerosols and ozone. In the Danish weather model runs, the integrated water vapour variations were available and have been accounted for. In the comparisons involving MÉRA, we assumed an atmosphere with an integrated water vapour load of 2.5 g·cm⁻², which is a typical load for a mid-latitude location like Ireland. Variable atmospheric water vapour was of course included in the MÉRA run, but these data are not readily available.

Figure 10 shows a time series of weighted daily CSI for MÉRA, calculated using hourly CSI weighted by hourly GHI and summed over the day for hours when the sun was above the horizon. Note that GHI for each hour was calculated using the accumulated GHI for consecutive hours extracted from the 0 to 33 h forecasts run at 00 Z . CSI was also averaged over all stations, and a 30-day running mean was applied to aid visualisation. CSI trends based on observations and on MÉRA are shown in the figure, as well as the bias in MÉRA relative to observations. On average, there is very good agreement between MÉRA and observations, but monthly mean biases in daily weighted CSI of ±0.2 can be seen. Figure 11 shows a frequency distribution of weighted daily CSI for the TUCSON stations. Overall, there is good agreement, but a small shift towards higher CSIs in MÉRA for most of the stations. MÉRA has significantly more cases for the lowest CSI bin, i.e., for clear sky indices from 0.0 to 0.1, than the observations.

Similarly, the frequency distribution of the Danish daily weighted CSIs are shown in Figure 12. In this case, the modelled CSIs never exceed 1.0, but the observed CSIs do so in approximately 3% of the cases. In the lower part of the frequency distribution, it can be seen that the NEA HARMONIE-AROME runs have significantly more cases in the lowest CSI bin, i.e., for clear sky indices from 0.0 to 0.1, than the observations, as was also seen for the MÉRA data (Figure 11).

3.3. Variability Index

The daily variability index (VI) of [8] is essentially the relative length of the daily curve of actual measured GHI to the daily curve of GHI in clear sky conditions. See Equation (1) below, where CSGHI is the GHI under clear sky conditions, k is the time index (here, hour of the day) and $\Delta$t is 60 min in this case, as we only have hourly GHI available from the model runs. Using this definition, the relative variability of morning and evening periods with low GHI are not weighted disproportionately to the periods with high GHI during midday. The advantage of using the VI is that a high VI indicates situations with subgrid-scale variability that is unresolved in the weather model runs. It is interesting to test whether the models have errors due to sub-grid cloud variability. In general, the VI provides a method to classify cloud situations.

$$\mathrm{VI}=\frac{{\sum }_{k=2}^n\sqrt{{({\mathrm{GHI}}_k-{\mathrm{GHI}}_{k-1})}^2+\Delta t^2}}{{\sum }_{k=2}^n\sqrt{{({\mathrm{CSGHI}}_k-{\mathrm{CSGHI}}_{k-1})}^2+\Delta t^2}}·\sqrt{\Delta t}$$

(1)

Figure 13 shows a time series of daily VI averaged over all TUCSON stations with a 30-day running mean applied to aid visualisation. VI trends based on observations and on MÉRA are shown in the figure, as well as the bias in MÉRA relative to observations. On average, there is very good agreement between MÉRA and observations with monthly mean biases in the range ±2. A pronounced seasonal cycle in VI is also evident over the seven-year period with the greatest variability during spring/summer when convective events are more common in Ireland compared to lower variability in autumn/winter, where overcast days with low stratus cloud are more typical. Similar results were found for the HARMONIE-AROME NEA run (not shown here).

The final part of this analysis merges both the CSI and VI with the aim of testing the GHI forecasts from the HARMONIE-AROME weather model under specific cloud conditions in order to identify weaknesses in the model. Figure 13 shows that the VI for Ireland mostly ranges between five and 10. Figure 14 shows very good agreement between MÉRA and observations for weighted daily CSI for VI values in the range of 10 to 15 (i.e., the most variable days cloud-wise). Here, we have plotted the weighted daily CSI for three VI bin ranges (0 to 5, 5 to 10 and 10 to 15). The middle VI range (5 to 10) shows good agreement overall with a small shift towards higher CSIs in MÉRA. The lowest VI range (0 to 5) shows a negative bias, which is consistent with the results for the lowest CSI bin in Figure 11.

Another way of visualising the CSI and VI results is shown in Figure 15. Here, binned 2D histograms of the modelled and observed CSI values are shown for the HARMONIE-AROME NEA model compared with the Danish GHI stations. These are similar to contingency tables used for precipitation evaluation and are made for four VI bins: 0 to 5, 5 to 7.5, 7.5 to 10 and 10 to 15. Ideally, all the data should follow the diagonal at which the model CSI values are equal to the observed CSI values. Here, negative deviations in the model data are seen for the lowest CSI and VI values in the upper left panel. This is consistent with the frequency distribution seen for the lowest CSI bin in Figure 12. Another negative deviation is seen for the highest CSI values in the bin that contains most of the cases in the lower left panel of Figure 15.

4. Discussion

The high resolution MÉRA dataset captures the variation in global SW irradiance much better than the coarser resolution ERA-Interim dataset. The bias in ERA-Interim relative to observations is mostly positive and up to 30 Wm⁻² in summer. On the other hand, MÉRA has much lower biases on average (mostly within 10 Wm⁻² and positively biased during the first half of the year, but negatively biased during the last half of the summer months). A similar analysis was carried out using data from the Irish TUCSON network of 20 stations available from 2008 onwards where the data are archived at one-minute resolution. This dataset yielded similar results (not included) to those shown in Figure 3. The MÉRA dataset is thus better suited for solar resource assessment than ERA-Interim for the region covered, i.e., Ireland and the United Kingdom.

Analyses of both datasets as a function of the day of the year and month by month show the variations of the model bias during the course of the year. This could be related to the monthly climatological fields that are used in the model. Several of these are needed to simulate the surface variability such as the evapotranspiration of the assumed vegetation. In countries with a large crop fraction, the vegetation varies greatly, in particular as a result of harvest. This also affects the cloud cover and, in turn, the GHI. Monthly aerosol climatologies [24,25] are also assumed in HARMONIE-AROME, which affect the GHI [26,27].

In Figure 14 and Figure 15, there is a shift towards lower CSI for MÉRA compared to observations for the low VI range (0 to 5), which is consistent with cloud water loads that are too high. Tests against satellite-derived cloud water loads also show this [7]. This bias does not count much for the overall GHI statistics, but it could be important for the cloud microphysical processes in the model that initiate precipitation.

Although the differences between the results shown here for the NWP model output for Ireland and Denmark could be due to geographical differences, we find it likely that changes in the HARMONIE-AROME model from Version 38h1.2 to Version 40h1.1 are important. In particular, the change of the cloud inhomogeneity factor from 0.7 to 1.0, as discussed by [21], reduces the cloud transmittance. This means that the cloud transmittance is lower in Version 40h1.1 as run for HARMONIE-AROME NEA than in Version 38h1.2 used for MÉRA. In theory, cloud inhomogeneity effects are largest when the clouds are most variable. This fits with the negative bias seen in Figure 9 and in the lower left panel of Figure 15, which is not seen in the right-hand panel of Figure 14. We suggest a cloud inhomogeneity factor that is reduced when variable clouds are present to be tested using the next model version. Another possible explanation for this is too high aerosol optical depth assumed in the model climatologies over Denmark. This possibility should also be tested.

In order to quality-assure GHI measurements, the ideal radiation station includes a three-component measurement setup, where the direct, diffuse and global irradiance are measured simultaneously. Since the pyranometer calibration can drift over time, such a redundancy of irradiance measurements is recommended [28]. We do not have country-wide measurements of this type, which would have enabled us to draw more firm conclusions from these analyses.

When computing the CSI values for the Danish NEA data, we accounted for the variations in the integrated water vapour. This made the CSI value range more homogeneous across the year and is recommended for computing clear sky irradiance and CSI values. In Figure 11, the CSI frequency distribution plot of the MÉRA data displays CSI values greater than 1.0. These occur during winter, when the atmosphere is very dry and the modelled CSI therefore becomes too high. In Figure 12, CSI values larger than 1.0 do not occur for the HARMONIE-AROME NEA data, but they still occur for the observations. This could be due to too high aerosol absorption and scattering in the assumed clear sky model [6,7].

5. Conclusions

The analyses presented here demonstrate new aspects of how to evaluate and improve numerical weather models. These should be made a standard part of the overall weather model evaluation analyses. In particular, ground-based GHI stations such as those used here give information about the complex cloud overlap and 3D cloud effects that occur at the sub-grid scale level in modelled and satellite-derived data.

Specifically, we have shown the average GHI bias for individual days of the year for the Met Éireann reanalysis (MÉRA) with HARMONIE-AROME 38h1.2 to be within ±10 Wm² and to be significantly less than the corresponding bias of the ERA-Interim dataset. We have also shown that HARMONIE-AROME 40h1.1 operational runs for Denmark have a larger negative GHI bias during summer than the HARMONIE-AROME 38h1.2 runs for Ireland. In both cases, we find that on days with very low GHI, HARMONIE-AROME has too low GHI. This indicates that too much cloud water is present in the model in cases with thick stratiform clouds.

We suggest that the GHI evaluation be performed as a standard for HARMONIE-AROME and weather models in general and the model cloud physics be adjusted to improve the issues that we have demonstrated.