2.1. Study area

The SEA region is located at 36-38º N and covers 57 200 km², approximately 7.5% of the Turkish territory. The region presents a diverse topography. Its southern part is an almost homogeneous, flat area with a mean elevation of around 744 masl. This area extends around the lower Tigris River basin. On the other hand, the eastern part of the region presents a very complex topography, with steep elevation gradients reaching altitudes over 3800 m.

In terms of GHI, the SEA region is the most significant in Turkey with annual total solar radiation of 1460 kWh m^-2 and 2993 h of sunshine, while average values for Turkey are 1311 kWh m^-2 y^-1 and 2640 h, respectively (^{GDRE, 2018}). Kilis, Ceylanpınar, Bozova, Şırnak, and Mardin stations have GHI measurements in SEA. Table I href="?a5-143-158"> shows the geographical characteristics of these stations.

2.2. Cloudy conditions

Clouds and their accompanying weather patterns are among the most important atmospheric phenomena limiting solar radiation at the Earth’s surface. GHI has a nonlinear dependency on cloud cover due to its dynamic chaotic characteristics (^{Paulescu and Badescu, 2011}; ^{Chiacchio and Vitolo, 2012}). Therefore, clouds are an important source of uncertainty in the prediction of solar irradiance.

Cloudiness or cloud fraction is a meteorological measurement traditionally made by human observers. Cloud cover measurements are reported by meteorological services in SYNOP and METAR reports (^{González et al., 2012}). Because of their subjectivity, these measurements are not always accurate.

The daily sunshine index (s/S) was considered as a criterion, in order to determine the cloud conditions of the study periods. s/S was obtained by normalizing s, the bright sunshine duration, with S, the astronomical day length. The values of daily index greater than 0.75 were considered descriptive of clear-sky conditions and values between 0.25 and 0.75 were considered as partly cloudy conditions. Simulating solar radiation in clear sky conditions (the absence of clouds) is generally easier than in partly-cloudy ones.

2.3. Model setup

The WRF model (^{Skamarock et al., 2008}; ^{Wei et al., 2008}) used in the present study (v. 3.3) was integrated using three nested domains with grid spacing values of 36, 12, and 4 km. The SEA region was covered with 4-km grid resolution, as shown in Figure 1.

Fig. 1 Model domains and TSMS stations in the region (Ks: Kilis; B: Bozova; C: Ceylanpınar; M: Mardin, and S: Şırnak) 

The simulations were driven based on the National Center for Environmental Prediction (NCEP) Global Forecast System (GFS) data. The temporal resolution of the GFS analysis was 6 h, while the spatial resolution was 1º. Model configurations are given in Table II. The atmospheric column was decomposed into 35 vertical levels. Only the innermost domain was used in the evaluation. The WRF model ran for 72 h in each case. The evaluations of the 72-h GHI forecasts were based on the hourly values measured on August 20-22, 2011 for clear days and on May 8-10, 2011 for partly cloudy days for the given stations in SEA. These dates are chosen after examining observations of all stations during 2011. The model uses two-way nested boundary conditions. No cumulus parameterization was used for domain 3.

The WRF model offers numerous physical options such as PBL, surface physics, radiation, microphysics and cumulus parameterization. The WRF model was operated with four different configurations defined as WRF0, WRF1, WRF2 and WRF3. The physical parameterizations and scheme choices used in this study are presented in Table III.

Radiation is the main force that regulates the surface energy budget. An atmospheric RTM calculates the radiative transfer of electromagnetic radiation through the atmosphere. Both SW and longwave (LW) radiation parameterizations in an atmospheric model calculate radiation fluxes and heating rates (^{Price et al., 2013}). Specifically, radiation schemes resolve atmospheric heating from radiative flux divergence and surface downward LW and SW radiation. The Rapid Radiative Transfer Model (RRTM) for general circulation models (GCMs) (RRTMG) is the modified version of the RRTM for SW solar irradiance (^{Clough et al., 2005}; ^{Iacono et al., 2008}). In the present study, RRTMG (^{Iacono et al., 2008}), RRTM (^{Mlawer et al., 1997}), Fu-Liou-Gu (FLG) (^{Gu et al., 2011}), and Goddard radiation schemes (^{Chou and Suarez, 1994}) were used for SW and LW with different microphysics groups in four configurations. Table IV presents a comparison of these four radiation schemes.

For example, the number of bands for RRTM and RRTMG are equal to 14 for SW and 16 for LW, respectively. RRTM performs the calculations through the Discrete Ordinates Model (DISORT) for radiative transfer using four streams. As seen in Table IV, both models employ gaseous absorption coefficients extracted from the line-by-line radiation model and the computationally efficient correlated-k method for radiation transfer calculations (^{Goody et al., 1989}). Furthermore, the RRTMG scheme uses a two-stream algorithm for multiple scattering and then uses reduced set of g-intervals for integration over the extinction in each band (Table IV). Moreover, the RRTMG-SW includes the Monte-Carlo Independent Column Approximation (McICA) method to represent sub-grid cloud variability (^{Barker et al., 2006}, ²⁰⁰⁸). ^{Ruiz-Arias et al. (2013)} have shown that the RRTMG has high potential skill for long-term assessment of clear-sky GHI. Furthermore, ^{Mlawer et al. (1997)} showed that the RRTM-LW is an accurate and efficient correlated-k, LW RTM. FLG (^{Fu and Liou, 1993}) uses 6 SW and 12 LW spectral bands. FLG uses Delta 4-stream approximation for SW and Delta 2-stream approximation for LW.

The Goddard-SW solar radiation model has been developed at the NASA (^{Chou and Suarez, 1999}). Goddard-SW includes extinction by water vapor (H₂O), ozone (O₃), carbon dioxide (CO₂), nitrous oxide (N₂O), methane (CH₄), and CFC-11, CFC-12 (^{Chou et al., 2001}). The SW solar radiation model contains attenuation by H₂O, O₃, oxygen, CO₂, aerosols and clouds and the fluxes are calculated at 11 spectral bands, using the two-stream added approximation. Its accuracy is expected to be within a few W m^-2, whereas the atmospheric heating rate between 0.01 hPa and the surface is accurate to within 5%, relative to line-by-line calculations (^{Chou and Suarez, 1999}). ^{Iacono and Nehrkom (2010)} showed that the Goddard-SW model produces downward SW surface fluxes that are about 30 W m^-2 higher than RRTMG-SW. Moreover, the Goddard-LW uses 10 spectral bands. In this study, we employed the newest Goddard scheme under WRF v3.3.

A proper representation of clouds in solar irradiance prediction is rather problematic. ^{Mathiesen et al. (2013)} showed that the forecast accuracy is low in case of the regions with dynamic cloud systems. In the present study, the WRF Single-Moment 6-Class (WSM6) (^{Hong and Lim, 2006}) and the ^{Thompson (Thompson et al., 2004}) microphysics schemes were used for representing microphysical processes in clouds while only the Kain-Fritsch scheme (^{Kain and Fritsch, 1990}; ^{Kain, 2004}) was considered for cumulus parameterization.

In microphysics, the WSM6 scheme uses recent scientific developments regarding ice nucleation and includes a dependence on aerosol concentration (^{Skamarock et al., 2008}). Besides, the scheme includes a parameterization for calculating the direct radiative effects of several types of aerosol, including maritime, continental, urban, mineral dust, soot (black carbon), and sea salt. The WSM6 has been one of the options of microphysical processes in the WRF model since 2004. It is based on the revised ice-microphysics suggested by Hong, Dudhia and Chen (HDC) (^{Hong et al., 2004}; Hong and Lim, 2006). In addition to the characteristics mentioned above, the WSM6 behaves realistically in response to the appropriate grid resolvable forcing (^{Hong and Lim, 2006}).

The Thompson scheme, developed for use in WRF and other mesoscale models, is another scheme used for parameterization of microphysics. The Thompson microphysics scheme includes improvements to the earlier bulk scheme of ^{Reisner et al. (1998}). It was tested and compared with idealistic case studies, specifically for mid-latitudes (^{Hall et al., 2005}).

There are many different PBL parameterization schemes in the literature (^{Kleczek et al., 2014}). Many of them are designed to perform well during specific weather conditions. The variation of the predicted parameters such as temperature, mixing height and GHI with different parameterization schemes were examined by several investigators. For example, ^{Shin and Hong (2011)} compared five PBL schemes in the WRF model for a single day from a field program. They reported a large variation in the mixing height values computed by the five schemes both in daytime and nighttime. ^{Hu et al. (2010)} compared the 3-month mean diurnal variations of the simulated mixing heights, using different combinations of parameterization schemes in the WRF model with the observed values for Southeast Texas. Similarly, ^{Han et al. (2008)} evaluated five PBL schemes in the MM5 model for the East Asian domain and reported a large difference in the mixing heights predicted by the various combinations of parameterization schemes. ^{Ruiz-Arias et al. (2008)} conducted an evaluation study for the performance of different parameterizations of the model. In this study, the model’s performance was found highly dependent on the sky conditions. For example, solar radiation RMSE values were about one order of magnitude higher during broken-clouds and overcast conditions compared to clear-sky ones. More recently, ^{Zempila et al. (2016)} evaluated the shortwave radiation schemes such as Dudhia, updated RRTMG, updated Goddard and the GFDL schemes based on the WRF-model derived GHI predictions for cloudy and cloudless conditions. All SW radiation schemes presented better GHI predictions during the warm period due to the reduced presence of clouds.

The Monin-Obukhov similarity theory with Carlson-Boland viscous sublayer was also used for representing surface-layer parameterization; the Noah model was selected as land surface model for surface physics. The Yonsei University (YSU) scheme for PBL (^{Hong and Lim, 2006}) and the Kain-Fritsch scheme for cumulus parameterizations (^{Kain and Fritsch, 1990}) were used.

2.4. Error metrics

In the present study, we employed two error metrics for the evaluation of the model simulations according to the observations: mean bias error (MBE) and root mean square error (RMSE).

where x _e is the estimated data and x _m is the measured data. In addition to these metrics, normalized values rMBE and rRMSE are defined as relative MBE and relative RMSE, respectively. Nearest grid-point is used for verification with the station observations.

3.1. Partly cloudy case (May 8-10, 2011)

The period May 8-10, 2011 during which the s/S varied between 0.25 and 0.75, was selected as the partly-cloudy case. Seventy two-hour simulation results of the model based on the four different configurations named WRF0, WRF1, WRF2 and WRF3 for partly cloudy days at five stations are shown in Figure 2.

Fig. 2 Modelled and measured GHI values at the TSMS stations: (a) Bozova, (b) Ceylanpınar, (c) Kilis, (d) Mardin, and (e) Şırnak in the SEA region for May 8-10, 2011 period with partly cloudy conditions. 

The MBE and RMSE are computed for all stations and on an average throughout the SEA region according to Eqs. (1) and (2); the results are shown in Table V. The rMBE and rRMSE values at each station for all configurations in partly cloudy conditions are calculated and shown in Figures 3 and 4, respectively. MBEs are found positive indicating overestimation of GHI for the partly cloudy conditions at five stations and for all configurations (Table V). The findings are supported by the results of ^{Zempila et al. (2016)}. The overestimated WRF-derived hourly GHI values can also be clearly seen in Figure 2. As seen in Table V, MBE varies from 49.49 (WRF2) to 80.98 W m^-2 (WRF3) depending on the configuration and site of the region for the partly cloudy situation.

Fig. 3 rMBE values for the SEA region on May 8-10, 2011 for partly-cloudy conditions. 

Fig. 4 rRMSE values for the SEA region on May 8-10, 2011 for partly-cloudy conditions. 

Moreover, the partly-cloudy days on May 8-10, 2011 are linked with high RMSEs. For example, the regional average RMSE for the scheme configurations varied between 152.29 and 161.82 W m^-2 at SEA. These results are in agreement with the findings of ^{Mathiesen and Kleissl (2011)}, who also presented higher biases during partly-cloudy conditions in the USA. The MBE values reported for NAM forecasts under moderately-cloudy conditions in the USA are positive and vary between 10.8 and 62.6 W m^-2 while the RMSE values vary between 93.5 and 149 W m^-2.

Furthermore, from the relative view of the errors, the rMBE varies between 17 and 30% depending on the configuration and site while the rRMSE varies between 41 and 70%, respectively (Figs. 3 and 4). The difference between the forecasting values of WRF-derived hourly GHI and the observations shows that the WRF3 combination has the highest error among the four configurations in the region for partly-cloudy situation (Fig. 2). The lowest RMSE was found in Bozova (in the range of 139.49 and 151.40 W m^-2, depending on the configuration) while the highest was found at Ceylanpınar (177.23 and 191.17 W m^-2) and Mardin (173.52 and 179.68 W m^-2).

In summary, the WRF model with the four configurations used in the study did not adequately represent the fluctuations in the partly-cloudy case. This may be explained by the impact of clouds and aerosols (^{Sezen et al., 2013}). Although the impacts of the aerosols were not considered in this study, spring is the favorable season for dust transport from the southern desert areas to many regions of Turkey including SEA. Specifically Mardin is very well known as a dusty city. As a final remark for the partly-cloudy, case the WRF0 configuration had a better performance compared to other configurations.

3.2. Clear-sky case (August 20-22, 2011)

In this study, the August 20-22, 2011 period has been selected to represent the clear sky situations using s/S values higher than 0.75. The modelled and measured GHI values for 72 h at SEA sites are shown in Figure 5. As seen from Table VI, MBEs are positive, indicating overestimation of GHI for the clear-sky days at the five stations and for all configurations, similar to the partly cloudy case in May. The results of the WRF model at each site in the region, in terms of performance of each statistical metric, are also shown in Figures 6 and 7. The comparison of the four model configurations presents significant differences in the relative MBE and RMSE at the five sites in the region. The relative MBE is between 10 and 18% and the relative RMSE is between 20 and 39% throughout the region depending on the WRF model configuration (Figs. 6 and 7). Furthermore, from a regional perspective, the average MBE and RMSE values varied in the ranges of 36.73-48.26 W m^-2 and 71.24-81.13 W m^-2, respectively, depending on the configuration.

Fig. 5 Modeled and measured GHI values at the TSMS stations (a) Bozova, (b) Ceylanpınar, (c) Kilis, (d) Mardin, and (e) Şırnak in the SEA region for August 20-22, 2011, with clear-sky conditions. 

Fig. 6 rMBE values for the SEA region on August 20-22, 2011 for clear-sky conditions. 

Fig. 7 rRMSE values for the SEA region on August 20-22, 2011 for clear-sky conditions. 

According to the statistical metrics, the model performance is better for the clear sky case than the partly cloudy conditions, as expected.

GHI forecast by the WRF2 model configuration is close to the observations particularly at Şırnak and Mardin, with the lowest RMSE values of 55.42 and 57.43 W m^-2, respectively (Table VI). These sites are located at the rugged mountainous area with higher altitudes in SEA (Table I). RMSE values are supported by the findings of ^{Ohtake et al. (2015)} in Japan using the Japan Meteorological Agency (JMA) mesoscale model.

In the WRF2 model configuration, rMBE varies between 10 and 17% while rRMSE varies between 20 and 35% depending on the site within the study region. These results by the WRF2 configuration are in agreement with the results by ^{Diagne et al. (2014)}, who reported rMBE of 6.34% and rRMSE of 37.18% at Saint Pierre together with an rMBE of 16.08% and rRMSE of 57.3% at Tampan, two sites on Reunion Island in the Indian Ocean. Generally, the statistics of WRF2 indicated slightly better performance than WRF0 and WRF1 and significantly better performance than WRF3 for clear sky conditions. WRF3 presents poor results compared to the other configurations. Hence, we can emphasize that it is a combined effect of all parameterization schemes which translates into a good simulation of GHI by the WRF model except for the WRF3 combination.